36 template<
class BasicTurbulenceModel>
42 const scalar kMin = this->kMin_.value();
56 template<
class BasicTurbulenceModel>
70 if (isA<wallFvPatch>(curPatch))
74 const scalarField& nutw = this->nut_.boundaryField()[patchi];
78 this->U_.boundaryField()[patchi].snGrad()
82 = this->mesh_.Sf().boundaryField()[patchi];
85 = this->mesh_.magSf().boundaryField()[patchi];
91 = (faceAreas[facei]/magFaceAreas[facei])*snGradU[facei];
96 Rw[facei] = -nutw[facei]*2*
devSymm(gradUw);
103 template<
class BasicTurbulenceModel>
109 const label maxDiffs = 5;
117 if (maxDiffs < nDiffs)
119 Info<<
"More than " << maxDiffs <<
" times" 120 <<
" Reynolds-stress realizability checks failed." 121 <<
" Skipping further comparisons." <<
endl;
130 <<
"Reynolds stress " << r <<
" at cell " << celli
131 <<
" does not obey the constraint: Rxx >= 0" 136 if (r.xx()*r.yy() -
sqr(r.xy()) < 0)
139 <<
"Reynolds stress " << r <<
" at cell " << celli
140 <<
" does not obey the constraint: Rxx*Ryy - sqr(Rxy) >= 0" 148 <<
"Reynolds stress " << r <<
" at cell " << celli
149 <<
" does not obey the constraint: det(R) >= 0" 164 template<
class BasicTurbulenceModel>
167 const word& modelName,
168 const alphaField&
alpha,
174 const word& propertiesName
225 if (couplingFactor_.value() < 0.0 || couplingFactor_.value() > 1.0)
228 <<
"couplingFactor = " << couplingFactor_
229 <<
" is not in range 0 - 1" <<
nl 237 template<
class BasicTurbulenceModel>
244 template<
class BasicTurbulenceModel>
252 template<
class BasicTurbulenceModel>
257 tk.ref().rename(
"k");
262 template<
class BasicTurbulenceModel>
266 return devRhoReff(this->U_);
270 template<
class BasicTurbulenceModel>
283 IOobject::groupName(
"devRhoReff", this->alphaRhoPhi_.group()),
284 this->runTime_.timeName(),
289 this->alpha_*this->rho_*R_
290 - (this->alpha_*this->rho_*this->
nu())
297 template<
class BasicTurbulenceModel>
298 template<
class RhoFieldType>
302 const RhoFieldType&
rho,
306 if (couplingFactor_.value() > 0.0)
312 (1.0 - couplingFactor_)*this->alpha_*
rho*this->
nut(),
333 this->alpha_*
rho*this->
nut(),
345 template<
class BasicTurbulenceModel>
352 return DivDevRhoReff(this->rho_,
U);
356 template<
class BasicTurbulenceModel>
364 return DivDevRhoReff(
rho,
U);
368 template<
class BasicTurbulenceModel>
375 template<
class BasicTurbulenceModel>
void correctWallShearStress(volSymmTensorField &R) const
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
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)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
tmp< fvVectorMatrix > DivDevRhoReff(const RhoFieldType &rho, volVectorField &U) const
Return the source term for the momentum equation.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
constexpr char nl
The newline '\n' character (0x0a)
tmp< GeometricField< Type, faPatchField, areaMesh > > div(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
Ostream & endl(Ostream &os)
Add newline and flush stream.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
virtual tmp< volSymmTensorField > R() const
Return the Reynolds stress tensor.
SymmTensor< Cmpt > devSymm(const SymmTensor< Cmpt > &st)
Return the deviatoric part of the symmetric part of a SymmTensor.
virtual tmp< fvVectorMatrix > divDevRhoReff(volVectorField &U) const
Return the source term for the momentum equation.
Generic dimensioned Type class.
dimensionedScalar det(const dimensionedSphericalTensor &dt)
GeometricField< vector, fvPatchField, volMesh > volVectorField
#define forAll(list, i)
Loop across all elements in list.
constexpr const char *const group
Group name for atomic constants.
void checkRealizabilityConditions(const volSymmTensorField &R) const
virtual bool read()=0
Re-read model coefficients if they have changed.
SymmTensor< scalar > symmTensor
SymmTensor of scalars, i.e. SymmTensor<scalar>.
A class for handling words, derived from Foam::string.
tmp< GeometricField< Type, faPatchField, areaMesh > > laplacian(const GeometricField< Type, faPatchField, areaMesh > &vf, const word &name)
virtual tmp< volSymmTensorField > devRhoReff() const
Return the effective stress tensor.
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
virtual void validate()
Validate the turbulence fields after construction.
void boundNormalStress(volSymmTensorField &R) const
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil) - psi) *pSat, rhoMin);# 1 "/home/chef2/andy/OpenFOAM/release/v2406/OpenFOAM-v2406/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=clamp((rho - rholSat)/(rhovSat - rholSat), zero_one{});alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
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.
dimensionedSymmTensor dev2(const dimensionedSymmTensor &dt)
#define R(A, B, C, D, E, F, K, M)
#define WarningInFunction
Report a warning using Foam::Warning.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
Base-class for all transport models used by the incompressible turbulence models. ...
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
const polyBoundaryMesh & patches
messageStream Info
Information stream (stdout output on master, null elsewhere)
A class for managing temporary objects.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Tensor of scalars, i.e. Tensor<scalar>.
Defines the attributes of an object for which implicit objectRegistry management is supported...
tmp< GeometricField< typename outerProduct< vector, Type >::type, faPatchField, areaMesh >> grad(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
SymmTensor< Cmpt > devTwoSymm(const SymmTensor< Cmpt > &st)
Return the deviatoric part of twice the symmetric part of a SymmTensor.