31 template<
class EquationOfState>
34 const EquationOfState& st,
38 const typename janafThermo<EquationOfState>::coeffArray& highCpCoeffs,
39 const typename janafThermo<EquationOfState>::coeffArray& lowCpCoeffs,
40 const bool convertCoeffs
52 highCpCoeffs_[coefLabel] = highCpCoeffs[coefLabel]*this->
R();
53 lowCpCoeffs_[coefLabel] = lowCpCoeffs[coefLabel]*this->
R();
60 highCpCoeffs_[coefLabel] = highCpCoeffs[coefLabel];
61 lowCpCoeffs_[coefLabel] = lowCpCoeffs[coefLabel];
67 template<
class EquationOfState>
87 template<
class EquationOfState>
94 EquationOfState(name, jt),
101 highCpCoeffs_[coefLabel] = jt.highCpCoeffs_[coefLabel];
102 lowCpCoeffs_[coefLabel] = jt.lowCpCoeffs_[coefLabel];
109 template<
class EquationOfState>
115 if (T < Tlow_ || T > Thigh_)
118 <<
"attempt to use janafThermo<EquationOfState>" 119 " out of temperature range " 120 << Tlow_ <<
" -> " << Thigh_ <<
"; T = " << T
123 return min(
max(T, Tlow_), Thigh_);
132 template<
class EquationOfState>
139 template<
class EquationOfState>
146 template<
class EquationOfState>
153 template<
class EquationOfState>
157 return highCpCoeffs_;
161 template<
class EquationOfState>
169 template<
class EquationOfState>
178 ((((a[4]*T + a[3])*T + a[2])*T + a[1])*T + a[0])
183 template<
class EquationOfState>
193 ((((a[4]/5.0*T + a[3]/4.0)*T + a[2]/3.0)*T + a[1]/2.0)*T + a[0])*T
195 ) + EquationOfState::H(p, T);
199 template<
class EquationOfState>
206 return Ha(p, T) -
Hf();
210 template<
class EquationOfState>
224 template<
class EquationOfState>
234 (((a[4]/4.0*T + a[3]/3.0)*T + a[2]/2.0)*T + a[1])*T + a[0]*
log(T)
240 template<
class EquationOfState>
251 - (((a[4]/20.0*T + a[3]/12.0)*T + a[2]/6.0)*T + a[1]/2.0)*T
259 template<
class EquationOfState>
267 return (((4*a[4]*T + 3*a[3])*T + 2*a[2])*T + a[1]);
273 template<
class EquationOfState>
274 inline void Foam::janafThermo<EquationOfState>::operator+=
279 scalar Y1 = this->
Y();
281 EquationOfState::operator+=(jt);
283 if (
mag(this->
Y()) > small)
286 const scalar Y2 = jt.Y()/this->
Y();
288 Tlow_ =
max(Tlow_, jt.Tlow_);
289 Thigh_ =
min(Thigh_, jt.Thigh_);
298 <<
"Tcommon " << Tcommon_ <<
" for " 299 << (this->
name().size() ? this->
name() :
"others")
300 <<
" != " << jt.Tcommon_ <<
" for " 301 << (jt.name().size() ? jt.name() :
"others")
308 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
312 highCpCoeffs_[coefLabel] =
313 Y1*highCpCoeffs_[coefLabel]
314 + Y2*jt.highCpCoeffs_[coefLabel];
316 lowCpCoeffs_[coefLabel] =
317 Y1*lowCpCoeffs_[coefLabel]
318 + Y2*jt.lowCpCoeffs_[coefLabel];
326 template<
class EquationOfState>
333 EquationOfState eofs = jt1;
336 if (
mag(eofs.Y()) < small)
350 const scalar Y1 = jt1.Y()/eofs.Y();
351 const scalar Y2 = jt2.Y()/eofs.Y();
359 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
363 highCpCoeffs[coefLabel] =
364 Y1*jt1.highCpCoeffs_[coefLabel]
365 + Y2*jt2.highCpCoeffs_[coefLabel];
367 lowCpCoeffs[coefLabel] =
368 Y1*jt1.lowCpCoeffs_[coefLabel]
369 + Y2*jt2.lowCpCoeffs_[coefLabel];
375 &&
notEqual(jt1.Tcommon_, jt2.Tcommon_)
379 <<
"Tcommon " << jt1.Tcommon_ <<
" for " 380 << (jt1.name().size() ? jt1.name() :
"others")
381 <<
" != " << jt2.Tcommon_ <<
" for " 382 << (jt2.name().size() ? jt2.name() :
"others")
389 max(jt1.Tlow_, jt2.Tlow_),
390 min(jt1.Thigh_, jt2.Thigh_),
399 template<
class EquationOfState>
408 s*
static_cast<const EquationOfState&
>(jt),
418 template<
class EquationOfState>
427 static_cast<const EquationOfState&>(jt1)
428 == static_cast<const EquationOfState&>(jt2)
431 const scalar Y1 = jt2.Y()/eofs.Y();
432 const scalar Y2 = jt1.Y()/eofs.Y();
440 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
444 highCpCoeffs[coefLabel] =
445 Y1*jt2.highCpCoeffs_[coefLabel]
446 - Y2*jt1.highCpCoeffs_[coefLabel];
448 lowCpCoeffs[coefLabel] =
449 Y1*jt2.lowCpCoeffs_[coefLabel]
450 - Y2*jt1.lowCpCoeffs_[coefLabel];
456 &&
notEqual(jt2.Tcommon_, jt1.Tcommon_)
460 <<
"Tcommon " << jt2.Tcommon_ <<
" for " 461 << (jt2.name().size() ? jt2.name() :
"others")
462 <<
" != " << jt1.Tcommon_ <<
" for " 463 << (jt1.name().size() ? jt1.name() :
"others")
470 max(jt2.Tlow_, jt1.Tlow_),
471 min(jt2.Thigh_, jt1.Thigh_),
bool notEqual(const Scalar s1, const Scalar s2)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
static const int nCoeffs_
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
scalar Thigh() const
Return const access to the high temperature limit.
const dimensionedScalar & Tstd
Standard temperature.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const coeffArray & highCpCoeffs() const
Return const access to the high temperature poly coefficients.
scalar Gstd(const scalar T) const
Gibbs free energy of the mixture in the standard state [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 dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
A class for handling words, derived from string.
scalar S(const scalar p, const scalar T) const
Entropy [J/kg/K].
scalar Tlow() const
Return const access to the low temperature limit.
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
janafThermo(const EquationOfState &st, const scalar Tlow, const scalar Thigh, const scalar Tcommon, const coeffArray &highCpCoeffs, const coeffArray &lowCpCoeffs, const bool convertCoeffs=false)
Construct from components.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
scalar Hf() const
Enthalpy of formation [J/kg].
#define R(A, B, C, D, E, F, K, M)
#define WarningInFunction
Report a warning using Foam::Warning.
PtrList< volScalarField > & Y
const coeffArray & lowCpCoeffs() const
Return const access to the low temperature poly coefficients.
JANAF tables based thermodynamics package templated into the equation of state.
scalar Cp(const scalar p, const scalar T) const
scalar Ha(const scalar p, const scalar T) const
Absolute enthalpy [J/kg].
dimensioned< scalar > mag(const dimensioned< Type > &)
scalar Cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/kg/K].
scalar Tcommon() const
Return const access to the common temperature.