31 template<
class EquationOfState>
34 const EquationOfState& st,
38 const typename janafThermo<EquationOfState>::coeffArray& highCpCoeffs,
39 const typename janafThermo<EquationOfState>::coeffArray& lowCpCoeffs
49 highCpCoeffs_[coefLabel] = highCpCoeffs[coefLabel];
50 lowCpCoeffs_[coefLabel] = lowCpCoeffs[coefLabel];
55 template<
class EquationOfState>
75 template<
class EquationOfState>
82 EquationOfState(name, jt),
89 highCpCoeffs_[coefLabel] = jt.highCpCoeffs_[coefLabel];
90 lowCpCoeffs_[coefLabel] = jt.lowCpCoeffs_[coefLabel];
97 template<
class EquationOfState>
103 if (T < Tlow_ || T > Thigh_)
107 "janafThermo<EquationOfState>::limit(const scalar T) const" 108 ) <<
"attempt to use janafThermo<EquationOfState>" 109 " out of temperature range " 110 << Tlow_ <<
" -> " << Thigh_ <<
"; T = " << T
113 return min(
max(T, Tlow_), Thigh_);
122 template<
class EquationOfState>
129 template<
class EquationOfState>
136 template<
class EquationOfState>
143 template<
class EquationOfState>
147 return highCpCoeffs_;
151 template<
class EquationOfState>
159 template<
class EquationOfState>
167 return RR*((((a[4]*T + a[3])*T + a[2])*T + a[1])*T + a[0]);
171 template<
class EquationOfState>
181 ((((a[4]/5.0*T + a[3]/4.0)*T + a[2]/3.0)*T + a[1]/2.0)*T + a[0])*T
187 template<
class EquationOfState>
194 return ha(p, T) -
hc();
198 template<
class EquationOfState>
212 template<
class EquationOfState>
223 (((a[4]/4.0*T + a[3]/3.0)*T + a[2]/2.0)*T + a[1])*T + a[0]*
log(T)
232 template<
class EquationOfState>
233 inline void Foam::janafThermo<EquationOfState>::operator+=
238 scalar molr1 = this->nMoles();
240 EquationOfState::operator+=(jt);
242 molr1 /= this->nMoles();
243 scalar molr2 = jt.nMoles()/this->nMoles();
245 Tlow_ =
max(Tlow_, jt.Tlow_);
246 Thigh_ =
min(Thigh_, jt.Thigh_);
252 "janafThermo<EquationOfState>::operator+=" 253 "(const janafThermo<EquationOfState>& jt) const" 254 ) <<
"Tcommon " << Tcommon_ <<
" for " 255 << (this->
name().size() ? this->
name() :
"others")
256 <<
" != " << jt.Tcommon_ <<
" for " 257 << (jt.name().size() ? jt.name() :
"others")
264 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
268 highCpCoeffs_[coefLabel] =
269 molr1*highCpCoeffs_[coefLabel]
270 + molr2*jt.highCpCoeffs_[coefLabel];
272 lowCpCoeffs_[coefLabel] =
273 molr1*lowCpCoeffs_[coefLabel]
274 + molr2*jt.lowCpCoeffs_[coefLabel];
279 template<
class EquationOfState>
280 inline void Foam::janafThermo<EquationOfState>::operator-=
285 scalar molr1 = this->nMoles();
287 EquationOfState::operator-=(jt);
289 molr1 /= this->nMoles();
290 scalar molr2 = jt.nMoles()/this->nMoles();
292 Tlow_ =
max(Tlow_, jt.Tlow_);
293 Thigh_ =
min(Thigh_, jt.Thigh_);
299 "janafThermo<EquationOfState>::operator-=" 300 "(const janafThermo<EquationOfState>& jt) const" 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 molr1*highCpCoeffs_[coefLabel]
317 - molr2*jt.highCpCoeffs_[coefLabel];
319 lowCpCoeffs_[coefLabel] =
320 molr1*lowCpCoeffs_[coefLabel]
321 - molr2*jt.lowCpCoeffs_[coefLabel];
328 template<
class EquationOfState>
335 EquationOfState eofs = jt1;
338 scalar molr1 = jt1.nMoles()/eofs.nMoles();
339 scalar molr2 = jt2.nMoles()/eofs.nMoles();
347 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
351 highCpCoeffs[coefLabel] =
352 molr1*jt1.highCpCoeffs_[coefLabel]
353 + molr2*jt2.highCpCoeffs_[coefLabel];
355 lowCpCoeffs[coefLabel] =
356 molr1*jt1.lowCpCoeffs_[coefLabel]
357 + molr2*jt2.lowCpCoeffs_[coefLabel];
363 &&
notEqual(jt1.Tcommon_, jt2.Tcommon_)
369 "(const janafThermo<EquationOfState>& jt1," 370 " const janafThermo<EquationOfState>& jt2)" 371 ) <<
"Tcommon " << jt1.Tcommon_ <<
" for " 372 << (jt1.name().size() ? jt1.name() :
"others")
373 <<
" != " << jt2.Tcommon_ <<
" for " 374 << (jt2.name().size() ? jt2.name() :
"others")
381 max(jt1.Tlow_, jt2.Tlow_),
382 min(jt1.Thigh_, jt2.Thigh_),
390 template<
class EquationOfState>
397 EquationOfState eofs = jt1;
400 scalar molr1 = jt1.nMoles()/eofs.nMoles();
401 scalar molr2 = jt2.nMoles()/eofs.nMoles();
409 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
413 highCpCoeffs[coefLabel] =
414 molr1*jt1.highCpCoeffs_[coefLabel]
415 - molr2*jt2.highCpCoeffs_[coefLabel];
417 lowCpCoeffs[coefLabel] =
418 molr1*jt1.lowCpCoeffs_[coefLabel]
419 - molr2*jt2.lowCpCoeffs_[coefLabel];
425 &&
notEqual(jt1.Tcommon_, jt2.Tcommon_)
431 "(const janafThermo<EquationOfState>& jt1," 432 " const janafThermo<EquationOfState>& jt2)" 433 ) <<
"Tcommon " << jt1.Tcommon_ <<
" for " 434 << (jt1.name().size() ? jt1.name() :
"others")
435 <<
" != " << jt2.Tcommon_ <<
" for " 436 << (jt2.name().size() ? jt2.name() :
"others")
443 max(jt1.Tlow_, jt2.Tlow_),
444 min(jt1.Thigh_, jt2.Thigh_),
452 template<
class EquationOfState>
461 s*
static_cast<const EquationOfState&
>(jt),
471 template<
class EquationOfState>
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 ))
word name(const complex &)
Return a string representation of a complex.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
A class for handling words, derived from string.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
errorManipArg< error, int > exit(error &err, const int errNo=1)
scalar Thigh() const
Return const access to the high temperature limit.
scalar cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kmol K)].
scalar Tcommon() const
Return const access to the common temperature.
dimensionedScalar log(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
#define WarningIn(functionName)
Report a warning using Foam::Warning.
scalar hc() const
Chemical enthalpy [J/kmol].
const dimensionedScalar Tstd
Standard temperature.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
scalar s(const scalar p, const scalar T) const
Entropy [J/(kmol K)].
janafThermo(const EquationOfState &st, const scalar Tlow, const scalar Thigh, const scalar Tcommon, const coeffArray &highCpCoeffs, const coeffArray &lowCpCoeffs)
Construct from components.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
const coeffArray & lowCpCoeffs() const
Return const access to the low temperature poly coefficients.
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
const scalar RR
Universal gas constant (default in [J/(kmol K)])
bool notEqual(const Scalar s1, const Scalar s2)
scalar Tlow() const
Return const access to the low temperature limit.
const coeffArray & highCpCoeffs() const
Return const access to the high temperature poly coefficients.
static const int nCoeffs_
scalar hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kmol].
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
JANAF tables based thermodynamics package templated into the equation of state.
scalar ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kmol].