48 template<
class BasicMomentumTransportModel>
57 Cmu_*
sqr(k_)/(this->nutMaxCoeff_*Cc2*this->nu())
62 template<
class BasicMomentumTransportModel>
65 this->nut_ = Cmu_*
sqr(k_)/epsilon_;
66 this->nut_.correctBoundaryConditions();
73 template<
class BasicMomentumTransportModel>
96 gasTurbulencePtr_(nullptr),
175 this->groupName(
"k"),
176 this->runTime_.
name(),
187 this->groupName(
"epsilon"),
188 this->runTime_.
name(),
198 if (
type == typeName)
200 this->printCoeffs(
type);
203 const phaseModel& phase = refCast<const phaseModel>(this->properties());
206 if (phase.
index() == 1)
215 this->runTime_.name(),
231 this->runTime_.name(),
247 this->runTime_.name(),
264 this->runTime_.name(),
277 template<
class BasicMomentumTransportModel>
280 static bool initialised =
false;
293 template<
class BasicMomentumTransportModel>
298 Cmu_.readIfPresent(this->coeffDict());
299 C1_.readIfPresent(this->coeffDict());
300 C2_.readIfPresent(this->coeffDict());
301 C3_.readIfPresent(this->coeffDict());
302 Cp_.readIfPresent(this->coeffDict());
303 sigmak_.readIfPresent(this->coeffDict());
304 sigmaEps_.readIfPresent(this->coeffDict());
315 template<
class BasicMomentumTransportModel>
319 if (!gasTurbulencePtr_)
324 refCast<const phaseModel>(this->properties());
338 momentumTransportModel::typeName,
345 return *gasTurbulencePtr_;
349 template<
class BasicMomentumTransportModel>
353 this->gasTurbulence();
368 (6*this->Cmu_/(4*
sqrt(3.0/2.0)))
373 volScalarField fAlphad((180 + (-4.71e3 + 4.26e4*alphag)*alphag)*alphag);
375 return sqr(1 + (Ct0 - 1)*
exp(-fAlphad));
379 template<
class BasicMomentumTransportModel>
388 template<
class BasicMomentumTransportModel>
405 template<
class BasicMomentumTransportModel>
411 return alphal*rholEff() + alphag*rhogEff();
415 template<
class BasicMomentumTransportModel>
425 return (alphal*rholEff()*fc + alphag*rhogEff()*fd)/rhom_();
429 template<
class BasicMomentumTransportModel>
440 (alphal*rholEff()*fc + alphag*rhogEff()*Ct2_()*fd)
441 /(alphal*rholEff() + alphag*rhogEff()*Ct2_());
445 template<
class BasicMomentumTransportModel>
467 template<
class BasicMomentumTransportModel>
472 this->gasTurbulence();
508 template<
class BasicMomentumTransportModel>
512 return fvm::Su(bubbleG()/rhom_(), km_());
516 template<
class BasicMomentumTransportModel>
520 return fvm::Su(C3_*epsilonm_()*bubbleG()/(rhom_()*km_()), epsilonm_());
524 template<
class BasicMomentumTransportModel>
527 const phaseModel& phase = refCast<const phaseModel>(this->properties());
530 if (phase.
index() == 0)
534 this->gasTurbulence();
543 if (!this->turbulence_)
558 this->gasTurbulence();
643 epsilonm == mix(epsilonl, epsilong);
654 -
fvm::SuSp(((2.0/3.0)*C1_)*divUm, epsilonm)
655 -
fvm::Sp(C2_*epsilonm/km, epsilonm)
660 epsEqn.
ref().relax();
666 const volScalarField Cc2(rhom/(alphal*rholEff() + alphag*rhogEff()*Ct2_()));
688 bound(km, this->kMin_);
694 epsilonl = Cc2*epsilonm;
701 epsilong = Ct2_()*epsilonl;
703 nutg = Ct2_()*(this->nu()/gasTurbulence.nu())*nutl;
Bound the given scalar field where it is below the specified minimum.
void updateCoeffs()
Update the boundary condition coefficients.
Generic GeometricField class.
void correctBoundaryConditions()
Correct boundary field.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
static word groupName(Name name, const word &group)
Templated abstract base class for RAS turbulence models.
Mixture k-epsilon turbulence model for two-phase gas-liquid systems.
virtual tmp< fvScalarMatrix > epsilonSource() const
tmp< volScalarField > Ct2() const
tmp< volScalarField > rhogEff() const
tmp< volScalarField > mixU(const volScalarField &fc, const volScalarField &fd) const
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField > rhom() const
tmp< volScalarField > bubbleG() const
autoPtr< volScalarField > rhom_
autoPtr< volScalarField > km_
autoPtr< volScalarField > Ct2_
virtual void correctNut()
Correct the eddy-viscosity nut.
virtual tmp< fvScalarMatrix > kSource() const
tmp< surfaceScalarField > mixFlux(const surfaceScalarField &fc, const surfaceScalarField &fd) const
tmp< volScalarField > rholEff() const
autoPtr< volScalarField > epsilonm_
tmp< volScalarField > mix(const volScalarField &fc, const volScalarField &fd) const
void boundEpsilonm(const volScalarField &Cc2)
Bound epsilonm.
virtual bool read()
Re-read model coefficients if they have changed.
mixtureKEpsilon(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity, const word &type=typeName)
Construct from components.
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Generic dimensioned Type class.
Class to represent a interface between phases where one phase is considered dispersed within the othe...
Eddy viscosity turbulence model base class.
Finite volume constraints.
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
tmp< fvMatrix< Type > > source(const VolField< Type > &field) const
Return source for an equation.
BasicMomentumTransportModel::alphaField alphaField
BasicMomentumTransportModel::rhoField rhoField
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
scalar rho(scalar p, scalar T) const
Liquid density [kg/m^3].
tmp< volScalarField > d() const
Return the Sauter-mean diameter.
label index() const
Return the index of the phase.
const word & name() const
Return the name of this phase.
virtual const volScalarField & rho() const =0
Return the density field.
Class to represent a system of phases and model interfacial transfers between them.
const phaseModelList & phases() const
Return the phase models.
const ModelType & lookupInterfacialModel(const phaseInterface &interface) const
Return a sub model for an interface.
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.
virtual tmp< volScalarField > Cvm() const =0
Return the virtual mass coefficient.
Abstract base class for all fluid physical properties.
A class for handling words, derived from string.
Foam::fvConstraints & fvConstraints(Foam::fvConstraints::New(mesh))
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
Calculate the matrix for implicit and explicit sources.
Info<< "Reading strained laminar flame speed field Su\n"<< endl;volScalarField Su(IOobject("Su", runTime.name(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field betav\n"<< endl;volScalarField betav(IOobject("betav", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Lobs\n"<< endl;volScalarField Lobs(IOobject("Lobs", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field CT\n"<< endl;volSymmTensorField CT(IOobject("CT", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Nv\n"<< endl;volScalarField Nv(IOobject("Nv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field nsv\n"<< endl;volSymmTensorField nsv(IOobject("nsv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);IOdictionary PDRProperties(IOobject("PDRProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));autoPtr< PDRDragModel > drag
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)
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
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 > > Su(const DimensionedField< Type, volMesh > &, const VolField< Type > &)
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, 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 > > ddt(const VolField< Type > &vf)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedScalar pos(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
bool read(const char *, int32_t &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
word name(const bool)
Return a word representation of a bool.
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
const dimensionSet dimless
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimDensity
dimensioned< scalar > mag(const dimensioned< Type > &)
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.
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.