• Author / Uploaded
• imaghost13

Citation preview

1

EQUATIONS OF STATE When the effects of complex phase behavior cannot be accurately calculated using simple

approaches, it is often desirable to utilize some sort of equation of state (EOS). An EOS approach is often recommended when dealing with volatile oils and retrograde condensate gases. The two most common equations of state utilized in petroleum engineering applications are the Peng-Robinson (PR) and Soave-Redlich-Kwong (SRK) equations, which historically were derived from van der Waals’ (vdW) equation. These three equations are called “cubic” because they result in a cubic representation for the molar volume. The basic equations are the following: Ideal Gas pvm = RT van der Waals  ac   p +  ( vm − b ) = RT vm   Soave-Redlich-Kwong  acα ( T )  p+  ( vm − b ) = RT vm ( vm + b )   Peng-Robinson   acα ( T ) p+  ( vm − b ) = RT vm ( vm + b ) + b ( vm − b )   The parameters ac , α ( T ) , and b are empirically determined from experimental data (for pure components, critical temperature and pressure, and a specified point on the vapor pressure curve), α ( T ) being a function of temperature and having a value of one at the critical temperature. Note that the parameters do not have the same value in each equation. Equations of State (MAM 04.Feb.2000)

1

We will focus on the two equations of most interest, the SRK and PR EOSs, first writing the equations in the following common form.   acα ( T ) p+  ( v − b ) = RT ( vm + f1b ) ( vm + f 2b )  m  Using the following substitutions z=

A=

B=

pvm (definition of gas deviation factor) RT acα ( T ) p R 2T 2

(dimensionless)

bp (dimensionless) RT

Equation can be algebraically rewritten in the following form: z 3 − a2 z 2 + a1 z − a0 = 0 where a0 = B  A + wB ( B + 1) 

a1 = A − B  B ( u − w ) + u  a2 = 1 − B ( u − 1) f1 = f2 = u = f1 + f 2 w = f1 f 2

σ = u 2 − 4w

SRK 1 0 1 0 1

PR 1+ 2 1− 2 2 -1 2 2

Note that even though Eq. is written in terms of z , Eq. can always be used to determine the molar volume if desired. Equations of State (MAM 04.Feb.2000)

2

The values of ac and b are determined by noting that at the critical point, ∂p ∂vm

=0 Tc

∂2 p =0 ∂vm2 T c

Since α ( Tc ) = 1 , evaluating these expressions at the critical point results in the following expressions ac = Ω a

b = Ωb

R 2Tc2 pc

RTc pc

With the following values determined from the various EOSs.

Ωa =

SRK 0.427480

PR 0.457236

Ωb =

0.086640

0.077796

zc =

0.333333

0.307401

Note that the value of the z-factor at the critical point is also given. The numbers above can be expressed analytically, as given at the end of this document. For both the SRK and PR equations, the temperature-dependent parameter is expressed as

(

)

α ( T ) = 1 + m 1 − Tr   

2

where, m = mo + m1ω − m2ω 2 Equations of State (MAM 04.Feb.2000)

3

Values of the coefficients in Eq. can be obtained from the following table.

m0 =

SRK 0.480

PR 0.37464

m1 = m2 =

1.574 0.176

1.54226 0.26992

The parameter ω is called the Pitzer accentric factor and is empirically determined from the actual vapor pressure curve for a substance according to the following relationship:

(

)

ω = − log pvr T =0.7 + 1 r

where, pvr T =0.7 = r

pv ( T = 0.7Tc ) pc

and pv is vapor pressure. The values of ω are typically determined from published tables, along with critical temperature and pressure. It should be noted that these cubic forms of the equation of state can apply to liquids as well as gases. Even though z has traditionally been used for gases, there is nothing restricting this to be the case, since all of the quantities in Eq. are defined for liquids as well as gases. Of course the values of z for liquids will be smaller than those for gases since liquids are denser. To address the more complex question of “phase”, the EOSs must be extended to include the concept of chemical potential, sometimes called the Gibbs molar free energy ( G ), changes in which are defined by the following relationship. dG = vm dp For an ideal gas, it can be shown that

Equations of State (MAM 04.Feb.2000)

4

dG = RT d ( ln p ) For non-ideal fluids, a property called “fugacity” ( f )is defined having the following two properties: dG = RT d ( ln f ) lim f = p p →0

Note that the fugacity can be thought of as the potential for transfer between phases. This means that molecules will tend to move phase-wise, so as to minimize the fugacity. Further defining a fugacity coefficient, ϕ as  p z −1  f ϕ = = exp  ∫ dp  p  0 p  results in the following expression 1 ln ϕ = z − 1 − ln z + RT

vm

 RT  − p   dvm ∫ v  ∞ m

The terms in Eq. can, of course, be evaluated using an EOS, resulting in the following: u +σ  z+B A  2 ln ϕ = z − 1 − ln ( z − B ) − ln  u − σ Bσ  z + B  2

2

    

SOLUTIONS In general an EOS can be solved for either one of the three parameters z (or vm ), p , or

T . First we’ll deal with what’s called a “flash” calculation, solving for z when pressure and temperature are known.

Equations of State (MAM 04.Feb.2000)

5

Equation is a cubic equation. In general the following possibilities are possible for the roots (solutions) of this equation: 1. Three real roots, some or all possibly repeated. 2. One real root and two imaginary roots. For pure components, away from the vapor pressure line, there is generally one real root and two imaginary roots. The value of z will be relatively large for gases and relatively small for liquids. Near the vapor pressure line, however, there will be three real roots. It is in this region that fugacities are needed to determine phase. When there are three real roots, the largest root represents the z-factor of the gaseous phase, while the smallest root represents the z-factor of the liquid phase. The middle root is nonphysical. A determination of which phase is present can be made by comparing the fugacities for both phases. The phase with the lowest fugacity is the one present. On the vapor pressure line, of course, there are two phases simultaneously present. When two phases are present concurrently, this means that their fugacities must be equal, since molecules can be in either phase. In fact, the way vapor pressure is usually determined at a given pressure is to iterate on the above equations, seeking the pressure that results in equal liquid and gas fugacities. For pure components, then, the procedure for determining the z-factor (and thus the specific molar volume and/or density) is the following: 1. For a given pressure and temperature, determine the coefficients in Eq. . 2. Solve this equation for the largest and smallest real roots using a cubic root solver. 3. Calculate both liquid and gaseous fugacity coefficients (if there are two) using Eq. . Use the z-factor from the lowest fugacity. 4. Calculate densities, ρ =

Equations of State (MAM 04.Feb.2000)

zRT pM , and/or molar volumes, vm = , as desired. p zRT

6

If it is desired to find a vapor pressure at a given temperature, iterate on the above procedure, until a pressure is selected so that the fugacities calculated in Step 3 are equal. The Excel Solver works quite well for this.

3

MIXTURES For mixtures an additional consideration must be taken into account, composition. With

so-called flash calculations, the composition of the total system (mixture) is usually taken as a “known”. The mole fractions of each component are given the symbols z j (not to be confused with the z-factor), where j is the component index, spanning the total number of components in the system. Likewise the mole fractions of the liquid phase are given the symbols x j and the mole fractions of the gaseous phase y j . By definition,

∑x =∑y =∑z j

j

j

j

j

j

= 1.

Often compositions of the liquid and gaseous phases are expressed in terms of equilibrium rations, K j , defined by

Kj =

yj xj

Note that K j will be small (but never zero) for components that prefer to be in the liquid phase, and much greater than one (but never infinite) for components that prefer to be in the gaseous phase. The values of K j can be determined from correlations, but the more modern approach is to determine them through EOSs. The final composition variable needed is the total mole fraction of the mixture that is in the gaseous phase, n% g . This variable is, of course, not defined outside the two-phase envelope, and has a value ranging from zero (at a bubble point) to one (at a dew point) inside the envelope. The z j , K j , and n% g then fully define the composition of a two-phase mixture. Liquid and gas compositions can be determined by molar balances using the following relationships. Equations of State (MAM 04.Feb.2000)

7

xj =

zj

1 + n% g ( K j − 1)

yj = K jxj =

K jzj

1 + n% g ( K j − 1)

Note that a bubble points, x j = z j , while at dew points, y j = z j . Dealing with mixtures also requires that chemical potentials and fugacities must be calculated for each component, defined by the following relationships. dG j = RT d ( ln f j ) lim f j = y j p p →0

Equation states that fugacities must approach partial pressures (ideal behavior) as pressure approaches zero. At equilibrium, each component’s fugacity must be the same in both phases, i.e., f gj = f lj for all components. We can also define the fugacity coefficient for each component as

ϕj =

fj yj p

With these definitions, it then turns out that Kj =

ϕlj ϕ gj

To apply the EOSs for mixtures, effective mixture coefficients for Eq. must be determined using some sort of “mixing rules”. For both the PR and SRK EOSs, the following mixing rules are used to obtain the effective mixture values of b and acα .

Equations of State (MAM 04.Feb.2000)

8

acα = ∑ yi aciα i ∑ y j acjα j ( 1 − δ ij ) i

j

b = ∑ y jbj j

Note that for a single component, these two equations will yield one components values of acα and b . Also, evn though these two equations are based on gas mixtures, they can also be applied to liquid mixtures by substituting the x j for the y j . Equation is a double summation, with both indexes i and j going over the entire number of components in the system. The δ ij are called “binary interaction coefficients” and are empirical measures of the attractive and repulsive forces between molecules of unlike size. Note that δ ii = δ jj = 0 and

δ ij = δ ji = 0 . There are many ways that binary interaction coefficients are characterized (Ahmed, 1989), but in general they increase as the relative difference between molecular weights increase. When no data is available, values of zero are sometimes use. Often binary interaction coefficients are used to “history match” EOS calculations against actual PVT experiments. EOS calculations of the fugacity coefficients for each phase are then made using the following equation. u +σ  z+B 2 A ln ϕ j = ( z − 1) B′j − ln ( z − B ) − ( A′j − B′j ) ln  u −σ Bσ z+B  2

    

where,

A′j =

B′j =

2 acjα j ∑ yi aciα i ( 1 − δ ij ) i

acα

bj b

Equations of State (MAM 04.Feb.2000)

9

Note that Eq. reverts to Eq. where there is only one component present. For two-phase mixtures, Eq. must be solved twice, once for the liquid phase and once for the gas phase. When two phases exist, both EOS calculations will yield three real roots. The smallest root of the liquid equation should be taken as the z-factor for the liquid phase, and the largest root of the gas equation should be taken as the gas z-factor. The procedure for doing a flash calculation on a mixture the is the following. 1. Calculate the acjα j and b j for each component. These are functions of pressure and temperature, but not composition. 2. Guess values for the K j . Various correlations are available for estimating these values (e.g., McCain, 1990). 3. Find a value of n% g which maintains molar balances, i.e., that ensures

∑x =∑y j

j

j

j

= 1 (more on this later).

4. Use the x j in the mixing equations to find the coefficients and solve the EOS for the liquid. Select the smallest root. Likewise use the y j to solve the EOS for the gas phase. Select the largest root. 5. Calculate fugacity coefficients for each component in each phase using Eq. . Use these to determine the K j from Eq. . 6. Repeat from Step 2, using the calculated K j as the new trial values. Stop when the calculated values are near the trial values. McCain (1990) suggests the following possible error calculation to determine convergence. Error = ∑

(K

T j

− K Cj

)

2

K Tj K Cj

Iterations should be stopped when the Eq. results in a value less than some specified tolerance. Step 3 is usually done with a Newton-Raphson type iteration to find the value of n% g . The following procedure is typically used.

Equations of State (MAM 04.Feb.2000)

10

When the z j and K j are known (as in Step 3 above), the correct value of n% g is the one that ensures molar balances. Although either the x j or y j equation may be used for this purpose, here we will focus on the x j equation, Eq. . Using Newton-Raphson iteration to solve this equation results in the following. f ( n% g ) = 1− ∑ j

f ′ ( n% g) =∑ j

zj

1 + n% g ( K j − 1) z j ( K j − 1)

1 + n%  g ( K j − 1)  

2

Recall that Newton-Raphson iterative involves successive guessing, with the “new” guess calculated from the “old” one by

k +1 g

n% = n% − k g

( ) f ′ ( n%) k f n% g k g

One of the roots of Eq. is always n% g = 0 , so the initial guess should start well away from this value. Outside the two-phase envelope, the value of n% g is undefined and may take on nonphysical values. The following relationships define how to determine whether the calculation is being done inside the two-phase envelope or not. Liquid phase

∑z K

Gaseous phase

∑K

j

zj

j

Bubble point

j

1

j

The above procedure can also be used to calculate bubble and dew point pressures, by simply iterating on pressure until the appropriate Eq. or is true. Again, the Excel Solver is a good way to do this type of calculation.

4

EXACT VALUES OF EOS PARAMETERS

van der Waals ac = b=

27 R 2Tc2 R 2Tc2 = 0.421875 64 pc pc

1 RTc RT = 0.125 c 8 pc pc

zc =

3 = 0.375 8

Soave-Redlich-Kwong  1 + 21 3 + 22 3  R 2Tc2 R 2Tc2 ac =  = 0.427480  9 pc2   pc  21 3 − 1  RTc RT b= = 0.086640 c  pc  3  pc zc =

1 = 0.333333 3

Peng-Robinson

Equations of State (MAM 04.Feb.2000)

12

 1715  −112 + 80 2 − 16 2 − 13   1 2205 1134 ac = + + 23 13 1536  16 2 − 13 16 2 − 13  13  − 162 16 2 − 13 + 45 16 2 − 13 

(

)

(

(

)

)

(

343  16 2 − 35 + 16 2 − 13 1  b= 1008 1536  + + 144 16 2 − 13 13  16 2 − 13 

(

)

 1  7 zc = 11 − 32  16 2 − 13 

(

5

(

)

13

(

+ 16 2 − 13

)

13

)

)

    R 2Tc2 R 2Tc2 = 0.457236  pc  pc  23  

   RTc = 0.077796 RTc 13  p pc  c 

  = 0.307401  

REFERENCES

Ahmed, T.: Hydrocarbon Phase Behavior, Gulf Publishing Co., Houston (1989). McCain, W.D., Jr.: The Properties of Petroleum Fluids, PennWell Publishing Co., Tulsa (1990).

Equations of State (MAM 04.Feb.2000)

13