Go to the documentation of this file.
32 #include "twoPhaseSystem.H"
33 #include "dragModel.H"
34 #include "virtualMassModel.H"
48 template<
class BasicTurbulenceModel>
49 mixtureKEpsilon<BasicTurbulenceModel>::mixtureKEpsilon
57 const word& propertiesName,
73 liquidTurbulencePtr_(
nullptr),
144 this->runTime_.timeName(),
156 this->runTime_.timeName(),
164 bound(k_, this->kMin_);
165 bound(epsilon_, this->epsilonMin_);
167 if (
type == typeName)
169 this->printCoeffs(
type);
174 template<
class BasicTurbulenceModel>
180 const volScalarField::Boundary& ebf =
epsilon.boundaryField();
186 if (isA<fixedValueFvPatchScalarField>(ebf[patchi]))
188 ebt[patchi] = fixedValueFvPatchScalarField::typeName;
196 template<
class BasicTurbulenceModel>
204 const volScalarField::Boundary& refBf =
211 isA<inletOutletFvPatchScalarField>(bf[patchi])
212 && isA<inletOutletFvPatchScalarField>(refBf[patchi])
215 refCast<inletOutletFvPatchScalarField>
216 (bf[patchi]).refValue() =
217 refCast<const inletOutletFvPatchScalarField>
218 (refBf[patchi]).refValue();
224 template<
class BasicTurbulenceModel>
240 this->runTime_.timeName(this->runTime_.startTime().value())
291 correctInletOutlet(km_(), kl);
305 mix(epsilonl, epsilong),
306 epsilonBoundaryTypes(epsilonl)
309 correctInletOutlet(epsilonm_(), epsilonl);
315 template<
class BasicTurbulenceModel>
320 Cmu_.readIfPresent(this->coeffDict());
321 C1_.readIfPresent(this->coeffDict());
322 C2_.readIfPresent(this->coeffDict());
323 C3_.readIfPresent(this->coeffDict());
324 Cp_.readIfPresent(this->coeffDict());
325 sigmak_.readIfPresent(this->coeffDict());
326 sigmaEps_.readIfPresent(this->coeffDict());
335 template<
class BasicTurbulenceModel>
338 this->nut_ = Cmu_*
sqr(k_)/epsilon_;
339 this->nut_.correctBoundaryConditions();
342 BasicTurbulenceModel::correctNut();
346 template<
class BasicTurbulenceModel>
350 if (!liquidTurbulencePtr_)
356 refCast<const twoPhaseSystem>(gas.fluid());
359 liquidTurbulencePtr_ =
373 return *liquidTurbulencePtr_;
377 template<
class BasicTurbulenceModel>
381 this->liquidTurbulence();
393 (6*this->Cmu_/(4*
sqrt(3.0/2.0)))
395 *(liquidTurbulence.
k_/liquidTurbulence.
epsilon_)
398 volScalarField fAlphad((180 + (-4.71e3 + 4.26e4*alphag)*alphag)*alphag);
400 return sqr(1 + (Ct0 - 1)*
exp(-fAlphad));
404 template<
class BasicTurbulenceModel>
409 return fluid.otherPhase(gas).rho();
413 template<
class BasicTurbulenceModel>
422 + virtualMass.
Cvm()*
fluid.otherPhase(gas).rho();
426 template<
class BasicTurbulenceModel>
432 return alphal*rholEff() + alphag*rhogEff();
436 template<
class BasicTurbulenceModel>
446 return (
alphal*rholEff()*fc + alphag*rhogEff()*fd)/rhom_();
450 template<
class BasicTurbulenceModel>
461 (
alphal*rholEff()*fc + alphag*rhogEff()*Ct2_()*fd)
462 /(
alphal*rholEff() + alphag*rhogEff()*Ct2_());
466 template<
class BasicTurbulenceModel>
488 template<
class BasicTurbulenceModel>
492 this->liquidTurbulence();
526 template<
class BasicTurbulenceModel>
529 return fvm::Su(bubbleG()/rhom_(), km_());
533 template<
class BasicTurbulenceModel>
536 return fvm::Su(C3_*epsilonm_()*bubbleG()/(rhom_()*km_()), epsilonm_());
540 template<
class BasicTurbulenceModel>
547 if (&gas != &
fluid.phase1())
551 this->liquidTurbulence();
556 if (!this->turbulence_)
574 this->liquidTurbulence();
655 bound(km, this->kMin_);
656 epsilonm == mix(epsilonl, epsilong);
657 bound(epsilonm, this->epsilonMin_);
668 -
fvm::SuSp(((2.0/3.0)*C1_)*divUm, epsilonm)
669 -
fvm::Sp(C2_*epsilonm/km, epsilonm)
674 epsEqn.
ref().relax();
679 bound(epsilonm, this->epsilonMin_);
701 bound(km, this->kMin_);
707 epsilonl = Cc2*epsilonm;
714 epsilong = Ct2_()*epsilonl;
716 nutg = Ct2_()*(liquidTurbulence.nu()/this->
nu())*nutl;
tmp< surfaceScalarField > mixFlux(const surfaceScalarField &fc, const surfaceScalarField &fd) const
tmp< volScalarField > mix(const volScalarField &fc, const volScalarField &fd) const
Info<< "Reading strained laminar flame speed field Su\n"<< endl;volScalarField Su(IOobject("Su", runTime.timeName(), 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
Defines the attributes of an object for which implicit objectRegistry management is supported,...
BasicTurbulenceModel::rhoField rhoField
tmp< volScalarField > mixU(const volScalarField &fc, const volScalarField &fd) const
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
A class for handling words, derived from Foam::string.
Class which solves the volume fraction equations for two phases.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
void clear() const noexcept
A class for managing temporary objects.
Templated abstract base class for RAS turbulence models.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
void correctInletOutlet(volScalarField &vsf, const volScalarField &refVsf) const
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
bool read(const char *buf, int32_t &val)
Same as readInt32.
wordList epsilonBoundaryTypes(const volScalarField &epsilon) const
tmp< volScalarField > bubbleG() const
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
static const word propertiesName
Default name of the turbulence properties dictionary.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
virtual void correctNut()
dimensionedScalar exp(const dimensionedScalar &ds)
BasicTurbulenceModel::transportModel transportModel
tmp< fvMatrix< Type > > Su(const DimensionedField< Type, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual tmp< fvScalarMatrix > epsilonSource() const
#define forAll(list, i)
Loop across all elements in list.
tmp< volScalarField > rhom() const
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
dimensionedScalar pow3(const dimensionedScalar &ds)
virtual tmp< volScalarField > Cvm() const =0
Return the virtual mass coefficient.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Bound the given scalar field if it has gone unbounded.
Mixture k-epsilon turbulence model for two-phase gas-liquid systems.
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
tmp< volScalarField > rholEff() const
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Generic dimensioned Type class.
tmp< volScalarField > Ct2() const
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
virtual bool read()
Re-read model coefficients if they have changed.
Calculate the matrix for implicit and explicit sources.
Base-class for all transport models used by the incompressible turbulence models.
void correctBoundaryConditions()
Correct boundary field.
scalar rho(scalar p, scalar T) const
Liquid density [kg/m^3].
tmp< volScalarField > rhogEff() const
surfaceScalarField phig("phig", -rhorAUf *ghf *fvc::snGrad(rho) *mesh.magSf())
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
virtual tmp< fvScalarMatrix > kSource() const
Eddy viscosity turbulence model base class.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
BasicTurbulenceModel::alphaField alphaField
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)