Go to the documentation of this file.
42 template<
class BasicTurbulenceModel>
64 template<
class BasicTurbulenceModel>
73 pow(k_, 1.5)/epsilon_,
89 template<
class BasicTurbulenceModel>
92 this->nut_ =
min(CmuKEps_*
sqr(k_)/epsilon_, this->Cmu_*v2_*Ts());
93 this->nut_.correctBoundaryConditions();
96 BasicTurbulenceModel::correctNut();
102 template<
class BasicTurbulenceModel>
111 const word& propertiesName,
224 this->runTime_.timeName(),
236 this->runTime_.timeName(),
248 this->runTime_.timeName(),
260 this->runTime_.timeName(),
270 bound(k_, this->kMin_);
271 bound(epsilon_, this->epsilonMin_);
275 if (
type == typeName)
277 this->printCoeffs(
type);
284 template<
class BasicTurbulenceModel>
289 Cmu_.readIfPresent(this->coeffDict());
290 CmuKEps_.readIfPresent(this->coeffDict());
291 C1_.readIfPresent(this->coeffDict());
292 C2_.readIfPresent(this->coeffDict());
293 CL_.readIfPresent(this->coeffDict());
294 Ceta_.readIfPresent(this->coeffDict());
295 Ceps2_.readIfPresent(this->coeffDict());
296 Ceps3_.readIfPresent(this->coeffDict());
297 sigmaK_.readIfPresent(this->coeffDict());
298 sigmaEps_.readIfPresent(this->coeffDict());
307 template<
class BasicTurbulenceModel>
310 if (!this->turbulence_)
339 1.0/Ts*((C1_ -
N)*v2_ - 2.0/3.0*k_*(C1_ - 1.0))
345 1.4*(1.0 + 0.05*
min(
sqrt(k_/v2_), scalar(100)))
349 epsilon_.boundaryFieldRef().updateCoeffs();
364 epsEqn.
ref().relax();
366 epsEqn.
ref().boundaryManipulate(epsilon_.boundaryFieldRef());
369 bound(epsilon_, this->epsilonMin_);
389 bound(k_, this->kMin_);
398 - 1.0/L2/k_*(v2fAlpha - C2_*
G)
Defines the attributes of an object for which implicit objectRegistry management is supported,...
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
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.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
const dimensionedScalar G
Newtonian constant of gravitation.
static constexpr const zero Zero
Global zero.
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.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
virtual bool read()
Re-read model coefficients if they have changed.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
tmp< volScalarField > Ts() const
Return time scale, Ts.
dimensionedScalar pow025(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
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)
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.
v2f(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
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.
BasicTurbulenceModel::alphaField alphaField
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
BasicTurbulenceModel::transportModel transportModel
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
static word groupName(StringType name, const word &group)
Create dot-delimited name.group.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
virtual void correctNut()
Eddy viscosity turbulence model base class.
const Vector< label > N(dict.get< Vector< label >>("N"))
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Abstract base-class for v2-f models to provide BCs access to the v2 and f fields.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
tmp< volScalarField > Ls() const
Return length scale, Ls.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)