31 template<
class EquationOfState>
40 <<
"attempt to evaluate hPowerThermo<EquationOfState>" 41 " for negative temperature " << T
47 template<
class EquationOfState>
54 EquationOfState(name, jt),
64 template<
class EquationOfState>
67 const EquationOfState& st,
82 template<
class EquationOfState>
93 template<
class EquationOfState>
106 template<
class EquationOfState>
116 template<
class EquationOfState>
119 const scalar p,
const scalar T
122 return c0_*
pow(T/Tref_, n0_) + EquationOfState::Cp(p, T);
126 template<
class EquationOfState>
129 const scalar p,
const scalar T
132 return Hs(p, T) + Hc();
136 template<
class EquationOfState>
139 const scalar p,
const scalar T
143 c0_*(
pow(T, n0_ + 1) -
pow(
Tstd, n0_ + 1))/(
pow(Tref_, n0_)*(n0_ + 1))
144 + EquationOfState::H(p, T);
148 template<
class EquationOfState>
155 template<
class EquationOfState>
158 const scalar p,
const scalar T
163 + EquationOfState::S(p, T);
167 template<
class EquationOfState>
170 const scalar p,
const scalar T
179 template<
class EquationOfState>
182 const scalar p,
const scalar T
193 template<
class EquationOfState>
194 inline void Foam::hPowerThermo<EquationOfState>::operator+=
199 scalar Y1 = this->
Y();
201 EquationOfState::operator+=(ct);
203 if (
mag(this->
Y()) > small)
206 const scalar Y2 = ct.Y()/this->
Y();
208 Hf_ = Y1*Hf_ + Y2*ct.Hf_;
209 c0_ = Y1*c0_ + Y2*ct.c0_;
210 n0_ = Y1*n0_ + Y2*ct.n0_;
211 Tref_ = Y1*Tref_ + Y2*ct.Tref_;
218 template<
class EquationOfState>
227 static_cast<const EquationOfState&>(ct1)
228 + static_cast<const EquationOfState&>(ct2)
231 if (
mag(eofs.Y()) < small)
247 ct1.Y()/eofs.Y()*ct1.c0_
248 + ct2.Y()/eofs.Y()*ct2.c0_,
249 ct1.Y()/eofs.Y()*ct1.n0_
250 + ct2.Y()/eofs.Y()*ct2.n0_,
251 ct1.Y()/eofs.Y()*ct1.Tref_
252 + ct2.Y()/eofs.Y()*ct2.Tref_,
253 ct1.Y()/eofs.Y()*ct1.Hf_
254 + ct2.Y()/eofs.Y()*ct2.Hf_
260 template<
class EquationOfState>
269 s*
static_cast<const EquationOfState&
>(ct),
278 template<
class EquationOfState>
287 static_cast<const EquationOfState&>(ct1)
288 == static_cast<const EquationOfState&>(ct2)
294 ct2.Y()/eofs.Y()*ct2.c0_
295 - ct1.Y()/eofs.Y()*ct1.c0_,
296 ct2.Y()/eofs.Y()*ct2.n0_
297 - ct1.Y()/eofs.Y()*ct1.n0_,
298 ct2.Y()/eofs.Y()*ct2.Tref_
299 - ct1.Y()/eofs.Y()*ct1.Tref_,
300 ct2.Y()/eofs.Y()*ct2.Hf_
301 - ct1.Y()/eofs.Y()*ct1.Hf_
A list of keyword definitions, which are a keyword followed by any number of values (e...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
scalar Ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kg].
scalar S(const scalar p, const scalar T) const
Entropy [J/(kg K)].
autoPtr< hPowerThermo > clone() const
Construct and return a clone.
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 dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
A class for handling words, derived from string.
errorManip< error > abort(error &err)
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
const dimensionedScalar Tstd
Standard temperature.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
scalar Hc() const
Chemical enthalpy [J/kg].
static autoPtr< hPowerThermo > New(const dictionary &dict)
Selector from dictionary.
PtrList< volScalarField > & Y
scalar Cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kg K)].
Power-function based thermodynamics package templated on EquationOfState.
scalar dGdT(const scalar p, const scalar T) const
Derivative of Gibbs free energy w.r.t. temperature.
dimensioned< scalar > mag(const dimensioned< Type > &)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.