Ind. Eng. Chem. Res. 1999, 38, 2485-2490
2485
GENERAL RESEARCH A Method To Predict Equilibrium Conditions of Gas Hydrate Formation in Porous Media Matthew A. Clarke, Mehran Pooladi-Darvish, and P. Raj Bishnoi* Department of Chemical & Petroleum Engineering, University of Calgary, 2500 University Drive NW, Calgary, Alberta, Canada T2N 1N4
A predictive model has been developed to determine the incipient hydrate formation conditions in porous media. The only additional information that is needed to determine the incipient hydrate formation conditions is the pore radius, surface energy per unit area, and wetting angle. It was found that the model performed well in predicting the experimental data of Handa and Stupin (J. Phys. Chem. 1992, 96, 8599). Introduction Gas hydrates are crystalline compounds of water that belong to a group of solids called clathrates. They are formed from mixtures of water and low molar mass gases at high pressures and low temperatures. Through hydrogen bonding, water molecules form a framework containing relatively large cavities that can be occupied by certain gas molecules, which stabilize the structure due to van der Waals forces. The hydrate-forming gases include light alkanes (methane to isobutane), carbon dioxide, hydrogen sulfide, nitrogen, oxygen, xenon, krypton, and chlorine. Gas hydrates are known to occur in one of three crystal structures: structures I, II,2,3 and H.4,5 Structures I and II consist of two types of cavities, and structure H consists of three types of cavities. Both methane and ethane form hydrates of structure I, a unit crystal structure being body-centered cubic, consisting of two types of cavities, pentagonal dodecahedron (small cavity) and tetrakaidecahedron (large cavity). In the petroleum industry, it is desirable to avoid the formation of gas hydrates. When gas hydrates form, they tend to agglomerate and block pipelines and process equipment. However, naturally occurring gas hydrates that form in the permafrost region or in deep oceans represent a vast untouched natural gas reserve.5,6 Although the exact amount of gas in the hydrate form is not known, it is believed to be comparable to the known amount of gas in the free state. Numerous methods for the recovery of natural gas from hydrate fields have been proposed.3,7-9 These techniques include thermal decomposition, depressuriation, and chemcial injection. To fully exploit hydrate reserves, it will be necessary to know the decomposition/formation conditions of the gas hydrate in porous media. The equilibrium conditions for hydrate formation in free water are well established, as are the methods to calculate them. Sir Humphry Davis made the first experimental observation of gas hydrates in 1811. Carson and Katz10 wrote the first publication stressing * To whom correspondence should be addressed. Tel: (403) 220-6695. Fax: (403) 284-4852. E-mail:
[email protected].
that hydrates behave as solid solutions. They were also able to study the four-phase equilibrium of gas mixtures in the presence of gas hydrates and of liquids rich in hydrocarbons. Katz was also the first person to present graphical and tabulated methods to predict the temperature and pressure at which hydrates form. Ng and Robinson11 investigated the conditions for hydrate formation in liquid hydrocarbon-water systems. Englezos and Bishnoi12 studied the effect of electrolytes and alcohols on the equilibrium formation conditions of gas hydrates. Experimental data from various sources are confirmed and reported by Sloan.5 The statistical thermodynamic model of van der Waals and Platteeuw13 is most commonly used to describe the chemical potential of water in the hydrate phase. The statistical model is based on the threedimensional generalization of localized adsorption. van der Waals derived relationships for the establishment of equilibrium conditions. The theory of Lennard-Jones and Devonshire made possible the determination of a partition function for the molecules in the cavities. Hence, it became possible to determine expressions for the equilibrium vapor pressure, the chemical potential, and the hydration number. Parrish and Prausnitz14 presented the first computational procedure for hydrate equilibrium that could be implemented in a computer program. Their model used the theory of van der Waals and Platteeuw13 along with the Kihara spherical potential function. This paper illustrated the usefulness of statistical thermodynamics in a practical industrial situation. Holder et al.15 simplified the model of Parrish and Prausnitz14 by eliminating the reference hydrate. Englezos and Bishnoi12 presented a predictive model to determine the incipient hydrate formation conditions in the presence of electrolytes. While there are plenty of data available for gas hydrate formation in free water, there are only a handful of studies that have dealt with gas hydrate equilibria in porous media. Handa and Stupin1 presented experimentally determined conditions of methane and propane hydrate equilibrium in 70-Å-radius silica gel pores. For both methane and propane, it was
10.1021/ie980625u CCC: $18.00 © 1999 American Chemical Society Published on Web 05/06/1999
2486 Ind. Eng. Chem. Res., Vol. 38, No. 6, 1999
observed that the presence of the porous silica gel caused hydrates to form at pressures between 20 and 70% higher than what are needed to form hydrates in free water. Yousif and Sloan9 related the permeability of Berea sandstone to the additional pressure for the onset of hydrates. Melnikov and Nesterov16 derived and applied a model to predict the equilibrium conditions in porous media. However, Nesterov and Melnikov 16 did not present predictions for hydrate formation in submicronsized pores. The conditions of hydrate decomposition in porous media will depend on the rock and fluid properties, such as the wetting angle and the pore radius. When gas hydrates form in porous media, it is necessary to account for the additional forces that result from interactions with the media, mainly the capillary forces. The effect of the capillary forces is to lower the activity of water in the pore. This, in turn, causes a depression of the freezing point of water in the pore. In the present study, a predictive method that incorporates the properties of the porous media with the model of van der Waals and Platteeuw13 to describe the chemical potential of water in the hydrate phase and an appropriate equation of state to describe the vapor phase is presented. Mathematical Model Hydrate Formation at a Plane Surface. For threephase vapor-liquid-hydrate equilibria, the basic equations for the equilibrium condition in gas hydrates are
Lo MT-Lo µMT w - µw ) ∆µw
The term µLwo is the chemical potential of pure water at the system (p, T), with the pressure being that abovethe flat surface. The right-hand side of eq 5 is commonly represented by15 o ∆µow ∆µMT-L w ) + RTf RTo
∆µow µMT w ) + RTf RTo
(2)
RTf
MT-Lo p∆νw
∫p
o
)
∆µow RTo
+
∫p
p
RTf
o ∆νMT-L w
o
dp -
RTf
dp -
µLwo
-
RTf
In eq 3, fj is the fugacity of hydrate former j, other than water, in the hydrate phase. The isofugacity criteria for a hydrate former in the vapor phase and in the hydrate phase, implied by eq 2, are usually incorporated by taking the fugacity in eq 3 as being that calculated from the equation of state for the vapor phase. Because the concentration of water in the vapor phase is usually negligible, the equilibrium relation for water becomes
(4)
The chemical potential difference of water in the empty hydrate lattice and that in the pure liquid state at the (p, T) of the system is
o
RT2
MT-Lo Tf∆hw
∫T
o
RT2
∫T
o ∆hMT-L w
Tf o
dT
µLwo dT + RTf
RT2
dT +
νm ln(1 + ∑Cmjfj) ∑ m j
Lo µsol w ) µw + RTf ln aw
(8)
(9)
where the definition of the activity is
aw ) fw/fow (3)
∫T
The chemical potential of water in the liquid solution phase is
NH
m)1
L µH w ) µw
dp -
Equations 3 and 7 can be combined to give the chemical potential of water in the hydrate phase.
V (i ) 1, NH) µH i ) µi
Cmjfj) ∑ νm ln(1 + ∑ j)1
RTf
o
(7)
µH w
MT µH w ) µw - RTf
∫p
MT-Lo Tf∆hw
The first term on the right-hand side is the chemical potential difference at the reference point (To, po). The reference point is usually taken as the melting temperature of ice and zero pressure. The second term on the right-hand side represents the pressure correction from the reference pressure to the formation pressure, at the formation temperature, and the third term is the correction from the reference temperature to the formation temperature, at the reference pressure. Equations 5 and 6 can be combined to give an expression for the chemical potential of water in the empty lattice
(1)
2
MT-Lo p∆νw
(6)
µLi ) µVi (i ) 1, N)
where N ) total number of components and NH ) number of hydrate-forming components. The fugacity, or chemical potential, of a component in the vapor or liquid phase may be calculated using a suitable equation of state. In the present study, the Trebble-Bishnoi equation of state17 is used. The model of van der Waals and Platteeuw13 is generally used to calculate the chemical potential (or fugacity) of water in the hydrate phase. The model is given by
(5)
(10)
The term fow is the fugacity of pure water in the standard state. The standard state used in eq 10 has to be consistent with that used in eqs 8 and 9 for the chemical potential of pure water in the standard state, µLwo. At equilibrium, the chemical potential of water in the hydrate has to be equal to that in the solution phase. Thus, at equilibrium, eqs 8 and 9 give
∆µow RTo
+
∫p
o ∆νMT-L w
p o
RTf
dp -
∫T
Tf
o ∆hMT-L w
dT RT2 νm ln(1 + Cmjfj) ) ln aw (11)
∑ m
o
∑j
Hydrate Formation in Porous Media. When hydrates are formed in free liquid water, or free ice, then it is possible to neglect the surface effects on the
Ind. Eng. Chem. Res., Vol. 38, No. 6, 1999 2487
the cross-sectional area of the cylindrical pore.
∆p ) pg - pl )
Figure 1. Fluid in a capillary (adapted from Adamson18).
equilibrium conditions. However, when hydrates are forming in small capillaries, it is not possible to neglect the surface effects on the equilibrium conditions. When surface effects are not neglected, the differential of the Gibbs free energy becomes
dG ) -S dT + V dP + σ dAs +
∑i µi dni
σ ) (dG/dAs)T,P,ni
[
∆Fgh ) σ
]
y′ y′′ + (1 + y′)3/2 x(1 + y′2)1/2
[
]
p dx + ∫ (1 +x dp p2)3/2 (1 + p2)1/2
[
W ) 2πσ
]
xp (1 + p2)1/2
x)r,p)tan φ
fw
fow
)
νl (p - pg) RT l
(20)
Substituting eq 19 into eq 20 gives the activity of water as
ln aw ) -
2σνl cos θ rRT
(21)
If the surface energy between the hydrate and the water is neglected, then the hydrate has to form at the same pressure as the liquid. This may be justifiable because we are seeking the incipient condition, i.e., the point at which the amount of hydrate is zero. Thus, in porous media, the upper limit of integration in eq 11 should be the liquid side pressure, pl. Also, because we are calculating the incipient point, the pressure used to evaluate the fugacity of the hydrate former in the hydrate phase can be taken as being the same as that in the continuous gas phase. This pressure is what would be measured experimentally. Equation 21 is substituted into eq 11 to give the equilibrium relationship for hydrate formation in porous media.
∆µow RTo
+
∫p
pl
o ∆νMT-L w
RTf
o
∑ m
dp -
νm ln(1 +
∫T
∑j
Tf o
o ∆hMT-L w
RT2
dT -
Cmjφjyjpg) )
-2σνl rRTf
cos θ (22)
(14)
The total weight W of the column of liquid in the capillary may be obtained exactly from the above equation. When the substitutions p ) y′ and y′′ ) pp′ are made, the above equation may be written as follows:18
W ) 2πσ
()
ln aw ) ln
(13)
The major implication of eq 12 is that at equilibrium the pressure in the liquid phase will not be equal to that in the vapor phase. In a capillary system, the pressure will be greater on the side of the interface that contains the center of curvature of the curved surface. The exact treatment of capillary rise must take into account the deviation of the meniscus from sphericity. That is, the curvature must correspond to the ∆p ) ∆Fgy at each point on the meniscus, where y is the distance above the flat liquid surface. Figure 1, which is adapted from Adamson,18 shows a fluid in a capillary of radius r. The angle between the capillary wall and the meniscus is the wetting angle, θ. The formal statement of the condition is obtained by writing the Young and Laplace equation for a general point (x, y) on the meniscus. It is still assumed that the capillary is circular in cross section so that the meniscus shape is that of a figure of revolution. The differential equation becomes (the prime and double prime indicate the first and second derivatives, respectively)
(19)
When the wetting angle is 90°, eq 19 shows that the pressure is the same in both phases. This corresponds to a flat surface. For the present study, it is assumed that the wetting angle is 0° and the surface energy per unit area is taken to be the same as that of water on silica gel, 72 mJ/m2.19,20 The standard state used in eqs 8 and 9 is the formation temperature and the pressure above the flat surface. Therefore, the activity of water can be written as
(12)
where σ, the surface energy per unit area (or surface tension), is defined as
2σ cos θ r
Because the pressure in the liquid phase is related to that in the vapor phase by eq 19, eq 22 only contains one unknown, the vapor-phase pressure. The only additional information that is required to determine the incipient conditions in porous media is the mean pore size, the surface energy, and the wetting angle. Results and Discussion
(15) (16)
x)0,p)0
φ ) π/2 - θ
(17)
W ) 2πrσ cos θ
(18)
The pressure difference between the gas and liquid phases is the weight of the column of liquid divided by
The computational procedure used in the current work is outlined in Figure 2. The pressure is updated by using the secant method to determine the pressure that satisfies eq 22. In the cases that the secant method fails, the bisection method is used. The program was written in C++. The results of the calculations are shown for methane and propane in Figures 3 and 4, respectively. The calculated values are compared to the experimental values reported by Handa and Stupin.1
2488 Ind. Eng. Chem. Res., Vol. 38, No. 6, 1999
Figure 2. Computational flow diagram.
The data of Handa and Stupin1 show a gradual change in slope in the p-T curve. This is in contrast to the data for hydrate formation in free water, where the slope of the p-T curve changes abruptly at the melting point of water. The gradual change in slope of the experimental data is due to the fact that in porous materials water melts over a temperature range rather than at a definite point. This is due to the pore size distribution. One implication of melting in a porous medium is that the quadruple point, the point at which ice, liquid water, hydrate, and vapor are in equilibrium, is no longer a unique point. In the present work, it was assumed that the water remained in the liquid state at temperatures above the quadruple point. The quadruple point was taken to be the average melting point, as
Figure 3. Methane hydrate formation conditions.
determined by Handa and Stupin.1 This temperature, which is 267.5 K, corresponds to the melting temperature in the average size pore.1 The effect of this is to move to the left the point at which the slope of the p-T curve changes slope. This is analogous to what was done by Englezos and Bishnoi12 for electrolytes. The predictions for methane are seen to match the data well. The maximum deviation between the experimental pressure and the calculated pressure is 16%. This point occurs at a temperature of 268.2 K. At a temperature of 276.2 K, the percent difference between the experimental pressure and the calculated pressure is only 0.95%. The larger deviation in the intermediate temperature range is due to the fact that the pore ice melts over a temperature range. As a result, at a given temperature between approximately 265 and 268 K, there could be vapor-liquid water-hydrate equilibrium in the smaller pores and vapor-ice-hydrate equilibrium in the larger pores. The predictions for the propane hydrate only match the data qualitatively, with the match being best at the extreme ends. For propane, the maximum percent deviation between the experimental and the calculated pressure is 29%. This occurs at a temperature of 268.8 K. The temperatures at which the deviation is greatest fall in the range over which the pore ice is melting. In this range it is difficult to characterize the equilibrium as being hydrate-ice-gas or hydrate-liquid watergas. This is because at a given temperature it is possible to have ice in larger pores while still having liquid water in the smaller pores. It is noted that the model predictions cross the experimental data at approximately 271 K. The predictions that are presented are based upon the average pore size in the silica gel, which is 70 Å. Thus, in both Figures 3 and 4, the model predicts a sharp change in slope at the melting point of water in the 70-Å pores, which is approximately 267.5 K.1 This will only occur in porous media if all of the pores are of uniform
Ind. Eng. Chem. Res., Vol. 38, No. 6, 1999 2489
Figure 4. Propane hydrate formation conditions. Table 1. Comparison of Computed ∆p Values Using the Proposed Model and the Model of Yousif and Sloan9 ∆p (kPa) core sample A B C
porosity
permeability (md)
Yousif and Sloan
proposed model
0.1977 0.188 0.223
83.88 96.89 394.7
64 70 38
79 88 68
size and shape. In reality, each pore will contribute to the shape of the p-T curve. Thus, the total behavior should account for the pore size distribution. Yousif and Sloan9 assumed a bundle of equal diameter tubes for the porous media to relate the rock permeability to an equivalent pore diameter and derived an expression for the increase in the incipient hydrate formation pressure due to porous media.
∆p ) 3.2639σxφ/k
(23)
r ) x8k/φ
(24)
Equation 23 requires permeability in md (milli-darcy) and surface energy per unit area in dyn/cm, whereas eq 24 requires the permeability in m2. In their work, the additional pressure to form hydrates in three different Berea sandstone cores was calculated. It was experimentally determined that an additional pressure of approximately 78 kPa was required to form methane hydrates in the sandstone at 273.7 K. The proposed model was also used to calculate the additional pressure required to form hydrates in the sandstone cores. Table 1 compares the additional pressure calculated from eq 23 to that calculated using the proposed model. The correlation of Yousif and Sloan was also used to calculate the additional pressure required to form hydrates in silica gel. Handa and Stupin1 reported the pore volume to be 1.11 cm3/g. When the pore volume and the skeletal density are known,21 the porosity can
be calculated as 0.6875.22 The permeability is calculated from eq 24 and is found to be 0.004 21 md. Thus, eq 23 would predict a constant additional pressure of 20.5 MPa to form hydrates. This is much larger than what is actually required to form hydrates in silica gel. The model of Yousif and Sloan9 appears to be suited only to relatively large pores (pore radius > 10-7 m). Calculations using the proposed model show that the contribution of the large pores to the equilibrium pressure is small. However, as the pore sizes decrease into the submicron size, their contribution becomes larger. Conclusion A model has been developed to predict the incipient hydrate formation conditions in porous media. The only additional information that is needed is the surface energy, mean pore radius, and wetting angle. The proposed model was able to predict the experimental data of Handa and Stupin1 reasonably well, with maximum deviations of 15% for methane and 29% for propane. The match between the experimental data and the predicted data may be improved upon if the pore size distribution is taken into account. Acknowledgment Financial support was provided by Natural Sciences and Engineering Research Council of Canada (NSERC). Notation As ) interfacial contact area (m2) Cij ) Langmuir constant of species i in cavity of type j (Pa-1) fi ) fugacity of species i (Pa) G ) Gibbs free energy (J/mol) h ) enthalpy (J/mol) k ) permeability (md or m2) p ) pressure (Pa)
2490 Ind. Eng. Chem. Res., Vol. 38, No. 6, 1999 ∆p ) increase in the incipient pressure due to porous media (Pa) r ) pore radius (m) R ) universal gas constant (J/mol‚K) T ) temperature (K) v ) volume (m3) Vm ) molar volume (m3/mol) yi ) vapor-phase mole fraction Greek Letters θ ) wetting angle (rad) µ ) chemical potential (J/mol) σ ) surface energy per unit area (J/m2) φ ) porosity φ ) fugacity coefficient ν ) fractional occupancy Subscripts and Superscripts f ) formation conditions H ) hydrate l ) liquid MT ) empty lattice o ) reference conditions Sol ) solution w ) water
Literature Cited (1) Handa, Y. P.; Stupin, D. Thermodynamic Properties and Dissociation Characteristics of Methane and Propane Hydrates in 70-Å-Radius Silica Gel Pores. J. Phys. Chem. 1992, 96, 8599. (2) von Stackelberg, M.; Mu¨ller, H. R. Feste Gas Hydrate II. Z. Elektrochem. 1954, 58, 25. (3) Ripmeester, J. A.; Ratcliffe, C. I.; Tse, J. S. The Nuclear Magnetic Resonance of 129Xe Trapped in Clathrates and Some Other Solids. J. Chem. Soc. Faraday Trans. 1988, 84, 3731. (4) Makogan, Y. F. Hydrates of Natural Gas, Translated from Russian by W. J. Cieslewicz; Penn Well Books: Tulsa, OK, 1981. (5) Sloan, E. D. Clathrate Hydrates of Natural Gases, 2nd ed.; Marcel Dekker Inc.: New York, 1998. (6) Holder, G. D.; Malone, R. D.; Lawson, W. F. Effects of Gas Composition and Geothermal Properties on the Thickness and Depth of the Natural Hydrate Zones. JPT, J. Pet. Technol. 1987, 39, 1147. (7) Holder, G. D.; Kamath, V. A.; Godbole, S. P. Potential of Natural Gas Hydrates as an Energy Resource. Annu. Rev. Energy 1984, 9, 427.
(8) Jamaluddin, A. K. M.; Kalogerakis, N.; Bishnoi, P. R. Modelling of Decomposition of a Synthetic Core of Methane Gas Hydrate by Coupling Intrinsic Kinetics with Heat Transfer Rates. Can. J. Chem. Eng. 1989, 67, 948. (9) Yousif, M. H.; Sloan, E. D. Experimental investigation of Hydrates Formation and Dissociation in Consolidated Porous Media. SPE Reservoir Eng. 1991, 25, 452. (10) Carson, D. B.; Katz, D. M. Gas Hydrates. Petr. Tr. AIME 1942, 146, 150. (11) Ng, H. J.; Robinson, D. B. The Measurement and Prediction of Hydrate Formation in Liquid Hydrocarbon-Water System. Ind. Eng. Chem. Fundam. 1976, 15, 59. (12) Englezos, P.; Bishnoi, P. R. Prediction of Gas Hydrate Formation Conditions in Aqueous Electrolyte Solutions. AIChE J. 1998, 34, 1718. (13) van der Waals, J. H.; Platteeuw, J. C. Clathrate Solutions. Adv. Chem. Phys. 1959, 2 (1), 1. (14) Parrish, W. R.; Prausnitz, J. M. Dissociation Pressures of Gas Hydrates formed by Gas Mixtures. Ind. Eng. Chem. Process Des. Dev. 1972, 11, 26. (15) Holder, G. D.; Corbin, G.; Papadopoulos, K. D. Thermodynamic and Molecular Properties of Gas Hydrates Containing Methane, Argon and Krypton. Ind. Eng. Chem. Fundam. 1980, 19, 282. (16) Melnikov, V.; Nesterov, A. Modelling of Gas Hydrates Formation in Porous Media. Proceedings of the 2nd International Conference on Natural Gas Hydrates, Toulouse, France, 1996; INPENSIGC: France, 1996; p 541. (17) Trebble, M. A.; Bishnoi, P. R. Extension of the TrebbleBishnoi Equation of State to Fluid Mixtures. Fluid Phase Equilibr. 1987, 40, 1. (18) Adamson, A. The Physical Chemistry of Surfaces, 2nd ed.; John Wiley and Sons: New York, 1967. (19) Brinker, C. J.; Scherer, G. W. The Physics and Chemistry of Sol-Gel Processing; Academic Press: New York, 1990. (20) Vigil, G.; Xu, Z.; Steinberg, S.; Israelachvili, J. Interactions of Silica Surfaces. J. Colloid. Interface Sci. 1994, 165, 367. (21) Desphande, R.; Hua, D.; Smith, D. M.; Brinker, C. J. Pore structure Evaluation in Silica Gel During Aging/Drying III. J. NonCryst. Solids 1992, 144, 32. (22) Smith, J. M. Chemical Engineering Kinetics; McGrawHill: New York, 1981.
Received for review September 30, 1998 Revised manuscript received March 1, 1999 Accepted March 4, 1999 IE980625U