32 template<
class BasicPsiThermo,
class MixtureType>
40 scalarField& TuCells = this->Tu_.primitiveFieldRef();
41 scalarField& CpCells = this->Cp_.primitiveFieldRef();
42 scalarField& CvCells = this->Cv_.primitiveFieldRef();
43 scalarField& psiCells = this->psi_.primitiveFieldRef();
44 scalarField& muCells = this->mu_.primitiveFieldRef();
45 scalarField& kappaCells = this->kappa_.primitiveFieldRef();
49 const typename MixtureType::thermoMixtureType& thermoMixture =
50 this->cellThermoMixture(celli);
52 const typename MixtureType::transportMixtureType& transportMixture =
53 this->cellTransportMixture(celli, thermoMixture);
55 TCells[celli] = thermoMixture.THE
62 CpCells[celli] = thermoMixture.Cp(pCells[celli], TCells[celli]);
63 CvCells[celli] = thermoMixture.Cv(pCells[celli], TCells[celli]);
64 psiCells[celli] = thermoMixture.psi(pCells[celli], TCells[celli]);
66 muCells[celli] = transportMixture.mu(pCells[celli], TCells[celli]);
68 transportMixture.kappa(pCells[celli], TCells[celli]);
70 TuCells[celli] = this->cellReactants(celli).THE
78 volScalarField::Boundary& pBf =
79 this->p_.boundaryFieldRef();
81 volScalarField::Boundary& TBf =
82 this->T_.boundaryFieldRef();
84 volScalarField::Boundary& TuBf =
85 this->Tu_.boundaryFieldRef();
87 volScalarField::Boundary& CpBf =
88 this->Cp_.boundaryFieldRef();
90 volScalarField::Boundary& CvBf =
91 this->Cv_.boundaryFieldRef();
93 volScalarField::Boundary& psiBf =
94 this->psi_.boundaryFieldRef();
96 volScalarField::Boundary& heBf =
97 this->
he().boundaryFieldRef();
99 volScalarField::Boundary& heuBf =
100 this->heu().boundaryFieldRef();
102 volScalarField::Boundary& muBf =
103 this->mu_.boundaryFieldRef();
105 volScalarField::Boundary& kappaBf =
106 this->kappa_.boundaryFieldRef();
125 const typename MixtureType::thermoMixtureType&
126 thermoMixture = this->patchFaceThermoMixture(
patchi, facei);
128 const typename MixtureType::transportMixtureType&
130 this->patchFaceTransportMixture
131 (
patchi, facei, thermoMixture);
133 phe[facei] = thermoMixture.HE(pp[facei], pT[facei]);
135 pCp[facei] = thermoMixture.Cp(pp[facei], pT[facei]);
136 pCv[facei] = thermoMixture.Cv(pp[facei], pT[facei]);
137 ppsi[facei] = thermoMixture.psi(pp[facei], pT[facei]);
138 pmu[facei] = transportMixture.mu(pp[facei], pT[facei]);
139 pkappa[facei] = transportMixture.kappa(pp[facei], pT[facei]);
146 const typename MixtureType::thermoMixtureType&
147 thermoMixture = this->patchFaceThermoMixture(
patchi, facei);
149 const typename MixtureType::transportMixtureType&
151 this->patchFaceTransportMixture
152 (
patchi, facei, thermoMixture);
154 pT[facei] = thermoMixture.THE(phe[facei], pp[facei], pT[facei]);
156 pCp[facei] = thermoMixture.Cp(pp[facei], pT[facei]);
157 pCv[facei] = thermoMixture.Cv(pp[facei], pT[facei]);
158 ppsi[facei] = thermoMixture.psi(pp[facei], pT[facei]);
159 pmu[facei] = transportMixture.mu(pp[facei], pT[facei]);
160 pkappa[facei] = transportMixture.kappa(pp[facei], pT[facei]);
163 this->patchFaceReactants(
patchi, facei)
164 .THE(pheu[facei], pp[facei], pTu[facei]);
173 template<
class BasicPsiThermo,
class MixtureType>
177 const word& phaseName
180 heThermo<BasicPsiThermo, MixtureType>(mesh, phaseName),
197 MixtureType::thermoType::heName() +
'u',
203 this->volScalarFieldProperty
205 MixtureType::thermoType::heName() +
'u',
207 &MixtureType::cellReactants,
208 &MixtureType::patchFaceReactants,
209 &MixtureType::thermoMixtureType::HE,
213 this->heuBoundaryTypes()
216 this->heuBoundaryCorrection(this->heu_);
220 this->psi_.oldTime();
226 template<
class BasicPsiThermo,
class MixtureType>
233 template<
class BasicPsiThermo,
class MixtureType>
242 this->psi_.oldTime();
253 template<
class BasicPsiThermo,
class MixtureType>
261 return this->cellSetProperty
263 &MixtureType::cellReactants,
264 &MixtureType::thermoMixtureType::HE,
272 template<
class BasicPsiThermo,
class MixtureType>
280 return this->patchFieldProperty
282 &MixtureType::patchFaceReactants,
283 &MixtureType::thermoMixtureType::HE,
285 this->p_.boundaryField()[
patchi],
291 template<
class BasicPsiThermo,
class MixtureType>
295 return this->volScalarFieldProperty
299 &MixtureType::cellProducts,
300 &MixtureType::patchFaceProducts,
301 &MixtureType::thermoMixtureType::THE,
309 template<
class BasicPsiThermo,
class MixtureType>
313 return this->volScalarFieldProperty
316 this->psi_.dimensions(),
317 &MixtureType::cellReactants,
318 &MixtureType::patchFaceReactants,
326 template<
class BasicPsiThermo,
class MixtureType>
332 return this->volScalarFieldProperty
335 this->psi_.dimensions(),
336 &MixtureType::cellProducts,
337 &MixtureType::patchFaceProducts,
345 template<
class BasicPsiThermo,
class MixtureType>
349 return this->volScalarFieldProperty
353 &MixtureType::cellReactants,
354 &MixtureType::patchFaceReactants,
362 template<
class BasicPsiThermo,
class MixtureType>
368 return this->volScalarFieldProperty
372 &MixtureType::cellProducts,
373 &MixtureType::patchFaceProducts,
#define forAll(list, i)
Loop across all elements in list.
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Mesh data needed to do the Finite Volume discretisation.
Enthalpy/Internal energy for a mixture.
virtual tmp< volScalarField > mub() const
Dynamic viscosity of burnt gas [kg/m/s].
virtual volScalarField & heu()
Unburnt gas enthalpy [J/kg].
virtual void correct()
Update properties.
virtual tmp< volScalarField > muu() const
Dynamic viscosity of unburnt gas [kg/m/s].
heheuPsiThermo(const fvMesh &, const word &phaseName)
Construct from mesh and phase name.
virtual tmp< volScalarField > Tb() const
Burnt gas temperature [K].
virtual tmp< volScalarField > psib() const
Burnt gas compressibility [s^2/m^2].
virtual tmp< volScalarField > psiu() const
Unburnt gas compressibility [s^2/m^2].
virtual ~heheuPsiThermo()
Destructor.
A class for managing temporary objects.
A class for handling words, derived from string.
volScalarField scalarField(fieldObject, mesh)
const volScalarField & psi
#define InfoInFunction
Report an information message using Foam::Info.
const dimensionedScalar mu
Atomic mass unit.
label calculate(const fvMesh &mesh, const labelHashSet &patchIDs, const scalar minFaceFraction, GeometricField< scalar, PatchField, GeoMesh > &distance)
Calculate distance data from patches.
const dimensionSet dimDynamicViscosity
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
const dimensionSet dimTemperature
const dimensionSet dimMass
word name(const complex &)
Return a string representation of a complex.
fvPatchField< scalar > fvPatchScalarField