createFluidFields.H 7.13 KB
Newer Older
1
// Initialise fluid field pointer lists
2
PtrList<rhoReactionThermo> thermoFluid(fluidRegions.size());
3
4
5
PtrList<volScalarField> rhoFluid(fluidRegions.size());
PtrList<volVectorField> UFluid(fluidRegions.size());
PtrList<surfaceScalarField> phiFluid(fluidRegions.size());
6
PtrList<uniformDimensionedScalarField> hRefFluid(fluidRegions.size());
7
8
PtrList<volScalarField> ghFluid(fluidRegions.size());
PtrList<surfaceScalarField> ghfFluid(fluidRegions.size());
9
10
PtrList<compressible::turbulenceModel> turbulenceFluid(fluidRegions.size());
PtrList<CombustionModel<rhoReactionThermo>> reactionFluid(fluidRegions.size());
11
PtrList<volScalarField> p_rghFluid(fluidRegions.size());
12
13
14
PtrList<radiation::radiationModel> radiation(fluidRegions.size());
PtrList<volScalarField> KFluid(fluidRegions.size());
PtrList<volScalarField> dpdtFluid(fluidRegions.size());
15
16
17
PtrList<multivariateSurfaceInterpolationScheme<scalar>::fieldTable>
    fieldsFluid(fluidRegions.size());
PtrList<volScalarField> QdotFluid(fluidRegions.size());
18
19
20
21
22

List<scalar> initialMassFluid(fluidRegions.size());
List<bool> frozenFlowFluid(fluidRegions.size(), false);

PtrList<IOMRFZoneList> MRFfluid(fluidRegions.size());
23
PtrList<fv::options> fluidFvOptions(fluidRegions.size());
24

25
26
List<label> pRefCellFluid(fluidRegions.size());
List<scalar> pRefValueFluid(fluidRegions.size());
27

28
29
const uniformDimensionedVectorField& g = meshObjects::gravity::New(runTime);

30
31
32
33
34
35
36
// Populate fluid field pointer lists
forAll(fluidRegions, i)
{
    Info<< "*** Reading fluid mesh thermophysical properties for region "
        << fluidRegions[i].name() << nl << endl;

    Info<< "    Adding to thermoFluid\n" << endl;
37
    thermoFluid.set(i, rhoReactionThermo::New(fluidRegions[i]).ptr());
38
39
40
41
42
43

    Info<< "    Adding to rhoFluid\n" << endl;
    rhoFluid.set
    (
        i,
        new volScalarField
Henry's avatar
Henry committed
44
        (
45
            IOobject
Henry's avatar
Henry committed
46
            (
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
                "rho",
                runTime.timeName(),
                fluidRegions[i],
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
            ),
            thermoFluid[i].rho()
        )
    );

    Info<< "    Adding to UFluid\n" << endl;
    UFluid.set
    (
        i,
        new volVectorField
Henry's avatar
Henry committed
62
        (
63
            IOobject
Henry's avatar
Henry committed
64
            (
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
                "U",
                runTime.timeName(),
                fluidRegions[i],
                IOobject::MUST_READ,
                IOobject::AUTO_WRITE
            ),
            fluidRegions[i]
        )
    );

    Info<< "    Adding to phiFluid\n" << endl;
    phiFluid.set
    (
        i,
        new surfaceScalarField
Henry's avatar
Henry committed
80
        (
81
            IOobject
Henry's avatar
Henry committed
82
            (
83
84
85
86
87
88
89
90
91
92
93
                "phi",
                runTime.timeName(),
                fluidRegions[i],
                IOobject::READ_IF_PRESENT,
                IOobject::AUTO_WRITE
            ),
            linearInterpolate(rhoFluid[i]*UFluid[i])
                & fluidRegions[i].Sf()
        )
    );

94
95
    Info<< "    Adding to hRefFluid\n" << endl;
    hRefFluid.set
96
97
    (
        i,
98
        new uniformDimensionedScalarField
Henry's avatar
Henry committed
99
        (
100
101
102
103
104
105
106
107
            IOobject
            (
                "hRef",
                runTime.constant(),
                fluidRegions[i],
                IOobject::READ_IF_PRESENT,
                IOobject::NO_WRITE
            ),
108
            dimensionedScalar("hRef", dimLength, Zero) // uses name
109
110
111
112
113
        )
    );

    dimensionedScalar ghRef
    (
114
115
116
        mag(g.value()) > SMALL
      ? g & (cmptMag(g.value())/mag(g.value()))*hRefFluid[i]
      : dimensionedScalar("ghRef", g.dimensions()*dimLength, 0)
117
118
119
120
121
122
    );

    Info<< "    Adding to ghFluid\n" << endl;
    ghFluid.set
    (
        i,
123
124
125
        new volScalarField
        (
            "gh",
126
            (g & fluidRegions[i].C()) - ghRef
127
        )
128
129
130
131
132
133
    );

    Info<< "    Adding to ghfFluid\n" << endl;
    ghfFluid.set
    (
        i,
134
135
136
        new surfaceScalarField
        (
            "ghf",
137
            (g & fluidRegions[i].Cf()) - ghRef
138
139
140
        )
    );

141
142
    Info<< "    Adding to turbulenceFluid\n" << endl;
    turbulenceFluid.set
143
144
145
146
147
148
149
150
151
    (
        i,
        compressible::turbulenceModel::New
        (
            rhoFluid[i],
            UFluid[i],
            phiFluid[i],
            thermoFluid[i]
        ).ptr()
152
153
    );

154
155
156
157
158
159
160
161
162
163
164
    Info<< "    Adding to reactionFluid\n" << endl;
    reactionFluid.set
    (
        i,
        CombustionModel<rhoReactionThermo>::New
        (
            thermoFluid[i],
            turbulenceFluid[i]
        )
    );

165
166
167
168
    p_rghFluid.set
    (
        i,
        new volScalarField
Henry's avatar
Henry committed
169
        (
170
            IOobject
Henry's avatar
Henry committed
171
            (
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
                "p_rgh",
                runTime.timeName(),
                fluidRegions[i],
                IOobject::MUST_READ,
                IOobject::AUTO_WRITE
            ),
            fluidRegions[i]
        )
    );

    // Force p_rgh to be consistent with p
    p_rghFluid[i] = thermoFluid[i].p() - rhoFluid[i]*ghFluid[i];

    fluidRegions[i].setFluxRequired(p_rghFluid[i].name());

187
    Info<< "    Adding to radiationFluid\n" << endl;
188
189
190
191
192
193
194
195
196
197
198
199
200
    radiation.set
    (
        i,
        radiation::radiationModel::New(thermoFluid[i].T())
    );

    initialMassFluid[i] = fvc::domainIntegrate(rhoFluid[i]).value();

    Info<< "    Adding to KFluid\n" << endl;
    KFluid.set
    (
        i,
        new volScalarField
Henry's avatar
Henry committed
201
        (
202
203
204
205
206
207
208
209
210
211
            "K",
            0.5*magSqr(UFluid[i])
        )
    );

    Info<< "    Adding to dpdtFluid\n" << endl;
    dpdtFluid.set
    (
        i,
        new volScalarField
Henry's avatar
Henry committed
212
        (
213
            IOobject
Henry's avatar
Henry committed
214
            (
215
216
217
218
219
                "dpdt",
                runTime.timeName(),
                fluidRegions[i]
            ),
            fluidRegions[i],
220
            dimensionedScalar(thermoFluid[i].p().dimensions()/dimTime, Zero)
221
222
223
        )
    );

224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
    Info<< "    Adding to fieldsFluid\n" << endl;
    fieldsFluid.set
    (
        i,
        new multivariateSurfaceInterpolationScheme<scalar>::fieldTable
    );
    forAll(thermoFluid[i].composition().Y(), j)
    {
        fieldsFluid[i].add(thermoFluid[i].composition().Y()[j]);
    }
    fieldsFluid[i].add(thermoFluid[i].he());

    Info<< "    Adding to QdotFluid\n" << endl;
    QdotFluid.set
    (
        i,
        new volScalarField
        (
            IOobject
            (
                "Qdot",
                runTime.timeName(),
                fluidRegions[i],
                IOobject::READ_IF_PRESENT,
                IOobject::AUTO_WRITE
            ),
            fluidRegions[i],
251
            dimensionedScalar(dimEnergy/dimVolume/dimTime, Zero)
252
253
254
        )
    );

255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
    const dictionary& pimpleDict =
        fluidRegions[i].solutionDict().subDict("PIMPLE");
    pimpleDict.readIfPresent("frozenFlow", frozenFlowFluid[i]);

    Info<< "    Adding MRF\n" << endl;
    MRFfluid.set
    (
        i,
        new IOMRFZoneList(fluidRegions[i])
    );

    Info<< "    Adding fvOptions\n" << endl;
    fluidFvOptions.set
    (
        i,
270
        new fv::options(fluidRegions[i])
271
    );
272

273
    turbulenceFluid[i].validate();
274

275
    pRefCellFluid[i] = -1;
276
    pRefValueFluid[i] = 0.0;
277
278
279
280
281
282
283
284

    if (p_rghFluid[i].needReference())
    {
        setRefCell
        (
            thermoFluid[i].p(),
            p_rghFluid[i],
            pimpleDict,
285
286
            pRefCellFluid[i],
            pRefValueFluid[i]
287
288
289
        );
    }

290
}