37namespace incompressible
55 return exp((-scalar(2.5) + Rt/20.0)/
pow3(scalar(1) + Rt/130.0));
60 exp(-6.0/
sqr(scalar(1) + Rt/50.0))
61 *(scalar(1) + 3.0*
exp(-Rt/10.0));
69 return scalar(1) - 0.3*
exp(-
sqr(Rt));
90 const word& propertiesName,
152 qMin_(
"qMin",
sqrt(kMin_)),
153 zetaMin_(
"zetaMin", epsilonMin_/(2*qMin_)),
159 IOobject::groupName(
"k", alphaRhoPhi.group()),
172 IOobject::groupName(
"epsilon", alphaRhoPhi.group()),
185 IOobject::groupName(
"q", alphaRhoPhi.group()),
192 k_.boundaryField().types()
199 IOobject::groupName(
"zeta", alphaRhoPhi.group()),
205 bound(epsilon_, epsilonMin_)/(2.0*q_),
206 epsilon_.boundaryField().types()
211 if (
type == typeName)
264 zetaEqn.
ref().relax();
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Bound the given scalar field if it has gone unbounded.
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
bool readIfPresent(const word &key, const dictionary &dict)
Update the value of the Switch if it is found in the dictionary.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Generic dimensioned Type class.
bool readIfPresent(const dictionary &dict)
Eddy viscosity turbulence model base class.
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
Gibson and Dafa'Alla's q-zeta two-equation low-Re turbulence model for incompressible flows.
tmp< volScalarField > f2() const
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
dimensionedScalar sigmaZeta_
dimensionedScalar qMin_
Lower limit of q.
tmp< volScalarField > DqEff() const
Return the effective diffusivity for q.
tmp< volScalarField > DzetaEff() const
Return the effective diffusivity for epsilon.
virtual void correctNut()
dimensionedScalar zetaMin_
Lower limit of zeta.
tmp< volScalarField > fMu() const
virtual bool read()
Re-read model coefficients if they have changed.
BasicTurbulenceModel::transportModel transportModel
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< volScalarField > magSqrGradGrad(const GeometricField< Type, fvPatchField, volMesh > &vf)
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 SuSp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
RASModel< turbulenceModel > RASModel
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)