Go to the documentation of this file.
43 template<
class BasicTurbulenceModel>
46 return nuTilda_/this->
nu();
50 template<
class BasicTurbulenceModel>
57 return chi3/(chi3 +
pow3(Cv1_));
61 template<
class BasicTurbulenceModel>
68 return 1.0 - chi/(1.0 + chi*fv1);
72 template<
class BasicTurbulenceModel>
86 + fv2(chi, fv1)*nuTilda_/
sqr(kappa_*y_),
93 template<
class BasicTurbulenceModel>
123 template<
class BasicTurbulenceModel>
129 this->nut_ = nuTilda_*fv1;
133 BasicTurbulenceModel::correctNut();
137 template<
class BasicTurbulenceModel>
140 correctNut(fv1(this->chi()));
146 template<
class BasicTurbulenceModel>
155 const word& propertiesName,
208 Cw1_(Cb1_/
sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_),
251 this->runTime_.timeName(),
261 if (
type == typeName)
263 this->printCoeffs(
type);
270 template<
class BasicTurbulenceModel>
275 sigmaNut_.readIfPresent(this->coeffDict());
276 kappa_.readIfPresent(this->coeffDict());
278 Cb1_.readIfPresent(this->coeffDict());
279 Cb2_.readIfPresent(this->coeffDict());
280 Cw1_ = Cb1_/
sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
282 Cw3_.readIfPresent(this->coeffDict());
283 Cv1_.readIfPresent(this->coeffDict());
284 Cs_.readIfPresent(this->coeffDict());
293 template<
class BasicTurbulenceModel>
303 template<
class BasicTurbulenceModel>
307 <<
"Turbulence kinetic energy not defined for "
308 <<
"Spalart-Allmaras model. Returning zero field"
316 this->runTime_.timeName(),
326 template<
class BasicTurbulenceModel>
330 <<
"Turbulence kinetic energy dissipation rate not defined for "
331 <<
"Spalart-Allmaras model. Returning zero field"
339 this->runTime_.timeName(),
349 template<
class BasicTurbulenceModel>
353 <<
"Specific dissipation rate not defined for "
354 <<
"Spalart-Allmaras model. Returning zero field"
362 this->runTime_.timeName(),
371 template<
class BasicTurbulenceModel>
374 if (!this->turbulence_)
405 nuTildaEqn.
ref().relax();
410 nuTilda_.correctBoundaryConditions();
tmp< volScalarField > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
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.
virtual void correctNut()
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
dimensionedTensor skew(const dimensionedTensor &dt)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
tmp< volScalarField > fv1(const volScalarField &chi) const
Templated abstract base class for RAS turbulence models.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
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.
bool read(const char *buf, int32_t &val)
Same as readInt32.
static const wallDist & New(const fvMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
Ostream & endl(Ostream &os)
Add newline and flush stream.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
virtual tmp< volScalarField > omega() const
Return the (estimated) specific dissipation rate.
This boundary condition applies a zero-gradient condition from the patch internal field onto the patc...
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
tmp< volScalarField > fv2(const volScalarField &chi, const volScalarField &fv1) const
tmp< volScalarField > Stilda(const volScalarField &chi, const volScalarField &fv1) const
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
tmp< volScalarField > fw(const volScalarField &Stilda) const
dimensionedScalar pow6(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
bool readIfPresent(const dictionary &dict)
dimensionedScalar pow3(const dimensionedScalar &ds)
virtual bool read()
Re-read model coefficients if they have changed.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Bound the given scalar field if it has gone unbounded.
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
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.
BasicTurbulenceModel::rhoField rhoField
Generic dimensioned Type class.
const uniformDimensionedVectorField & g
BasicTurbulenceModel::alphaField alphaField
void correctBoundaryConditions()
Correct boundary field.
BasicTurbulenceModel::transportModel transportModel
Spalart-Allmaras one-eqn mixing-length model for incompressible and compressible external flows.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
dimensionedScalar sqrt(const dimensionedScalar &ds)
tmp< volScalarField > chi() const
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.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
Eddy viscosity turbulence model base class.
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 tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.