37 template<
class BasicTurbulenceModel>
50 template<
class BasicTurbulenceModel>
60 template<
class BasicTurbulenceModel>
68 min(
max(gammaIntEff_, scalar(0.1)), scalar(1))
73 template<
class BasicTurbulenceModel>
98 (1 -
sqr((gammaInt_() - 1.0/ce2_)/(1 - 1.0/ce2_)))
107 template<
class BasicTurbulenceModel>
118 this->runTime_.timeName(),
129 const scalar ReThetat = ReThetat_[celli];
136 + 120.656e-4*ReThetat
137 - 868.230e-6*
sqr(ReThetat)
138 + 696.506e-9*
pow3(ReThetat)
139 - 174.105e-12*
pow4(ReThetat)
141 ReThetat - 593.11 - 0.482*(ReThetat - 1870);
148 template<
class BasicTurbulenceModel>
161 this->runTime_.timeName(),
175 const scalar ReThetat = ReThetat_[celli];
181 - 119.270e-4*ReThetat
182 - 132.567e-6*
sqr(ReThetat);
184 else if (ReThetat < 596)
188 - 123.939e-2*ReThetat
189 + 194.548e-5*
sqr(ReThetat)
190 - 101.695e-8*
pow3(ReThetat);
192 else if (ReThetat < 1200)
194 Flength[celli] = 0.5 - 3
e-4*(ReThetat - 596);
198 Flength[celli] = 0.3188;
201 const scalar Fsublayer =
202 exp(-
sqr(
sqr(y[celli])*omega[celli]/(200*nu[celli])));
204 Flength[celli] = Flength[celli]*(1 - Fsublayer) + 40*Fsublayer;
211 template<
class BasicTurbulenceModel>
226 this->runTime_.timeName(),
243 max(100*
sqrt((2.0/3.0)*k[celli])/Us[celli], scalar(0.027))
258 const scalar lambda0 =
lambda;
262 const scalar Flambda =
269 - 405.689*
pow3(lambda)
273 + 0.275*(1 -
exp(-35*lambda))
277 (1173.51 - 589.428*Tu + 0.2196/
sqr(Tu))
283 const scalar Flambda =
290 -405.689*
pow3(lambda)
294 + 0.275*(1 -
exp(-35*lambda))
298 331.50*
pow((Tu - 0.5658), -0.671)
299 *Flambda*nu[celli]/Us[celli];
302 lambda =
sqr(thetat)/nu[celli]*dUsds[celli];
303 lambda =
max(
min(lambda, 0.1), -0.1);
305 lambdaErr =
mag(lambda - lambda0);
307 maxIter =
max(maxIter, ++iter);
309 }
while (lambdaErr > lambdaErr_);
311 ReThetat0[celli] =
max(thetat*Us[celli]/nu[celli], scalar(20));
314 if (maxIter > maxLambdaIter_)
317 <<
"Number of lambda iterations exceeds maxLambdaIter(" 318 << maxLambdaIter_ <<
')'<<
endl;
325 template<
class BasicTurbulenceModel>
347 max(Fonset2 - Fonset3, scalar(0))
355 template<
class BasicTurbulenceModel>
364 const word& propertiesName,
435 this->coeffDict_.lookupOrDefault(
"lambdaErr", 1
e-6)
439 this->coeffDict_.lookupOrDefault(
"maxLambdaIter", 10)
448 this->runTime_.timeName(),
461 this->runTime_.timeName(),
474 this->runTime_.timeName(),
485 template<
class BasicTurbulenceModel>
490 ca1_.readIfPresent(this->coeffDict());
491 ca2_.readIfPresent(this->coeffDict());
492 ce1_.readIfPresent(this->coeffDict());
493 ce2_.readIfPresent(this->coeffDict());
494 sigmaThetat_.readIfPresent(this->coeffDict());
495 cThetat_.readIfPresent(this->coeffDict());
496 this->coeffDict().readIfPresent(
"lambdaErr", lambdaErr_);
497 this->coeffDict().readIfPresent(
"maxLambdaIter", maxLambdaIter_);
508 template<
class BasicTurbulenceModel>
537 alpha()*
rho()*(cThetat_/t)*(1 - Fthetat)
547 Pthetat*ReThetat0(Us, dUsds, nu) -
fvm::Sp(Pthetat, ReThetat_)
551 ReThetatEqn.
ref().relax();
566 *ca1_*Flength(nu)*S*
sqrt(gammaInt_()*Fonset(Rev, ReThetac, RT))
573 alpha()*
rho()*ca2_*Omega*Fturb*gammaInt_()
583 Pgamma -
fvm::Sp(ce1_*Pgamma, gammaInt_)
584 + Egamma -
fvm::Sp(ce2_*Egamma, gammaInt_)
588 gammaIntEqn.
ref().relax();
598 min(2*
max(Rev/(3.235*ReThetac) - 1, scalar(0))*Freattach, scalar(2))
602 gammaIntEff_ =
max(gammaInt_(), gammaSep);
606 template<
class BasicTurbulenceModel>
609 if (!this->turbulence_)
615 correctReThetatGammaInt();
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
BasicTurbulenceModel::transportModel transportModel
#define forAll(list, i)
Loop across all elements in list.
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.
void clear() const
If object pointer points to valid object:
BasicTurbulenceModel::rhoField rhoField
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.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedTensor skew(const dimensionedTensor &dt)
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
T & ref() const
Return non-const reference or generate a fatal error.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
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.
tmp< volScalarField::Internal > ReThetac() const
Empirical correlation for critical Reynolds number where the.
dimensionedScalar exp(const dimensionedScalar &ds)
A class for handling words, derived from string.
SolverPerformance< Type > solve(fvMatrix< Type > &, const dictionary &)
Solve returning the solution statistics given convergence tolerance.
static word groupName(Name name, const word &group)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Specialisation for RAS of the generic kOmegaSSTBase base class. For more information, see Description of kOmegaSSTBase.H.
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual bool read()
Re-read model coefficients if they have changed.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< volScalarField::Internal > Fonset(const volScalarField::Internal &Rev, const volScalarField::Internal &ReThetac, const volScalarField::Internal &RT) const
Transition onset location control function.
tmp< volScalarField::Internal > Flength(const volScalarField::Internal &nu) const
Empirical correlation that controls the length of the.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Return k production rate.
tensor Ry(const scalar &omega)
Rotational transformation tensor about the y-axis by omega radians.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
#define WarningInFunction
Report a warning using Foam::Warning.
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)
Langtry-Menter 4-equation transitional SST model based on the k-omega-SST RAS model.
dimensionedScalar pow4(const dimensionedScalar &ds)
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
dimensioned< scalar > mag(const dimensioned< Type > &)
void correctReThetatGammaInt()
Solve the turbulence equations and correct the turbulence viscosity.
const doubleScalar e
Elementary charge.
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.
A class for managing temporary objects.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
BasicTurbulenceModel::alphaField alphaField
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
virtual tmp< volScalarField > F1(const volScalarField &CDkOmega) const
Modified form of the k-omega SST F1 function.
const dimensionSet dimVelocity
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Modified form of the k-omega SST k production rate.