Published on PNG 520: Phase Behavior of Natural Gas and Condensate Fluids (https://www.e-education.psu.edu/png520)

Home > Course Outline > Cubic EOS and Their Behavior (III)

Cubic EOS and Their Behavior (III)

Module Goals

Module Goal: To demonstrate thermodynamic quantification using modern cubic EOS.

Module Objective: To assess the relative merit of applying the most common EOS.

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+ αa v ˜ 2 +2b v ˜ − b 2 )( v ˜ −b )=RT This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(11.1a)
 

or, explicitly in pressure,

P= RT v ˜ −b − αa v ˜ 2 +2b v ˜ − b 2 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(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− T r ) ⌋ 2 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(11.1c)
 
a=0.45724 R 2 T c 2 P c This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(11.2a)
 

The PR cubic expression in Z becomes:

b=0.07780 R T c P c This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(11.2b)
 

where:

Z 3 −( 1−B ) Z 2 +( A−2B−3 B 2 )Z−( AB− B 2 − B 3 )=0 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(11.3a)
 
A= αaP R 2 T 2 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(11.3b)
 
B= bP RT This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(11.3c)
 

Similar to SRK, the PR mixing rules are:

( αa ) m = ∑ ​ ∑ ​ y i y j ( αa ) ij  ;  ( αa ) ij = ( αa ) i ( αa ) j ( 1− k ij ) This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(11.4a)
 
b m = ∑ i y i b i This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
(11.4b)
 

where binary interaction parameters (kij) again play the important empirical role of helping to better fit experimental data. Due to the empirical character of these interaction parameters, kij’s calculated for PR EOS are unlikely to be the same as the kij’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 @ Pr < Tr/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 Zc for Various Substances

CO2 = 0.2744
CH4 = 0.2862
C2H6 = 0.2793
nC5 = 0.2693
nC6 = 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 (Pr = Tr = 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 Zc (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” Zc 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” Zc 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 (kij) 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: x3 + ax2 + bx + c = 0, the first step is to calculate the parameters:

Q= a 2 −3b 9  and R= 2 a 3 −9ab+27c 54 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.5)

Now let M = R2 – Q3 be the discriminant. We then consider the following cases:

  1. If M < 0 (R2 < Q3), the polynomial has three real roots. For this case, compute θ=arccos( R Q 3 ) This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. and calculate the three distinct real roots as:

    x 1 =−( 2 Q cos θ 3 )− a 3 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.6a)

    x 2 =−( 2 Q cos θ+2π 3 )− a 3 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.6b)

    x 3 =−( 2 Q cos θ−2π 3 )− a 3 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.6c)

    Note that x1, x2, x3 are not given in any special order, and that θ This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. has to be calculated in radians.

  2. If M > 0 (R2 > Q3), the polynomial has only one real root. Compute:

    S= −R+ M 3 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.7a)

    T= −R− M 3 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.7b)

    and calculate the real root as follows:

    x 1 =S+T− a 3 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (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, 2nd Edition, Cambridge Univ. Press, (1992), p. 179.
  • Spiegel, M., Liu, J., Mathematical handbook of formulas and tables, 2nd 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'=−sign( R ) abs( R )+ M 3 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. T'=Q/S' This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (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'− a 3 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.

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

x 1 + x 2 + x 3 =−a This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.8a)

x 1 x 2 + x 2 x 3 + x 3 x 1 =+b This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.8b)

x 1 x 2 x 3 =−c This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.8c)

Self-Check

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

With general cubic expressions,

x 3 −6 x 2 +11x−6=0⇒ x 1 =1 x 2 =2 x 3 =3 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
x 3 +7 x 2 +49x+343=0 x 1 =−7 x 2 =0+7i x 3 =0−7i This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.
x 3 +2 x 2 +3x+4=0⇒ x 1 =−1.65063 x 2 =−0.174685+1.54687i x 3 =−0.174685−1.54687i This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.

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

v 3 −7.8693 v 2 +13.3771v−6.5354=0 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.

=> One real root (one phase),

v 1 =5.7357

Another example,

v 3 −15.6368 v 2 +30.315v−14.8104=0

=> Two possible phases (three real roots),

v 1 = 0.807582 ( liquid phase ) v 2 = 1.36174  (rejected) v 3 = 13.4675( gas phase )

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

z 3 −1.0595 z 2 +0.2215z−0.01317=0 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.

=> One phase,

z=0.8045

Another example,

z 3 − z 2 +0.089z−0.0013=0

=> Two possible phases,

z liquid =0.0183012 z x =0.0786609 (useless) z vapor =0.903038 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers.

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 new = x old − f( x old ) f'( x old ) This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.9)

Considering a cubic equation f( x )= x 3 +a x 2 +bx+c=0 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. , the previous equation takes the form:

x new = x old − x old 3 +a x old 2 +b x old +c 3 x old 2 +2a x old +b This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.10)

The iterations continue until no significant improvement for “xnew” is achieved, i.e., | xnew – xold | < 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 ‘x1’ as the root already known, the other two roots are calculated by solving the system of equations:

x 2 + x 3 =−a− x 1 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.11a)

x 2 x 3 =−c/ x 1 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.11b)

which leads to a quadratic expression.

This procedure can be reduced to the following steps:

  1. Let x 3 +a x 2 +bx+c=0 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. be the original cubic polynomial and “E” the root which is already known (x1 = E). Then, we may factorize such a cubic expression as:

    ( x−E )( x 2 +Fx+G )=0 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.12a)

    where:
     

    F=a+E This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.12b)

    G=−c/E This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.12c)

  2. Solve for x2, x3 by using the quadratic expression formulae,

    x 1 =E This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.13a)

    x 2 = −F+ F 2 −4G 2 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (11.13b)

    x 3 = −F− F 2 −4G 2 This equation is not rendering properly due to an incompatible browser. See Technical Requirements in the Orientation for a list of compatible browsers. (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

  1. What is the meaning of Zc? (critical compressibility factor). Why do you think it is used as a measure of the performance of an EOS?

  2. What do you think should be the value of “Z” at the critical point? Is it “one” (1)? Why? Is this to be expected?


Source URL: https://www.e-education.psu.edu/png520/m11.html