Go to the documentation of this file.
32 template<
class Thermo,
template<
class>
class Type>
42 template<
class Thermo,
template<
class>
class Type>
57 <<
"Negative initial temperature T0: " <<
T0
63 scalar Ttol =
T0*tol_;
71 (Test - ((this->*F)(
p, Test) -
f)/(this->*dFdT)(
p, Test));
73 if (iter++ > maxIter_)
76 <<
"Maximum number of iterations exceeded: " << maxIter_
80 }
while (
mag(Tnew - Test) > Ttol);
88 template<
class Thermo,
template<
class>
class Type>
101 template<
class Thermo,
template<
class>
class Type>
105 return Type<thermo<Thermo, Type>>::energyName();
109 template<
class Thermo,
template<
class>
class Type>
113 return Type<thermo<Thermo, Type>>::Cpv(*
this,
p,
T);
117 template<
class Thermo,
template<
class>
class Type>
122 volatile const scalar
Cp = this->
Cp(p,
T);
124 const scalar
Cp = this->
Cp(p,
T);
127 return Cp/(
Cp - this->CpMCv(
p,
T));
131 template<
class Thermo,
template<
class>
class Type>
139 return Type<thermo<Thermo, Type>>::CpByCpv(*
this,
p,
T);
143 template<
class Thermo,
template<
class>
class Type>
147 return Type<thermo<Thermo, Type>>::HE(*
this,
p,
T);
151 template<
class Thermo,
template<
class>
class Type>
155 return this->
Ha(p,
T) -
T*this->S(
p,
T);
159 template<
class Thermo,
template<
class>
class Type>
163 return this->
Ea(p,
T) -
T*this->S(
p,
T);
167 template<
class Thermo,
template<
class>
class Type>
171 return this->
Cp(p,
T)*this->W();
175 template<
class Thermo,
template<
class>
class Type>
179 return this->
Ha(p,
T)*this->W();
183 template<
class Thermo,
template<
class>
class Type>
187 return this->
Hs(p,
T)*this->W();
191 template<
class Thermo,
template<
class>
class Type>
195 return this->Hc()*this->W();
199 template<
class Thermo,
template<
class>
class Type>
203 return this->S(
p,
T)*this->W();
207 template<
class Thermo,
template<
class>
class Type>
211 return this->HE(
p,
T)*this->W();
215 template<
class Thermo,
template<
class>
class Type>
219 return this->
Cv(p,
T)*this->W();
223 template<
class Thermo,
template<
class>
class Type>
227 return this->
Es(p,
T)*this->W();
231 template<
class Thermo,
template<
class>
class Type>
235 return this->
Ea(p,
T)*this->W();
239 template<
class Thermo,
template<
class>
class Type>
243 return this->
G(p,
T)*this->W();
247 template<
class Thermo,
template<
class>
class Type>
251 return this->
A(p,
T)*this->W();
255 template<
class Thermo,
template<
class>
class Type>
259 scalar arg = -this->
Y()*this->
G(
Pstd, T)/(
RR*
T);
272 template<
class Thermo,
template<
class>
class Type>
280 template<
class Thermo,
template<
class>
class Type>
284 const scalar nm = this->
Y()/this->W();
286 if (
equal(nm, SMALL))
297 template<
class Thermo,
template<
class>
class Type>
304 const scalar nm = this->
Y()/this->W();
306 if (
equal(nm, SMALL))
317 template<
class Thermo,
template<
class>
class Type>
325 const scalar nm = this->
Y()/this->W();
327 if (
equal(nm, SMALL))
338 template<
class Thermo,
template<
class>
class Type>
346 return Type<thermo<Thermo, Type>>::THE(*
this,
he,
p,
T0);
350 template<
class Thermo,
template<
class>
class Type>
370 template<
class Thermo,
template<
class>
class Type>
390 template<
class Thermo,
template<
class>
class Type>
410 template<
class Thermo,
template<
class>
class Type>
430 template<
class Thermo,
template<
class>
class Type>
438 const scalar nm = this->
Y()/this->W();
440 if (
equal(nm, SMALL))
442 return -this->dGdT(
Pstd,
T)*this->
Y()/
RR;
446 return -(nm/
T + this->dGdT(
Pstd,
T)*this->
Y()/
RR);
451 template<
class Thermo,
template<
class>
class Type>
455 return this->dCpdT(
p,
T)*this->W();
460 template<
class Thermo,
template<
class>
class Type>
466 Thermo::operator+=(st);
470 template<
class Thermo,
template<
class>
class Type>
473 Thermo::operator*=(
s);
479 template<
class Thermo,
template<
class>
class Type>
488 static_cast<const Thermo&
>(st1) +
static_cast<const Thermo&
>(st2)
493 template<
class Thermo,
template<
class>
class Type>
502 s*
static_cast<const Thermo&
>(st)
507 template<
class Thermo,
template<
class>
class Type>
516 static_cast<const Thermo&
>(st1) ==
static_cast<const Thermo&
>(st2)
const scalar RR
Universal gas constant: default in [J/(kmol K)].
scalar ea(const scalar p, const scalar T) const
Absolute internal energy [J/kmol].
A class for handling words, derived from Foam::string.
void limit(const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phi, surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su, const PsiMaxType &psiMax, const PsiMinType &psiMin, const bool returnCorr)
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))
scalar Cpv(const scalar p, const scalar T) const
Heat capacity at constant pressure/volume [J/(kg K)].
const dimensionedScalar G
Newtonian constant of gravitation.
scalar a(const scalar p, const scalar T) const
Helmholtz free energy [J/kmol].
scalar Hs(const scalar p, const scalar T) const
const dimensionedScalar Pstd
Standard pressure.
scalar he(const scalar p, const scalar T) const
Enthalpy/Internal energy [J/kmol].
const volScalarField & Cv
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
volVectorField F(fluid.F())
thermo(const Thermo &sp)
Construct from components.
dimensionedScalar exp(const dimensionedScalar &ds)
scalar THs(const scalar Hs, const scalar p, const scalar T0) const
Temperature from sensible enthalpy given an initial T0.
scalar Kp(const scalar p, const scalar T) const
Equilibrium constant [] i.t.o. partial pressures.
scalar Ea(const scalar p, const scalar T) const
void operator*=(const scalar)
CGAL::Exact_predicates_exact_constructions_kernel K
scalar TEa(const scalar E, const scalar p, const scalar T0) const
Temperature from absolute internal energy.
static word heName()
Name of Enthalpy/Internal energy.
scalar THE(const scalar H, const scalar p, const scalar T0) const
Temperature from enthalpy or internal energy.
scalar K(const scalar p, const scalar T) const
Equilibrium constant [] i.t.o fugacities.
scalar hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kmol].
word name(const complex &c)
Return string representation of complex.
scalar dcpdT(const scalar p, const scalar T) const
Derivative of cp w.r.t. temperature.
scalar hc() const
Chemical enthalpy [J/kmol].
scalar s(const scalar p, const scalar T) const
Entropy [J/(kmol K)].
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
scalar es(const scalar p, const scalar T) const
Sensible internal energy [J/kmol].
scalar CpByCpv(const scalar p, const scalar T) const
Ratio of heat capacity at constant pressure to that at.
errorManip< error > abort(error &err)
PtrList< volScalarField > & Y
scalar G(const scalar p, const scalar T) const
Gibbs free energy [J/kg].
scalar ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kmol].
scalar Kn(const scalar p, const scalar T, const scalar n) const
Equilibrium constant [] i.t.o. number of moles.
scalar cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kmol K)].
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
scalar dKcdTbyKc(const scalar p, const scalar T) const
Derivative of B (acooding to Niemeyer et al.) w.r.t. temperature.
scalar Es(const scalar p, const scalar T) const
scalar A(const scalar p, const scalar T) const
Helmholtz free energy [J/kg].
scalar TEs(const scalar E, const scalar p, const scalar T0) const
Temperature from sensible internal energy.
scalar Kx(const scalar p, const scalar T) const
Equilibrium constant [] i.t.o. mole-fractions.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const volScalarField & Cp
scalar Ha(const scalar p, const scalar T) const
scalar gamma(const scalar p, const scalar T) const
Gamma = Cp/Cv [].
scalar g(const scalar p, const scalar T) const
Gibbs free energy [J/kmol].
scalar HE(const scalar p, const scalar T) const
Enthalpy/Internal energy [J/kg].
scalar Kc(const scalar p, const scalar T) const
Equilibrium constant i.t.o. molar concentration.
scalar THa(const scalar H, const scalar p, const scalar T0) const
Temperature from absolute enthalpy.
scalar cv(const scalar p, const scalar T) const
Heat capacity at constant volume [J/(kmol K)].