28 #include "twoPhaseSystem.H" 39 template<
class BasicTurbulenceModel>
40 NicenoKEqn<BasicTurbulenceModel>::NicenoKEqn
48 const word& propertiesName,
64 gasTurbulencePtr_(
nullptr),
98 this->printCoeffs(type);
105 template<
class BasicTurbulenceModel>
110 alphaInversion_.readIfPresent(this->coeffDict());
111 Cp_.readIfPresent(this->coeffDict());
112 Cmub_.readIfPresent(this->coeffDict());
123 template<
class BasicTurbulenceModel>
126 typename BasicTurbulenceModel::transportModel
130 if (!gasTurbulencePtr_)
136 refCast<const twoPhaseSystem>(liquid.fluid());
151 return *gasTurbulencePtr_;
155 template<
class BasicTurbulenceModel>
159 this->gasTurbulence();
164 *(
mag(this->U_ - gasTurbulence.
U()));
166 this->nut_.correctBoundaryConditions();
169 BasicTurbulenceModel::correctNut();
173 template<
class BasicTurbulenceModel>
177 this->gasTurbulence();
181 refCast<const twoPhaseSystem>(liquid.fluid());
188 Cp_*
sqr(magUr)*fluid.
drag(gas).
K()/liquid.rho()
195 template<
class BasicTurbulenceModel>
207 max(alphaInversion_ - alpha, scalar(0))
211 this->Ce_*
sqrt(gasTurbulence.
k())/this->
delta(),
218 template<
class BasicTurbulenceModel>
225 this->gasTurbulence();
227 const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
231 + phaseTransferCoeff*gasTurbulence.
k()
232 -
fvm::Sp(phaseTransferCoeff, this->k_);
BasicTurbulenceModel::alphaField alphaField
const volVectorField & U() const
Access function to velocity field.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Generic dimensioned Type class.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
BasicTurbulenceModel::rhoField rhoField
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Abstract base class for turbulence models (RAS, LES and laminar).
One equation eddy-viscosity model.
Class which solves the volume fraction equations for two phases.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
const transportModel & transport() const
Access function to incompressible transport model.
static const word propertiesName
Default name of the turbulence properties dictionary.
virtual bool read()
Read model coefficients if they have changed.
A class for handling words, derived from string.
static word groupName(Name name, const word &group)
tmp< volScalarField > bubbleG() const
tmp< volScalarField > phaseTransferCoeff() const
const alphaField & alpha() const
Access function to phase fraction.
virtual tmp< fvScalarMatrix > kSource() const
Templated abstract base class for multiphase compressible turbulence models.
const phaseModel & otherPhase(const phaseModel &phase) const
Constant access the phase not given as an argument.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
One-equation SGS model for the continuous phase in a two-phase system including bubble-generated turb...
const Time & time() const
Return time.
const dragModel & drag(const phaseModel &phase) const
Return the drag model for the given phase.
dimensioned< scalar > mag(const dimensioned< Type > &)
A class for managing temporary objects.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
const objectRegistry & db() const
Return the local objectRegistry.
BasicTurbulenceModel::transportModel transportModel
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual tmp< volScalarField > K(const volScalarField &Ur) const =0
The dragfunction K used in the momentum eq.
virtual void correctNut()
dimensionedScalar deltaT() const
Return time step.