28 template<
class EquationOfState>
31 const EquationOfState& st,
46 template<
class EquationOfState>
53 EquationOfState(
name, ct),
61 template<
class EquationOfState>
74 template<
class EquationOfState>
84 template<
class EquationOfState>
95 template<
class EquationOfState>
102 return Cv_*(
T - Tref_) + Esref_ + EquationOfState::E(
p,
T);
106 template<
class EquationOfState>
113 return Es(
p,
T) + Hf();
117 template<
class EquationOfState>
124 template<
class EquationOfState>
135 template<
class EquationOfState>
147 template<
class EquationOfState>
161 template<
class EquationOfState>
167 scalar Y1 = this->
Y();
169 EquationOfState::operator+=(ct);
171 if (
mag(this->
Y()) > small)
180 <<
"Tref " << Tref_ <<
" for "
181 << (this->
name().size() ? this->
name() :
"others")
182 <<
" != " << ct.Tref_ <<
" for "
183 << (ct.name().size() ? ct.name() :
"others")
188 const scalar Y2 = ct.Y()/this->
Y();
190 Cv_ = Y1*Cv_ + Y2*ct.Cv_;
191 Hf_ = Y1*Hf_ + Y2*ct.Hf_;
192 Esref_ = Y1*Esref_ + Y2*ct.Esref_;
199 template<
class EquationOfState>
208 static_cast<const EquationOfState&
>(ct1)
209 +
static_cast<const EquationOfState&
>(ct2)
212 if (
mag(eofs.Y()) < small)
227 eConstThermo<EquationOfState>::debug
232 <<
"Tref " << ct1.Tref_ <<
" for "
233 << (ct1.name().size() ? ct1.name() :
"others")
234 <<
" != " << ct2.Tref_ <<
" for "
235 << (ct2.name().size() ? ct2.name() :
"others")
239 return eConstThermo<EquationOfState>
242 ct1.Y()/eofs.Y()*ct1.Cv_
243 + ct2.Y()/eofs.Y()*ct2.Cv_,
244 ct1.Y()/eofs.Y()*ct1.Hf_
245 + ct2.Y()/eofs.Y()*ct2.Hf_,
247 ct1.Y()/eofs.Y()*ct1.Esref_
248 + ct2.Y()/eofs.Y()*ct2.Esref_
254 template<
class EquationOfState>
258 const eConstThermo<EquationOfState>& ct
261 return eConstThermo<EquationOfState>
263 s*
static_cast<const EquationOfState&
>(ct),
272 template<
class EquationOfState>
275 const eConstThermo<EquationOfState>& ct1,
276 const eConstThermo<EquationOfState>& ct2
281 static_cast<const EquationOfState&
>(ct1)
282 ==
static_cast<const EquationOfState&
>(ct2)
287 eConstThermo<EquationOfState>::debug
292 <<
"Tref " << ct1.Tref_ <<
" for "
293 << (ct1.name().size() ? ct1.name() :
"others")
294 <<
" != " << ct2.Tref_ <<
" for "
295 << (ct2.name().size() ? ct2.name() :
"others")
299 return eConstThermo<EquationOfState>
302 ct2.Y()/eofs.Y()*ct2.Cv_
303 - ct1.Y()/eofs.Y()*ct1.Cv_,
304 ct2.Y()/eofs.Y()*ct2.Hf_
305 - ct1.Y()/eofs.Y()*ct1.Hf_,
307 ct2.Y()/eofs.Y()*ct2.Esref_
308 - ct1.Y()/eofs.Y()*ct1.Esref_
scalar Cp(const scalar p, const scalar T) const
scalar Cv(const scalar p, const scalar T) const
scalar Es(const scalar p, const scalar T) const
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Internal energy based thermodynamics package using a constant heat capacity at constant volume:
scalar Hf() const
Enthalpy of formation [J/kg].
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
scalar S(const scalar p, const scalar T) const
Entropy [J/kg/K].
scalar Cv(const scalar p, const scalar T) const
Heat capacity at constant volume [J/kg/K].
scalar Gstd(const scalar T) const
Gibbs free energy of the mixture in the standard state [J/kg].
autoPtr< eConstThermo > clone() const
Construct and return a clone.
scalar Es(const scalar p, const scalar T) const
Sensible internal energy [J/kg].
eConstThermo(const EquationOfState &st, const scalar Cv, const scalar Hf, const scalar Tref, const scalar Esref)
Construct from components.
scalar Ea(const scalar p, const scalar T) const
Absolute internal energy [J/kg].
A class for handling words, derived from string.
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
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.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const dimensionedScalar Pstd
Standard pressure.
const dimensionedScalar Tstd
Standard temperature.
tmp< fvMatrix< Type > > S(const Pair< tmp< volScalarField::Internal >> &, const VolField< Type > &)
errorManipArg< error, int > exit(error &err, const int errNo=1)
bool notEqual(const Scalar s1, const Scalar s2)
dimensionedScalar log(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
PtrList< volScalarField > & Y