Go to the documentation of this file.
30 template<
class EquationOfState>
33 const EquationOfState& st,
46 template<
class EquationOfState>
50 const eConstThermo<EquationOfState>& ct
53 EquationOfState(
name, ct),
59 template<
class EquationOfState>
67 template<
class EquationOfState>
77 template<
class EquationOfState>
87 template<
class EquationOfState>
98 template<
class EquationOfState>
105 return Cv_*
T + EquationOfState::E(
p,
T);
109 template<
class EquationOfState>
116 return Es(
p,
T) + Hc();
120 template<
class EquationOfState>
127 template<
class EquationOfState>
138 template<
class EquationOfState>
149 template<
class EquationOfState>
161 template<
class EquationOfState>
167 scalar Y1 = this->
Y();
169 EquationOfState::operator+=(ct);
171 if (
mag(this->
Y()) > SMALL)
174 const scalar Y2 = ct.Y()/this->
Y();
176 Cv_ = Y1*Cv_ + Y2*ct.Cv_;
177 Hf_ = Y1*Hf_ + Y2*ct.Hf_;
184 template<
class EquationOfState>
193 static_cast<const EquationOfState&
>(ct1)
194 +
static_cast<const EquationOfState&
>(ct2)
197 if (
mag(eofs.Y()) < SMALL)
208 return eConstThermo<EquationOfState>
211 ct1.Y()/eofs.Y()*ct1.Cv_
212 + ct2.Y()/eofs.Y()*ct2.Cv_,
213 ct1.Y()/eofs.Y()*ct1.Hf_
214 + ct2.Y()/eofs.Y()*ct2.Hf_
220 template<
class EquationOfState>
224 const eConstThermo<EquationOfState>& ct
227 return eConstThermo<EquationOfState>
229 s*
static_cast<const EquationOfState&
>(ct),
236 template<
class EquationOfState>
239 const eConstThermo<EquationOfState>& ct1,
240 const eConstThermo<EquationOfState>& ct2
245 static_cast<const EquationOfState&
>(ct1)
246 ==
static_cast<const EquationOfState&
>(ct2)
249 return eConstThermo<EquationOfState>
252 ct2.Y()/eofs.Y()*ct2.Cv_
253 - ct1.Y()/eofs.Y()*ct1.Cv_,
254 ct2.Y()/eofs.Y()*ct2.Hf_
255 - ct1.Y()/eofs.Y()*ct1.Hf_
scalar Cv(const scalar p, const scalar T) const
Heat capacity at constant volume [J/(kg K)].
autoPtr< eConstThermo > clone() const
Construct and return a clone.
static autoPtr< eConstThermo > New(const dictionary &dict)
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))
Constant properties thermodynamics package templated on an equation of state.
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
scalar dGdT(const scalar p, const scalar T) const
Derivative of Gibbs free energy w.r.t. temperature.
const volScalarField & Cv
const dimensionedScalar Tstd
Standard temperature.
scalar Hc() const
Chemical enthalpy [J/kg].
scalar Ea(const scalar p, const scalar T) const
Absolute internal energy [J/kg].
word name(const complex &c)
Return string representation of complex.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar Es(const scalar p, const scalar T) const
Sensible internal energy [J/kg].
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
dimensionedScalar log(const dimensionedScalar &ds)
PtrList< volScalarField > & Y
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
scalar S(const scalar p, const scalar T) const
Entropy [J/(kg K)].
scalar Es(const scalar p, const scalar T) const
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const volScalarField & Cp
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.