33 template<
class Specie>
55 template<
class Specie>
71 template<
class Specie>
79 template<
class Specie>
92 template<
class Specie>
99 const scalar Z = this->Z(
p,
T);
100 return p/(Z*this->
R()*
T);
104 template<
class Specie>
107 const scalar
Pr =
p/Pc_;
108 const scalar Tr =
T/Tc_;
109 const scalar
B = 0.07780*
Pr/Tr;
110 const scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
113 const scalar Z = this->Z(
p,
T);
121 *
log((Z + 2.414*
B)/(Z - 0.414*
B))
126 template<
class Specie>
129 const scalar Tr =
T/Tc_;
130 const scalar a = 0.45724*
sqr(
RR*Tc_)/Pc_;
131 const scalar
b = 0.07780*
RR*Tc_/Pc_;
132 const scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
136 const scalar
B =
b*
p/(
RR*
T);
138 const scalar Z = this->Z(
p,
T);
143 const scalar
M = (
sqr(Z) + 2*
B*Z -
sqr(
B))/(Z -
B);
144 const scalar
N = ap*
B/(
b*
RR);
146 const scalar root2 =
sqrt(2.0);
150 app*(
T/(2*root2*
b))*
log((Z + (root2 + 1)*
B)/(Z - (root2 - 1)*
B))
157 template<
class Specie>
160 const scalar
Pr =
p/Pc_;
161 const scalar Tr =
T/Tc_;
162 const scalar
B = 0.07780*
Pr/Tr;
163 const scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
166 const scalar Z = this->Z(
p,
T);
173 *
log((Z + 2.414*
B)/(Z - 0.414*
B))
178 template<
class Specie>
181 const scalar a = 0.45724*
sqr(
RR*Tc_)/Pc_;
182 const scalar
b = 0.07780*
RR*Tc_/Pc_;
183 const scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
185 const scalar
B =
b*
p/(
RR*
T);
187 const scalar Z = this->Z(
p,
T);
191 const scalar root2 =
sqrt(2.0);
195 app*(
T/(2*root2*
b))*
log((Z + (root2 + 1)*
B)/(Z - (root2 - 1)*
B))
201 template<
class Specie>
208 const scalar
Pr =
p/Pc_;
209 const scalar Tr =
T/Tc_;
210 const scalar
B = 0.07780*
Pr/Tr;
211 const scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
213 const scalar Z = this->Z(
p,
T);
222 *
log((Z + 2.414*
B)/(Z - 0.414*
B))
228 template<
class Specie>
235 const scalar Z = this->Z(
p,
T);
237 return 1.0/(Z*this->
R()*
T);
241 template<
class Specie>
248 const scalar Tr =
T/Tc_;
249 const scalar a = 0.45724*
sqr(
RR*Tc_)/Pc_;
250 const scalar
b = 0.07780*
RR*Tc_/Pc_;
251 const scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
255 const scalar
B =
b*
p/(
RR*
T);
257 const scalar a2 =
B - 1;
258 const scalar a1 =
A - 2*
B - 3*
sqr(
B);
261 const scalar Q = (3*a1 - a2*a2)/9.0;
262 const scalar Rl = (9*a2*a1 - 27*
a0 - 2*a2*a2*a2)/54.0;
264 const scalar Q3 = Q*Q*Q;
265 const scalar
D = Q3 + Rl*Rl;
272 const scalar qm = 2*
sqrt(-Q);
273 const scalar r1 = qm*
cos(th/3.0) - a2/3.0;
279 root =
max(r1,
max(r2, r3));
284 const scalar D05 =
sqrt(
D);
285 const scalar S =
cbrt(Rl + D05);
296 root = S + Tl - a2/3.0;
303 template<
class Specie>
310 const scalar Tr =
T/Tc_;
311 const scalar a = 0.45724*
sqr(
RR*Tc_)/Pc_;
312 const scalar
b = 0.07780*
RR*Tc_/Pc_;
313 const scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
317 const scalar
B =
b*
p/(
RR*
T);
319 const scalar Z = this->Z(
p,
T);
322 const scalar
M = (
sqr(Z) + 2*
B*Z -
sqr(
B))/(Z -
B);
323 const scalar
N = ap*
B/(
b*
RR);
325 return this->
R()*
sqr(M -
N)/(
sqr(
M) - 2*
A*(Z +
B));
331 template<
class Specie>
337 scalar Y1 = this->
Y();
338 Specie::operator+=(pg);
340 if (
mag(this->
Y()) > SMALL)
343 const scalar Y2 = pg.Y()/this->
Y();
345 Tc_ = Y1*Tc_ + Y2*pg.Tc_;
346 Vc_ = Y1*Vc_ + Y2*pg.Vc_;
347 Zc_ = Y1*Zc_ + Y2*pg.Zc_;
348 Pc_ =
RR*Zc_*Tc_/Vc_;
349 omega_ = Y1*omega_ + Y2*pg.omega_;
354 template<
class Specie>
357 Specie::operator*=(
s);
364 template<
class Specie>
367 const PengRobinsonGas<Specie>& pg1,
368 const PengRobinsonGas<Specie>& pg2
373 static_cast<const Specie&>(pg1)
374 + static_cast<const Specie&>(pg2)
377 if (
mag(sp.Y()) < SMALL)
379 return PengRobinsonGas<Specie>
391 const scalar Y1 = pg1.Y()/sp.Y();
392 const scalar Y2 = pg2.Y()/sp.Y();
394 const scalar Tc = Y1*pg1.Tc_ + Y2*pg2.Tc_;
395 const scalar Vc = Y1*pg1.Vc_ + Y2*pg2.Vc_;
396 const scalar Zc = Y1*pg1.Zc_ + Y2*pg2.Zc_;
398 return PengRobinsonGas<Specie>
405 Y1*pg1.omega_ + Y2*pg2.omega_
411 template<
class Specie>
415 const PengRobinsonGas<Specie>& pg
418 return PengRobinsonGas<Specie>
420 s*
static_cast<const Specie&
>(pg),
430 template<
class Specie>
433 const PengRobinsonGas<Specie>& pg1,
434 const PengRobinsonGas<Specie>& pg2
439 static_cast<const Specie&>(pg1)
440 == static_cast<const Specie&>(pg2)
443 const scalar Y1 = pg1.Y()/sp.Y();
444 const scalar Y2 = pg2.Y()/sp.Y();
446 const scalar Tc = Y2*pg2.Tc_ - Y1*pg1.Tc_;
447 const scalar Vc = Y2*pg2.Vc_ - Y1*pg1.Vc_;
448 const scalar Zc = Y2*pg2.Zc_ - Y1*pg1.Zc_;
450 return PengRobinsonGas<Specie>
457 Y2*pg2.omega_ - Y1*pg1.omega_
dimensionedScalar Pr("Pr", dimless, laminarTransport)
dimensionedScalar acos(const dimensionedScalar &ds)
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...
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
scalar psi(scalar p, scalar T) const
Return compressibility rho/p [s^2/m^2].
dimensionedScalar sqrt(const dimensionedScalar &ds)
scalar Z(scalar p, scalar T) const
Return compression factor [].
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.
void operator*=(const scalar)
scalar H(const scalar p, const scalar T) const
Return enthalpy departure [J/kg].
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
scalar Cv(scalar p, scalar T) const
Return Cv departure [J/(kg K].
autoPtr< PengRobinsonGas > clone() const
Construct and return a clone.
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
scalar Cp(scalar p, scalar T) const
Return Cp departure [J/(kg K].
dimensionedScalar cos(const dimensionedScalar &ds)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
A class for handling words, derived from Foam::string.
dimensionedScalar cbrt(const dimensionedScalar &ds)
static autoPtr< PengRobinsonGas > New(const dictionary &dict)
constexpr scalar pi(M_PI)
PengRobinsonGas gas equation of state.
const dimensionedScalar Pstd
Standard pressure.
PengRobinsonGas(const Specie &sp, const scalar &Tc, const scalar &Vc, const scalar &Zc, const scalar &Pc, const scalar &omega)
Construct from components.
scalar S(const scalar p, const scalar T) const
Return entropy [J/(kg K)].
const Vector< label > N(dict.get< Vector< label >>("N"))
dimensionedScalar pow3(const dimensionedScalar &ds)
#define R(A, B, C, D, E, F, K, M)
PtrList< volScalarField > & Y
const dimensionedScalar a0
Bohr radius: default SI units: [m].
const scalar RR
Universal gas constant: default in [J/(kmol K)].
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
const dimensionedScalar & D
scalar E(const scalar p, const scalar T) const
Return internal energy departure [J/kg].
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 alpha
Fine-structure constant: default SI units: [].
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)