41template<
class BasicTurbulenceModel>
44 return nuTilda_/this->
nu();
48template<
class BasicTurbulenceModel>
55 return chi3/(chi3 +
pow3(Cv1_));
59template<
class BasicTurbulenceModel>
66 return 1.0 - chi/(1.0 + chi*fv1);
70template<
class BasicTurbulenceModel>
76 return Ct3_*
exp(-Ct4_*
sqr(chi));
80template<
class BasicTurbulenceModel>
90template<
class BasicTurbulenceModel>
104 + fv2(chi, fv1)*nuTilda_/
sqr(kappa_*dTilda),
111template<
class BasicTurbulenceModel>
135 tr.ref().boundaryFieldRef() == 0.0;
141template<
class BasicTurbulenceModel>
155template<
class BasicTurbulenceModel>
179 if (lowReCorrection_)
192 (1 - Cb1_/(Cw1_*
sqr(kappa_)*fwStar_)*(ft2 + (1 - ft2)*fv2))
193 /
max(SMALL, (fv1*
max(scalar(1
e-10), 1 - ft2)))
202template<
class BasicTurbulenceModel>
211 min(tdTilda.
ref().ref(), tdTilda(), y_);
216template<
class BasicTurbulenceModel>
222 this->nut_ = nuTilda_*fv1;
226 BasicTurbulenceModel::correctNut();
230template<
class BasicTurbulenceModel>
233 correctNut(fv1(this->chi()));
239template<
class BasicTurbulenceModel>
248 const word& propertiesName,
300 Cw1_(Cb1_/
sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_),
407 if (
type == typeName)
409 this->printCoeffs(
type);
416template<
class BasicTurbulenceModel>
421 sigmaNut_.readIfPresent(this->coeffDict());
422 kappa_.readIfPresent(*
this);
424 Cb1_.readIfPresent(this->coeffDict());
425 Cb2_.readIfPresent(this->coeffDict());
426 Cw1_ = Cb1_/
sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
428 Cw3_.readIfPresent(this->coeffDict());
429 Cv1_.readIfPresent(this->coeffDict());
430 Cs_.readIfPresent(this->coeffDict());
432 CDES_.readIfPresent(this->coeffDict());
433 ck_.readIfPresent(this->coeffDict());
435 lowReCorrection_.readIfPresent(
"lowReCorrection", this->coeffDict());
436 Ct3_.readIfPresent(this->coeffDict());
437 Ct4_.readIfPresent(this->coeffDict());
438 fwStar_.readIfPresent(this->coeffDict());
447template<
class BasicTurbulenceModel>
458template<
class BasicTurbulenceModel>
471 this->mesh_.time().timeName(),
482 dTildaF = dTilda(chi, fv1,
fvc::grad(this->U_));
485 return sqr(this->
nut()/ck_/dTildaF);
489template<
class BasicTurbulenceModel>
502 this->mesh_.time().timeName(),
515template<
class BasicTurbulenceModel>
518 if (!this->turbulence_)
539 const volScalarField Stilda(this->Stilda(chi, fv1, Omega, dTilda));
548 Cb1_*
alpha*
rho*Stilda*nuTilda_*(scalar(1) - ft2)
551 (Cw1_*fw(Stilda, dTilda) - Cb1_/
sqr(kappa_)*ft2)
558 nuTildaEqn.
ref().relax();
563 nuTilda_.correctBoundaryConditions();
const uniformDimensionedVectorField & g
const dimensionSet & dimensions() const
Return dimensions.
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Templated abstract base class for DES models.
SpalartAllmarasDES DES turbulence model for incompressible and compressible flows.
BasicTurbulenceModel::alphaField alphaField
tmp< volScalarField > chi() const
BasicTurbulenceModel::rhoField rhoField
tmp< volScalarField > fw(const volScalarField &Omega, const volScalarField &dTilda) const
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.
virtual void correct()
Correct nuTilda and related properties.
tmp< volScalarField > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
tmp< volScalarField > fv1(const volScalarField &chi) const
virtual tmp< volScalarField > LESRegion() const
Return the LES field indicator.
virtual void correctNut()
BasicTurbulenceModel::transportModel transportModel
tmp< volScalarField > ft2(const volScalarField &chi) const
tmp< volScalarField > fv2(const volScalarField &chi, const volScalarField &fv1) const
virtual tmp< volScalarField > dTilda(const volScalarField &chi, const volScalarField &fv1, const volTensorField &gradU) const
Length scale.
virtual bool read()
Re-read model coefficients if they have changed.
vector Omega() const
Return the current Omega vector.
tmp< volScalarField::Internal > Stilda() const
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Generic dimensioned Type class.
bool readIfPresent(const dictionary &dict)
const quaternion & r() const
A class for managing temporary objects.
Interface to run-time selectable methods to calculate the distance-to-wall and normal-to-wall fields.
A class for handling words, derived from Foam::string.
This boundary condition applies a zero-gradient condition from the patch internal field onto the patc...
const volScalarField & psi
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.
dimensionedScalar pow6(const dimensionedScalar &ds)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedScalar exp(const dimensionedScalar &ds)
const dimensionSet dimless
Dimensionless.
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
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.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
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)
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 neg(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
dimensionedTensor skew(const dimensionedTensor &dt)