32 template<
class BasicThermo,
class MixtureType>
36 class PatchFaceMixture,
45 CellMixture cellMixture,
46 PatchFaceMixture patchFaceMixture,
55 IOobject::groupName(psiName, this->
group()),
65 psi[celli] = ((this->*cellMixture)(celli).*psiMethod)(
args[celli] ...);
77 ((this->*patchFaceMixture)(
patchi, facei).*psiMethod)
88 template<
class BasicThermo,
class MixtureType>
89 template<
class CellMixture,
class Method,
class ... Args>
93 CellMixture cellMixture,
108 ((this->*cellMixture)(
cells[celli]).*psiMethod)(
args[celli] ...);
115 template<
class BasicThermo,
class MixtureType>
116 template<
class PatchFaceMixture,
class Method,
class ... Args>
120 PatchFaceMixture patchFaceMixture,
135 ((this->*patchFaceMixture)(
patchi, facei).*psiMethod)
145 template<
class BasicThermo,
class MixtureType>
157 template<
class BasicThermo,
class MixtureType>
165 return psi.primitiveField();
169 template<
class BasicThermo,
class MixtureType>
177 if (isA<gradientEnergyFvPatchScalarField>(hBf[
patchi]))
179 refCast<gradientEnergyFvPatchScalarField>(hBf[
patchi]).gradient()
180 = hBf[
patchi].fvPatchField::snGrad();
182 else if (isA<mixedEnergyFvPatchScalarField>(hBf[
patchi]))
184 refCast<mixedEnergyFvPatchScalarField>(hBf[
patchi]).refGrad()
185 = hBf[
patchi].fvPatchField::snGrad();
193 template<
class BasicThermo,
class MixtureType>
197 const word& phaseName
200 BasicThermo(mesh, phaseName),
201 MixtureType(*this, mesh, phaseName),
207 BasicThermo::phasePropertyName
209 MixtureType::thermoType::heName(),
217 volScalarFieldProperty
221 &MixtureType::cellThermoMixture,
222 &MixtureType::patchFaceThermoMixture,
223 &MixtureType::thermoMixtureType::HE,
227 this->heBoundaryTypes(),
228 this->heBoundaryBaseTypes()
235 BasicThermo::phasePropertyName(
"Cp", phaseName),
247 BasicThermo::phasePropertyName(
"Cv", phaseName),
255 heBoundaryCorrection(he_);
261 template<
class BasicThermo,
class MixtureType>
268 template<
class BasicThermo,
class MixtureType>
272 if (MixtureType::thermoType::enthalpy())
283 template<
class BasicThermo,
class MixtureType>
290 return volScalarFieldProperty
294 &MixtureType::cellThermoMixture,
295 &MixtureType::patchFaceThermoMixture,
296 &MixtureType::thermoMixtureType::HE,
303 template<
class BasicThermo,
class MixtureType>
310 return cellSetProperty
312 &MixtureType::cellThermoMixture,
313 &MixtureType::thermoMixtureType::HE,
315 cellSetScalarList(this->p_,
cells),
321 template<
class BasicThermo,
class MixtureType>
328 return patchFieldProperty
330 &MixtureType::patchFaceThermoMixture,
331 &MixtureType::thermoMixtureType::HE,
333 this->p_.boundaryField()[
patchi],
339 template<
class BasicThermo,
class MixtureType>
343 return volScalarFieldProperty
347 &MixtureType::cellThermoMixture,
348 &MixtureType::patchFaceThermoMixture,
356 template<
class BasicThermo,
class MixtureType>
363 return volScalarFieldProperty
367 &MixtureType::cellThermoMixture,
368 &MixtureType::patchFaceThermoMixture,
376 template<
class BasicThermo,
class MixtureType>
383 return cellSetProperty
385 &MixtureType::cellThermoMixture,
388 cellSetScalarList(this->p_,
cells),
394 template<
class BasicThermo,
class MixtureType>
401 return patchFieldProperty
403 &MixtureType::patchFaceThermoMixture,
406 this->p_.boundaryField()[
patchi],
412 template<
class BasicThermo,
class MixtureType>
416 return volScalarFieldProperty
420 &MixtureType::cellThermoMixture,
421 &MixtureType::patchFaceThermoMixture,
429 template<
class BasicThermo,
class MixtureType>
436 return volScalarFieldProperty
440 &MixtureType::cellThermoMixture,
441 &MixtureType::patchFaceThermoMixture,
449 template<
class BasicThermo,
class MixtureType>
456 return cellSetProperty
458 &MixtureType::cellThermoMixture,
461 cellSetScalarList(this->p_,
cells),
467 template<
class BasicThermo,
class MixtureType>
474 return patchFieldProperty
476 &MixtureType::patchFaceThermoMixture,
479 this->p_.boundaryField()[
patchi],
485 template<
class BasicThermo,
class MixtureType>
489 return volScalarFieldProperty
493 &MixtureType::cellThermoMixture,
494 &MixtureType::patchFaceThermoMixture,
495 &MixtureType::thermoMixtureType::Hf
500 template<
class BasicThermo,
class MixtureType>
507 return patchFieldProperty
509 &MixtureType::patchFaceThermoMixture,
512 this->p_.boundaryField()[
patchi],
518 template<
class BasicThermo,
class MixtureType>
526 return patchFieldProperty
528 &MixtureType::patchFaceThermoMixture,
531 this->p_.boundaryField()[
patchi],
537 template<
class BasicThermo,
class MixtureType>
544 return patchFieldProperty
546 &MixtureType::patchFaceThermoMixture,
547 &MixtureType::thermoMixtureType::gamma,
549 this->p_.boundaryField()[
patchi],
555 template<
class BasicThermo,
class MixtureType>
563 template<
class BasicThermo,
class MixtureType>
570 if (MixtureType::thermoType::enthalpy())
581 template<
class BasicThermo,
class MixtureType>
589 return volScalarFieldProperty
593 &MixtureType::cellThermoMixture,
594 &MixtureType::patchFaceThermoMixture,
595 &MixtureType::thermoMixtureType::THE,
603 template<
class BasicThermo,
class MixtureType>
611 return cellSetProperty
613 &MixtureType::cellThermoMixture,
614 &MixtureType::thermoMixtureType::THE,
617 cellSetScalarList(this->p_,
cells),
623 template<
class BasicThermo,
class MixtureType>
631 return patchFieldProperty
633 &MixtureType::patchFaceThermoMixture,
634 &MixtureType::thermoMixtureType::THE,
637 this->p_.boundaryField()[
patchi],
643 template<
class BasicThermo,
class MixtureType>
647 return volScalarFieldProperty
651 &MixtureType::cellThermoMixture,
652 &MixtureType::patchFaceThermoMixture,
658 template<
class BasicThermo,
class MixtureType>
664 return patchFieldProperty
666 &MixtureType::patchFaceThermoMixture,
673 template<
class BasicThermo,
class MixtureType>
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.
Generic GeometricBoundaryField class.
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
A List with indirect addressing.
Dimension set for the base types.
Mesh data needed to do the Finite Volume discretisation.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
heThermo(const fvMesh &, const word &phaseName)
Construct from mesh.
virtual const volScalarField & Cv() const
Heat capacity at constant volume [J/kg/K].
virtual tmp< volScalarField > W() const
Molecular weight [kg/kmol].
virtual const volScalarField & Cp() const
Heat capacity at constant pressure [J/kg/K].
virtual tmp< volScalarField > ha() const
Absolute enthalpy [J/kg/K].
static UIndirectList< scalar > cellSetScalarList(const volScalarField &psi, const labelList &cells)
Return an indirect list of a field for the given set of cells.
tmp< scalarField > cellSetProperty(CellMixture cellMixture, Method psiMethod, const labelList &cells, const Args &... args) const
Return a scalarField of the given property on a cell set.
virtual ~heThermo()
Destructor.
virtual tmp< volScalarField > hc() const
Enthalpy of formation [J/kg].
void heBoundaryCorrection(volScalarField &he)
Correct the enthalpy/internal energy field boundaries.
tmp< scalarField > patchFieldProperty(PatchFaceMixture patchFaceMixture, Method psiMethod, const label patchi, const Args &... args) const
Return a scalarField of the given property on a patch.
virtual tmp< volScalarField > gamma() const
Gamma = Cp/Cv [].
virtual const volScalarField & Cpv() const
Heat capacity at constant pressure/volume [J/kg/K].
virtual tmp< volScalarField > hs() const
Sensible enthalpy [J/kg/K].
virtual volScalarField & he()
Enthalpy/Internal energy [J/kg].
virtual tmp< volScalarField > THE(const volScalarField &h, const volScalarField &p, const volScalarField &T0) const
Temperature from enthalpy/internal energy.
tmp< volScalarField > volScalarFieldProperty(const word &psiName, const dimensionSet &psiDim, CellMixture cellMixture, PatchFaceMixture patchFaceMixture, Method psiMethod, const Args &... args) const
Return a volScalarField of the given property.
virtual bool read()
Read thermophysical properties dictionary.
A class for managing temporary objects.
T & ref() const
Return non-const reference or generate a fatal error.
A class for handling words, derived from string.
volScalarField scalarField(fieldObject, mesh)
const volScalarField & psi
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
const char *const group
Group name for atomic constants.
const dimensionedScalar h
Planck constant.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const dimensionSet dimEnergy
const dimensionSet dimTemperature
const dimensionSet dimMoles
const dimensionSet dimMass
word name(const complex &)
Return a string representation of a complex.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const scalarList W(::W(thermo))
Foam::argList args(argc, argv)