30 template<
class Thermo,
template<
class>
class Type>
40 template<
class Thermo,
template<
class>
class Type>
55 <<
"Negative initial temperature T0: " << T0
61 scalar Ttol = T0*tol_;
69 (Test - ((this->*F)(p, Test) - f)/(this->*dFdT)(p, Test));
71 if (iter++ > maxIter_)
74 <<
"Maximum number of iterations exceeded: " << maxIter_
78 }
while (
mag(Tnew - Test) > Ttol);
86 template<
class Thermo,
template<
class>
class Type>
99 template<
class Thermo,
template<
class>
class Type>
103 return Type<thermo<Thermo, Type>>::energyName();
107 template<
class Thermo,
template<
class>
class Type>
111 return this->Cp(p, T) - this->CpMCv(p, T);
115 template<
class Thermo,
template<
class>
class Type>
119 return Type<thermo<Thermo, Type>>::Cpv(*
this, p, T);
123 template<
class Thermo,
template<
class>
class Type>
127 const scalar Cp = this->Cp(p, T);
128 return Cp/(Cp - this->CpMCv(p, T));
132 template<
class Thermo,
template<
class>
class Type>
140 return Type<thermo<Thermo, Type>>::CpByCpv(*
this, p, T);
144 template<
class Thermo,
template<
class>
class Type>
148 return Type<thermo<Thermo, Type>>::HE(*
this, p, T);
152 template<
class Thermo,
template<
class>
class Type>
156 return this->Hs(p, T) - p/this->
rho(p, T);
160 template<
class Thermo,
template<
class>
class Type>
164 return this->Ha(p, T) - p/this->
rho(p, T);
168 template<
class Thermo,
template<
class>
class Type>
172 return this->Ha(p, T) - T*this->S(p, T);
176 template<
class Thermo,
template<
class>
class Type>
180 return this->Ea(p, T) - T*this->S(p, T);
184 template<
class Thermo,
template<
class>
class Type>
188 return this->Cp(p, T)*this->
W();
192 template<
class Thermo,
template<
class>
class Type>
196 return this->Ha(p, T)*this->
W();
200 template<
class Thermo,
template<
class>
class Type>
204 return this->Hs(p, T)*this->
W();
208 template<
class Thermo,
template<
class>
class Type>
212 return this->Hc()*this->
W();
216 template<
class Thermo,
template<
class>
class Type>
220 return this->S(p, T)*this->
W();
224 template<
class Thermo,
template<
class>
class Type>
228 return this->HE(p, T)*this->
W();
232 template<
class Thermo,
template<
class>
class Type>
236 return this->Cv(p, T)*this->
W();
240 template<
class Thermo,
template<
class>
class Type>
244 return this->Es(p, T)*this->
W();
248 template<
class Thermo,
template<
class>
class Type>
252 return this->Ea(p, T)*this->
W();
256 template<
class Thermo,
template<
class>
class Type>
260 return this->
G(p, T)*this->
W();
264 template<
class Thermo,
template<
class>
class Type>
268 return this->A(p, T)*this->
W();
272 template<
class Thermo,
template<
class>
class Type>
276 scalar arg = -this->
Y()*this->
G(
Pstd, T)/(
RR*
T);
289 template<
class Thermo,
template<
class>
class Type>
297 template<
class Thermo,
template<
class>
class Type>
301 const scalar nm = this->
Y()/this->
W();
303 if (
equal(nm, small))
314 template<
class Thermo,
template<
class>
class Type>
321 const scalar nm = this->
Y()/this->
W();
323 if (
equal(nm, small))
329 return Kp(p, T)*
pow(
Pstd/p, nm);
334 template<
class Thermo,
template<
class>
class Type>
342 const scalar nm = this->
Y()/this->
W();
344 if (
equal(nm, small))
350 return Kp(p, T)*
pow(n*
Pstd/p, nm);
355 template<
class Thermo,
template<
class>
class Type>
363 return Type<thermo<Thermo, Type>>::THE(*
this, he, p, T0);
367 template<
class Thermo,
template<
class>
class Type>
387 template<
class Thermo,
template<
class>
class Type>
407 template<
class Thermo,
template<
class>
class Type>
427 template<
class Thermo,
template<
class>
class Type>
447 template<
class Thermo,
template<
class>
class Type>
455 const scalar nm = this->
Y()/this->
W();
457 if (
equal(nm, small))
459 return -this->dGdT(
Pstd, T)*this->
Y()/
RR;
463 return -(nm/T + this->dGdT(
Pstd, T)*this->
Y()/
RR);
468 template<
class Thermo,
template<
class>
class Type>
472 return this->dCpdT(p, T)*this->
W();;
478 template<
class Thermo,
template<
class>
class Type>
479 inline void Foam::species::thermo<Thermo, Type>::operator+=
484 Thermo::operator+=(st);
488 template<
class Thermo,
template<
class>
class Type>
491 Thermo::operator*=(s);
497 template<
class Thermo,
template<
class>
class Type>
506 static_cast<const Thermo&
>(st1) + static_cast<const Thermo&>(st2)
511 template<
class Thermo,
template<
class>
class Type>
520 s*
static_cast<const Thermo&
>(st)
525 template<
class Thermo,
template<
class>
class Type>
534 static_cast<const Thermo&
>(st1) == static_cast<const Thermo&>(st2)
scalar TEs(const scalar E, const scalar p, const scalar T0) const
Temperature from sensible internal energy.
scalar Cv(const scalar p, const scalar T) const
Heat capacity at constant volume [J/(kg K)].
scalar dKcdTbyKc(const scalar p, const scalar T) const
Derivative of B (acooding to Niemeyer et al.) w.r.t. temperature.
scalar G(const scalar p, const scalar T) const
Gibbs free energy [J/kg].
scalar A(const scalar p, const scalar T) const
Helmholtz free energy [J/kg].
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
scalar ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kmol].
scalar he(const scalar p, const scalar T) const
Enthalpy/Internal energy [J/kmol].
const dimensionedScalar G
Newtonian constant of gravitation.
scalar Kn(const scalar p, const scalar T, const scalar n) const
Equilibrium constant [] i.t.o. number of moles.
scalar dcpdT(const scalar p, const scalar T) const
Derivative of cp w.r.t. temperature.
scalar s(const scalar p, const scalar T) const
Entropy [J/(kmol K)].
scalar Kx(const scalar p, const scalar T) const
Equilibrium constant [] i.t.o. mole-fractions.
scalar Es(const scalar p, const scalar T) const
Sensible internal energy [J/kg].
scalar g(const scalar p, const scalar T) const
Gibbs free energy [J/kmol].
CGAL::Exact_predicates_exact_constructions_kernel K
static word heName()
Name of Enthalpy/Internal energy.
scalar CpByCpv(const scalar p, const scalar T) const
Ratio of heat capacity at constant pressure to that at.
scalar THa(const scalar H, const scalar p, const scalar T0) const
Temperature from absolute enthalpy.
scalar HE(const scalar p, const scalar T) const
Enthalpy/Internal energy [J/kg].
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().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))
scalar Ea(const scalar p, const scalar T) const
Absolute internal energy [J/kg].
dimensionedScalar exp(const dimensionedScalar &ds)
A class for handling words, derived from string.
scalar hc() const
Chemical enthalpy [J/kmol].
scalar Cpv(const scalar p, const scalar T) const
Heat capacity at constant pressure/volume [J/(kg K)].
scalar gamma(const scalar p, const scalar T) const
Gamma = Cp/Cv [].
const dimensionedScalar Pstd
Standard pressure.
errorManip< error > abort(error &err)
scalar a(const scalar p, const scalar T) const
Helmholtz free energy [J/kmol].
scalar cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kmol K)].
scalar ea(const scalar p, const scalar T) const
Absolute internal energy [J/kmol].
scalar TEa(const scalar E, const scalar p, const scalar T0) const
Temperature from absolute internal energy.
scalar hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kmol].
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
bool equal(const T &s1, const T &s2)
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
const scalarList W(::W(thermo))
PtrList< volScalarField > & Y
void operator*=(const scalar)
thermo(const Thermo &sp)
Construct from components.
scalar Kc(const scalar p, const scalar T) const
Equilibrium constant i.t.o. molar concentration.
scalar THs(const scalar Hs, const scalar p, const scalar T0) const
Temperature from sensible enthalpy given an initial T0.
dimensioned< scalar > mag(const dimensioned< Type > &)
const scalar RR
Universal gas constant (default in [J/(kmol K)])
scalar cv(const scalar p, const scalar T) const
Heat capacity at constant volume [J/(kmol K)].
scalar es(const scalar p, const scalar T) const
Sensible internal energy [J/kmol].
scalar K(const scalar p, const scalar T) const
Equilibrium constant [] i.t.o fugacities.
scalar THE(const scalar H, const scalar p, const scalar T0) const
Temperature from enthalpy or internal energy.
scalar Kp(const scalar p, const scalar T) const
Equilibrium constant [] i.t.o. partial pressures.