31 template<
class BasePhaseModel>
35 const tmp<volScalarField>& pressureWork
40 scalar pressureWorkAlphaLimit =
41 this->thermo_->properties()
42 .lookupOrDefault(
"pressureWorkAlphaLimit", 0.0);
44 if (pressureWorkAlphaLimit > 0)
48 max(alpha - pressureWorkAlphaLimit, scalar(0))
49 /
max(alpha - pressureWorkAlphaLimit, pressureWorkAlphaLimit)
61 template<
class BasePhaseModel>
64 const phaseSystem&
fluid,
65 const word& phaseName,
66 const bool referencePhase,
70 BasePhaseModel(fluid, phaseName, referencePhase, index)
76 template<
class BasePhaseModel>
83 template<
class BasePhaseModel>
88 this->thermo_->correct();
92 template<
class BasePhaseModel>
99 template<
class BasePhaseModel>
106 const tmp<volVectorField> tU(this->
U());
109 const tmp<surfaceScalarField> talphaRhoPhi(this->alphaRhoPhi());
112 const tmp<volScalarField> tcontErr(this->continuityError());
115 tmp<volScalarField> tK(this->
K());
120 tmp<fvScalarMatrix> tEEqn
135 if (he.name() == this->thermo_->phasePropertyName(
"e"))
137 tEEqn.ref() += filterPressureWork
147 else if (this->thermo_->dpdt())
149 tEEqn.ref() -= filterPressureWork(alpha*this->
fluid().
dpdt());
virtual void correctThermo()
Correct the thermodynamics.
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
fluidReactionThermo & thermo
virtual bool isothermal() const
Return whether the phase is isothermal.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Class which represents a phase for which the temperature (strictly energy) varies. Returns the energy equation and corrects the thermodynamic model.
tmp< volScalarField > rho() const
Average density.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
GeometricField< vector, fvPatchField, volMesh > volVectorField
CGAL::Exact_predicates_exact_constructions_kernel K
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
AnisothermalPhaseModel(const phaseSystem &fluid, const word &phaseName, const bool referencePhase, const label index)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
const phaseSystem & fluid() const
Return the phase system.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
virtual ~AnisothermalPhaseModel()
Destructor.
A class for managing temporary objects.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
volScalarField::Internal contErr((fvc::ddt(rho)+fvc::div(rhoPhi) -(fvModels.source(alpha1, rho1)&rho1) -(fvModels.source(alpha2, rho2)&rho2))())
virtual tmp< fvScalarMatrix > heEqn()
Return the enthalpy equation.