39 template<
class BasicTurbulenceModel>
42 this->nut_ = this->Cmu_*
sqr(k_)/epsilon_;
43 this->nut_.correctBoundaryConditions();
46 BasicTurbulenceModel::correctNut();
52 template<
class BasicTurbulenceModel>
61 const word& propertiesName,
183 this->runTime_.timeName(),
195 this->runTime_.timeName(),
203 if (type == typeName)
205 this->printCoeffs(type);
207 this->boundNormalStress(this->R_);
208 bound(epsilon_, this->epsilonMin_);
209 k_ = 0.5*
tr(this->R_);
216 template<
class BasicTurbulenceModel>
221 Cmu_.readIfPresent(this->coeffDict());
222 C1_.readIfPresent(this->coeffDict());
223 C2_.readIfPresent(this->coeffDict());
224 Ceps1_.readIfPresent(this->coeffDict());
225 Ceps2_.readIfPresent(this->coeffDict());
226 Cs_.readIfPresent(this->coeffDict());
227 Ceps_.readIfPresent(this->coeffDict());
229 wallReflection_.readIfPresent(
"wallReflection", this->coeffDict());
230 kappa_.readIfPresent(this->coeffDict());
231 Cref1_.readIfPresent(this->coeffDict());
232 Cref2_.readIfPresent(this->coeffDict());
243 template<
class BasicTurbulenceModel>
249 (Cs_*(this->k_/this->epsilon_))*this->R_ +
I*this->
nu()
254 template<
class BasicTurbulenceModel>
260 (Ceps_*(this->k_/this->epsilon_))*this->R_ +
I*this->
nu()
265 template<
class BasicTurbulenceModel>
268 if (!this->turbulence_)
290 epsilon_.boundaryFieldRef().updateCoeffs();
299 Ceps1_*alpha*rho*G*epsilon_/k_
300 -
fvm::Sp(Ceps2_*alpha*rho*epsilon_/k_, epsilon_)
304 epsEqn.
ref().relax();
306 epsEqn.
ref().boundaryManipulate(epsilon_.boundaryFieldRef());
309 bound(epsilon_, this->epsilonMin_);
320 if (isA<wallFvPatch>(curPatch))
327 G[celli]/(0.5*
mag(
tr(P[celli])) + small),
340 +
fvm::Sp(C1_*alpha*rho*epsilon_/k_, R)
343 - (2.0/3.0*(1 - C1_)*
I)*alpha*rho*epsilon_
344 - C2_*alpha*rho*
dev(P)
356 Cref1_*R - ((Cref2_*C2_)*(k_/epsilon_))*
dev(P)
360 ((3*
pow(Cmu_, 0.75)/kappa_)*(alpha*rho*
sqrt(k_)/y_))
369 this->boundNormalStress(R);
376 this->correctWallShearStress(R);
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
tmp< volSymmTensorField > DREff() const
Return the effective diffusivity for R.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
const dimensionedScalar G
Newtonian constant of gravitation.
T & ref() const
Return non-const reference or generate a fatal error.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Generic dimensioned Type class.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual bool read()
Read model coefficients if they have changed.
static Switch lookupOrAddToDict(const word &, dictionary &, const Switch &defaultValue=false)
Construct from dictionary, supplying default value so that if the.
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
tmp< volSymmTensorField > DepsilonEff() const
Return the effective diffusivity for epsilon.
static const wallDist & New(const fvMesh &mesh)
static tmp< GeometricField< symmTensor, fvPatchField, volMesh > > New(const word &name, const Mesh &, const dimensionSet &, const word &patchFieldType=fvPatchField< symmTensor >::calculatedType())
Return a temporary field constructed from name, mesh, dimensionSet.
Templated abstract base class for RAS turbulence models.
bool read(const char *, int32_t &)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
static const Identity< scalar > I
virtual void correctNut()
Update the eddy-viscosity.
A class for handling words, derived from string.
virtual const labelUList & faceCells() const
Return faceCells.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
BasicTurbulenceModel::transportModel transportModel
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
LRR(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-7/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 > &)
BasicTurbulenceModel::alphaField alphaField
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
virtual void correct()
Solve the turbulence equations and correct eddy-Viscosity and.
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
#define R(A, B, C, D, E, F, K, M)
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
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.
dimensioned< scalar > mag(const dimensioned< Type > &)
BasicTurbulenceModel::rhoField rhoField
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
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...
Reynolds-stress turbulence model base class.