43template<
class BasicTurbulenceModel>
46 return nuTilda_/this->
nu();
50template<
class BasicTurbulenceModel>
58 return chi3/(chi3 +
pow3(Cv1_));
62template<
class BasicTurbulenceModel>
69 return scalar(1) - chi/(scalar(1) + chi*fv1);
73template<
class BasicTurbulenceModel>
90 Omega + fv2(chi(), fv1())*nuTilda_()/
sqr(kappa_*y_),
97template<
class BasicTurbulenceModel>
131template<
class BasicTurbulenceModel>
134 this->nut_ = nuTilda_*this->fv1(this->chi());
135 this->nut_.correctBoundaryConditions();
138 BasicTurbulenceModel::correctNut();
144template<
class BasicTurbulenceModel>
153 const word& propertiesName,
206 Cw1_(Cb1_/
sqr(kappa_) + (scalar(1) + Cb2_)/sigmaNut_),
259 if (
type == typeName)
261 this->printCoeffs(
type);
268template<
class BasicTurbulenceModel>
273 sigmaNut_.readIfPresent(this->coeffDict());
274 kappa_.readIfPresent(this->coeffDict());
276 Cb1_.readIfPresent(this->coeffDict());
277 Cb2_.readIfPresent(this->coeffDict());
278 Cw1_ = Cb1_/
sqr(kappa_) + (scalar(1) + Cb2_)/sigmaNut_;
280 Cw3_.readIfPresent(this->coeffDict());
281 Cv1_.readIfPresent(this->coeffDict());
282 Cs_.readIfPresent(this->coeffDict());
291template<
class BasicTurbulenceModel>
297 (nuTilda_ + this->
nu())/sigmaNut_
302template<
class BasicTurbulenceModel>
306 const scalar Cmu = 0.09;
313 this->runTime_.timeName(),
316 cbrt(this->fv1(this->chi()))
320 this->nut_.boundaryField().types()
325template<
class BasicTurbulenceModel>
329 const scalar Cmu = 0.09;
337 this->runTime_.timeName(),
340 pow(this->fv1(this->chi()), 0.5)
342 /(nuTilda_ + this->nut_ + nutSMALL),
343 this->nut_.boundaryField().types()
348template<
class BasicTurbulenceModel>
352 const scalar betaStar = 0.09;
360 this->runTime_.timeName(),
363 this->epsilon()/(betaStar*(this->k() + k0)),
364 this->nut_.boundaryField().types()
369template<
class BasicTurbulenceModel>
372 if (!this->turbulence_)
395 Cb1_*
alpha()*
rho()*Stilda*nuTilda_()
400 nuTildaEqn.
ref().relax();
405 nuTilda_.correctBoundaryConditions();
Bound the given scalar field if it has gone unbounded.
const uniformDimensionedVectorField & g
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet & dimensions() const
Return dimensions.
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.
Templated abstract base class for RAS turbulence models.
BasicTurbulenceModel::alphaField alphaField
tmp< volScalarField > chi() const
BasicTurbulenceModel::rhoField rhoField
virtual tmp< volScalarField > k() const
Return the (estimated) turbulent kinetic energy.
virtual void correct()
Solve the turbulence equations and correct the turbulent viscosity.
tmp< volScalarField > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
tmp< volScalarField > fv1(const volScalarField &chi) const
tmp< volScalarField::Internal > fv2(const volScalarField::Internal &chi, const volScalarField::Internal &fv1) const
virtual tmp< volScalarField > epsilon() const
Return the (estimated) turbulent kinetic energy dissipation rate.
virtual void correctNut()
Update nut with the latest available nuTilda.
BasicTurbulenceModel::transportModel transportModel
tmp< volScalarField::Internal > fw(const volScalarField::Internal &Stilda) const
virtual tmp< volScalarField > omega() const
Return the (estimated) specific dissipation rate.
tmp< volScalarField::Internal > Stilda() const
virtual bool read()
Re-read model coefficients if they have changed.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Generic dimensioned Type class.
bool readIfPresent(const dictionary &dict)
Eddy viscosity turbulence model base class.
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.
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.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
bool read(const char *buf, int32_t &val)
Same as readInt32.
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.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
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.
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.
dimensionedScalar cbrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
dimensionedTensor skew(const dimensionedTensor &dt)