41template<
class BasicTurbulenceModel>
42tmp<volScalarField> EBRSM<BasicTurbulenceModel>::calcL()
const
48 pow(k_, 1.5)/epsilon_,
64template<
class BasicTurbulenceModel>
65tmp<volVectorField> EBRSM<BasicTurbulenceModel>::calcN()
const
74template<
class BasicTurbulenceModel>
75tmp<volScalarField> EBRSM<BasicTurbulenceModel>::calcTau()
const
94template<
class BasicTurbulenceModel>
102 return (Cmu_/sigma*tau)*this->R_ + this->
nu()*
I;
106template<
class BasicTurbulenceModel>
113 return this->nut_/
sigma + this->
nu();
117template<
class BasicTurbulenceModel>
118tmp<volScalarField> EBRSM<BasicTurbulenceModel>::Ceps1Prime
124 return Ceps1_*(scalar(1) + A1_*(scalar(1) -
pow3(f_))*G/this->epsilon_);
128template<
class BasicTurbulenceModel>
129void EBRSM<BasicTurbulenceModel>::correctNut()
131 this->nut_ = Cmu_*k_*calcTau();
132 this->nut_.correctBoundaryConditions();
135 BasicTurbulenceModel::correctNut();
141template<
class BasicTurbulenceModel>
150 const word& propertiesName,
166 simpleGradientDiffusion_
170 "simpleGradientDiffusion",
325 IOobject::groupName(
"f", alphaRhoPhi.group()),
337 IOobject::groupName(
"k", alphaRhoPhi.group()),
349 IOobject::groupName(
"epsilon", alphaRhoPhi.group()),
362 if (type == typeName)
371template<
class BasicTurbulenceModel>
376 simpleGradientDiffusion_.readIfPresent
378 "simpleGradientDiffusion",
381 g1_.readIfPresent(this->coeffDict());
382 g1star_.readIfPresent(this->coeffDict());
383 g3_.readIfPresent(this->coeffDict());
384 g3star_.readIfPresent(this->coeffDict());
385 g4_.readIfPresent(this->coeffDict());
386 g5_.readIfPresent(this->coeffDict());
387 Cmu_.readIfPresent(this->coeffDict());
388 Ceps1_.readIfPresent(this->coeffDict());
389 Ceps2_.readIfPresent(this->coeffDict());
390 sigmaK_.readIfPresent(this->coeffDict());
391 sigmaEps_.readIfPresent(this->coeffDict());
392 A1_.readIfPresent(this->coeffDict());
393 Ct_.readIfPresent(this->coeffDict());
394 Cl_.readIfPresent(this->coeffDict());
395 Ceta_.readIfPresent(this->coeffDict());
396 Cstability_.readIfPresent(this->coeffDict());
405template<
class BasicTurbulenceModel>
408 if (!this->turbulence_)
476 epsilon_.boundaryFieldRef().updateCoeffs();
484 simpleGradientDiffusion_
488 +
fvm::Sp(Cstability_/k_, epsilon_)
490 alpha()*
rho()*Ceps1Prime()*G()/tau()
492 + Cstability_*epsilon_()/k_()
498 epsEqn.
ref().relax();
500 epsEqn.
ref().boundaryManipulate(epsilon_.boundaryFieldRef());
504 bound(epsilon_, this->epsilonMin_);
526 (g3_ - g3star_*
mag(
B))*S
542 - scalar(5)*epsilon_/k_*
545 - 0.5*(
R && nn)*(nn +
I)
557 (scalar(1) - fCube)*tPhiW + fCube*tPhiH
563 (scalar(1) - fCube)*epsilon_/k_
575 fCube*(g1_*epsilon_ + g1star_*G)/(scalar(2)*k_)
581 fCube*(g1_*epsilon_ + g1star_*G)*
oneThirdI
591 simpleGradientDiffusion_
618 this->boundNormalStress(
R);
621 this->checkRealizabilityConditions(
R);
625 bound(k_, this->kMin_);
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
#define R(A, B, C, D, E, F, K, M)
Bound the given scalar field if it has gone unbounded.
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Templated abstract base class for RAS turbulence models.
virtual void printCoeffs(const word &type)
Print model coefficients.
dimensionedScalar epsilonMin_
Lower limit of epsilon.
dimensionedScalar kMin_
Lower limit of k.
Manceau (2015)'s elliptic-blending Reynolds-stress turbulence model for incompressible and compressib...
BasicTurbulenceModel::alphaField alphaField
BasicTurbulenceModel::rhoField rhoField
virtual void correct()
Solve the transport equations and correct the turbulent viscosity.
BasicTurbulenceModel::transportModel transportModel
virtual bool read()
Re-read model coefficients if they have changed.
Reynolds-stress turbulence model base class.
void boundNormalStress(volSymmTensorField &R) const
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Generic dimensioned Type class.
A class for managing temporary objects.
void clear() const noexcept
A class for handling words, derived from Foam::string.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
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.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
GeometricField< vector, fvPatchField, volMesh > volVectorField
const dimensionSet dimless
Dimensionless.
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
bool read(const char *buf, int32_t &val)
Same as readInt32.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensionedScalar pow3(const dimensionedScalar &ds)
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
static const Identity< scalar > I
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
static const sphericalTensor twoThirdsI(2.0/3.0)
static const sphericalTensor oneThirdI(1.0/3.0)
dimensionedScalar pow025(const dimensionedScalar &ds)
const dimensionedScalar & D