36template<
class BasicTurbulenceModel>
42 scalar kMin = this->kMin_.value();
61template<
class BasicTurbulenceModel>
75 if (isA<wallFvPatch>(curPatch))
79 const scalarField& nutw = this->nut_.boundaryField()[patchi];
83 this->U_.boundaryField()[patchi].snGrad()
87 = this->mesh_.Sf().boundaryField()[patchi];
90 = this->mesh_.magSf().boundaryField()[patchi];
96 = (faceAreas[facei]/magFaceAreas[facei])*snGradU[facei];
101 Rw[facei] = -nutw[facei]*2*
dev(
symm(gradUw));
108template<
class BasicTurbulenceModel>
114 const label maxDiffs = 5;
122 if (maxDiffs < nDiffs)
124 Info<<
"More than " << maxDiffs <<
" times"
125 <<
" Reynolds-stress realizability checks failed."
126 <<
" Skipping further comparisons." <<
endl;
135 <<
"Reynolds stress " << r <<
" at cell " << celli
136 <<
" does not obey the constraint: Rxx >= 0"
144 <<
"Reynolds stress " << r <<
" at cell " << celli
145 <<
" does not obey the constraint: Rxx*Ryy - sqr(Rxy) >= 0"
153 <<
"Reynolds stress " << r <<
" at cell " << celli
154 <<
" does not obey the constraint: det(R) >= 0"
169template<
class BasicTurbulenceModel>
172 const word& modelName,
179 const word& propertiesName
208 IOobject::groupName(
"R", alphaRhoPhi.group()),
221 IOobject::groupName(
"nut", alphaRhoPhi.group()),
230 if (couplingFactor_.value() < 0.0 || couplingFactor_.value() > 1.0)
233 <<
"couplingFactor = " << couplingFactor_
234 <<
" is not in range 0 - 1" <<
nl
242template<
class BasicTurbulenceModel>
249template<
class BasicTurbulenceModel>
257template<
class BasicTurbulenceModel>
262 tk.
ref().rename(
"k");
267template<
class BasicTurbulenceModel>
271 return devRhoReff(this->U_);
275template<
class BasicTurbulenceModel>
288 IOobject::groupName(
"devRhoReff", this->alphaRhoPhi_.group()),
289 this->runTime_.timeName(),
294 this->alpha_*this->rho_*R_
295 - (this->alpha_*this->rho_*this->nu())
302template<
class BasicTurbulenceModel>
303template<
class RhoFieldType>
307 const RhoFieldType&
rho,
311 if (couplingFactor_.value() > 0.0)
317 (1.0 - couplingFactor_)*this->alpha_*
rho*this->
nut(),
325 *this->alpha_*
rho*this->
nut()*fvc::grad(
U),
328 - fvc::div(this->alpha_*
rho*this->
nu()*
dev2(
T(fvc::grad(
U))))
329 - fvm::laplacian(this->alpha_*
rho*this->nuEff(),
U)
338 this->alpha_*
rho*this->
nut(),
342 + fvc::div(this->alpha_*
rho*R_)
343 - fvc::div(this->alpha_*
rho*this->
nu()*
dev2(
T(fvc::grad(
U))))
344 - fvm::laplacian(this->alpha_*
rho*this->nuEff(),
U)
350template<
class BasicTurbulenceModel>
357 return DivDevRhoReff(this->rho_,
U);
361template<
class BasicTurbulenceModel>
369 return DivDevRhoReff(
rho,
U);
373template<
class BasicTurbulenceModel>
380template<
class BasicTurbulenceModel>
#define R(A, B, C, D, E, F, K, M)
Defines the attributes of an object for which implicit objectRegistry management is supported,...
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
virtual bool read()
Re-read model coefficients if they have changed.
Reynolds-stress turbulence model base class.
BasicTurbulenceModel::alphaField alphaField
virtual tmp< volSymmTensorField > R() const
Return the Reynolds stress tensor.
void boundNormalStress(volSymmTensorField &R) const
BasicTurbulenceModel::rhoField rhoField
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
virtual void validate()
Validate the turbulence fields after construction.
virtual bool read()=0
Re-read model coefficients if they have changed.
tmp< fvVectorMatrix > DivDevRhoReff(const RhoFieldType &rho, volVectorField &U) const
Return the source term for the momentum equation.
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
BasicTurbulenceModel::transportModel transportModel
virtual tmp< volSymmTensorField > devRhoReff() const
Return the effective stress tensor.
void correctWallShearStress(volSymmTensorField &R) const
void checkRealizabilityConditions(const volSymmTensorField &R) const
Generic dimensioned Type class.
tmp< volVectorField > divDevRhoReff()
Return the effective viscous stress (laminar + turbulent).
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
const polyBoundaryMesh & patches
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define WarningInFunction
Report a warning using Foam::Warning.
dimensionedSymmTensor dev2(const dimensionedSymmTensor &dt)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedScalar det(const dimensionedSphericalTensor &dt)
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
errorManipArg< error, int > exit(error &err, const int errNo=1)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.