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
30Info<< "\nReading field U\n" << endl;
31volVectorField U
32(
33 IOobject
34 (
35 "U",
36 runTime.timeName(),
37 mesh,
38 IOobject::MUST_READ,
39 IOobject::AUTO_WRITE
40 ),
41 mesh
42);
43
45
46mesh.setFluxRequired(p.name());
47
48Info<< "Creating turbulence model\n" << endl;
49autoPtr<compressible::RASModel> turbulence
50(
51 compressible::New<compressible::RASModel>
52 (
53 rho,
54 U,
55 phi,
56 thermo
57 )
58);
59
60
61Info<< "Creating the unstrained laminar flame speed\n" << endl;
62autoPtr<laminarFlameSpeed> unstrainedLaminarFlameSpeed
63(
64 laminarFlameSpeed::New(thermo)
65);
66
67
68Info<< "Reading strained laminar flame speed field Su\n" << endl;
69volScalarField Su
70(
71 IOobject
72 (
73 "Su",
74 runTime.timeName(),
75 mesh,
76 IOobject::MUST_READ,
77 IOobject::AUTO_WRITE
78 ),
79 mesh
80);
81
82Info<< "Reading field betav\n" << endl;
83volScalarField betav
84(
85 IOobject
86 (
87 "betav",
88 mesh.facesInstance(),
89 mesh,
90 IOobject::MUST_READ,
91 IOobject::NO_WRITE
92 ),
93 mesh
94);
95
96Info<< "Reading field Lobs\n" << endl;
97volScalarField Lobs
98(
99 IOobject
100 (
101 "Lobs",
102 mesh.facesInstance(),
103 mesh,
104 IOobject::MUST_READ,
105 IOobject::NO_WRITE
106 ),
107 mesh
108);
109
110Info<< "Reading field CT\n" << endl;
111volSymmTensorField CT
112(
113 IOobject
114 (
115 "CT",
116 mesh.facesInstance(),
117 mesh,
118 IOobject::MUST_READ,
119 IOobject::NO_WRITE
120 ),
121 mesh
122);
123
124Info<< "Reading field Nv\n" << endl;
125volScalarField Nv
126(
127 IOobject
128 (
129 "Nv",
130 mesh.facesInstance(),
131 mesh,
132 IOobject::MUST_READ,
133 IOobject::NO_WRITE
134 ),
135 mesh
136);
137
138Info<< "Reading field nsv\n" << endl;
139volSymmTensorField nsv
140(
141 IOobject
142 (
143 "nsv",
144 mesh.facesInstance(),
145 mesh,
146 IOobject::MUST_READ,
147 IOobject::NO_WRITE
148 ),
149 mesh
150);
151
152IOdictionary PDRProperties
153(
154 IOobject
155 (
156 "PDRProperties",
157 runTime.constant(),
158 mesh,
159 IOobject::MUST_READ_IF_MODIFIED,
160 IOobject::NO_WRITE
161 )
162);
163
164//- Create the drag model
165autoPtr<PDRDragModel> drag = PDRDragModel::New
166(
167 PDRProperties,
168 *turbulence,
169 rho,
170 U,
171 phi
172);
173
174//- Create the flame-wrinkling model
175autoPtr<XiModel> flameWrinkling = XiModel::New
176(
177 PDRProperties,
178 thermo,
179 *turbulence,
180 Su,
181 rho,
182 b,
183 phi
184);
185
186Info<< "Calculating turbulent flame speed field St\n" << endl;
187volScalarField St
188(
189 IOobject
190 (
191 "St",
192 runTime.timeName(),
193 mesh,
194 IOobject::NO_READ,
195 IOobject::AUTO_WRITE
196 ),
197 flameWrinkling->Xi()*Su
198);
199
200
201multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
202
203if (composition.contains("ft"))
204{
205 fields.add(composition.Y("ft"));
206}
207
209fields.add(thermo.he());
210fields.add(thermo.heu());
212
213#include "createDpdt.H"
214
215#include "createK.H"
216
217#include "createMRF.H"
218#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
const volScalarField & betav
Foam::argList args(argc, argv)
autoPtr< XiModel > flameWrinkling
Create the flame-wrinkling model.
Definition: createFields.H:175
Info<< "Creating the unstrained laminar flame speed\n"<< endl;autoPtr< laminarFlameSpeed > unstrainedLaminarFlameSpeed(laminarFlameSpeed::New(thermo))
Info<< "Reading strained laminar flame speed field Su\n"<< endl;volScalarField Su(IOobject("Su", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field betav\n"<< endl;volScalarField betav(IOobject("betav", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Lobs\n"<< endl;volScalarField Lobs(IOobject("Lobs", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field CT\n"<< endl;volSymmTensorField CT(IOobject("CT", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Nv\n"<< endl;volScalarField Nv(IOobject("Nv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field nsv\n"<< endl;volSymmTensorField nsv(IOobject("nsv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);IOdictionary PDRProperties(IOobject("PDRProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));autoPtr< PDRDragModel > drag
Definition: createFields.H:165
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