createFields.H
Go to the documentation of this file.
1Info<< "Reading thermophysical properties\n" << endl;
2
3autoPtr<psiuReactionThermo> pThermo
4(
5 psiuReactionThermo::New(mesh)
6);
7psiuReactionThermo& thermo = pThermo();
8thermo.validate(args.executable(), "ha", "ea");
9
10basicSpecieMixture& composition = thermo.composition();
11
12volScalarField rho
13(
14 IOobject
15 (
16 "rho",
17 runTime.timeName(),
18 mesh,
19 IOobject::NO_READ,
20 IOobject::AUTO_WRITE
21 ),
22 thermo.rho()
23);
24
25volScalarField& p = thermo.p();
26
27volScalarField& b = composition.Y("b");
28Info<< "min(b) = " << min(b).value() << endl;
29
30
31Info<< "\nReading field U\n" << endl;
32volVectorField U
33(
34 IOobject
35 (
36 "U",
37 runTime.timeName(),
38 mesh,
39 IOobject::MUST_READ,
40 IOobject::AUTO_WRITE
41 ),
42 mesh
43);
44
46
47mesh.setFluxRequired(p.name());
48
49Info<< "Creating turbulence model\n" << endl;
50autoPtr<compressible::turbulenceModel> turbulence
51(
52 compressible::turbulenceModel::New
53 (
54 rho,
55 U,
56 phi,
57 thermo
58 )
59);
60
61
62Info<< "Creating field Xi\n" << endl;
63volScalarField Xi
64(
65 IOobject
66 (
67 "Xi",
68 runTime.timeName(),
69 mesh,
70 IOobject::MUST_READ,
71 IOobject::AUTO_WRITE
72 ),
73 mesh
74);
75
76
77Info<< "Creating the unstrained laminar flame speed\n" << endl;
78autoPtr<laminarFlameSpeed> unstrainedLaminarFlameSpeed
79(
80 laminarFlameSpeed::New(thermo)
81);
82
83
84Info<< "Reading strained laminar flame speed field Su\n" << endl;
85volScalarField Su
86(
87 IOobject
88 (
89 "Su",
90 runTime.timeName(),
91 mesh,
92 IOobject::MUST_READ,
93 IOobject::AUTO_WRITE
94 ),
95 mesh
96);
97
98dimensionedScalar SuMin = 0.01*Su.average();
99dimensionedScalar SuMax = 4*Su.average();
100
101Info<< "Calculating turbulent flame speed field St\n" << endl;
102volScalarField St
103(
104 IOobject
105 (
106 "St",
107 runTime.timeName(),
108 mesh,
109 IOobject::NO_READ,
110 IOobject::AUTO_WRITE
111 ),
112 Xi*Su
113);
114
115
116multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
117
118if (composition.contains("ft"))
119{
120 fields.add(composition.Y("ft"));
121}
122
124fields.add(thermo.he());
125fields.add(thermo.heu());
126
127#include "createDpdt.H"
128
129#include "createK.H"
130
131#include "createMRF.H"
132#include "createFvOptions.H"
surfaceScalarField & phi
const word & executable() const noexcept
Name of executable without the path.
Definition: argListI.H:51
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
U
Definition: pEqn.H:72
volScalarField & p
Creates and initialises the face-flux field phi.
dynamicFvMesh & mesh
engineTime & runTime
compressible::turbulenceModel & turbulence
zeroField Su
Definition: alphaSuSp.H:1
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
Foam::argList args(argc, argv)
Info<< "Creating the unstrained laminar flame speed\n"<< endl;autoPtr< laminarFlameSpeed > unstrainedLaminarFlameSpeed(laminarFlameSpeed::New(thermo))
volScalarField & b
Definition: createFields.H:27
Info<< "Reading thermophysical properties\n"<< endl;autoPtr< psiReactionThermo > pThermo(psiReactionThermo::New(mesh))
basicSpecieMixture & composition
Definition: createFields.H:8
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:97