Go to the documentation of this file.
32 template<
class Specie>
52 template<
class Specie>
67 template<
class Specie>
75 template<
class Specie>
88 template<
class Specie>
95 return rho0_*
pow((
p + B_)/(p0_ + B_), 1.0/gamma_);
99 template<
class Specie>
110 template<
class Specie>
121 template<
class Specie>
132 template<
class Specie>
143 template<
class Specie>
150 scalar
n = 1 - 1.0/gamma_;
157 template<
class Specie>
165 (rho0_/(gamma_*(p0_ + B_)))
166 *
pow((
p + B_)/(p0_ + B_), 1.0/gamma_ - 1.0);
170 template<
class Specie>
177 template<
class Specie>
190 template<
class Specie>
196 scalar Y1 = this->
Y();
197 Specie::operator+=(pf);
199 if (
mag(this->
Y()) > SMALL)
202 const scalar Y2 = pf.Y()/this->
Y();
204 p0_ = Y1*p0_ + Y2*pf.p0_;
205 rho0_ = Y1*rho0_ + Y2*pf.rho0_;
206 gamma_ = Y1*gamma_ + Y2*pf.gamma_;
207 B_ = Y1*B_ + Y2*pf.B_;
212 template<
class Specie>
215 Specie::operator*=(
s);
221 template<
class Specie>
230 static_cast<const Specie&
>(pf1)
231 +
static_cast<const Specie&
>(pf2)
234 if (
mag(sp.Y()) < SMALL)
247 const scalar Y1 = pf1.Y()/sp.Y();
248 const scalar Y2 = pf2.Y()/sp.Y();
250 return adiabaticPerfectFluid<Specie>
253 Y1*pf1.p0_ + Y2*pf2.p0_,
254 Y1*pf1.rho0_ + Y2*pf2.rho0_,
255 Y1*pf1.gamma_ + Y2*pf2.gamma_,
256 Y1*pf1.B_ + Y2*pf2.B_
262 template<
class Specie>
266 const adiabaticPerfectFluid<Specie>& pf
269 return adiabaticPerfectFluid<Specie>
271 s*
static_cast<const Specie&
>(pf),
280 template<
class Specie>
283 const adiabaticPerfectFluid<Specie>& pf1,
284 const adiabaticPerfectFluid<Specie>& pf2
289 static_cast<const Specie&
>(pf1)
290 ==
static_cast<const Specie&
>(pf2)
293 const scalar Y1 = pf1.Y()/sp.Y();
294 const scalar Y2 = pf2.Y()/sp.Y();
296 return adiabaticPerfectFluid<Specie>
299 Y2*pf2.p0_ - Y1*pf1.p0_,
300 Y2*pf2.rho0_ - Y1*pf1.rho0_,
301 Y2*pf2.gamma_ - Y1*pf1.gamma_,
302 Y2*pf2.B_ - Y1*pf1.B_
scalar Cv(scalar p, scalar T) const
Return Cv departure [J/(kg K].
scalar psi(scalar p, scalar T) const
Return compressibility rho/p [s^2/m^2].
A class for handling words, derived from Foam::string.
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))
const dimensionedScalar Pstd
Standard pressure.
scalar S(const scalar p, const scalar T) const
Return entropy [J/(kg K)].
scalar Z(scalar p, scalar T) const
Return compression factor [].
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
Adiabatic perfect fluid equation of state.
scalar E(const scalar p, const scalar T) const
Return internal energy departure [J/kg].
scalar H(const scalar p, const scalar T) const
Return enthalpy departure [J/kg].
word name(const complex &c)
Return string representation of complex.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
scalar Cp(scalar p, scalar T) const
Return Cp departure [J/(kg K].
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.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
static autoPtr< adiabaticPerfectFluid > New(const dictionary &dict)
autoPtr< adiabaticPerfectFluid > clone() const
Construct and return a clone.
const volScalarField & p0
adiabaticPerfectFluid(const Specie &sp, const scalar p0, const scalar rho0, const scalar gamma, const scalar B)
Construct from components.
void operator*=(const scalar)