Go to the documentation of this file.
33 template<
class EquationOfState>
36 const EquationOfState& st,
40 const typename janafThermo<EquationOfState>::coeffArray& highCpCoeffs,
41 const typename janafThermo<EquationOfState>::coeffArray& lowCpCoeffs,
42 const bool convertCoeffs
52 for (label coefLabel=0; coefLabel<
nCoeffs_; coefLabel++)
55 lowCpCoeffs_[coefLabel] =
lowCpCoeffs[coefLabel]*this->
R();
60 for (label coefLabel=0; coefLabel<
nCoeffs_; coefLabel++)
69 template<
class EquationOfState>
89 template<
class EquationOfState>
96 EquationOfState(
name, jt),
101 for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
103 highCpCoeffs_[coefLabel] = jt.highCpCoeffs_[coefLabel];
104 lowCpCoeffs_[coefLabel] = jt.lowCpCoeffs_[coefLabel];
111 template<
class EquationOfState>
117 if (T < Tlow_ || T > Thigh_)
120 <<
"attempt to use janafThermo<EquationOfState>"
121 " out of temperature range "
122 << Tlow_ <<
" -> " << Thigh_ <<
"; T = " <<
T
125 return min(
max(
T, Tlow_), Thigh_);
134 template<
class EquationOfState>
141 template<
class EquationOfState>
148 template<
class EquationOfState>
155 template<
class EquationOfState>
159 return highCpCoeffs_;
163 template<
class EquationOfState>
171 template<
class EquationOfState>
180 ((((a[4]*
T + a[3])*
T + a[2])*
T + a[1])*
T + a[0])
185 template<
class EquationOfState>
195 ((((a[4]/5.0*
T + a[3]/4.0)*
T + a[2]/3.0)*
T + a[1]/2.0)*
T + a[0])*
T
201 template<
class EquationOfState>
208 return Ha(
p,
T) - Hc();
212 template<
class EquationOfState>
226 template<
class EquationOfState>
236 (((a[4]/4.0*
T + a[3]/3.0)*
T + a[2]/2.0)*
T + a[1])*
T + a[0]*
log(
T)
238 ) + EquationOfState::S(
p,
T);
242 template<
class EquationOfState>
250 return -((a[0] + a[5]/
T)/
T + a[1]/2 +
T*(a[2]/3 +
T*(a[3]/4 +
T*a[4]/5)));
254 template<
class EquationOfState>
263 (((4*a[4]*
T + 3*a[3])*
T + 2*a[2])*
T + a[1]);
268 template<
class EquationOfState>
274 scalar Y1 = this->
Y();
276 EquationOfState::operator+=(jt);
278 if (
mag(this->
Y()) > SMALL)
281 const scalar Y2 = jt.Y()/this->
Y();
283 Tlow_ =
max(Tlow_, jt.Tlow_);
284 Thigh_ =
min(Thigh_, jt.Thigh_);
293 <<
"Tcommon " << Tcommon_ <<
" for "
294 << (this->
name().size() ? this->
name() :
"others")
295 <<
" != " << jt.Tcommon_ <<
" for "
296 << (jt.name().size() ? jt.name() :
"others")
303 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
307 highCpCoeffs_[coefLabel] =
308 Y1*highCpCoeffs_[coefLabel]
309 + Y2*jt.highCpCoeffs_[coefLabel];
311 lowCpCoeffs_[coefLabel] =
312 Y1*lowCpCoeffs_[coefLabel]
313 + Y2*jt.lowCpCoeffs_[coefLabel];
321 template<
class EquationOfState>
328 EquationOfState eofs = jt1;
331 if (
mag(eofs.Y()) < SMALL)
345 const scalar Y1 = jt1.Y()/eofs.Y();
346 const scalar Y2 = jt2.Y()/eofs.Y();
354 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
358 highCpCoeffs[coefLabel] =
359 Y1*jt1.highCpCoeffs_[coefLabel]
360 + Y2*jt2.highCpCoeffs_[coefLabel];
362 lowCpCoeffs[coefLabel] =
363 Y1*jt1.lowCpCoeffs_[coefLabel]
364 + Y2*jt2.lowCpCoeffs_[coefLabel];
370 &&
notEqual(jt1.Tcommon_, jt2.Tcommon_)
374 <<
"Tcommon " << jt1.Tcommon_ <<
" for "
375 << (jt1.name().size() ? jt1.name() :
"others")
376 <<
" != " << jt2.Tcommon_ <<
" for "
377 << (jt2.name().size() ? jt2.name() :
"others")
381 return janafThermo<EquationOfState>
384 max(jt1.Tlow_, jt2.Tlow_),
385 min(jt1.Thigh_, jt2.Thigh_),
394 template<
class EquationOfState>
398 const janafThermo<EquationOfState>& jt
401 return janafThermo<EquationOfState>
403 s*
static_cast<const EquationOfState&
>(jt),
413 template<
class EquationOfState>
416 const janafThermo<EquationOfState>& jt1,
417 const janafThermo<EquationOfState>& jt2
422 static_cast<const EquationOfState&
>(jt1)
423 ==
static_cast<const EquationOfState&
>(jt2)
426 const scalar Y1 = jt2.Y()/eofs.Y();
427 const scalar Y2 = jt1.Y()/eofs.Y();
435 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
439 highCpCoeffs[coefLabel] =
440 Y1*jt2.highCpCoeffs_[coefLabel]
441 - Y2*jt1.highCpCoeffs_[coefLabel];
443 lowCpCoeffs[coefLabel] =
444 Y1*jt2.lowCpCoeffs_[coefLabel]
445 - Y2*jt1.lowCpCoeffs_[coefLabel];
451 &&
notEqual(jt2.Tcommon_, jt1.Tcommon_)
455 <<
"Tcommon " << jt2.Tcommon_ <<
" for "
456 << (jt2.name().size() ? jt2.name() :
"others")
457 <<
" != " << jt1.Tcommon_ <<
" for "
458 << (jt1.name().size() ? jt1.name() :
"others")
462 return janafThermo<EquationOfState>
465 max(jt2.Tlow_, jt1.Tlow_),
466 min(jt2.Thigh_, jt1.Thigh_),
int debug
Static debugging option.
const coeffArray & highCpCoeffs() const
Return const access to the high temperature poly coefficients.
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))
scalar Tlow() const
Return const access to the low temperature limit.
static const int nCoeffs_
JANAF tables based thermodynamics package templated into the equation of state.
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
Ostream & endl(Ostream &os)
Add newline and flush stream.
const coeffArray & lowCpCoeffs() const
Return const access to the low temperature poly coefficients.
const dimensionedScalar Tstd
Standard temperature.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
scalar Cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kg K)].
#define R(A, B, C, D, E, F, K, M)
FixedList< scalar, nCoeffs_ > coeffArray
word name(const complex &c)
Return string representation of complex.
scalar Ha(const scalar p, const scalar T) const
Absolute 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.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedScalar log(const dimensionedScalar &ds)
PtrList< volScalarField > & Y
errorManipArg< error, int > exit(error &err, const int errNo=1)
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
scalar Tcommon() const
Return const access to the common temperature.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const volScalarField & Cp
scalar Ha(const scalar p, const scalar T) const
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.
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
scalar Thigh() const
Return const access to the high temperature limit.
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
bool notEqual(const Scalar s1, const Scalar s2)
#define WarningInFunction
Report a warning using Foam::Warning.
scalar Hc() const
Chemical enthalpy [J/kg].
scalar S(const scalar p, const scalar T) const
Entropy [J/(kg K)].