kOmegaSSTLM.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2016-2021 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "kOmegaSSTLM.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 namespace RASModels
33 {
34 
35 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
36 
37 template<class BasicMomentumTransportModel>
39 (
40  const volScalarField& CDkOmega
41 ) const
42 {
43  const volScalarField Ry(this->y_*sqrt(this->k_)/this->nu());
44  const volScalarField F3(exp(-pow(Ry/120.0, 8)));
45 
46  return max(kOmegaSST<BasicMomentumTransportModel>::F1(CDkOmega), F3);
47 }
48 
49 
50 template<class BasicMomentumTransportModel>
52 (
54 ) const
55 {
56  return gammaIntEff_*kOmegaSST<BasicMomentumTransportModel>::Pk(G);
57 }
58 
59 
60 template<class BasicMomentumTransportModel>
63 (
64  const volScalarField::Internal& F1,
65  const volScalarField::Internal& F2
66 ) const
67 {
68  return
69  min(max(gammaIntEff_, scalar(0.1)), scalar(1))
71 }
72 
73 
74 template<class BasicMomentumTransportModel>
76 (
77  const volScalarField::Internal& Us,
78  const volScalarField::Internal& Omega,
79  const volScalarField::Internal& nu
80 ) const
81 {
82  const volScalarField::Internal& omega = this->omega_();
83  const volScalarField::Internal& y = this->y_();
84 
85  const volScalarField::Internal delta(375*Omega*nu*ReThetat_()*y/sqr(Us));
86  const volScalarField::Internal ReOmega(sqr(y)*omega/nu);
87  const volScalarField::Internal Fwake(exp(-sqr(ReOmega/1e5)));
88 
90  (
91  IOobject::groupName("Fthetat", this->alphaRhoPhi_.group()),
92  min
93  (
94  max
95  (
96  Fwake*exp(-pow4((y/delta))),
97  (1 - sqr((gammaInt_() - 1.0/ce2_)/(1 - 1.0/ce2_)))
98  ),
99  scalar(1)
100  )
101  );
102 }
103 
104 
105 template<class BasicMomentumTransportModel>
108 {
110  (
112  (
113  IOobject::groupName("ReThetac", this->alphaRhoPhi_.group()),
114  this->mesh_,
115  dimless
116  )
117  );
118  volScalarField::Internal& ReThetac = tReThetac.ref();
119 
120  forAll(ReThetac, celli)
121  {
122  const scalar ReThetat = ReThetat_[celli];
123 
124  ReThetac[celli] =
125  ReThetat <= 1870
126  ?
127  ReThetat
128  - 396.035e-2
129  + 120.656e-4*ReThetat
130  - 868.230e-6*sqr(ReThetat)
131  + 696.506e-9*pow3(ReThetat)
132  - 174.105e-12*pow4(ReThetat)
133  :
134  ReThetat - 593.11 - 0.482*(ReThetat - 1870);
135  }
136 
137  return tReThetac;
138 }
139 
140 
141 template<class BasicMomentumTransportModel>
143 (
144  const volScalarField::Internal& nu
145 ) const
146 {
148  (
150  (
151  IOobject::groupName("Flength", this->alphaRhoPhi_.group()),
152  this->mesh_,
153  dimless
154  )
155  );
156  volScalarField::Internal& Flength = tFlength.ref();
157 
158  const volScalarField::Internal& omega = this->omega_();
159  const volScalarField::Internal& y = this->y_();
160 
161  forAll(ReThetat_, celli)
162  {
163  const scalar ReThetat = ReThetat_[celli];
164 
165  if (ReThetat < 400)
166  {
167  Flength[celli] =
168  398.189e-1
169  - 119.270e-4*ReThetat
170  - 132.567e-6*sqr(ReThetat);
171  }
172  else if (ReThetat < 596)
173  {
174  Flength[celli] =
175  263.404
176  - 123.939e-2*ReThetat
177  + 194.548e-5*sqr(ReThetat)
178  - 101.695e-8*pow3(ReThetat);
179  }
180  else if (ReThetat < 1200)
181  {
182  Flength[celli] = 0.5 - 3e-4*(ReThetat - 596);
183  }
184  else
185  {
186  Flength[celli] = 0.3188;
187  }
188 
189  const scalar Fsublayer =
190  exp(-sqr(sqr(y[celli])*omega[celli]/(200*nu[celli])));
191 
192  Flength[celli] = Flength[celli]*(1 - Fsublayer) + 40*Fsublayer;
193  }
194 
195  return tFlength;
196 }
197 
198 
199 template<class BasicMomentumTransportModel>
202 (
203  const volScalarField::Internal& Us,
204  const volScalarField::Internal& dUsds,
205  const volScalarField::Internal& nu
206 ) const
207 {
209  (
211  (
212  IOobject::groupName("ReThetat0", this->alphaRhoPhi_.group()),
213  this->mesh_,
214  dimless
215  )
216  );
217  volScalarField::Internal& ReThetat0 = tReThetat0.ref();
218 
219  const volScalarField& k = this->k_;
220 
221  label maxIter = 0;
222 
223  forAll(ReThetat0, celli)
224  {
225  const scalar Tu
226  (
227  max(100*sqrt((2.0/3.0)*k[celli])/Us[celli], scalar(0.027))
228  );
229 
230  // Initialise lambda to zero.
231  // If lambda were cached between time-steps convergence would be faster
232  // starting from the previous time-step value.
233  scalar lambda = 0;
234 
235  scalar lambdaErr;
236  scalar thetat;
237  label iter = 0;
238 
239  do
240  {
241  // Previous iteration lambda for convergence test
242  const scalar lambda0 = lambda;
243 
244  if (Tu <= 1.3)
245  {
246  const scalar Flambda =
247  dUsds[celli] <= 0
248  ?
249  1
250  - (
251  - 12.986*lambda
252  - 123.66*sqr(lambda)
253  - 405.689*pow3(lambda)
254  )*exp(-pow(Tu/1.5, 1.5))
255  :
256  1
257  + 0.275*(1 - exp(-35*lambda))
258  *exp(-Tu/0.5);
259 
260  thetat =
261  (1173.51 - 589.428*Tu + 0.2196/sqr(Tu))
262  *Flambda*nu[celli]
263  /Us[celli];
264  }
265  else
266  {
267  const scalar Flambda =
268  dUsds[celli] <= 0
269  ?
270  1
271  - (
272  -12.986*lambda
273  -123.66*sqr(lambda)
274  -405.689*pow3(lambda)
275  )*exp(-pow(Tu/1.5, 1.5))
276  :
277  1
278  + 0.275*(1 - exp(-35*lambda))
279  *exp(-2*Tu);
280 
281  thetat =
282  331.50*pow((Tu - 0.5658), -0.671)
283  *Flambda*nu[celli]/Us[celli];
284  }
285 
286  lambda = sqr(thetat)/nu[celli]*dUsds[celli];
287  lambda = max(min(lambda, 0.1), -0.1);
288 
289  lambdaErr = mag(lambda - lambda0);
290 
291  maxIter = max(maxIter, ++iter);
292 
293  } while (lambdaErr > lambdaErr_);
294 
295  ReThetat0[celli] = max(thetat*Us[celli]/nu[celli], scalar(20));
296  }
297 
298  if (maxIter > maxLambdaIter_)
299  {
301  << "Number of lambda iterations exceeds maxLambdaIter("
302  << maxLambdaIter_ << ')'<< endl;
303  }
304 
305  return tReThetat0;
306 }
307 
308 
309 template<class BasicMomentumTransportModel>
311 (
312  const volScalarField::Internal& Rev,
313  const volScalarField::Internal& ReThetac,
314  const volScalarField::Internal& RT
315 ) const
316 {
317  const volScalarField::Internal Fonset1(Rev/(2.193*ReThetac));
318 
319  const volScalarField::Internal Fonset2
320  (
321  min(max(Fonset1, pow4(Fonset1)), scalar(2))
322  );
323 
324  const volScalarField::Internal Fonset3(max(1 - pow3(RT/2.5), scalar(0)));
325 
327  (
328  IOobject::groupName("Fonset", this->alphaRhoPhi_.group()),
329  max(Fonset2 - Fonset3, scalar(0))
330  );
331 }
332 
333 
334 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
335 
336 template<class BasicMomentumTransportModel>
338 (
339  const alphaField& alpha,
340  const rhoField& rho,
341  const volVectorField& U,
342  const surfaceScalarField& alphaRhoPhi,
343  const surfaceScalarField& phi,
344  const transportModel& transport,
345  const word& type
346 )
347 :
349  (
350  alpha,
351  rho,
352  U,
353  alphaRhoPhi,
354  phi,
355  transport
356  ),
357 
358  ca1_
359  (
361  (
362  "ca1",
363  this->coeffDict_,
364  2
365  )
366  ),
367  ca2_
368  (
370  (
371  "ca2",
372  this->coeffDict_,
373  0.06
374  )
375  ),
376  ce1_
377  (
379  (
380  "ce1",
381  this->coeffDict_,
382  1
383  )
384  ),
385  ce2_
386  (
388  (
389  "ce2",
390  this->coeffDict_,
391  50
392  )
393  ),
394  cThetat_
395  (
397  (
398  "cThetat",
399  this->coeffDict_,
400  0.03
401  )
402  ),
403  sigmaThetat_
404  (
406  (
407  "sigmaThetat",
408  this->coeffDict_,
409  2
410  )
411  ),
412  lambdaErr_
413  (
414  this->coeffDict_.lookupOrDefault("lambdaErr", 1e-6)
415  ),
416  maxLambdaIter_
417  (
418  this->coeffDict_.lookupOrDefault("maxLambdaIter", 10)
419  ),
420  deltaU_("deltaU", dimVelocity, small),
421 
422  ReThetat_
423  (
424  IOobject
425  (
426  IOobject::groupName("ReThetat", alphaRhoPhi.group()),
427  this->runTime_.timeName(),
428  this->mesh_,
431  ),
432  this->mesh_
433  ),
434 
435  gammaInt_
436  (
437  IOobject
438  (
439  IOobject::groupName("gammaInt", alphaRhoPhi.group()),
440  this->runTime_.timeName(),
441  this->mesh_,
444  ),
445  this->mesh_
446  ),
447 
448  gammaIntEff_
449  (
450  IOobject
451  (
452  IOobject::groupName("gammaIntEff", alphaRhoPhi.group()),
453  this->runTime_.timeName(),
454  this->mesh_
455  ),
456  this->mesh_,
458  )
459 {}
460 
461 
462 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
463 
464 template<class BasicMomentumTransportModel>
466 {
468  {
469  ca1_.readIfPresent(this->coeffDict());
470  ca2_.readIfPresent(this->coeffDict());
471  ce1_.readIfPresent(this->coeffDict());
472  ce2_.readIfPresent(this->coeffDict());
473  sigmaThetat_.readIfPresent(this->coeffDict());
474  cThetat_.readIfPresent(this->coeffDict());
475  this->coeffDict().readIfPresent("lambdaErr", lambdaErr_);
476  this->coeffDict().readIfPresent("maxLambdaIter", maxLambdaIter_);
477 
478  return true;
479  }
480  else
481  {
482  return false;
483  }
484 }
485 
486 
487 template<class BasicMomentumTransportModel>
489 {
490  // Local references
491  const alphaField& alpha = this->alpha_;
492  const rhoField& rho = this->rho_;
493  const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_;
494  const volVectorField& U = this->U_;
495  const volScalarField& k = this->k_;
496  const volScalarField& omega = this->omega_;
497  const tmp<volScalarField> tnu = this->nu();
498  const volScalarField::Internal& nu = tnu()();
499  const volScalarField::Internal& y = this->y_();
500  const Foam::fvModels& fvModels(Foam::fvModels::New(this->mesh_));
502  (
503  Foam::fvConstraints::New(this->mesh_)
504  );
505 
506  // Fields derived from the velocity gradient
507  tmp<volTensorField> tgradU = fvc::grad(U);
508  const volScalarField::Internal Omega(sqrt(2*magSqr(skew(tgradU()()))));
509  const volScalarField::Internal S(sqrt(2*magSqr(symm(tgradU()()))));
510  const volScalarField::Internal Us(max(mag(U()), deltaU_));
511  const volScalarField::Internal dUsds((U() & (U() & tgradU()()))/sqr(Us));
512  tgradU.clear();
513 
514  const volScalarField::Internal Fthetat(this->Fthetat(Us, Omega, nu));
515 
516  {
517  const volScalarField::Internal t(500*nu/sqr(Us));
518  const volScalarField::Internal Pthetat
519  (
520  alpha()*rho()*(cThetat_/t)*(1 - Fthetat)
521  );
522 
523  // Transition onset momentum-thickness Reynolds number equation
524  tmp<fvScalarMatrix> ReThetatEqn
525  (
526  fvm::ddt(alpha, rho, ReThetat_)
527  + fvm::div(alphaRhoPhi, ReThetat_)
528  - fvm::laplacian(alpha*rho*DReThetatEff(), ReThetat_)
529  ==
530  Pthetat*ReThetat0(Us, dUsds, nu) - fvm::Sp(Pthetat, ReThetat_)
531  + fvModels.source(alpha, rho, ReThetat_)
532  );
533 
534  ReThetatEqn.ref().relax();
535  fvConstraints.constrain(ReThetatEqn.ref());
536  solve(ReThetatEqn);
537  fvConstraints.constrain(ReThetat_);
538  bound(ReThetat_, 0);
539  }
540 
541  const volScalarField::Internal ReThetac(this->ReThetac());
542  const volScalarField::Internal Rev(sqr(y)*S/nu);
543  const volScalarField::Internal RT(k()/(nu*omega()));
544 
545  {
546  const volScalarField::Internal Pgamma
547  (
548  alpha()*rho()
549  *ca1_*Flength(nu)*S*sqrt(gammaInt_()*Fonset(Rev, ReThetac, RT))
550  );
551 
552  const volScalarField::Internal Fturb(exp(-pow4(0.25*RT)));
553 
554  const volScalarField::Internal Egamma
555  (
556  alpha()*rho()*ca2_*Omega*Fturb*gammaInt_()
557  );
558 
559  // Intermittency equation
560  tmp<fvScalarMatrix> gammaIntEqn
561  (
562  fvm::ddt(alpha, rho, gammaInt_)
563  + fvm::div(alphaRhoPhi, gammaInt_)
564  - fvm::laplacian(alpha*rho*DgammaIntEff(), gammaInt_)
565  ==
566  Pgamma - fvm::Sp(ce1_*Pgamma, gammaInt_)
567  + Egamma - fvm::Sp(ce2_*Egamma, gammaInt_)
568  + fvModels.source(alpha, rho, gammaInt_)
569  );
570 
571  gammaIntEqn.ref().relax();
572  fvConstraints.constrain(gammaIntEqn.ref());
573  solve(gammaIntEqn);
574  fvConstraints.constrain(gammaInt_);
575  bound(gammaInt_, 0);
576  }
577 
578  const volScalarField::Internal Freattach(exp(-pow4(RT/20.0)));
579  const volScalarField::Internal gammaSep
580  (
581  min(2*max(Rev/(3.235*ReThetac) - 1, scalar(0))*Freattach, scalar(2))
582  *Fthetat
583  );
584 
585  gammaIntEff_ = max(gammaInt_(), gammaSep);
586 }
587 
588 
589 template<class BasicMomentumTransportModel>
591 {
592  if (!this->turbulence_)
593  {
594  return;
595  }
596 
597  // Correct ReThetat and gammaInt
598  correctReThetatGammaInt();
599 
600  // Correct k and omega
602 }
603 
604 
605 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
606 
607 } // End namespace RASModels
608 } // End namespace Foam
609 
610 // ************************************************************************* //
scalar delta
BasicMomentumTransportModel::alphaField alphaField
Definition: kOmegaSST.H:68
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
static dimensioned< Type > lookupOrAddToDict(const word &, dictionary &, const dimensionSet &dims=dimless, const Type &defaultValue=pTraits< Type >::zero)
Construct from dictionary, with default value.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:237
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &)
Return a temporary field constructed from name, mesh.
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField::Internal > ReThetat0(const volScalarField::Internal &Us, const volScalarField::Internal &dUsds, const volScalarField::Internal &nu) const
Return the transition onset momentum-thickness Reynolds number.
Definition: kOmegaSSTLM.C:202
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedTensor skew(const dimensionedTensor &dt)
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
dimensionedSymmTensor sqr(const dimensionedVector &dv)
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
BasicMomentumTransportModel::transportModel transportModel
Definition: kOmegaSST.H:70
const dimensionSet dimless
label k
Boltzmann constant.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField::Internal &F1, const volScalarField::Internal &F2) const
Modified form of the k-omega SST epsilon/k.
Definition: kOmegaSSTLM.C:63
kOmegaSSTLM(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &type=typeName)
Construct from components.
Definition: kOmegaSSTLM.C:338
tmp< volScalarField::Internal > ReThetac() const
Empirical correlation for critical Reynolds number where the.
Definition: kOmegaSSTLM.C:107
scalar y
dimensionedScalar exp(const dimensionedScalar &ds)
Foam::fvConstraints & fvConstraints
#define F3(B, C, D)
Definition: SHA1.C:175
A class for handling words, derived from string.
Definition: word.H:59
tmp< fvMatrix< Type > > S(const Pair< tmp< volScalarField::Internal >> &, const GeometricField< Type, fvPatchField, volMesh > &)
static word groupName(Name name, const word &group)
Finite volume models.
Definition: fvModels.H:60
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kOmegaSSTLM.C:590
phi
Definition: correctPhi.H:3
Specialisation for RAS of the generic kOmegaSSTBase base class. For more information, see Description of kOmegaSSTBase.H.
Definition: kOmegaSST.H:57
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual bool read()
Re-read model coefficients if they have changed.
Definition: kOmegaSSTLM.C:465
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
const dimensionSet dimVelocity
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:96
tmp< volScalarField::Internal > Fonset(const volScalarField::Internal &Rev, const volScalarField::Internal &ReThetac, const volScalarField::Internal &RT) const
Transition onset location control function.
Definition: kOmegaSSTLM.C:311
tmp< volScalarField::Internal > Flength(const volScalarField::Internal &nu) const
Empirical correlation that controls the length of the.
Definition: kOmegaSSTLM.C:143
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Return k production rate.
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
tensor Ry(const scalar &omega)
Rotational transformation tensor about the y-axis by omega radians.
Definition: transform.H:106
Foam::fvModels & fvModels
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
U
Definition: pEqn.H:72
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:33
#define WarningInFunction
Report a warning using Foam::Warning.
Finite volume constraints.
Definition: fvConstraints.H:57
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
dimensionedScalar pow4(const dimensionedScalar &ds)
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > mag(const dimensioned< Type > &)
void correctReThetatGammaInt()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kOmegaSSTLM.C:488
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField::Internal &F1, const volScalarField::Internal &F2) const
Return epsilon/k which for standard RAS is betaStar*omega.
tmp< volScalarField::Internal > Fthetat(const volScalarField::Internal &Us, const volScalarField::Internal &Omega, const volScalarField::Internal &nu) const
Freestream blending-function.
Definition: kOmegaSSTLM.C:76
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
virtual tmp< volScalarField > F1(const volScalarField &CDkOmega) const
Modified form of the k-omega SST F1 function.
Definition: kOmegaSSTLM.C:39
Namespace for OpenFOAM.
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Modified form of the k-omega SST k production rate.
Definition: kOmegaSSTLM.C:52
BasicMomentumTransportModel::rhoField rhoField
Definition: kOmegaSST.H:69