Go to the documentation of this file.
35 template<
class BasicTurbulenceModel>
41 scalar kMin = this->kMin_.value();
60 template<
class BasicTurbulenceModel>
68 volSymmTensorField::Boundary& RBf =
R.boundaryFieldRef();
74 if (isA<wallFvPatch>(curPatch))
78 const scalarField& nutw = this->nut_.boundaryField()[patchi];
82 this->U_.boundaryField()[patchi].snGrad()
86 = this->mesh_.Sf().boundaryField()[patchi];
89 = this->mesh_.magSf().boundaryField()[patchi];
95 = (faceAreas[facei]/magFaceAreas[facei])*snGradU[facei];
100 Rw[facei] = -nutw[facei]*2*
dev(
symm(gradUw));
109 template<
class BasicTurbulenceModel>
112 const word& modelName,
119 const word& propertiesName
148 IOobject::groupName(
"R", alphaRhoPhi.group()),
149 this->runTime_.timeName(),
161 IOobject::groupName(
"nut", alphaRhoPhi.group()),
162 this->runTime_.timeName(),
170 if (couplingFactor_.value() < 0.0 || couplingFactor_.value() > 1.0)
173 <<
"couplingFactor = " << couplingFactor_
174 <<
" is not in range 0 - 1" <<
nl
182 template<
class BasicTurbulenceModel>
189 template<
class BasicTurbulenceModel>
197 template<
class BasicTurbulenceModel>
202 tk.
ref().rename(
"k");
207 template<
class BasicTurbulenceModel>
217 IOobject::groupName(
"devRhoReff", this->alphaRhoPhi_.group()),
218 this->runTime_.timeName(),
223 this->alpha_*this->rho_*R_
224 - (this->alpha_*this->rho_*this->
nu())
231 template<
class BasicTurbulenceModel>
232 template<
class RhoFieldType>
236 const RhoFieldType&
rho,
240 if (couplingFactor_.value() > 0.0)
246 (1.0 - couplingFactor_)*this->alpha_*
rho*this->
nut(),
267 this->alpha_*
rho*this->
nut(),
279 template<
class BasicTurbulenceModel>
286 return DivDevRhoReff(this->rho_,
U);
290 template<
class BasicTurbulenceModel>
298 return DivDevRhoReff(
rho,
U);
302 template<
class BasicTurbulenceModel>
309 template<
class BasicTurbulenceModel>
Defines the attributes of an object for which implicit objectRegistry management is supported,...
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
virtual bool read()=0
Re-read model coefficients if they have changed.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
void correctWallShearStress(volSymmTensorField &R) const
tmp< GeometricField< typename outerProduct< vector, Type >::type, faPatchField, areaMesh >> grad(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
dimensionedSymmTensor dev2(const dimensionedSymmTensor &dt)
tmp< GeometricField< Type, faPatchField, areaMesh > > div(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
#define forAll(list, i)
Loop across all elements in list.
#define R(A, B, C, D, E, F, K, M)
virtual tmp< fvVectorMatrix > divDevRhoReff(volVectorField &U) const
Return the source term for the momentum equation.
ReynoldsStress(const word &modelName, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName)
Construct from components.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
tmp< fvVectorMatrix > DivDevRhoReff(const RhoFieldType &rho, volVectorField &U) const
Return the source term for the momentum equation.
SymmTensor< scalar > symmTensor
SymmTensor of scalars.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
BasicTurbulenceModel::rhoField rhoField
Generic dimensioned Type class.
virtual void validate()
Validate the turbulence fields after construction.
BasicTurbulenceModel::alphaField alphaField
errorManipArg< error, int > exit(error &err, const int errNo=1)
virtual tmp< volSymmTensorField > devRhoReff() const
Return the effective stress tensor.
BasicTurbulenceModel::transportModel transportModel
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
virtual tmp< volSymmTensorField > R() const
Return the Reynolds stress tensor.
const polyBoundaryMesh & patches
void boundNormalStress(volSymmTensorField &R) const
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
tmp< GeometricField< Type, faPatchField, areaMesh > > laplacian(const GeometricField< Type, faPatchField, areaMesh > &vf, const word &name)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)