Go to the documentation of this file.
41 template<
class BasicTurbulenceModel>
45 this->nut_ = Cmu_*phit_*k_*T_;
46 this->nut_.correctBoundaryConditions();
49 BasicTurbulenceModel::correctNut();
53 template<
class BasicTurbulenceModel>
73 template<
class BasicTurbulenceModel>
80 pow(k_, 1.5)/epsilon_,
98 template<
class BasicTurbulenceModel>
107 const word& propertiesName,
247 this->runTime_.timeName(),
259 this->runTime_.timeName(),
271 this->runTime_.timeName(),
283 this->runTime_.timeName(),
295 this->runTime_.timeName(),
310 bound(k_, this->kMin_);
311 bound(epsilon_, this->epsilonMin_);
312 bound(phit_, phitMin_);
315 if (
type == typeName)
317 this->printCoeffs(
type);
322 mag(sigmaK_.value()) < VSMALL
323 ||
mag(sigmaEps_.value()) < VSMALL
324 ||
mag(sigmaPhit_.value()) < VSMALL
328 <<
"Non-zero values are required for the model constants:" <<
nl
329 <<
"sigmaK = " << sigmaK_ <<
nl
330 <<
"sigmaEps = " << sigmaEps_ <<
nl
331 <<
"sigmaPhit = " << sigmaPhit_ <<
nl
339 template<
class BasicTurbulenceModel>
344 Cmu_.readIfPresent(this->coeffDict());
345 Ceps1a_.readIfPresent(this->coeffDict());
346 Ceps1b_.readIfPresent(this->coeffDict());
347 Ceps1c_.readIfPresent(this->coeffDict());
348 Ceps2_.readIfPresent(this->coeffDict());
349 Cf1_.readIfPresent(this->coeffDict());
350 Cf2_.readIfPresent(this->coeffDict());
351 CL_.readIfPresent(this->coeffDict());
352 Ceta_.readIfPresent(this->coeffDict());
353 CT_.readIfPresent(this->coeffDict());
354 sigmaK_.readIfPresent(this->coeffDict());
355 sigmaEps_.readIfPresent(this->coeffDict());
356 sigmaPhit_.readIfPresent(this->coeffDict());
365 template<
class BasicTurbulenceModel>
368 if (!this->turbulence_)
400 Ceps1a_*(Ceps1b_ + Ceps1c_*
sqrt(1.0/phit_()))
405 epsilon_.boundaryFieldRef().updateCoeffs();
425 epsEqn.
ref().relax();
427 epsEqn.
ref().boundaryManipulate(epsilon_.boundaryFieldRef());
430 bound(epsilon_, this->epsilonMin_);
451 bound(k_, this->kMin_);
462 (Cf1_ - 1.0)*(phit_() - 2.0/3.0)/T_()
464 +(Cf2_*(2.0/3.0)*
divU)
491 /(k_()*sigmaPhit_*phit_())
498 phitEqn.
ref().relax();
502 bound(phit_, phitMin_);
Defines the attributes of an object for which implicit objectRegistry management is supported,...
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
virtual void correct()
Solve the transport equations and correct the turbulent viscosity.
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 dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
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.
BasicTurbulenceModel::alphaField alphaField
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)
BasicTurbulenceModel::rhoField rhoField
dimensionedScalar pow025(const dimensionedScalar &ds)
tmp< volScalarField > Ls() const
Return the turbulent length scale, L.
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
BasicTurbulenceModel::transportModel transportModel
virtual void correctNut()
Update nut with the latest available k, phit, and T.
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.
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.
Generic dimensioned Type class.
Abstract base-class for the k-epsilon-phit-f model to provide boundary condition access to the phit a...
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
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.
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.
Eddy viscosity turbulence model base class.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
The k-epsilon-phit-f turbulence closure model for incompressible and compressible flows.
virtual bool read()
Re-read model coefficients if they have changed.
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...
tmp< volScalarField > Ts() const
Return the turbulent time scale, T.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)