Benedict-Webb-Rubin Equation of state

The Benedict-Webb-Rubin equation of state (EOS) is an early fluid model used for describing real gas behavior¹ developed by M. Benedict, G.B. Webb and L.C. Rubin in 1940. The original BWR EOS has been modified to more accurately predict real gas behavior by, among others, Starling², Jacobsen et al.³ and was also modified to describe multi-phase fluid systems by Starling in 1966⁴. Unlike the cubic EOS models like the Peng-Robinson (PR) and Soave-Redlich-Kwong (SRK) EOS models, the BWR models do not solve a polynomial representation of the Z-factor to describe the phase behavior of the fluid. The various BWR EOS models are so-called Virial EOS models.

Original BWR EOS

The original BWR EOS model can be described by

$\begin{equation}\label{eq:original_bwr} p = \rho RT + (B_0RT-A_0-\frac{C_0}{T^2})\rho^2 +(bRT-a)\rho^3 + \alpha a \rho^6 +\frac{c\rho^3}{T^2}(1+\gamma\rho^2)\exp{(-\gamma\rho^2)} \end{equation}$

where the parameters $a$, $b$, $c$, $A$, $B$, $C$, $\alpha$, $\gamma$ are all empirical coefficients. The model parameters behave similar to those in a cubic EOS model and have their own mixing rules. These mixing rules are given in the mixing rules section below.

Modified BWR EOS

The modified BWR model by Starling in 1970 introduced three new parameters ($D_0$, $E_0$ and $d$) and is given by

$\begin{equation}\label{eq:starling_bwr} p = p_1 + p_2 + p_3 + p_4 + p_5 \end{equation}$

where

$\begin{equation}\label{eq:starling_bwr_1} p_1 = \rho RT \end{equation}$

$\begin{equation}\label{eq:starling_bwr_2} p_2 = (B_0RT-A_0-\frac{C_0}{T^2} + \frac{D_0}{T^3} - \frac{E_0}{T^4})\rho^2 \end{equation}$

$\begin{equation}\label{eq:starling_bwr_3} p_3 = (bRT-a-\frac{d}{T})\rho^3 \end{equation}$

$\begin{equation}\label{eq:starling_bwr_4} p_4 = \alpha (a + \frac{d}{T}) \rho^6 \end{equation}$

$\begin{equation}\label{eq:starling_bwr_5} p_5 = \frac{c\rho^3}{T^2}(1+\gamma\rho^2)\exp{(-\gamma\rho^2)} \end{equation}$

Standing provides the parameter values for several light hydrocarbon (normal paraffins) and also gives a correlation for estimating the parameters for light hydrocarbons using the acentric factor.

McCarty Modification of BWR EOS

In 1974 McCarty expanded the BWR EOS model by introducing even more coefficients to the general structure of the model⁵. The McCarty modified BWR EOS is given by

$\begin{equation}\label{eq:mccarty_bwr} p = \sum_{n=1}^9a_n(T)\rho^n + \sum_{n=10}^{15}a_n(T)\rho^{2n-17}\exp{\{-\gamma\rho^2\}} \end{equation}$

where $a_n$ are functions of temperature (similar to $p_1$ to $p_3$) and $\gamma$ is an empirical constant.

The list of $a_n$ can be found in McCarty's paper or in Pedersen's book Phase Behavior of Petroleum Reservoir Fluids⁶.

Mixing Rules

The mixing rules for the parameters in the original and modified BWR EOS are given by

$\begin{equation}\label{eq:A_mix} A_0 = \sum_{i=1}^N \sum_{i=1}^N u_i u_j \sqrt{A_{0i}A_{0j}}(1-k_{ij}) \end{equation}$

$\begin{equation}\label{eq:B_mix} B_0 = \sum_{i=1}^N u_i B_{0i} \end{equation}$

$\begin{equation}\label{eq:C_mix} C_0 = \sum_{i=1}^N \sum_{i=1}^N u_i u_j \sqrt{C_{0i}C_{0j}}(1-k_{ij}) \end{equation}$

$\begin{equation}\label{eq:D_mix} D_0 = \sum_{i=1}^N \sum_{i=1}^N u_i u_j \sqrt{D_{0i}D_{0j}}(1-k_{ij}) \end{equation}$

$\begin{equation}\label{eq:E_mix} E_0 = \sum_{i=1}^N \sum_{i=1}^N u_i u_j \sqrt{E_{0i}E_{0j}}(1-k_{ij}) \end{equation}$

$\begin{equation}\label{eq:a_mix} a = [\sum_{i=1}^N u_i a_i^{1/3}]^3 \end{equation}$

$\begin{equation}\label{eq:b_mix} b = [\sum_{i=1}^N u_i b_i^{1/3}]^3 \end{equation}$

$\begin{equation}\label{eq:c_mix} c = [\sum_{i=1}^N u_i c_i^{1/3}]^3 \end{equation}$

$\begin{equation}\label{eq:d_mix} d = [\sum_{i=1}^N u_i d_i^{1/3}]^3 \end{equation}$

$\begin{equation}\label{eq:alpha_mix} \alpha = [\sum_{i=1}^N u_i \alpha_i^{1/3}]^3 \end{equation}$

$\begin{equation}\label{eq:gamma_mix} \gamma = [\sum_{i=1}^N u_i \gamma_i^{1/3}]^3 \end{equation}$

Soave's Modified BWR EOS

In 1996 Soave published a modification to the BWR (SBWR) EOS model⁷ shown by

$\begin{equation}\label{eq:sbwr_96} Z=\frac{p}{RT\rho}=1 + B \cdot \rho+C\cdot \rho^2 + D \cdot \rho^4 + E \cdot \rho^2 \cdot (1+F\cdot \rho^2)\cdot \exp{(-F\cdot \rho^2)} \end{equation}$

In 1999 Soave⁸ further simplified the Eq. \eqref{eq:sbwr_96}, by removing the $C\cdot \rho^2$ term resulting in

$\begin{equation}\label{eq:sbwr_99} Z=\frac{p}{RT\rho}=1 + B \cdot \rho + D \cdot \rho^4 + E \cdot \rho^2 \cdot (1+F\cdot \rho^2)\cdot \exp{(-F\cdot \rho^2)} \end{equation}$

Even though Eq. \eqref{eq:sbwr_99} has one less term, Soave shows that the updated EOS is a more accurate model. Soave also gives all the relevant equations needed to calculate the fugacity, using the SBWR EOS model.

One advantageous feature of the SBWR that is that the EOS coefficients (B, D, E and F) are calculated using a composition ($u_i$) and mixing rules for the component properties, namely the critical pressure ($p_c$), critical temperature ($T_c$), critical Z-factor ($Z_c$), and acentric factor ($\omega$). The advantage of this approach is that it does not need any additional measured component specific properties.

The mixing rules for the component properties are given by

$\begin{equation}\label{eq:mixing_rule_a} Z_{cm} = \frac{\sum_i u_i Z_{ci} T_{ci} / p_{ci}}{\sum_i u_i T_{ci} / p_{ci}} \end{equation}$

$\begin{equation}\label{eq:mixing_rule_b} T_{cm} = \frac{S_1}{(\sqrt{S_2}+\sqrt{S_3})^2} \end{equation}$

$\begin{equation}\label{eq:mixing_rule_c} p_{cm} = T_{cm}/S_3 \end{equation}$

$\begin{equation}\label{eq:mixing_rule_d} \omega_{m} = \sqrt{\frac{S_2}{S_3}} / 1.25 \end{equation}$

where

$\begin{equation}\label{eq:summa_1} S_1 = \sum_i \sum_j u_i u_j (1-k_{ij}) \frac{T_{ci}T_{cj}}{\sqrt{p_{ci}p_{cj}}}(1+m_i)(1+m_j) \end{equation}$

$\begin{equation}\label{eq:summa_2} S_2 = \sum_i \sum_j u_i u_j (1-k_{ij}) \sqrt{\frac{T_{ci}T_{cj}}{p_{ci}p_{cj}}} m_i m_j \end{equation}$

$\begin{equation}\label{eq:summa_3} S_3 = \sum_i u_i \frac{T_{ci}}{p_{ci}} \end{equation}$

and

$\begin{equation}\label{eq:m_eq} m_i = 1.25 \cdot \omega_i \end{equation}$

With the mixture properties calculated using Eq. \eqref{eq:mixing_rule_a} to \eqref{eq:mixing_rule_d}, the EOS properties can be calculated by

$\begin{equation}\label{eq:SBWR_B} B = (\frac{RT_{cm}}{p_{cm}}) \beta \end{equation}$

$\begin{equation}\label{eq:SBWR_D} D = (\frac{RT_{cm}}{p_{cm}})^4 \delta \end{equation}$

$\begin{equation}\label{eq:SBWR_E} E = (\frac{RT_{cm}}{p_{cm}})^2 \epsilon \end{equation}$

$\begin{equation}\label{eq:SBWR_F} F = (\frac{RT_{cm}}{p_{cm}})^2 \phi \end{equation}$

where

$\begin{equation}\label{eq:SBWR_beta} \beta = \beta_c + 0.422(1-T_r^{-1.6}) + 0.234 \omega_{m} (1-T_r^{-3}) \end{equation}$

$\begin{equation}\label{eq:SBWR_delta} \delta = \delta_c [1 + d_1 (T_r^{-1}-1) + d_2 (T_r^{-1}-1)^2] \end{equation}$

$\begin{equation}\label{eq:SBWR_epsilon} \epsilon = \epsilon_c + e_1 (T_r^{-1}-1) + e_2 (T_r^{-1}-1)^2 + e_3 (T_r^{-1}-1)^3 \end{equation}$

$\begin{equation}\label{eq:SBWR_phi} \phi = f \cdot Z_{cm}^2 \end{equation}$

and

$\begin{equation}\label{eq:d_1} d_1 = 0.4912 + 0.6478 \cdot \omega_m \end{equation}$

$\begin{equation}\label{eq:d_2} d_2 = 0.300 + 0.3619 \cdot \omega_m \end{equation}$

$\begin{equation}\label{eq:e_1} e_1 = 0.0841 + 0.1318 \cdot \omega_m + 0.0018 \cdot \omega_m^2 \end{equation}$

$\begin{equation}\label{eq:e_2} e_2 = 0.0750 + 0.2408 \cdot \omega_m - 0.0140 \cdot \omega_m^2 \end{equation}$

$\begin{equation}\label{eq:e_3} e_3 = - 0.0065 + 0.1798 \cdot \omega_m - 0.0078 \cdot \omega_m^2 \end{equation}$

$\begin{equation}\label{eq:beta_c} \beta_c = b_c \cdot Z_{cm} \end{equation}$

$\begin{equation}\label{eq:delta_c} \delta_c = d_c \cdot Z_{cm}^4 \end{equation}$

$\begin{equation}\label{eq:epsilon_c} \epsilon_c = e_c \cdot Z_{cm}^2 \end{equation}$

Finally, the equations for $b_c$, $d_c$, $e_c$, and $f_c=f$ are calculated from

$\begin{equation}\label{eq:f} f = 0.77 \end{equation}$

$\begin{equation}\label{eq:e_c} e_c = \frac{2 - 5 Z_{cm}}{(1+f+3f^2-2f^3)\exp{(-f)}} \end{equation}$

$\begin{equation}\label{eq:d_c} d_c = \frac{1-2 Z_{cm} - e_c (1+f-2f^2)\exp{(-f)}}{3} \end{equation}$

$\begin{equation}\label{eq:b_c} b_c = Z_{cm} - 1 - d_c - e_c (1+f)\exp{(-f)} \end{equation}$

M. Benedict, G. B. Webb, and L. C. Rubin. An empirical equation for thermodynamic properties of light hydrocarbons and their mixtures i. methane, ethane, propane and n-butane. The Journal of Chemical Physics, 8:334–345, 1940. doi:https://doi.org/10.1063/1.1750658. ↩
K. E. Starling and J. E. Powers. Enthalpy of mixtures by modified bwr equation. Industrial & Engineering Chemistry Fundamentals, 9:531–537, 1970. ↩
R. T. Jacobsen and R. B. Stewart. Thermodynamic properties of nitrogen including liquid and vapor phases from 63k to 2000k with pressures to 10,000 bar. Journal of Physical and Chemical Reference Data, 2:757–922, 1973. doi:https://doi.org/10.1063/1.3253132. ↩
K.E. Starling. A new approach for determining equation-of-state parameters using phase equilibria data. Society of Petroleum Engineers Journal, 6:paper SPE–1481–PA, 1966. doi:https://doi.org/10.2118/1481-PA. ↩
R. D. McCarty. A modified benedict-webb-rubin equation of state for methane using recent experimental data. Cryogenics, 14:276–280, 1974. doi:https://doi.org/10.1016/0011-2275 $74$ 90228-8. ↩
K. S. Pedersen, P. L. Christensen, and J. A. Shaikh. Phase behavior of petroleum reservoir fluids. CRC press, 2014. ↩
G. S. Soave. A noncubic equation of state for the treatment of hydrocarbon fluids at reservoir conditions. Industrial & engineering chemistry research, 34:3981–3994, 1995. ↩
G. S. Soave. An effective modification of the benedict–webb–rubin equation of state. Fluid Phase Equilibria, 164:157–172, 1999. ↩