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 clamp(
T, Tlow_, Thigh_);
133 template<
class EquationOfState>
140 template<
class EquationOfState>
147 template<
class EquationOfState>
154 template<
class EquationOfState>
158 return highCpCoeffs_;
162 template<
class EquationOfState>
170 template<
class EquationOfState>
177 const coeffArray& a = coeffs(
T);
179 ((((a[4]*
T + a[3])*
T + a[2])*
T + a[1])*
T + a[0])
184 template<
class EquationOfState>
191 const coeffArray& a = coeffs(
T);
194 ((((a[4]/5.0*
T + a[3]/4.0)*
T + a[2]/3.0)*
T + a[1]/2.0)*
T + a[0])*
T 200 template<
class EquationOfState>
207 return Ha(
p,
T) - Hc();
211 template<
class EquationOfState>
214 const coeffArray& a = lowCpCoeffs_;
225 template<
class EquationOfState>
232 const coeffArray& a = coeffs(
T);
235 (((a[4]/4.0*
T + a[3]/3.0)*
T + a[2]/2.0)*
T + a[1])*
T + a[0]*
log(
T)
237 ) + EquationOfState::S(
p,
T);
241 template<
class EquationOfState>
247 const coeffArray& a = coeffs(
T);
252 - (((a[4]/20.0*
T + a[3]/12.0)*
T + a[2]/6.0)*
T + a[1]/2.0)*
T 260 template<
class EquationOfState>
267 const coeffArray& a = coeffs(
T);
269 (((4*a[4]*
T + 3*a[3])*
T + 2*a[2])*
T + a[1]);
274 template<
class EquationOfState>
280 scalar Y1 = this->
Y();
282 EquationOfState::operator+=(jt);
284 if (
mag(this->
Y()) > SMALL)
287 const scalar Y2 = jt.Y()/this->
Y();
289 Tlow_ =
max(Tlow_, jt.Tlow_);
290 Thigh_ =
min(Thigh_, jt.Thigh_);
295 && !
equal(Tcommon_, jt.Tcommon_)
299 <<
"Tcommon " << Tcommon_ <<
" for " 300 << (this->
name().size() ? this->
name() :
"others")
301 <<
" != " << jt.Tcommon_ <<
" for " 302 << (jt.name().size() ? jt.name() :
"others")
309 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
313 highCpCoeffs_[coefLabel] =
314 Y1*highCpCoeffs_[coefLabel]
315 + Y2*jt.highCpCoeffs_[coefLabel];
317 lowCpCoeffs_[coefLabel] =
318 Y1*lowCpCoeffs_[coefLabel]
319 + Y2*jt.lowCpCoeffs_[coefLabel];
327 template<
class EquationOfState>
330 const janafThermo<EquationOfState>& jt1,
331 const janafThermo<EquationOfState>& jt2
334 EquationOfState eofs = jt1;
337 if (
mag(eofs.Y()) < SMALL)
339 return janafThermo<EquationOfState>
351 const scalar Y1 = jt1.Y()/eofs.Y();
352 const scalar Y2 = jt2.Y()/eofs.Y();
360 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
364 highCpCoeffs[coefLabel] =
365 Y1*jt1.highCpCoeffs_[coefLabel]
366 + Y2*jt2.highCpCoeffs_[coefLabel];
368 lowCpCoeffs[coefLabel] =
369 Y1*jt1.lowCpCoeffs_[coefLabel]
370 + Y2*jt2.lowCpCoeffs_[coefLabel];
376 && !
equal(jt1.Tcommon_, jt2.Tcommon_)
380 <<
"Tcommon " << jt1.Tcommon_ <<
" for " 381 << (jt1.name().size() ? jt1.name() :
"others")
382 <<
" != " << jt2.Tcommon_ <<
" for " 383 << (jt2.name().size() ? jt2.name() :
"others")
387 return janafThermo<EquationOfState>
390 max(jt1.Tlow_, jt2.Tlow_),
391 min(jt1.Thigh_, jt2.Thigh_),
400 template<
class EquationOfState>
404 const janafThermo<EquationOfState>& jt
407 return janafThermo<EquationOfState>
409 s*
static_cast<const EquationOfState&
>(jt),
419 template<
class EquationOfState>
422 const janafThermo<EquationOfState>& jt1,
423 const janafThermo<EquationOfState>& jt2
428 static_cast<const EquationOfState&>(jt1)
429 == static_cast<const EquationOfState&>(jt2)
432 const scalar Y1 = jt2.Y()/eofs.Y();
433 const scalar Y2 = jt1.Y()/eofs.Y();
441 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
445 highCpCoeffs[coefLabel] =
446 Y1*jt2.highCpCoeffs_[coefLabel]
447 - Y2*jt1.highCpCoeffs_[coefLabel];
449 lowCpCoeffs[coefLabel] =
450 Y1*jt2.lowCpCoeffs_[coefLabel]
451 - Y2*jt1.lowCpCoeffs_[coefLabel];
457 && !
equal(jt2.Tcommon_, jt1.Tcommon_)
461 <<
"Tcommon " << jt2.Tcommon_ <<
" for " 462 << (jt2.name().size() ? jt2.name() :
"others")
463 <<
" != " << jt1.Tcommon_ <<
" for " 464 << (jt1.name().size() ? jt1.name() :
"others")
468 return janafThermo<EquationOfState>
471 max(jt2.Tlow_, jt1.Tlow_),
472 min(jt2.Thigh_, jt1.Thigh_),
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
bool equal(const T &a, const T &b)
Compare two values for equality.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
scalar Thigh() const
Return const access to the high temperature limit.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
scalar Hc() const
Chemical enthalpy [J/kg].
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].
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
A class for handling words, derived from Foam::string.
scalar S(const scalar p, const scalar T) const
Entropy [J/(kg K)].
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
scalar Tlow() const
Return const access to the low temperature limit.
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
const dimensionedScalar Tstd
Standard temperature.
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.
const volScalarField & Cp
int debug
Static debugging option.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
#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 Ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kg].
static constexpr int nCoeffs_
scalar Ha(const scalar p, const scalar T) const
FixedList< scalar, nCoeffs_ > coeffArray
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 Cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kg K)].
volScalarField H(IOobject("H", runTime.timeName(), mesh.thisDb(), IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
dimensionSet clamp(const dimensionSet &a, const dimensionSet &range)
scalar Tcommon() const
Return const access to the common temperature.