Peng-Robinson EOS (1976)

The Peng-Robinson EOS has become the most popular equation of state for natural gas systems in the petroleum industry. During the decade of the 1970’s, D. Peng was a PhD student of Prof. D.B. Robinson at the University of University of Alberta (Edmonton, Canada). The Canadian Energy Board sponsored them to develop an EOS specifically focused on natural gas systems. When you compare the performance of the PR EOS and the SRK EOS, they are pretty close to a tie; they are “neck to neck,” except for a slightly better behavior by the PR EOS at the critical point. A slightly better performance around critical conditions makes the PR EOS somewhat better suited to gas/condensate systems.

Peng and Robinson introduced the following modified vdW EOS:

(P + \frac{α a}{{\tilde{v}}^{2} + 2 b \tilde{v} - b^{2}}) (\tilde{v} - b) = R T

(11.1a)

or, explicitly in pressure,

P = \frac{R T}{\tilde{v} - b} - \frac{α a}{{\tilde{v}}^{2} + 2 b \tilde{v} - b^{2}}

(11.1b)

where:

Peng and Robinson conserved the temperature dependency of the attractive term and the acentric factor introduced by Soave. However, they presented different fitting parameters to describe this dependency (equation 4.11c), and further manipulated the denominator of the pressure correction (attractive) term. As we have seen before, coefficients “a” and “b” are made functions of the critical properties by imposing the criticality conditions. This yields:

α = {⌊ 1 + (0.37464 + 1.54226 ω - 0.26992 ω^{2}) (1 - \sqrt{T_{r}}) ⌋}^{2}

(11.1c)

a = 0.45724 \frac{R^{2} T_{c}^{2}}{P_{c}}

(11.2a)

The PR cubic expression in Z becomes:

b = 0.07780 \frac{R T_{c}}{P_{c}}

(11.2b)

where:

Z^{3} - (1 - B) Z^{2} + (A - 2 B - 3 B^{2}) Z - (A B - B^{2} - B^{3}) = 0

(11.3a)

A = \frac{α a P}{R^{2} T^{2}}

(11.3b)

B = \frac{b P}{R T}

(11.3c)

Similar to SRK, the PR mixing rules are:

{(α a)}_{m} = \sum^{​} \sum^{​} y_{i} y_{j} {(α a)}_{i j}; {(α a)}_{i j} = \sqrt{{(α a)}_{i} {(α a)}_{j}} (1 - k_{i j})

(11.4a)

b_{m} = \sum_{i} y_{i} b_{i}

(11.4b)

where binary interaction parameters (k_ij) again play the important empirical role of helping to better fit experimental data. Due to the empirical character of these interaction parameters, k_ij’s calculated for PR EOS are unlikely to be the same as the k_ij’s calculated for SRK EOS for the same pair of molecules.

Comparative Assessment of RK, SRK, and PR EOS

Over the years, these EOS have been tested, and some comparisons can be made. As an engineer, you have to be able to decide which EOS best fits your purposes.

Redlich Kwong EOS:

Generally good for gas phase properties.
Poor for liquid phase properties.
Better when used in conjunction with a correlation for liquid phase behavior.
Satisfactory for gas phase fugacity calculation @ P_r < T_r/3.
Satisfactory for enthalpy departure and entropy departure calculations.

Soave Redlich Kwong & Peng Robinson EOS

Serve similar functions as the Redlich Kwong EOS but require more parameters.
PR obtains better liquid densities than SRK.
Overall, PR does a better job (slightly) for gas and condensate systems than SRK. However, for polar systems, SRK always makes a better prediction, but in the petroleum engineering business we do not usually deal with those.

Zc as a Measure of Goodness of an EOS

Some values for critical compressibility factors are shown below:

Experimental Z_c for Various Substances

CO₂ = 0.2744
CH₄ = 0.2862
C₂H₆ = 0.2793
nC₅ = 0.2693
nC₆ = 0.2659

The values of critical compressibility factors shown here are relatively close to each other, but, in actuality, they are different. They are, in fact, substance-dependent. This is a striking finding if we recall our highly-praised Principle of Corresponding States. Didn’t we say that at the same reduced conditions, all substances “must” have, at the very least, the same compressibility factor, Z? Here is a case where we have different substances at the same corresponding states (P_r = T_r = 1, right at the critical point) but different “Z” values. After all the Corresponding States Principle is not infallible (as was stated by Pitzer). As we recall, he proposed the introduction of a third parameter (acentric factor) into the corresponding state definition to alleviate these kinds of “problems.”

At the very least, we may say that the values of Z_c (compressibility factor at the critical point) of different substances are “close enough” among themselves. That is, they are not “grossly different,” so as to say that the application of the two-parameter corresponding states principle would be outrageous at the critical point. The fact of the matter is, as a consequence of the Corresponding States Principle, all cubic EOS predict a “unique” and “universal” value of Z at the critical point, regardless of the substance. The list below tells us how they perform.

“Universal” Z_c predicted by different EOS

Ideal EOS = 1.000
vdW EOS = 0.375
RK EOS = 0.333
SRK EOS = 0.333
PR EOS = 0.301

From the list above, we would have been expecting some kind of “average” Z_c of 0.27 or so. But none of the equations of state we have studied is capable of predicting a value that low. The “best” job is done by PR EOS, which provides the “closest” match to the real values observed for most substances. This illustrates why the PR EOS performs somewhat better near critical conditions.

Advantages of Using Cubic Equations of State

All cubic equations of state have their foundation in vdW EOS. The use of cubic equations of state has become widespread because of their advantages:

Simplicity of application
Only a few parameters need to be determined
Low computational overhead is required to implement them. This was a critical issue, particularly in the early stages of computers; it is not really anymore. Nevertheless, this feature is still a “plus.”

The engineer using cubic equations of state must also be aware of the disadvantages that they all share. The most important ones are the limited accuracy that they can provide, particularly for complex systems. In these cases, the procurement of empirical adjustments through the use of the binary interaction parameters (k_ij) is essential.

Solution Techniques for Cubic Expressions and Root Finding

Hopefully, we have convinced you that the use of cubic equations of state can represent a very meaningful and advantageous way of modeling the PVT behavior of petroleum fluids. What we need now are the tools that will allow us to get the information that we want out of them. Even though cubic equations of state are explicit in pressure, pressure is not the common unknown to be calculated in the typical problem. In the most common problem, pressure and temperature are known and we want either molar volume (or its reciprocal, molar density) or compressibility factor (the most likely case). Therefore, we are faced very often with the need to solve for the roots of a cubic expression. Here we present a number of approaches that may be followed.

Analytical scheme

Given the cubic polynomial with real coefficients: x³ + ax² + bx + c = 0, the first step is to calculate the parameters:

$Q = \frac{a^{2} - 3 b}{9} and R= \frac{2 a^{3} - 9 a b + 27 c}{54}$ (11.5)

Now let M = R² – Q³ be the discriminant. We then consider the following cases:

If M < 0 (R² < Q³), the polynomial has three real roots. For this case, compute $θ = \arccos (\frac{R}{\sqrt{Q^{3}}})$ and calculate the three distinct real roots as:

$x_{1} = - (2 \sqrt{Q} \cos \frac{θ}{3}) - \frac{a}{3}$ (11.6a)

$x_{2} = - (2 \sqrt{Q} \cos \frac{θ + 2 π}{3}) - \frac{a}{3}$ (11.6b)

$x_{3} = - (2 \sqrt{Q} \cos \frac{θ - 2 π}{3}) - \frac{a}{3}$ (11.6c)

Note that x₁, x₂, x₃ are not given in any special order, and that $θ$ has to be calculated in radians.
If M > 0 (R² > Q³), the polynomial has only one real root. Compute:

$S = \sqrt[3]{- R + \sqrt{M}}$ (11.7a)

$T = \sqrt[3]{- R - \sqrt{M}}$ (11.7b)
and calculate the real root as follows:
$x_{1} = S + T - \frac{a}{3}$ (11.7c)

Two complex roots (complex conjugates) may be found as well. However, they are of no interest for our purposes, and thus no formulas are provided. Such formulas can be found in the following suggested readings:

W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flannery, Numerical Recipes in Fortran, 2^nd Edition, Cambridge Univ. Press, (1992), p. 179.
Spiegel, M., Liu, J., Mathematical handbook of formulas and tables, 2^nd Edition, Schaum’s Outline Series, McGraw Hill, p.10.

Sometimes, the equations for S and T listed above cause problems while programming. This usually happens whenever the computer/calculator performs the cubic root of a negative quantity. If you want to avoid such a situation, you may compute S’ and T’ instead:

$S' = - s i g n (R) \sqrt[3]{a b s (R) + \sqrt{M}}$ $T' = Q / S'$ (making T’=0 when S’=0)

where:

abs(R) = Absolute value of R

and

sign(R) = (+1) or (– 1) if R is positive or negative respectively.

It may be defined as:

sign(R) = R/Abs(R)

and then the real root is:

$x_{1} = S' + T' - \frac{a}{3}$

Keep in mind the following useful relationships among the roots of any cubic expression:

$x_{1} + x_{2} + x_{3} = - a$ (11.8a)

$x_{1} x_{2} + x_{2} x_{3} + x_{3} x_{1} = + b$ (11.8b)

$x_{1} x_{2} x_{3} = - c$ (11.8c)

Self-Check

Test your understanding of cubic root calculations by analytical means by solving the following examples.

With general cubic expressions,

\begin{array}{l} x^{3} - 6 x^{2} + 11 x - 6 = 0 \Rightarrow \\ x_{1} = 1 \\ x_{2} = 2 \\ x_{3} = 3 \end{array}

\begin{array}{l} x^{3} + 7 x^{2} + 49 x + 343 = 0 \\ x_{1} = - 7 \\ x_{2} = 0 + 7 i \\ x_{3} = 0 - 7 i \end{array}

\begin{array}{l} x^{3} + 2 x^{2} + 3 x + 4 = 0 \Rightarrow \\ x_{1} = - 1.65063 \\ x_{2} = - 0.174685 + 1.54687 i \\ x_{3} = - 0.174685 - 1.54687 i \end{array}

with EOS in terms of volume (v), for a pure component,

$v^{3} - 7.8693 v^{2} + 13.3771 v - 6.5354 = 0$

=> One real root (one phase),

$v_{1} = 5.7357$

Another example,

$v^{3} - 15.6368 v^{2} + 30.315 v - 14.8104 = 0$

=> Two possible phases (three real roots),

$\begin{array}{l} v_{1} = 0.807582 (liquid phase) \\ v_{2} = 1.36174 (rejected) \\ v_{3} = 13.4675 (gas phase) \end{array}$

with EOS in terms of the compressibility factor (z), for a pure component,

$z^{3} - 1.0595 z^{2} + 0.2215 z - 0.01317 = 0$

=> One phase,

$z = 0.8045$

Another example,

$z^{3} - z^{2} + 0.089 z - 0.0013 = 0$

=> Two possible phases,

$\begin{array}{l} z_{l i q u i d} = 0.0183012 \\ z_{x} = 0.0786609 (useless) \\ z_{vapor} = 0.903038 \end{array}$

Numerical Scheme

The Newton-Raphson method provides a useful scheme for solving for a non-explicit variable from any form of equation (not only cubic ones). Newton-Raphson is an iterative procedure with a fast convergence, although it is not always capable of providing an answer — because a first guess close enough to the actual answer must be provided.

In solving for “x” in any equation of the type f(x)=0, the method provides a new estimate (new guess) closer to the actual answer, based on the previous estimate (old guess). This is written as follows:

$x_{n e w} = x_{o l d} - \frac{f (x_{o l d})}{f' (x_{o l d})}$ (11.9)

Considering a cubic equation $f (x) = x^{3} + a x^{2} + b x + c = 0$ , the previous equation takes the form:

$x_{n e w} = x_{o l d} - \frac{x_{o l d}^{3} + a x_{o l d}^{2} + b x_{o l d} + c}{3 x_{o l d}^{2} + 2 a x_{o l d} + b}$ (11.10)

The iterations continue until no significant improvement for “x_new” is achieved, i.e., | x_new – x_old | < tolerance. An educated guess must be provided as the starting value for the iterations. If you are solving a cubic equation in Z (compressibility factor), it is usually recommend to take Z = bP/RT as the starting guess for the compressibility of the liquid phase and Z = 1 for the vapor root.

Semi-analytical Scheme

If you used the previous numerical approach to calculate one of the roots of the cubic expression, the semi-analytical scheme can give you the other two real roots (when they exist). By using the relationships given before, with the value ‘x₁’ as the root already known, the other two roots are calculated by solving the system of equations:

$x_{2} + x_{3} = - a - x_{1}$ (11.11a)

$x_{2} x_{3} = - c / x_{1}$ (11.11b)

which leads to a quadratic expression.

This procedure can be reduced to the following steps:

Let $x^{3} + a x^{2} + b x + c = 0$ be the original cubic polynomial and “E” the root which is already known (x₁ = E). Then, we may factorize such a cubic expression as:

$(x - E) (x^{2} + F x + G) = 0$ (11.12a)
where:

$F = a + E$ (11.12b)

$G = - c / E$ (11.12c)
Solve for x₂, x₃ by using the quadratic expression formulae,

$x_{1} = E$ (11.13a)

$x_{2} = \frac{- F + \sqrt{F^{2} - 4 G}}{2}$ (11.13b)

$x_{3} = \frac{- F - \sqrt{F^{2} - 4 G}}{2}$ (11.13c)

Action Item

Answer the following problems, and submit your answers to the drop box in Canvas that has been created for this module.

Please note:

Your answers must be submitted in the form of a Microsoft Word document.
Include your Penn State Access Account user ID in the name of your file (for example, "module2_abc123.doc").
The due date for this assignment will be sent to the class by e-mail in Canvas.
Your grade for the assignment will appear in the drop box approximately one week after the due date.
You can access the drop box for this module in Canvas by clicking on the Lessons tab, and then locating the drop box on the list that appears.

Problem Set

What is the meaning of Z_c? (critical compressibility factor). Why do you think it is used as a measure of the performance of an EOS?
What do you think should be the value of “Z” at the critical point? Is it “one” (1)? Why? Is this to be expected?

Cubic EOS and Their Behavior (III)

Module Goals

Peng-Robinson EOS (1976)

Comparative Assessment of RK, SRK, and PR EOS

Zc as a Measure of Goodness of an EOS

Some values for critical compressibility factors are shown below:

Experimental Z_c for Various Substances

“Universal” Z_c predicted by different EOS

Advantages of Using Cubic Equations of State

Solution Techniques for Cubic Expressions and Root Finding

Analytical scheme

Self-Check

Numerical Scheme

Action Item

Problem Set

Cubic EOS and Their Behavior (III)

Module Goals

Peng-Robinson EOS (1976)

Comparative Assessment of RK, SRK, and PR EOS

Zc as a Measure of Goodness of an EOS

Some values for critical compressibility factors are shown below:

Experimental Zc for Various Substances

“Universal” Zc predicted by different EOS

Advantages of Using Cubic Equations of State

Solution Techniques for Cubic Expressions and Root Finding

Analytical scheme

Self-Check

Numerical Scheme

Action Item

Problem Set

Experimental Z_c for Various Substances

“Universal” Z_c predicted by different EOS