createFluidFields.H
Go to the documentation of this file.
1 // Initialise fluid field pointer lists
2 PtrList<rhoReactionThermo> thermoFluid(fluidRegions.size());
3 PtrList<volScalarField> rhoFluid(fluidRegions.size());
4 PtrList<volVectorField> UFluid(fluidRegions.size());
5 PtrList<surfaceScalarField> phiFluid(fluidRegions.size());
6 PtrList<uniformDimensionedVectorField> gFluid(fluidRegions.size());
7 PtrList<uniformDimensionedScalarField> hRefFluid(fluidRegions.size());
8 PtrList<uniformDimensionedScalarField> pRefFluid(fluidRegions.size());
9 PtrList<volScalarField> ghFluid(fluidRegions.size());
10 PtrList<surfaceScalarField> ghfFluid(fluidRegions.size());
11 PtrList<compressible::momentumTransportModel>
13 PtrList<rhoReactionThermophysicalTransportModel>
15 PtrList<CombustionModel<rhoReactionThermo>> reactionFluid(fluidRegions.size());
16 PtrList<volScalarField> p_rghFluid(fluidRegions.size());
17 PtrList<radiationModel> radiation(fluidRegions.size());
18 PtrList<volScalarField> KFluid(fluidRegions.size());
19 PtrList<volScalarField> dpdtFluid(fluidRegions.size());
20 PtrList<multivariateSurfaceInterpolationScheme<scalar>::fieldTable>
21  fieldsFluid(fluidRegions.size());
22 
23 List<scalar> initialMassFluid(fluidRegions.size());
24 
25 PtrList<IOMRFZoneList> MRFfluid(fluidRegions.size());
26 PtrList<fv::options> fluidFvOptions(fluidRegions.size());
27 
28 // Populate fluid field pointer lists
30 {
31  Info<< "*** Reading fluid mesh thermophysical properties for region "
32  << fluidRegions[i].name() << nl << endl;
33 
34  Info<< " Adding to thermoFluid\n" << endl;
36 
37  Info<< " Adding to rhoFluid\n" << endl;
38  rhoFluid.set
39  (
40  i,
41  new volScalarField
42  (
43  IOobject
44  (
45  "rho",
46  runTime.timeName(),
47  fluidRegions[i],
48  IOobject::NO_READ,
49  IOobject::AUTO_WRITE
50  ),
51  thermoFluid[i].rho()
52  )
53  );
54 
55  Info<< " Adding to UFluid\n" << endl;
56  UFluid.set
57  (
58  i,
59  new volVectorField
60  (
61  IOobject
62  (
63  "U",
64  runTime.timeName(),
65  fluidRegions[i],
66  IOobject::MUST_READ,
67  IOobject::AUTO_WRITE
68  ),
69  fluidRegions[i]
70  )
71  );
72 
73  Info<< " Adding to phiFluid\n" << endl;
74  phiFluid.set
75  (
76  i,
78  (
79  IOobject
80  (
81  "phi",
82  runTime.timeName(),
83  fluidRegions[i],
84  IOobject::READ_IF_PRESENT,
85  IOobject::AUTO_WRITE
86  ),
88  & fluidRegions[i].Sf()
89  )
90  );
91 
92  Info<< " Adding to gFluid\n" << endl;
93  gFluid.set
94  (
95  i,
97  (
98  IOobject
99  (
100  "g",
101  runTime.constant(),
102  fluidRegions[i],
103  IOobject::MUST_READ,
104  IOobject::NO_WRITE
105  )
106  )
107  );
108 
109  Info<< " Adding to hRefFluid\n" << endl;
110  hRefFluid.set
111  (
112  i,
114  (
115  IOobject
116  (
117  "hRef",
118  runTime.constant(),
119  fluidRegions[i],
120  IOobject::READ_IF_PRESENT,
121  IOobject::NO_WRITE
122  ),
124  )
125  );
126 
127  Info<< " Adding to pRefFluid\n" << endl;
128  pRefFluid.set
129  (
130  i,
132  (
133  IOobject
134  (
135  "pRef",
136  runTime.constant(),
137  fluidRegions[i],
138  IOobject::READ_IF_PRESENT,
139  IOobject::NO_WRITE
140  ),
142  )
143  );
144 
145  dimensionedScalar ghRef(- mag(gFluid[i])*hRefFluid[i]);
146 
147  Info<< " Adding to ghFluid\n" << endl;
148  ghFluid.set
149  (
150  i,
151  new volScalarField
152  (
153  "gh",
154  (gFluid[i] & fluidRegions[i].C()) - ghRef
155  )
156  );
157 
158  Info<< " Adding to ghfFluid\n" << endl;
159  ghfFluid.set
160  (
161  i,
163  (
164  "ghf",
165  (gFluid[i] & fluidRegions[i].Cf()) - ghRef
166  )
167  );
168 
169  Info<< " Adding to turbulenceFluid\n" << endl;
170  turbulenceFluid.set
171  (
172  i,
174  (
175  rhoFluid[i],
176  UFluid[i],
177  phiFluid[i],
178  thermoFluid[i]
179  ).ptr()
180  );
181 
182  Info<< " Adding to thermophysicalTransport\n" << endl;
184  (
185  i,
187  (
188  turbulenceFluid[i],
189  thermoFluid[i]
190  ).ptr()
191  );
192 
193  Info<< " Adding to reactionFluid\n" << endl;
194  reactionFluid.set
195  (
196  i,
198  (
199  thermoFluid[i],
200  turbulenceFluid[i]
201  )
202  );
203 
204  p_rghFluid.set
205  (
206  i,
207  new volScalarField
208  (
209  IOobject
210  (
211  "p_rgh",
212  runTime.timeName(),
213  fluidRegions[i],
214  IOobject::MUST_READ,
215  IOobject::AUTO_WRITE
216  ),
217  fluidRegions[i]
218  )
219  );
220 
221  // Force p_rgh to be consistent with p
222  p_rghFluid[i] = thermoFluid[i].p() - rhoFluid[i]*ghFluid[i] - pRefFluid[i];
223 
224  fluidRegions[i].setFluxRequired(p_rghFluid[i].name());
225 
226  Info<< " Adding to radiationFluid\n" << endl;
227  radiation.set
228  (
229  i,
231  );
232 
234 
235  Info<< " Adding to KFluid\n" << endl;
236  KFluid.set
237  (
238  i,
239  new volScalarField
240  (
241  "K",
242  0.5*magSqr(UFluid[i])
243  )
244  );
245 
246  Info<< " Adding to dpdtFluid\n" << endl;
247  dpdtFluid.set
248  (
249  i,
250  new volScalarField
251  (
252  IOobject
253  (
254  "dpdt",
255  runTime.timeName(),
256  fluidRegions[i]
257  ),
258  fluidRegions[i],
260  (
261  thermoFluid[i].p().dimensions()/dimTime,
262  0
263  )
264  )
265  );
266 
267  Info<< " Adding to fieldsFluid\n" << endl;
268  fieldsFluid.set
269  (
270  i,
271  new multivariateSurfaceInterpolationScheme<scalar>::fieldTable
272  );
273  forAll(thermoFluid[i].composition().Y(), j)
274  {
275  fieldsFluid[i].add(thermoFluid[i].composition().Y()[j]);
276  }
277  fieldsFluid[i].add(thermoFluid[i].he());
278 
279  Info<< " Adding MRF\n" << endl;
280  MRFfluid.set
281  (
282  i,
283  new IOMRFZoneList(fluidRegions[i])
284  );
285 
286  Info<< " Adding fvOptions\n" << endl;
287  fluidFvOptions.set
288  (
289  i,
290  new fv::options(fluidRegions[i])
291  );
292 
293  turbulenceFluid[i].validate();
294 }
PtrList< uniformDimensionedScalarField > pRefFluid(fluidRegions.size())
volScalarField & he
Definition: YEEqn.H:50
PtrList< uniformDimensionedVectorField > gFluid(fluidRegions.size())
PtrList< volVectorField > UFluid(fluidRegions.size())
PtrList< fv::options > fluidFvOptions(fluidRegions.size())
PtrList< IOMRFZoneList > MRFfluid(fluidRegions.size())
UniformDimensionedField< vector > uniformDimensionedVectorField
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > linearInterpolate(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: linear.H:108
basicSpecieMixture & composition
PtrList< volScalarField > p_rghFluid(fluidRegions.size())
PtrList< uniformDimensionedScalarField > hRefFluid(fluidRegions.size())
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:58
UniformDimensionedField< scalar > uniformDimensionedScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
List< scalar > initialMassFluid(fluidRegions.size())
PtrList< volScalarField > ghFluid(fluidRegions.size())
PtrList< rhoReactionThermo > thermoFluid(fluidRegions.size())
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
PtrList< volScalarField > KFluid(fluidRegions.size())
PtrList< volScalarField > dpdtFluid(fluidRegions.size())
PtrList< multivariateSurfaceInterpolationScheme< scalar >::fieldTable > fieldsFluid(fluidRegions.size())
const dimensionSet dimPressure
PtrList< fvMesh > fluidRegions(fluidNames.size())
PtrList< compressible::momentumTransportModel > turbulenceFluid(fluidRegions.size())
dimensioned< scalar > magSqr(const dimensioned< Type > &)
PtrList< surfaceScalarField > ghfFluid(fluidRegions.size())
static const char nl
Definition: Ostream.H:260
forAll(fluidRegions, i)
const volScalarField & T
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
PtrList< volScalarField > & Y
PtrList< rhoReactionThermophysicalTransportModel > thermophysicalTransportFluid(fluidRegions.size())
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
PtrList< radiationModel > radiation(fluidRegions.size())
volScalarField & p
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
PtrList< surfaceScalarField > phiFluid(fluidRegions.size())
PtrList< CombustionModel< rhoReactionThermo > > reactionFluid(fluidRegions.size())
PtrList< volScalarField > rhoFluid(fluidRegions.size())