Go to the documentation of this file.
41 template<
class BasicEddyViscosityModel>
59 (scalar(1)/betaStar_)*
sqrt(k_)/(omega_*y_),
60 scalar(500)*(this->
mu()/this->rho_)/(
sqr(y_)*omega_)
62 (4*alphaOmega2_)*k_/(CDkOmegaPlus*
sqr(y_))
71 template<
class BasicEddyViscosityModel>
78 (scalar(2)/betaStar_)*
sqrt(k_)/(omega_*y_),
79 scalar(500)*(this->
mu()/this->rho_)/(
sqr(y_)*omega_)
88 template<
class BasicEddyViscosityModel>
93 150*(this->
mu()/this->rho_)/(omega_*
sqr(y_)),
101 template<
class BasicEddyViscosityModel>
115 template<
class BasicEddyViscosityModel>
122 this->nut_ = a1_*k_/
max(a1_*omega_, b1_*F23()*
sqrt(S2));
123 this->nut_.correctBoundaryConditions();
128 template<
class BasicEddyViscosityModel>
135 template<
class BasicEddyViscosityModel>
141 return min(
G, (c1_*betaStar_)*this->k_()*this->omega_());
145 template<
class BasicEddyViscosityModel>
153 return betaStar_*omega_();
157 template<
class BasicEddyViscosityModel>
168 (c1_/a1_)*betaStar_*omega_()
174 template<
class BasicEddyViscosityModel>
188 template<
class BasicEddyViscosityModel>
202 template<
class BasicEddyViscosityModel>
223 template<
class BasicEddyViscosityModel>
233 const word& propertiesName
236 BasicEddyViscosityModel
373 this->runTime_.timeName(),
385 this->runTime_.timeName(),
422 bound(k_, this->kMin_);
423 bound(omega_, this->omegaMin_);
425 setDecayControl(this->coeffDict_);
431 template<
class BasicEddyViscosityModel>
437 decayControl_.readIfPresent(
"decayControl",
dict);
442 omegaInf_.read(
dict);
444 Info<<
" Employing decay control with kInf:" << kInf_
445 <<
" and omegaInf:" << omegaInf_ <<
endl;
450 omegaInf_.value() = 0;
455 template<
class BasicEddyViscosityModel>
460 alphaK1_.readIfPresent(this->coeffDict());
461 alphaK2_.readIfPresent(this->coeffDict());
462 alphaOmega1_.readIfPresent(this->coeffDict());
463 alphaOmega2_.readIfPresent(this->coeffDict());
464 gamma1_.readIfPresent(this->coeffDict());
465 gamma2_.readIfPresent(this->coeffDict());
466 beta1_.readIfPresent(this->coeffDict());
467 beta2_.readIfPresent(this->coeffDict());
468 betaStar_.readIfPresent(this->coeffDict());
469 a1_.readIfPresent(this->coeffDict());
470 b1_.readIfPresent(this->coeffDict());
471 c1_.readIfPresent(this->coeffDict());
472 F3_.readIfPresent(
"F3", this->coeffDict());
474 setDecayControl(this->coeffDict());
483 template<
class BasicEddyViscosityModel>
486 if (!this->turbulence_)
507 this->
type() +
":GbyNu",
513 omega_.boundaryFieldRef().updateCoeffs();
527 GbyNu0 = GbyNu(GbyNu0, F23(), S2());
541 alpha()*
rho()*(
F1() - scalar(1))*CDkOmega()/omega_(),
555 bound(omega_, this->omegaMin_);
568 +
alpha()*
rho()*betaStar_*omegaInf_*kInf_
579 bound(k_, this->kMin_);
virtual tmp< fvScalarMatrix > Qsas(const volScalarField::Internal &S2, const volScalarField::Internal &gamma, const volScalarField::Internal &beta) const
virtual tmp< volScalarField > F2() const
Defines the attributes of an object for which implicit objectRegistry management is supported,...
void boundaryManipulate(typename GeometricField< Type, fvPatchField, volMesh >::Boundary &values)
Manipulate based on a boundary field.
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.
void clear() const noexcept
const dimensionedScalar mu
Atomic mass unit.
const dimensionedScalar G
Newtonian constant of gravitation.
A class for managing temporary objects.
virtual tmp< fvScalarMatrix > kSource() const
virtual tmp< volScalarField > F1(const volScalarField &CDkOmega) const
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
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 > omegaSource() const
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField &F1, const volTensorField &gradU) const
Return epsilon/k which for standard RAS is betaStar*omega.
void setDecayControl(const dictionary &dict)
static const wallDist & New(const fvMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Return k production rate.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
virtual void correctNut()
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
dimensionedScalar pow4(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual tmp< volScalarField::Internal > GbyNu(const volScalarField::Internal &GbyNu0, const volScalarField::Internal &F2, const volScalarField::Internal &S2) const
Return G/nu.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
dimensionedScalar tanh(const dimensionedScalar &ds)
messageStream Info
Information stream (uses stdout - output is on the master only)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Bound the given scalar field if it has gone unbounded.
virtual bool read()
Re-read model coefficients if they have changed.
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
BasicTurbulenceModel::rhoField rhoField
Generic dimensioned Type class.
BasicTurbulenceModel::alphaField alphaField
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
virtual tmp< volScalarField > F23() const
BasicTurbulenceModel::transportModel transportModel
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const dimensionedScalar e
Elementary charge.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
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.
virtual tmp< volScalarField > F3() const
Base class implementation of the k-omega-SST turbulence model for incompressible and compressible flo...
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
static Switch getOrAddToDict(const word &key, dictionary &dict, const Switch deflt=switchType::FALSE)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)