42template<
class BasicTurbulenceModel>
55template<
class BasicTurbulenceModel>
65template<
class BasicTurbulenceModel>
73 min(
max(gammaIntEff_, scalar(0.1)), scalar(1))
78template<
class BasicTurbulenceModel>
108 (1 -
sqr((gammaInt_() - 1.0/ce2_)/(1 - 1.0/ce2_)))
117template<
class BasicTurbulenceModel>
128 this->runTime_.timeName(),
139 const scalar ReThetat = ReThetat_[celli];
146 + 120.656e-4*ReThetat
147 - 868.230e-6*
sqr(ReThetat)
148 + 696.506e-9*
pow3(ReThetat)
149 - 174.105e-12*
pow4(ReThetat)
151 ReThetat - 593.11 - 0.482*(ReThetat - 1870);
158template<
class BasicTurbulenceModel>
171 this->runTime_.timeName(),
185 const scalar ReThetat = ReThetat_[celli];
191 - 119.270e-4*ReThetat
192 - 132.567e-6*
sqr(ReThetat);
194 else if (ReThetat < 596)
198 - 123.939e-2*ReThetat
199 + 194.548e-5*
sqr(ReThetat)
200 - 101.695e-8*
pow3(ReThetat);
202 else if (ReThetat < 1200)
204 Flength[celli] = 0.5 - 3
e-4*(ReThetat - 596);
208 Flength[celli] = 0.3188;
211 const scalar Fsublayer =
214 Flength[celli] = Flength[celli]*(1 - Fsublayer) + 40*Fsublayer;
221template<
class BasicTurbulenceModel>
236 this->runTime_.timeName(),
253 max(100*
sqrt((2.0/3.0)*
k[celli])/
Us[celli], scalar(0.027))
268 const scalar lambda0 =
lambda;
272 const scalar Flambda =
287 (1173.51 - 589.428*Tu + 0.2196/
sqr(Tu))
293 const scalar Flambda =
308 331.50*
pow((Tu - 0.5658), -0.671)
309 *Flambda*
nu[celli]/
Us[celli];
317 maxIter =
max(maxIter, ++iter);
319 }
while (lambdaErr > lambdaErr_);
321 ReThetat0[celli] =
max(thetat*
Us[celli]/
nu[celli], scalar(20));
324 if (maxIter > maxLambdaIter_)
327 <<
"Number of lambda iterations exceeds maxLambdaIter("
328 << maxLambdaIter_ <<
')'<<
endl;
335template<
class BasicTurbulenceModel>
357 max(Fonset2 - Fonset3, scalar(0))
365template<
class BasicTurbulenceModel>
374 const word& propertiesName,
446 this->coeffDict_.template getOrDefault<scalar>(
"lambdaErr", 1
e-6)
450 this->coeffDict_.template getOrDefault<label>(
"maxLambdaIter", 10)
458 IOobject::groupName(
"ReThetat", alphaRhoPhi.group()),
471 IOobject::groupName(
"gammaInt", alphaRhoPhi.group()),
484 IOobject::groupName(
"gammaIntEff", alphaRhoPhi.group()),
492 if (
type == typeName)
494 this->printCoeffs(
type);
501template<
class BasicTurbulenceModel>
506 ca1_.readIfPresent(this->coeffDict());
507 ca2_.readIfPresent(this->coeffDict());
508 ce1_.readIfPresent(this->coeffDict());
509 ce2_.readIfPresent(this->coeffDict());
510 sigmaThetat_.readIfPresent(this->coeffDict());
511 cThetat_.readIfPresent(this->coeffDict());
512 this->coeffDict().readIfPresent(
"lambdaErr", lambdaErr_);
513 this->coeffDict().readIfPresent(
"maxLambdaIter", maxLambdaIter_);
522template<
class BasicTurbulenceModel>
551 alpha()*
rho()*(cThetat_/t)*(1 - Fthetat)
561 Pthetat*ReThetat0(
Us, dUsds,
nu) -
fvm::Sp(Pthetat, ReThetat_)
565 ReThetatEqn.
ref().relax();
580 *ca1_*Flength(
nu)*S*
sqrt(gammaInt_()*Fonset(Rev, ReThetac, RT))
587 alpha()*
rho()*ca2_*Omega*Fturb*gammaInt_()
597 Pgamma -
fvm::Sp(ce1_*Pgamma, gammaInt_)
598 + Egamma -
fvm::Sp(ce2_*Egamma, gammaInt_)
602 gammaIntEqn.
ref().relax();
612 min(2*
max(Rev/(3.235*ReThetac) - 1, scalar(0))*Freattach, scalar(2))
616 gammaIntEff_ =
max(gammaInt_(), gammaSep);
620template<
class BasicTurbulenceModel>
623 if (!this->turbulence_)
632 correctReThetatGammaInt();
Bound the given scalar field if it has gone unbounded.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Langtry-Menter 4-equation transitional SST model based on the k-omega-SST RAS model.
BasicTurbulenceModel::alphaField alphaField
BasicTurbulenceModel::rhoField rhoField
tmp< volScalarField::Internal > Flength(const volScalarField::Internal &nu) const
Empirical correlation that controls the length of the.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField::Internal > Fthetat(const volScalarField::Internal &Us, const volScalarField::Internal &Omega, const volScalarField::Internal &nu) const
Freestream blending-function.
void correctReThetatGammaInt()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField::Internal > ReThetac() const
Empirical correlation for critical Reynolds number where the.
BasicTurbulenceModel::transportModel transportModel
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField &F1, const volTensorField &gradU) const
Modified form of the k-omega SST epsilon/k.
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Modified form of the k-omega SST k production rate.
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 > ReThetat0(const volScalarField::Internal &Us, const volScalarField::Internal &dUsds, const volScalarField::Internal &nu) const
Return the transition onset momentum-thickness Reynolds number.
virtual bool read()
Re-read model coefficients if they have changed.
Implementation of the k-omega-SST turbulence model for incompressible and compressible flows.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField &F1, const volTensorField &gradU) const
Return epsilon/k which for standard RAS is betaStar*omega.
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Return k production rate.
A class for managing temporary objects.
void clear() const noexcept
A class for handling words, derived from Foam::string.
#define WarningInFunction
Report a warning using Foam::Warning.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedScalar exp(const dimensionedScalar &ds)
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensionedScalar pow3(const dimensionedScalar &ds)
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
tensor Ry(const scalar omega)
Rotational transformation tensor about the y-axis by omega radians.
const dimensionSet dimVelocity
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
dimensionedScalar pow4(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
dimensionedTensor skew(const dimensionedTensor &dt)
dimensionedScalar lambda("lambda", dimTime/sqr(dimLength), laminarTransport)
#define forAll(list, i)
Loop across all elements in list.