65 #ifndef adiabaticPerfectFluid_H
66 #define adiabaticPerfectFluid_H
77 template<
class Specie>
class adiabaticPerfectFluid;
79 template<
class Specie>
80 inline adiabaticPerfectFluid<Specie>
operator+
82 const adiabaticPerfectFluid<Specie>&,
83 const adiabaticPerfectFluid<Specie>&
86 template<
class Specie>
87 inline adiabaticPerfectFluid<Specie>
operator*
93 template<
class Specie>
100 template<
class Specie>
112 template<
class Specie>
165 return "adiabaticPerfectFluid<" +
word(Specie::typeName_()) +
'>';
178 inline scalar
rho(scalar
p, scalar
T)
const;
181 inline scalar
H(
const scalar
p,
const scalar
T)
const;
184 inline scalar
Cp(scalar
p, scalar
T)
const;
187 inline scalar
E(
const scalar
p,
const scalar
T)
const;
190 inline scalar
Cv(scalar
p, scalar
T)
const;
193 inline scalar
Sp(
const scalar
p,
const scalar
T)
const;
196 inline scalar
Sv(
const scalar
p,
const scalar
T)
const;
199 inline scalar
psi(scalar
p, scalar
T)
const;
202 inline scalar
Z(scalar
p, scalar
T)
const;
205 inline scalar
CpMCv(scalar
p, scalar
T)
const;
208 inline scalar
alphav(
const scalar
p,
const scalar
T)
const;
246 friend Ostream& operator<< <Specie>
static const Foam::dimensionedScalar B("B", Foam::dimless, 18.678)
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Adiabatic perfect fluid equation of state for liquids:
scalar Cv(scalar p, scalar T) const
Return Cv contribution [J/(kg K].
adiabaticPerfectFluid(const Specie &sp, const scalar p0, const scalar rho0, const scalar gamma, const scalar B)
Construct from components.
scalar Sv(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cv/T [J/kg/K].
scalar E(const scalar p, const scalar T) const
Return internal energy contribution [J/kg].
scalar psi(scalar p, scalar T) const
Return compressibility [s^2/m^2].
scalar H(const scalar p, const scalar T) const
Return enthalpy contribution [J/kg].
scalar alphav(const scalar p, const scalar T) const
Return volumetric coefficient of thermal expansion [1/T].
static word typeName()
Return the instantiated type name.
void write(Ostream &os) const
Write to Ostream.
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
scalar Cp(scalar p, scalar T) const
Return Cp contribution [J/(kg K].
scalar Sp(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cp/T [J/kg/K].
autoPtr< adiabaticPerfectFluid > clone() const
Construct and return a clone.
static const bool isochoric
Is the equation of state is isochoric i.e. rho = const.
void operator+=(const adiabaticPerfectFluid &)
static const bool incompressible
Is the equation of state is incompressible i.e. rho != f(p)
scalar Z(scalar p, scalar T) const
Return compression factor [].
void operator*=(const scalar)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
A list of keyword definitions, which are a keyword followed by any number of values (e....
A class for handling words, derived from string.
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))
word name(const complex &)
Return a string representation of a complex.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)