32 template<
class BasicThermo,
class MixtureType>
40 if (isA<gradientEnergyFvPatchScalarField>(hBf[patchi]))
42 refCast<gradientEnergyFvPatchScalarField>(hBf[
patchi]).gradient()
43 = hBf[
patchi].fvPatchField::snGrad();
45 else if (isA<mixedEnergyFvPatchScalarField>(hBf[patchi]))
47 refCast<mixedEnergyFvPatchScalarField>(hBf[
patchi]).refGrad()
48 = hBf[
patchi].fvPatchField::snGrad();
54 template<
class BasicThermo,
class MixtureType>
64 this->cellMixture(celli).HE(pCells[celli], TCells[celli]);
67 volScalarField::Boundary& heBf = he_.boundaryFieldRef();
73 this->p_.boundaryField()[
patchi],
74 this->T_.boundaryField()[
patchi],
79 this->heBoundaryCorrection(he_);
86 template<
class BasicThermo,
class MixtureType>
93 BasicThermo(mesh, phaseName),
94 MixtureType(*
this, mesh, phaseName),
100 BasicThermo::phasePropertyName
102 MixtureType::thermoType::heName()
111 this->heBoundaryTypes(),
112 this->heBoundaryBaseTypes()
119 template<
class BasicThermo,
class MixtureType>
124 const word& phaseName
127 BasicThermo(mesh, dict, phaseName),
128 MixtureType(*
this, mesh, phaseName),
134 BasicThermo::phasePropertyName
136 MixtureType::thermoType::heName()
145 this->heBoundaryTypes(),
146 this->heBoundaryBaseTypes()
155 template<
class BasicThermo,
class MixtureType>
162 template<
class BasicThermo,
class MixtureType>
169 const fvMesh& mesh = this->T_.mesh();
189 this->cellMixture(celli).HE(pCells[celli], TCells[celli]);
203 this->patchFaceMixture(patchi, facei).HE(pp[facei], Tp[facei]);
211 template<
class BasicThermo,
class MixtureType>
224 he[celli] = this->cellMixture(cells[celli]).HE(p[celli], T[celli]);
231 template<
class BasicThermo,
class MixtureType>
245 this->patchFaceMixture(patchi, facei).HE(p[facei], T[facei]);
252 template<
class BasicThermo,
class MixtureType>
256 const fvMesh& mesh = this->T_.mesh();
273 hcCells[celli] = this->cellMixture(celli).Hc();
284 hcp[facei] = this->patchFaceMixture(patchi, facei).Hc();
292 template<
class BasicThermo,
class MixtureType>
306 this->patchFaceMixture(patchi, facei).Cp(p[facei], T[facei]);
313 template<
class BasicThermo,
class MixtureType>
317 const fvMesh& mesh = this->T_.mesh();
334 this->cellMixture(celli).Cp(this->p_[celli], this->T_[celli]);
348 this->patchFaceMixture(patchi, facei).Cp(pp[facei], pT[facei]);
356 template<
class BasicThermo,
class MixtureType>
371 this->patchFaceMixture(patchi, facei).Cv(p[facei], T[facei]);
378 template<
class BasicThermo,
class MixtureType>
382 const fvMesh& mesh = this->T_.mesh();
399 this->cellMixture(celli).Cv(this->p_[celli], this->T_[celli]);
408 this->p_.boundaryField()[
patchi],
409 this->T_.boundaryField()[
patchi],
418 template<
class BasicThermo,
class MixtureType>
432 this->patchFaceMixture(patchi, facei).gamma(p[facei], T[facei]);
439 template<
class BasicThermo,
class MixtureType>
443 const fvMesh& mesh = this->T_.mesh();
460 this->cellMixture(celli).gamma(this->p_[celli], this->T_[celli]);
473 pgamma[facei] = this->patchFaceMixture(patchi, facei).gamma
485 template<
class BasicThermo,
class MixtureType>
499 this->patchFaceMixture(patchi, facei).Cpv(p[facei], T[facei]);
506 template<
class BasicThermo,
class MixtureType>
510 const fvMesh& mesh = this->T_.mesh();
527 this->cellMixture(celli).Cpv(this->p_[celli], this->T_[celli]);
541 this->patchFaceMixture(patchi, facei).Cpv(pp[facei], pT[facei]);
549 template<
class BasicThermo,
class MixtureType>
563 this->patchFaceMixture(patchi, facei).CpByCpv(p[facei], T[facei]);
570 template<
class BasicThermo,
class MixtureType>
574 const fvMesh& mesh = this->T_.mesh();
590 CpByCpv[celli] = this->cellMixture(celli).CpByCpv
597 volScalarField::Boundary& CpByCpvBf =
608 pCpByCpv[facei] = this->patchFaceMixture(patchi, facei).CpByCpv
620 template<
class BasicThermo,
class MixtureType>
635 this->cellMixture(cells[celli]).THE(h[celli], p[celli], T0[celli]);
642 template<
class BasicThermo,
class MixtureType>
656 T[facei] = this->patchFaceMixture
660 ).THE(h[facei], p[facei], T0[facei]);
667 template<
class BasicThermo,
class MixtureType>
672 const fvMesh& mesh = this->T_.mesh();
689 WCells[celli] = this->cellMixture(celli).W();
699 Wp[facei] = this->patchFaceMixture(patchi, facei).W();
707 template<
class BasicThermo,
class MixtureType>
713 const fvMesh& mesh = this->T_.mesh();
719 W[facei] = this->patchFaceMixture(patchi, facei).W();
726 template<
class BasicThermo,
class MixtureType>
736 template<
class BasicThermo,
class MixtureType>
745 this->p_.boundaryField()[
patchi],
746 this->T_.boundaryField()[
patchi],
748 )*this->alpha_.boundaryField()[
patchi];
752 template<
class BasicThermo,
class MixtureType>
762 template<
class BasicThermo,
class MixtureType>
769 this->p_.boundaryField()[
patchi],
770 this->T_.boundaryField()[
patchi],
773 *this->alpha_.boundaryField()[
patchi];
777 template<
class BasicThermo,
class MixtureType>
790 template<
class BasicThermo,
class MixtureType>
801 this->p_.boundaryField()[
patchi],
802 this->T_.boundaryField()[
patchi],
806 this->alpha_.boundaryField()[
patchi]
812 template<
class BasicThermo,
class MixtureType>
825 template<
class BasicThermo,
class MixtureType>
836 this->p_.boundaryField()[
patchi],
837 this->T_.boundaryField()[
patchi],
841 this->alpha_.boundaryField()[
patchi]
847 template<
class BasicThermo,
class MixtureType>
virtual tmp< volScalarField > kappa() const
Thermal diffusivity for temperature of mixture [W/m/K].
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual tmp< scalarField > THE(const scalarField &he, const scalarField &p, const scalarField &T0, const labelList &cells) const
Temperature from enthalpy/internal energy for cell-set.
scalar Cv(const scalar p, const scalar T) const
#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< volScalarField > gamma() const
Gamma = Cp/Cv [].
virtual tmp< volScalarField > CpByCpv() const
Heat capacity ratio [].
A list of keyword definitions, which are a keyword followed by any number of values (e...
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
virtual ~heThermo()
Destructor.
void size(const label)
Override size to be inconsistent with allocated storage.
virtual tmp< volScalarField > Cpv() const
Heat capacity at constant pressure/volume [J/kg/K].
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
const Time & time() const
Return the top-level database.
virtual volScalarField & he()
Enthalpy/Internal energy [J/kg].
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
virtual tmp< volScalarField > hc() const
Chemical enthalpy [J/kg].
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
virtual tmp< volScalarField > Cp() const
Heat capacity at constant pressure [J/kg/K].
A class for handling words, derived from string.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
virtual tmp< volScalarField > alphahe() const
Thermal diffusivity for energy of mixture [kg/m/s].
volScalarField scalarField(fieldObject, mesh)
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
virtual tmp< volScalarField > alphaEff(const volScalarField &alphat) const
Effective thermal turbulent diffusivity of mixture [kg/m/s].
const volScalarField & cp
const dimensionSet dimMoles(0, 0, 0, 0, 1, 0, 0)
virtual void rename(const word &newName)
Rename.
const dimensionSet dimEnergy
label size() const
Return the number of elements in the UPtrList.
virtual bool read()
Read thermophysical properties dictionary.
virtual tmp< volScalarField > kappaEff(const volScalarField &) const
Effective thermal turbulent diffusivity for temperature.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
virtual tmp< volScalarField > Cv() const
Heat capacity at constant volume [J/kg/K].
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
Enthalpy/Internal energy for a mixture.
const scalarList W(::W(thermo))
Mesh data needed to do the Finite Volume discretisation.
scalar Cp(const scalar p, const scalar T) const
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
A class for managing temporary objects.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
void heBoundaryCorrection(volScalarField &he)
Correct the enthalpy/internal energy field boundaries.
virtual tmp< volScalarField > W() const
Molecular weight [kg/kmol].