30 template<
class ThermoType>
44 thermoMixture_(this->specieThermos()),
45 transportMixture_(this->specieThermos())
51 template<
class ThermoType>
62 template<
class ThermoType>
63 template<
class Method,
class ... Args>
75 psi += Y_[i]*(specieThermos_[i].*psiMethod)(
args ...);
82 template<
class ThermoType>
83 template<
class Method,
class ... Args>
96 rPsi += Y_[i]/(specieThermos_[i].*psiMethod)(
args ...);
103 template<
class ThermoType>
104 template<
class Method,
class ... Args>
116 psi += X_[i]*(specieThermos_[i].*psiMethod)(
args ...);
123 template<
class ThermoType>
131 template<
class ThermoType>
142 template<
class ThermoType>
150 scalar psiByRho2 = 0;
154 const scalar rhoi = specieThermos_[i].rho(
p,
T);
155 const scalar psii = specieThermos_[i].psi(
p,
T);
157 oneByRho += Y_[i]/rhoi;
161 psiByRho2 += Y_[i]*psii/
sqr(rhoi);
165 return psiByRho2/
sqr(oneByRho);
169 template<
class ThermoType>
173 return massWeighted(&ThermoType::Hf);
177 #define thermoMixtureFunction(Func) \
178 template<class ThermoType> \
180 Foam::valueMulticomponentMixture<ThermoType>::thermoMixture::Func \
186 return massWeighted(&ThermoType::Func, p, T); \
198 template<class ThermoType>
219 template<
class ThermoType>
231 template<
class ThermoType>
243 template<
class ThermoType>
255 Y[i] = this->
Y()[i][celli];
258 return thermoMixture_;
262 template<
class ThermoType>
275 Y[i] = this->
Y()[i].boundaryField()[
patchi][facei];
278 return thermoMixture_;
282 template<
class ThermoType>
296 X[i] = this->
Y()[i][celli]/this->specieThermos()[i].W();
305 return transportMixture_;
309 template<
class ThermoType>
325 this->
Y()[i].boundaryField()[
patchi][facei]
326 /this->specieThermos()[i].W();
335 return transportMixture_;
scalar Ha(const scalar p, const scalar T) const
scalar Hs(const scalar p, const scalar T) const
scalar Cp(const scalar p, const scalar T) const
scalar Cv(const scalar p, const scalar T) const
#define forAll(list, i)
Loop across all elements in list.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Mesh data needed to do the Finite Volume discretisation.
Foam::multicomponentMixture.
scalar psi(scalar p, scalar T) const
Return compressibility [s^2/m^2].
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
scalar W() const
Molecular weight [kg/kmol].
scalar mu(const scalar p, const scalar T) const
Dynamic viscosity [kg/m/s].
scalar kappa(const scalar p, const scalar T) const
Thermal conductivity [W/m/K].
Thermophysical properties mixing class which applies mass-fraction weighted mixing to thermodynamic p...
const transportMixtureType & patchFaceTransportMixture(const label patchi, const label facei) const
const thermoMixtureType & cellThermoMixture(const label celli) const
const transportMixtureType & cellTransportMixture(const label celli) const
valueMulticomponentMixture(const dictionary &, const fvMesh &, const word &)
Construct from dictionary, mesh and phase name.
const thermoMixtureType & patchFaceThermoMixture(const label patchi, const label facei) const
A class for handling words, derived from string.
const volScalarField & psi
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)
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
const dimensionedScalar mu
Atomic mass unit.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const scalarList W(::W(thermo))
Foam::argList args(argc, argv)
PtrList< volScalarField > & Y
#define thermoMixtureFunction(Func)