Go to the documentation of this file.
42 template<
class BasicTurbulenceModel>
55 template<
class BasicTurbulenceModel>
65 template<
class BasicTurbulenceModel>
73 min(
max(gammaIntEff_, scalar(0.1)), scalar(1))
78 template<
class BasicTurbulenceModel>
108 (1 -
sqr((gammaInt_() - 1.0/ce2_)/(1 - 1.0/ce2_)))
117 template<
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);
158 template<
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;
221 template<
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;
335 template<
class BasicTurbulenceModel>
357 max(Fonset2 - Fonset3, scalar(0))
365 template<
class BasicTurbulenceModel>
374 const word& propertiesName,
446 this->coeffDict_.template getOrDefault<scalar>(
"lambdaErr", 1
e-6)
450 this->coeffDict_.template getOrDefault<label>(
"maxLambdaIter", 10)
459 this->runTime_.timeName(),
472 this->runTime_.timeName(),
485 this->runTime_.timeName(),
492 if (
type == typeName)
494 this->printCoeffs(
type);
501 template<
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_);
522 template<
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);
620 template<
class BasicTurbulenceModel>
623 if (!this->turbulence_)
632 correctReThetatGammaInt();
Defines the attributes of an object for which implicit objectRegistry management is supported,...
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Modified form of the k-omega SST k production rate.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
A class for handling words, derived from Foam::string.
void clear() const noexcept
dimensionedTensor skew(const dimensionedTensor &dt)
Implementation of the k-omega-SST turbulence model for incompressible and compressible flows.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
const dimensionedScalar G
Newtonian constant of gravitation.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimVelocity
tmp< volScalarField::Internal > Fthetat(const volScalarField::Internal &Us, const volScalarField::Internal &Omega, const volScalarField::Internal &nu) const
Freestream blending-function.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
void correctReThetatGammaInt()
Solve the turbulence equations and correct the turbulence viscosity.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField &F1, const volTensorField &gradU) const
Return epsilon/k which for standard RAS is betaStar*omega.
tensor Ry(const scalar &omega)
Rotational transformation tensor about the y-axis by omega radians.
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.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Return k production rate.
dimensionedScalar exp(const dimensionedScalar &ds)
BasicTurbulenceModel::alphaField alphaField
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
BasicTurbulenceModel::rhoField rhoField
dimensionedScalar pow4(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
dimensionedScalar pow3(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
Bound the given scalar field if it has gone unbounded.
Langtry-Menter 4-equation transitional SST model based on the k-omega-SST RAS model.
virtual tmp< volScalarField > F1(const volScalarField &CDkOmega) const
Modified form of the k-omega SST F1 function.
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedScalar lambda("lambda", dimTime/sqr(dimLength), laminarTransport)
static dimensioned< Type > getOrAddToDict(const word &name, dictionary &dict, const dimensionSet &dims=dimless, const Type &deflt=Type(Zero))
Construct dimensioned from dictionary, with default value.
BasicTurbulenceModel::transportModel transportModel
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField::Internal > ReThetac() const
Empirical correlation for critical Reynolds number where the.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField &F1, const volTensorField &gradU) const
Modified form of the k-omega SST epsilon/k.
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.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
label k
Boltzmann constant.
const dimensionedScalar e
Elementary charge.
virtual bool read()
Re-read model coefficients if they have changed.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
#define WarningInFunction
Report a warning using Foam::Warning.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
const dimensionSet dimless
Dimensionless.