39 template<
class BasicMomentumTransportModel>
42 return exp(-3.4/
sqr(scalar(1) +
sqr(k_)/(this->nu()*epsilon_)/50.0));
46 template<
class BasicMomentumTransportModel>
51 - 0.3*
exp(-
min(
sqr(
sqr(k_)/(this->nu()*epsilon_)), scalar(50.0)));
55 template<
class BasicMomentumTransportModel>
58 this->nut_ = Cmu_*fMu()*
sqr(k_)/epsilon_;
59 this->nut_.correctBoundaryConditions();
64 template<
class BasicMomentumTransportModel>
73 dimVolume*this->rho_.dimensions()*k_.dimensions()
80 template<
class BasicMomentumTransportModel>
89 dimVolume*this->rho_.dimensions()*epsilon_.dimensions()
98 template<
class BasicMomentumTransportModel>
181 this->runTime_.timeName(),
194 this->runTime_.timeName(),
202 bound(k_, this->kMin_);
203 bound(epsilon_, this->epsilonMin_);
205 if (type == typeName)
207 this->printCoeffs(type);
214 template<
class BasicMomentumTransportModel>
219 Cmu_.readIfPresent(this->coeffDict());
220 C1_.readIfPresent(this->coeffDict());
221 C2_.readIfPresent(this->coeffDict());
222 C3_.readIfPresent(this->coeffDict());
223 sigmak_.readIfPresent(this->coeffDict());
224 sigmaEps_.readIfPresent(this->coeffDict());
235 template<
class BasicMomentumTransportModel>
238 if (!this->turbulence_)
273 C1_*alpha*rho*G*epsilon_/k_
274 -
fvm::SuSp(((2.0/3.0)*C1_ - C3_)*alpha*rho*divU, epsilon_)
275 -
fvm::Sp(C2_*f2()*alpha*rho*epsilon_/k_, epsilon_)
281 epsEqn.
ref().relax();
283 epsEqn.
ref().boundaryManipulate(epsilon_.boundaryFieldRef());
286 bound(epsilon_, this->epsilonMin_);
296 alpha*rho*G -
fvm::SuSp(2.0/3.0*alpha*rho*divU, k_)
297 -
fvm::Sp(alpha*rho*(epsilon_ + D)/k_, k_)
306 bound(k_, this->kMin_);
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
BasicMomentumTransportModel::alphaField alphaField
fvMatrix< scalar > fvScalarMatrix
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
void clear() const
If object pointer points to valid object:
virtual bool read()
Re-read model coefficients if they have changed.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
T & ref() const
Return non-const reference or generate a fatal error.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
dimensionedScalar sqrt(const dimensionedScalar &ds)
tmp< volScalarField > magSqrGradGrad(const GeometricField< Type, fvPatchField, volMesh > &vf)
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Generic dimensioned Type class.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Eddy viscosity turbulence model base class.
virtual tmp< fvScalarMatrix > kSource() const
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
Templated abstract base class for RAS turbulence models.
BasicMomentumTransportModel::rhoField rhoField
const dimensionSet dimVolume(pow3(dimLength))
bool read(const char *, int32_t &)
virtual tmp< fvScalarMatrix > epsilonSource() const
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedScalar exp(const dimensionedScalar &ds)
A class for handling words, derived from string.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
tmp< volScalarField > f2() const
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Bound the given scalar field if it has gone unbounded.
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-8/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
LaunderSharmaKE(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &type=typeName)
Construct from components.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
virtual void correctNut()
tmp< volScalarField > fMu() const
A class for managing temporary objects.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
BasicMomentumTransportModel::transportModel transportModel
const dimensionedScalar & G
Newtonian constant of gravitation.