Go to the documentation of this file.
31 #include "phaseSystem.H"
43 template<
class BasePhaseModel>
47 const word& phaseName,
51 BasePhaseModel(
fluid, phaseName, index),
56 fluid.subDict(phaseName)
62 fluid.mesh().solverDict(
"Yi")
68 this->thermo_->getOrDefault(
"inertSpecie", word::null)
73 inertIndex_ = this->thermo_->composition().species()[
inertSpecie];
80 if (i != inertIndex_ && this->thermo_->composition().active(i))
82 const label j = YActive_.size();
83 YActive_.resize(j + 1);
84 YActive_.set(j, &
Y[i]);
92 template<
class BasePhaseModel>
99 template<
class BasePhaseModel>
106 IOobject::groupName(
"Yt", this->
name()),
118 if (i != inertIndex_)
124 if (inertIndex_ != -1)
126 Yi[inertIndex_] = scalar(1) -
Yt;
127 Yi[inertIndex_].max(0);
142 template<
class BasePhaseModel>
149 template<
class BasePhaseModel>
160 +
fvm::div(alphaRhoPhi, Yi,
"div(" + alphaRhoPhi.name() +
",Yi)")
177 template<
class BasePhaseModel>
181 return this->thermo_->composition().Y();
185 template<
class BasePhaseModel>
189 return this->thermo_->composition().Y(
name);
193 template<
class BasePhaseModel>
197 return this->thermo_->composition().Y();
201 template<
class BasePhaseModel>
209 template<
class BasePhaseModel>
Defines the attributes of an object for which implicit objectRegistry management is supported,...
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
psiReactionThermo & thermo
Calculate the divergence of the given field.
virtual ~MultiComponentPhaseModel()=default
Destructor.
tmp< GeometricField< Type, faPatchField, areaMesh > > div(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
virtual const PtrList< volScalarField > & Y() const
Constant access the species mass fractions.
Calculate the matrix for the divergence of the given field and flux.
#define forAll(list, i)
Loop across all elements in list.
virtual PtrList< volScalarField > & YRef()
Access the species mass fractions.
MultiComponentPhaseModel(const phaseSystem &fluid, const word &phaseName)
#define R(A, B, C, D, E, F, K, M)
volScalarField Yt(0.0 *Y[0])
Calculate the matrix for the laplacian of the field.
virtual bool pure() const
Return whether the phase is pure (i.e., not multi-component)
A list of pointers to objects of type <T>, without allocation/deallocation management of the pointers...
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
const word inertSpecie(thermo.get< word >("inertSpecie"))
PtrList< volScalarField > & Y
Calculate the matrix for implicit and explicit sources.
virtual const UPtrList< volScalarField > & YActive() const
Return the active species mass fractions.
virtual tmp< fvScalarMatrix > YiEqn(volScalarField &Yi)
Return the species fraction equation.
tmp< GeometricField< Type, faPatchField, areaMesh > > ddt(const dimensioned< Type > dt, const faMesh &mesh)
virtual void correctThermo()
Correct the thermodynamics.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
virtual UPtrList< volScalarField > & YActiveRef()
Access the active species mass fractions.
Calculate the first temporal derivative.
Class to represent a system of phases and model interfacial transfers between them.
tmp< GeometricField< Type, faPatchField, areaMesh > > laplacian(const GeometricField< Type, faPatchField, areaMesh > &vf, const word &name)
static tmp< GeometricField< Type, faePatchField, edgeMesh > > interpolate(const GeometricField< Type, faPatchField, areaMesh > &tvf, const edgeScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Calculate the matrix for the first temporal derivative.
const dimensionSet dimless
Dimensionless.