Go to the documentation of this file.
30 template<
class EquationOfState>
33 const EquationOfState& st,
50 template<
class EquationOfState>
57 EquationOfState(
name, ct),
65 template<
class EquationOfState>
73 template<
class EquationOfState>
83 template<
class EquationOfState>
93 template<
class EquationOfState>
104 template<
class EquationOfState>
107 const scalar
p,
const scalar
T
110 return Hs(
p,
T) + Hc();
114 template<
class EquationOfState>
117 const scalar
p,
const scalar
T
124 template<
class EquationOfState>
131 template<
class EquationOfState>
134 const scalar
p,
const scalar
T
137 return Cp_*
log(
T/
Tstd) + EquationOfState::S(
p,
T);
141 template<
class EquationOfState>
144 const scalar
p,
const scalar
T
151 template<
class EquationOfState>
154 const scalar
p,
const scalar
T
163 template<
class EquationOfState>
169 scalar Y1 = this->
Y();
171 EquationOfState::operator+=(ct);
173 if (
mag(this->
Y()) > SMALL)
176 const scalar Y2 = ct.Y()/this->
Y();
178 Cp_ = Y1*Cp_ + Y2*ct.Cp_;
179 Hf_ = Y1*Hf_ + Y2*ct.Hf_;
186 template<
class EquationOfState>
195 static_cast<const EquationOfState&
>(ct1)
196 +
static_cast<const EquationOfState&
>(ct2)
199 if (
mag(eofs.Y()) < SMALL)
212 return hRefConstThermo<EquationOfState>
215 ct1.Y()/eofs.Y()*ct1.Cp_
216 + ct2.Y()/eofs.Y()*ct2.Cp_,
217 ct1.Y()/eofs.Y()*ct1.Hf_
218 + ct2.Y()/eofs.Y()*ct2.Hf_,
219 ct1.Y()/eofs.Y()*ct1.Tref_
220 + ct2.Y()/eofs.Y()*ct2.Tref_,
221 ct1.Y()/eofs.Y()*ct1.Href_
222 + ct2.Y()/eofs.Y()*ct2.Href_
228 template<
class EquationOfState>
232 const hRefConstThermo<EquationOfState>& ct
235 return hRefConstThermo<EquationOfState>
237 s*
static_cast<const EquationOfState&
>(ct),
246 template<
class EquationOfState>
249 const hRefConstThermo<EquationOfState>& ct1,
250 const hRefConstThermo<EquationOfState>& ct2
255 static_cast<const EquationOfState&
>(ct1)
256 ==
static_cast<const EquationOfState&
>(ct2)
259 return hRefConstThermo<EquationOfState>
262 ct2.Y()/eofs.Y()*ct2.Cp_
263 - ct1.Y()/eofs.Y()*ct1.Cp_,
264 ct2.Y()/eofs.Y()*ct2.Hf_
265 - ct1.Y()/eofs.Y()*ct1.Hf_
scalar Ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kg].
A class for handling words, derived from Foam::string.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const volScalarField & cp
scalar Hs(const scalar p, const scalar T) const
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
autoPtr< hRefConstThermo > clone() const
Construct and return a clone.
scalar S(const scalar p, const scalar T) const
Entropy [J/(kg K)].
const dimensionedScalar Tstd
Standard temperature.
scalar Cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kg K)].
word name(const complex &c)
Return string representation of complex.
scalar Hc() const
Chemical enthalpy [J/kg].
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar dGdT(const scalar p, const scalar T) const
Derivative of Gibbs free energy w.r.t. temperature.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
dimensionedScalar log(const dimensionedScalar &ds)
static autoPtr< hRefConstThermo > New(const dictionary &dict)
Selector from dictionary.
PtrList< volScalarField > & Y
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const volScalarField & Cp
Constant properties thermodynamics package templated into the EquationOfState.