42 template<
class BasicEddyViscosityModel>
60 (scalar(1)/betaStar_)*
sqrt(k_)/(omega_*y_),
61 scalar(500)*(this->
mu()/this->rho_)/(
sqr(y_)*omega_)
63 (4*alphaOmega2_)*k_/(CDkOmegaPlus*
sqr(y_))
72 template<
class BasicEddyViscosityModel>
79 (scalar(2)/betaStar_)*
sqrt(k_)/(omega_*y_),
80 scalar(500)*(this->
mu()/this->rho_)/(
sqr(y_)*omega_)
89 template<
class BasicEddyViscosityModel>
94 150*(this->
mu()/this->rho_)/(omega_*
sqr(y_)),
102 template<
class BasicEddyViscosityModel>
116 template<
class BasicEddyViscosityModel>
123 this->nut_ = a1_*k_/
max(a1_*omega_, b1_*F23()*
sqrt(S2));
124 this->nut_.correctBoundaryConditions();
129 template<
class BasicEddyViscosityModel>
136 template<
class BasicEddyViscosityModel>
146 template<
class BasicEddyViscosityModel>
152 return min(
G, (c1_*betaStar_)*this->k_()*this->omega_());
156 template<
class BasicEddyViscosityModel>
163 return betaStar_*omega_();
167 template<
class BasicEddyViscosityModel>
182 template<
class BasicEddyViscosityModel>
193 (c1_/a1_)*betaStar_*omega_()*
max(a1_*omega_(), b1_*
F2*
sqrt(S2))
198 template<
class BasicEddyViscosityModel>
209 template<
class BasicEddyViscosityModel>
220 template<
class BasicEddyViscosityModel>
238 template<
class BasicEddyViscosityModel>
242 const alphaField&
alpha,
248 const word& propertiesName
251 BasicEddyViscosityModel
438 bound(k_, this->kMin_);
439 bound(omega_, this->omegaMin_);
441 setDecayControl(this->coeffDict_);
447 template<
class BasicEddyViscosityModel>
453 decayControl_.readIfPresent(
"decayControl",
dict);
458 omegaInf_.read(
dict);
460 Info<<
" Employing decay control with kInf:" << kInf_
461 <<
" and omegaInf:" << omegaInf_ <<
endl;
466 omegaInf_.value() = 0;
471 template<
class BasicEddyViscosityModel>
476 alphaK1_.readIfPresent(this->coeffDict());
477 alphaK2_.readIfPresent(this->coeffDict());
478 alphaOmega1_.readIfPresent(this->coeffDict());
479 alphaOmega2_.readIfPresent(this->coeffDict());
480 gamma1_.readIfPresent(this->coeffDict());
481 gamma2_.readIfPresent(this->coeffDict());
482 beta1_.readIfPresent(this->coeffDict());
483 beta2_.readIfPresent(this->coeffDict());
484 betaStar_.readIfPresent(this->coeffDict());
485 a1_.readIfPresent(this->coeffDict());
486 b1_.readIfPresent(this->coeffDict());
487 c1_.readIfPresent(this->coeffDict());
488 F3_.readIfPresent(
"F3", this->coeffDict());
490 setDecayControl(this->coeffDict());
499 template<
class BasicEddyViscosityModel>
502 if (!this->turbulence_)
508 const alphaField&
alpha = this->alpha_;
509 const rhoField&
rho = this->rho_;
517 const volScalarField::Internal
divU 524 volScalarField::Internal GbyNu0(this->GbyNu0(tgradU(), S2));
525 volScalarField::Internal
G(this->GName(),
nut*GbyNu0);
541 omega_.boundaryFieldRef().updateCoeffs();
546 omega_.boundaryFieldRef().template evaluateCoupled<coupledFvPatch>();
565 const volScalarField::Internal
beta(this->
beta(
F1));
567 GbyNu0 = GbyNu(GbyNu0, F23(), S2());
570 tmp<fvScalarMatrix> omegaEqn
581 alpha()*
rho()*(
F1() - scalar(1))*CDkOmega()/omega_(),
590 omegaEqn.ref().relax();
592 omegaEqn.ref().boundaryManipulate(omega_.boundaryFieldRef());
595 bound(omega_, this->omegaMin_);
600 tmp<fvScalarMatrix> kEqn
609 +
alpha()*
rho()*betaStar_*omegaInf_*kInf_
620 bound(k_, this->kMin_);
void correct(GeometricField< Type, PatchField, GeoMesh > &field)
Apply correction to field.
virtual tmp< volScalarField > F23() const
dimensionedScalar tanh(const dimensionedScalar &ds)
virtual tmp< volScalarField::Internal > GbyNu(const volScalarField::Internal &GbyNu0, const volScalarField::Internal &F2, const volScalarField::Internal &S2) const
Return G/nu.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField &F1, const volTensorField &gradU) const
Return epsilon/k which for standard RAS is betaStar*omega.
virtual tmp< volScalarField::Internal > GbyNu0(const volTensorField &gradU, const volScalarField &S2) const
Return (G/nu)_0.
A list of keyword definitions, which are a keyword followed by a number of values (eg...
const dimensionedScalar G
Newtonian constant of gravitation.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual void correctNut()
tmp< GeometricField< Type, faPatchField, areaMesh > > div(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, any/none. Also accepts 0/1 as a string and shortcuts t/f, y/n.
Base class implementation of the k-omega-SST turbulence model for incompressible and compressible flo...
Generic dimensioned Type class.
tmp< GeometricField< Type, faPatchField, areaMesh > > ddt(const dimensioned< Type > dt, const faMesh &mesh)
const dimensionSet dimless
Dimensionless.
zeroField SuSp(const Foam::zero, const GeometricField< Type, faPatchField, areaMesh > &)
A no-op source.
virtual tmp< fvScalarMatrix > Qsas(const volScalarField::Internal &S2, const volScalarField::Internal &gamma, const volScalarField::Internal &beta) const
GeometricField< vector, fvPatchField, volMesh > volVectorField
virtual tmp< fvScalarMatrix > omegaSource() const
constexpr const char *const group
Group name for atomic constants.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
kOmegaSSTBase(const kOmegaSSTBase &)=delete
No copy construct.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const dimensionSet dimVolume(pow3(dimLength))
static word scopedName(const std::string &scope, const word &name)
Create scope:name or scope_name string.
const dimensionedScalar e
Elementary charge.
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< volScalarField > S2(const volTensorField &gradU) const
Return square of strain rate.
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
virtual tmp< volScalarField > F2() const
Bound the given scalar field if it has gone unbounded.
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()
virtual tmp< fvScalarMatrix > kSource() const
const dimensionedScalar mu
Atomic mass unit.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Return k production rate.
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Base-class for all transport models used by the incompressible turbulence models. ...
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar pow4(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
virtual bool read()
Re-read model coefficients if they have changed.
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
Interface to run-time selectable methods to calculate the distance-to-wall and normal-to-wall fields...
virtual tmp< volScalarField > F1(const volScalarField &CDkOmega) const
virtual tmp< volScalarField > F3() const
A class for managing temporary objects.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Defines the attributes of an object for which implicit objectRegistry management is supported...
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
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.