39 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
46 tmp<volScalarField> CDkOmegaPlus =
max
52 tmp<volScalarField> arg1 =
min
58 (scalar(1)/betaStar_)*
sqrt(k_)/(omega_*y_),
59 scalar(500)*this->nu()/(
sqr(y_)*omega_)
61 (4*alphaOmega2_)*k_/(CDkOmegaPlus*
sqr(y_))
69 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
74 tmp<volScalarField> arg2 =
min
78 (scalar(2)/betaStar_)*
sqrt(k_)/(omega_*y_),
79 scalar(500)*this->nu()/(
sqr(y_)*omega_)
87 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
92 tmp<volScalarField> arg3 =
min
94 150*this->nu()/(omega_*
sqr(y_)),
101 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
103 kOmegaSST<MomentumTransportModel, BasicMomentumTransportModel>::kOmegaSST::
106 tmp<volScalarField> f23(
F2());
117 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
124 this->nut_ = a1_*k_/
max(a1_*omega_, b1_*
F2*
sqrt(S2));
125 this->nut_.correctBoundaryConditions();
132 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
140 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
147 return min(
G, (c1_*betaStar_)*this->k_()*this->omega_());
151 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
159 return betaStar_*omega_();
163 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
178 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
194 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
216 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
228 MomentumTransportModel
363 this->groupName(
"k"),
364 this->runTime_.
name(),
375 this->groupName(
"omega"),
376 this->runTime_.
name(),
384 bound(k_, this->kMin_);
385 bound(omega_, this->omegaMin_);
391 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
396 alphaK1_.readIfPresent(this->coeffDict());
397 alphaK2_.readIfPresent(this->coeffDict());
398 alphaOmega1_.readIfPresent(this->coeffDict());
399 alphaOmega2_.readIfPresent(this->coeffDict());
400 gamma1_.readIfPresent(this->coeffDict());
401 gamma2_.readIfPresent(this->coeffDict());
402 beta1_.readIfPresent(this->coeffDict());
403 beta2_.readIfPresent(this->coeffDict());
404 betaStar_.readIfPresent(this->coeffDict());
405 a1_.readIfPresent(this->coeffDict());
406 b1_.readIfPresent(this->coeffDict());
407 c1_.readIfPresent(this->coeffDict());
408 F3_.readIfPresent(
"F3", this->coeffDict());
419 template<
class MomentumTransportModel,
class BasicMomentumTransportModel>
422 if (!this->turbulence_)
453 omega_.boundaryFieldRef().updateCoeffs();
478 (c1_/a1_)*betaStar_*omega_()
479 *
max(a1_*omega_(), b1_*F23()*
sqrt(S2()))
485 alpha()*
rho()*(
F1() - scalar(1))*CDkOmega()/omega_(),
488 + Qsas(S2(), gamma, beta)
493 omegaEqn.
ref().relax();
495 omegaEqn.
ref().boundaryManipulate(omega_.boundaryFieldRef());
498 bound(omega_, this->omegaMin_);
519 bound(k_, this->kMin_);
Bound the given scalar field where it is below the specified minimum.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Generic dimensioned Type class.
Finite volume constraints.
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
BasicMomentumTransportModel::alphaField alphaField
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField::Internal &F1, const volScalarField::Internal &F2) const
Return epsilon/k which for standard RAS is betaStar*omega.
kOmegaSST(const word &type, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity)
Construct from components.
virtual tmp< fvScalarMatrix > Qsas(const volScalarField::Internal &S2, const volScalarField::Internal &gamma, const volScalarField::Internal &beta) const
virtual tmp< fvScalarMatrix > omegaSource() const
virtual void correctNut()
virtual tmp< fvScalarMatrix > kSource() const
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Return k production rate.
virtual bool read()
Re-read model coefficients if they have changed.
BasicMomentumTransportModel::rhoField rhoField
A class for managing temporary objects.
void clear() const
If object pointer points to valid object:
T & ref() const
Return non-const reference or generate a fatal error.
Abstract base class for all fluid physical properties.
Interface to run-time selectable methods to calculate the distance-to-wall and normal-to-wall fields.
A class for handling words, derived from string.
Foam::fvConstraints & fvConstraints(Foam::fvConstraints::New(mesh))
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
void correct(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su)
const dimensionedScalar G
Newtonian constant of gravitation.
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const VolField< Type > &)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
const dimensionSet dimless
dimensionedScalar tanh(const dimensionedScalar &ds)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimTime
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimVolume
VolField< scalar > volScalarField
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
dimensionedScalar pow4(const dimensionedScalar &ds)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
bool bound(volScalarField &, const dimensionedScalar &min)
Bound the given scalar field where it is below the specified min value.
word name(const complex &)
Return a string representation of a complex.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.