32 template<
class Specie>
48 template<
class Specie>
61 template<
class Specie>
69 template<
class Specie>
82 template<
class Specie>
85 return rho0_ + psi_*
p;
89 template<
class Specie>
96 template<
class Specie>
102 template<
class Specie>
109 template<
class Specie>
116 template<
class Specie>
119 return -
log((rho0_ + psi_*
p)/(rho0_ + psi_*
Pstd))/(
T*psi_);
123 template<
class Specie>
130 template<
class Specie>
137 template<
class Specie>
146 template<
class Specie>
152 scalar Y1 = this->
Y();
153 Specie::operator+=(pf);
155 if (
mag(this->
Y()) > SMALL)
158 const scalar Y2 = pf.Y()/this->
Y();
160 psi_ = Y1*psi_ + Y2*pf.psi_;
161 rho0_ = Y1*rho0_ + Y2*pf.rho0_;
166 template<
class Specie>
169 Specie::operator*=(
s);
175 template<
class Specie>
178 const linear<Specie>& pf1,
179 const linear<Specie>& pf2
184 static_cast<const Specie&>(pf1)
185 + static_cast<const Specie&>(pf2)
188 if (
mag(sp.Y()) < SMALL)
190 return linear<Specie>
199 const scalar Y1 = pf1.Y()/sp.Y();
200 const scalar Y2 = pf2.Y()/sp.Y();
202 return linear<Specie>
205 Y1*pf1.psi_ + Y2*pf2.psi_,
206 Y1*pf1.rho0_ + Y2*pf2.rho0_
212 template<
class Specie>
216 const linear<Specie>& pf
219 return linear<Specie>
221 s*
static_cast<const Specie&
>(pf),
228 template<
class Specie>
231 const linear<Specie>& pf1,
232 const linear<Specie>& pf2
237 static_cast<const Specie&>(pf1)
238 == static_cast<const Specie&>(pf2)
241 const scalar Y1 = pf1.Y()/sp.Y();
242 const scalar Y2 = pf2.Y()/sp.Y();
244 return linear<Specie>
247 Y2*pf2.psi_ - Y1*pf1.psi_,
248 Y2*pf2.rho0_ - Y1*pf1.rho0_
scalar H(const scalar p, const scalar T) const
Return enthalpy departure [J/kg].
dimensionedScalar log(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
A list of keyword definitions, which are a keyword followed by a number of values (eg...
static autoPtr< linear > New(const dictionary &dict)
Central-differencing interpolation scheme class.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
scalar psi(scalar p, scalar T) const
Return compressibility rho/p [s^2/m^2].
scalar S(const scalar p, const scalar T) const
Return entropy [J/(kg K)].
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
A class for handling words, derived from Foam::string.
scalar Cv(scalar p, scalar T) const
Return Cv departure [J/(kg K].
void operator*=(const scalar)
const dimensionedScalar Pstd
Standard pressure.
linear(const fvMesh &mesh)
Construct from mesh.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
PtrList< volScalarField > & Y
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
const volScalarField & psi
autoPtr< linear > clone() const
Construct and return a clone.
scalar Cp(scalar p, scalar T) const
Return Cp departure [J/(kg K].
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))
scalar E(const scalar p, const scalar T) const
Return internal energy departure [J/kg].
scalar Z(scalar p, scalar T) const
Return compression factor [].