40 template<
class BasePhaseModel>
43 const phaseSystem& fluid,
44 const word& phaseName,
48 BasePhaseModel(fluid, phaseName, index),
53 fluid.subDict(phaseName)
59 fluid.
mesh().solverDict(
"Yi")
65 this->thermo_->lookupOrDefault(
"inertSpecie",
word::null)
70 inertIndex_ = this->thermo_->composition().species()[
inertSpecie];
77 template<
class BasePhaseModel>
84 template<
class BasePhaseModel>
99 PtrList<volScalarField>& Yi =
Y();
103 if (i != inertIndex_)
109 if (inertIndex_ != -1)
111 Yi[inertIndex_] = scalar(1) -
Yt;
112 Yi[inertIndex_].
max(0);
127 template<
class BasePhaseModel>
142 this->thermo_->composition().species()[inertIndex_],
147 !this->thermo_->composition().active
149 this->thermo_->composition().species()[Yi.member()]
155 return tmp<fvScalarMatrix>();
165 +
fvm::div(alphaRhoPhi, Yi,
"div(" + alphaRhoPhi.name() +
",Yi)")
166 -
fvm::Sp(this->continuityError(), Yi)
183 template<
class BasePhaseModel>
187 return this->thermo_->composition().Y();
191 template<
class BasePhaseModel>
195 return this->thermo_->composition().Y();
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
virtual tmp< fvScalarMatrix > YiEqn(volScalarField &Yi)
Return the species fraction equation.
Calculate the matrix for the laplacian of the field.
virtual void correctThermo()
Correct the thermodynamics.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Calculate the first temporal derivative.
virtual const PtrList< volScalarField > & Y() const
Constant access the species mass fractions.
static word groupName(Name name, const word &group)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Calulate the matrix for the first temporal derivative.
static const word null
An empty word.
Calculate the divergence of the given field.
MultiComponentPhaseModel(const phaseSystem &fluid, const word &phaseName, const label index)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Info<< "Reading field U\"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
word name(const complex &)
Return a string representation of a complex.
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 ~MultiComponentPhaseModel()
Destructor.
Calculate the matrix for the divergence of the given field and flux.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
#define R(A, B, C, D, E, F, K, M)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
void max(const dimensioned< Type > &)
PtrList< volScalarField > & Y
A class for managing temporary objects.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const word inertSpecie(thermo.lookup("inertSpecie"))
Calculate the matrix for implicit and explicit sources.
volScalarField Yt(0.0 *Y[0])