Equations of State and Phase Equilibria

Chapter

5

Equations of State and Phase Equilibria A phase is the part of a system that is uniform in physical and chemical properties, homogeneous in composition, and separated from other, coexisting phases by definite boundary surfaces. The most important phases occurring in petroleum production are the hydrocarbon liquid phase and gas phase. Water is also commonly present as an additional liquid phase. These can coexist in equilibrium when the variables describing change in the entire system remain constant in time and position. The chief variables that determine the state of equilibrium are system temperature, system pressure, and composition. The conditions under which these different phases can exist are a matter of considerable practical importance in designing surface separation facilities and developing compositional models. These types of calculations are based on the following two concepts, equilibrium ratios and flash calculations, which are discussed next.

EQUILIBRIUM RATIOS As indicated in Chapter 1, a system that contains only one component is considered the simplest type of hydrocarbon system. The qualitative understanding of the relationship that exists between temperature, T, pressure, p, and volume, V, of pure components can provide an excellent basis for understanding the phase behavior of complex hydrocarbon mixtures. In a multicomponent system, the partitioning and the tendency of components to escape between the liquid and gas phase are described by the equilibrium ratio “Ki” of a given component. The equilibrium ratio is defined as the ratio of the mole fraction of the component in the gas phase “yi” to the mole fraction of the component in the liquid phase “xi”. Mathematically, the relationship is expressed as Ki ¼

yi xi

Equations of State and PVT Analysis. http://dx.doi.org/10.1016/B978-0-12-801570-4.00005-2 Copyright # 2016 Elsevier Inc. All rights reserved.

(5.1)

467

468 CHAPTER 5 Equations of State and Phase Equilibria

where Ki ¼ equilibrium ratio of component i yi ¼ mole fraction of component i in the gas phase xi ¼ mole fraction of component i in the liquid phase The above expression indicates that when the Ki-value for a given component is greater than 1, it signifies that the component tends to concentrate in the gas phase. Phase equilibria and other phase behavior calculations require accurate estimation equilibrium ratio values. At pressures below 100 psia, Raoult’s and Dalton’s laws for ideal solutions provide a simplified means of predicting equilibrium ratios. Raoult’s law states that the partial pressure, pi, of a component in a multicomponent system is the product of its mole fraction in the liquid phase, xi, and the vapor pressure of the component, pvi: pi ¼ xi pvi

(5.2)

where pi ¼ partial pressure of a component i, psia pvi ¼ vapor pressure of component i, psia xi ¼ mole fraction of component i in the liquid phase Dalton’s law states that the partial pressure of a component is the product of its mole fraction in the gas phase “yi” and system pressure “p”; that is, p i ¼ yi p

(5.3)

where p ¼ system pressure, psia. At equilibrium and in accordance with the previously cited laws, the partial pressure exerted by a component in the gas phase must be equal to the partial pressure exerted by the same component in the liquid phase. Therefore, equating the equations describing the two laws, gives xi pvi ¼ yi p

By rearranging the preceding relationship and introducing the concept of the equilibrium ratio gives yi pvi ¼ ¼ Ki xi p

(5.4)

Eq. (5.4) suggests that for an “ideal solution” with an overall composition of “zi” that exist at low pressure and temperature exhibits the following two properties: n

n

The value of equilibrium ratio “Ki” for any component is independent of the overall composition. As shown in Fig. 5.1, the vapor pressure “pvi” is only a function of temperature, which implies that the K-value is only a function of system pressure “p” and temperature “T”.

Equilibrium Ratios 469

n FIGURE 5.1 Vapor pressure chart for hydrocarbon components. Courtesy of the Gas Processors Suppliers Association. Published in the GPSA Engineering Data Book, Tenth Edition, 1987.

470 CHAPTER 5 Equations of State and Phase Equilibria

Taking the logarithm of both sides of Eq. (5.4) and rearranging gives log ðKi Þ ¼ log ðpvi Þ  log ðpÞ

This relationship indicates that for a constant temperature, T (implying a constant pvi), the component equilibrium ratio “Ki” is a linear function of pressure with a slope of 1 when plotted on a log-log scale. It is appropriate at this stage to introduce and define the following nomenclatures: zi ¼ mole fraction of component i in the entire hydrocarbon mixture nt ¼ total number of moles of the hydrocarbon mixture, lb-mol nL ¼ total number of moles in the liquid phase nυ ¼ total number of moles in the vapor (gas) phase By definition, nt ¼ nL + nυ

(5.5)

Eq. (5.5) indicates that the total number of moles in the system is equal to the total number of moles in the liquid phase plus the total number of moles in the vapor phase. A material balance on the ith component results in zi nt ¼ xi nL + yi nv

(5.6)

where zint ¼ total number of moles of component i in the system xinL ¼ total number of moles of component i in the liquid phase yinv ¼ total number of moles of component i in the vapor phase Also, by the definition of the total mole fraction in a hydrocarbon system, we may write X zi ¼ 1

(5.7)

i

X xi ¼ 1

(5.8)

i

X yi ¼ 1

(5.9)

i

It is convenient to perform all the phase equilibria calculations on the basis of 1 mol of the hydrocarbon mixture; that is, nt ¼ 1. That assumption reduces Eqs. (5.5), (5.6) to nL + nυ ¼ 1

(5.10)

xi nL + yi nυ ¼ zi

(5.11)

Combining Eqs. (5.4), (5.11) to eliminate yi from Eq. (5.11) gives xi nL + ðxi Ki Þnυ ¼ z

Flash Calculations 471

Solving for xi yields zi nL + nυ Ki

xi ¼

(5.12)

Eq. (5.11) can also be solved for yi by combining it with Eq. (5.4) to eliminate xi, which gives yi ¼

zi Ki ¼ xi Ki nL + nυ Ki

(5.13)

Combining Eq. (5.12) with Eq. (5.8) and Eq. (5.13) with Eq. (5.9) results in X X xi ¼ i

X

i

yi ¼

i

X i

zi ¼1 nL + nυ Ki

(5.14)

zi Ki ¼1 nL + nυ Ki

(5.15)

As X i

yi 

X xi ¼ 0 i

then X i

X zi Ki zi  ¼0 nL + nυ Ki n + L nυ Ki i

Rearranging the above expression to give X zi ðKi  1Þ i

nL + nυ Ki

¼0

Replacing nL in the above relationship with (1  nυ) yields f ðnυ Þ ¼

X zi ðKi  1Þ ¼0 nυ ðKi  1Þ + 1 i

(5.16)

Eqs. (5.12)–(5.16) provide the necessary relationships to perform volumetric and compositional calculations on a hydrocarbon system. These types of calculation are referred to as flash calculations. Details of performing flash calculations are described next.

FLASH CALCULATIONS Flash calculations are an integral part of all reservoir and process engineering calculations. They are required whenever it is desirable to know the amounts (in moles) of hydrocarbon liquid and gas coexisting in a reservoir

472 CHAPTER 5 Equations of State and Phase Equilibria

or a vessel at a given pressure and temperature. These calculations also are performed to determine the composition of the existing hydrocarbon phases. Given the overall composition of a hydrocarbon system at a specified pressure and temperature, flash calculations are performed to determine the moles of the gas phase, nυ, moles of the liquid phase, nL, composition of the liquid phase, xi, and composition of the gas phase, yi. The computational steps for determining nL, nv, yi, and xi of a hydrocarbon mixture with a known overall composition of zi and characterized by a set of equilibrium ratios, Ki, are summarized in the following steps: Step 1. Calculate number of moles of gas phase “nυ”: Eq. (5.16) can be solved for the number of moles of the vapor phase nυ by using the Newton-Raphson iteration technique. In applying this iterative technique, assume any arbitrary value of nυ between 0 and 1, such as nυ ¼ 0.5. A good assumed value may be calculated from the following relationship: nυ ¼ A=ðA + BÞ

where A¼

X ½zi ðKi  1Þ i

 X  1 zi 1 B¼ Ki i

Using the assumed value of nv, evaluate the function f(nυ) as given by Eq. (5.16); that is, f ðnυ Þ ¼

X zi ðKi  1Þ nv ðKi  1Þ + 1 i

If the absolute value of the function f(nυ) is smaller than a preset tolerance, such as 106, then the assumed value of nυ is the desired solution. If the absolute value of f(nυ) is greater than the preset tolerance, then a new value of (nυ)new is calculated from the following expression: ðnυ Þnew ¼ nυ 

f ðnυ Þ f 0 ðnυ Þ

with the derivative f 0 (nυ) as given by f 0 ðnυ Þ ¼ 

( X i

zi ðKi  1Þ2 ½nυ ðKi  1Þ + 12

)

Flash Calculations 473

where (nυ)new is the new value of nυ to be used for the next iteration. This procedure is repeated with the new value of nυ until convergence is achieved; that is, when jf ðnυ Þj  eps or alternatively jðnυ Þnew  ðnυ Þj  eps

in which eps is a selected tolerance, such as eps ¼ 106. Step 2. Calculate number of moles of liquid phase “nL”: The number of moles of the liquid phase “nL” can be determined from Eq. (5.10) to give nL ¼ 1  nυ

Step 3. Calculate the composition of the liquid phase “xi”: Given nv and nL, calculate the composition of the liquid phase by applying Eq. (5.12): xi ¼

zi nL + nυ Ki

Step 4. Calculate the composition of the gas phase “yi”: Determine the composition of the gas phase from Eq. (5.13): yi ¼

zi Ki ¼ xi Ki nL + nυ Ki

EXAMPLE 5.1 A hydrocarbon mixture with the following overall composition is flashed in a separator at 50 psia and 100°F: Component

zi

C3 i-C4 n-C4 i-C5 n-C5 C6

0.20 0.10 0.10 0.20 0.20 0.20

Assuming an ideal solution behavior, perform flash calculations. Solution Step 1. Determine the vapor pressure pvi from the Cox chart (Fig. 5.1) and calculate the equilibrium ratios using Eq. (5.4). The results are shown below.

474 CHAPTER 5 Equations of State and Phase Equilibria

Component

zi

pvi at 100°F

Ki 5 pvi/50

C3 i-C4 n-C4 i-C5 n-C5 C6

0.20 0.10 0.10 0.20 0.20 0.20

190 72.2 51.6 20.44 15.57 4.956

3.80 1.444 1.032 0.4088 0.3114 0.09912

Step 2. Solve Eq. (5.16) for nυ using the Newton-Raphson method: ðnυ Þn ¼ nυ 

f ðnυ Þ f 0 ðnυ Þ

Iteration

nυ

f(nv)

0 1 2 3 4

0.08196579 0.1079687 0.1086363 0.1086368 0.1086368

3.073(102) 8.894(104) 7.60(107) 1.49(108) 0.0

to give nυ ¼ 0.1086368. Step 3. Solve for nL: nL ¼ 1  nυ nL ¼ 1  0:1086368 ¼ 0:8913631 Step 4. Solve for xi and yi to yield xi ¼

zi nL + nυ Ki

yi ¼ xi Ki The component results are shown below. Component

zi

Ki

xi 5 zi/(0.8914 + 0.1086 Ki)

yi 5 xi Ki

C3 i-C4 n-C4 i-C5 n-C5 C6

0.20 0.10 0.10 0.20 0.20 0.20

3.80 1.444 1.032 0.4088 0.3114 0.09912

0.1534 0.0954 0.0997 0.2137 0.2162 0.2216

0.5829 0.1378 0.1029 0.0874 0.0673 0.0220

Flash Calculations 475

Binary System Note that, for a binary system (ie, a two-component system), flash calculations can be performed without restoring to the preceding iterative technique. Flash calculations can be performed by applying the following steps. Step 1. Solve for the composition of the liquid phase, xi. For a twocomponent system, Eqs. (5.8), (5.9) can be expanded as X xi ¼ x1 + x2 ¼ 1 i

X yi ¼ y1 + y2 ¼ K1 x1 + K2 x2 ¼ 1 i

Solving these expressions for the liquid composition, x1 and x2, gives x1 ¼

1  K2 K1  K2

and x2 ¼ 1  x1

where x1 ¼ mole fraction of the first component in the liquid phase x2 ¼ mole fraction of the second component in the liquid phase K1 ¼ equilibrium ratio of the first component K2 ¼ equilibrium ratio of the second component Step 2. Solve for the composition of the gas phase, yi. From the definition of the equilibrium ratio, calculate the composition of the liquid as follows: y1 ¼ x1 K1 y2 ¼ x2 K2 ¼ 1  y1

Step 3. Solve for the number of moles of the vapor phase, nυ, and liquid phase, nL. Arrange Eq. (5.12) to solve for nυ by using the mole fraction and K-value of one of the two components to give nυ ¼

z1  x 1 x1 ðK1  1Þ

and nL ¼ 1  nυ

Exact results will be obtained if selecting the second component; that is, nυ ¼

z2  x 2 x2 ðK2  1Þ

476 CHAPTER 5 Equations of State and Phase Equilibria

and nL ¼ 1  nυ

where z1 ¼ mole fraction of the first component in the binary system x1 ¼ mole fraction of the first component in the liquid phase z2 ¼ mole fraction of the second component in the binary system x2 ¼ mole fraction of the second component in the liquid phase K1 ¼ equilibrium ratio of the first component K2 ¼ equilibrium ratio of the second component The equilibrium ratio, as calculated by Eq. (5.4) in terms of vapor pressure and system pressure, proved inadequate. Eq. (5.4) was developed based on the following assumptions: n n

The vapor phase is an ideal gas as described by Dalton’s law. The liquid phase is an ideal solution as described by Raoult’s law.

The combination of the above two assumptions is unrealistic and results in inaccurate predictions of equilibrium ratios at high pressures.

EQUILIBRIUM RATIOS FOR REAL SOLUTIONS For a real solution, the equilibrium ratios are no longer a function of the pressure and temperature alone but also the composition of the hydrocarbon mixture. This observation can be stated mathematically as Ki ¼ K ðp, T, zi Þ

Numerous methods have been proposed for predicting the equilibrium ratios of hydrocarbon mixtures. These correlations range from a simple mathematical expression to a complicated expression containing several compositional dependent variables. The following methods are presented: Wilson’s correlation, Standing’s correlation, the convergence pressure method, and Whitson and Torp’s correlation.

Wilson’s Correlation Wilson (1968) proposed a simplified thermodynamic expression for estimating K-values. The proposed expression has the following form: Ki ¼

   pci Tci exp 5:37ð1 + ωi Þ 1  p T

(5.17)

Equilibrium Ratios for Real Solutions 477

where pci ¼ critical pressure of component i, psia p ¼ system pressure, psia Tci ¼ critical temperature of component i, °R T ¼ system temperature, °R This relationship generates reasonable values for the equilibrium ratio when applied at low pressures. The main advantage of using Wilson’s expression to estimate the equilibrium ratios is in its ability to provide a consistent set of K-values at a specified pressure and temperature; that is, KC1 > KC2 > KC3 >, and so on. As detailed later in this chapter, equations of state rely on iterative technique to calculate the K-value with stating values based on Wilson’s equation.

Standing’s Correlation Several authors, including Hoffmann et al. (1953), Brinkman and Sicking (1960), Kehn (1964), and Dykstra and Mueller (1965), suggested that any hydrocarbon or nonhydrocarbon component could be uniquely characterized by combining its boiling point temperature, critical temperature, and critical pressure into a characterization parameter, which is defined by the following expression: Fi ¼ bi ½1=Tbi  1=T 

(5.18)

log ðpci =14:7Þ ½1=Tbi  1=Tci 

(5.19)

with bi ¼

where Fi ¼ component characterization factor Tbi ¼ normal boiling point of component i, °R. Based on the characterization parameter “Fi,” Standing (1979) derived a set of equations that fit the graphical presentations of the experimentally determined equilibrium ratio of Katz and Hachmuth (1937). The proposed form of the correlation is based on an observation that plots of log(Ki p) versus Fi at a given pressure often form straight lines with a slope of “c” and intercept of “a”. The basic equation of the straight-line relationship is given by log ðKi pÞ ¼ a + cFi

Solving for the equilibrium ratio, Ki, gives

478 CHAPTER 5 Equations of State and Phase Equilibria

1 Ki ¼ 10ða + cFi Þ p

(5.20)

where the coefficients “a” and “c” in the relationship are the intercept and the slope of the line, respectively. From six isobar plots of log(Ki p) versus Fi for 18 sets of equilibrium ratio values, Standing correlated the coefficients “a” and “c” with the pressure, to give   a ¼ 1:2 + 0:00045p + 15 108 p2   c ¼ 0:89  0:00017p  3:5 108 p2

(5.21) (5.22)

Standing pointed out that the predicted values of the equilibrium ratios of N2, CO2, H2S, and C1–C6 can be improved considerably by changing n n

the correlating parameter “bi” of Eq. (5.19) the boiling point “Tbi” of these components.

The author proposed the following optimized modified values: Component

bi

Tbi (°R)

N2 CO2 H2S C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 n-C6 n-C7 n-C8 n-C9 n-C10

470 652 1136 300 1145 1799 2037 2153 2368 2480 2738 2780 3068 3335 3590 3828

109 194 331 94 303 416 471 491 542 557 610 616 616 718 763 805

K-Value of the Plus Fraction When making flash calculations, the equilibrium ratio for each component including the plus fraction (eg, C7+) must be calculated at the prevailing pressure and temperature. Katz and Hachmuth (1937) proposed a

Equilibrium Ratios for Real Solutions 479

“rule-of-thumb” that might be applied for light hydrocarbon mixtures that suggest the K-value of C7+ can be taken as 15% of the K of C7; that is, KC7 + ¼ 0:15KC7

Standing offered an alternative approach for determining the K-value of the heptanes and heavier fractions. By imposing experimental equilibrium ratio values for C7+ on Eq. (5.20), Standing calculated the corresponding characterization factors, Fi, for the plus fraction; that is, 

log FC7 + ¼

KC7 + p c



a

The calculated FC7 + values were used to specify or FIND the pure normal paraffin hydrocarbon having the K-value equivalent to that of the C7+ fraction. Standing suggested the following computational steps for determining the parameters b and Tb of the heptanes-plus fraction. Step 1. Determine from the following relationship the number of carbon atoms “n” of the normal paraffin hydrocarbon that have the K-value of the C7+ fraction: n ¼ 7:30 + 0:0075ðT  460Þ + 0:0016p

(5.23)

Step 2. Calculate the correlating parameter “b” and the boiling point “Tb” of that normal paraffin component from the following expression: b ¼ 1013 + 324n  4:256n2

(5.24)

Tb ¼ 301 + 59:85n  0:971n2

(5.25)

Values of “b” and “Tb” can then be used to calculate the characterization parameter “Fi” of that particular normal paraffin hydrocarbon and be used to calculate its K-value, which is assigned to the C7+. It is interesting to note that numerous experimental phase equilibria data suggest that the equilibrium ratio for carbon dioxide can be closely approximated by the following relationship: KCO2 ¼

pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi KC1 KC2

where KCO2 ¼ equilibrium ratio of CO2 at system pressure p and temperature T KC1 ¼ equilibrium ratio of methane at p and T KC2 ¼ equilibrium ratio of ethane at p and T

480 CHAPTER 5 Equations of State and Phase Equilibria

It should be pointed out that Standing’s proposed correlations are valid for pressures and temperatures that are less than 1000 psia and 200°F, respectively. This limitation on the use of Standing method lends itself appropriately for performing separation calculations, which usually operate within the range of Standing’s limit of applications.

A Simplified K-Value Method For light hydrocarbon mixtures that exist at pressure below 1000 psia, the K-value for the methane fraction can be estimated from the following relationship:   B exp A  T KC1 ¼ p

with

  A ¼ 2:0 107 p2  0:0005p + 9:5073 B ¼ 0:0001p2  0:456p + 856

where p ¼ pressure, psia T ¼ temperature, °R. It is possible to correlate the equilibrium ratios of C2–C7 with that of the Kvalue of C1 by the following relationship: Ki ¼

KC1 Ri lnðpKC1 Þ

with the component characterization parameter Ri for component “i” as given by Ri ¼ ai T  bi

The temperature, T, is in °R. Values of the coefficients ai and bi for C2–C7 are tabulated below: Component

ai

bi

C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7

0.00665 0.00431 0.00327 0.00245 0.00203 0.00147 0.00094 0.00091

1.29114 1.269876 1.206566 1.1002 0.910736 0.737 0.528 0.516263

Equilibrium Ratios for Real Solutions 481

This correlation provides a very rough estimate of the K-values for the above listed components at pressures below 1000 psia. The K-value for the C7+ can be estimated from the “rule-of-thumb” equation (ie, KC7 + ¼ 0:15KC7 ).

EXAMPLE 5.2 A hydrocarbon mixture with the following composition is flashed at 1000 psia and 150°F. Component

zi

CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.009 0.003 0.535 0.115 0.088 0.023 0.023 0.015 0.015 0.015 0.159

If the molecular weight and specific gravity of C7+ are 150.0 and 0.78, respectively, calculate the equilibrium ratios using: n n n

Wilson’s correlation Standing’s correlation Simplified method

Solution Using Wilson’s method Step 1. Calculate the critical pressure, critical temperature, and acentric factor of C7+ by using the characterization method of Riazi and Daubert discussed in Chapter 2 and Example 2.1, to give Tc ¼ 1139:4°R pc ¼ 320:3 psia ω ¼ 0:5067 Step 2. Apply Eq. (5.17) to get the results shown below: Ki 5

   pci Tci exp 5:37ð1 + ωi Þ 1 1000 610

482 CHAPTER 5 Equations of State and Phase Equilibria

Component

P0 c (psia)

T0 c (°R)

ω

Ki

CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

1071 493 667.8 707.8 616.3 529.1 550.7 490.4 488.6 436.9 320.3

547.9 227.6 343.37 550.09 666.01 734.98 765.65 829.1 845.7 913.7 1139.4

0.225 0.040 0.0104 0.0986 0.1524 0.1848 0.2010 0.2223 0.2539 0.3007 0.5069

2.0923 16.343 7.155 1.263 0.349 0.144 0.106 0.046 0.036 0.013 0.00029

Using Standing’s method Step 1. Calculate the coefficients a and c from Eqs. (5.21), (5.22) to give   a ¼ 1:2 + 0:00045p + 15 108 p2   a ¼ 1:2 + 0:00045 ð1000Þ + 15 108 ð1000Þ2 ¼ 1:80   c ¼ 0:89  0:00017p  3:5 108 p2   c ¼ 0:89  0:00017 ð1000Þ  3:5 108 ð1000Þ2 ¼ 0:685 Step 2. Calculate the number of carbon atoms, n, from Eq. (5.23) to give n ¼ 7:30 + 0:0075 ðT  460Þ + 0:0016p n ¼ 7:3 + 0:0075 ð150Þ + 0:0016 ð1000Þ ¼ 10:025 Step 3. Determine the parameter b and the boiling point, Tb, for the hydrocarbon component with n carbon atoms by using Eqs. (5.24), (5.25): b ¼ 1013 + 324n  4:256n2 b ¼ 1013 + 324ð10:025Þ  4:256ð10:025Þ2 ¼ 3833:369 Tb ¼ 301 + 59:85n  0:971n2 Tb ¼ 301 + 59:85ð10:025Þ  0:971ð10:025Þ2 ¼ 803:41°R Step 4. Apply Eqs. (5.18), (5.20) to give the results shown below: Fi ¼ bi ½1=Tbi  1=T  1 Ki ¼ 10ða + cFi Þ p

Equilibrium Ratios for Real Solutions 483

Component

bi

Tbi

Fi

Ki

CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

652 470 300 1145 1799 2037 2153 2368 2480 2738 3833.369

194 109 94 303 416 471 491 542 557 610 803.41

2.292 3.541 2.700 1.902 1.375 0.985 0.855 0.487 0.387 0 1.513

2.344 16.811 4.462 1.267 0.552 0.298 0.243 0.136 0.116 0.063 0.0058

The Convergence Pressure Method Early high-pressure phase equilibria studies revealed that, when a hydrocarbon mixture of a fixed overall composition is held at a constant temperature as the pressure increases, the equilibrium values of all components converge toward a common value of unity at certain pressure. This pressure is termed the convergence pressure, pk, of the hydrocarbon mixture. The convergence pressure essentially is used to correlate the effect of the composition on equilibrium ratios. The concept of the convergence pressure can be better explained by examining Fig. 5.2. The figure is a schematic diagram of a typical set of equilibrium ratios plotted versus pressure on log-log paper for a hydrocarbon mixture held at a constant temperature. The illustration shows a tendency of the equilibrium ratios to converge isothermally to a value of Ki ¼ 1 for all components at a specific pressure; that is, convergence pressure. A different hydrocarbon mixture may exhibit a different convergence pressure. Standing (1977) suggested that the convergence pressure could be roughly correlated linearly with the molecular weight of the heptanes-plus fraction. Whitson and Torp (1981) expressed this relationship by the following equation: pk ¼ 60MC7 +  4200

where MC7 + is the molecular weight of the heptanes-plus fraction.

(5.26)

484 CHAPTER 5 Equations of State and Phase Equilibria

Constant temperature

Me

tha

ne

10 Eth

an

Pro

e

Equilibrium ratio “K-value”

pa

ne

Bu

tan

1.0

es

Pen

tan

es

He

xan

es

0.1 Hep

tan

es

1000 Pressure

100

PK

10,000

n FIGURE 5.2 Equilibrium ratios for a hydrocarbon system.

Whitson and Torp reformulated Wilson’s equation [Eq. (5.17)] to yield more accurate results at higher pressures by incorporating the convergence pressure into the correlation, to give Ki ¼

 A1      pci pci Tci exp 5:37Að1 + ωi Þ 1  pk p T

(5.27)

with the exponent “A” roughly approximated by the following expression: A¼ 1

 0:7 p pk

where p ¼ system pressure, psig pk ¼ convergence pressure, psig T ¼ system temperature, °R ωi ¼ acentric factor of component i

Equilibrium Ratios for the Plus Fractions 485

EXAMPLE 5.3 Rework Example 5.2 and calculate the equilibrium ratios by using Whitson and Torp’s method. Solution Step 1. Determine the convergence pressure from Eq. (5.26) to give pk ¼ 60MC7 +  4200 ¼ 9474 Step 2. Calculate the coefficient A:  0:7 p A¼1 pk  A¼1

 1000 0:7 ¼ 0:793 9474

Step 3. Calculate the equilibrium ratios from Eq. (5.27) to give the results shown below:    p 0:7931 p Tci ci ci Ki ¼ exp 5:37Að1 + ωi Þ 1  9474 1000 610

Component

pc

T0 c (°R)

ω

Ki

CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

1071 493 667.8 707.8 616.3 529.1 550.7 490.4 488.6 436.9 320.3

547.9 227.6 343.37 550.09 666.01 734.98 765.65 829.1 845.7 913.7 1139.4

0.225 0.040 0.0104 0.0986 0.1524 0.1848 0.2010 0.2223 0.2539 0.3007 0.5069

2.9 14.6 7.6 2.1 0.7 0.42 0.332 0.1794 0.150 0.0719 0.683(103)

EQUILIBRIUM RATIOS FOR THE PLUS FRACTIONS The equilibrium ratios of the plus fractions often behave in a manner different from the other components of a system. This is because the plus fraction, in itself, is a mixture of components. Several techniques have been proposed for estimating the K-value of the plus fractions. Three methods are presented next.

486 CHAPTER 5 Equations of State and Phase Equilibria

n n n

Campbell’s method Winn’s method Katz’s method

Campbell’s Method Campbell (1976) proposed that the plot of log(Ki) versus (Tci)2 for each component is a linear relationship for any hydrocarbon system. Campbell suggested that, by drawing the best straight line through the points for propane through hexane components, the resulting line can be extrapolated to obtain the K-value of the plus fraction. Alternatively, a plot of log (Ki) versus (1=Tbi ) of each fraction in the mixture could also produce a straight line that can be extrapolated to obtain the equilibrium ratio of the plus fraction from the reciprocal of its average boiling point.

Winn’s Method Winn (1954) proposed the following expression for determining the equilibrium ratio of heavy fractions with a boiling point above 210°F: KC + ¼

KC7 ðKC2 =KC7 Þb

(5.28)

where KC + ¼value of the plus fraction KC7 ¼ K-value of n-heptane at system pressure, temperature, and convergence pressure KC2 ¼ K-value of ethane b ¼ volatility exponent Winn correlated, graphically, the volatility component, b, of the heavy fraction, with the atmosphere boiling point, as shown in Fig. 5.3. This graphical correlation can be expressed mathematically by the following equation: b ¼ a1 + a2 ðTb  460Þ + a3 ðTb  460Þ2 + a4 ðTb  460Þ3 + a5 =ðT  460Þ (5.29)

where Tb ¼ boiling point, °R a1–a5 ¼ coefficients with the following values: a1 ¼ 1:6744337 a2 ¼ 3:4563079  103 a3 ¼ 6:1764103  106 a4 ¼ 2:4406839  109 a5 ¼ 2:9289623  102

Equilibrium Ratios for the Plus Fractions 487

n FIGURE 5.3 Volatility exponent “b”.

Katz’s Method Katz et al. (1959) suggested that a factor of 0.15 times the equilibrium ratio for the heptane component gives a reasonably close approximation to the equilibrium ratio for heptanes and heavier. This suggestion is expressed mathematically by the following equation: KC7 + ¼ 0:15KC7

(5.30)

K-values for Nonhydrocarbon Components Lohrenze et al. (1963) developed the following set of correlations that express the K-values of H2S, N2, and CO2 as a function of pressure, temperature, and the convergence pressure, pk. For H2 S :

  " p 0:8 1399:2204 6:3992127 + lnðKH2 S Þ ¼ 1  pk T #   18:215052 1112446:2  ln ðpÞ 0:76885112 +  T T2

488 CHAPTER 5 Equations of State and Phase Equilibria

For N2 :       p 0:4 1184:2409 11:294748   0:90459907 ln ðpÞ ln KN 2 ¼ 1  pk T For CO2 :      p 0:6 152:7291 7:0201913   ln ðpÞ ln KCO 2 ¼ 1  pk T  # 1719:2856 644740:69 1:8896974  + T T2

where T ¼ temperature, °R p ¼ pressure, psia pk ¼ convergence pressure, psia

VAPOR-LIQUID EQUILIBRIUM CALCULATIONS The vast amount of experimental and theoretical work performed on equilibrium ratio studies indicates their importance in solving phase equilibrium problems in reservoir and process engineering. The basic phase equilibrium calculations include: n n n

the dew point pressure “Pd” bubble point pressure “Pb” separator calculations.

The use of the equilibrium ratio in performing these applications is discussed next.

Dew Point Pressure The dew point pressure, pd, of a hydrocarbon system is defined as the pressure at which an infinitesimal quantity of liquid is in equilibrium with a large quantity of gas. For a total of 1 lb-mol of a hydrocarbon mixture (ie, nt ¼ 1), the following conditions are applied at the dew point pressure: nL  0 nv  1 yi ¼ zi

Imposing the above constraints on Eq. (5.14), gives X i

X X zi zi zi ¼ ¼ ¼1 nL + nυ Ki 0 + ð1:0ÞKi Ki i i

where zi is the overall composition of the hydrocarbon system.

(5.31)

Vapor-Liquid Equilibrium Calculations 489

The solution of Eq. (5.31) for the dew point pressure, pd, involves a trialand-error process, as summarized in the following steps: Step 1. Assume a trial value of pd. A good starting value can be obtained by combining Wilson’s equation (ie, Eq. (5.17)) with Eq. (5.31) to give X zi i

Ki

¼

8 > < X> i

9 > > =

zi    ¼ 1 > > p T > : ci exp 5:37ð1 + ωi Þ 1  ci > ; pd T

(5.32)

Solving the above expression for pd to give an initial estimate of pd as pd ¼

1

8 > < X> i

9 > > = zi    > > T > :pci exp 5:37ð1 + ωi Þ 1  ci > ; T

(5.33)

Step 2. Using the assumed dew point pressure, calculate the equilibrium ratio, Ki, for each component at the system temperature. Step 3. Calculate the summation of Eq. (5.31) using the K-values of step 2 X

zi =Ki

i

Step 4. The correct value of the dew point pressure is obtained when the sum  1; however, if n the sum ≪ 1, repeat steps 2 and 3 at a higher initial value of pressure n the sum ≫ 1, repeat the calculations with a lower initial value of pd. The above outlined procedure using Wilson’s correlation is intended to only provide approximations of the K-values and Pd. These approximations are used as initial values in equations of state application. However, using Eq. (5.33) will generate a very low value for the dew point pressure; perhaps close to the lower dew point value. The main reason is that when applying Eq. (5.33), the term “zi/[Pc exp(5.37(1 + ω)(1 Tc/T))]” for the plus fraction is much greater than the nearest component value by over 40 times, resulting in a large value of the denominator of Eq. (5.33). Because it’s only an approximation, the plus fraction should be excluded.

490 CHAPTER 5 Equations of State and Phase Equilibria

EXAMPLE 5.4 A natural gas reservoir has the following composition: Component

zi

C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.7778 0.0858 0.0394 0.0083 0.016 0.0064 0.0068 0.0078 0.0517

The heptanes-plus fraction is characterized by a molecular weight and specific gravity of 150 and 0.78, respectively. The reservoir exists at 250°F; estimate the dew point pressure of the gas mixture. Solution Step 1. Calculate the critical properties of the C7+ using the Riazi-Daubert equation: θ ¼ a ðMÞb γ c exp½d M + eγ + f γ M  Tc ¼ 544:2ð150Þ0:2998 ð0:78Þ1:0555 exp 1:3478  104 ð150Þ  0:61641ð0:78Þ + 0 ¼ 1139:4°R  pc ¼ 4:5203  104 ð150Þ0:8063 ð0:78Þ1:6015 exp 1:8078  103 ð150Þ  0:3084ð0:78Þ + 0 ¼ 320:3 psia  Tb ¼ 6:77857ð150Þ0:401673 ð0:78Þ1:58262 exp 3:77409  103 ð150Þ  2:984036ð0:78Þ 4:25288  103 ð150Þð0:78Þ ¼ 825:26°R

Use Edmister’s Eq. (2.21) to estimate the acentric factor: 3½log ðpc =14:70Þ 1 7½Tc =Tb  1 3½log ð320:3=14:7Þ ω¼  1 ¼ 0:5067 7½1139:4=825:26  1 ω¼

Step 2. Construct the following working table and approximate the dew point pressure by applying Eq. (5.33): pd ¼

8 > < X>

1

9 > > =

zi    > T > :pci exp 5:37ð1 + ωi Þ 1  ci > ; T

i>

Vapor-Liquid Equilibrium Calculations 491

Component

zi

Pci

Tci

ωi

Sum

C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.7778 0.0858 0.0394 0.0083 0.016 0.0064 0.0068 0.0078 0.0517

667.8 707.8 616.3 529.1 550.7 490.4 488.6 436.9 320.3

343.37 550.09 666.01 734.98 765.65 829.1 845.7 913.7 1139.4

0.0104 0.0986 0.1524 0.1848 0.201 0.2223 0.2539 0.3007 0.5067

7.0699E-05 3.2101E-05 4.357E-05 1.9623E-05 4.8166E-05 3.9247E-05 5.0404E-05 0.00013244 0.02153184

n n

Total sum of the last column ¼ 0.02197 to give a dew point pressure ¼ 1/0.02197 ¼ 45 psi Excluding the C7+, sum ¼ 0.000436 to give a dew point pressure ¼ 1/0.000436 ¼ 2292 psi

Bubble Point Pressure At the bubble point pb, the hydrocarbon system is essentially liquid, except for an infinitesimal amount of vapor. For a total of 1 lb-mol of the hydrocarbon mixture, the following conditions are applied at the bubble point pressure: nL  1 nv  0 xi ¼ zi

Obviously, under these conditions, xi ¼ zi. Imposing the above constraints to Eq. (5.15) yields X i

X zi K i X zi K i ¼ ¼ ðzi K i Þ ¼ 1 nL + nv Ki 1 + ð0ÞKi i i

(5.34)

Following the procedure outlined for the dew point pressure determination, Eq. (5.34) is solved for the bubble point pressure, pb, by assuming various pressures and determining the pressure that produces K-values satisfying Eq. (5.34). During the iterative process, if X ðzi Ki Þ < 1, then the assumed pressure is high i X ðzi Ki Þ > 1, then the assumed pressure is low i

Wilson’s equation can be used to give an estimate of the bubble point pressure for equations of state application:

492 CHAPTER 5 Equations of State and Phase Equilibria

   X pci Tci ¼1 zi exp 5:37ð1 + ωÞ 1  pb T i

Solving for the bubble point pressure gives pb 

   X Tci zi pci exp 5:37ð1 + ωÞ 1  T i

(5.35)

EXAMPLE 5.5 A crude oil reservoir exits at 200°F with a composition and associated component properties as tabulated below: Component

zi

Pci

Tci

ωi

C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.42 0.05 0.05 0.03 0.02 0.01 0.01 0.01 0.40

667.8 707.8 616.3 529.1 550.7 490.4 488.6 436.9 230.4

343.37 550.09 666.01 734.98 765.65 829.1 845.7 913.7 1279.8

0.0104 0.0986 0.1524 0.1848 0.2010 0.2223 0.2539 0.3007 1.0653

The heptanes-plus fraction is characterized by the following measured properties: ðMÞC7 + ¼ 216:0 ðγ ÞC7 + ¼ 0:8605 ðTb ÞC7 + ¼ 977°R Estimate the bubble point pressure. Solution Step 1. Calculate the critical pressure and temperature of the C7+ by the RiaziDaubert equation, Eq. (2.4), to give pc ¼ 3:12281  109 ð977Þ2:3125 ð0:8605Þ2:3201 ¼ 230:4psi Tc ¼ 24:27870ð977Þ0:58848 ð0:8605Þ0:3596 ¼ 1279:8°R Step 2. Calculate the acentric factor by employing the Edmister correlation [Eq. (2.21)] to yield ω¼

3½ log ðpc =14:70Þ ¼ 1:0653 7½Tc =Tb  1

Vapor-Liquid Equilibrium Calculations 493

Step 3. Construct the following table: Component

zi

Pci

Tci

ωi

C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.42 0.05 0.05 0.03 0.02 0.01 0.01 0.01 0.40

667.8 707.8 616.3 529.1 550.7 490.4 488.6 436.9 230.4

343.37 550.09 666.01 734.98 765.65 829.1 845.7 913.7 1279.8

0.0104 0.0986 0.1524 0.1848 0.2010 0.2223 0.2539 0.3007 1.0653

Step 4. Estimate the bubble point pressure from Eq. (5.35) as follows:    X Tci zi pci exp 5:37ð1 + ωÞ 1  pb  T i

Component

zi

Pci

Tci

ωi

Eq. (5.35)

C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C P7+

0.42 0.05 0.05 0.03 0.02 0.01 0.01 0.01 0.4

667.8 707.8 616.3 529.1 550.7 490.4 488.6 436.9 320.3

343.37 550.09 666.01 734.98 765.65 829.1 845.7 913.7 1279.8

0.0104 0.0986 0.1524 0.1848 0.201 0.2223 0.2539 0.3007 1.0633

4620.64266 133.639009 45.2146052 12.6894565 6.6436649 1.630693 1.34909342 0.58895645 0.01761403 4822.41575

The estimated bubble point pressure is 4822 psia; notice the impact of the methane on the estimated saturation pressure.

Separator Calculations Produced reservoir fluids are complex mixtures of different physical characteristics. As a well stream flows from the high-temperature, high-pressure petroleum reservoir, it experiences pressure and temperature reductions. Gases evolve from the liquids and the well stream changes in character. The physical separation of these phases is by far the most common of all field-processing operations and one of the most critical. The manner in which the hydrocarbon phases are separated at the surface influences the stock-tank oil recovery. The principal means of surface separation of gas and oil is the conventional stage separation.

494 CHAPTER 5 Equations of State and Phase Equilibria

Stage separation is a process in which gaseous and liquid hydrocarbons are flashed (separated) into vapor and liquid phases by two or more separators. These separators usually are operated in a series at consecutively lower pressures. Each condition of pressure and temperature at which hydrocarbon phases are flashed is called a stage of separation. Traditionally, the stock tank is considered a separate stage of separation. Example of two-stage separation processes are shown in Fig. 5.4. To describe the liquid-gas separation process, it is convenient to group the composition of a hydrocarbon mixture into three groups of components: n

n

n

Volatile (light) components: The group includes nitrogen, methane, and ethane. Intermediate components: This group comprises components with intermediate volatility such as propane through hexane and CO2. Heavy components: This group contains components of less volatility such as heptane and heavier components.

Mechanically, there are two types of gas-oil separation a. differential separation b. flash separation.

n FIGURE 5.4 Schematic drawing of two-stage separation processes.

Vapor-Liquid Equilibrium Calculations 495

Differential Separation Process In this process, the pressure is reduced in stages and the liberated gas is removed from contact with the oil as soon as the pressure is reduced. Initially, the liberated gas is composed mainly of lighter components with intermediate and heavy components remaining in the liquid. When the gas is separated in this manner, the maximum amount of heavy and intermediate components remain in the liquid resulting in a minimum shrinkage of the oil. A greater stock-tank oil recovery will be obtained from the differential liberation process due to the fact that the lighter components liberated earlier at higher pressures and is not present at lower pressures to attract the intermediate and heavy components and pull them into the gas phase.

Flash Separation In this flash liberation process, the liberated gas remains in contact with the oil until its instantaneous removal at the final separation pressure. A maximum proportion of intermediate and heavy components are attracted into the gas phase by this process resulting in a maximum oil shrinkage and, therefore, lower oil recovery. It should be pointed out that the separation process of the solution gas taking place in the reservoir as the pressure drops below the bubble point pressure can be described experimentally by the gas differential liberation process. The main reason for this classification is that due to the high mobility of the liberated gas, the gas moves faster than the oil toward the producing well, a process similarly simulated in the laboratory through a differential liberation test by removing the liberated gas as soon as it forms. However, there is a brief period when the liberated gas saturation is less that its critical saturation; that is, gas is immobile and remains in contact with the oil, which can be described by flash liberation. The flow from the wellbore through surface facilities can be described by a series of flash separations. However, it is believed that the produced hydrocarbon system undergoes a flash liberation in the primary separator followed by a differential process as the actual separation occurs where gas and liquid are removed from contact to a subsequent stage of separations. As the number of stages increases, the differential aspect of the overall separation becomes greater. The purpose of stage separation then is to reduce the pressure on the produced oil in steps to maximize the stock-tank oil recovery results. Separator calculations are basically performed to determine

496 CHAPTER 5 Equations of State and Phase Equilibria

a. Optimum separation conditions in terms of separator pressure and temperature. b. Composition of the separated gas and oil phases. c. Oil formation volume factor. d. Producing gas-oil ratio. e. API gravity of the stock-tank oil. Stage separation of oil and gas is accomplished with a series of separator stages operating at sequentially reduced pressures when the liquid is discharged from a higher pressure separator into a lower pressure separator. The objective is to select the optimum set of working separators pressures that yield the maximum liquid recovery of liquid at stock-tank conditions. The optimization problem is constrained because the pressure of the first stage (primary separator) must be less and also impacted by the wellhead pressure and the pressure loss occurring between the two nodes. Most of optimization is performed on the subsequent separation stages, particularly the second stage. The optimization procedure is based on and supported by the following field and phase behavior observations: n

n

n

When the separator pressure is high, large amounts of light components remain in the liquid phase; however, they are liberated at the stock-tank conditions along with other valuable components (intermediate and heavy components) that attracted to the to light-end molecules. This process, when it occurs, will result in a large amount of the stock-tank oil in the gas phase at the stock tank. If the pressure is too low, large amounts of light components are separated from the liquid, and they attract substantial quantities of intermediate and heavier components resulting in a substantial stocktank oil recovery. An intermediate pressure, called the optimum separator pressure, should be selected to maximize the oil volume accumulation in the stock tank. Selecting the optimum pressure will yield 1. A maximum in the stock-tank API gravity. 2. A minimum in the oil formation volume factor (ie, less oil shrinkage). 3. A minimum in the producing gas-oil ratio (gas solubility).

This separator pressure optimization can be performed by conducting a separator test. Typical results of the test are tabulated below with the concept of determining the optimum separator pressure is shown in Fig. 5.5 by graphically plotting the API gravity, Bo, and Rs as a function of pressure, which indicates an optimum separator pressure of 100 psig.

Vapor-Liquid Equilibrium Calculations 497

Rs API

Rs

Bo

Optimum separator pressure

API

First stage separator pressure

n FIGURE 5.5 Effect of separator pressure on API, Bo, and GOR.

Psep (psig)

Tsep (°F)

Rsfb

°API

Bofb

50 to 0 Total 100 to 0 Total 200 to 0 Total 300 to 0 Total

75 75

737 41 778 676 92 768 602 178 780 549 246 795

40.5

1.481

40.7

1.474

40.4

1.483

40.1

1.495

75 75 75 75 75 75

The computational steps of the “separator calculations” are described by the following steps and in conjunction with Fig. 5.6, which schematically shows a bubble point reservoir flowing into a surface separation unit consisting of n-stages operating at successively lower pressures. Step 1. Calculate the volume of oil occupied by 1 lb-mol of the crude oil with a given overall composition “zi” at the reservoir pressure

498 CHAPTER 5 Equations of State and Phase Equilibria

Gas removed 1

Gas removed 2

(nv)1 yi

(nv)2 yi xi (nL)1

zi Separator 1

Stage 1

xi (nL)3 Gas removed

Gas removed 3 (nv)3 yi xi (nL)2

Separator 2

n

(nL)st =

Stage 2

(nL)i

i =1

Separator 3

Stage 3

Stock-tank

Well stream “zi” Pb

1 mol Vo, ft3 n FIGURE 5.6 Schematic illustration of n-separation stages.

(or bubble point pressure) and temperature. This volume, denoted as Vo, is calculated by recalling and applying the equation that defines the number of moles; that is, n¼

m ρo Vo ¼ ¼1 Ma Ma

(5.36)

Solving for the oil volume gives X Vo ¼

Ma ¼ ρo

ðzi M i Þ

i

ρo

(5.37)

where m ¼ total weight of 1 lb-mol of the crude oil, lb/mol Vo ¼ volume of 1 lb-mol of the crude oil at reservoir conditions, ft3/mol Ma ¼ apparent molecular weight Mi ¼ molecular weight of component i ρo ¼ density of the reservoir oil at the bubble point pressure, lb/ft3 It should be pointed out that the density can be calculated from the crude oil composition using the Standing-Katz or Alani-Kennedy method.

Vapor-Liquid Equilibrium Calculations 499

Step 2. Given the composition of the hydrocarbon stream “zi” to the first separator and the operating conditions of the separator, calculate a set of equilibrium ratios that describes the hydrocarbon mixture. Step 3. Assuming a total of 1 mol of the feed (ie, nt ¼ 1) entering the first separator, perform flash calculations to obtain n compositions of the separated gas “yi” and liquid “xi” n number of moles of the gas “nv” and the liquid “nL”. Designating these moles as (nL)1 and (nυ)1, the actual number of moles of the gas and the liquid leaving the first separation stage are ½nυ1 a ¼ ðnt Þðnυ Þ1 ¼ ð1Þðnυ Þ1 ½nL1 a ¼ ðnt ÞðnL Þ1 ¼ ð1ÞðnL Þ1

where [nv1]a ¼ actual number of moles of vapor leaving the first separator [nL1]a ¼ actual number of moles of liquid leaving the first separator. Step 4. Using the composition of the liquid “xi” leaving the first separator as the feed for the second separator “zi ¼ xi” and calculate the equilibrium ratios of the hydrocarbon mixture at the prevailing pressure and temperature of the second-stage separator. Step 5. Assume 1 mol of the feed from first-stage separator to the second stage, perform flash calculations to determine the compositions and quantities of the gas and liquid leaving the second separation stage. The actual number of moles of the two phases then is calculated from ½nυ2 a ¼ ½nL1 a ðnυ Þ2 ¼ ð1ÞðnL Þ1 ðnυ Þ2 ½nL2 a ¼ ½nL1 a ðnL Þ2 ¼ ð1ÞðnL Þ1 ðnL Þ2

where [nv2]a, [nL2]a ¼ actual moles of gas and liquid leaving separator 2 (nv)2, (nL)2 ¼ moles of gas and liquid as determined from flash calculations Step 6. The previously outlined procedure is repeated for each separation stage, including the stock-tank storage, and the calculated moles and compositions are recorded. The total number of moles of liberated solution gas from all stages are calculated as ðnν Þt ¼

n X ðnva Þi ¼ ðnν Þ1 + ðnL Þ1 ðnν Þ2 + ðnL Þ1 ðnL Þ2 ðnν Þ3 + ⋯ + ðnL Þ1 …ðnL Þn1 ðnν Þn i¼1

500 CHAPTER 5 Equations of State and Phase Equilibria

In a more compact form, this expression can be written as ðnν Þt ¼ ðnν Þ1 +

" # n i1 Y X ðnν Þi ðnL Þj i¼2

(5.38)

j¼1

where (nv)t ¼ total moles of liberated gas from all stages, lb-mol/mol of feed n ¼ number of separation stages. The total moles of liquid remaining in the stock tank can also be calculated as ðnL Þst ¼ nL1 nL2 …nLn

or ðnL Þst ¼

n Y

ðnL Þi

(5.39)

i¼1

where (nL)st ¼ number of moles of liquid remaining in the stock tank. Step 7. Calculate the volume, in scf, of all the liberated solution gas from Vg ¼ 379:4ðnυ Þt

(5.40)

where Vg ¼ total volume of the liberated solution gas scf/mol of feed. Step 8. Determine the volume of stock-tank oil occupied by (nL)st moles of liquid from ðVo Þst ¼

ðnL Þst ðMa Þst ðρo Þst

(5.41)

where (Vo)st ¼ volume of stock-tank oil, ft3/mol of feed (Ma)st ¼ apparent molecular weight of the stock-tank oil (ρo)st ¼ density of the stock-tank oil, lb/ft3. This desnsity is calculated from the oil composition at stock-tank Step 9. Calculate the specific gravity and the API gravity of the stock-tank oil by applying the expression γo ¼

ðρo Þst 62:4

(5.42)

Vapor-Liquid Equilibrium Calculations 501

Step 10. Calculate the total gas-oil ratio (gas solubility) Rs: GOR ¼

Vg ð5:615Þð379:4Þðnv Þt ¼ ðVo Þst =5:615 ðnL Þst ðMÞst =ðρo Þst GOR ¼

2130:331ðnv Þt ðρo Þst ðnL Þst ðMÞst

(5.43)

where GOR ¼ gas-oil ratio, scf/STB. Step 11. Calculate the oil formation volume factor from the relationship Bo ¼

Vo ðVo Þst

Substituting Eqs. (5.36), (5.41) with equivalent term in the above expression gives Bo ¼

Ma ðρo Þst ρo ðnL Þst ðMa Þst

(5.44)

where Bo ¼ oil formation volume factor, bbl/STB Ma ¼ apparent molecular weight of the feed (Ma)st ¼ apparent molecular weight of the stock-tank oil ρo ¼ density of crude oil at reservoir conditions, lb/ft3 The separator pressure can be optimized by calculating the API gravity, GOR, and Bo in the manner just outlined at different assumed pressures. The optimum pressure corresponds to a maximum in the API gravity and a minimum in gas-oil ratio and oil formation volume factor.

EXAMPLE 5.6 A crude oil, with the composition as follows, exists at its bubble point pressure of 1708.7 psia and a temperature of 131°F. The surface separation facilities include two separators and stock tank with the following operating conditions: Stage #

Pressure (psia)

Temperature (°F)

1 2 Stock tank

400 350 14.7

72 72 60

502 CHAPTER 5 Equations of State and Phase Equilibria

The composition of the crude oil “feed” to the primary separator is tabulated below: Mi Component zi CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.0008 0.0164 0.2840 0.0716 0.1048 0.0420 0.0420 0.0191 0.1912 0.0405 0.3597

44.0 28.0 164.0 30.0 44.0 584.0 58.0 72.0 72.0 86.0 252

Given that the molecular weight and specific gravity of C7+ are 252 and 0.8429, calculate n n n

Oil formation volume factor at the bubble point pressure Bofb Gas solubility at the bubble point pressure Rsfb. Stock-tank density and the API gravity of the hydrocarbon system.

Solution Step 1. Calculate the apparent molecular weight of the crude oil to give X Ma ¼ zi Mi ¼ 113:5102 Step 2. Calculate the oil density at bubble point pressure by using the Standing and Katz correlations to yield ρob ¼ 44:794 lb=ft3 Step 3. Using the pressure and temperature of the primary separator (ie, 400 psia and 72°F) calculate the K-value for each component by applying the Standing’s K-value, to give 1 Ki ¼ 10ða + cFi Þ p Component

zi

Ki

CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.0008 0.0164 0.2840 0.0716 0.1048 0.0420 0.0420 0.0191 0.0191 0.0405 0.3597

3.509 39.90 8.850 1.349 0.373 0.161 0.120 0.054 0.043 0.018 0.0021

Vapor-Liquid Equilibrium Calculations 503

Step 4. Perform flash calculations on the hydrocarbon stream to the first separator to give nL ¼ 0:7209 nv ¼ 0:29791 Step 5. Solve for xi and yi to give xi ¼

zi nL + nv Ki

yi ¼ xi Ki

Component

zi

Ki

xi

yi

Mi

CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.0008 0.0164 0.2840 0.0716 0.1048 0.0420 0.0420 0.0191 0.0191 0.0405 0.3597

3.509 39.90 8.850 1.349 0.373 0.161 0.120 0.054 0.043 0.018 0.0021

0.0005 0.0014 0.089 0.0652 0.1270 0.0548 0.0557 0.0259 0.0261 0.0558 0.4986

0.0018 0.0552 0.7877 0.0880 0.0474 0.0088 0.0067 0.0014 0.0011 0.0010 0.0009

44.0 28.0 16.0 30.0 44.0 58.0 58.0 72.0 72.0 86.0 252

Step 6. Use the calculated liquid composition as the feed for the second separator or flash the composition at the operating condition of the separator. The results are shown in the table below, with nL ¼ 0.9851 and nυ ¼ 0.0149. Component

zi

Ki

xi

yi

Mi

CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.0005 0.0014 0.089 0.0652 0.1270 0.0548 0.0557 0.0259 0.0261 0.0558 0.4986

3.944 46.18 10.06 1.499 0.4082 0.1744 0.1291 0.0581 0.0456 0.0194 0.00228

0.0005 0.0008 0.0786 0.0648 0.1282 0.0555 0.0564 0.0263 0.0264 0.0566 0.5061

0.0018 0.0382 0.7877 0.0971 0.0523 0.0097 0.0072 0.0015 0.0012 0.0011 0.0012

44.0 28.0 16.0 30.0 44.0 58.0 58.0 72.0 72.0 86.0 252

Step 7. Repeat the calculation for the stock-tank stage, to give the results in the following table, with nL ¼ 0.6837 and nυ ¼ 0.3163.

504 CHAPTER 5 Equations of State and Phase Equilibria

Component

zi

Ki

xi

yi

CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.0005 0.0008 0.0784 0.0648 0.1282 0.0555 0.0564 0.0263 0.0264 0.0566 0.5061

81.14 1159 229 27.47 6.411 2.518 1.805 0.7504 0.573 0.2238 0.03613

0000 0000 0.0011 0.0069 0.0473 0.0375 0.0450 0.0286 0.02306 0.0750 0.7281

0.0014 0.026 0.2455 0.1898 0.3030 0.0945 0.0812 0.0214 0.0175 0.0168 0.0263

Step 8. Calculate the actual number of moles of the liquid phase at the stocktank conditions from Eq. (5.39): ðnL Þst ¼

n Y

ðnL Þi

i¼1

ðnL Þst ¼ ð1Þð0:7209Þð0:9851Þð0:6837Þ ¼ 0:48554 Step 9. Calculate the total number of moles of the liberated gas from the entire surface separation system: nυ ¼ 1  ðnL Þst ¼ 1  0:48554 ¼ 0:51446 Step 10. Calculate apparent molecular weight of the stock-tank oil from its composition, to give ðMa Þst ¼ ΣMi xi ¼ 200:6 Step 11. Using the composition of the stock-tank oil, calculate the density of the stock-tank oil by using Standing’s correlation, to give ðρo Þst ¼ 50:920 Step 12. Calculate the API gravity of the stock-tank oil from the specific gravity: γ ¼ ρo =62:4 ¼ 50:920=62:4 ¼ 0:81660°=60° API ¼ ð141:5=γ Þ  131:5 ¼ ð141:5=0:816Þ  131:5 ¼ 41:9 Step 13. Calculate the gas solubility from Eq. (5.43) to give Rs ¼

2130:331ðnv Þt ðρo Þst ðnL Þst ðMÞst

Rs ¼

2130:331ð0:51446Þð50:92Þ ¼ 573:0 scf=STB 0:48554ð200:6Þ

Step 14. Calculate Bo from Eq. (5.44) to give

Equations of State 505

Bo ¼

Ma ðρo Þst ρo ðnL Þst ðMa Þst

Bo ¼

ð113:5102Þð50:92Þ ¼ 1:325bbl=STB ð44:794Þð0:48554Þð200:6Þ

To optimize the operating pressure of the separator, these steps should be repeated several times under different assumed pressures and the results, in terms of API, Bo, and Rs, should be expressed graphically and used to determine the optimum pressure.

EQUATIONS OF STATE An equation of state (EOS) is an analytical expression relating the pressure, p, to the temperature, T, and the volume, V. A proper description of this PVT relationship for real hydrocarbon fluids is essential in determining the volumetric and phase behavior of petroleum reservoir fluids and predicting the performance of surface separation facilities; these can be described accurately by equations of state. In general, most equations of state require only the critical properties and acentric factor of individual components. The main advantage of using an EOS is that the same equation can be used to model the behavior of all phases, thereby assuring consistency when performing phase equilibria calculations. The best known and the simplest example of an equation of state is the ideal gas equation, expressed mathematically by the expression p¼

RT V

(5.45)

where V ¼ gas volume in ft3 per 1 mol of gas, ft3/mol. The PVT relationship is used to describe the volumetric behavior of hydrocarbon gases at pressures close to the atmospheric pressure for which it was experimentally derived. However, the extreme limitations of the applicability of Eq. (5.45) prompted numerous attempts to develop an equation of state suitable for describing the behavior of real fluids at extended ranges of pressures and temperatures. A review of recent developments and advances in the field of empirical cubic equations of state are presented next, along with samples of their applications in petroleum engineering. Four equations of state are discussed here: n n n n

van der Waals Redlich-Kwong Soave-Redlich-Kwong Peng-Robinson

506 CHAPTER 5 Equations of State and Phase Equilibria

The Van der Waals Equation of State “VdW EOS” In developing the ideal gas EOS, two assumptions were made: n

n

The volume of the gas molecules is insignificant compared to the total volume and distance between the molecules. There are no attractive or repulsive forces between the molecules.

Van der Waals (1873) attempted to eliminate these two assumptions in developing an empirical equation of state for real gases. In his attempt to eliminate the first assumption, van der Waals pointed out that the gas molecules occupy a significant fraction of the volume at higher pressures and proposed that the volume of the molecules, denoted by the parameter b, be subtracted from the actual molar volume, V, in Eq. (5.45), to give p¼

RT V b

where the parameter b is known as the covolume and considered to reflect the volume of molecules. The variable V represents the actual volume in ft3 per 1 mol of gas. To eliminate the second assumption, van der Waals subtracted a corrective term, denoted by a/V2, from this equation to account for the attractive forces between molecules. In a mathematical form, van der Waals proposed the following expression: p¼

RT a  V  b V2

(5.46)

where p ¼ system pressure, psia T ¼ system temperature, °R R ¼ gas constant, 10.73 psi-ft3/lb-mol, °R V ¼ volume, ft3/mol a ¼ “attraction” parameter b ¼ “repulsion” parameter The two parameters, a and b, are constants characterizing the molecular properties of the individual components. The symbol a is considered a measure of the intermolecular attractive forces between the molecules. Eq. (5.46) shows the following important characteristics: 1. At low pressures, the volume of the gas phase is large in comparison with the volume of the molecules. The parameter b becomes negligible in comparison with V and the attractive forces term, a/V2, becomes insignificant; therefore, the van der Waals equation reduces to the ideal gas equation [Eq. (5.45)].

Equations of State 507

2. At high pressure (ie, p ! 1) the volume, V, becomes very small and approaches the value b, which is the actual molecular volume, as mathematically given by lim V ðpÞ ¼ b

p!1

The van der Waals or any other equation of state can be expressed in a more generalized form as follows: p ¼ prepulsion  pattraction

(5.47)

where the repulsion pressure term, prepulsion, is represented by the term RT/(V  b), and the attraction pressure term, pattraction, is described by a/V2. In determining the values of the two constants, a and b, for any pure substance, van der Waals observed that the critical isotherm has a horizontal slope and an inflection point at the critical point, as shown in Fig. 5.7. This observation can be expressed mathematically as follows: 

 @p ¼0 @V Tc, pc  2  @ p ¼0 @V 2 Tc, pc

(5.48)

Differentiating Eq. (5.46) with respect to the volume at the critical point results in 

@p @V

 Tc , pc

Liquid

¼

RTc ðVc  bÞ

+

2a ¼0 Vc3

(5.49)

Gas

H

VdW EOS critical point condition

[ Pressure

2

∂p ∂2p ] T , p ,V = 0 & [ 2 ] T , p ,V = 0 c c c ∂V ∂V c c c T1 < T2 < T3 < Tc < T4

C

Pc

Liquid – Vapor

G F

B Vc

A

Volume V

n FIGURE 5.7 An idealized pressure-volume relationship for a pure component.

T4 Tc T3 T E 2 T1

508 CHAPTER 5 Equations of State and Phase Equilibria



@2p @V 2

 Tc , pc

¼

2RTc ðVc  bÞ3

+

6a ¼0 Vc4

(5.50)

Solving Eqs. (5.49), (5.50) simultaneously for the parameters a and b gives   1 Vc 3   8 a¼ RTc Vc 9 b¼

(5.51)

(5.52)

Eq. (5.51) suggests that the covolume “b” is approximately 0.333 of the critical volume “Vc” of the substance. Numerous experimental studies revealed that the co-volume “b” is in the range 0.24–0.28 of the critical volume and pure component. Applying Eq. (5.48) to the critical point (i.e., by setting T ¼ Tc, p ¼ pc, and V ¼ Vc) and combining with Eqs. (5.51), (5.52) yields pc Vc ¼ ð0:375ÞRTc

(5.53)

Eq. (5.53) indicates that regardless of the type of the substance, the Van der Waals EOS gives a universal critical gas compressibility factor “Zc” of 0.375. Experimental studies show that Zc values for substances are ranging between 0.23 and 0.31. Eq. (5.53) can be combined with Eqs. (5.51), (5.52) to give more convenient and traditional expressions for calculating the parameters a and b, to yield a ¼ Ωa

R2 Tc2 pc

(5.54)

RTc pc

(5.55)

b ¼ Ωb

where R ¼ gas constant, 10.73 psia-ft3/lb-mol-°R pc ¼ critical pressure, psia Tc ¼ critical temperature, °R Ωa ¼ 0.421875 Ωb ¼ 0.125 Eq. (5.46) can also be expressed in a cubic form in terms of the volume “V” as follows: p¼

RT a  2 V b V

Equations of State 509

Rearranging gives       RT 2 a ab V  b+ V + V ¼0 p p p 3

(5.56)

Eq. (5.56) is referred to as the van der Waals two-parameter cubic equation of state. The term two-parameter refers to parameters a and b. The term cubic equation of state implies that the equation have three possible values for the volume “V,” at least one of which is real. Perhaps the most significant features of Eq. (5.56) are that it describes the liquid-condensation phenomenon and the passage from the gas to the liquid phase as the gas is compressed. These important features of the van der Waals EOS are discussed next in conjunction with Fig. 5.8. Consider a pure substance with a p-V behavior as shown in Fig. 5.8. Assume that the substance is kept at a constant temperature, T, below its critical temperature. At this temperature, Eq. (5.56) has three real roots (volumes) for each specified pressure, p. A typical solution of Eq. (5.56) at constant temperature, T, is shown graphically by the dashed isotherm, the constant temperature curve DWEZB, in Fig. 5.8. The three values of V are the intersections B, E, and D on the horizontal line, corresponding to a fixed value of the pressure. This dashed calculated line (DWEZB) then appears to give a continuous transition from the gaseous phase to the liquid phase, but in reality, the transition is abrupt and discontinuous, with both liquid and vapor existing along the straight horizontal line DB. Examining the graphical solution of Eq. (5.56) shows that the largest root (volume), indicated by point D, corresponds to the volume of the saturated vapor, while the smallest positive volume, indicated by point B, corresponds to the volume of the C

Pressure

Pure component isotherm as calculated from VdW EOS W D

P B

E Z

Volume

n FIGURE 5.8 Pressure-volume diagram for a pure component.

T

510 CHAPTER 5 Equations of State and Phase Equilibria

saturated liquid. The third root, point E, has no physical meaning. Note that these values become identical as the temperature approaches the critical temperature, Tc, of the substance. Eq. (5.56) can be expressed in a more practical form in terms of the compressibility factor, Z. Replacing the molar volume, V, in Eq. (5.56) with ZRT/p gives       RT 2 a ab V3 ¼ b + V + V ¼0 p p p   2      RT ZRT a ZRT ab +  ¼0 V3  b + p p p p p

or Z 3  ð1 + BÞZ2 + AZ  AB ¼ 0

(5.57)

where A¼

ap R2 T 2

(5.58)

bp RT

(5.59)

B¼

Z ¼ compressibility factor p ¼ system pressure, psia T ¼ system temperature, °R Eq. (5.57) yields one real root1 in the one-phase region and three real roots in the two-phase region (where system pressure equals the vapor pressure of the substance). In the two-phase region, the largest positive root corresponds to the compressibility factor of the vapor phase, Zv, while the smallest positive root corresponds to that of the liquid phase, ZL. An important practical application of Eq. (5.57) is density calculations, as illustrated in the following example.

EXAMPLE 5.7 A pure propane is held in a closed container at 100°F. Both gas and liquid are present. Calculate, using the van der Waals EOS, the density of the gas and liquid phases.

1

In some super critical regions, Eq. (5.57) can yield three real roots for Z. From the three real roots, the largest root is the value of the compressibility with physical meaning.

Equations of State 511

Solution Step 1. Determine the vapor pressure, pv, of the propane from the Cox chart. This is the only pressure at which two phases can exist at the specified temperature: pv ¼ 185 psi Step 2. Calculate parameters a and b from Eqs. (5.54), (5.55), respectively: a ¼ Ωa

R2 Tc2 ð10:73Þ2 ð666Þ2 ¼ 0:421875 ¼ 34,957:4 pc 616:3

and b ¼ Ωb

RTc 10:73ð666Þ ¼ 0:125 ¼ 1:4494 616:3 pc

Step 3. Compute the coefficients A and B by applying Eqs. (5.58), (5.59), respectively: A¼

ap ð34, 957:4Þð185Þ ¼ ¼ 0:179122 R2 T 2 ð10:73Þ2 ð560Þ2

B¼

bp ð1:4494Þð185Þ ¼ ¼ 0:044625 RT ð10:73Þð560Þ

Step 4. Substitute the values of A and B into Eq. (5.57), to give Z 3  ð1 + BÞZ2 + AZ  AB ¼ 0 Z 3  1:044625Z2 + 0:179122Z  0:007993 ¼ 0 Step 5. Solve the above third-degree polynomial by extracting the largest and smallest roots of the polynomial using the appropriate direct or iterative method, to give Zv ¼ 0:72365 ZL ¼ 0:07534 Step 6. Solve for the density of the gas and liquid phases using Eq. (2.17): ρg ¼

pM ð185Þð44:0Þ ¼ ¼ 1:87 lb=ft3 Z v RT ð0:72365Þð10:73Þð560Þ

and ρL ¼

pM ð185Þð44Þ ¼ ¼ 17:98 lb=ft3 ZL RT ð0:07534Þð10:73Þð560Þ

The van der Waals equation of state, despite its simplicity, provides a correct description, at least qualitatively, of the PVT behavior of substances in the liquid and gaseous states. Yet, it is not accurate enough to be suitable for design purposes. With the rapid development of computers, the equation-of-state approach for calculating physical properties and phase equilibria proved to be a powerful tool, and much energy was devoted to the development of new and accurate equations

512 CHAPTER 5 Equations of State and Phase Equilibria

of state. Many of these equations are a modification of the van der Waals EOS and range in complexity from simple expressions containing 2 or 3 parameters to complicated forms containing more than 50 parameters. Although the complexity of any equation of state presents no computational problem, most authors prefer to retain the simplicity found in the van der Waals cubic equation while improving its accuracy through modifications. All equations of state generally are developed for pure fluids first and then extended to mixtures through the use of mixing rules.

Redlich-Kwong’s Equation of State “RK EOS” Redlich and Kwong (1949) observed that the van der Waals a/V2 term does not contain the system temperature to account for its impact on the intermolecular attractive force between the molecules. Redlich and Kwong demonstrated that by a simple adjustment of the van der Waals’ a/V2 term to explicitly include the system temperature could considerably improve the prediction of the volumetric and physical properties of the vapor phase. Redlich and Kwong replaced the attraction pressure term with a generalized temperature dependence term, as given in the following form: p¼

RT a pffiffiffi  V  b V ðV + bÞ T

(5.60)

in which T is the system temperature, in °R. Redlich and Kwong noted that, as the system pressure becomes very large (ie, p ! 1) the molar volume, V, of the substance shrinks to about 26% of its critical volume, Vc, regardless of the system temperature. Accordingly, they constructed Eq. (5.60) to satisfy the following condition: b ¼ 0:26Vc

(5.61)

Imposing the critical point conditions (as expressed by Eq. (5.48)) on Eq. (5.60) and solving the resulting two equations simultaneously, gives a ¼ Ωa

R2 Tc2:5 pc

(5.62)

RTc pc

(5.63)

b ¼ Ωb

where Ωa ¼ 0.42747 Ωb ¼ 0.08664.

Equations of State 513

It should be noted that by equating Eq. (5.63) with (5.61) gives 0:26Vc ¼ Ωb

RTc pc

Arranging the above expression: pc Vc ¼ 0:33RTc

(5.64)

Eq. (5.64) shows that the RK EOS produces a universal critical compressibility factor (Zc) of 0.333 for all substances. As indicated earlier, the critical gas compressibility ranges from 0.23 to 0.31 for most of the substances, with an average value of 0.27. Replacing the molar volume, V, in Eq. (5.60) with ZRT/p and rearranging gives   Z3  Z2 + A  B  B2 Z  AB ¼ 0

(5.65)

where A¼

ap R2 T 2:5

(5.66)

bp RT

(5.67)

B¼

As in the van der Waals EOS, Eq. (5.65) yields one real root in the one-phase region (gas-phase region or liquid-phase region) and three real roots in the two-phase region. In the latter case, the largest root corresponds to the compressibility factor of the gas phase, Zv, while the smallest positive root corresponds to that of the liquid, ZL. EXAMPLE 5.8 Rework Example 5.7 using the Redlich-Kwong equation of state. Solution Step 1. Calculate the parameters a, b, A, and B. a ¼ Ωa

R2 Tc2:5 ð10:73Þ2 ð666Þ2:5 ¼ 0:4247 ¼ 914, 1101:1 pc 616:3

b ¼ Ωb

RTc ð10:73Þð666Þ ¼ 0:08664 ¼ 1:0046 pc 616:3

A¼

ð914, 110:1Þð185Þ ap ¼ ¼ 0:197925 R2 T 2:5 ð10:73Þ2 ð560Þ2:5

B¼

bp ð1:0046Þð185Þ ¼ ¼ 0:03093 RT ð10:73Þð560Þ

Step 2. Substitute the parameters A and B into Eq. (5.65) and extract the largest and the smallest roots, to give

514 CHAPTER 5 Equations of State and Phase Equilibria

  Z3  Z2 + A  B  B2 Z  AB ¼ 0 Z 3  Z2 + 0:1660384Z  0:0061218 ¼ 0 Largest root : Z v ¼ 0:802641 Smallest root : Z L ¼ 0:0527377 Step 3. Solve for the density of the liquid phase and gas phase: ρL ¼

pM ð185Þð44Þ ¼ ¼ 25:7 lb=ft3 ZL RT ð0:0527377Þð10:73Þð560Þ

ρv ¼

pM ð185Þð44Þ ¼ ¼ 1:688 lb=ft3 Zv RT ð0:0802641Þð10:73Þð560Þ

Redlich and Kwong extended the application of their equation to hydrocarbon liquid and gas mixtures by employing the following mixing rules: " #2 n X pffiffiffiffi am ¼ xi ai (5.68) i¼1

bm ¼

n X ½xi bi 

(5.69)

i¼1

where n ¼ number of components in the mixture ai ¼ Redlich-Kwong a parameter for the ith component as given by Eq. (5.62) bi ¼ Redlich-Kwong b parameter for the ith component as given by Eq. (5.63) am ¼ parameter a for mixture bm ¼ parameter b for mixture xi ¼ mole fraction of component i in the liquid phase To calculate am and bm for a hydrocarbon gas mixture with a composition of yi, use Eqs. (5.68), (5.69) and replace xi with yi: " #2 n X pffiffiffiffi am ¼ yi ai i¼1 n X bm ¼ ½yi bi  i¼1

Eq. (5.65) gives the compressibility factor of the gas phase or the liquid with the coefficients A and B as defined by Eqs. (5.66), (5.67) as A¼

am p R2 T 2:5

B¼

bm p RT

The application of the Redlich and Kwong equation of state for hydrocarbon mixtures can be best illustrated through the following two examples (Examples 5.9 and 5.10).

Equations of State 515

EXAMPLE 5.9 Calculate the density of a crude oil with the composition at 4000 psia and 160°F given in the table below. Use the Redlich-Kwong EOS. Component

xi

M

pc

Tc

C1 C2 C3 n-C4 n-C5 C6 C7+

0.45 0.05 0.05 0.03 0.01 0.01 0.40

16.043 30.070 44.097 58.123 72.150 84.00 215

666.4 706.5 616.0 527.9 488.6 453 285

343.33 549.92 666.06 765.62 845.8 923 1287

Solution Step 1. Determine the parameters ai and bi for each component using Eqs. (5.62), (5.63): ai ¼ 0:47274

R2 Tc2:5 pc

bi ¼ 0:08664

RTc pc

The results are shown in the table below. Component xi C1 C2 C3 n-C4 n-C5 C6 C7+

M

pc

0.45 16.043 0.05 30.070 0.05 44.097 0.03 58.123 0.01 72.150 0.01 84.00 0.40 215

Tc

666.4 343.33 706.5 549.92 616.0 666.06 527.9 765.62 488.6 845.8 453 923 285 1287

ai

bi

161,044.3 493,582.7 914,314.8 1,449,929 2,095,431 2,845,191 1.022348 (107)

0.4780514 0.7225732 1.004725 1.292629 1.609242 1.945712 4.191958

Step 2. Calculate the mixture parameters, am and bm, from Eqs. (5.68), (5.69), to give " #2 n X pffiffiffiffi am ¼ xi ai ¼ 2,591, 967 i¼1

and bm ¼

n X ½xi bi  ¼ 2:0526 i¼1

Step 3. Compute the coefficients A and B by using Eqs. (5.66), (5.67), to produce

516 CHAPTER 5 Equations of State and Phase Equilibria

A¼

am p 2,591, 967ð4000Þ ¼ ¼ 9:406539 R2 T 2:5 10:732 ð620Þ2:5

B¼

bm p 2:0526ð4000Þ ¼ ¼ 1:234049 RT 10:73ð620Þ

Step 4. Solve Eq. (5.65) for the largest positive root, to yield Z3  Z2 + 6:934845Z  11:60813 ¼ 0 ZL ¼ 1:548126 Step 5. Calculate the apparent molecular weight of the crude oil: X xi Mi ¼ 110:2547 Ma ¼ Step 6. Solve for the density of the crude oil: ρL ¼

pMa ð4000Þð100:2547Þ ¼ ¼ 38:93 lb=ft3 ZL RT ð10:73Þð620Þð1:548120Þ

Note that calculating the liquid density as by Standing’s method gives a value of 46.23 lb/ft3.

EXAMPLE 5.10 Using the the RK EOS, calculate the density of a gas phase with the composition at 4000 psia and 160°F shown in the following table. Component

yi

M

pc

Tc

C1 C2 C3 C4 C5 C6 C7+

0.86 0.05 0.05 0.02 0.01 0.005 0.005

16.043 30.070 44.097 58.123 72.150 84.00 215

666.4 706.5 616.0 527.9 488.6 453 285

343.33 549.92 666.06 765.62 845.8 923 1287

Solution Step 1. Determine the parameters ai and bi for each component using Eqs. (5.62), (5.63): ai ¼ 0:47274

R2 Tc2:5 pc

bi ¼ 0:08664

RTc pc

The results are shown in the table below.

Equations of State 517

Component yi C1 C2 C3 C4 C5 C6 C7+

M

pc

0.86 16.043 0.05 30.070 0.05 44.097 0.02 58.123 0.01 72.150 0.005 84.00 0.005 215

Tc

666.4 343.33 706.5 549.92 616.0 666.06 527.9 765.62 488.6 845.8 453 923 285 1287

ai

bi

161,044.3 493,582.7 914,314.8 1,449,929 2,095,431 2,845,191 1.022348(107)

0.4780514 0.7225732 1.004725 1.292629 1.609242 1.945712 4.191958

Step 2. Calculate am and bm using Eqs. (5.68), (5.69), to give " #2 n X pffiffiffiffi am ¼ yi ai ¼ 241,118 i¼1

X bi xi ¼ 0:5701225 bm ¼ Step 3. Calculate the coefficients A and B by applying Eqs. (5.66), (5.67), to yield A¼

am p 241,118ð4000Þ ¼ ¼ 0:8750 R2 T 2:5 10:733 ð620Þ2:5

B¼

bm p 0:5701225ð4000Þ ¼ ¼ 0:3428 RT 10:73ð620Þ

Step 4. Solve Eq. (5.65) for Zv, to give Z 3  Z 2 + 0:414688Z  0:29995 ¼ 0 Z v ¼ 0:907 Step 5. Calculate the apparent molecular weight and density of the gas mixture: P Ma ¼ yi Mi ¼ 20:89 ρv ¼

pMa ð4000Þð20:89Þ ¼ ¼ 13:85 lb=ft3 Zv RT ð10:73Þð620Þð0:907Þ

Soave-Redlich-Kwong Equation of State “SRK EOS” One of the most significant milestones in the development of cubic equations of state is accredited to Soave (1972) with his proposed modification of the RK EOS attraction pressure term. Soave replaced the explicit temperature term (a/T0.5) in Eq. (5.60) with a more generalized temperaturedependent term, denoted by [a α (T)]; that is, p¼

RT aαðT Þ  V  b V ðV + bÞ

(5.70)

518 CHAPTER 5 Equations of State and Phase Equilibria

where α(T) is a dimensionless parameter that is defined by the following expression: h pffiffiffiffiffiffiffiffiffiffi i2 αðT Þ ¼ 1 + m 1  T=Tc

Or equivalently in terms of reduced temperature “Tr”:   pffiffiffiffiffi 2 α ðT Þ ¼ 1 + m 1  T r

(5.71)

Soave imposed the following conditions on the a dimensionless parameter α(T): n

n

α(T) ¼ 1 when the reduced temperature Tr, ¼ 1; that is, α(T ¼ Tc) ¼ 1 when T/Tc ¼ 1 At temperatures other than the critical temperature, Soave regressed on the vapor pressure of pure components to develop a temperature correction parameter “m” that correlated with the acentric factor “ω” by the following relationship: m ¼ 0:480 + 1:74ω  0:176ω2

(5.72)

where Tr ¼ reduced temperature, T/Tc ω ¼ acentric factor of the substance T ¼ system temperature, °R For simplicity and convenience, the temperature-dependent term “α(T)” is replaced by the symbol “α” for the remainder of this chapter. For any pure component, the constants “a” and “b” in Eq. (5.70) are found by imposing the classical van der Waals’s critical point constraints on Eq. (5.70) and solving the resulting two equations, to give a ¼ Ωa

R2 Tc2 pc

(5.73)

RTc pc

(5.74)

b ¼ Ωb

where Ωa ¼ 0:42747 Ωb ¼ 0:08664

Edmister and Lee (1986) pointed out that the following constrain satisfies the critical isotherm as given by   ðV  Vc Þ3 ¼ V 3  ½3Vc V 2 + 3Vc2 V  Vc3 ¼ 0

(5.75)

Equations of State 519

Eq. (5.70) also can be expressed into a cubic form, to give V3 

      RT 2 aα bRT ðaαÞb V +   b2 V  ¼0 p p p p

(5.76)

At the critical point, Eqs. (5.75), (5.76) are identical with α ¼ 1. Equating the corresponding coefficients in both equations gives RTc pc

(5.77)

a bRTc   b2 pc pc

(5.78)

ab pc

(5.79)

3Vc ¼ 3Vc2 ¼

Vc3 ¼

Solving these equations simultaneously for parameters a and b will give relationships identical to those given by Eqs. (5.73), (5.74), in addition: n

rearranging Eq. (5.77) gives 1 pc Vc ¼ RTc 3

which indicates that the SRK EOS gives a universal critical gas compressibility factor “Zc” of 0.333. n

combining Eq. (5.74) with (5.77), gives the expected value of the covolume of 26% of the critical volume; that is, b ¼ 0:26Vc

The compressibility factor “Z” can be introduced into Eq. (5.76) by replacing the molar volume “V” with (ZRT/p) and rearranging to give   Z3  Z2 + A  B  B2 Z  AB ¼ 0

(5.80)

with A¼

ðaαÞp ðRT Þ2

B¼

bp RT

(5.81)

(5.82)

520 CHAPTER 5 Equations of State and Phase Equilibria

where p ¼ system pressure, psia T ¼ system temperature, °R R ¼ 10.730 psia ft3/lb-mol °R EXAMPLE 5.11 Rework Example 5.7 and solve for the density of the two phases using the SRK EOS. Solution Step 1. Determine the critical pressure, critical temperature, and acentric factor from Table 1.1 of Chapter 1, to give Tc ¼ 666:01°R pc ¼ 616:3 psia ω ¼ 0:1524 Step 2. Calculate the reduced temperature: Tr ¼ T=Tc ¼ 560=666:01 ¼ 0:8408 Step 3. Calculate the parameter m by applying Eq. (5.72), to yield m ¼ 0:480 + 1:574ω  0:176ω2 m ¼ 0:480 + 1:574ð0:1524Þ  0:176ð1:524Þ2 Step 4. Solve for parameter α using Eq. (5.71), to give   pffiffiffiffiffi 2 α ¼ 1 + m 1  Tr pffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2   α ¼ 1 + 0:7052 1  0:8408 ¼ 1:120518 Step 5. Compute the coefficients a and b by applying Eqs. (5.73), (5.74), to yield α ¼ Ωa

R2 Tc2 10:732 ð666:01Þ2 ¼ 35, 427:6 ¼ 0:42747 pc 616:3

α ¼ Ωb

RTc 10:73ð666:01Þ ¼ 0:08664 ¼ 1:00471 616:3 pc

Step 6. Calculate the coefficients A and B from Eqs. (5.81), (5.82) to produce A¼

ðaαÞp ð35, 427:6Þð1:120518Þ185 ¼ ¼ 0:203365 R2 T 2 10:732 ð560Þ2

B¼

bp ð1:00471Þð185Þ ¼ ¼ 0:034658 RT ð10:73Þð560Þ

Step 7. Solve Eq. (5.80) for ZL and Zv:   Z 3  Z2 + A  B  B2 Z + AB ¼ 0   Z 3  Z2 + 0:203365  0:034658  0:0346582 Z + ð0:203365Þð0:034658Þ ¼ 0 Solving this third-degree polynomial gives ZL ¼ 0:06729 Zv ¼ 0:80212

Equations of State 521

Step 8. Calculate the gas and liquid density, to give ρv ¼

pM ð185Þð44:0Þ ¼ 1:6887 lb=ft3 ¼ Z v RT ð0:802121Þð10:73Þð560Þ

ρL ¼

pM ð185Þð44:0Þ ¼ ¼ 20:13 lb=ft3 L Z RT ð0:06729Þð10:73Þð560Þ

The Mixture Mixing Rule To use Eq. (5.80) with mixtures, mixing rules are required to determine the terms (aα) and b for the mixtures. Soave adopted the following mixing rules: ðaαÞm ¼

XX pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  xi xj ai aj αi αj 1  kij i

(5.83)

j

bm ¼

X ½xi bi 

(5.84)

i

with A¼

ðaαÞm p ðRT Þ2

(5.85)

and B¼

bm p RT

(5.86)

The parameter kij is an empirically determined correction factor called the binary interaction coefficient (BIC), which is included to characterize any binary system formed by components i and j in the hydrocarbon mixture. These binary interaction coefficients are used to model the intermolecular interaction through empirical adjustment of the (aα)m term as represented mathematically by Eq. (5.83). They are dependent on the difference in molecular size of components in a binary system and they are characterized by the following properties: n

The interaction between hydrocarbon components increases as the relative difference between their molecular weights increases: ki, j + 1 > ki, j

n

Hydrocarbon components with the same molecular weight have a binary interaction coefficient of zero: ki, j ¼ 0

522 CHAPTER 5 Equations of State and Phase Equilibria

n

The binary interaction coefficient matrix is symmetric: ki, j ¼ kj, i

Slot-Petersen (1987) and Vidal and Daubert (1978) presented theoretical background to the meaning of the interaction coefficient and techniques for determining their values. Graboski and Daubert (1978) and Soave (1972) suggested that no binary interaction coefficients are required for hydrocarbon-hydrocarbon pairs (eg, kij ¼ 0), except between C1 and C7+. Whiston and Brule (2000) proposed the following expression for calculating the binary interaction coefficient between methane and the heavy fractions: kC1 C7 + ¼ 0:18  h

16:668 vci 1:1311 + ðvci Þ1=3

i6

where vci is the heavy fraction critical volume, which can be estimated from vci ¼ 0:4804 + 0:06011 Mi + 0:00001076 ðMi Þ2

where vci ¼ critical volume of the C7+, or its lumped component, ft3/lbm Mi ¼ molecular weight If nonhydrocarbons are present in the mixture, including binary interaction parameters can greatly improve the volumetric and phase behavior predictions of the mixture by the SRK EOS. Reid et al. (1987) proposed the following tabulated binary interaction coefficients for applications with the SRK EOS: Component

N2

CO2

H2S

N2 CO2 H2S C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7 +

0 0 0.12 0.02 0.06 0.08 0.08 0.08 0.08 0.08 0.08 0.08

0 0 0.12 0.12 0.15 0.15 0.15 0.15 0.15 0.15 0.15 0.15

0 0 0 0.08 0.07 0.07 0.06 0.06 0.06 0.06 0.05 0.03

In applying Eq. (5.75) for the compressibility factor of the liquid phase “ZL” the composition of the liquid “xi” is used to calculate the coefficients A and B

Equations of State 523

of Eqs. (5.85), (5.86) through the use of the mixing rules as described by Eqs. (5.83), (5.84). For determining the compressibility factor of the gas phase, Zv, the previously outlined procedure is used with composition of the gas phase, yi, replacing xi. EXAMPLE 5.12 A two-phase hydrocarbon system exists in equilibrium at 4000 psia and 160°F. The system has the composition shown in the following table. Component

xi

yi

pc

Tc

ωi

C1 C2 C3 C4 C5 C6 C7+

0.45 0.05 0.05 0.03 0.01 0.01 0.40

0.86 0.05 0.05 0.02 0.01 0.005 0.0005

666.4 706.5 616.0 527.9 488.6 453 285

343.33 549.92 666.06 765.62 845.8 923 1160

0.0104 0.0979 0.1522 0.1852 0.2280 0.2500 0.5200

The heptanes-plus fraction has the following properties: m ¼ 215 Pc ¼ 285 psia Tc ¼ 700°F ω ¼ 0:52 Assuming kij ¼ 0, calculate the density of each phase by using the SRK EOS. Solution Step 1. Calculate the parameters α, a, and b by applying Eqs. (5.71), (5.73), (5.74):   pffiffiffiffiffi 2 α ¼ 1 + m 1  Tr R2 Tc2 pc RTc b ¼ 0:08664 pc a ¼ 0:42747

The results are shown in the table below. Component

pc

Tc

ωi

αi

ai

bi

C1 C2 C3 C4 C5 C6 C7+

666.4 706.5 616.0 527.9 488.6 453 285

343.33 549.92 666.06 765.62 845.8 923 1160

0.0104 0.0979 0.1522 0.1852 0.2280 0.2500 0.5200

0.6869 0.9248 1.0502 1.1616 1.2639 1.3547 1.7859

8689.3 21,040.8 35,422.1 52,390.3 72,041.7 94,108.4 232,367.9

0.4780 0.7225 1.0046 1.2925 1.6091 1.9455 3.7838

524 CHAPTER 5 Equations of State and Phase Equilibria

Step 2. Calculate the mixture parameters (aα)m and bm for the gas phase and liquid phase by applying Eqs. (5.83), (5.84), to give the following. For the gas phase, using yi: XX pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  yi yj ai aj αi αj 1  kij ¼ 9219:3 ðaαÞm ¼ X i j bm ¼ ½yi bi  ¼ 0:5680 i

For the liquid phase, using xi: XX pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  xi xj ai aj αi αj 1  kij ¼ 104,362:9 ðaαÞm ¼ X i j bm ¼ ½xi bi  ¼ 1:8893 i

Step 3. Calculate the coefficients A and B for each phase by applying Eqs. (5.85), (5.86), to yield the following. For the gas phase, A¼

ðaαÞm p ð9219:3Þð4000Þ ¼ ¼ 0:8332 R2 T 2 ð10:73Þ2 ð620Þ2

B¼

bm p ð0:5680Þð4000Þ ¼ ¼ 0:3415 RT ð10:73Þð620Þ

For the liquid phase, A¼

ðaαÞm p ð104, 362:9Þð4000Þ ¼ ¼ 9:4324 R2 T 2 ð10:73Þ2 ð620Þ2

B¼

bm p ð1:8893Þð4000Þ ¼ ¼ 1:136 RT ð10:73Þð620Þ

Step 4. Solve Eq. (5.80) for the compressibility factor of the gas phase, to produce   Z3  Z2 + A  B  B2 Z + AB ¼ 0 Z 3  Z2 + ð0:8332  0:3415  0:34152ÞZ + ð0:8332Þð0:3415Þ ¼ 0 Solving this polynomial for the largest root gives Zv ¼ 0:9267 Step 5. Solve Eq. (5.80) for the compressibility factor of the liquid phase, to produce   Z3  Z 2 + A  B  B2 Z + AB ¼ 0 Z 3  Z 2 + ð9:4324  1:136  1:1362ÞZ + ð9:4324Þð1:136Þ ¼ 0 Solving the polynomial for the smallest root gives ZL ¼ 1:4121 Step 6. Calculate the apparent molecular weight of the gas phase and liquid phase from their composition, to yield the following:

Equations of State 525

For the gas phase: Ma ¼

X

For the liquid phase: Ma ¼

yi Mi ¼ 20:89

X xi Mi ¼ 100:25

Step 7. Calculate the density of each phase. ρ¼

pMa RTZ

For the gas phase: ρv ¼

ð4000Þð20:89Þ ¼ 13:556 lb=ft3 ð10:73Þð620Þð0:9267Þ

For the liquid phase: ρL ¼

ð4000Þð100:25Þ ¼ 42:68 lb=ft3 ð10:73Þð620Þð1:4121Þ

It is appropriate at this time to introduce and define the concept of the fugacity and the fugacity coefficient of the component. The fugacity “f” is a measure of the molar Gibbs energy of a real gas. It is evident from the definition that the fugacity has the units of pressure; in fact, the fugacity may be looked on as a vapor pressure modified to correctly represent the tendency of the molecules from one phase to escape into the other. In a mathematical form, the fugacity of a pure component is defined by the following expression: f o ¼ p exp

Z p    Z1 dp p o

(5.87)

where fo ¼ fugacity of a pure component, psia p ¼ pressure, psia Z ¼ compressibility factor The ratio of the fugacity to the pressure, f/p, is called the fugacity coefficient “Φ”. This dimensionless property is calculated from Eq. (5.87) as fo ¼ Φ ¼ exp p

Z p    Z1 dp p o

Soave applied this generalized thermodynamic relationship to Eq. (5.70) to determine the fugacity coefficient of a pure component, to give  o   f A Z+B ¼ ln ðΦÞ ¼ Z  1  ln ðZ  BÞ  ln ln p B Z

(5.88)

For a hydrocarbon multicomponent mixture that exits in the two-phase region, the component fugacity in each phase is introduced to develop a criterion for thermodynamic equilibrium. Physically, the fugacity of a

526 CHAPTER 5 Equations of State and Phase Equilibria

component i in one phase with respect to the fugacity of the component in a second phase is a measure of the potential for transfer of the component between phases. The phase with the lower component fugacity accepts the component from the phase with a higher component fugacity. Equal fugacities of a component in the two phases results in a zero net transfer. A zero transfer for all components implies that the hydrocarbon system is in thermodynamic equilibrium. Therefore, the condition of the thermodynamic equilibrium can be expressed mathematically by fiv ¼ fiL , 1  i  n

(5.89)

where f vi ¼ fugacity of component i in the gas phase, psi f Li ¼ fugacity of component i in the liquid phase, psi n ¼ number of components in the system The fugacity coefficient of component i in a hydrocarbon liquid mixture or hydrocarbon gas mixture is a function of the system pressure, mole fraction, and fugacity of the component. The fugacity coefficient is defined as For a component i in the gas phase: Φvi ¼

fiv yi p

For a component i in the liquid phase: ΦLi ¼

fiL xi p

(5.90)

(5.91)

where Φvi ¼ fugacity coefficient of component i in the vapor phase ΦLi ¼ fugacity coefficient of component i in the liquid phase.

Equilibrium Ratio “Ki” The equilibrium ratio is defined mathematically by Eq. (5.1) as Ki ¼ yi =xi

When two phases exist in equilibrium, it suggests that fiL ¼ fiv , which leads to redefining the K-value in terms of the fugacity of components as  L fi =ðxi pÞ ΦL ¼ iv Ki ¼  v Φi f i ð y i pÞ

(5.92)

Reid et al. (1987) defined the fugacity coefficient of component i in a hydrocarbon mixture by the following generalized thermodynamic relationship: 

1 lnðΦi Þ ¼ RT

Z V

1

@p RT  @ni V

 dV  ln ðZ Þ

(5.93)

Equations of State 527

where V ¼ total volume of n moles of the mixture ni ¼ number of moles of component i Z ¼ compressibility factor of the hydrocarbon mixture By combining the above thermodynamic definition of the fugacity with the SRK EOS [Eq. (5.70)], Soave proposed the following expression for the fugacity coefficient of component i in the liquid phase:  ln

          bi ðZL  1Þ fiL A 2Ψ i bi B ¼ ln ΦLi ln 1 + L  ln ZL  B   xi p bm B ðaαÞm bm Z (5.94)

where Ψi ¼

X pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  xj ai aj αi αj 1  kij

(5.95)

j

ðaαÞm ¼

XX pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  xi xj ai aj αi αj 1  kij i

(5.96)

j

Eq. (5.94) is also used to determine the fugacity coefficient of component in the gas phase, Φvi , by using the composition of the gas phase, yi, in calculating A, B, Zv, and other composition-dependent terms; or,        bi ðZv  1Þ A 2Ψ i bi B ln 1 + v ln Φvi ¼  lnðZv  BÞ   bm B ðaαÞm bm Z

where X pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  yj ai aj αi αj 1  kij j XX  pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  ðaαÞm ¼ yi yj ai aj αi αj 1  kij

Ψi ¼

i

j

Gibbs Energy and Selecting the Correct Z-Factor In solving the cubic Z-factor equation for hydrocarbon mixtures, one or three real roots may exist. When three roots exist with different values designated as n n n

Largest root: ZLargest Middle root: Zmiddle Smallest root: ZSmallest,

The middle root ZMiddle, always is discarded as it has nonphysical value or meaning. However, selecting of the correct root between the remaining two

528 CHAPTER 5 Equations of State and Phase Equilibria

roots is important to avoid phase equilibria convergence problems, which will lead the EOS to generate unrealistic results. The correct root is chosen as the one with the lowest normalized Gibbs energy function, g*. The normalized Gibbs energy function is defined in terms of the composition of the hydrocarbon phase and the individual component fugacity in that phase. Mathematically, the normalized Gibbs energy function is defined by the following expressions: Normalized Gibbs energy for gas phase: ggas ¼

n X   yi ln fiv i¼1

Normalized Gibbs energy for liquid phase:

gLiquid

¼

n X   xi ln fiL i¼1

Therefore, assume that a liquid phase with a composition of xi has three Z-factor roots; the middle root is discarded automatically and the remaining two are designated ZL1 and ZL2 . To select the correct root, calculate the normalized Gibbs energy function using the two remaining roots:

1

n X   xi ln fiL

2

n X   ¼ xi ln fiL

gZL ¼ gZL

i¼1

i¼1

in which the fugacity fLi is determined from Eq. (5.94) using the two remaining roots, ZL1 and ZL2 , as

      L bi ðZL1  1Þ A 2Ψ i bi B ln 1 + fi zL1 ¼ xi p + exp  lnðZL1  BÞ   bm B ðaαÞm bm ZL1  

    L bi ðZL2  1Þ A 2Ψ i bi B ln 1 +  lnðZL2  BÞ   fi zL2 ¼ xi + exp bm B ðaαÞm bm ZL2

If gZL1
n X   yi ln fiv i¼1

gZg2

n X   ¼ yi ln fiv i¼1

Equations of State 529

in which the fugacity f vi is determined by using the two remaining roots, Zg1 and Zg2, as 

       v   bi Zg1  1 A 2Ψ i bi B ln 1 + fi zg1 ¼ xi p + exp  ln Zg1  B   bm B ðaαÞm bm Zg1 

       v   bi Zg2  1 A 2Ψ i bi B ln 1 + fi zg2 ¼ xi + exp  ln Zg2  B   bm B ðaαÞm bm Zg2

If gZg1 < gZg2 then Zg1 is chosen as the correct root; otherwise, Zg2 is chosen.

Chemical Potential The chemical potential “μi” is a quantity that measures how much energy a component brings to a mixture as the component transfers from one phase to another. In a mixture, the chemical potential is defined as the slope of the free energy “G” of the system with respect to a change in the number of moles of that specific component. Thus, it is the partial derivative of the free energy with respect to the amount of that specific component, which is given the name chemical potential: 

@G μi ¼ @ni

 p, T , nj

In a mixture, the normalized Gibbs free energy “g” is not the same in both phases because their compositions are different. Instead, the criterion for equilibrium is defined mathematically by i @G h gas ¼0 ¼ μi  μliquid i @ni

This criterion implies that, at equilibrium, the chemical potential of each component μi has the same value in each phase: liquid μgas i ¼ μi

The chemical potential of component i is related to the activity a^i by μi ¼ μoi + RT lnða^i Þ

with a^i ¼

fi fio

530 CHAPTER 5 Equations of State and Phase Equilibria

where a^i ¼ activity of component i in mixture at p and T μoi ¼ chemical potential of pure of component i at p and T f oi ¼ fugacity of pure component i at p and T fi ¼ fugacity of component i in mixture at p and T

MODIFICATIONS OF THE SRK EOS To improve the pure component vapor pressure predictions by the SRK equation of state, Groboski and Daubert (1978) proposed a new expression for calculating parameter m of Eq. (5.72). The proposed relationship, which originated from analyzing extensive experimental data for pure hydrocarbons, has the following form: m ¼ 0:48508 + 1:55171ω  0:15613ω2

(5.97)

Sim and Daubert (1980) pointed out that, because the coefficients of Eq. (5.97) were determined by analyzing vapor pressure data of lowmolecular-weight hydrocarbons, it is unlikely that Eq. (5.97) will suffice for high-molecular-weight petroleum fractions. Realizing that the acentric factors for the heavy petroleum fractions are calculated from an equation such as the Edmister correlation or the Lee and Kesler correlation, the authors proposed the following expressions for determining parameter m: n

If the acentric factor is determined using the Edmister correlation, then m ¼ 0:431 + 1:57ωi  0:161ðωÞ2i

n

(5.98)

If the acentric factor is determined using the Lee and Kesler correction, then m ¼ 0:315 + 1:60ωi  0:166ω2i

(5.99)

Elliot and Daubert (1985) stated that the optimal binary interaction coefficient, kij, would minimize the error in the representation of all the thermodynamic properties of a mixture. Properties of particular interest in phase equilibrium calculations include bubble point pressure, dew point pressure, and equilibrium ratios. The authors proposed a set of relationships for determining interaction coefficients for asymmetric mixtures2 that contain 2

An asymmetric mixture is defined as one in which two of the components are considerably different in their chemical behavior. Mixtures of methane with hydrocarbons of 10 or more carbon atoms can be considered asymmetric. Mixtures containing gases such as nitrogen or hydrogen are asymmetric.

Modifications of the SRK EOS 531

methane, N2, CO2, and H2S. Referring to the principal component as i and other fractions as j, Elliot and Daubert proposed the following expressions: Binary Interaction between N2 – Components kij ¼ 0:107089 + 2:9776kij1

(5.100)

Binary Interaction between CO2 – Components

kij ¼ 0:08058  0:77215kij1  1:8404 kij1

(5.101)

Binary Interaction between H2S – Components kij ¼ 0:07654 + 0:017921kij1

(5.102)

Methane with components > C9 (Nonane) 2 kij ¼ 0:17985 + 2:6958kij1 + 10:853 kij1

(5.103)

where kij1 ¼

 2  εi  εj 2εi εj

(5.104)

and εi ¼

pffiffiffiffi 0:480453 ai bi

(5.105)

The two parameters, ai and bi, of Eq. (5.105) were defined previously by Eqs. (5.73), (5.74); that is, a ¼ Ωa

R2 Tc2 pc

b ¼ Ωb

RTc pc

Volume Translation (Volume Shift) Parameter The major drawback in the SRK EOS is that the critical compressibility factor takes on the unrealistic universal critical compressibility of 0.333 for all substances. Consequently, the molar volumes typically are overestimated; that is, densities are underestimated. Peneloux et al. (1982) developed a procedure for improving the volumetric predictions of the SRK EOS by introducing a volume correction parameter, ci, into the equation. This third parameter does not change the vapor/liquid

532 CHAPTER 5 Equations of State and Phase Equilibria

equilibrium conditions determined by the unmodified SRK equation (ie, the equilibrium ratio Ki) but it modifies the liquid and gas volumes. The proposed methodology, known as the volume translation method, uses the following expressions: L Vcorr ¼ VL 

X ðxi ci Þ

(5.106)

i v Vcorr ¼ Vv 

X ðy i c i Þ

(5.107)

i

where VL ¼ uncorrected liquid molar volume (ie, VL ¼ ZLRT/p), ft3/mol Vv ¼ uncorrected gas molar volume Vv ¼ ZvRT/p, ft3/mol VLcorr ¼ corrected liquid molar volume, ft3/mol Vvcorr ¼ corrected gas molar volume, ft3/mol xi ¼ mole fraction of component i in the liquid phase yi ¼ mole fraction of component i in the gas phase The authors proposed six schemes for calculating the correction factor, ci, for each component. For petroleum fluids and heavy hydrocarbons, Peneloux and coworkers suggested that the best correlating parameter for the volume correction factor ci is the Rackett compressibility factor, ZRA. The correction factor then is defined mathematically by the following relationship: ci ¼ 4:437978ð0:29441  ZRA ÞTci =pci

(5.108)

where ci ¼ volume shift coefficient for component i, ft3/lb-mol Tci ¼ critical temperature of component i, °R pci ¼ critical pressure of component i, psia The parameter ZRA is a unique constant for each compound. The values of ZRA, in general, are not much different from those of the critical compressibility factors, Zc. If their values are not available, Peneloux et al. proposed the following correlation for calculating ci:  ci ¼ ð0:0115831168 + 0:411844152ωi Þ

where ωi ¼ acentric factor of component i.

Tci pci

 (5.109)

Modifications of the SRK EOS 533

EXAMPLE 5.13 Rework Example 5.12 by incorporating the Peneloux et al. volume correction approach in the solution. Key information from Example 5.12 includes: For gas: Zv ¼ 0:9267, Ma ¼ 20:89 For liquid: Z L ¼ 1:4212, Ma ¼ 100:25 T ¼ 160°F, and p ¼ 4000psi Solution Step 1. Calculate the correction factor ci using Eq. (5.109):   Tci ci ¼ ð0:0115831168 + 0:411844152ωi Þ pci The results are shown in the table below. Component xi

pc

C1 C2 C3 C4 C5 C6 C P7+

666.4 343.33 706.5 549.92 616.0 666.06 527.9 765.62 488.6 845.8 453 923 285 1160

0.45 0.05 0.05 0.03 0.01 0.01 0.40

Tc

ωi

ci

cixi

yi

ciyi

0.0104 0.0979 0.1522 0.1852 0.2280 0.2500 0.5200

0.00839 0.03807 0.07729 0.1265 0.19897 0.2791 0.91881

0.003776 0.001903 0.003861 0.00379 0.001989 0.00279 0.36752 0.38564

0.86 0.05 0.05 0.02 0.01 0.005 0.005

0.00722 0.00190 0.00386 0.00253 0.00198 0.00139 0.00459 0.02349

Step 2. Calculate the uncorrected volume of the gas and liquid phase using the compressibility factors as calculated in Example 5.12. Vv ¼

RTZv ð10:73Þð620Þð0:9267Þ ¼ ¼ 1:54119ft3 =mol p 4000

VL ¼

RTZL ð10:73Þð620Þð1:4121Þ ¼ ¼ 2:3485ft3 =mol p 4000

Step 3. Calculate the corrected gas and liquid volumes by applying Eqs. (5.106), (5.107). X L ¼ VL  ðxi ci Þ ¼ 2:3485  0:38564 ¼ 1:962927ft3 =mol Vcorr i

v Vcorr

X ¼V  ðyi ci Þ ¼ 1:54119  0:02394 ¼ 1:5177 ft3 =mol v

i

Step 4. Calculate the corrected compressibility factors. v Zcorr ¼

ð4000Þð1:5177Þ ¼ 0:91254 ð10:73Þð620Þ

L ¼ Zcorr

ð4000Þð1:962927Þ ¼ 1:18025 ð10:73Þð620Þ

534 CHAPTER 5 Equations of State and Phase Equilibria

Step 5. Determine the corrected densities of both phases. ρ¼

pMa RTZ

ð4000Þð20:89Þ ¼ 13:767lb=ft3 ð10:73Þð620Þð0:91254Þ ð4000Þð100:25Þ ρL ¼ ¼ 51:07lb=ft3 ð10:73Þð620Þð1:18025Þ ρv ¼

As pointed out by Whitson and Brule (2000), when the volume shift is introduced to the EOS for mixtures, the resulting expressions for fugacity are h  L   pi fi modified ¼ fiL original exp ci RT

and h  v   pi fi modified ¼ fiv original exp ci RT

This implies that the fugacity ratios are unchanged by the volume shift:  L  L fi f  iv  modified ¼  v  original fi modified fi modified

The Peng-Robinson Equation of State Peng and Robinson (1976a) conducted a comprehensive study to evaluate the use of the SRK equation of state for predicting the behavior of naturally occurring hydrocarbon systems. They illustrated the need for an improvement in the ability of the equation of state to predict liquid densities and other fluid properties, particularly in the vicinity of the critical region. As a basis for creating an improved model, Peng and Robinson (PR) proposed the following expression: p¼

RT aα  V  b ðV + bÞ2  cb2

where a, b, and α have the same significance as they have in the SRK model, and parameter c is a whole number optimized by analyzing the values of the terms Zc and b/Vc as obtained from the equation. It is generally accepted that Zc should be close to 0.28 and b/Vc should be approximately 0.26. An optimized value of c ¼ 2 gives Zc ¼ 0.307 and (b/Vc) ¼ 0.253. Based on this value

Modifications of the SRK EOS 535

of c, Peng and Robinson proposed the following equation of state (commonly referred to as PR EOS): p¼

RT aα  V  b V ðV + bÞ + bðV  bÞ

(5.110)

Imposing the classical critical point conditions [Eq. (5.48)] on Eq. (5.110) and solving for parameters a and b yields a ¼ Ωa

R2 Tc2 pc

(5.111)

RTc pc

(5.112)

b ¼ Ωb

where Ωa ¼ 0:45724 Ωb ¼ 0:07780

This equation predicts a universal critical gas compressibility factor, Zc, of 0.307, compared to 0.333 for the SRK model. Peng and Robinson also adopted Soave’s approach for calculating parameter α:   pffiffiffiffiffi 2 α ¼ 1 + m 1  Tr

(5.113)

where m ¼ 0.3796 + 1.54226ω  0.2699ω2. Peng and Robinson (1978) proposed the following modified expression for m that is recommended for heavier components with acentric values ω > 0.49: m ¼ 0:379642 + 1:48503ω  0:1644ω2 + 0:016667ω3

(5.114)

Rearranging Eq. (5.110) into the compressibility factor form, gives     Z3 + ðB  1ÞZ 2 + A  3B2  2B Z  AB  B2  B3 ¼ 0

(5.115)

where A and B are given by Eqs. (5.81), (5.82) for pure component and by Eqs. (5.85), (5.86) for mixtures; that is, A¼ B¼

with

ðaαÞm p ðRT Þ2 bm p RT

XX pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  xi xj ai aj αi αj 1  kij i j X bm ¼ ½xi bi 

ðaαÞm ¼

i

536 CHAPTER 5 Equations of State and Phase Equilibria

EXAMPLE 5.14 Using the composition given in Example 5.12, calculate the density of the gas phase and liquid phase by using the Peng-Robinson EOS. Assume kij ¼ 0. The components are described in the following table. Component

xi

yi

pc

Tc

ωi

C1 C2 C3 C4 C5 C6 C7+

0.45 0.05 0.05 0.03 0.01 0.01 0.40

0.86 0.05 0.05 0.02 0.01 0.005 0.0005

666.4 706.5 616.0 527.9 488.6 453 285

343.33 549.92 666.06 765.62 845.8 923 1160

0.0104 0.0979 0.1522 0.1852 0.2280 0.2500 0.5200

The heptanes-plus fraction has the following properties: M ¼ 215 pc ¼ 285 psia Tc ¼ 700°F ω ¼ 0:582 Solution Step 1. Calculate the parameters a, b, and α for each component in the system by applying Eqs. (5.111)–(5.113): R2 Tc2 pc RTc b ¼ 0:07780 pc   pffiffiffiffiffi 2 α ¼ 1 + m 1  Tr a ¼ 0:45724

The results are shown in the table below. Component pc C1 C2 C3 C4 C5 C6 C7+

Tc

666.4 343.33 706.5 549.92 616.0 666.06 527.9 765.62 488.6 845.8 453 923 285 1160

ωi

m

αi

ai

0.0104 0.0979 0.1522 0.1852 0.2280 0.2500 0.5200

1.8058 1.1274 0.9308 0.80980 0.73303 0.67172 0.53448

0.7465 9294.4 0.9358 22,506.1 1.0433 37,889.0 1.1357 56,038.9 1.2170 77,058.9 1.2882 100,662.3 1.6851 248,550.5

bi 0.4347 0.6571 0.9136 1.1755 1.6091 1.4635 3.4414

Step 2. Calculate the mixture parameters (aα)m and bm for the gas and liquid phase, to give the following.

Modifications of the SRK EOS 537

For the gas phase, XX pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  yi yj ai aj αi αj 1  kij ¼ 10,423:54 ðaαÞm ¼ i j X bm ¼ ðyi bi Þ ¼ 0:862528 i

For the liquid phase, XX pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  xi xj ai aj αi αj 1  kij ¼ 107,325:4 ðaαÞm ¼ i j X bm ¼ ðxi bi Þ ¼ 1:696543 i

Step 3. Calculate the coefficients A and B, to give the following. For the gas phase, A¼

ðaαÞm p ð10, 423:54Þð4000Þ ¼ ¼ 0:94209 R2 T 2 ð10:73Þ2 ð620Þ2

B¼

bm p ð0:862528Þð4000Þ ¼ ¼ 0:30669 RT ð10:73Þð620Þ

For the liquid phase, A¼

ðaαÞm p ð107, 325:4Þð4000Þ ¼ ¼ 9:700183 R2 T 2 ð10:73Þ2 ð620Þ2

B¼

bm p ð1:696543Þð4000Þ ¼ ¼ 1:020078 RT ð10:73Þð620Þ

Step 4. Solve Eq. (5.115) for the compressibility factor of the gas and liquid phases to give     Z3 + ðB  1ÞZ 2 + A  3B2  2B Z  AB  B2  B3 ¼ 0 For the gas phase: Substituting for A ¼ 0.94209 and B ¼ 0.30669 in the above equation gives     Z 3 + ðB  1ÞZ2 + A  3B2  2B2  2B Z  AB  B2  B3 ¼ 0 h i Z 3 ¼ ð0:30669  1ÞZ 2 + 0:94209  3ð0:30669Þ2  2ð0:30669Þ Z h i  ð0:94209Þð0:30669Þ  ð0:30669Þ2  ð0:30669Þ3 ¼ 0 Solving for Z gives Z v ¼ 0:8625 For the liquid phase: Substituting A ¼ 9.700183 and B ¼ 1.020078 in Eq. (5.115) gives     Z + ðB  1ÞZ 2 + A  3B2  2B Z  AB  B2  B3 ¼ 0 h i Z3 ¼ ð1:020078  1ÞZ2 + ð9:700183Þ  3ð1:020078Þ2  2ð1:020078Þ Z h i  ð9:700183Þð1:020078Þ  ð1:020078Þ2  ð1:020078Þ3 ¼ 0 3

538 CHAPTER 5 Equations of State and Phase Equilibria

Solving for Z gives ZL ¼ 1:2645 Step 5. Calculate the density of both phases. Density of the vapor phase: pMa Zv RT ð4000Þð20:89Þ ρv ¼ ¼ 14:566lb=ft3 ð10:73Þð620Þð0:8625Þ ρv ¼

Density of the liquid phase: pMa Zv RT ð4000Þð100:25Þ ¼ 47:67lb=ft3 ρL ¼ ð10:73Þð620Þð1:2645Þ

ρL ¼

Applying the thermodynamic relationship, as given by Eq. (5.88), to Eq. (5.111) yields the following expression for the fugacity of a pure component:  pffiffiffi #   "   Z+ 1+ 2 B f A pffiffiffi  ln ¼ lnðΦÞ ¼ Z  1  lnðZ  BÞ  pffiffiffiffiffiffi ln p Z 1 2 B 2 2B

(5.116)

The fugacity coefficient of component i in a hydrocarbon liquid mixture is calculated from the following expression:  L   bi ðZL  1Þ   f ln  ln ZL  B ¼ ln ΦLi ¼ xi p bm   " L  pffiffiffi  # Z + 1 + 2B A 2Ψ i bi pffiffiffi   ln   pffiffiffi ð aα Þ b Z L + 1  2B 2 2B m m

(5.117)

where the mixture parameters bm, B, A, Ψ i, and (aα)m are defined previously. Eq. (5.117) also is used to determine the fugacity coefficient of any component in the gas phase by replacing the composition of the liquid phase, xi, with the composition of the gas phase, yi, in calculating the compositiondependent terms of the equation:      bi ðZ v  1Þ fv A 2Ψ i bi  lnðZ v  BÞ  pffiffiffi  ¼ ln ΦLi ¼ yi p bm 2 2B ðaαÞm bm "  pffiffiffi  # v Z + 1 + 2B pffiffiffi  ln  Z v + 1  2B

 ln

Modifications of the SRK EOS 539

The table below shows the set of binary interaction coefficient, kij, traditionally used when predicting the volumetric behavior of hydrocarbon mixture with the Peng-Robinson equation of state.

Component CO2 N2 H2S

C1

C2

C3

i-C4

n-C4

i-C5

n-C5

C6

C7

C8

C9

C10

CO2 N2 H2S C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7 C8 C9 C10

0.105 0.025 0.070 0

0.130 0.010 0.085 0.005 0

0.125 0.090 0.080 0.010 0.005 0

0.120 0.095 0.075 0.035 0.005 0.000 0

0.115 0.095 0.075 0.025 0.010 0.000 0.005 0

0.115 0.100 0.070 0.050 0.020 0.015 0.005 0.005 0

0.115 0.100 0.070 0.030 0.020 0.015 0.005 0.005 0.000 0

0.115 0.110 0.070 0.030 0.020 0.010 0.005 0.005 0.000 0.000 0

0.115 0.115 0.060 0.035 0.020 0.005 0.005 0.005 0.000 0.000 0.000 0

0.115 0.120 0.060 0.040 0.020 0.005 0.005 0.005 0.000 0.000 0.000 0.000 0

0.115 0.120 0.060 0.040 0.020 0.005 0.005 0.005 0.000 0.000 0.000 0.000 0.000 0

0.115 0.125 0.055 0.045 0.020 0.005 0.005 0.005 0.000 0.000 0.000 0.000 0.000 0.00 0

0

0 0

0.135 0.130 0

Note that kij ¼ kij.

To improve the predictive capability of the PR EOS when applying for mixtures containing nonhydrocarbon components (eg, N2, CO2, and CH4), Nikos et al. (1986) proposed a generalized correlation for generating the binary interaction coefficient, kij. The authors correlated these coefficients with system pressure, temperature, and the acentric factor. These generalized correlations originated with all the binary experimental data available in the literature. The authors proposed the following generalized form for kij: kij ¼ λ2 Trj2 + λ1 Trj + λ0

(5.118)

where i refers to the principal component, N2, CO2, or CH4, and j refers to the other hydrocarbon component of the binary. The acentric-factordependent coefficients, λ0, λ1, and λ2 are determined for each set of binaries by applying the following expressions: For nitrogen-hydrocarbons,      2 λ0 ¼ 0:1751787  0:7043 log ωj  0:862066 log ωj

(5.119)

     2 λ1 ¼ 0:584474  1:328 log ωj + 2:035767 log ωj

(5.120)

540 CHAPTER 5 Equations of State and Phase Equilibria      2 λ2 ¼ 2:257079 + 7:869765 log ωj + 13:50466 log ωj + 8:3864½ log ðωÞ3 (5.121)

Nikos et al. also suggested the following pressure correction:   kij0 ¼ kij 1:04  4:2  105 p

(5.122)

where p is the pressure in psi. For methane-hydrocarbons,      2 λ0 ¼ 0:01664  0:37283 log ωj + 1:31757 log ωj

(5.123)

λ1 ¼ 0:48147 + 3:35342 log ðωi Þ  1:0783½ log ðωi Þ2

(5.124)

λ2 ¼ 0:4114  3:5072 log ðωi Þ  0:78798½ log ðωi Þ2

(5.125)

For CO2-hydrocarbons,   λ0 ¼ 0:4025636 + 0:1748927 log ωj   λ1 ¼ 0:94812  0:6009864 log ωj

(5.127)

λ2 ¼ 0:741843368 + 0:441775 log ðωi Þ

(5.128)

(5.126)

For the CO2 interaction parameters, the following pressure correction is suggested:   kij0 ¼ kij 1:044269  4:375  105 p

(5.129)

Stryjek and Vera (1986) proposed an improvement in the reproduction of vapor pressures of a pure component by the PR EOS in the reduced temperature range from 0.7 to 1.0, by replacing the m term in Eq. (5.113) with the following expression: mo ¼ 0:378893 + 1:4897153  0:17131848ω2 + 0:0196554ω3

(5.130)

To reproduce vapor pressures at reduced temperatures below 0.7, Stryjek and Vera further modified the m parameter in the PR equation by introducing an adjustable parameter, m1, characteristic of each compound to Eq. (5.113). They proposed the following generalized relationship for the parameter m:   pffiffiffiffiffi m ¼ mo + m1 1 + Tr ð0:7  Tr Þ

where Tr ¼ reduced temperature of the pure component m0 ¼ defined by Eq. (5.130) m1 ¼ adjustable parameter

(5.131)

Modifications of the SRK EOS 541

For all components with a reduced temperature above 0.7, Stryjek and Vera recommended setting m1 ¼ 0. For components with a reduced temperature greater than 0.7, the optimum values of m1 for compounds of industrial interest are tabulated below, Component

m1

Nitrogen Carbon dioxide Water Methane Ethane Propane Butane Pentane Hexane Heptane Octane Nonane Decane Undecane Dodecane Tridecane Tetradecane Pentadecane Hexadecane

0.01996 0.04285 0.0664 0.0016 0.02669 0.03136 0.03443 0.03946 0.05104 0.04648 0.04464 0.04104 0.0451 0.02919 0.05426 0.04157 0.02686 0.01892 0.02665

Due to the totally empirical nature of the parameter m1, Stryjek and Vera could not find a generalized correlation for m1 in terms of pure component parameters. They pointed out that these values of m1 should be used without changes. Jhaveri and Youngren (1984) pointed out that, when applying the PengRobinson equation of state to reservoir fluids, the error associated with the equation in the prediction of gas-phase Z-factors ranged from 3% to 5% and the error in the liquid density predictions ranged from 6% to 12%. Following the procedure proposed by Peneloux and coworkers (see the SRK EOS), Jhaveri and Youngren introduced the volume correction parameter, ci, to the PR EOS. This third parameter has the same units as the second parameter, bi, of the unmodified PR equation and is defined by the following relationship: ci ¼ Si bi

(5.132)

where Si is a dimensionless parameter, called the shift parameter, and bi is the Peng-Robinson covolume, as given by Eq. (5.112).

542 CHAPTER 5 Equations of State and Phase Equilibria

As in the SRK EOS, the volume correction parameter, ci, does not change the vapor/liquid equilibrium conditions, equilibrium ratio Ki. The corrected hydrocarbon phase volumes are given by the following expressions: L Vcorr ¼ VL 

X ðxi ci Þ

v ¼ Vv  Vcorr

X ðy i c i Þ

i¼1

i¼1

where VL, Vv ¼ volumes of the liquid phase and gas phase as calculated by the unmodified PR EOS, ft3/mol VLcorr, Vvcorr ¼ corrected volumes of the liquid and gas phase, ft3/mol Whitson and Brule (2000) point out that the volume translation (correction) concept can be applied to any two-constant cubic equation, thereby eliminating the volumetric deficiency associated with the application of EOS. Whitson and Brule extended the work of Jhaveri and Youngren and tabulated the shift parameters, Si, for a selected number of pure components. These tabulated values, which follow, are used in Eq. (5.132) to calculate the volume correction parameter, ci, for the Peng-Robinson and SRK equations of state. Component

PR EOS

N2 CO2 H2S C1 C2 C3 i-C4 n-C4 i-C5 n-C5 n-C6 n-C7 n-C8 n-C9 n-C10

0.1927 0.0817 0.1288 0.1595 0.1134 0.0863 0.0844 0.0675 0.0608 0.0390 0.0080 0.0033 0.0314 0.0408 0.0655

SRK EOS 0.0079 0.0833 0.0466 0.0234 0.0605 0.0825 0.0830 0.0975 0.1022 0.1209 0.1467 0.1554 0.1794 0.1868 0.2080

Jhaveri and Youngren proposed the following expression for calculating the shift parameter for the C7+: SC7 + ¼ 1 

d ð M Þe

Modifications of the SRK EOS 543

where M ¼ molecular weight of the heptanes-plus fraction d, e ¼ positive correlation coefficients The authors proposed that, in the absence of the experimental information needed for calculating e and d, the power coefficient e could be set equal to 0.2051 and the coefficient d adjusted to match the C7+ density with the values of d ranging from 2.2 to 3.2. In general, the following values may be used for C7+ fractions, by hydrocarbon family: Hydrocarbon family

d

e

Paraffins Naphthenes Aromatics

2.258 3.044 2.516

0.1823 0.2324 0.2008

To use the Peng-Robinson equation of state to predict the phase and volumetric behavior of mixtures, one must be able to provide the critical pressure, the critical temperature, and the acentric factor for each component in the mixture. For pure compounds, the required properties are well defined and known. Nearly all naturally occurring petroleum fluids contain a quantity of heavy fractions that are not well defined. These heavy fractions often are lumped together as a heptanes-plus fraction. The problem of how to adequately characterize the C7+ fractions in terms of their critical properties and acentric factors has been long recognized in the petroleum industry. Changing the characterization of C7+ fractions present in even small amounts can have a profound effect on the PVT properties and the phase equilibria of a hydrocarbon system as predicted by the Peng-Robinson equation of state. The usual approach for such situations is to “tune” the parameters in the EOS in an attempt to improve the accuracy of prediction. During the tuning process, the critical properties of the heptanes-plus fraction and the binary interaction coefficients are adjusted to obtain a reasonable match with experimental data available on the hydrocarbon mixture. Recognizing that the inadequacy of the predictive capability of the PR EOS lies with the improper procedure for calculating the parameters a, b, and α of the equation for the C7+ fraction, Ahmed (1991) devised an approach for determining these parameters from the following two readily measured physical properties of C7+: the molecular weight, M7+, and the specific gravity, γ 7+. The approach is based on generating 49 density values for the C7+ by applying the Riazi and Daubert correlation. These values were subsequently subjected to 10 temperature and 10 pressure values in the range of 60–300°F and 14.7–7000 psia, respectively. The Peng-Robinson EOS was then

544 CHAPTER 5 Equations of State and Phase Equilibria

applied to match the 4900 generated density values by optimizing the parameters a, b, and α using a nonlinear regression model. The optimized parameters for the heptanes-plus fraction are given by the following expressions. For the parameter α of C7+, "

rffiffiffiffiffiffiffiffi!#2 520 α¼ 1+m 1 T

(5.133)

with m as defined by m¼

D A4 A7 + A2 M7 + + A3 M72 + + + A5 γ 7 + + A6 γ 27 + + M7 + γ7 + A0 + A1 D

(5.134)

with parameter D as defined by the ratio of the molecular weight to the specific gravity of the heptanes-plus fraction: D¼

M7 + γ7 +

where M7+ ¼ molecular weight of C7+ γ 7+ ¼ specific gravity of C7+ A0-A7 ¼ coefficients as given in the table below Coefficient a A0 A1 A2 A3 A4 A5 A6 A7

b

m

2.433525  10 6.8453198 36.91776 8.3201587  103 1.730243  102 5.2393763  102 2 6 0.18444102  10 6.2055064  10 1.7316235  102 2 9 3.6003101  10 9.0910383  10 1.3743308  105 7 3.4992796  10 13.378898 12.718844 2.838756  107 7.9492922 10.246122 1.1325365  107 3.1779077 7.6697942 6.418828  106 1.7190311 2.6078099 7

For parameters a and b of C7+, the following generalized correlation is proposed: " a or b ¼

3  X i¼0

Ai D

i



#

" # 6  X  A4 A7 i4 + + Ai γ 7 + + D γ 7+ i¼5

(5.135)

The coefficients A0–A7 are included in the table above. To further improve the predictive capability of the Peng-Robinson EOS, Ahmed (1991)

Modifications of the SRK EOS 545

optimized the coefficients a, b, and m for nitrogen, CO2, and methane by matching 100 Z-factor values for each of these components. Using a nonlinear regression model, the optimized values in the table below are recommended. Component

a

CO2 N2 C1

1.499914  10 4.5693589  103 7.709708  103 4

b

m (Eq. (5.133))

0.41503575 0.4682582 0.46749727

0.73605717 0.97962859 0.549765

To provide the modified PR EOS with a consistent procedure for determining the binary interaction coefficient, kij, the following computational steps are proposed. Step 1. Calculate the binary interaction coefficient between methane and the heptanes-plus fraction from kC1 C7 + ¼ 0:00189T  1:167059

where the temperature, T, is in °R. Step 2. Set kCO2 -N2 ¼ 0:12 KCO2 -hydrocarbon ¼ 0:10 kN2 -hydrocarbon ¼ 0:10

Step 3. Adopting the procedure recommended by Pedersen, Thomassen, and Fredenslund (1989), calculate the binary interaction coefficients between components heavier than methane (ie, C2, C3, etc.) and the heptanes-plus fraction from kCn C7 + ¼ 0:8kCðn1Þ C7 +

where n is the number of carbon atoms of component Cn; for example, Binary interaction coefficient between C2 and C7 + : kC2 C7 + ¼ 0:8kC1 C7 + Binary interaction coefficient between C3 and C7 + : kC3 C7 + ¼ 0:8kC2 C7 + Step 4. Determine the remaining kij from "  kij  kiC7 +

Mj

5

 ðM i Þ5

ðMC7 + Þ5  ðMi Þ5

#

546 CHAPTER 5 Equations of State and Phase Equilibria

where M is the molecular weight of any specified component. For example, the binary interaction coefficient between propane, C3, and butane, C4, is " kC3 C4 ¼ kC3 C7 +

ðM C 4 Þ5  ðM C 3 Þ5

#

ðMC7 + Þ5  ðMC3 Þ5

Summary of Equations of State A summary of the different equations of state discussed in this chapter are listed in the table below in the form of p ¼ prepulsion  pattraction EOS

prepulsion

pattraction

a

b

Ideal

RT V RT V b RT V b RT V b RT V b

0

0

0

a V2

R2 Tc2 pc R2 Tc2:5 Ωa pc R2 Tc2 Ωa pc R2 Tc2 Ωa pc

Ωb

vdW RK SRK PR

a pffiffiffi V ð V + bÞ T aαðT Þ V ð V + bÞ aαðT Þ V ðV + bÞ + bðV  bÞ

Ωa

RTc pc RTc Ωb pc RTc Ωb pc RTc Ωb pc

EQUATION-OF-STATE APPLICATIONS Most petroleum engineering applications rely on the use of the equation of state due to its simplicity, consistency, and accuracy when properly tuned. Cubic equations of state have found widespread acceptance as tools that permit the convenient and flexible calculation of the complex phase behavior of reservoir fluids. Some of these applications are presented in this section, including: n n n n n n

Equilibrium ratio “K-value” Dew point pressure “pd” Bubble point pressure “pb” Three-phase flash calculations Vapor pressure “pv” Simulating PVT laboratory experiments.

Equilibrium Ratio A flow diagram is presented in Fig. 5.9 to illustrate the procedure of determining equilibrium ratios of a hydrocarbon mixture. For this type of calculation, the system temperature, T, the system pressure, p, and the overall

Equation-of-State Applications 547

Given: zi, P, T

KA i =

Assume KAi Ki

yi, nV

ZV, FVi

No, set

Xi, nL

Perform flash calculations

Z L, F Li

Calculate Ki

e > 0.0001 KA i =

Tci pci )] Ki = p exp[5.37(1 + wi)(1− T

n

Ki

Test for convergence

e ≤ 0.0001

Ki

−1 ≤ e

KA i where e ≈ 0.0001 i=1

Yes

Solution gives: Ki, xi, yi, nL, nV, Z L, Z V

n FIGURE 5.9 Flow diagram of equilibrium ratio determination by EOS.

composition of the mixture, zi, must be known. The procedure is summarized in the following steps in conjunction with Fig. 5.9. Step 1. Assume an initial value of the equilibrium ratio for each component in the mixture at the specified system pressure and temperature. Wilson’s equation can provide the starting values of Ki: KiA ¼

pci exp½5:37ð1 + ωi Þð1  Tci =T Þ p

where KA i is an assumed equilibrium ratio for component i. Step 2. Using the overall composition and the assumed K-values, perform flash calculations to determine xi, yi, nL, and nυ. Step 3. Using the calculated composition of the liquid phase, xi, determine the fugacity coefficient ΦLi for each component in the liquid phase by applying Eq. (5.117): ln



ΦLi



  " L  pffiffiffi #  L  Z + 1+ 2 B bi ðZL  1Þ A 2Ψ i bi pffiffiffi  ln  ln Z  B  pffiffiffi  ¼ bm ZL  1  2 B 2 2B ðaαÞm bm

548 CHAPTER 5 Equations of State and Phase Equilibria

Step 4. Repeat step 3 using the calculated composition of the gas phase, yi, to determine Φvi .   " v  pffiffiffi #  L  bi ðZv  1Þ Z + 1+ 2 B A 2Ψ i bi v pffiffiffi  ln ln Φi ¼  lnðZ  BÞ  pffiffiffi  bm 2 2B ðaαÞm bm Zv  1  2 B

Step 5. Calculate the new set of equilibrium ratios from Ki ¼

ΦLi Φvi

Step 6. Check for the solution by applying the following constraint: n  X 2 Ki =KiA  1  ε i¼1

where ε ¼ preset error tolerance, such as 0.0001 n ¼ number of components in the system. If these conditions are satisfied, then the solution has been reached. If not, steps 1–6 are repeated using the calculated equilibrium ratios as initial values.

Dew Point Pressure A saturated vapor exists for a given temperature at the pressure at which an infinitesimal amount of liquid first appears. This pressure is referred to as the dew point pressure, pd. The dew point pressure of a mixture is described mathematically by the following two conditions: yi ¼ zi , for 1  i  n and nv ¼ 1 n   X zi ¼1 K i i1

(5.136) (5.137)

Applying the definition of Ki in terms of the fugacity coefficient to Eq. (5.137) gives n   X zi i¼1

ki

" n X  ¼ i¼1

# X zi  f v  zi i  ¼1 ¼ ΦLi zi pd ΦLi =Φvi

or pd ¼

n  v X f i

i¼1

ΦLi

Equation-of-State Applications 549

This above equation is arranged to give f ðpd Þ ¼

n  v X f i

i¼1

ΦLi

 pd ¼ 0

(5.138)

where pd ¼ dew point pressure, psia f vi ¼ fugacity of component i in the vapor phase, psia ΦLi ¼ fugacity coefficient of component i in the liquid phase Eq. (5.35) can be solved for the dew point pressure by using the NewtonRaphson iterative method. To use the iterative method, the derivative of Eq. (5.138) with respect to the dew point pressure, pd, is required. This derivative is given by the following expression: "    # n ΦLi @fiv =@pd  fiv @ΦLi =@pd @f X 1 ¼  L 2 @pd i¼1 Φ

(5.139)

i

The two derivatives in this equation can be approximated numerically as follows:  v  @f v f ðpd + Δpd Þ  fiv ðpd  Δpd Þ ¼ i @pd 2Δpd

(5.140)

 L  @ΦLi Φ ðpd + Δpd Þ  ΦLi ðpd  Δpd Þ ¼ i @pd 2Δpd

(5.141)

and

where Δpd ¼ pressure increment (eg, 5 psia) fiv ðpd + Δpd Þ ¼ fugacity of component i at (pd + Δpd) f vi (pd–Δpd) ¼ fugacity of component i at (pd  Δpd) ΦLi ðpd + Δpd Þ ¼ fugacity coefficient of component i at (pd + Δpd) ΦLi (pd–Δpd) ¼ fugacity coefficient of component i at (pd  Δpd) ΦLi ¼ fugacity coefficient of component i at pd The computational procedure of determining pd is summarized in the following steps. Step 1. Assume an initial value for the dew point pressure, pA d. Step 2. Using the assumed value of pA d , calculate a set of equilibrium ratios for the mixture by using any of the previous correlations, such as Wilson’s correlation.

550 CHAPTER 5 Equations of State and Phase Equilibria

Step 3. Calculate the composition of the liquid phase, that is, the composition of the droplets of liquid, by applying the mathematical definition of Ki, to give xi ¼

zi Ki

It should be noted that yi ¼ zi. L Step 4. Calculate pA d using the composition of the gas phase, zi and Φi , using the composition of liquid phase, xi, at the following three pressures: n n n

pA d pA d + Δpd pA d –Δpd

where pA d is the assumed dew point pressure and Δpd is a selected pressure increment of 5–10 psi. Step 5. Evaluate the function f(pd) [ie, Eq. (5.138)] and its derivative using Eqs. (5.139)–(5.141). Step 6. Using the values of the function f(pd) and the derivative @f/@pd as determined in step 5, calculate a new dew point pressure by applying the Newton-Raphson formula: pd ¼ pA d  f ðpd Þ=½@f =@pd 

(5.142)

Step 7. The calculated value of pd is checked numerically against the assumed value by applying the following condition: jpd  pA dj5

If this condition is met, then the correct dew point pressure, pd, has been found. Otherwise, steps 3–6 are repeated using the calculated pd as the new value for the next iteration. A set of equilibrium ratios must be calculated at the new assumed dew point pressure from ki ¼

ΦLi Φvi

Bubble Point Pressure The bubble point pressure, pb, is defined as the pressure at which the first bubble of gas is formed. Accordingly, the bubble point pressure is defined mathematically by the following equations: yi ¼ zi and nL ¼ 1:0

(5.143)

n X ½zi Ki  ¼ 1

(5.144)

i¼1

Equation-of-State Applications 551

Introducing the concept of the fugacity coefficient into Eq. (5.144) gives 2  L 3 fi   n n 6 X X 7 ΦLi z 6zi i pb 7 ¼ 1 zi v ¼ v 5 4 Φi Φi i¼1 i¼1

Rearranging, pb ¼

n  L X f i

Φvi

i¼1

or f ðpb Þ ¼

n  L X f i

i¼1

  p b ¼0 v

Φi

(5.145)

The iteration sequence for calculation of pb from the preceding function is similar to that of the de-point pressure, which requires differentiating that function with respect to the bubble point pressure: "    # n Φvi @fiL =@pb  fiL @Φvi =@pb @f X 1 ¼  v 2 @pb i¼1 Φ

(5.146)

i

The phase envelope can be constructed by either calculating the saturation pressures, pb and pd, as a function of temperature or calculating the bubble point and dew point temperatures at selected pressures. The selection to calculate a pressure or temperature depends on how steep or flat the phase envelope is at a given point. In terms of the slope of the envelope at a point, the criteria are n

Calculate the bubble point or dew point pressure if      @ lnðpÞ   lnðp2 Þ  lnðp1 Þ     @ lnðT Þ lnðT Þ  lnðT Þ < 20 2 1

n

Calculate the bubble point or dew point temperature if      @ lnðpÞ   lnðp2 Þ  lnðp1 Þ     @ lnðT Þ lnðT Þ  lnðT Þ > 20 2 1

Three-Phase Equilibrium Calculations Two- and three-phase equilibria occur frequently during the processing of hydrocarbon and related systems. Peng and Robinson (1976b) proposed a three-phase equilibrium calculation scheme of systems that exhibit a water-rich liquid phase, a hydrocarbon-rich liquid phase, and a vapor phase.

552 CHAPTER 5 Equations of State and Phase Equilibria

In applying the principle of mass conversation to 1 mol of a waterhydrocarbon in a three-phase state of thermodynamic equilibrium at a fixed temperature, T, and pressure, p, gives nL + nw + nv ¼ 1

(5.147)

nL xi + nw xwi + nv yi ¼ zi

(5.148)

n n n n X X X X xi ¼ xwi ¼ yi ¼ zi ¼ 1

(5.149)

i

i¼1

i¼1

i¼1

where nL, nw, nυ ¼ number of moles of the hydrocarbon-rich liquid, the waterrich liquid, and the vapor, respectively xi, xwi, yi ¼ mole fraction of component i in the hydrocarbon-rich liquid, the water rich liquid, and the vapor, respectively The equilibrium relations between the compositions of each phase are defined by the following expressions: Ki ¼

yi ΦLi ¼ xi Φvi

(5.150)

Kwi ¼

yi Φw ¼ i xwi Φvi

(5.151)

and

where Ki ¼ equilibrium ratio of component i between vapor and hydrocarbonrich liquid Kwi ¼ equilibrium ratio of component i between the vapor and waterrich liquid ΦLi ¼ fugacity coefficient of component i in the hydrocarbon-rich liquid Φvi ¼ fugacity coefficient of component i in the vapor phase Φw i ¼ fugacity coefficient of component i in the water-rich liquid Combining Eqs. (5.147)–(5.151) gives the following conventional nonlinear equations: 3

2 X X6 6 xi ¼ 4 i¼1

i¼1

7 7¼1  5 Ki  Ki + Ki Kwi

zi  nL ð1  Ki Þ + nw

(5.152)

Equation-of-State Applications 553

3

2

X X6 7 zi ðKi =Kwi Þ 7¼1 6   xwi ¼ 5 4 K i i¼1 i¼1 n ð1  K Þ + n  Ki + Ki L i w Kwi 3 2 X

yi ¼

i¼1

X6 6 4 i¼1

7 7¼1  5 Ki  Ki + Ki Kwi

zi Ki  nL ð1  Ki Þ + nw

(5.153)

(5.154)

Assuming that the equilibrium ratios between phases can be calculated, these equations are combined to solve for the two unknowns, nL and nυ, and hence xi, xwi, and yi. The nature of the specific equilibrium calculation is what determines the appropriate combination of Eqs. (5.152)–(5.154). The combination of these three expressions then can be used to determine the phase and volumetric properties of the three-phase system. There are essentially three types of phase behavior calculations for the threephase system: bubble point prediction, dew point prediction, and flash calculation. Peng and Robinson (1980) proposed the following combination schemes of Eqs. (5.152)–(5.154). For the bubble point pressure determination, X X xi  xwi ¼ 0 i

i

and "

X

# yi  1 ¼ 0

i

Substituting Eqs. (5.152)–(5.154) in this relationships gives f ðnL , nw Þ ¼

X

 zi ð1  Ki =Kwi Þ ¼0 nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki

(5.155)

 zi K i 1¼0 nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki

(5.156)

i

and gðnL , nw Þ ¼

X i

For the dew point pressure, X X xwi  yi ¼ 0 "i # i X xi  1 ¼ 0 i

554 CHAPTER 5 Equations of State and Phase Equilibria

Combining with Eqs. (5.152)–(5.154) yields f ðnL , nw Þ ¼

X

 zi Ki ð1=Kwi  1Þ ¼0 nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki

(5.157)

 zi 1¼0 nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki

(5.158)

i

and gðnL , nw Þ ¼

X i

For flash calculations, X X xi  yi ¼ 0 i

i

" # X xwi  1 ¼ 0 i

or f ðnL , nw Þ ¼

X

 zi ð1  Ki Þ ¼0 nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki

(5.159)

 zi Ki =Kwi  1:0 ¼ 0 nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki

(5.160)

i

and gðnL , nw Þ ¼

X i

Note that, in performing any of these property predictions, we always have two unknown variables, nL and nw, and between them, two equations. Providing that the equilibrium ratios and the overall composition are known, the equations can be solved simultaneously using the appropriate iterative technique, such as the Newton-Raphson method. The application of this iterative technique for solving Eqs. (5.159), (5.160) is summarized in the following steps. Step 1. Assume initial values for the unknown variables nL and nw. Step 2. Calculate new values of nL and nw by solving the following two linear equations: 

nL nw

new

   1  f ðnL , nw Þ @f =@nL @f =@nw nL  ¼ gðnL , nw Þ nw @g=@nL @g=@nw 

where f(nL, nw) ¼ value of the function f(nL, nw) as expressed by Eq. (5.159) and g(nL, nw) ¼ value of the function g(nL, nw) as expressed by Eq. (5.160). The first derivative of these functions with respect to nL and nw is given by the following expressions:

Equation-of-State Applications 555

ð@f =@nL Þ ¼

( X

zi ð1  Ki Þ2

)

½nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki 2 ( ) X zi ð1  Ki ÞðKi =Kwi  Ki Þ i¼1

ð@f =@nw Þ ¼

i¼1

ð@g=@nL Þ ¼

( X

½nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki 2 zi ðKi =Kwi Þð1  Ki Þ

)

½nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki 2 ( ) X zi ðKi Kwi ÞðKi =Kwi  Ki Þ i¼1

ð@g=@nw Þ ¼

½nL ð1  Ki Þ + nw ðKi =Kwi  Ki Þ + Ki 2

i¼1

Step 3. The new calculated values of nL and nw then are compared with the initial values. If no changes in the values are observed, then the correct values of nL and nw have been obtained. Otherwise, the steps are repeated with the new calculated values used as initial values. Peng and Robinson (1980) proposed two modifications when using their equation of state for three-phase equilibrium calculations: n

The first modification concerns the use of the parameter α as expressed by Eq. (5.113) for the water compound. Peng and Robinson suggested that, when the reduced temperature of this compound is less than 0.85, the following equation should be applied:    0:5 2 αw ¼ 1:0085677 + 0:82154 1  Trw

(5.161)

where

n

Trw is the reduced temperature of the water component; that is, Trw ¼ T/ Tcw. The second important modification of the PR EOS is the introduction of alternative mixing rules and new interaction coefficients for the parameter (aα). A temperature-dependent binary interaction coefficient was introduced into the equation, to give ðaαÞwater ¼

XXh i

 0:5  i xwi xwj ai aj αi αj 1  τij

(5.162)

j

where τij ¼ temperature-dependent binary interaction coefficient and xwi ¼ mole fraction of component i in the water phase. Peng and Robinson proposed graphical correlations for determining this parameter for each aqueous binary pair. Lim and Ahmed (1984) expressed these graphical correlations mathematically by the following generalized equation:

556 CHAPTER 5 Equations of State and Phase Equilibria

 τij ¼ a1

T Tci

2 

pci pcj

2

 + a2

T Tci

  pci + a3 pcj

(5.163)

where T ¼ system temperature, °R Tci ¼ critical temperature of the component of interest, °R pci ¼ critical pressure of the component of interest, psia pcj ¼ critical pressure of the water compound, psia Values of the coefficients a1, a2, and a3 of this polynomial are given in the table below for selected binaries. Component i

a1

a2

a3

C1 C2 C3 n-C4 n-C6

0 0 18.032 0 49.472

1.659 2.109 9.441 2.800 5.783

0.761 0.607 1.208 0.488 0.152

For selected nonhydrocarbon components, values of interaction parameters are given by the following expressions: For N2-H2O binary, τij ¼ 0:402ðT=Tci Þ  1:586

(5.164)

where τij ¼ binary parameter between nitrogen and the water compound Tci ¼ critical temperature of nitrogen, °R. For CO2-H2O binary,  2   T T  0:503 τij ¼ 0:074 + 0:478 Tci Tci

(5.165)

where Tci is the critical temperature of CO2. In the course of making phase equilibrium calculations, it is always desirable to provide initial values for the equilibrium ratios so the iterative procedure can proceed as reliably and rapidly as possible. Peng and Robinson adopted Wilson’s equilibrium ratio correlation to provide initial K-values for the hydrocarbon/vapor phase: Ki ¼

   pci Tci exp 5:3727ð1 + ωi Þ 1  p T

Equation-of-State Applications 557

While for the water-vapor phase, Peng and Robinson proposed the following expression for Kwi: Kwi ¼ 106 ½pci T=ðTci pÞ

Vapor Pressure The calculation of the vapor pressure of a pure component through an equation of state usually is made by the same trial-and-error algorithms used to calculate vapor-liquid equilibria of mixtures. Soave (1972) suggests that the van der Waals (vdW), Soave-Redlich-Kwong (SRK), and the PengRobinson (PR) equations of state can be written in the following generalized form: p¼

RT aα  V  b V 2 + μVb + wb2

(5.166)

with R2 Tc2 pc RTc b ¼ Ωb pc a ¼ Ωa

where the values of u, w, Ωa, and Ωb for three different equations of state are given in the table below. EOS

u

w

Ωa

Ωb

vdW SRK PR

0 1 2

0 0 1

0.421875 0.42748 0.45724

0.125 0.08664 0.07780

Soave introduced the reduced pressure, pr, and reduced temperature, Tr, to these equations, to give A¼

aαp αpr ¼ Ωa Tr R2 T 2

(5.167)

bp pr ¼ Ωb Tr RT

(5.168)

B¼

and

  A Ωa α ¼ B Ωb Tr

where pr ¼ p=pc Tr ¼ T=Tc

(5.169)

558 CHAPTER 5 Equations of State and Phase Equilibria

In the cubic form and in terms of the Z-factor, the three equations of state can be written as vdW : Z3  Z2 ð1 +BÞ + ZA  AB  ¼0 SRK : Z 3  Z2 + Z A  B  B2  AB ¼ 0     PR : Z 3  Z2 ð1  BÞ + Z A  3B2  2B  AB  B2  B2 ¼ 0

And the pure component fugacity coefficient is given by  o f A vdW : ln ¼ Z  1  lnðZ  BÞ  p Z  o     f A B ¼ Z  1  lnðZ  BÞ  SRK : ln ln 1 + p B Z  pffiffiffi #  o   " Z+ 1+ 2 B f A pffiffiffi  ¼ Z  1  lnðZ  BÞ  pffiffiffi ln PR : ln p Z 1 2 B 2 2B

A typical iterative procedure for the calculation of pure component vapor pressure at any temperature, T, through one of these EOS is summarized next. Calculate the reduced temperature; that is, Tr ¼ T/Tc. Calculate the ratio A/B from Eq. (5.169). Assume a value for B. Solve the selected equation of state, such as SRK EOS, to obtain ZL and Z (ie, smallest and largest roots) for both phases. Step 5. Substitute ZL and Zv into the pure component fugacity coefficient and obtain ln(fo/p) for both phases. Step 6. Compare the two values of f/p. If the isofugacity condition is not satisfied, assume a new value of B and repeat steps 3–6. Step 7. From the final value of B, obtain the vapor pressure from Eq. (5.168): Step 1. Step 2. Step 3. Step 4.

B ¼ Ωb

ðpv =pc Þ Tr

Solving for pv gives pv ¼

BTr pc Ωb

SIMULATING OF LABORATORY PVT DATA BY EQUATIONS OF STATE Before attempting to use an EOS-based compositional model in a full-field simulation study, it is essential the selected equation of state is capable of achieving a satisfactory match between EOS results and all the available

Simulating of Laboratory PVT Data by Equations of State 559

PVT test data. It should be pointed out that an equation of state generally is not predictive without tuning its parameters to match the relevant experimental data. Several PVT laboratory tests were presented in Chapter 4 and are repeated here for convenience to illustrate the procedure of applying an EOS to simulate these experiments, including n n n n n n n n

Constant volume depletion. Constant composition expansion. Differential liberation. Flash separation tests. Composite liberation. Swelling tests. Compositional gradients. Minimum miscibility pressure.

Constant Volume Depletion Test A reliable prediction of the pressure depletion performance of a gascondensate reservoir is necessary in determining reserves and evaluating field-separation methods. The predicted performance is also used in planning future operations and studying the economics of projects for increasing liquid recovery by gas cycling. Such predictions can be performed with aid of the experimental data collected by conducting constant volume depletion (CVD) tests on gas condensates. These tests are performed on a reservoir fluid sample in such a manner as to simulate depletion of the actual reservoir, assuming that retrograde liquid appearing during production remains immobile in the reservoir. The CVD test provides five important laboratory measurements that can be used in a variety of reservoir engineering predictions: n n n n n

Dew point pressure. Composition changes of the gas phase with pressure depletion. Compressibility factor at reservoir pressure and temperature. Recovery of original in-place hydrocarbons at any pressure. Retrograde condensate accumulation; that is, liquid saturation.

The laboratory procedure of the CVD test (with immobile condensate), as shown schematically in Fig. 5.10, is summarized in the following steps. Step 1. A measured amount, m, of a representative sample of the original reservoir fluid with a known overall composition of zi is charged to a visual PVT cell at the dew point pressure, pd (section A). The temperature of the PVT cell is maintained at the reservoir

560 CHAPTER 5 Equations of State and Phase Equilibria

Constant volume depletion (CVD) Gas Pd

Vi

Original gas @ Pd

P

P

Gas

Gas

Vt

Vi

Condensate Condensate Hg

Hg

Hg

A

B

C

n FIGURE 5.10 Schematic illustration of the CVD test.

temperature, T, throughout the experiment. The initial volume, Vi, of the saturated fluid is used as a reference volume. Step 2. The initial gas compressibility factor is calculated from the real gas equation: Zdew ¼

pd Vi ni RT

(5.170)

with ni ¼

m Ma

where pd ¼ dew point pressure, psia Vi ¼ initial gas volume, ft3 ni ¼ initial number of moles of the gas Ma ¼ apparent molecular weight, lb/lb-mol m ¼ mass of the initial gas in place, lb R ¼ gas constant, 10.73 T ¼ temperature, °R zd ¼ compressibility factor at dew-point pressure with the gas initially in place as expressed in standard units, scf, given by G ¼ 379:4ni

(5.171)

Simulating of Laboratory PVT Data by Equations of State 561

Step 3. The cell pressure is reduced from the saturation pressure to a predetermined level, p. This can be achieved by withdrawing mercury from the cell, as illustrated in Fig. 5.10, section B. During the process, a second phase (retrograde liquid) is formed. The fluid in the cell is brought to equilibrium and the total volume, Vt, and volume of the retrograde liquid, VL, are visually measured. This retrograde volume is reported as a percent of the initial volume, Vi: 

SL ¼

 VL 100 Vi

This value can be translated into condensate saturation in the presence of water saturation by So ¼ ð1  Sw ÞSL =100

Step 4. Mercury is reinjected into the PVT cell at constant pressure p while an equivalent volume of gas is removed. When the initial volume, Vi, is reached, mercury injection ceases, as illustrated in Fig. 5.10, section C. The volume of the removed gas is measured at the cell conditions and recorded as (Vgp)p,T. This step simulates a reservoir producing only gas, with retrograde liquid remaining immobile in the reservoir. Step 5. The removed gas is charged to analytical equipment, where its composition, yi, is determined and its volume is measured at standard conditions and recorded as (Vgp)sc. The corresponding moles of gas produced can be calculated from the expression 

 Vgp sc np ¼ 379:4

(5.172)

where np ¼ moles of gas produced (Vgp)sc ¼ volume of gas produced measured at standard conditions, scf. Step 6. The compressibility factor of the gas phase at cell pressure and temperature is calculated by   Vgp p, T Vactual Z¼ ¼ Videal Videal

where Videal ¼

RTnp p

(5.173)

562 CHAPTER 5 Equations of State and Phase Equilibria

Combining Eqs. (5.172), (5.173) and solving for the compressibility factor Z gives   379:4p Vgp p, T   Z¼ RT Vgp sc

(5.174)

Another property, the two-phase compressibility factor, also is calculated. The two-phase Z-factor represents the total compressibility of all the remaining fluid (gas and retrograde liquid) in the cell and is computed from the real gas law as Ztwo-phase ¼ 

pVi  ni  np ðRT Þ

(5.175)

where (ni  np) ¼ represents the remaining moles of fluid in the cell ni ¼ initial moles in the cell np ¼ cumulative moles of gas removed Step 7. The volume of gas produced as a percentage of gas initially in place is calculated by dividing the cumulative volume of the produced gas by the gas initially in place, both at standard conditions: "X %Gp ¼

Vgp

 # sc

G

100

(5.176)

or equivalently as " X # np 100 %Gp ¼ ðni Þoriginal

This experimental procedure is repeated several times, until a minimum test pressure is reached, after which the quantity and composition of the gas and retrograde liquid remaining in the cell are determined. This test procedure also can be conducted on a volatile oil sample. In this case, the PVT cell initially contains liquid, instead of gas, at its bubble point pressure.

Simulating of the CVD Test by EOS In the absence of the CVD test data on a specific gas-condensate system, predictions of pressure-depletion behavior can be obtained by using any of the well-established equations of state to compute the phase behavior when the composition of the total gas-condensate system is known. The stepwise computational procedure using the Peng-Robinson EOS as a

Simulating of Laboratory PVT Data by Equations of State 563

representative equation of state now is summarized in conjunction with the flow diagram shown in Fig. 5.11. Step 1. Assume that the original hydrocarbon system with a total composition zi occupies an initial volume of 1 ft3 at the dew point pressure pd and system temperature T: Vi ¼ 1

Step 2. Using the initial gas composition, calculate the gas compressibility factor Zdew at the dew point pressure from Eq. (5.115):     Z3 + ðB  1ÞZ 2 + A  3B2  2B Z  AB  B2  B3 ¼ 0

Given: zi, Pd, T

Calculate Zdew from: Z3 + (B − 1) Z2 + (A − 3 B2 − 2B)Z − (AB − B2 − B3) = 0

Assume Vi =1 ft3 Vp Calculate ni

ni =

ZRT

=

(1) (pd) Zdew TR

Assume depletion pressure “P”

Calculate Ki for zi @ P see Fig. 15-11

VL =

(nL)actual ZL RT P

Solution; gives: Ki, xi, yi, nL, nV, ZL, ZV

Calculate volume of condensate VL

Calculate volume of gas Vg Calculate TOTAL volume Vt = Vg + VL

Determine moles of gas removed nP & remaining nr xi (nL)actual + yi nr zi =

(nL)actual + nr

n FIGURE 5.11 Flow diagram for EOS simulating CVD test.

Combine the remaining phase and determine new Zi

Vg =

(nV)actual ZV RT

(Vg)P = Vt − 1 np =

P(Vg)P

ZL RT nr = (nV)actual − np

P

564 CHAPTER 5 Equations of State and Phase Equilibria

Step 3. Calculate the initial moles in place by applying the real gas law: ni ¼

ð1Þðpd Þ Zdew RT

(5.177)

where pd ¼ dew point pressure, psi Zdew ¼ gas compressibility factor at the dew point pressure. Step 4. Reduce the pressure to a predetermined value p. At this pressure level, the equilibrium ratios (K-values) are calculated from EOS. Results of the calculation at this stage include n n n n n n n

Equilibrium ratios (K-values). Composition of the liquid phase (ie, retrograde liquid), xi. Moles of the liquid phase, nL. Composition of the gas phase, yi. Moles of the gas phase, nυ. Compressibility factor of the liquid phase, ZL. Compressibility factor of the gas phase, Zv.

The composition of the gas phase, yi, should reasonably match the experimental gas composition at pressure, p. Step 5. Because flash calculations, as part of the K-value results, are usually performed assuming the total moles are equal to 1, calculate the actual moles of the liquid and gas phases from ðnL Þactual ¼ ni nL

(5.178)

ðnυ Þactual ¼ ni nυ

(5.179)

where nL and nυ are moles of liquid and gas, respectively, as determined from flash calculations. Step 6. Calculate the volume of each hydrocarbon phase by applying the expression VL ¼

ðnL Þactual ZL RT p

(5.180)

Vg ¼

ðnv Þactual Zv RT p

(5.181)

where VL ¼ volume of the retrograde liquid, ft3/ft3 Vg ¼ volume of the gas phase, ft3/ft3

Simulating of Laboratory PVT Data by Equations of State 565

ZL, Zv ¼ compressibility factors of the liquid and gas phase T ¼ cell (reservoir) temperature, °R As Vi ¼ 1, then SL ¼ ðVL Þ100

This value should match the experimental value if the equation of state is properly tuned. Step 7. Calculate the total volume of fluid in the cell: Vt ¼ VL + Vg

(5.182)

where Vt ¼ total volume of the fluid, ft3. Step 8. As the volume of the cell is constant at 1 ft3, remove the following excess gas volume from the cell:   Vgp p, T ¼ Vt  1

(5.183)

Step 9. Calculate the number of moles of gas removed: np ¼

  p Vgp p, T

(5.184)

Zv RT

Step 10. Calculate cumulative gas produced as a percentage of gas initially P in place by dividing cumulative moles of gas removed, np, by the original moles in place:

%Gp ¼

" X  # np ðni Þoriginal

100

Step 11. Calculate the two-phase gas deviation factor from the relationship: Ztwo-phase ¼ 

ðpÞð1Þ  ni  np RT

(5.185)

Step 12. Calculate the remaining moles of gas (nυ)r by subtracting the moles produced (np) from the actual number of moles of the gas phase (nυ)actual: ðnυ Þr ¼ ðnυ Þactual  np

(5.186)

566 CHAPTER 5 Equations of State and Phase Equilibria

Step 13. Calculate the new total moles and new composition remaining in the cell by applying the molal and component balances, respectively: ni ¼ ðnL Þactual + ðnυ Þr ; zi ¼

(5.187)

xi ðnL Þactual + yi ðnυ Þr ni

(5.188)

Step 14. Consider a new lower pressure and repeat steps 4–13.

Constant Composition Expansion Test Constant composition expansion experiments, commonly called pressure/ volume tests, are performed on gas condensates or crude oil to simulate the pressure/volume relations of these hydrocarbon systems. The test is conducted for the purposes of determining Saturation pressure (bubble point or dew-point pressure). Isothermal compressibility coefficients of the single-phase fluid in excess of saturation pressure. Compressibility factors of the gas phase. Total hydrocarbon volume as a function of pressure.

n n

n n

The experimental procedure, shown schematically in Fig. 5.12, involves placing a hydrocarbon fluid sample (oil or gas) in a visual PVT cell at

Voil = 100%

P2 > Pb

P1 >> Pb

100% Oil

Voil ≈ 100%

Voil = 100%

Vt1

Oil

Voil < 100%

P3 = Pb

Vt2

Oil

Voil << 100%

P4 << Pb

P5 <<< Pb

Gas

Gas

Vsat

Oil

Vt4

Oil

Hg

(a)

(b)

n FIGURE 5.12 Constant composition expansion test.

(c)

(d)

(e)

Vt5

Simulating of Laboratory PVT Data by Equations of State 567

reservoir temperature and a pressure in excess of the initial reservoir pressure (Fig. 5.12, section A). The pressure is reduced in steps at constant temperature by removing mercury from the cell, and the change in the hydrocarbon volume, V, is measured for each pressure increment. The saturation pressure (bubble point or dew point pressure) and the corresponding volume are observed and recorded (Fig. 5.12, section B). The volume at the saturation pressure is used as a reference volume. At pressure levels higher than the saturation pressure, the volume of the hydrocarbon system is recorded as a ratio of the reference volume. This volume, commonly termed the relative volume, is expressed mathematically by the following equation: Vrel ¼

V Vsat

(5.189)

where Vrel ¼ relative volume V ¼ volume of the hydrocarbon system Vsat ¼ volume at the saturation pressure Also, above the saturation pressure, the isothermal compressibility coefficient of the single-phase fluid usually is determined from the expression c¼

  1 @Vrel Vrel @p T

(5.190)

where c ¼ isothermal compressibility coefficient, psi1. For gas-condensate systems, the gas compressibility factor, Z, is determined in addition to the preceding experimentally determined properties. Below the saturation pressure, the two-phase volume, Vt, is measured relative to the volume at saturation pressure and expressed as Relativetotalvolume ¼

Vt Vsat

(5.191)

where Vt ¼ total hydrocarbon volume. Note that no hydrocarbon material is removed from the cell; therefore, the composition of the total hydrocarbon mixture in the cell remains fixed at the original composition.

568 CHAPTER 5 Equations of State and Phase Equilibria

Simulating of the CCE Test by EOS The simulation procedure utilizing the Peng-Robinson equation of state is illustrated in the following steps. Step 1. Given the total composition of the hydrocarbon system, zi, and saturation pressure (pb for oil systems, pd for gas systems), calculate the total volume occupied by 1 mol of the system. This volume corresponds to the reference volume Vsat (volume at the saturation pressure). Mathematically, the volume is calculated from the relationship Vsat ¼

ð1ÞZRT psat

(5.192)

where Vsat ¼ volume of saturation pressure, ft3/mol psat ¼ saturation pressure (dew point or bubble point pressure), psia T ¼ system temperature, °R Z ¼ compressibility factor, ZL or Zv depending on the type of system Step 2. The pressure is increased in steps above the saturation pressure, where the single phase still exists. At each pressure, the compressibility factor, ZL or Zv, is calculated by solving Eq. (5.115) and used to determine the fluid volume. V¼

ð1ÞZRT p

where V ¼ compressed liquid or gas volume at the pressure level p, ft3/mol Z ¼ compressibility factor of the compressed liquid or gas, ZL or Zv p ¼ system pressure, p > psat, psia R ¼ gas constant; 10.73 psi ft3/lb-mol °R The corresponding relative phase volume, Vrel, is calculated from the expression Vrel ¼

V Vsat

Step 3. The pressure is then reduced in steps below the saturation pressure psat, where two phases are formed. The equilibrium ratios are calculated and flash calculations are performed at each pressure level. Results of the calculation at each pressure level include Ki, xi, yi, nL, nv, Zv, and ZL. Because no hydrocarbon material is removed during pressure depletion, the original moles (ni ¼ 1) and

Simulating of Laboratory PVT Data by Equations of State 569

composition, zi, remain constant. The volumes of the liquid and gas phases can then be calculated from the expressions VL ¼

ð1ÞðnL ÞZL RT p

(5.193)

Vg ¼

ð1Þðnυ ÞZ v RT p

(5.194)

and Vt ¼ VL + Vg

where nL, nυ ¼ moles of liquid and gas as calculated from flash calculations ZL, Zv ¼ liquid and gas compressibility factors Vt ¼ total volume of the hydrocarbon system Step 4. Calculate the relative total volume from the following expression: Relativetotalvolume ¼

Vt Vsat

Differential Liberation Test The test is carried out on reservoir oil samples and involves charging a visual PVT cell with a liquid sample at the bubble point pressure and reservoir temperature, as described in detail in Chapter 4. The pressure is reduced in steps, usually 10–15 pressure levels, and all the liberated gas is removed and its volume, Gp, is measured under standard conditions. The volume of oil remaining, VL, also is measured at each pressure level. Note that the remaining oil is subjected to continual compositional changes as it becomes progressively richer in the heavier components. This procedure is continued to atmospheric pressure, where the volume of the residual (remaining) oil is measured and converted to a volume at 60°F, Vsc. The differential oil formation volume factors Bod (commonly called the relative oil volume factors) at all the various pressure levels are calculated by dividing the recorded oil volumes, VL, by the volume of residual oil, Vsc: Bod ¼

VL Vsc

(5.195)

The differential solution gas-oil ratio Rsd is calculated by dividing the volume of gas in the solution by the residual oil volume. Typical table and in Fig. 5.13.

570 CHAPTER 5 Equations of State and Phase Equilibria

Bod

Rsd

B od

R sd

Pressure

n FIGURE 5.13 Bod and Rsd versus pressure.

Pressure (psig)

Rsda

Bodb

Shrinkage (bbl/bbl)

2620 2350 2100 1850 1600 1350 1100 850 600 350 159 0 0

854 763 684 612 544 479 416 354 292 223 157 0 0

1.6 1.554 1.515 1.479 1.445 1.412 1.382 1.351 1.32 1.283 1.244 1.075 at 60°F ¼ 1.000

1.0000 0.9713 0.9469 0.9244 0.9031 0.8825 0.8638 0.8444 0.8250 0.8019 0.7775 0.6250

a b

Cubic feet of gas at 14.65 psia and 60°F per barrel of residual oil at 60°F. Barrels of oil at indicated pressure and temperature per barrel of residual oil at 60°F.

It should be pointed out that reporting the experimental data relative to the residual oil volume at 60°F gives the relative oil volume curve the appearance of the formation volume factor curve, leading to its misuse in reservoir calculations. A better method of reporting these data is in the form of a shrinkage curve, as shown in Fig. 5.14. The relative-oil volume data in Fig. 5.13 can be converted to a shrinkage curve by dividing each

Sod

Simulating of Laboratory PVT Data by Equations of State 571

Pressure

n FIGURE 5.14 Oil shrinkage curve versus pressure.

relative-oil volume factor Bod by the relative-oil volume factor at the bubble point, Bodb: Sod ¼

Bod Bodb

(5.196)

where Bod ¼ differential relative-oil volume factor at pressure p, bbl/STB Bodb ¼ differential relative-oil volume factor at the bubble point pressure, pb, psia, bbl/STB Sod ¼ differential oil shrinkage factor, bbl/bbl, of bubble point oil The shrinkage curve has a value of 1 at the bubble point and a value less than 1 at subsequent pressures below pb. In this suggested form, the shrinkage describes the changes in the volume of the bubble point oil as reservoir pressure declines. The results of the differential liberation test, when combined with the results of the flash separation test, provide a means of calculating the oil formation volume factors and gas solubility as a function of reservoir pressure. These calculations are outlined in the next section.

Simulating the Differential Liberation Test by EOS The simulation procedure of the test by the Peng-Robinson EOS is summarized in the following steps and in conjunction with the flow diagram shown in Fig. 5.15.

572 CHAPTER 5 Equations of State and Phase Equilibria

Given: zi, Pb, T

Calculate ZL from: Z3 + (B − 1) Z2 + (A − 3B2 − 2B)Z − (AB − B2 − B3) = 0

Assume ni = 1 mole

Calculate: Vsat

(1) Z LRT

nZRT Vsat =

P

=

Pb

Assume depletion pressure P < Pb

Calculate Ki for zi @ P See Fig. 15-11

VL =

(nL)actual ZLRT P

Calculate volume of condensate VL

Solution gives: Ki, xi, yi, nL, nV, ZL, ZV

Calculate volume of gas Vg

(nV)actual ZVRT

P GP = 379.4 (nV)actual (GP)total =

Set: zi = xi and ni = nL

No

Vg =

ΣG

p

i=1

Last pressure ?

Calculate Bod & Rsd

n FIGURE 5.15 Flow diagram for EOS simulating DE test.

Step 1. Starting with the saturation pressure, psat, and reservoir temperature, T, calculate the volume occupied by a total of 1 mol (ie, ni ¼ 1) of the hydrocarbon system with an overall composition of zi. This volume, Vsat, is calculated by applying Eq. (5.192): Vsat ¼

ð1ÞZRT psat

Step 2. Reduce the pressure to a predetermined value of p at which the equilibrium ratios are calculated and used in performing flash calculations. The actual number of moles of the liquid phase, with a composition of xi, and the actual number of moles of the gas phase, with a composition of yi, then are calculated from the expressions

Simulating of Laboratory PVT Data by Equations of State 573

ðnL Þactual ¼ ni nL ðnυ Þactual ¼ ni nυ

where (nL)actual ¼ actual number of moles of the liquid phase (nυ)actual ¼ actual number of moles of the gas phase nL, nυ ¼ number of moles of the liquid and gas phases as computed by performing flash calculations Step 3. Determine the volume of the liquid and gas phase from VL ¼

ZL RT ðnL Þactual p

(5.197)

Vg ¼

Zv RT ðnυ Þactual p

(5.198)

where VL, Vg ¼ volumes of the liquid and gas phases, ft3 ZL, Zv ¼ compressibility factors of the liquid and gas phases. The volume of the produced (liberated) solution gas as measured at standard conditions is determined from the relationship Gp ¼ 379:4ðnυ Þactual

(5.199)

where Gp ¼ gas produced during depletion pressure p, scf. The total cumulative gas produced at any depletion pressure, p, is the cumulative gas liberated from the crude oil sample during the pressure reduction process (sum of all the liberated gas from previous pressures and current pressure) as calculated from the expression 

Gp



¼ p

p X Gp psat

where (Gp)p is the cumulative produced solution gas at pressure level p. Step 4. Assume that all the equilibrium gas is removed from contact with the oil. This can be achieved mathematically by setting the overall composition, zi, equal to the composition of the liquid phase, xi, and setting the total moles equal to the liquid phase: zi ¼ x i ni ¼ ðnL Þactual

Step 5. Using the new overall composition and total moles, steps 2–4 are repeated. When the depletion pressure reaches the atmospheric pressure, the temperature is changed to 60°F and the residual-oil

574 CHAPTER 5 Equations of State and Phase Equilibria

volume is calculated. This residual-oil volume, calculated by Eq. (5.197), is referred to as Vsc. The total volume of gas that evolved from the oil and produced (Gp)Total is the sum of the all-liberated gas including that at atmospheric pressure: 

Gp



¼ Total

14:7 X Gp psat

Step 6. The calculated volumes of the oil and removed gas then are divided by the residual-oil volume to calculate the relative-oil volumes (Bod) and the solution GOR at all the selected pressure levels from Bod ¼

VL Vsc

h    i ð5:615Þðremaining gas in solutionÞ ð5:615Þ Gp Total  Gp p Rsd ¼ ¼ Vsc Vsc (5.200)

Flash Separation Tests Flash separation tests, commonly called separator tests, are conducted to determine the changes in the volumetric behavior of the reservoir fluid as the fluid passes through the separator (or separators) and then into the stock tank. The resulting volumetric behavior is influenced to a large extent by the operating conditions; that is, pressures and temperatures of the surface separation facilities. The primary objective of conducting separator tests, therefore, is to provide the essential laboratory information necessary for determining the optimum surface separation conditions, which in turn maximize the stock-tank oil production. In addition, the results of the test, when appropriately combined with the differential liberation test data, provide a means of obtaining the PVT parameters (Bo, Rs, and Bt) required for petroleum engineering calculations. The laboratory procedure of the flash separation test involves the following steps. Step 1. The oil sample is charged into a PVT cell at reservoir temperature and its bubble point pressure. Step 2. A small amount of the oil sample, with a volume of (Vo)pb, is removed from the PVT cell at constant pressure and flashed through a multistage separator system, with each separator at a fixed pressure and temperature. The gas liberated from each stage (separator) is removed and its volume and specific gravity measured. The volume of the remaining oil in the last stage (stock-tank condition) is recorded as (Vo)st.

Simulating of Laboratory PVT Data by Equations of State 575

Step 3. The total solution gas-oil ratio and the oil formation volume factor at the bubble point pressure are then calculated: Bofb ¼ ðVo Þpb =ðVo Þst   Rsfb ¼ Vg sc =ðVo Þst

(5.201) (5.202)

where Bofb ¼ bubble point oil formation volume factor, as measured by flash liberation, bbl of the bubble point oil/STB Rsfb ¼ bubble point solution gas-oil ratio, as measured by flash liberation, scf/STB (Vg)sc ¼ total volume of gas removed from separators, scf Step 4. Repeat steps 1–3. A typical example of a set of separator tests for a two-stage separation is shown in the table below.

Separator pressure (psig)

Temperature (°F)

Rsfba

°API

Bofbb

50 to 0

75 75

40.5

1.481

100 to 0

75 75

40.7

1.474

200 to 0

75 75

40.4

1.483

300 to 0

75 75

737 41 P ¼ 778 676 92 P ¼ 768 602 178 P ¼ 780 549 246 P ¼ 795

40.1

1.495

a

GOR in cubic feet of gas at 14.65 psia and 60°F per barrel of stock-tank oil at 60°F. FVF in barrels of saturated oil at 2620 psig and 220°F per barrel of stock-tank oil at 60°F.

b

The differential liberation data must be adjusted to account for the separation process taking place at the separator condition. The procedure, as discussed in Chapter 4, is based on multiplying the flash oil formation factor at the bubble point, Bofb, by the differential oil shrinkage factor Sod (as defined by Eq. (5.196)) at various reservoir pressures. Mathematically, this relationship is expressed as follows: Bo ¼ Bofb Sod

(5.203)

576 CHAPTER 5 Equations of State and Phase Equilibria

where Bo ¼ oil formation volume factor, bbl/STB Bofb ¼ bubble point oil formation volume factor, bbl of the bubble point oil/STB Sod ¼ differential oil shrinkage factor, bbl/bbl of bubble point oil The adjusted differential gas solubility data are given by Rs ¼ Rsfb  ðRsdb  Rsd Þ

where

Bofb Bodb

(5.204)

Rs ¼ gas solubility, scf/STB Rsfb ¼ bubble point solution gas-oil ratio as defined by Eq. (5.202), scf/ STB Rsdb ¼ solution gas-oil ratio at the bubble point pressure as measured by the differential liberation test, scf/STB Rsd ¼ solution gas-oil ratio at various pressure levels as measured by the differential liberation test, scf/STB These adjustments typically produce lower formation volume factors and gas solubilities than the differential liberation data. This adjustment procedure is illustrated graphically in Figs. 5.16 and 5.17.

n FIGURE 5.16 Adjustment of oil-volume curve to separate conditions.

Simulating of Laboratory PVT Data by Equations of State 577

n FIGURE 5.17 Adjustment of gas-in-solution curve to separator conditions.

Composite Liberation Test The laboratory procedures for conducting the differential liberation and flash liberation tests for the purpose of generating Bo and Rs versus p relationships are considered approximations to the actual relationships. Another type of test, called a composite liberation test, has been suggested by Dodson et al. (1953) and represents a combination of differential and flash liberation processes. The laboratory test, commonly called the Dodson test, provides a better means of describing the PVT relationships. The experimental procedure for this composite liberation, as proposed by Dodson et al., is summarized in the following steps. Step 1. A large representative fluid sample is placed into a cell at a pressure higher than the bubble point pressure of the sample. The temperature of the cell then is raised to reservoir temperature. Step 2. The pressure is reduced in steps by removing mercury from the cell, and the change in the oil volume is recorded. The process is repeated until the bubble point of the hydrocarbon is reached.

578 CHAPTER 5 Equations of State and Phase Equilibria

Step 3. A carefully measured small volume of oil is removed from the cell at constant pressure and flashed at temperatures and pressures equal to those in the surface separators and stock tank. The liberated gas volume and stock-tank oil volume are measured. The oil formation volume factor, Bo, and gas solubility, Rs, then are calculated from the measured volumes: Bo ¼ ðVo Þp, T =ðVo Þst   Rs ¼ Vg sc =ðVo Þst

(5.205) (5.206)

where (Vo)p,T ¼ volume of oil removed from the cell at constant pressure p (Vg)sc ¼ total volume of the liberated gas as measured under standard conditions (Vo)st ¼ volume of the stock-tank oil Step 4. The volume of the oil remaining in the cell is allowed to expand through a pressure decrement and the gas evolved is removed, as in the differential liberation. Step 5. Following the gas removal, step 3 is repeated and Bo and Rs calculated. Step 6. Steps 3–5 are repeated at several progressively lower reservoir pressures to secure a complete PVT relationship. Note that this type of test, while more accurately representing the PVT behavior of complex hydrocarbon systems, is more difficult and costly to perform than other liberation tests. Consequently, these experiments usually are not included in a routine fluid property analysis.

Simulating of the Composite Liberation Test by EOS The stepwise computational procedure for simulating the composite liberation test using the Peng-Robinson EOS is summarized next in conjunction with the flow chart shown in Fig. 5.18. Step 1. Assume a large reservoir fluid sample with an overall composition of zi is placed in a cell at the bubble point pressure, pb, and reservoir temperature, T. Step 2 Remove 1 lb-mol (ie, ni ¼ 1) of the liquid from the cell and calculate the corresponding volume at pb and T by applying the PengRobinson EOS to estimate the liquid compressibility factor, ZL: ðVo Þp, T ¼

ð1ÞZL RT pb

(5.207)

Simulating of Laboratory PVT Data by Equations of State 579

Given: zi, Pb, T

Assume n = 1 mole

Calculate: VP,T

Flash zi at surface separation conditions. See Fig.15-11

Calculate Z L from: Z3 + (B − 1) Z2 + (A − 3B2 − 2B)Z − (AB − B2 − B3) = 0 (1) ZLRT

nZRT VP,T =

P

=

P

Solution gives: Ki, xi, yi, nL, nV, ZL, ZV

Calculate B0 and Rs

Assume P < Pb and Perform flash calculations

Set: zi = xi and ni = 1

P No Last pressure ?

Yes Results Bo & Rs

n FIGURE 5.18 Flow diagram for EOS simulating composite liberation test.

where (Vo)p,T ¼ volume of 1 mol of the oil, ft3. Step 3. Flash the resulting volume (step 2) at temperatures and pressures equal to those in the surface separators and stock tank. Designating the resulting total liberated gas volume (Vg)sc and the volume of the stock-tank oil (Vo)st, calculate the Bo and Rs by applying Eqs. (5.205), (5.206), respectively: Bo ¼ ðVo Þp, T =ðVo Þst   Rs ¼ Vg sc =ðVo Þst

Step 4. Set the cell pressure to a lower pressure level, p, and using the original composition, zi, generate the Ki values through the application of the PR EOS.

580 CHAPTER 5 Equations of State and Phase Equilibria

Step 5. Perform flash calculations based on zi and the calculated ki values. Step 6. To simulate the differential liberation step of removing the equilibrium liberated gas from the cell at constant pressure, simply set the overall composition, zi, equal to the mole fraction of the equilibrium liquid, x, and the bubble point pressure, pb, equal to the new pressure level, p. Mathematically, this step is summarized by the following relationships: zi ¼ x i pb ¼ p

Step 7. Repeat steps 2–6.

Swelling Tests Swelling tests should be performed if a reservoir is to be depleted under gas injection or a dry gas cycling scheme. The swelling test can be conducted on gas-condensate or crude oil samples. The purpose of this laboratory experiment is to determine the degree to which the proposed injection gas will dissolve in the hydrocarbon mixture. The data that can be obtained during a test of this type include n n

The relationship of saturation pressure and volume of gas injected. The volume of the saturated fluid mixture compared to the volume of the original saturated fluid.

The test procedure, as shown schematically in Fig. 5.19, is summarized in the following steps. Step 1. A representative sample of the hydrocarbon mixture with a known overall composition, zi, is charged to a visual PVT cell at saturation pressure and reservoir temperature. The volume at the saturation pressure is recorded and designated (Vsat)orig. Step 2. A predetermined volume of gas with a composition similar to the proposed injection gas is added to the hydrocarbon mixture. The cell is pressured up until only one phase is present. The new saturation pressure and volume are recorded and designated ps and Vsat, respectively. The original saturation volume is used as a reference value and the results are presented as relative total volumes: Vrel ¼ Vsat =ðVsat Þorig

where Vrel ¼ relative total volume (Vsat)orig ¼ original saturation volume.

Simulating of Laboratory PVT Data by Equations of State 581

Swelling test Injection gas Psat

(Vsat)orig

Original fluid (gas or oil) @ saturation pressure Psat

(Psat)new ie Pb or Pd

P

Vt

Hg

Original fluid + injection gas

B

n FIGURE 5.19 Schematic illustration of the Swelling test.

Step 3. Repeat step 2 until the mole percent of the injected gas in the fluid sample reaches a preset value (around 80%). Typical laboratory results of a gas-condensate swelling test with lean gas of a composition are shown Figs. 5.20 and 5.21. These two figures show the gradual increase in the saturation pressure and the swollen volume with each addition of the injection lean gas. As shown in Figs. 5.20 and 5.21, the dew point increased from 3428 to 4880 psig and relative total volume from 1 to 2.5043 after the final addition.

Simulating of the Swelling Test by EOS The simulation procedure of the swelling test by the PR EOS is summarized in the following steps. Step 1. Assume the hydrocarbon system occupies a total volume of 1 bbl at the saturation pressure, ps, and reservoir temperature, T. Calculate the initial moles of the hydrocarbon system from the following expression: ni ¼

5:615psat ZRT

Original fluid + injection gas

Hg

Hg A

(Vsat)new

C

582 CHAPTER 5 Equations of State and Phase Equilibria

5000

Injected gas

4800

New dewpoint pressure

4600

Component

yi

CO2

0

N2

0

4400

C1

0.9468

4200

C2

0.0527

C3

0.0005

4000 3800

Cum. gas inj, scf/bbl

3600

0

1

3428

1.1224

3635

572

1.3542

4015

1523

1.9248

4610

2467

2.5043

4880

3200

0

500

1000

1500

2000

2500

3000

Dewpoint pressure, psig

190

3400

3000

Swollen volume, bbl/bbl

scf/bbl of dewpoint gas

n FIGURE 5.20 Dew point pressure during swelling of reservoir gas with lean gas at 200°F. Relative volume vs. volume of gas injected 2.6

Injected gas

Swelling ratio, (Vsat)new/(Vsat)orig

2.4 2.2 2 1.8

Component

yi

CO2

0

N2

0

C1

0.9468

C2

0.0527

C3

0.0005

1.6 Cum. gas inj, scf/bbl

1.4 1.2 1 0

500

1000

1500

2000

scf/bbl of the dewpoint gas

2500

3000

Swollen volume, bbl/bbl

Dewpoint pressure, psig

0

1

3428

190

1.1224

3635

572

1.3542

4015

1523

1.9248

4610

2467

2.5043

4880

n FIGURE 5.21 Reservoir gas swelling with lean gas at 200°F as expressed in terms.

where ni ¼ initial moles of the hydrocarbon system psat ¼ saturation pressure “pb or pd” depending on the type of the hydrocarbon system, psia Z ¼ gas or liquid compressibility factor Step 2. Given the composition yinj of the proposed injection gas, add a predetermined volume (as measured in scf) of the injection gas to

Simulating of Laboratory PVT Data by Equations of State 583

the original hydrocarbon system and calculate the new overall composition by applying the molal and component balances, respectively: nt ¼ ni + ninj , zi ¼

yinj ninj + ðzsat Þi ni nt

with ninj ¼

Vinj 379:4

where ninj ¼ total moles of the injection gas Vinj ¼ volume of the injection gas, scf (zsat)i ¼ mole fraction of component i in the saturated hydrocarbon system Step 3. Using the new composition, zi, calculate the new saturation pressure, pb or pd, using the Peng-Robinson EOS as outlined in this chapter. Step 4. Calculate the relative total volume (swollen volume) by applying the relationship Vrel ¼

Znt RT 5:615psat

where Vrel ¼ swollen volume psat ¼ saturation pressure. Step 5. Steps 2–4 are repeated until the last predetermined volume of the lean gas is combined with the original hydrocarbon system.

Compositional Gradients Vertical compositional variation within a reservoir is expected during early reservoir life. One might expect the reservoir fluids to have attained equilibrium at maturity due to molecular diffusion and mixing over geological time. However, the diffusive mixing may require many tens of million years to eliminate compositional heterogeneities. When the reservoir is considered mature, it often is assumed that fluids are at equilibrium at uniform temperature and pressure and the thermal and gravitational equilibrium are established, which lead to the uniformity of each component fugacity throughout all the coexisting phases. For a single phase, the uniformity of fugacity is equivalent to the uniformity of concentration.

584 CHAPTER 5 Equations of State and Phase Equilibria

The pressure and temperature, however, are not uniform throughout a reservoir. The temperature increases with depth with a gradient of about 1°F/100 ft. The pressure also changes according to the hydrostatic head of fluid in the formation. Therefore, compositional variations within a reservoir, particularly those with a tall column of fluid, should be expected. Table 5.1 shows the fluid compositional variation of a North Sea reservoir at different depths. Note that the methane concentration decreases from 72.30 to 54.92 mol% over a depth interval of only 81 m, indicating a major change in composition that must be considered when estimating reserves and production planning. In general, the hydrocarbon mixture is expected to get richer in heavier compounds and contains less of light components (such as methane) with increasing depth. The variations in the composition and temperature of fluid with depth will result in changes of saturation pressure as a function of depth. In the oil column, the bubble point pressure generally decreases with

Table 5.1 Variations in Fluid Composition with Depth in a Reservoir Fluid

D, Well 1

C, Well 2

B, Well 2

A, Well 2

Depth (meters subsea) Nitrogen Carbon dioxide Methane Ethane Propane i-Butane n-Butane n-Pentane Hexanes Heptanes Octanes Nonanes Decanes Undecanes-plus Molecular weight Undecanes-plus characteristics Molecular weight Specific gravity

3136 0.65 2.56 72.3 8.79 4.83 0.61 1.79 0.75 0.86 1.13 0.92 0.54 0.28 3.49 33.1

3156 0.59 2.48 64.18 8.85 5.6 0.68 2.07 0.94 1.24 2.14 2.18 1.51 0.91 6 43.6

3181 0.6 2.46 59.12 8.18 5.5 0.66 2.09 1.09 1.49 3.18 2.75 1.88 1.08 9.25 55.4

3217 0.53 2.44 54.92 9.02 6.04 0.74 2.47 1.33 1.71 3.15 2.96 2.03 1.22 10.62 61

260 0.848

267 0.8625

285 0.8722

290 0.8768

Simulating of Laboratory PVT Data by Equations of State 585

7500.

De

wp

oin

o Reserv

tp

7600.

re

ss

ur

e

ir press

7700.

p

g

144

g

7800.

dient

Depth/ft

; psi / ft

ure gra

D

7900.

Sat P

GOC

Res P pd = pb = pi

8000.

8100. 3300

3400

p

o

D

144

3500

o

; psi / ft

3600

t oin -p e ble sur b s Bu pre

3700

3800

3900

Pressure/psia

n FIGURE 5.22 Determination of GOC from pressure gradients.

decreasing methane concentration, whereas the gas condensate dew point pressure increases with increasing heavy fractions, as shown in Fig. 5.22. As illustrated in Fig. 5.22, the GOC is defined as the depth at which the fluid system changes from a mixture with a dew point to a mixture with a bubble point. This may occur at a saturated condition, in which the GOC “gas” is in thermodynamic equilibrium with the GOC “oil” at the gas-oil contact reservoir pressure. By this definition, the following equality applies: pb ¼ pd ¼ pGOC

where pb ¼ bubble point pressure pd ¼ dew point pressure pGOC ¼ reservoir pressure at GOC Fig. 5.23 shows an undersaturated GOC that describes the transition from dew point gas to bubble point oil through a “critical” or “new-critical” mixture with critical temperatures that equal the reservoir temperature. However, the critical pressure of this critical GOC mixture is lower than reservoir pressure.

586 CHAPTER 5 Equations of State and Phase Equilibria

Pressure & temperature in Dewpo t press

Retrograde gas-cap

ure

Depth

re

pres oint blep Bub

on

ure

carb

dro

press

peratu oir tem

f hy

TCo

rvoir

Reserv

Oil rim

Rese

sure

A compositional transition zone that contains a near critical hydrocarbon mixture that separate gas-cap gas from oil

n FIGURE 5.23 Concept of undersaturated GOC; near critical mixture separates gas-cap gas from oil rim.

Minimum Miscibility Pressure Several schemes for predicting MMP based on the equation-of-state approach have been proposed. Benmekki and Mansoori (1988) applied the PR EOS with different mixing rules that were joined with a newly formulated expression for the unlike three-body interactions between the injection gas and the reservoir fluid. Several authors used upward pressure increments to locate the MMP, taking advantage through extrapolation of the fact the K-values approach unity as the limiting coincident tie line is approached. A simplified scheme proposed by Ahmed (1997) is based on the fact that the MMP for a solvent/oil mixture generally is considered to correspond to the critical point at which the K-values for all components converge to unity. The following function, miscibility function, is found to provide accurate miscibility criteria as the solvent/hydrocarbon mixture approaches the critical point. This miscibility function, as defined next, approaches 0 or a very small value as the pressure reaches the MMP: FM ¼ 

" !# n X ½Φðyi Þvi 0 zi 1  ½Φðzi ÞLi i¼1

where FM ¼ miscibility function and Φ ¼ fugacity coefficient. This expression shows that, as the overall composition, zi, changes by gas injection and reaches the critical composition, the miscibility function is

Tuning EOS Parameters 587

monotonically decreasing and approaching 0 or a negative value. The specific steps of applying the miscibility function follow. Step 1. Select a convenient reference volume (eg, one barrel) of the crude oil with a specified initial composition and temperature. Step 2. Combine an incremental volume of the injection gas with the crude oil system and determine the overall composition, zi. Step 3. Calculate the bubble point pressure, pb, of the new composition by applying the PR EOS (as given by Eq. (5.145)) Step 4. Evaluate the miscibility function, FM, at this pressure using the calculated values Φ(zi), Φ(yi), and the overall composition at the bubble point pressure. Step 5. If the value of the miscibility function is small, then the calculated pb is approximately equal to the MMP. Otherwise, repeat steps 2–7.

TUNING EOS PARAMETERS Equations-of-state calculations are frequently burdened by the large number of fractions necessary to describe the hydrocarbon mixture for accurate phase behavior simulation. However, the cost and computer resources required for compositional reservoir simulation increase considerably with the number of components used in describing the reservoir fluid. For example, solving a two-phase compositional simulation problem of a reservoir system that is described by Nc components requires the simultaneous solution of 2Nc + 4 equations for each grid block. Therefore, reducing the number of components through lumping has been a common industry practice in compositional simulation to significantly speed up the simulation process. As discussed previously in this chapter, calculating EOS parameters (ie, ai, bi, and αi) requires component properties such as Tc, pc, and ω. These properties generally are well defined for pure components; however, determining these properties for the heavy fractions and lumped components rely on empirical correlations and the use of mixing rules. Inevitably, these empirical correlations and the selected mixing rules introduce uncertainties and errors in equation-of-state predictions. Tuning the selected equation of state is an important prerequisite to provide meaningful and reliable predictions. Tuning an equation of state refers to adjusting the parameters of the selected EOS to achieve a satisfactory match between the laboratory fluid PVT data and EOS results. The experimental data used should be closely relevant to the reservoir fluid and other recovery processes implemented in the field. The main objective of tuning EOS is to match all the available laboratory data, including

588 CHAPTER 5 Equations of State and Phase Equilibria

n n n n n n n

Saturation pressures; dew point and bubblepoint pressure DE data (Bod, Rsd, Bgd, ρ, etc.) CCE data (Y-function, Co, μ, ρ, etc.) CVD data (LDO, RF, Z-factor, etc.) Swelling tests data Separator tests data MMP

Manual adjustments through trial and error or by an automatic nonlinear regression algorithm are used to adjust the EOS parameters to achieve a match between laboratory and EOS results. The regression variables are based on selecting a number of EOS parameters that may be adjusted or tuned to achieve a match between the experimental and EOS predictions. The following EOS parameters are commonly used as regression variables: n

n

n

Properties of the undefined fractions that include Tc, pc, and ω or alternatively through the adjustments of Ωa and Ωb of these fractions. The Ω modifications should be interpreted as modifications of the critical properties, as they are related by the expressions a ¼ Ωa

R2 Tc2 pc

b ¼ Ωb

RTc pc

Binary interaction coefficient, kij, between the methane and C7+ fractions. When an injection gas contains significant amounts of nonhydrocarbon components, CO2 and N2, the binary interaction, kij, between these fractions and methane also should be modified.

The regression is a nonlinear mathematical model that places global upper and lower limits on each regression variable, ΩaC1 , ΩbC1 ; and so forth. The nonlinear mathematical model determines the optimum values of these regression variables that minimize the objective function, defined as   exp Nexp E  Ecal  X  j j  F¼ Wj     Eexp j j

where Nexp ¼ total number of experimental data points ¼ experimental (laboratory) value of observation j, such as pd, pb, Eexp j or ρob

Tuning EOS Parameters 589

Ecal j ¼ EOC-calculated value of observation j, such as pd, pb, or ρob Wj ¼ weight factor on observation j The values of the weight factors in the objective function are internally set with default or user override values. The default values are 1.0 with the exceptions of values of 40 and 20 for saturation pressure and density, respectively. Note that the calculated value of any observation, such as Ecal j , is a function of the all-regression variables: Ecal j ¼ f ðpC1 ,TC1 ,…Þ

The stepwise regression procedure in tuning EOS parameters is summarized next. Step 1. Split the plus fraction to C35 + or C45 + Step 2. Lump the extended compositional analysis into an optimum number of pseudocomponents; for example, F1, F2, F3: Comp Mol% CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 F1 F2 F3

MW

0.25 44.01 0.88 28.01 23.94 16.04 11.67 30.07 9.36 44.1 1.39 58.12 4.61 58.12 1.5 72.15 2.48 72.15 3.26 84 19.8198 111.46 13.5902 161.57 6.63966 241.13

Sp.gr

Pc

0.818 1071 0.8094 493.1 0.3 666.4 0.3562 706.5 0.507 616 0.5629 527.9 0.584 550.6 0.6247 490.4 0.6312 489.6 0.6781 457.1 0.7542 409.4 0.8061 308 0.8526 229.5

Tc

Vc

ω

Tb

547.9 227.49 343.33 549.92 666.06 734.46 765.62 829.1 845.8 910.1 1047.55 1199.96 1363.62

0.0342 0.0514 0.0991 0.0788 0.0737 0.0724 0.0702 0.0679 0.0675 0.064 0.06272 0.06299 0.06355

0.231 0.0372 0.0105 0.0992 0.1523 0.1852 0.1995 0.228 0.2514 0.2806 0.3653 0.4718 0.6199

351 140 201 333 416 471 491 542 557 607 714.42 864.2 1036.9

Step 3. Match PVT data by regressing on the following EOS parameters:

• • • • •

BIC: “kij” Critical pressure: Pc Critical temperature: Tc Shift volume: C Acentric factor: ω

A summary of the first is described in the following tabulated form with the EOS regression parameters represented by “x”:

590 CHAPTER 5 Equations of State and Phase Equilibria

Component N2 CO2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 F1 F2 F3

kij

pc

Tc

c

ω

x

x

x

x

x

x x x

x x x

x x x

x x x

x x x

Step 4. Having achieved a satisfactory match with the experimental data and obtained the “optimum” values of EOS parameters for C1, F1, F2, F3, and F4, further reduce the number of fractions representing the system by groping N2 with C1 and CO2 with C2. Match the available PVT data by only regression on the parameters of the newly combined groups, as shown below by “x”. Compare results with those obtained in step 2. If satisfactory, proceed to match with further grouping as shown in step 4.

Component

kij

Pc

Tc

c

ω

N2 + C1 CO2 + C2 C3 i-C4 n-C4 i-C5 n-C5 C6 F1 F2 F3

x

x x

x x

x x

x x

x x x

Step 5. Continue the regression process on EOS parameters on each newly lumped groups without adjusting EOS parameters formed in the previous step. This step is summarized in the following two tabulated forms:

Tuning EOS Parameters 591

Component

kij

N2 + C1 CO2 + C2 C3 i-C4 n-C4 i-C5 n-C5 C6 F1 F2 F3

x

Component

kij

N2 + C1 CO2 + C2 C3+ i-C4 n-C4 i-C5 n-C5+ C6 F1 F2 F3

x

Pc

Tc

c

ω

x x

x x

x x

x x

Pc

Tc

c

ω

x x

x x

x x

x x

x x x

x x x

As pointed out in Chapter 4, the Lohrenz, Bray, and Clark (LBC) viscosity correlation has become a standard in compositional reservoir simulation. When observed viscosity data are available, values of the coefficients a1– a5 in Eq. (4.86) and the critical volume of C7+ are adjusted and used as tuning parameters until a match with the experimental data is achieved. These adjustments are performed independent of the process of tuning equation-of-state parameters to match other PVT data. The key relationships of the LBC viscosity model follow: μob ¼ μo + 0

 4 a1 + a2 ρr + a3 ρ2r + a4 ρ3r + a5 ρ4r  0:0001 ξm 1

n C B X C B ½ðxi Mi Vci Þ + xC7 + VC7 + Cρo B A @ i¼1 i62C7 + ρr ¼ Ma

where a1 ¼ 0.1023 a2 ¼ 0.023364 a3 ¼ 0.058533 a4 ¼ 0.040758 a5 ¼ 0.0093324

592 CHAPTER 5 Equations of State and Phase Equilibria

PROBLEMS 1. A hydrocarbon system has the following overall composition: Component

zi

C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.30 0.10 0.05 0.03 0.03 0.02 0.02 0.05 0.40

given System pressure ¼ 2100 psia System temperature ¼ 150°F Specific gravity of C7+ ¼ 0.80 Molecular weight of C7+ ¼ 140 Calculate the equilibrium ratios of this system assuming both an ideal solution and real solution behavior. 2. A well is producing oil and gas with the following compositions at a gasoil ratio of 500 scf/STB:

Component

xi

yi

C1 C2 C3 n-C4 n-C5 C6 C7+

0.35 0.08 0.07 0.06 0.05 0.05 0.34

0.60 0.10 0.10 0.07 0.05 0.05 0.03

given Current reservoir pressure ¼ 3000 psia Bubble point pressure ¼ 2800 psia Reservoir temperature ¼ 120°F Molecular weight of C7+ ¼ 160 Specific gravity of C7+ ¼ 0.823

Problems 593

Calculate the overall composition of the system; that is, zi. 3. The hydrocarbon mixture with the following composition exists in a reservoir at 234°F and 3500 psig: Component

zi

C1 C2 C3 C4 C5 C6 C7+

0.3805 0.0933 0.0885 0.0600 0.0378 0.0356 0.3043

The C7+ has a molecular weight of 200 and a specific gravity of 0.8366. Calculate a. The bubble point pressure of the mixture. b. The compositions of the two phases if the mixture is flashed at 500 psia and 150°F. c. The density of the resulting liquid phase. d. The density of the resulting gas phase. e. The compositions of the two phases if the liquid from the first separator is further flashed at 14.7 psia and 60°F. f. The oil formation volume factor at the bubble point pressure. g. The original gas solubility. h. The oil viscosity at the bubble point pressure. 4. A crude oil exists in a reservoir at its bubble point pressure of 2520 psig and a temperature of 180°F. The oil has the following composition: Component

xi

CO2 N2 C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

0.0044 0.0045 0.3505 0.0464 0.0246 0.0683 0.0083 0.0080 0.0080 0.0546 0.4224

The molecular weight and specific gravity of C7+ are 225 and 0.8364. The reservoir initially contains 122 MMbbl of oil. The surface facilities

594 CHAPTER 5 Equations of State and Phase Equilibria

consist of two separation stages connected in series. The first separation stage operates at 500 psig and 100°F. The second stage operates under standard conditions. a. Characterize C7+ in terms of its critical properties, boiling point, and acentric factor. b. Calculate the initial oil in place in STB. c. Calculate the standard cubic feet of gas initially in solution. d. Calculate the composition of the free gas and the composition of the remaining oil at 2495 psig, assuming the overall composition of the system remains constant. 5. A pure n-butane exists in the two-phase region at 120°F. Calculate the density of the coexisting phase using the following equations of state: a. van der Waals. b. Redlich-Kwong. c. Soave-Redlich-Kwong. d. Peng-Robinson. 6. A crude oil system with the following composition exists at its bubble point pressure of 3250 psia and 155°F.

Component

xi

C1 C2 C3 C4 C5 C6 C7+

0.42 0.08 0.06 0.02 0.01 0.04 0.37

If the molecular weight and specific gravity of the heptanes-plus fraction are 225 and 0.823, calculate the density of the crude using the a. Standing-Katz density correlation. b. Alani-Kennedy density correlation. c. van der Waals correlation. d. Redlich-Kwong EOS. e. Soave-Redlich-Kwong correlation. f. Soave-Redlich-Kwong correlation with the shift parameter. g. Peneloux volume correction. h. Peng-Robinson EOS.

References 595

7. The following reservoir oil composition is available: Component

Mol% Oil

C1 C2 C3 i-C4 n-C4 i-C5 n-C5 C6 C7+

40.0 8.0 5.0 1.0 3.0 1.0 1.5 6.0 34.5

The system exists at 150°F. Estimate the MMP if the oil is to be displaced by a. Methane. b. Nitrogen. c. CO2. d. 90% CO2 and 10% N2. e. 90% CO2 and 10% C1. 8. The compositional gradients of the Nameless field are as follows: C1 + N2

C7+

C2–C6

TVD, ft

69.57 67.44 65.61 63.04 62.11 61.23 58.87 54.08

6.89 8.58 10.16 12.54 13.44 14.32 16.75 22.02

23.54 23.98 24.23 24.42 24.45 24.45 24.38 23.9

4500 4640 4700 4740 4750 4760 4800 5000

Estimate the location of the GOC.

REFERENCES Ahmed, T., 1991. A practical equation of state. (translation). SPE Reserv. Eng. 6 (1), 137–146. Ahmed, T., 1997. A generalized methodology for minimum miscibility pressure. In: Paper SPE 39034 presented at the Latin American Petroleum Conference, Rio de Janeiro, Brazil, August 30-September 3, 1997.

596 CHAPTER 5 Equations of State and Phase Equilibria

Benmekki, E., Mansoori, G., 1988. Minimum miscibility pressure prediction with EOS. SPE Reserv. Eng. Brinkman, F.H., Sicking, J.N., 1960. Equilibrium ratios for reservoir studies. Trans. AIME 219, 313–319. Campbell, J.M., 1976. Gas conditioning and processing. In: Campbell Petroleum Series, vol. 1. Campbell, Norman, OK. Dodson, C., Goodwill, D., Mayer, E., 1953. Application of laboratory PVT data to engineering problems. J. Pet. Technol. 5 (12), 287–298. Dykstra, H., Mueller, T.D., 1965. Calculation of phase composition and properties for lean-or enriched-gas drive. Soc. Petrol. Eng. J. 5 (3), 239–246. Edmister, W., Lee, B., 1986. Applied Hydrocarbon Thermodynamics, In: second ed. vol. 1. Gulf, Houston, TX, p. pp. 52. Elliot, J., Daubert, T., 1985. Revised procedure for phase equilibrium calculations with soave equation of state. Ind. Eng. Chem. Process. Des. Dev. 23, 743–748. Graboski, M.S., Daubert, T.E., 1978. A modified soave equation of state for phase equilibrium calculations, 1. Hydrocarbon system. Ind. Eng. Chem. Process. Des. Dev. 17, 443–448. Hoffmann, A.E., Crump, J.S., Hocott, R.C., 1953. Equilibrium constants for a gascondensate system. Trans. AIME 198, 1–10. Jhaveri, B.S., Youngren, G.K., 1984. Three-parameter modification of the Peng-Robinson equation of sate to improve volumetric predictions. In: Paper SPE 13118, presented at the SPE Annual Technical Conference, Houston, September 16–19, 1984. Katz, D.L., Hachmuth, K.H., 1937. Vaporization equilibrium constants in a crude oilnatural gas system. Ind. Eng. Chem. 29, 1072. Katz, D., et al., 1959. Handbook of Natural Gas Engineering. McGraw-Hill, New York. Kehn, D.M., 1964. Rapid analysis of condensate systems by chromatography. J. Pet. Technol. 16 (4), 435–440. Lim, D., Ahmed, T., 1984. Calculation of liquid dropout for systems containing water. In: Paper SPE 13094, presented at the 59th Annual Technical Conference of the SPE, Houston, September 16–19, 1984. Lohrenze, J., Clark, G., Francis, R., 1963. A compositional material balance for combination drive reservoirs. J. Pet. Technol. Nikos, V., et al., 1986. Phase behavior of systems comprising north sea reservoir fluids and injection gases. J. Pet. Technol. 38 (11), 1221–1233. Pedersen, K., Thomassen, P., Fredenslund, A., 1989. Characterization of gas condensate mixtures. In: Advances in Thermodynamics. Taylor and Francis, New York. Peneloux, A., Rauzy, E., Freze, R., 1982. A consistent correlation for Redlich-KwongSoave volumes. Fluid Phase Equilib. 8, 7–23. Peng, D., Robinson, D., 1976a. A new two constant equation of state. Ind. Eng. Chem. Fundam. 15 (1), 59–64. Peng, D., Robinson, D., 1976b. Two and three phase equilibrium calculations for systems containing water. Can. J. Chem. Eng. 54, 595–598. Peng, D., Robinson, D., 1978. The Characterization of the Heptanes and Their Fractions: Research Report 28. Gas Producers Association, Tulsa, OK. Peng, D., Robinson, D., 1980. Two and three phase equilibrium calculations for coal gasification and related processes. In: Thermodynamics of Aqueous Systems with Industrial Applications. In: ACS Symposium Series no. 133.

References 597

Redlich, O., Kwong, J., 1949. On the thermodynamics of solutions. An equation of state. Fugacities of gaseous solutions. Chem. Rev. 44, 233–247. Reid, R., Prausnitz, J.M., Sherwood, T., 1987. The Properties of Gases and Liquids, fourth ed. McGraw-Hill, New York. Sim, W.J., Daubert, T.E., 1980. Prediction of vapor-liquid equilibria of undefined mixtures. Ind. Eng. Chem. Process. Des. Dev. 19 (3), 380–393. Slot-Petersen, C., 1987. A systematic and consistent approach to determine binary interaction coefficients for the Peng-Robinson equation of state. In: Paper SPE 16941, presented at the 62nd Annual Technical Conference of the SPE, Dallas, September 27–30, 1987. Soave, G., 1972. Equilibrium constants from a modified Redlich-Kwong equation of state. Chem. Eng. Sci. 27, 1197–1203. Standing, M.B., 1977. Volumetric and Phase Behavior of Oil Field Hydrocarbon Systems. Society of Petroleum Engineers of AIME, Dallas, TX. Standing, M.B., 1979. A set of equations for computing equilibrium ratios of a crude oil/ natural gas system at pressures below 1,000 psia. J. Pet. Technol. 31 (9), 1193–1195. Stryjek, R., Vera, J.H., 1986. PRSV: an improvement Peng-Robinson equation of state for pure compounds and mixtures. Can. J. Chem. Eng. 64, 323–333. Van der Waals, J.D., 1873. On the Continuity of the Liquid and Gaseous State. Sijthoff, Leiden, The Netherlands (Ph.D. dissertation). Vidal, J., Daubert, T., 1978. Equations of state—reworking the old forms. Chem. Eng. Sci. 33, 787–791. Whiston, C., Brule´, M., 2000. “Phase Behavior.” Monograph volume 20, Society of Petroleum Engineers. Richardson, TX. Whitson, C.H., Brule, M.R., 2000. Phase Behavior. SPE, Richardson, TX. Whitson, C.H., Torp, S.B., 1981. Evaluating constant volume depletion data. In: Paper SPE 10067, presented at the SPE 56th Annual Fall Technical Conference, San Antonio, October 5–7, 1981. Wilson, G., 1968. A modified Redlich-Kwong EOS, application to general physical data calculations. In: Paper 15C, presented at the Annual AIChE National Meeting, Cleveland, May 4–7, 1968. Winn, F.W., 1954. Simplified nomographic presentation, hydrocarbon vapor-liquid equilibria. Chem. Eng. Prog. Symp. Ser. 33 (6), 131–135.