Go to the documentation of this file.
42 template<
class BasicTurbulenceModel>
45 this->nut_ = Cmu_*
sqr(k_)/epsilon_;
46 this->nut_.correctBoundaryConditions();
49 BasicTurbulenceModel::correctNut();
53 template<
class BasicTurbulenceModel>
61 dimVolume*this->rho_.dimensions()*k_.dimensions()
68 template<
class BasicTurbulenceModel>
76 dimVolume*this->rho_.dimensions()*epsilon_.dimensions()
85 template<
class BasicTurbulenceModel>
94 const word& propertiesName,
188 this->runTime_.timeName(),
200 this->runTime_.timeName(),
208 bound(k_, this->kMin_);
209 bound(epsilon_, this->epsilonMin_);
211 if (
type == typeName)
213 this->printCoeffs(
type);
220 template<
class BasicTurbulenceModel>
225 Cmu_.readIfPresent(this->coeffDict());
226 C1_.readIfPresent(this->coeffDict());
227 C2_.readIfPresent(this->coeffDict());
228 C3_.readIfPresent(this->coeffDict());
229 sigmak_.readIfPresent(this->coeffDict());
230 sigmaEps_.readIfPresent(this->coeffDict());
231 eta0_.readIfPresent(this->coeffDict());
232 beta_.readIfPresent(this->coeffDict());
241 template<
class BasicTurbulenceModel>
244 if (!this->turbulence_)
278 ((eta*(-eta/eta0_ + scalar(1)))/(beta_*eta3 + scalar(1)))
282 epsilon_.boundaryFieldRef().updateCoeffs();
298 epsEqn.
ref().relax();
300 epsEqn.
ref().boundaryManipulate(epsilon_.boundaryFieldRef());
303 bound(epsilon_, this->epsilonMin_);
325 bound(k_, this->kMin_);
Defines the attributes of an object for which implicit objectRegistry management is supported,...
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 correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
void clear() const noexcept
const dimensionedScalar G
Newtonian constant of gravitation.
A class for managing temporary objects.
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: [].
virtual tmp< fvScalarMatrix > kSource() const
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.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
virtual tmp< fvScalarMatrix > epsilonSource() const
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Renormalization group k-epsilon turbulence model for incompressible and compressible flows.
#define R(A, B, C, D, E, F, K, M)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
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.
BasicTurbulenceModel::rhoField rhoField
Generic dimensioned Type class.
BasicTurbulenceModel::alphaField alphaField
virtual void correctNut()
BasicTurbulenceModel::transportModel transportModel
virtual bool read()
Re-read model coefficients if they have changed.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
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.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Eddy viscosity turbulence model base class.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
const dimensionSet dimVolume(pow3(dimLength))
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)