- Email: [email protected]
Nucleation of CO2-hydrate in a porous medium
Fluid Phase Equilibria 200 (2002) 263–275
Nucleation of CO2 -hydrate in a porous medium O.Ye. Zatsepina, B.A. Buffett∗ Department of Earth and Ocean Sciences, University of British Columbia, Vancouver, BC, Canada V6T 1Z4 Accepted 19 April 2001
Abstract Experiments designed to crystallize gas hydrate from dissolved CO2 in natural porous media are used to study nucleation under varying thermodynamic conditions. We recover quantitative information from these experiments using a stochastic model for the nucleation process. Estimates of the model parameters are used to determine the average time for nucleation as a function of temperature and composition. © 2002 Elsevier Science B.V. All rights reserved. Keywords: Gas hydrate; Model; Nucleation time; Interfacial tension; Electrical conductivity; Porous media
1. Introduction Experiments on hydrate formation from aqueous solution (e.g. [1]) show that hydrate begins forming at thermodynamic conditions which are different from the conditions of equilibria. This suggests that the kinetic process of nucleation initiates hydrate formation and that equilibrium conditions do not resolve the question whether hydrate can form at a given value of pressure, temperature, and solute concentration. Nucleation of a hydrate crystal requires an excess energy to create a nucleus surface. Since a driving force is necessary to overcome a nucleation barrier, hydrate nucleates from the solution which is cooled below the equilibrium temperature. The thermodynamic force that drives nucleation can be related to the difference in the chemical potential of the hydrate components (e.g. water and hydrate-forming gas) in the liquid and hydrate phases at the existing thermodynamic conditions. Once hydrate begins to form and draw solute from the surrounding liquid, the difference in the chemical potential of the components in the two phases decrease. Hydrate nuclei stop growing when the driving force disappears. In the other words, the solute is drawn from solution during hydrate formation until its concentration is lowered to the equilibrium value. Fig. 1 shows how the difference in the chemical potentials of water in the liquid and hydrate phases varies with temperature when carbon dioxide is the hydrate-forming gas. The chemical potentials µl of water in the liquid phase are calculated using the Trebble–Bishnoi equation of state [2], while the chemical ∗
Corresponding author. Tel.: +1-604-822-2267; fax: +1-604-822-6047. E-mail address: [email protected] (B.A. Buffett). 0378-3812/02/$ – see front matter © 2002 Elsevier Science B.V. All rights reserved. PII: S 0 3 7 8 - 3 8 1 2 ( 0 2 ) 0 0 0 3 2 - 8
264
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 1. Chemical potential µ/RT of water in CO2 -hydrate (dashed line) and in liquid solution (solid line). The standard state corresponds to pure liquid water at pressure and temperature of interest. Pressure is 3 MPa and mole fraction of dissolved CO2 is X = 0.020. The equilibrium temperature Teq (P , X) is defined by the intersection of the two curves. The driving force for nucleation is proportional to the difference µhl = µh − µl . Stable nuclei of hydrate appear when µhl is sufficient to overcome the energy surface of the nuclei.
potentials µh of water in the hydrate phase are described by the thermodynamic model of van der Waals and Platteeuw[3]. The solid line shows µl /RT, while the dashed line indicates µh /RT (both potentials are referred to the value for pure liquid water at the pressure and temperature of interest). As Fig. 1 shows, the energy difference becomes large enough to overcome the surface energy barrier and launch the process of nucleation at sufficiently low temperature. In this paper, we use experiments to study the timescale, over which mixtures of carbon dioxide and water evolve toward their equilibrium state. Quantifying the rate of hydrate nucleation is a difficult experimental problem. First, we must be able to detect the appearance of hydrate nuclei. Second, we need to observe a large number of nucleation events in order to obtain a meaningful statistical average for the nucleation rate. Finally, we need to control the thermodynamic conditions of the experiments, which include pressure, temperature, and concentration of a solute. An experimental setup that deals with each of these challenges is described in [4]. The initial goal of this experiment was to crystallize gas hydrate from aqueous solution in the absence of bulk gas [5]. To achieve this goal, hydrate formation occurs in a water-saturated porous medium. The medium effectively isolates the aqueous solution in it from the bulk gas phase, which is used to maintain pressure in the experiment (pressurized hydrate-forming gas is supplied to the top of the pressure cell). The time required to equilibrate the pore liquid with the vapour phase depends on thickness of the porous medium because transport of a solute through the porous medium occurs solely by diffusion. When the porous medium is more than a few centimeters thick, the concentration of the solute in the interior of the porous medium remains nearly constant over a 2 days experiment.
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
265
In this experimental setup, hydrate formation can not be visually detected. Instead, we use electrical resistance measurements to monitor nucleation and growth of hydrate. We choose CO2 as the hydrate-forming gas mainly because it supplies ions when dissolved in water to alter electrical conductivity of the solution. The liquid conductivity, determined by concentration and mobility of these ions, can be expressed in terms of the solute concentration X and temperature T [4]. The conductivity increases with X because the abundance of the solute determines the number of ions, while temperature primarily alters the ion mobility. The measurements of electrical resistance R are inversely proportional to the pore liquid conductivity σ (X, T ) and can be represented by R(X, T ) =
A σ (X, T )
(1)
where A is a factor, which depends on porosity and geometric factors of the electrode configuration [4]. As hydrate nuclei appear in the pore fluid, they consume the solute, which alters the conductivity of solution. We infer these changes in the solute concentration by measuring electrical resistance. Since the resistance measurements are made over a large volume of the porous medium, we average over many nucleation sites and obtain a more reliable estimate of the average nucleation rate. Data collected in the experiment provide quantitative information about nucleation. To process this information and recover nucleation characteristics, we develop a simple theoretical model of hydrate formation in porous media. Estimates of the nucleation time are extracted from our experimental data and the results are extrapolated to other thermodynamic conditions.
2. Experimental results Fig. 2 presents measurements of resistance and temperature as functions of time collected in an experiment. The measurements show that changes in resistance reflect changes in temperature very closely. However, resistance continues to increase with time, even when temperature reaches a nearly constant value of T ≈ −1.5 ◦ C. (Freezing of water occurs at −1.7 ◦ C due to the presence of the solute.) This transition in the behavior of resistance versus temperature becomes more evident when resistance is plotted as a function of temperature (see Fig. 3). Fig. 3 demonstrates the steady increase in measurements of resistance with cooling from 15 ◦ C to approximately 0 ◦ C. Prior to hydrate formation, the solute concentration is constant and the changes in resistance reflect the effects of temperature on ion mobility. We find that these measurements are well approximated by assuming that the conductivity σ (X, T ) of the pore liquid varies linearly with T . This allows us to recover the temperature dependence of ion mobility using measurements of resistance and temperature. When hydrate forms in the pore liquid, the solute is transferred from solution into growing hydrate crystals. Increases in resistance due to depletion of solution with the solute concentration can be identified in Fig. 3 when temperature decreases from 0 to −1.5 ◦ C. The resistance increases sharply and continues to increase even when temperature is held constant at −1.5 ◦ C. We do not cool the apparatus below −1.5 ◦ C to prevent freezing of the water. We subsequently increase temperature and hydrate starts dissociating in order to re-establish higher equilibrium values of the solute concentration [6]. Measurements from two subsequent cycles of warming and cooling show that changes in resistance with temperature follow the same path with almost no offset [4]. This suggests that equilibrium states in our system are established on a timescale which is less than that of the controlled changes
266
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 2. Experimental measurements of resistance (A) and temperature (B) as a function of time from a typical experiment at P = 2.8 MPa. Measurements collected during cooling are indicated by the line, whereas the points represent individual measurements during warming. A small change in the slope of R(T ) near T = 4 ◦ C coincides with the disappearance of hydrate from solution during warming.
in temperature. On the basis of these results, we assume that changes in resistance with temperature during warming mirror changes in the solubility as the system adjusts through a sequence of equilibrium states. Correspondingly, the offset between cooling and warming paths in Fig. 3 is representative of kinetics of the process. The data suggest that hydrate starts forming at T ≈ 0 ◦ C. On the other hand, measurements of resistance at temperature above 4 ◦ C nearly coincide with measurements of resistance at the same temperature during the initial cooling. This implies that hydrate totally disappears from solution at
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
267
Fig. 3. Resistance as a function of temperature for the experimental data shown in Fig. 2.
T = 4 ◦ C and that the equilibrium temperature for aqueous solution and hydrate phases at the initial solute concentration is Teq = 4 ◦ C. It follows that the solution was cooled almost 4 ◦ C below the equilibrium temperature to initiate nucleation of hydrate. We infer nucleation characteristics from the time dependence of resistance measurements during the interval when temperature is nearly constant. Fig. 4 reproduces the measurements in Fig. 2, but corrects the values of resistance for the influence of temperature on mobility of ions. This procedure is not necessary because we are going to analyze values of resistance at constant temperature anyway. However, the corrected values of resistance reflect the effect of changes in the solute concentration and more clearly represent the kinetic aspects of the process. Fig. 4 shows that the resistance does not change much within the first 5 h of the experiment because the solute concentration X is constant prior to hydrate formation. As hydrate crystallizes from solution and X decreases, values of the resistance increase. This increase in the resistance is well fit to an exponential function of time t/τ , where the timescale τ is 4 h. We have several other experiments conducted at slightly different thermodynamic conditions. Changes in the resistance with time show the same pattern as that shown in Fig. 4. The exponential fit yields different values of timescale τ for different thermodynamic conditions. To understand these results and to extrapolate them to other thermodynamic situations, we develop a theoretical model for hydrate formation in porous media. An essential element of this model is the relationship between the number of pores with hydrate nuclei and the average electrical resistance of the medium.
268
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 4. Experimental measurements of resistance from Fig. 2, corrected for the effects of temperature. Resistance is nearly constant prior to the start of hydrate formation. Subsequent increases in resistance are attributed to the nucleation and growth of hydrate in the porous medium. A fit of (7) to the resistance measurements yields a nucleation time of τ = 4 h.
3. Theoretical modeling We consider the appearance of hydrate nuclei in a large number of pores. Nucleation does not occur in all pores at once as temperature decreases. Instead, the nuclei appear randomly. Under uniform thermodynamic conditions in identical pores, there is an equal probability of forming a nucleus in each of them. We also assume that the pores are chemically isolated. This means that once a nucleus forms, it equilibrates with the pore liquid without affecting the neighbouring pores. Under this approximation, nucleation of hydrate in pores is analogous to the freezing of water drops [7]. However, gas hydrate involves two components, so the details of the kinetic process as well as the thermodynamics of the bulk phases are different. Despite these differences, the freezing of water drops still serves as a useful conceptual model to illustrate the stochastic nature of nucleation. The probability of forming a nucleus in a pore of volume V in a time increment dt is VJ dt, where J is rate of nucleation. If N(t) denotes the number of pores without nuclei at time t, then dN = −N(t)VJ dt
(2)
For the special case of constant J , the number Nh of pores with hydrate nuclei is given by Nh (t) = N0 (1 − e−t/τ )
(3)
where N0 is the initial number of empty pores, τ = (JV)−1 is the timescale of nucleation, determined by the nucleation rate J and the pore volume V . The nucleation rate J can be generally represented by an Arrhenius equation − G J = J0 exp (4) kT where J0 is the pre-exponential factor and G is the activation energy.
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
An equation of the same form − G∗ J = J0 exp kT
269
(5)
holds for classical nucleation. In the theory of classical nucleation, G∗ is the increase in the Gibbs free energy due to nucleation of a stable embryo, which can be represented as [8]
G∗ =
16πvh2 γhl3 3 µ2hl
(6)
where vh is the molar volume of hydrate, µhl = µh − µl is the difference in chemical potential of water in the hydrate and liquid phase at a given thermodynamic condition. Surface tension between hydrate and water γhl can be considered as a lumped parameter that embodies the influences of pore walls and other possible complications on nucleation [8,9]. In classical nucleation, the rate of nucleation depends on local thermodynamic conditions through the value of µhl in G∗ . Our calculations show that G∗ is strongly dependent on T . This suggests that the activation energy G in (4) may have a similar dependence. Because G is expected to vary with temperature, we can not treat it as a constant parameter in the numerical model. Instead, we adopt classical nucleation and treat J0 and γhl as constant parameters to recover from the experimental data. The activation energy at experiment conditions can then be inferred from (6). 4. Numerical simulation We simulate the electrical resistance of the porous medium using a network of resistors. The value of resistance for each resistor is determined by the conductivity of the pore liquid, while the network defines the pathways for electrical current. Consequently, the resistance of the pore liquid changes when a nucleus forms and this alters the flow of current through the network. Since the electrical current is inversely proportional to bulk resistance, changes in bulk resistance can, in principle, be related to changes in the number of pores with nuclei. We establish this relationship below. The porous medium is considered as a large number of pores distributed in a three-dimensional cubic network. The whole network is subdivided into elementary cubes. The relationship between the porous medium and the network is shown in Fig. 5. Resistors are associated with the pores, where the value of resistance is inversely proportional to the conductivity of the pore liquid. The pores are connected by pore throats which intersect at the vertices of the cube. When voltage is applied, the network bonds (elementary cube sides) serve as pathways for electrical current. The flow of current through the circuit is affected by the distribution of resistances in the network. As a hydrate nucleus appears in a pore, the conductivity of the pore liquid changes. (We assume that the solute concentration evolves to its equilibrium value at once as hydrate nuclei appear in the pore liquid.) This alters the resistor distribution and causes a change in bulk resistance. We compute the bulk resistance by defining a system of algebraic equations based on conservation of charge at the vertices of the cube [10]. The algebraic equations are solved for the current along individual pathways using a relaxation techinque. The bulk resistance is evaluated using the net current across the network and the applied voltage. Each pore has an identical value of resistance at the start of the simulation. Prior to hydrate nucleation, concentration X0 of the solute is constant. Cooling causes a uniform increase in resistance due to the
270
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 5. (A) Schematic illustration of pore and pore throats in the porous medium. (B) Correspondence between the porous medium and the network resistor model.
temperature dependence of the liquid conductivity. Once hydrate nuclei begin to form in the pores, their subsequent growth depletes the pore liquid to the equilibrium concentration Xeq . All other pores remain at their initial concentration X0 , under the assumption that the individual pores are chemically isolated. A further complication is introduced by the fact that Xeq depends on temperature. We account for this in the simulation by assuming that pores with nuclei maintain the equilibrium concentration Xeq (T ) as temperature changes. The probability of forming a nucleus in a volume V within a time increment t is given by JV t. J is the rate defined in (5), where J0 and γhl are free parameters in the model. We use a coarse-grained sand as the porous medium to avoid influence of the pore size on thermodynamic conditions of hydrate formation and adopt a nominal pore volume of V = 10−9 m−3 [11]. The stochastic nature of nucleation is simulated using a sequence of random numbers from a uniform distribution between (0, 1). The probability of nucleation in a pore within a given time interval t is equal
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
271
Fig. 6. Predictions of hydrate nucleation in a porous medium at constant temperature. The stochastic model is illustrated using three different values of temperature (0, 1, and 1.5 ◦ C). (A) Predicted number Nh of pores with nuclei over time in a three-dimensional network of 400 pores. (B) Predicted change in the bulk resistance R due to the formation of hydrate nuclei in the pores. Initial and final resistances in the three simulations are different because of the temperature dependence of R.
to the probability that the random number lies between 0 and JV t. When the random number falls in this range, we assume that nucleation has occurred. This procedure is repeated at each pore in the network for each successive increment of time. The resistances of pores with nuclei are adjusted to represent local equilibrium and the resulting change in the bulk resistance is computed using the network model. Calculations confirm that changes in the bulk resistance R do reflect changes in the number of pores with nuclei Nh . Fig. 6 shows several examples of R and Nh as a function of time from simulations at constant temperature. The curves are slightly rough because of the stochastic nature of the simulation, but there is very good agreement in the trends. In fact, a plot of R versus Nh reveals a nearly linear dependence. Thus, R is expected to follow the trend in Nh (t) very closely. Fig. 7 shows how closely the change in Nh (t) from the simulation at constant temperature follows the simple exponential form in (3). In this particular example, we set the nucleation time τ using a representative estimate of the nucleation rate J . We can also recover an estimate of τ by fitting (3) to the output of the simulation. The good agreement between the input and recovered τ suggests that the
272
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 7. A comparison of the number of pores with nuclei between the simulation at constant temperature and the analytical prediction in (3). The nucleation time used in the simulation is recovered using a fit of (3) to the simulation output.
nucleation time can be estimated from the experimental data by assuming that the predicted change in bulk resistance can be approximated by R = R(X0 , T ) + [R(Xeq , T ) − R(X0 , T )][1 − e−t/τ ]
(7)
where R(X0 , T ) is the initial resistance before any nuclei form and R(Xeq , T ) the resistance that would be measured if all of the pore were filled with nuclei and the concencentration of dissolved CO2 was everywhere equal to the equilibrium value Xeq . The timescale is given by τ = (VJ)−1 . We use this approach in analyzing the experimental data. 5. Recovering nucleation characteristics Fig. 4 shows that the temperature-corrected resistance measurements are nearly constant during the first 5 h of the experiment. After this time the resistance increases sharply due to the onset of hydrate nucleation and growth. After 7 h, the temperature reaches a nearly steady value, while the resistance continues to increase. A fit of the exponential form in (7) to the data after 7 h yields the smooth curve shown in Fig. 4. The nucleation time τ recovered from the fit is 4 h at T = −1.5 ◦ C. In order to interpret τ we need to evaluate the thermodynamic conditions that prevailed prior to hydrate nucleation. Temperature is measured, while pressure is maintained at P = 2.8 MPa and the mole fraction X 0 of the solute in the porous medium is constant. We inferred the equilibrium temperature Teq = 4 ◦ C for hydrate stability from our experimental data (see discussion of Fig. 3). The solubility Xeq (T ) within the hydrate stability zone is known from previous experiments [4], so we take Xeq (T = 4 ◦ C) = X0 ≈ 0.02. Because we estimate τ at T = −1.5 ◦ C, while T eq = 4 ◦ C, the cooling T below the equilibrium temperature is T = 5.5 ◦ C. Similar results are obtained from two other experiments at slightly different pressure and temperature conditions. The solute mole fraction is nearly constant at X0 = 0.02 in all three experiments. (Small differences in X0 are inferred from measurable differences in T eq in each experiment.) A summary of results is listed in Table 1. Differences in the thermodynamic conditions for each experiment have been
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
273
Table 1 Summary of experiments P (MPa)
T (◦ C)
X (mole fraction)
Teq (◦ C)
T (◦ C)
τ (h)
2.8 2.7 2.6
−1.5 0.8 1.2
0.0197 0.0202 0.0201
4.0 5.5 5.2
5.5 4.7 4.0
4.0 ± 0.5 4.8 ± 0.5 5.9 ± 0.5
characterized in terms of T . The recovered values of τ vary systematically with T , although the range of T is not large. We explore the range of values for J0 and γhl that yields predictions for τ which are consistent with the experiments. To do so, we fit experimental estimates of τ to the prediction τ = (JV)−1 , where J is given by (5). Thermodynamic calculations, based on the technique described in [6], for µhl in (6) are carried out for the experimental conditions listed in Table 1. Using the three estimates of τ from Table 1, we obtain three equations for J0 and γhl . The best fitting values are J0 = 9 × 104 m3 s−1 and γhl = 0.0015 J m−2 . Fig. 8 shows a good agreement between the theoretical prediction and the experimental estimates of τ .
Fig. 8. Experimental estimates of τ (circles) and theoretical prediction (line) calculated using J0 = 9 × 104 m3 s−1 and σ = 0.0015 J m−2 .
274
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 9. Nucleation time τ for CO2 -hydrate as a function of temperature T . Pressure is fixed at P = 2.8 MPa and two values of the solute mole fraction are considered: X1 = 0.020 and X2 = 0.018. The equilibrium temperatures for the two concentrations are equal to Teq (P , X1 ) = 4 ◦ C and Teq (P , X2 ) = 1.5 ◦ C. These values are based on an experimentally-derived phase diagram [4].
Fig. 9 presents the temperature dependence of nucleation time τ predicted with the recovered values of J0 and γhl . To demonstrate the influence of the solute concentration on nucleation time, estimates of τ are given at X = 0.020 and 0.018 mole fraction. The plots in Fig. 9 indicate that τ increases with temperature and becomes infinite when T approaches Teq (X, P ). The plots also show that at lower temperatures, the nucleation time becomes relatively insensitive to T and more dependent on X. The larger nucleation time corresponds to the lower solute concentration. The dependence of τ on X arises in several ways. First, the value of X is important in the calculation of µhl . Second, we expect X to influence J0 by altering the number of elementary hydrate building units per unit volume of solution [12–15]. For simplicity, we assume that J0 is linearly dependent on X and adopt J0 = 8.1 × 104 m3 s−1 at X = 0.018. This approach is used to extrapolate our estimates of τ to other solute concentrations. 6. Conclusions We design an experiment to crystallize gas hydrate from dissolved CO2 in natural porous media. Experimental data collected in a series of experiments provide quantitative information on nucleation of CO2 -hydrate. We use a stochastic model for nucleation of hydrate in porous media to analyze these data and estimate the nucleation timescale under different thermodynamic conditions. The nucleation
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
275
time is extrapolated to other values of temperature and concentration of the solute using this simple model. Fitting the model to the experimental data yields estimates of the surface tension and the kinetic prefactor. The surface tension γhl , which we recover from the experimental data, is about 18 times smaller than the surface tension between ice and water (γhl = 0.0276 J m−2 ). The latter is often taken as a reasonable approximation for the surface tension between hydrate and water. The lower value recovered in our experiments implies that the effective surface energy of the nuclei is smaller than that predicted for homogeneous nucleation. This is consistent with our expectations that the pore walls and other foreign particles serve as nucleation sites and promote nucleation. In fact, calculations using γhl = 0.0267 yields the nucleation time that exceeds the age of the Earth when T = 5 ◦ C. Clearly, the porous medium plays an important role in hydrate nucleation. Acknowledgements We thank two anonymous reviewers for helpful comments and suggestions. This work is supported by the Institute of Applied Energy and by NSERC. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15]
E.D. Sloan, Clathrate Hydrates of Natural Gases, 2nd Edition, Marcel Dekker, New York, 1998. M.A. Trebble, P.R. Bishnoi, Fluid Phase Equilib. 40 (1988) 1–21. J.H. van der Waals, J.C. Platteeuw, Adv. Chem. Phys. 2 (1959) 1–57. O.Ye. Zatsepina, B.A. Buffett, Fluid Phase Equilib. 192 (2001) 85–102. B.A. Buffett, O.Ye. Zatsepina, Marine Geology 164 (2000) 66–77. O.Ye. Zatsepina, B.A. Buffett, Geophys. Res. Lett. 24 (1997) 1567–1570. T.W. Barlow, A.D.J. Haymet, Natural gas hydrates, Ann. New York Acad. Sci. 715 (1994) 549–551. P.V. Hobbs, Ice Physics, Clarendon Press, Oxford, 1974. K.F. Kelton, Solid State Phys. 45 (1991) 75–177. R.J. Suman, R.J. Knight, Geophysics 62 (1997) 1151–1162. V. Melnikov, V. Nesterov, in: Proceedings of the 2nd International Conference on Natural Gas Hydrates, Toulouse, France, 1996, pp. 541–549. D. Turnbull, J.C. Fisher, J. Chem. Phys. 17 (1949) 71–73. F.F. Abraham, Homogeneous Nucleation Theory, Academic Press, New York, 1974. P.M. Rodger, J. Phys. Chem. 94 (1990) 6080–6089. R.L. Christensen, E.D. Sloan, Natural gas hydrates, Ann. New York Acad. Sci. 715 (1994) 283–305.
Nucleation of CO2 -hydrate in a porous medium O.Ye. Zatsepina, B.A. Buffett∗ Department of Earth and Ocean Sciences, University of British Columbia, Vancouver, BC, Canada V6T 1Z4 Accepted 19 April 2001
Abstract Experiments designed to crystallize gas hydrate from dissolved CO2 in natural porous media are used to study nucleation under varying thermodynamic conditions. We recover quantitative information from these experiments using a stochastic model for the nucleation process. Estimates of the model parameters are used to determine the average time for nucleation as a function of temperature and composition. © 2002 Elsevier Science B.V. All rights reserved. Keywords: Gas hydrate; Model; Nucleation time; Interfacial tension; Electrical conductivity; Porous media
1. Introduction Experiments on hydrate formation from aqueous solution (e.g. [1]) show that hydrate begins forming at thermodynamic conditions which are different from the conditions of equilibria. This suggests that the kinetic process of nucleation initiates hydrate formation and that equilibrium conditions do not resolve the question whether hydrate can form at a given value of pressure, temperature, and solute concentration. Nucleation of a hydrate crystal requires an excess energy to create a nucleus surface. Since a driving force is necessary to overcome a nucleation barrier, hydrate nucleates from the solution which is cooled below the equilibrium temperature. The thermodynamic force that drives nucleation can be related to the difference in the chemical potential of the hydrate components (e.g. water and hydrate-forming gas) in the liquid and hydrate phases at the existing thermodynamic conditions. Once hydrate begins to form and draw solute from the surrounding liquid, the difference in the chemical potential of the components in the two phases decrease. Hydrate nuclei stop growing when the driving force disappears. In the other words, the solute is drawn from solution during hydrate formation until its concentration is lowered to the equilibrium value. Fig. 1 shows how the difference in the chemical potentials of water in the liquid and hydrate phases varies with temperature when carbon dioxide is the hydrate-forming gas. The chemical potentials µl of water in the liquid phase are calculated using the Trebble–Bishnoi equation of state [2], while the chemical ∗
Corresponding author. Tel.: +1-604-822-2267; fax: +1-604-822-6047. E-mail address: [email protected] (B.A. Buffett). 0378-3812/02/$ – see front matter © 2002 Elsevier Science B.V. All rights reserved. PII: S 0 3 7 8 - 3 8 1 2 ( 0 2 ) 0 0 0 3 2 - 8
264
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 1. Chemical potential µ/RT of water in CO2 -hydrate (dashed line) and in liquid solution (solid line). The standard state corresponds to pure liquid water at pressure and temperature of interest. Pressure is 3 MPa and mole fraction of dissolved CO2 is X = 0.020. The equilibrium temperature Teq (P , X) is defined by the intersection of the two curves. The driving force for nucleation is proportional to the difference µhl = µh − µl . Stable nuclei of hydrate appear when µhl is sufficient to overcome the energy surface of the nuclei.
potentials µh of water in the hydrate phase are described by the thermodynamic model of van der Waals and Platteeuw[3]. The solid line shows µl /RT, while the dashed line indicates µh /RT (both potentials are referred to the value for pure liquid water at the pressure and temperature of interest). As Fig. 1 shows, the energy difference becomes large enough to overcome the surface energy barrier and launch the process of nucleation at sufficiently low temperature. In this paper, we use experiments to study the timescale, over which mixtures of carbon dioxide and water evolve toward their equilibrium state. Quantifying the rate of hydrate nucleation is a difficult experimental problem. First, we must be able to detect the appearance of hydrate nuclei. Second, we need to observe a large number of nucleation events in order to obtain a meaningful statistical average for the nucleation rate. Finally, we need to control the thermodynamic conditions of the experiments, which include pressure, temperature, and concentration of a solute. An experimental setup that deals with each of these challenges is described in [4]. The initial goal of this experiment was to crystallize gas hydrate from aqueous solution in the absence of bulk gas [5]. To achieve this goal, hydrate formation occurs in a water-saturated porous medium. The medium effectively isolates the aqueous solution in it from the bulk gas phase, which is used to maintain pressure in the experiment (pressurized hydrate-forming gas is supplied to the top of the pressure cell). The time required to equilibrate the pore liquid with the vapour phase depends on thickness of the porous medium because transport of a solute through the porous medium occurs solely by diffusion. When the porous medium is more than a few centimeters thick, the concentration of the solute in the interior of the porous medium remains nearly constant over a 2 days experiment.
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
265
In this experimental setup, hydrate formation can not be visually detected. Instead, we use electrical resistance measurements to monitor nucleation and growth of hydrate. We choose CO2 as the hydrate-forming gas mainly because it supplies ions when dissolved in water to alter electrical conductivity of the solution. The liquid conductivity, determined by concentration and mobility of these ions, can be expressed in terms of the solute concentration X and temperature T [4]. The conductivity increases with X because the abundance of the solute determines the number of ions, while temperature primarily alters the ion mobility. The measurements of electrical resistance R are inversely proportional to the pore liquid conductivity σ (X, T ) and can be represented by R(X, T ) =
A σ (X, T )
(1)
where A is a factor, which depends on porosity and geometric factors of the electrode configuration [4]. As hydrate nuclei appear in the pore fluid, they consume the solute, which alters the conductivity of solution. We infer these changes in the solute concentration by measuring electrical resistance. Since the resistance measurements are made over a large volume of the porous medium, we average over many nucleation sites and obtain a more reliable estimate of the average nucleation rate. Data collected in the experiment provide quantitative information about nucleation. To process this information and recover nucleation characteristics, we develop a simple theoretical model of hydrate formation in porous media. Estimates of the nucleation time are extracted from our experimental data and the results are extrapolated to other thermodynamic conditions.
2. Experimental results Fig. 2 presents measurements of resistance and temperature as functions of time collected in an experiment. The measurements show that changes in resistance reflect changes in temperature very closely. However, resistance continues to increase with time, even when temperature reaches a nearly constant value of T ≈ −1.5 ◦ C. (Freezing of water occurs at −1.7 ◦ C due to the presence of the solute.) This transition in the behavior of resistance versus temperature becomes more evident when resistance is plotted as a function of temperature (see Fig. 3). Fig. 3 demonstrates the steady increase in measurements of resistance with cooling from 15 ◦ C to approximately 0 ◦ C. Prior to hydrate formation, the solute concentration is constant and the changes in resistance reflect the effects of temperature on ion mobility. We find that these measurements are well approximated by assuming that the conductivity σ (X, T ) of the pore liquid varies linearly with T . This allows us to recover the temperature dependence of ion mobility using measurements of resistance and temperature. When hydrate forms in the pore liquid, the solute is transferred from solution into growing hydrate crystals. Increases in resistance due to depletion of solution with the solute concentration can be identified in Fig. 3 when temperature decreases from 0 to −1.5 ◦ C. The resistance increases sharply and continues to increase even when temperature is held constant at −1.5 ◦ C. We do not cool the apparatus below −1.5 ◦ C to prevent freezing of the water. We subsequently increase temperature and hydrate starts dissociating in order to re-establish higher equilibrium values of the solute concentration [6]. Measurements from two subsequent cycles of warming and cooling show that changes in resistance with temperature follow the same path with almost no offset [4]. This suggests that equilibrium states in our system are established on a timescale which is less than that of the controlled changes
266
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 2. Experimental measurements of resistance (A) and temperature (B) as a function of time from a typical experiment at P = 2.8 MPa. Measurements collected during cooling are indicated by the line, whereas the points represent individual measurements during warming. A small change in the slope of R(T ) near T = 4 ◦ C coincides with the disappearance of hydrate from solution during warming.
in temperature. On the basis of these results, we assume that changes in resistance with temperature during warming mirror changes in the solubility as the system adjusts through a sequence of equilibrium states. Correspondingly, the offset between cooling and warming paths in Fig. 3 is representative of kinetics of the process. The data suggest that hydrate starts forming at T ≈ 0 ◦ C. On the other hand, measurements of resistance at temperature above 4 ◦ C nearly coincide with measurements of resistance at the same temperature during the initial cooling. This implies that hydrate totally disappears from solution at
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
267
Fig. 3. Resistance as a function of temperature for the experimental data shown in Fig. 2.
T = 4 ◦ C and that the equilibrium temperature for aqueous solution and hydrate phases at the initial solute concentration is Teq = 4 ◦ C. It follows that the solution was cooled almost 4 ◦ C below the equilibrium temperature to initiate nucleation of hydrate. We infer nucleation characteristics from the time dependence of resistance measurements during the interval when temperature is nearly constant. Fig. 4 reproduces the measurements in Fig. 2, but corrects the values of resistance for the influence of temperature on mobility of ions. This procedure is not necessary because we are going to analyze values of resistance at constant temperature anyway. However, the corrected values of resistance reflect the effect of changes in the solute concentration and more clearly represent the kinetic aspects of the process. Fig. 4 shows that the resistance does not change much within the first 5 h of the experiment because the solute concentration X is constant prior to hydrate formation. As hydrate crystallizes from solution and X decreases, values of the resistance increase. This increase in the resistance is well fit to an exponential function of time t/τ , where the timescale τ is 4 h. We have several other experiments conducted at slightly different thermodynamic conditions. Changes in the resistance with time show the same pattern as that shown in Fig. 4. The exponential fit yields different values of timescale τ for different thermodynamic conditions. To understand these results and to extrapolate them to other thermodynamic situations, we develop a theoretical model for hydrate formation in porous media. An essential element of this model is the relationship between the number of pores with hydrate nuclei and the average electrical resistance of the medium.
268
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 4. Experimental measurements of resistance from Fig. 2, corrected for the effects of temperature. Resistance is nearly constant prior to the start of hydrate formation. Subsequent increases in resistance are attributed to the nucleation and growth of hydrate in the porous medium. A fit of (7) to the resistance measurements yields a nucleation time of τ = 4 h.
3. Theoretical modeling We consider the appearance of hydrate nuclei in a large number of pores. Nucleation does not occur in all pores at once as temperature decreases. Instead, the nuclei appear randomly. Under uniform thermodynamic conditions in identical pores, there is an equal probability of forming a nucleus in each of them. We also assume that the pores are chemically isolated. This means that once a nucleus forms, it equilibrates with the pore liquid without affecting the neighbouring pores. Under this approximation, nucleation of hydrate in pores is analogous to the freezing of water drops [7]. However, gas hydrate involves two components, so the details of the kinetic process as well as the thermodynamics of the bulk phases are different. Despite these differences, the freezing of water drops still serves as a useful conceptual model to illustrate the stochastic nature of nucleation. The probability of forming a nucleus in a pore of volume V in a time increment dt is VJ dt, where J is rate of nucleation. If N(t) denotes the number of pores without nuclei at time t, then dN = −N(t)VJ dt
(2)
For the special case of constant J , the number Nh of pores with hydrate nuclei is given by Nh (t) = N0 (1 − e−t/τ )
(3)
where N0 is the initial number of empty pores, τ = (JV)−1 is the timescale of nucleation, determined by the nucleation rate J and the pore volume V . The nucleation rate J can be generally represented by an Arrhenius equation − G J = J0 exp (4) kT where J0 is the pre-exponential factor and G is the activation energy.
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
An equation of the same form − G∗ J = J0 exp kT
269
(5)
holds for classical nucleation. In the theory of classical nucleation, G∗ is the increase in the Gibbs free energy due to nucleation of a stable embryo, which can be represented as [8]
G∗ =
16πvh2 γhl3 3 µ2hl
(6)
where vh is the molar volume of hydrate, µhl = µh − µl is the difference in chemical potential of water in the hydrate and liquid phase at a given thermodynamic condition. Surface tension between hydrate and water γhl can be considered as a lumped parameter that embodies the influences of pore walls and other possible complications on nucleation [8,9]. In classical nucleation, the rate of nucleation depends on local thermodynamic conditions through the value of µhl in G∗ . Our calculations show that G∗ is strongly dependent on T . This suggests that the activation energy G in (4) may have a similar dependence. Because G is expected to vary with temperature, we can not treat it as a constant parameter in the numerical model. Instead, we adopt classical nucleation and treat J0 and γhl as constant parameters to recover from the experimental data. The activation energy at experiment conditions can then be inferred from (6). 4. Numerical simulation We simulate the electrical resistance of the porous medium using a network of resistors. The value of resistance for each resistor is determined by the conductivity of the pore liquid, while the network defines the pathways for electrical current. Consequently, the resistance of the pore liquid changes when a nucleus forms and this alters the flow of current through the network. Since the electrical current is inversely proportional to bulk resistance, changes in bulk resistance can, in principle, be related to changes in the number of pores with nuclei. We establish this relationship below. The porous medium is considered as a large number of pores distributed in a three-dimensional cubic network. The whole network is subdivided into elementary cubes. The relationship between the porous medium and the network is shown in Fig. 5. Resistors are associated with the pores, where the value of resistance is inversely proportional to the conductivity of the pore liquid. The pores are connected by pore throats which intersect at the vertices of the cube. When voltage is applied, the network bonds (elementary cube sides) serve as pathways for electrical current. The flow of current through the circuit is affected by the distribution of resistances in the network. As a hydrate nucleus appears in a pore, the conductivity of the pore liquid changes. (We assume that the solute concentration evolves to its equilibrium value at once as hydrate nuclei appear in the pore liquid.) This alters the resistor distribution and causes a change in bulk resistance. We compute the bulk resistance by defining a system of algebraic equations based on conservation of charge at the vertices of the cube [10]. The algebraic equations are solved for the current along individual pathways using a relaxation techinque. The bulk resistance is evaluated using the net current across the network and the applied voltage. Each pore has an identical value of resistance at the start of the simulation. Prior to hydrate nucleation, concentration X0 of the solute is constant. Cooling causes a uniform increase in resistance due to the
270
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 5. (A) Schematic illustration of pore and pore throats in the porous medium. (B) Correspondence between the porous medium and the network resistor model.
temperature dependence of the liquid conductivity. Once hydrate nuclei begin to form in the pores, their subsequent growth depletes the pore liquid to the equilibrium concentration Xeq . All other pores remain at their initial concentration X0 , under the assumption that the individual pores are chemically isolated. A further complication is introduced by the fact that Xeq depends on temperature. We account for this in the simulation by assuming that pores with nuclei maintain the equilibrium concentration Xeq (T ) as temperature changes. The probability of forming a nucleus in a volume V within a time increment t is given by JV t. J is the rate defined in (5), where J0 and γhl are free parameters in the model. We use a coarse-grained sand as the porous medium to avoid influence of the pore size on thermodynamic conditions of hydrate formation and adopt a nominal pore volume of V = 10−9 m−3 [11]. The stochastic nature of nucleation is simulated using a sequence of random numbers from a uniform distribution between (0, 1). The probability of nucleation in a pore within a given time interval t is equal
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
271
Fig. 6. Predictions of hydrate nucleation in a porous medium at constant temperature. The stochastic model is illustrated using three different values of temperature (0, 1, and 1.5 ◦ C). (A) Predicted number Nh of pores with nuclei over time in a three-dimensional network of 400 pores. (B) Predicted change in the bulk resistance R due to the formation of hydrate nuclei in the pores. Initial and final resistances in the three simulations are different because of the temperature dependence of R.
to the probability that the random number lies between 0 and JV t. When the random number falls in this range, we assume that nucleation has occurred. This procedure is repeated at each pore in the network for each successive increment of time. The resistances of pores with nuclei are adjusted to represent local equilibrium and the resulting change in the bulk resistance is computed using the network model. Calculations confirm that changes in the bulk resistance R do reflect changes in the number of pores with nuclei Nh . Fig. 6 shows several examples of R and Nh as a function of time from simulations at constant temperature. The curves are slightly rough because of the stochastic nature of the simulation, but there is very good agreement in the trends. In fact, a plot of R versus Nh reveals a nearly linear dependence. Thus, R is expected to follow the trend in Nh (t) very closely. Fig. 7 shows how closely the change in Nh (t) from the simulation at constant temperature follows the simple exponential form in (3). In this particular example, we set the nucleation time τ using a representative estimate of the nucleation rate J . We can also recover an estimate of τ by fitting (3) to the output of the simulation. The good agreement between the input and recovered τ suggests that the
272
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 7. A comparison of the number of pores with nuclei between the simulation at constant temperature and the analytical prediction in (3). The nucleation time used in the simulation is recovered using a fit of (3) to the simulation output.
nucleation time can be estimated from the experimental data by assuming that the predicted change in bulk resistance can be approximated by R = R(X0 , T ) + [R(Xeq , T ) − R(X0 , T )][1 − e−t/τ ]
(7)
where R(X0 , T ) is the initial resistance before any nuclei form and R(Xeq , T ) the resistance that would be measured if all of the pore were filled with nuclei and the concencentration of dissolved CO2 was everywhere equal to the equilibrium value Xeq . The timescale is given by τ = (VJ)−1 . We use this approach in analyzing the experimental data. 5. Recovering nucleation characteristics Fig. 4 shows that the temperature-corrected resistance measurements are nearly constant during the first 5 h of the experiment. After this time the resistance increases sharply due to the onset of hydrate nucleation and growth. After 7 h, the temperature reaches a nearly steady value, while the resistance continues to increase. A fit of the exponential form in (7) to the data after 7 h yields the smooth curve shown in Fig. 4. The nucleation time τ recovered from the fit is 4 h at T = −1.5 ◦ C. In order to interpret τ we need to evaluate the thermodynamic conditions that prevailed prior to hydrate nucleation. Temperature is measured, while pressure is maintained at P = 2.8 MPa and the mole fraction X 0 of the solute in the porous medium is constant. We inferred the equilibrium temperature Teq = 4 ◦ C for hydrate stability from our experimental data (see discussion of Fig. 3). The solubility Xeq (T ) within the hydrate stability zone is known from previous experiments [4], so we take Xeq (T = 4 ◦ C) = X0 ≈ 0.02. Because we estimate τ at T = −1.5 ◦ C, while T eq = 4 ◦ C, the cooling T below the equilibrium temperature is T = 5.5 ◦ C. Similar results are obtained from two other experiments at slightly different pressure and temperature conditions. The solute mole fraction is nearly constant at X0 = 0.02 in all three experiments. (Small differences in X0 are inferred from measurable differences in T eq in each experiment.) A summary of results is listed in Table 1. Differences in the thermodynamic conditions for each experiment have been
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
273
Table 1 Summary of experiments P (MPa)
T (◦ C)
X (mole fraction)
Teq (◦ C)
T (◦ C)
τ (h)
2.8 2.7 2.6
−1.5 0.8 1.2
0.0197 0.0202 0.0201
4.0 5.5 5.2
5.5 4.7 4.0
4.0 ± 0.5 4.8 ± 0.5 5.9 ± 0.5
characterized in terms of T . The recovered values of τ vary systematically with T , although the range of T is not large. We explore the range of values for J0 and γhl that yields predictions for τ which are consistent with the experiments. To do so, we fit experimental estimates of τ to the prediction τ = (JV)−1 , where J is given by (5). Thermodynamic calculations, based on the technique described in [6], for µhl in (6) are carried out for the experimental conditions listed in Table 1. Using the three estimates of τ from Table 1, we obtain three equations for J0 and γhl . The best fitting values are J0 = 9 × 104 m3 s−1 and γhl = 0.0015 J m−2 . Fig. 8 shows a good agreement between the theoretical prediction and the experimental estimates of τ .
Fig. 8. Experimental estimates of τ (circles) and theoretical prediction (line) calculated using J0 = 9 × 104 m3 s−1 and σ = 0.0015 J m−2 .
274
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
Fig. 9. Nucleation time τ for CO2 -hydrate as a function of temperature T . Pressure is fixed at P = 2.8 MPa and two values of the solute mole fraction are considered: X1 = 0.020 and X2 = 0.018. The equilibrium temperatures for the two concentrations are equal to Teq (P , X1 ) = 4 ◦ C and Teq (P , X2 ) = 1.5 ◦ C. These values are based on an experimentally-derived phase diagram [4].
Fig. 9 presents the temperature dependence of nucleation time τ predicted with the recovered values of J0 and γhl . To demonstrate the influence of the solute concentration on nucleation time, estimates of τ are given at X = 0.020 and 0.018 mole fraction. The plots in Fig. 9 indicate that τ increases with temperature and becomes infinite when T approaches Teq (X, P ). The plots also show that at lower temperatures, the nucleation time becomes relatively insensitive to T and more dependent on X. The larger nucleation time corresponds to the lower solute concentration. The dependence of τ on X arises in several ways. First, the value of X is important in the calculation of µhl . Second, we expect X to influence J0 by altering the number of elementary hydrate building units per unit volume of solution [12–15]. For simplicity, we assume that J0 is linearly dependent on X and adopt J0 = 8.1 × 104 m3 s−1 at X = 0.018. This approach is used to extrapolate our estimates of τ to other solute concentrations. 6. Conclusions We design an experiment to crystallize gas hydrate from dissolved CO2 in natural porous media. Experimental data collected in a series of experiments provide quantitative information on nucleation of CO2 -hydrate. We use a stochastic model for nucleation of hydrate in porous media to analyze these data and estimate the nucleation timescale under different thermodynamic conditions. The nucleation
O.Ye. Zatsepina, B.A. Buffett / Fluid Phase Equilibria 200 (2002) 263–275
275
time is extrapolated to other values of temperature and concentration of the solute using this simple model. Fitting the model to the experimental data yields estimates of the surface tension and the kinetic prefactor. The surface tension γhl , which we recover from the experimental data, is about 18 times smaller than the surface tension between ice and water (γhl = 0.0276 J m−2 ). The latter is often taken as a reasonable approximation for the surface tension between hydrate and water. The lower value recovered in our experiments implies that the effective surface energy of the nuclei is smaller than that predicted for homogeneous nucleation. This is consistent with our expectations that the pore walls and other foreign particles serve as nucleation sites and promote nucleation. In fact, calculations using γhl = 0.0267 yields the nucleation time that exceeds the age of the Earth when T = 5 ◦ C. Clearly, the porous medium plays an important role in hydrate nucleation. Acknowledgements We thank two anonymous reviewers for helpful comments and suggestions. This work is supported by the Institute of Applied Energy and by NSERC. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15]
E.D. Sloan, Clathrate Hydrates of Natural Gases, 2nd Edition, Marcel Dekker, New York, 1998. M.A. Trebble, P.R. Bishnoi, Fluid Phase Equilib. 40 (1988) 1–21. J.H. van der Waals, J.C. Platteeuw, Adv. Chem. Phys. 2 (1959) 1–57. O.Ye. Zatsepina, B.A. Buffett, Fluid Phase Equilib. 192 (2001) 85–102. B.A. Buffett, O.Ye. Zatsepina, Marine Geology 164 (2000) 66–77. O.Ye. Zatsepina, B.A. Buffett, Geophys. Res. Lett. 24 (1997) 1567–1570. T.W. Barlow, A.D.J. Haymet, Natural gas hydrates, Ann. New York Acad. Sci. 715 (1994) 549–551. P.V. Hobbs, Ice Physics, Clarendon Press, Oxford, 1974. K.F. Kelton, Solid State Phys. 45 (1991) 75–177. R.J. Suman, R.J. Knight, Geophysics 62 (1997) 1151–1162. V. Melnikov, V. Nesterov, in: Proceedings of the 2nd International Conference on Natural Gas Hydrates, Toulouse, France, 1996, pp. 541–549. D. Turnbull, J.C. Fisher, J. Chem. Phys. 17 (1949) 71–73. F.F. Abraham, Homogeneous Nucleation Theory, Academic Press, New York, 1974. P.M. Rodger, J. Phys. Chem. 94 (1990) 6080–6089. R.L. Christensen, E.D. Sloan, Natural gas hydrates, Ann. New York Acad. Sci. 715 (1994) 283–305.