41 template<
class BasicTurbulenceModel>
42 tmp<volScalarField> EBRSM<BasicTurbulenceModel>::calcL()
const 48 pow(k_, 1.5)/epsilon_,
64 template<
class BasicTurbulenceModel>
65 tmp<volVectorField> EBRSM<BasicTurbulenceModel>::calcN()
const 74 template<
class BasicTurbulenceModel>
75 tmp<volScalarField> EBRSM<BasicTurbulenceModel>::calcTau()
const 94 template<
class BasicTurbulenceModel>
95 tmp<volSymmTensorField> EBRSM<BasicTurbulenceModel>::D
102 return (Cmu_/
sigma*tau)*this->R_ + this->
nu()*
I;
106 template<
class BasicTurbulenceModel>
107 tmp<volScalarField> EBRSM<BasicTurbulenceModel>::D
113 return this->nut_/
sigma + this->
nu();
117 template<
class BasicTurbulenceModel>
118 tmp<volScalarField> EBRSM<BasicTurbulenceModel>::Ceps1Prime
124 return Ceps1_*(scalar(1) + A1_*(scalar(1) -
pow3(f_))*
G/this->epsilon_);
128 template<
class BasicTurbulenceModel>
129 void EBRSM<BasicTurbulenceModel>::correctNut()
131 this->nut_ = Cmu_*k_*calcTau();
132 this->nut_.correctBoundaryConditions();
135 BasicTurbulenceModel::correctNut();
141 template<
class BasicTurbulenceModel>
144 const alphaField&
alpha,
150 const word& propertiesName,
166 simpleGradientDiffusion_
170 "simpleGradientDiffusion",
362 if (type == typeName)
371 template<
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());
405 template<
class BasicTurbulenceModel>
408 if (!this->turbulence_)
414 const alphaField&
alpha = this->alpha_;
415 const rhoField&
rho = this->rho_;
421 ReynoldsStress<RASModel<BasicTurbulenceModel>>
::correct();
430 tmp<volSymmTensorField> tP = -
twoSymm(
R & gradU);
441 tmp<volScalarField> tinvLsqr = scalar(1)/
sqr(calcL());
445 tmp<fvScalarMatrix> fEqn
461 tmp<volVectorField> tn = calcN();
465 tmp<volScalarField> ttau = calcTau();
472 tmp<volScalarField> tCeps1Prime = Ceps1Prime(
G);
478 epsilon_.boundaryFieldRef().template evaluateCoupled<coupledFvPatch>();
481 tmp<fvScalarMatrix> epsEqn
486 simpleGradientDiffusion_
490 +
fvm::Sp(Cstability_/k_, epsilon_)
494 + Cstability_*epsilon_()/k_()
500 epsEqn.ref().relax();
502 epsEqn.ref().boundaryManipulate(epsilon_.boundaryFieldRef());
506 bound(epsilon_, this->epsilonMin_);
513 tmp<volSymmTensorField> tPhiH;
528 (g3_ - g3star_*
mag(
B))*S
536 tmp<volSymmTensorField> tPhiW;
538 tmp<volSymmTensorField> tnn =
symm(
n*
n);
544 - scalar(5)*epsilon_/k_*
547 - 0.5*(
R && nn)*(nn +
I)
552 tmp<fvSymmTensorMatrix> REqn;
559 (scalar(1) - fCube)*tPhiW + fCube*tPhiH
565 (scalar(1) - fCube)*epsilon_/k_
577 fCube*(g1_*epsilon_ + g1star_*
G)/(scalar(2)*k_)
593 simpleGradientDiffusion_
620 this->boundNormalStress(
R);
623 this->checkRealizabilityConditions(
R);
627 bound(k_, this->kMin_);
void correct(GeometricField< Type, PatchField, GeoMesh > &field)
Apply correction to field.
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
static const sphericalTensor twoThirdsI(2.0/3.0)
dimensionedScalar kMin_
Lower limit of k.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
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.
GeometricField< tensor, fvPatchField, volMesh > volTensorField
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar pow025(const dimensionedScalar &ds)
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.
Generic dimensioned Type class.
const dimensionSet dimless
Dimensionless.
GeometricField< vector, fvPatchField, volMesh > volVectorField
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
GeometricField< sphericalTensor, fvPatchField, volMesh > volSphericalTensorField
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
bool read(const char *buf, int32_t &val)
Same as readInt32.
Templated abstract base class for RAS turbulence models.
constexpr const char *const group
Group name for atomic constants.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
static const sphericalTensor oneThirdI(1.0/3.0)
static const Identity< scalar > I
void clear()
Clear the list, i.e. set size to zero.
A class for handling words, derived from Foam::string.
Manceau (2015)'s elliptic-blending Reynolds-stress turbulence model for incompressible and compressib...
tmp< GeometricField< Type, faPatchField, areaMesh > > laplacian(const GeometricField< Type, faPatchField, areaMesh > &vf, const word &name)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
virtual void printCoeffs(const word &type)
Print model coefficients.
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
SolverPerformance< Type > solve(faMatrix< Type > &, const dictionary &solverControls)
Solve returning the solution statistics given convergence tolerance.
void updateCoeffs()
Update the boundary condition coefficients.
void boundNormalStress(volSymmTensorField &R) const
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(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()
dimensioned< Type > T() const
Return transpose.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
virtual bool read()
Re-read model coefficients if they have changed.
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)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Base-class for all transport models used by the incompressible turbulence models. ...
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
virtual void correct()
Solve the transport equations and correct the turbulent viscosity.
const dimensionedScalar & D
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...
dimensionedScalar epsilonMin_
Lower limit of epsilon.
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
Reynolds-stress turbulence model base class.
SymmTensor< Cmpt > devTwoSymm(const SymmTensor< Cmpt > &st)
Return the deviatoric part of twice the symmetric part of a SymmTensor.
static constexpr const zero Zero
Global zero (0)