Go to the documentation of this file.
34 template<
class EquationOfState>
37 const EquationOfState& st,
41 const typename janafThermo<EquationOfState>::coeffArray& highCpCoeffs,
42 const typename janafThermo<EquationOfState>::coeffArray& lowCpCoeffs,
43 const bool convertCoeffs
53 for (label coefLabel=0; coefLabel<
nCoeffs_; coefLabel++)
56 lowCpCoeffs_[coefLabel] =
lowCpCoeffs[coefLabel]*this->
R();
61 for (label coefLabel=0; coefLabel<
nCoeffs_; coefLabel++)
70 template<
class EquationOfState>
90 template<
class EquationOfState>
97 EquationOfState(
name, jt),
100 Tcommon_(jt.Tcommon_)
102 for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
104 highCpCoeffs_[coefLabel] = jt.highCpCoeffs_[coefLabel];
105 lowCpCoeffs_[coefLabel] = jt.lowCpCoeffs_[coefLabel];
112 template<
class EquationOfState>
118 if (T < Tlow_ || T > Thigh_)
121 <<
"attempt to use janafThermo<EquationOfState>"
122 " out of temperature range "
123 << Tlow_ <<
" -> " << Thigh_ <<
"; T = " <<
T
126 return min(
max(
T, Tlow_), Thigh_);
135 template<
class EquationOfState>
142 template<
class EquationOfState>
149 template<
class EquationOfState>
156 template<
class EquationOfState>
160 return highCpCoeffs_;
164 template<
class EquationOfState>
172 template<
class EquationOfState>
181 ((((a[4]*
T + a[3])*
T + a[2])*
T + a[1])*
T + a[0])
186 template<
class EquationOfState>
196 ((((a[4]/5.0*
T + a[3]/4.0)*
T + a[2]/3.0)*
T + a[1]/2.0)*
T + a[0])*
T
202 template<
class EquationOfState>
209 return Ha(
p,
T) - Hc();
213 template<
class EquationOfState>
227 template<
class EquationOfState>
237 (((a[4]/4.0*
T + a[3]/3.0)*
T + a[2]/2.0)*
T + a[1])*
T + a[0]*
log(
T)
239 ) + EquationOfState::S(
p,
T);
243 template<
class EquationOfState>
254 - (((a[4]/20.0*
T + a[3]/12.0)*
T + a[2]/6.0)*
T + a[1]/2.0)*
T
262 template<
class EquationOfState>
271 (((4*a[4]*
T + 3*a[3])*
T + 2*a[2])*
T + a[1]);
276 template<
class EquationOfState>
282 scalar Y1 = this->
Y();
284 EquationOfState::operator+=(jt);
286 if (
mag(this->
Y()) > SMALL)
289 const scalar Y2 = jt.Y()/this->
Y();
291 Tlow_ =
max(Tlow_, jt.Tlow_);
292 Thigh_ =
min(Thigh_, jt.Thigh_);
301 <<
"Tcommon " << Tcommon_ <<
" for "
302 << (this->
name().size() ? this->
name() :
"others")
303 <<
" != " << jt.Tcommon_ <<
" for "
304 << (jt.name().size() ? jt.name() :
"others")
311 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
315 highCpCoeffs_[coefLabel] =
316 Y1*highCpCoeffs_[coefLabel]
317 + Y2*jt.highCpCoeffs_[coefLabel];
319 lowCpCoeffs_[coefLabel] =
320 Y1*lowCpCoeffs_[coefLabel]
321 + Y2*jt.lowCpCoeffs_[coefLabel];
329 template<
class EquationOfState>
336 EquationOfState eofs = jt1;
339 if (
mag(eofs.Y()) < SMALL)
353 const scalar Y1 = jt1.Y()/eofs.Y();
354 const scalar Y2 = jt2.Y()/eofs.Y();
362 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
366 highCpCoeffs[coefLabel] =
367 Y1*jt1.highCpCoeffs_[coefLabel]
368 + Y2*jt2.highCpCoeffs_[coefLabel];
370 lowCpCoeffs[coefLabel] =
371 Y1*jt1.lowCpCoeffs_[coefLabel]
372 + Y2*jt2.lowCpCoeffs_[coefLabel];
378 &&
notEqual(jt1.Tcommon_, jt2.Tcommon_)
382 <<
"Tcommon " << jt1.Tcommon_ <<
" for "
383 << (jt1.name().size() ? jt1.name() :
"others")
384 <<
" != " << jt2.Tcommon_ <<
" for "
385 << (jt2.name().size() ? jt2.name() :
"others")
389 return janafThermo<EquationOfState>
392 max(jt1.Tlow_, jt2.Tlow_),
393 min(jt1.Thigh_, jt2.Thigh_),
402 template<
class EquationOfState>
406 const janafThermo<EquationOfState>& jt
409 return janafThermo<EquationOfState>
411 s*
static_cast<const EquationOfState&
>(jt),
421 template<
class EquationOfState>
424 const janafThermo<EquationOfState>& jt1,
425 const janafThermo<EquationOfState>& jt2
430 static_cast<const EquationOfState&
>(jt1)
431 ==
static_cast<const EquationOfState&
>(jt2)
434 const scalar Y1 = jt2.Y()/eofs.Y();
435 const scalar Y2 = jt1.Y()/eofs.Y();
443 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
447 highCpCoeffs[coefLabel] =
448 Y1*jt2.highCpCoeffs_[coefLabel]
449 - Y2*jt1.highCpCoeffs_[coefLabel];
451 lowCpCoeffs[coefLabel] =
452 Y1*jt2.lowCpCoeffs_[coefLabel]
453 - Y2*jt1.lowCpCoeffs_[coefLabel];
459 &&
notEqual(jt2.Tcommon_, jt1.Tcommon_)
463 <<
"Tcommon " << jt2.Tcommon_ <<
" for "
464 << (jt2.name().size() ? jt2.name() :
"others")
465 <<
" != " << jt1.Tcommon_ <<
" for "
466 << (jt1.name().size() ? jt1.name() :
"others")
470 return janafThermo<EquationOfState>
473 max(jt2.Tlow_, jt1.Tlow_),
474 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.
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)].
static constexpr int nCoeffs_
#define R(A, B, C, D, E, F, K, M)
FixedList< scalar, nCoeffs_ > coeffArray
scalar Ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kg].
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
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
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
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 Gstd(const scalar T) const
Gibbs free energy of the mixture in the standard state [J/kg].
scalar Hc() const
Chemical enthalpy [J/kg].
scalar S(const scalar p, const scalar T) const
Entropy [J/(kg K)].