45template<
class BasePhaseModel,
class phaseThermo>
53 BasePhaseModel(
fluid, phaseName),
71 <<
" The selected thermo is pure. Use a multicomponent thermo."
80 thermoPtr_().template getOrDefault<bool>(
"addDiffusion",
false);
115template<
class BasePhaseModel,
class phaseThermo>
131 if (i != inertIndex_)
144 scalarField& xbf = X_[i].boundaryFieldRef()[patchi];
145 const scalarField& ybf =
Y()[i].boundaryField()[patchi];
149 xbf[facei] = Wbf[facei]*ybf[facei]/Wi.
value();
155 X_[inertIndex_] = 1.0 - Xtotal;
160template<
class BasePhaseModel,
class phaseThermo>
165 for(label i=1; i< species_.size(); i++)
167 W += X_[i]*
thermo().composition().W(i);
172 Y()[i] = X_[i]*
thermo().composition().W(i)/W;
174 Info<<
Y()[i].name() <<
" mass fraction = "
175 <<
" Min(Y) = " <<
min(
Y()[i]).value()
176 <<
" Max(Y) = " <<
max(
Y()[i]).value()
182template<
class BasePhaseModel,
class phaseThermo>
190template<
class BasePhaseModel,
class phaseThermo>
198template<
class BasePhaseModel,
class phaseThermo>
201 return thermo().correct();
205template<
class BasePhaseModel,
class phaseThermo>
218 scalar cAlpha(MULEScontrols.
get<scalar>(
"cYi"));
260 if (inertIndex_ != i)
269 "phi" + Yi.
name() +
"Corr",
274 "div(phi," + Yi.
name() +
')'
317 if (phi1p[facei] < 0)
319 phiYiCorrp[facei] = Yip[facei]*phi1p[facei];
350 if (inertIndex_ != i)
372 phiYiCorr += YiEqn.
flux();
374 if (nYiSubCycles > 1)
379 !(++YiSubCycle).end();
430 X_[inertIndex_] = scalar(1) -
Yt;
431 X_[inertIndex_].
max(0.0);
433 calculateMassFractions();
437template<
class BasePhaseModel,
class phaseThermo>
441 return thermoPtr_->composition().Y();
445template<
class BasePhaseModel,
class phaseThermo>
449 return thermoPtr_->composition().Y();
453template<
class BasePhaseModel,
class phaseThermo>
CMULES: Multidimensional universal limiter for explicit corrected implicit solution.
surfaceScalarField phic(mixture.cAlpha() *mag(alphaPhic/mesh.magSf()))
volScalarField Yt(0.0 *Y[0])
const volScalarField & alpha1
const volScalarField & alpha2
const Mesh & mesh() const
Return mesh.
void max(const dimensioned< Type > &dt)
Use the maximum of the field and specified value.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const word & name() const noexcept
Return the object name.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Class which represents a phase with multiple species. Returns the species' mass fractions,...
void calculateMassFractions()
Transfor volume fraction into mass fractions.
virtual void solveYi(PtrList< volScalarField::Internal > &, PtrList< volScalarField::Internal > &)
Solve species fraction equation.
virtual void correct()
Correct phase thermo.
label inertIndex() const
Return inert species index.
autoPtr< phaseThermo > thermoPtr_
Thermophysical model.
PtrList< volScalarField > X_
Ptr list of volumetric fractions for species.
scalar Sct_
Schmidt number.
virtual const phaseThermo & thermo() const
Access to thermo.
void calculateVolumeFractions()
Transfor mass fraction into volume fractions.
label inertIndex_
Inert species index.
bool addDiffusion_
Add diffusion transport on Yi's Eq.
hashedWordList species_
Species table.
virtual const PtrList< volScalarField > & Y() const
Constant access the species mass fractions.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
const T * set(const label i) const
dimensionedScalar deltaT() const
Return time step.
static word timeName(const scalar t, const int precision=precision_)
bool empty() const noexcept
Deprecated(2020-07) True if the managed pointer is null.
void reset(autoPtr< T > &&other) noexcept
Delete managed object and set to new given pointer.
static word phasePropertyName(const word &name, const word &phaseName)
static const word dictName
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
const Type & value() const
Return const reference to value.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
Mesh data needed to do the Finite Volume discretisation.
const Time & time() const
Return the top-level database.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
virtual bool fixesValue() const
Return true if this patch field fixes a value.
Basic first-order Euler implicit/explicit ddt using only the current and previous time-step values.
Basic second-order convection using face-gradients and Gauss' theorem.
virtual bool coupled() const
Return true if this patch field is coupled.
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
label find(const word &val) const
Find index of the value (searches the hash).
A class representing the concept of a field of 1 used to avoid unnecessary manipulations for objects ...
const dictionary & solver(const word &name) const
Return the solver controls dictionary for the given field.
const dictionary & solverDict(const word &name) const
Return the solver controls dictionary for the given field.
Perform a subCycleTime on a field.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
const surfaceScalarField & phi() const
Return the mixture flux.
const fvMesh & mesh() const
Return the mesh.
Upwind differencing scheme class.
A class for handling words, derived from Foam::string.
A class representing the concept of a field of 0 used to avoid unnecessary manipulations for objects ...
basicSpecieMixture & composition
PtrList< volScalarField > & Y
surfaceScalarField phir(fvc::flux(UdmModel.Udm()))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
compressible::turbulenceModel & turbulence
Calculate the substantive (total) derivative.
Calculate the first temporal derivative.
Calculate the divergence of the given field.
Calculate the face-flux of the given field.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for the laplacian of the field.
Calculate the finiteVolume matrix for implicit and explicit sources.
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
void limit(const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phi, surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su, const PsiMaxType &psiMax, const PsiMinType &psiMin, const bool returnCorr)
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
const dimensionSet dimMoles(0, 0, 0, 0, 1, 0, 0)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
multiphaseSystem::phaseModelList & phases
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object.