4. Equations of State¶
4.1. Single Component¶
4.1.1. Virial¶
Todo
merge docs with those in realgas.partial_molar_properties
.
There, the definitions of residual properties are displayed and here wee need only
write simplified forms for the specific equation of state.
4.1.1.1. Theory¶
The second order virial equation of state is [GP07]
Where the composition dependency of \(B\) is given by the exact mixing rule
where \(B_{ij}=B_{ji}\), and \(B_{ii}\) and \(B_{jj}\) are virial coefficients for the pure species
In this package, the useful correlation
or
So that, combining Equations (1) and (3), the compressibility can be calculated from dimensionless quantities as
where [SVanNessA05]
so that
so that the following derivatives can be computed as
Which allow the \(H^\text{R}\), \(S^\text{R}\), and \(G^\text{R}\) to be readily computed as follows [GP07]
The cross coefficients are calculated as
so that the cross derivatives can be computed as
or, equivalently,
4.1.1.2. Fugacity Coefficients¶
The Fugacity coefficients are calculated as
For the virial equation of state, this becomes [VanNessA82]
where both \(i\) and \(j\) indices run over all species
and
4.1.1.3. Residual Partial Molar Properties¶
The partial molar residual free energy of component k is
The partial molar residual enthalpy of component k is
where \(\frac{\partial \delta_{ij}}{\partial T_r}\) is given by (39) so that we obtain
Since
In terms of partial molar properties, then
By comparing Equation (16) and (17) it is observed that
where \(\frac{\partial \delta_{ij}}{\partial T_r}\) is given by (39)
The partial molar residual volume of component i is calculated as
which simplifies to
from which we obtain the intuitive result of
-
class
realgas.eos.virial.
Virial
(pow: callable = <ufunc 'power'>, exp: callable = <ufunc 'exp'>)[source]¶ - Parameters
pow (callable, optional) – function for computing power, defaults to numpy.power
exp (callable, optional) – function for computing logarithm, defaults to numpy.exp
-
class
realgas.eos.virial.
SecondVirial
(dippr_no: str = None, compound_name: str = None, cas_number: str = None, pow: callable = <ufunc 'power'>, **kwargs)[source]¶ Virial equation of state for one component. See [GP07][SVanNessA05]
-
calc_Z_from_dimensionless
(P_r, T_r)[source]¶ - Parameters
P_r – reduced pressure, dimensionless
T_r – reduced temperature, dimensionless
- Returns
Equation (4)
-
calc_Z_from_units
(P, T)[source]¶ - Parameters
P – pressure in Pa
T – temperature in K
- Returns
Equation (1)
-
ln_hat_phi_k_expr
(P, T)[source]¶ logarithm of fugacity coefficient
Note
single-component version
In this case, Equation (14) simplifies to
\[\ln\hat{\phi}_i = \frac{PB}{RT}\]- Parameters
P (float) – pressure in Pa
T (float) – temperature in K
-
plot_Z_vs_P
(T, P_min, P_max, symbol='o', ax=None, **kwargs)[source]¶ Plot compressibility as a function of pressure
- Parameters
T (float) – temperature [K]
P_min (float) – minimum pressure for plotting [Pa]
P_max (float) – maximum pressure for plotting [Pa]
phase (str) – phase type (liquid or vapor), defaults to vapor
symbol (str) – marker symbol, defaults to ‘o’
ax (plt.axis) – matplotlib axes for plotting, defaults to None
kwargs – keyword arguments for plotting
-
4.1.2. Cubic¶
4.1.2.1. Theory¶
The generic cubic equation of state is [GP07]
where \(\epsilon\) and \(\sigma\) are pure numbers (the same for all substances), and \(a(T)\) and \(b\) are given by the following equations
The compressibility factor can be calculated by solving the following equation
where
and
An iterative routine to calculate \(Z\) using Equation (21) is implemented, following [GP07], where
that is continued until the following is true
4.1.2.2. Residual Molar Properties¶
where the following functions have been defined
4.1.2.3. Fugacity Coefficients¶
The pure-component fugacity coefficient is defined as
So that, from Equation (27),
4.1.2.4. Mixtures¶
Warning
Not implemented yet!
Todo
Implement mixtures with cubic equations of state
-
class
realgas.eos.cubic.
Cubic
(sigma: float, epsilon: float, Omega: float, Psi: float, dippr_no: str = None, compound_name: str = None, cas_number: str = None, log: callable = <ufunc 'log'>, exp: callable = <ufunc 'exp'>, name: str = 'cubic', **kwargs)[source]¶ Generic Cubic Equation of State
- Parameters
sigma – Parameter defined by specific equation of state, \(\sigma\)
epsilon – Parameter defined by specific equation of state, \(\epsilon\)
Omega – Parameter defined by specific equation of state, \(\Omega\)
Psi – Parameter defined by specific equation of state, \(\Psi\)
tol – tolerance for iteration (see Equation (25)), set to 0.01
log – function for computing natural log, defaults to
numpy.log
exp – function for computing exponential, defaults to
numpy.exp
-
G_R_RT_expr
(P, V, T)[source]¶ Dimensionless residual gibbs
- Parameters
P – pressure in Pa
V – molar volume in m**3/mol
T – temperature in K
- Returns
\(\frac{G^\text{R}}{RT}\), see Equation (27)
-
H_R_RT_expr
(P, V, T)[source]¶ Dimensionless residual enthalpy
- Parameters
P – pressure in Pa
V – molar volume in m**3/mol
T – temperature in K
- Returns
\(\frac{H^\text{R}}{RT}\), see Equation (26)
-
I_expr
(P, V, T)[source]¶ - Parameters
P – pressure in Pa
V – molar volume in m**3/mol
T – temperature in K
- Returns
\(I\) (see Equation (29))
-
S_R_R_expr
(P, V, T)[source]¶ Dimensionless residual entropy
- Parameters
P – pressure in Pa
V – molar volume in m**3/mol
T – temperature in K
- Returns
\(\frac{S^\text{R}}{R}\), see Equation (28)
-
Z_right_hand_side
(Z, beta, q)[source]¶ Estimate of compressibility of vapor [GP07], used for iterative methods
- Returns
RHS of residual, see Equation (30)
-
alpha_expr
(T_r)[source]¶ An empirical expression, specific to a particular form of the equation of state
- Parameters
T_r – reduced temperature (T/Tc), dimensionless
- Returns
\(\alpha (T_r)\) see Table ??
-
beta_expr
(T, P)[source]¶ - Parameters
T – temperature in K
P – pressure in Pa
- Returns
\(\beta\) (see Equation (22))
-
cardano_constants
(T, P)[source]¶ - Parameters
T – temperature [T]
P – pressure [Pa]
- Returns
cardano constants p, q
- Return type
tuple
-
coefficients
(T, P)[source]¶ Polynomial oefficients for cubic equation of state
\[Z^3 c_0 + Z^2c_1 + Z*c_2 + c_3 = 0\]- Returns
(c_0, c_1, c_2, c_3)
-
d_ln_alpha_d_ln_Tr
(T_r)[source]¶ - Parameters
T_r – reduced temperature [dimensionless]
- Returns
Expression for \(\frac{\mathrm{d} \ln \alpha(T_\mathrm{r})}{\mathrm{d} \ln T_\mathrm{r}}\)
-
hat_phi_i_expr
(*args)[source]¶ expression for fugacity coefficient :returns: \(\exp\left({\ln{\hat{\phi}_i}}\right)\)
-
iterate_to_solve_Z
(T, P) → float[source]¶ Iterate to compute \(Z\) using Equation (24) util termination condition (Equation (25) is met)
- Parameters
T – temperature in K
P – pressure in Pa
- Returns
compressibility factor
-
ln_hat_phi_k_expr
(P, V, T)[source]¶ - Parameters
P – pressure in Pa
T – temperature in K
V – molar volume in m**3/mol
- Returns
\(\ln{\hat{\phi}_k\), see Equation \(ln_hat_phi_i\)
-
num_roots
(T, P)[source]¶ Find number of roots
- Parameters
T – temperature in K
P – pressure in Pa
- Returns
number of roots
-
plot_Z_vs_P
(T: float, P_min: float, P_max: float, symbol='o', ax: matplotlib.pyplot.axis = None, fig: matplotlib.pyplot.figure = None, **kwargs)[source]¶ Plot compressibility as a function of pressure
- Parameters
T – temperature [K]
P_min – minimum pressure for plotting [Pa]
P_max – maximum pressure for plotting [Pa]
symbol – marker symbol, defaults to ‘o’
ax – matplotlib axes for plotting, defaults to None
kwargs – keyword arguments for plotting
-
class
realgas.eos.cubic.
RedlichKwong
(**kwargs)[source]¶ Redlich-Kwong Equation of State [RK49]
This Equation of state has the following parameters [SVanNessA05]
Symbol
Value
\(\alpha(T_\text{r})\)
\(1/\sqrt{T_\text{r}}\)
\(\sigma\)
1
\(\epsilon\)
0
\(\Omega\)
0.08664
\(\Psi\)
0.42748
>>> from realgas.eos.cubic import RedlichKwong >>> model = RedlichKwong(compound_name='Propane') >>> model.sigma 1 >>> model.epsilon 0 >>> model.Omega 0.08664 >>> model.Psi 0.42748
-
class
realgas.eos.cubic.
SoaveRedlichKwong
(**kwargs)[source]¶ Soave Redlich-Kwong Equation of State [Soa72]
This equation of state has the following parameters
Symbol
Value
\(\alpha(T_\text{r})\)
\(\left[1 + f_w(1-\sqrt{T_\text{r}})\right]^2\)
\(\sigma\)
\(1\)
\(\epsilon\)
\(0\)
\(\Omega\)
\(0.08664\)
\(\Psi\)
\(0.42748\)
where
(32)¶\[f_w = 0.480 + 1.574\omega - 0.176\omega^2\]>>> from realgas.eos.cubic import SoaveRedlichKwong >>> model = SoaveRedlichKwong(compound_name='Water') >>> model.Omega 0.08664 >>> model.sigma 1 >>> model.epsilon 0 >>> model.Psi 0.42748 >>> model.f_w_rule(0.) 0.48 >>> model.f_w_rule(1.) 1.878
- Parameters
f_w (float, derived) – empirical expression used in \(\\alpha\) [dimensionless], see Equation (32)
-
alpha_expr
(T_r)[source]¶ An empirical expression, specific to a particular form of the equation of state
- Parameters
T_r – reduced temperature (T/Tc), dimensionless
- Returns
\(\alpha (T_r)\) see Table ??
-
class
realgas.eos.cubic.
PengRobinson
(**kwargs)[source]¶ Peng-Robinson Equation of State [PR76]
This equation of state has the following parameters
Symbol
Value
\(\alpha(T_\text{r})\)
\(\left[1 + f_w(1-\sqrt{T_\text{r}})\right]^2\)
\(\sigma\)
\(1 + \sqrt{2}\)
\(\epsilon\)
\(1 - \sqrt{2}\)
\(\Omega\)
\(0.07780\)
\(\Psi\)
\(0.45724\)
where
(33)¶\[f_w = 0.480 + 1.574\omega - 0.176\omega^2\]>>> from realgas.eos.cubic import PengRobinson >>> model = PengRobinson(compound_name='Water') >>> model.Omega 0.0778 >>> model.sigma 2.4142135 >>> model.epsilon -0.414213 >>> model.Psi 0.45724 >>> model.f_w_rule(0.) 0.37464 >>> model.f_w_rule(1.) 1.64698
- Parameters
f_w (float, derived) – empirical expression used in \(\\alpha\) [dimensionless?]
4.2. Multicomponent¶
4.2.1. Virial¶
-
class
realgas.eos.virial.
MixingRule
(pow: callable = <ufunc 'power'>, exp: callable = <ufunc 'exp'>)[source]¶ Van der Walls mixing rule
combining rules from [PLdeAzevedo86]
(34)¶\[\omega_{ij} = \frac{\omega_i + \omega_j}{2}\](35)¶\[T_{\text{c}ij} = \sqrt{T_{\text{c}i}T_{\text{c}j}}(1-k_{ij})\](36)¶\[P_{\text{c}ij} = \frac{Z_{\text{c}ij}RT_{\text{c}ij}}{V_{\text{c}ij}}\](37)¶\[Z_{\text{c}ij} = \frac{Z_{\text{c}i} + Z_{\text{c}j}}{2}\](38)¶\[V_{\text{c}ij} = \left(\frac{V_{\text{c}i}^{1/3} + V_{\text{c}j}^{1/3}}{2}\right)^3\]-
T_cij_rule
(T_ci, T_cj, k_ij)[source]¶ - Parameters
T_ci – critical temperature of component i [K]
T_cj – ** j [K]
k_ij – k_ij parameter
- Returns
Equation (35)
-
V_cij_rule
(V_ci, V_cj)[source]¶ - Parameters
V_ci – critical molar volume of component i [m**3/mol]
V_cj – critical molar volume of component j [m**3/mol]
- Returns
Equation eq:Vc_combine
-
-
class
realgas.eos.virial.
SecondVirialMixture
(pow: callable = <ufunc 'power'>, exp: callable = <ufunc 'exp'>, **kwargs)[source]¶ Second virial with mixing rule from
MixingRule
Note
can only input both custom critical properties or both from DIPPR–cant have mixed at the moment
- Parameters
pow – function to calculate power, defaults to numpy.power
exp – function to calculate exponential,d efaults to numpy.exp
kwargs – key-word arguments to pass to
RealMixture
-
B_ij_expr
(i: int, j: int, T)[source]¶ - Parameters
i – index of first component
j – index of second component
T – temperature [K]
- Returns
Equation (12)
-
B_mix_expr
(y_k: List[Union[float, Any]], T)[source]¶ - Parameters
y_k – mole fractions of each component \(k\)
T – temperature in K
- Returns
Equation (2)
-
M_R_dimensionless
(method: callable, ys: List[Union[float, Any]], P: float, T: float)[source]¶ Residual property of \(X\) for mixture.
Similar to Equation (3) but in dimensionless form
- Parameters
method (callable) – function to compute partial molar property of compound
-
bar_GiR_RT
(cas_k: str, ys: List[Union[float, Any]], P, T)[source]¶ Dimensionless residual partial molar free energy of component \(i\)
- Parameters
cas_k – cas number of component k
ys – mole fractions
P – pressure in Pa
T – temperature in K
- Returns
Equation (16)
-
bar_HiR_RT
(cas_k: str, ys: List[Union[float, Any]], P, T)[source]¶ Dimensionless residual partial molar enthalpy of component \(i\)
- Parameters
cas_k – cas number for component of interest
ys – mole fractions
P – pressure in Pa
T – temperature in K
-
bar_HiR_star
(T_star, cas_k: str, ys: List[Union[float, Any]], P, T)[source]¶ Returns the scaled entropy useful for computations \(\bar{H}_i^{\text{R},\star}\), as defined in Equation (10)
- Parameters
cas_k – cas number of component k
ys – mole fractions
P – pressure in Pa
T – temperature in K
-
bar_SiR_R
(cas_k: str, ys: List[Union[float, Any]], P, T)[source]¶ Dimensionless residual partial molar entropy of component \(i\)
- Parameters
cas_k – cas number for component of interest
ys – mole fractions
P – pressure in Pa
T – temperature in K
- Returns
Equation (18)
-
bar_ViR_RT
(cas_k: str, ys: List[Union[float, Any]], P, T)[source]¶ residual Partial molar volume for component i
Note
This expression does not depend on \(P\)
- Parameters
cas_k – cas number for component of interest
ys – mole fractions
P – pressure in Pa
T – temperature in K
- Returns
Equation (19)
-
calc_Z_from_units
(y_k: List, P, T)[source]¶ - Parameters
y_k – mole fractions of each component \(k\)
P – pressure in Pa
T – temperature in K
- Returns
Equation (1)
-
d_Bij_d_Trij
(i: int, j: int, T)[source]¶ - Parameters
i – index for component i
j – index for componen t j
T – temperature in K
- Returns
Equation (13)
-
d_dij_d_Tr
(i, j, T)[source]¶ - (39)¶\[\frac{\partial \delta_{ij}}{\partial T_{\text{r}ij}} = 2\frac{\partial B_{ij}}{\partial T_{\text{r}ij}}-\frac{\partial B_{ii}}{\partial T_{\text{r}ij}}-\frac{\partial B_{jj}}{\partial T_{\text{r}ij}}\]
Todo
test this with symbolic differentiation of d_ik expression
- Parameters
i – index for component i
j – index for componen t j
T – temperature [K]
-
d_ik_expr
(i: int, k: int, T)[source]¶ - Parameters
i – index of component i
k – index of component k
T – temperature [K]
- Returns
Equation (15)
-
fugacity_i_expr
(cas_i: str, ys: List[Union[float, Any]], P, T)[source]¶ Fugacity of component i in mixture \(f_i=\hat{\phi}_i y_i P\)
- Parameters
cas_i – cas number for component of interest
ys – mole fractions
P – pressure in Pa
T – temperature in K
-
get_w_Tc_Pc
(i: int, j=None)[source]¶ Returns critical constants for calculation based off of whetner i = j or not
- Returns
(\(w\), \(T_c\), \(P_c\))
- Return type
tuple
-
ln_hat_phi_k_expr
(k: int, ys: List[Union[float, Any]], P, T)[source]¶ logarithm of fugacity coefficient
Note
The order of
ys
corresponds to the order of components made during initialization- Parameters
k – index of component k
ys – mole fractions ordered as component order
P – pressure in Pa
T – temperature in K
- Returns
Equation (14)
-
plot_residual_HSG
(P, T, ax=None, fig=None) → Tuple[matplotlib.pyplot.figure, matplotlib.pyplot.subplot][source]¶ Plot dimensionless residual properties as a function of mole fraction
- Parameters
P – pressure in Pa
T – Temperature in K
ax – matplotlib ax, defaults to None
fig – matplotlib figure, defautls to None