31 template<
class Specie>
53 template<
class Specie>
69 template<
class Specie>
80 template<
class Specie>
91 template<
class Specie>
107 template<
class Specie>
114 scalar Z = this->Z(p, T);
115 return p/(Z*this->
R()*
T);
119 template<
class Specie>
124 scalar B = 0.07780*Pr/Tr;
125 scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
128 scalar Z = this->Z(p, T);
134 - 2.078*(1 + kappa)*
sqrt(alpha)
135 *
log((Z + 2.414*B)/(Z - 0.414*B))
140 template<
class Specie>
144 scalar a = 0.45724*
sqr(
RR*Tc_)/Pc_;
145 scalar
b = 0.07780*
RR*Tc_/Pc_;
146 scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
149 scalar A = a*alpha*p/
sqr(
RR*T);
150 scalar B = b*p/(
RR*
T);
152 scalar Z = this->Z(p, T);
154 scalar ap = kappa*a*(kappa/Tc_ - (1 +
kappa)/
sqrt(T*Tc_));
157 scalar
M = (
sqr(Z) + 2*B*Z -
sqr(B))/(Z - B);
158 scalar
N = ap*B/(b*
RR);
160 const scalar root2 =
sqrt(2.0);
163 app*(T/(2*root2*
b))*
log((Z + (root2 + 1)*B)/(Z - (root2 - 1)*B))
164 +
RR*
sqr(M - N)/(
sqr(M) - 2*A*(Z + B))
169 template<
class Specie>
178 scalar B = 0.07780*Pr/Tr;
179 scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
181 scalar Z = this->Z(p, T);
189 *
log((Z + 2.414*B)/(Z - 0.414*B))
194 template<
class Specie>
201 scalar Z = this->Z(p, T);
203 return 1.0/(Z*this->
R()*
T);
207 template<
class Specie>
215 scalar a = 0.45724*
sqr(
RR*Tc_)/Pc_;
216 scalar
b = 0.07780*
RR*Tc_/Pc_;
217 scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
220 scalar A = a*alpha*p/
sqr(
RR*T);
221 scalar B = b*p/(
RR*
T);
224 scalar a1 = A - 2*B - 3*
sqr(B);
227 scalar Q = (3*a1 - a2*a2)/9.0;
228 scalar Rl = (9*a2*a1 - 27*a0 - 2*a2*a2*a2)/54.0;
231 scalar D = Q3 + Rl*Rl;
238 scalar qm = 2*
sqrt(-Q);
239 scalar r1 = qm*
cos(th/3.0) - a2/3.0;
243 root =
max(r1,
max(r2, r3));
248 scalar D05 =
sqrt(D);
249 scalar S =
pow(Rl + D05, 1.0/3.0);
253 Tl = -
pow(
mag(Rl - D05), 1.0/3.0);
257 Tl =
pow(Rl - D05, 1.0/3.0);
260 root = S + Tl - a2/3.0;
267 template<
class Specie>
275 scalar a = 0.45724*
sqr(
RR*Tc_)/Pc_;
276 scalar
b = 0.07780*
RR*Tc_/Pc_;
277 scalar
kappa = 0.37464 + 1.54226*omega_ - 0.26992*
sqr(omega_);
280 scalar A = alpha*a*p/
sqr(
RR*T);
281 scalar B = b*p/(
RR*
T);
283 scalar Z = this->Z(p, T);
285 scalar ap = kappa*a*(kappa/Tc_ - (1 +
kappa)/
sqrt(T*Tc_));
286 scalar
M = (
sqr(Z) + 2*B*Z -
sqr(B))/(Z - B);
287 scalar
N = ap*B/(b*
RR);
289 return RR*
sqr(M - N)/(
sqr(M) - 2*A*(Z + B));
295 template<
class Specie>
296 inline void Foam::PengRobinsonGas<Specie>::operator+=
301 scalar molr1 = this->nMoles();
302 Specie::operator+=(pg);
304 molr1 /= this->nMoles();
305 scalar molr2 = pg.nMoles()/this->nMoles();
307 Tc_ = molr1*Tc_ + molr2*pg.Tc_;
308 Vc_ = molr1*Vc_ + molr2*pg.Vc_;
309 Zc_ = molr1*Zc_ + molr2*pg.Zc_;
310 Pc_ =
RR*Zc_*Tc_/Vc_;
311 omega_ = molr1*omega_ + molr2*pg.omega_;
315 template<
class Specie>
316 inline void Foam::PengRobinsonGas<Specie>::operator-=
321 scalar molr1 = this->nMoles();
323 Specie::operator-=(pg);
325 molr1 /= this->nMoles();
326 scalar molr2 = pg.nMoles()/this->nMoles();
328 Tc_ = molr1*Tc_ - molr2*pg.Tc_;
329 Vc_ = molr1*Vc_ - molr2*pg.Vc_;
330 Zc_ = molr1*Zc_ - molr2*pg.Zc_;
331 Pc_ =
RR*Zc_*Tc_/Vc_;
332 omega_ = molr1*omega_ - molr2*pg.omega_;
336 template<
class Specie>
339 Specie::operator*=(s);
346 template<
class Specie>
353 scalar nMoles = pg1.nMoles() + pg2.nMoles();
354 scalar molr1 = pg1.nMoles()/nMoles;
355 scalar molr2 = pg2.nMoles()/nMoles;
357 const scalar Tc = molr1*pg1.Tc_ + molr2*pg2.Tc_;
358 const scalar Vc = molr1*pg1.Vc_ + molr2*pg2.Vc_;
359 const scalar Zc = molr1*pg1.Zc_ + molr2*pg2.Zc_;
363 static_cast<const Specie&
>(pg1)
364 + static_cast<const Specie&>(pg2),
369 molr1*pg1.omega_ + molr2*pg2.omega_
374 template<
class Specie>
381 scalar nMoles = pg1.nMoles() + pg2.nMoles();
382 scalar molr1 = pg1.nMoles()/nMoles;
383 scalar molr2 = pg2.nMoles()/nMoles;
385 const scalar Tc = molr1*pg1.Tc_ + molr2*pg2.Tc_;
386 const scalar Vc = molr1*pg1.Vc_ + molr2*pg2.Vc_;
387 const scalar Zc = molr1*pg1.Zc_ + molr2*pg2.Zc_;
391 static_cast<const Specie&
>(pg1)
392 - static_cast<const Specie&>(pg2),
397 molr1*pg1.omega_ - molr2*pg2.omega_
402 template<
class Specie>
411 s*
static_cast<const Specie&
>(pg),
421 template<
class Specie>
scalar cpMcv(scalar p, scalar T) const
Return (cp - cv) [J/(kmol K].
static autoPtr< PengRobinsonGas > New(Istream &is)
dimensionedScalar Pr("Pr", dimless, laminarTransport)
dimensionedScalar acos(const dimensionedScalar &ds)
dimensionedScalar log(const dimensionedScalar &ds)
scalar cp(scalar p, scalar T) const
Return cp departure [J/(kmol K].
A list of keyword definitions, which are a keyword followed by any number of values (e...
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
scalar h(const scalar p, const scalar T) const
Return enthalpy departure [J/kmol].
dimensionedScalar sqrt(const dimensionedScalar &ds)
void operator*=(const scalar)
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
autoPtr< PengRobinsonGas > clone() const
Construct and return a clone.
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.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
scalar s(const scalar p, const scalar T) const
Return entropy [J/(kmol K)].
dimensionedScalar cos(const dimensionedScalar &ds)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
A class for handling words, derived from string.
scalar Z(scalar p, scalar T) const
Return compression factor [-].
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
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.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
#define R(A, B, C, D, E, F, K, M)
scalar psi(scalar p, scalar T) const
Return compressibility rho/p [s^2/m^2].
const dimensionedScalar a0
Bohr radius: default SI units: [m].
dimensioned< scalar > mag(const dimensioned< Type > &)
const scalar RR
Universal gas constant (default in [J/(kmol K)])
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].