Stages in the Dynamics of Hydrate Formation and Consequences for Design of Experiments for Hydrate Formation in Sediments

Abstract

Natural gas hydrates in sediments can never reach thermodynamic equilibrium. Every section of any hydrate-filled reservoir is unique and resides in a stationary balance that depends on many factors. Fluxes of hydrocarbons from below support formation of new hydrate, and inflow of water through fracture systems leads to hydrate dissociation. Mineral/fluid/hydrate interaction and geochemistry are some of the many other factors that determine local hydrate saturation in the pores. Even when using real sediments from coring it is impossible to reproduce in the laboratory a natural gas hydrate reservoir which has developed over geological time-scales. In this work we discuss the various stages of hydrate formation, with a focus on dynamic rate limiting processes which can lead to trapped pockets of gas and trapped liquid water inside hydrate. Heterogeneous hydrate nucleation on the interface between liquid water and the phase containing the hydrate former rapidly leads to mass transport limiting films of hydrate. These hydrate films can delay the onset of massive, and visible, hydrate growth by several hours. Heat transport in systems of liquid water and hydrate is orders of magnitude faster than mass transport. We demonstrate that a simple mass transport model is able to predict induction times for selective available experimental data for CO₂ hydrate formation and CH₄ hydrate formation. Another route to hydrate nucleation is towards mineral surfaces. CH₄ cannot adsorb directly but can get trapped in water structures as a secondary adsorption. H₂S has a significant dipole moment and can adsorb directly on mineral surfaces. The quadropole-moment in CO₂ also plays a significant role in adsorption on minerals. Hydrate that nucleates toward minerals cannot stick to the mineral surfaces so the role of these nucleation sites is to produce hydrate cores for further growth elsewhere in the system. Various ways to overcome these obstacles and create realistic hydrate saturation in laboratory sediment are also discussed.

1. Introduction

Natural gas hydrates are crystalline compounds which are mainly stabilized by hydrogen bonds that forms cavities which enclathrate small hydrocarbons. Inorganic gases, like for instance CO₂ and H₂S, also form hydrates. The extra stabilization of the water structure by the inclusion of these “guest” molecules makes hydrates stable also for temperatures above zero Celsius if the pressure is high enough. The formation of these ice-like compounds during processing and transport of hydrocarbons has motivated substantial amounts of hydrate research in the past. The enormous quantity of natural gas hydrates spread around the world stimulates creativity in development of new methods for commercial dissociation and use of these energy sources. Adding steam or hot water for thermal dissociation of the hydrate is very expensive. Substantial portions of the added heat are lost to minerals. Adding chemicals is also very expensive. The pressure reduction method is the approach that has drawn the most attention in the past. In this method the local pressure is reduced to below the stability limit pressure of the hydrate. One drawback of the pressure reduction method is that the heat necessary for dissociation of the hydrate still needs to be supplied. Studies conducted on pilot plant scale [1,2,3,4] have so far not been successful, and freezing down (from the heat drawn away by dissociation demand) has been only one of several problems. Production of sand and water are other drawbacks. In addition to the stability limit considerations, hydrate dissociation in sediment is limited by the transport processes across a thin interface (1.2 nm) between the hydrate and surrounding phases. These phase transition dynamics are implicitly coupled to the dynamics of the flow and all phase transitions in every pore. Finally, the coupled dynamics of all pores are connected to the reservoir flow. Pressure reduction attacks the top level of dynamics. The transfer of the impact all the way down to the phase transition dynamics will impose significant delays. Hence, pressure reduction is a slow method from a theoretical point of view. An experimental challenge for pressure reduction is to establish realistic boundary conditions that can supply heat from surrounding sediments. Many published experiments even report pressure reduction experiments at constant temperature as controlled by surrounding heat baths. Temperature changes, on the other hand, distribute very rapidly through liquid water and hydrate, making thermal approaches dynamically efficient.

The strong focus on reduction of CO₂ emissions during the latest two decades has also triggered rapid development of various ways to utilize CO₂, like for instance the use of CO₂ for enhanced oil recovery. The use of CO₂ for production of natural gas hydrate is another option. The exchange between CH₄ hydrate and CO₂ hydrate is possible through two primary mechanisms. The first mechanism is a solid-state exchange which has been verified experimentally for temperatures substantially inside the ice region for water [5,6] and also theoretically through molecular dynamics (MD) simulations [7]. For temperatures above freezing the water surrounding the hydrate in the pores is a combination of liquid water and pore bounded water. A typical hydrate saturation for Alaska is 75% or lower, like in the area used for the Ignik Sikimu pilot test [8,9]. Hydrate saturation is the percentage of available pore volume that is filled with hydrate. For temperatures in the liquid water region another mechanism [10,11] than solid state conversion is possible. Injected CO₂ will form new hydrate with free pore water. Released heat from this hydrate formation, and other factors, will dissociate the in situ hydrate. Similar to the formation of CH₄ hydrate the formation of a new CO₂ hydrate will also go through the various stages of the hydrate formation process. From a physical point of view any fluid/solid phase transition in the classical regime consists of two stages. The nucleation of hydrate requires that the free energy change of the formation be negative and the free energy difference has to be large enough to supply the work needed to push away the molecules in the original phases. The necessary hydrate particle volume needed to win the competition between benefits of the phase transition, and the penalty of the push work, is the critical size hydrate core. After this first nucleation stage the hydrate will grow steadily. but will compete with other hydrate cores for mass, and may diminish in competition with more stable cores. And also in the steady growth stage the heat released from the growth of neighboring hydrate cores can lead to dissociation. Also, as discussed further later, dilution of fluid phases will lead to dissociation. Even the two physically well-defined stages of hydrate growth are fairly complex in non-equilibrium systems. Rapid formation of hydrate films on the interface between hydrate former phase and water leads to a dramatic slow-down in supply of hydrate formers from typical liquid phase diffusivity (coefficients in the order of 10⁻⁷ m²/s to 10⁻⁸ m²/s) down to diffusivity for transport through hydrate (coefficient from 10⁻¹⁵ m²/s to 10⁻¹⁷ m²/s). Dramatic reduction in mass transport rates and limited shear forces (limited or no flow) lead to very long times, relative to the typical nanoseconds time scale for nucleation, before visible hydrate can be observed at millimeter scales. In addition to the nucleation stage and the growth stage the term induction frequently occurs in the open literature. This time can be defined as the time to reach visible hydrate or to reach the onset of massive hydrate growth rate (depends on monitoring method resolution).

We utilize classical nucleation theory (CNT) and residual thermodynamics to examine the dynamics of hydrate formation. This is used as the basis for considering various ways to break hydrate films and increase formation rate. The primary goal is to shed light on how to design laboratory experiments for hydrate formation between water and CH₄ or CO₂ in sediments. Reducing trapped gas pockets and trapped liquid water to a minimum will be an important step towards creating hydrates in sediments that in important ways can be compared to natural gas hydrates in nature.

2. Methodologies

Kvamme and Tanaka [12] used MD simulations to calculate the chemical potentials of water in ice and empty clathrates of structures I and II hydrate. They then used experimental data for heat of freezing at zero Celsius and experimental data for specific heat capacity of liquid water in order to get a continuous description of water chemical potential. Ideal gas chemical potential for the three phases (ice, empty clathrate of structure I and empty clathrate of structure II) were trivially sampled from momentum space during the MD simulations. The values calculated by Kvamme and Tanaka [12] are plotted in Figure 1 along with parameters listed in

Table 1. Coefficients a_i for chemical potential of water in various phases; ${\mu }_{H_{2}O}={\displaystyle \sum _{i=1}^{i=6}{a}_i}{\left(\frac{T_{C_{H_{2}O}}}{T}\right)}^{(i-1)}$ with $T_{C_{H_{2}O}}=647.096 \mathrm{K}$.

i	Str. I ${\mathit{a}}_{\mathit{i}}(\mathit{k}\mathit{J}/\mathit{m}\mathit{o}\mathit{l}\mathit{e})$	Str. II ${\mathit{a}}_{\mathit{i}}(\mathit{k}\mathit{J}/\mathit{m}\mathit{o}\mathit{l}\mathit{e})$	Liquid ${\mathit{a}}_{\mathit{i}}(\mathit{k}\mathit{J}/\mathit{m}\mathit{o}\mathit{l}\mathit{e})$	Ice ${\mathit{a}}_{\mathit{i}}(\mathit{k}\mathit{J}/\mathit{m}\mathit{o}\mathit{l}\mathit{e})$
1	−90.48377	−92.72076	−109.03432	−87.13485
2	42.46997	45.09842	58.33379	41.31529
3	−14.61482	−15.72658	−19.48522	−17.71143
4	1.76419	1.92639	2.32631	3.88149
5	0	0	0	−0.43451
6	0	0	0	0.01971

. Gas (or fluid phases) of CH₄ and CO₂ are modelled using residual thermodynamics in the manner that is similar for gas/liquid equilibrium in hydrocarbon systems. For this purpose we utilize the Soave-Redlich-Kwong (SRK) equation of state [13]. The only difference is that we use chemical potential rather than fugacity. For that reason, the gas phase needs model molecules for the ideal gas term for CH₄ and CO₂. CH₄ is modelled as a spherical molecule and the ideal gas chemical potential is an analytical integral over the Boltzmann factors of translational momentums in the x, y and z directions. The equations for this can be found in any textbook on physical chemistry or statistical mechanics and are not repeated here. The model we use for CO₂ is a linear molecule with only two independent rotational degrees of freedom with equal moments of inertia due to the symmetry. The translational part is given by the molecular weight and the center of mass. For completeness we also give the three rotational momentums for the water model utilized in our earlier study so that it will be possible to reproduce the results presented in this work in

Table 2. Moments of inertia for H₂O and CO₂ models utilized.

Axis	H2O Moment of Inertia (kg·m2)	CO2 Moment of Inertia (kg·m2)
x	2.938·10−47	0.7173·10−45
y	1.020·10−47	0.7173·10−45
z	1.918·10−47	0

below.

More advanced theories, like for instance phase field theories [14,15,16,17,18], may be more rigorous in the implicit couplings between the various contributions to the phase transition dynamics but are also more complex in terms of distinguishing between the various stages of hydrate growth. Classical nucleation theory (CNT) is still rigorous enough to illustrate nucleation times. In principle CNT can be used to determine two parameters related to volume and shape of a hydrate core [19] in transition to growth. For the main goal of this work, however, we stick to a spherical model. A critical size hydrate core is defined as the size when the hydrate core enters the stable growth period.

We also limit ourselves to some of the important hydrate phase transitions involved in hydrate formation in porous media. Formation of hydrate on the interface between a separate phase containing hydrate formers, and hydrate formation towards mineral surfaces, are the two fastest routes. When these two routes have created a hydrate film there will also be some hydrate formation from dissolved hydrate former in water. For methane this is limited by fairly low solubility, while the higher solubility of CO₂ in water makes this route far more significant.

In order for a hydrate particle to reach critical size, the hydrate former molecules must be transported from an outer boundary on the liquid water side and through an interface layer of gradually increased water structure. The more structured the water, the slower the transport. The necessary time for enough hydrate former molecules to reach surface of the hydrate core at critical size is the nucleation time. Further growth through the hydrate film is very slow as indicated above. Since heat transport is very fast, and maybe two to three orders of magnitude faster than mass transport in liquid water, it will be several orders of magnitude faster than mass transport through hydrate. Assuming a constant mass transport rate through the hydrate film in order to supply the hydrate film with hydrate formers we are then able to estimate hydrate film thickness as a function of time. This will provide a basis for evaluating the necessary dynamics needed to break hydrate films, and correspondingly ensure that measures can be taken to prevent blocking hydrate films.

2.1. Hydrate Formation from Water and a Separate Phase Containing Hydrate Formers

Hydrate phase transition along the equilibrium curve is reversible. At equilibrium the free energy changes of the hydrate formation from liquid water (or ice) and hydrate formers coming from gas, liquid or fluid state is zero. Outside of equilibrium the molar free energy change, as given by Equation (1) below, has to be negative.

$$\begin{array}{l}\Delta G^{(H_1)}=x_{H_{2}O}^{H_1}\left({\mu }_{H_{2}O}^{H_1}(T,P,{\overrightarrow{x}}^{H_1})-{\mu }_{H_{2}O}^{water}(T,P,\overrightarrow{x})\right)\\ +{\displaystyle \sum _i{x}_i^{H_1}\left({\mu }_i^{H_1}(T,P,{\overrightarrow{x}}^{H_1})-{\mu }_i^{gas}(T,P,{\overrightarrow{y}}^{gas})\right)}\end{array}$$

(1)

The superscript H₁ denotes this specific heterogeneous phase transition. T is temperature, P is pressure. x is mole-fraction in either liquid or hydrate (denoted with a superscript H). y is mole-fraction in gas (or liquid, or fluid) hydrate former phase; i is an index for hydrate formers. Superscript water denotes a water phase. Generally this is ice, liquid or adsorbed water on mineral surfaces. In this work we only consider liquid water; µ is chemical potential. x is mole-fraction in liquid water or hydrate (as given by superscripts). Vector sign denote mole-fractions of all components in the actual phase.

Symmetric excess formulation for liquid water chemical potential is given by:

$${\mu }_{H_{2}O}^{water}(T,P,\overrightarrow{x})={\mu }_{H_{2}O}^{pure,H_{2}O}(T,P)+RT\ln \left[x_{H_{2}O}{\gamma }_{H_{2}O}(T,P,\overrightarrow{x})\right]\approx {\mu }_{H_{2}O}^{pure,H_{2}O}(T,P)+RT\ln \left[x_{H_{2}O}\right]$$

(2)

$\lim ({\gamma }_{H_{2}O})=1.0$ when x_H2O approaches unity.

Water as superscript on the left-hand side distinguishes liquid water phase from water in the hydrate phase. A right-hand side approximation of Equation (2) is not necessary but good enough for the purpose of this work. The alternative would be to use a model for the activity coefficient or utilization of the Gibbs-Duhem relation. In our phase field theory (PFT) modelling of CO₂ hydrate phase transition dynamics studies [14,15] we used the latter approach. Our PFT models are fairly complex and a simpler kinetic model might be more useful in order to visualize the various stages of the hydrate formation.

Water chemical potential in the hydrate structure is given by [16]:

$${\mu }_{H_{2}O}^H={\mu }_{H_{2}O}^{O,H}-{\displaystyle \sum }_{k=1,2}RT{v}_k\ln \left(1+{\displaystyle \sum }_i{h}_{ij}\right)$$

(3)

H denotes hydrate and superscript O on first term on right hand side denotes empty clathrate. Calculated values for water chemical potentials in empty hydrates of structure I and II are readily available from model water (TIP4P) simulations [12] as discussed above. Cavities per water in structure I hydrate, ν_k is 1/23 for small cavities and 3/23 for large cavities. $h_{ki}$ is the canonical partitition function for s guest of type i in cavity type k. For a rigid water lattice the result is a Boltzmann integral over all possible water-guest and guest-guest interactions and a function of the free energy of the huest molecule [12]. This is the most common way to calculate $h_{ki}$ in various available codes for hydrate equilibrium. A different formulation of $h_{ki}$ utilize a perturbation approach in which the movements of the guest molecule, relative to energy minimum position in the cavity, is approximated by an harmonic oscillator. The advantage of this approach is that some frequencies of guest movements may interfere with the water lattice librational frequencies. As such we directly also get calculations for these effects, which are typically included as empirical correction, incorporated. For CO₂ a comparison with a rigid lattice calculations and the harmonic oscillator approach reveal a destabilization effect of 1 kJ/mole due to CO₂ movements in the large cavity of structure I [12]:

$$h_{ki}=e^{\beta \left[{\mu }_{ki}-\Delta g_{ki}\right]}$$

(4)

β is the inverse of the universal gas constant times temperature. μ_ki is the chemical potential of guest molecule i in hydrate cavity of type k. At equilibrium this chemical potential is equal to the chemical potential for the same molecule in the phase it comes from during the hydrate formation. For Equation (1) that means gas, liquid or fluid as a separate phase. In a non-equilibrium situation the guest chemical potentials are adjusted for distances from equilibrium through a Taylor expansion, as discussed later. The free energies of inclusion (latter term in the exponent) are given in

Table 3. Coefficients a_i for free energy of guest type j inclusion in cavity type k; $\Delta g_{kj}={\displaystyle \sum _{i=1}^{i=6}{a}_i}{\left(\frac{T_{C_j}}{T}\right)}^{(i-1)}$ with $T_{C_{C{H}_4}}=190.56 \mathrm{K}$ and $T_{C_{C{O}_2}}=304.13 \mathrm{K}$.

i	CH4		CO2 Large		CO2 Small
	Large	Small	T ≤ 283.14 K	T > 283.14 K
1	17.97150	−42.47683	41.52168	−17.87093	0.19929
2	−23.44013	119.24124	−41.96874	−17.89249	−28.28735
3	−161.81535	−183.19565	−70.72691	17.38136	−11.94528
4	45.20561	128.39252	−11.81084	−29.68940	−2.66250
5	36.67261	−54.98784	16.73045	−19.90321	3.85653
6	138.00217	−78.55671	21.91621	25.22112	3.21040

below.

Filling fractions in the various cavities, and mole-fractions in the hydrate are given by:

$${\theta }_{ki}=\frac{h_{ki}}{1+{\displaystyle \sum _j{h}_{ki}}}$$

(5)

${\theta }_{ki}$ is the filling fraction of component i in cavity type k:

$$x_i^H=\frac{{\theta }_{large,i}{\nu }_{large}+{\theta }_{small,i}{\nu }_{small}}{1+{\theta }_{large,i}{\nu }_{large}+{\theta }_{small,i}{\nu }_{small}}$$

(6)

ν is fraction of cavity per water. Corresponding mole-fraction water is then given by:

$$x_{H_{2}O}^H=1-{\displaystyle \sum _i{x}_i^H}$$

(7)

The associated hydrate free energy is then:

$$G^{(H)}=x_{H_{2}O}^H{\mu }_{H_{2}O}^H+{\displaystyle \sum _i{x}_i^H{\mu }_i^H}$$

(8)

where µ is chemical potential. H₂O subscripts denote water. i are hydrate formers. H superscripts denote hydrate. x is mole-fraction in hydrate (superscript H) and G is free energy.

Guest molecule i (in the case of this work either CO₂ or CH₄) chemical potential that enters Equations (4) and (8) at equilibrium is given by:

$${\mu }_i(T,P,\overrightarrow{y})={\mu }_i^{pure,idealgas}(T,P,\overrightarrow{y})+RT\ln \left[y_i{\varphi }_i(T,P,\overrightarrow{y})\right]$$

(9)

where y_i is the mole-fraction of component i in the gas (or liquid or fluid) mixture. ${\varphi }_i$ is the fugacity coefficient for component i. Chemical potential for a monoatomic model of methane in ideal gas state is trivial and analytical from statistical mechanics and the Boltzmann integral over translational momentums. For the rigid CO₂ model there are two additional rotational contributions. The necessary moments of inertia are given in Table 2. The resulting ideal gas chemical potential ends up as trivial functions of temperature and density. The SRK [13] equation of state is utilized for calculating the fugacity coefficient, and also the density needed for the ideal gas chemical potential calculations.

2.2. Temperature and Pressure as Driving Forces for Heterogeneous Hydrate Formation on Water/Hydrate Former Interface

With two components (one hydrate former and water) and three phases (water, hydrate former phase and hydrate), there are 12 independent thermodynamic variables, while the sum of conservation laws and conditions for thermodynamic equilibrium cover 11 variables. Hence, when both temperature and pressure are given locally in a reservoir, the system is over determined by 1, and there is no unique equilibrium state. But a pressure temperature equilibrium curve will still represent asymptotic limits at which hydrate formation ends (or enter a modus on infinite time to continue) at that two-dimensional projection of the thermodynamic variables.

For the heterogeneous case we therefore first calculate the equilibrium curve. For pure CH₄ or pure CO₂, and a defined temperature, the chemical potential for the guest is given by Equation (9) for a variable pressure. This chemical potential enters Equation (4) along with the free energy of inclusion, which is a function of temperature only for small and large cavities. Equations (2) and (3) are solved for the same chemical potential of water in liquid state and hydrate given the equilibrium pressure. For completeness the free energies of inclusion are given in Table 3 below. Calculated equilibrium curves for CH₄ and CO₂ are compared to experimental data in Figure 2. It is important to keep in mind that the parameters for free energies of inclusion have been calculated from MD simulations with the procedures described in Kvamme and Tanaka [12] for temperatures up to 280 K for CH₄. As such, all calculated results are predictions, and results for temperatures above 280 K are all extrapolations. For CO₂, on the other hand, a system of nine unit cells of hydrate structure I was surrounded on all sides with a 3 nm thick layer of CO₂ molecules. The system was run at 273.16 K and 283 K using SRK calculated densities in order to keep consistency with other calculations in this work that utilize the SRK equation of state. The system was then also run for 284 K and 290 K, which are both in the dense CO₂ region. The essential difference between the low temperature results and the high temperature results are the impact of the surrounding CO₂ molecules on the water lattice vibrations. The hydrate water structure is dominated by the 75% large cavities and the small cavity does not even benefit much from any stabilization of CO₂ due to the size and shape of the molecule. For the small cavities in structure I there is not statistical justification to distinguish between the effect of surround gas or liquid CO₂.

The models for extending thermodynamic properties of water and fluid phases are discussed in detail elsewhere [19,20,21] and will not be repeated here due to space limitations.

Two aspects regarding Figure 2a,b are worth mentioning. The equilibrium curve for CO₂ in Figure 2b does not have a discontinuity but there is a rapid change in density related to a phase transition. This leads to a rapid decrease in fugacity coefficients and a corresponding higher pressure needed to reach hydrate equilibrium. Some researchers actually presents Figure 2b as a continuous smooth curve without any rapid change. A second aspect is illustrated in Figure 2a. As can be seen, the equilibrium pressures for CO₂ hydrate are higher than the equilibrium pressures for CH₄ after the CO₂ density increase. This is often wrongly interpreted as a higher stability for the CH₄ hydrate above the temperature for the CO₂ phase transition. Pressure and temperature are independent thermodynamic variables. Free energies are the corresponding thermodynamic responses for the level of thermodynamic stability. In Figure 2b we plot free energies of CH₄ and CO₂ hydrate as function of temperature along the equilibrium pressures for the two types of hydrates.

Along the curves in Figure 2 the free energy changes in Equation (1) are trivially zero. Outside equilibrium all properties of fluids are continuous and can be calculated for any pressure and temperature. The necessary pressure correction for water chemical potential in Equation (2) is available from the molar volume of liquid water, which is almost constant and independent of temperature and pressure, time pressure minus equilibrium pressure for the actual temperature. Hydrate water chemical potential, on the other hand, is based on an equilibrium theory. The derivation [12] from the grand canonical ensemble in statistical mechanics leads to a Langmuir type of adsorption theory as expressed by Equation (3). For:

$$G_{\mathrm{Non}-\mathrm{equilibrium}}^H(T,P,\overrightarrow{x})=\begin{array}{c}{G}^{H,\mathrm{Eq}.}(T^{\mathrm{Eq}.},P^{\mathrm{Eq}.},{\overrightarrow{x}}^{\mathrm{Eq}.})+{\displaystyle \sum }_r{\frac{\partial G^H}{\partial x_r}|}_{P,T,x_{i\ne r}}\left(x_r-x_r^{\mathrm{Eq}.}\right)\\ +{\frac{\partial G^H}{\partial P}|}_{T,\overrightarrow{x}}\left(P-P^{\mathrm{Eq}.}\right)+{\frac{\partial G^H}{\partial T}|}_{P,\overrightarrow{x}}\left(T-T^{\mathrm{Eq}.}\right)\end{array}$$

(10)

we choose T to be equal to an equilibrium temperature, for which equilibrium pressure and compositions are calculated according to the equation and discussion above. The correction for pressure change is related to the partial molar volumes of water and hydrate formers in hydrate. Water partial molar volume is given by the structure of the hydrate while the occupation volume of the guest molecules can be calculated from Monte Carlo simulations according to the procedures described by Kvamme and Lund [26] and Kvamme and Førridahl [27]. The calculated volumes are 164.2 ų/molecule and 89.2 ų/molecule for CH₄ in large and small cavities, respectively. Corresponding values for CO₂ are 135.6 ų/molecule and 76.9 ų/molecule for large and small cavities, respectively. The value for CO₂ in small cavities is not practically interesting since the filling fractions of CO₂ in small cavities is practically zero.

Classical nucleation theory (CNT) for a spherical hydrate particle can be written as:

$$J=J_0{e}^{-\beta \Delta G^{Total}}$$

(11)

where J₀ is the mass transport flux supplying the hydrate growth. The phase transition in Equation (1) it will be the supply of CH₄ or CO₂ across an interface of gradually more structured water towards the hydrate core, as discussed in Kvamme et al. [19,20,21]. The units of J₀ will be moles/m²s for heterogeneous hydrate formation on the growing surface area of the hydrate crystal. β is the inverse of the gas constant times temperature and ΔG^Total is the molar free energy change of the phase transition. This molar free energy consists of two contributions (Equation (1) with Equation (10)) to correct for hydrate properties outside of equilibrium which are the free energy benefit of the phase transition and the second contribution is the penalty related to the work of pushing aside old phases. Even for a hydrate forming on the gas/water interface any hydrate core below critical size will be covered by water also on the side facing the gas due to capillary forces. Molar densities of liquid water and hydrate are reasonably close. It is therefore a fair approximation to multiply the molar free energy of the phase transition with molar density of hydrate times the volume of hydrate core. The penalty of the push work is the interface free energy times the surface area of the hydrate crystal. The total free energy change in extensive formulation (underlines indicate Joule units):

$$\Delta {\underset{\_}{G}}^{Total}=\Delta {\underset{\_}{G}}^{Phasetransition}+\Delta {\underset{\_}{G}}^{Pushwork}$$

(12)

The simplest possible geometry of a crystal is a sphere and for a core radius R the result is:

$$\Delta {\underset{\_}{G}}^{Total}=\frac{4}{3}\pi R^3{\rho }_N^H\Delta G^{Phasetransition}+4\pi R^2\gamma$$

(13)

where ${\rho }_N^H$ is the molar density of the hydrate and $\gamma$ is the interface free energy between the hydrate and the surrounding phase. There is no reliable value for interface free energy between hydrate and liquid water since the measurement of such a property is extremely difficult. We have not found any value in the open literature. Interface free energy between liquid water and ice is available [28] and the reported value is 29.1 mJ/m². We have used this value as an approximation for interface free energy between liquid water and hydrate in the last term of Equation (13).

Differentiation of (13) with respect to R gives the solution for maximum free energy radius (the critical core size):

$$R^{*}=-\frac{2\gamma }{{\rho }_N^H\Delta G^{Phasetransition}}$$

(14)

Superscript * denotes critical nuclei radius. Calculated critical radii for CH₄ hydrate at two different temperatures are given in Figure 3a. The associated nucleation times are based on integration of Fick’s law:

$$t(R)-t(o)={\displaystyle {\int }_{C_{C{H}_4}(R=0)}^{C_{C{H}_4}(R=12)}\frac{\partial C_{C{H}_4}(z)}{\left[-D_{C{H}_4}(z)\frac{{\partial }^2{C}_{C{H}_4}(z)}{\partial z^2}\right]}}$$

(15)

The concentration profiles for CH₄ as a function of distance from the liquid side of the interface (z = 0 Å) to the hydrate side of the interface (z = 12 Å) is available from Kvamme [19,20] and Kvamme et al. [21]. The diffusivity coefficient profile for CH₄ is assumed to be the same as for CO₂ since the diffusion of CH₄ and CO₂ through the liquid/water interface is dominated by the gradually increased structure of water towards the hydrate side of the interface. The profile is given by Figure 8a in Kvamme [19] and fitted to Equation (16) with the parameters in

Table 4. Parameters for Equation (16).

i	Parameter	i	Parameter	I	Parameter
1	0.979242	4	171.673	7	−9649.96
2	15.5427	5	6.76975	8	14,779.7
3	−88.5112	6	1939.55	9	−7496.15

below.

$$\left[\frac{D{(z)}_j}{D_{liquid,j}}\right]=9.5·{10}^{-{\displaystyle \sum _{i=1}^9{a}_i{\left(a\tan \left[(0.5z\frac{\pi }{2})/12\right]\right)}^{i-1}}}$$

(16)

Using D_liquid,j equal to 1.0·10⁻⁰⁸ m²/s and integrating Equation (15) to critical size. That is, for every supply of hydrate-former needed to grow a hydrate core, the transport has to cross the interface thickness at the mass transport penalty given by Equation (16). The number of transported molecules is then recalculated to provide a corresponding radius added to the core size. This latter recalculation involves the calculated filling fraction and the corresponding volume of hydrate water. With reference to the equilibrium curve in Figure 2a note that the equilibrium pressure for 283.16 K is around 80 bars and that is the reason for the exponential increase in critical size when the pressure approaches 80 bar. As expected, there is a substantial increase in nucleation time between these two temperatures. At 273.16 K and 150 bar the calculated critical radius is 12.0 Å and nucleation time is essentially instantaneous (less than one ns). Calculated critical radius at 283.16 K and 150 bars is 33.6 Å, and nucleation time is 26.5 ns.

For transport of CH₄ through the hydrate, in order to continue the growth of the film, it is now assumed that the diffusivity coefficient, and associated diffusion rate at the hydrate side of the interface, is equal to a stationary transport through the hydrate. This rate is expected to be slightly higher than regular diffusion through a block of hydrate without dynamics of couplings to heat transport and dynamic situations of partial local dissociation and reformation. See Kvamme et al. [21] for a brief summary of reported diffusivities of hydrate-formers through hydrate. In summary the diffusivities range from 10⁻¹⁵ m²/s to 10⁻¹⁷ m²/s. Quite a number of the studies are based on Monte Carlo studies of “cavity jumping” and rarely reflect any mechanism for how guest molecules are actually able to move around in the hydrate structure. Based on our MD studies [12] the mechanism seems more like a temporary destabilization of water hydrogen bonding structures between filled cavities and empty neighbor cavities. Nevertheless, more detailed studies are needed to verify these observations in a scientifically satisfactory manner. For our purpose in the context of this work we then might expect diffusivity coefficients which are higher than the range indicated above due to the couplings with heat transport dynamics rapid local dissociation/refreezing effects. The heat transport dynamics is implicitly coupled through the relationship:

$$\frac{\partial {\left[\frac{\Delta G^{Total}}{RT}\right]}_{P,\overrightarrow{N}}}{\partial T}=-\left[\frac{\Delta H^{Total}}{R{T}^2}\right]$$

(17)

which then connects Equations (12), (13) and (1) for any hydrate core since the number of moles in the actual radius of the core is simply the volume of the core times the molar density of the hydrate core. In the simple calculations below we disregard effects of heat transport dynamics. Based on our earlier calculations [15,16,17,18], and also in accordance with other sources in the literature, the heat transport rate through liquid water is 2–3 times faster than the corresponding mass transport rate. This ratio will be orders of magnitude higher for transport through hydrate. In accordance with this we expect that the mass transport rate is only affected by heat transport due to fast local phase transitions. Within the focus of this work we simply examine a couple of values for D_liq in (16) and plot the time needed to reach 1 mm hydrate film thickness. In Figure 4 we plot results for D_liq = 5·10⁻⁰⁸ m²/s in Equation (16) and the condition of temperature and pressure equal to that used in a magnetic resonance imaging experiment conducted at the ConocoPhillips research laboratory in Bartlesville (OK, USA). The experiment was conducted at a temperature of 276.25 K and 83 bars with initial equal volumes of liquid water and CH₄. See Kvamme et al. [29] for more details on the experiment. The container was made of plastic and as such was methane-wetting.

Nucleation studies for CO₂ have been published in a recent paper [19] and in order to save space we limit ourselves to estimates of growth as illustrated in Figure 5.

Based on the methods of Kvamme and Tanaka [12] it is easy to see that water molecules between a filled cavity and a neighboring empty cavity have higher librational movements and hydrogen bonding structures can temporarily break and let guest molecules move to the empty cavity. For CH₄ this can happen for all cavities while for the larger CO₂ only large cavity neighbors can participate in such transports. As such it is reasonable that the diffusivity coefficient for CO₂, through hydrate, is lower than the corresponding values for CH₄. And just as macroscopic heat transport is orders of magnitude faster than mass transport so is the local molecular heat transport involved in breaking and reforming hydrogen water hydrogen bonds during migration of a guest molecule between cavities.

The mechanical strength of hydrate films in a pore depends on the cross section area, thickness and strength of connection to the mineral surfaces. As discussed in several papers, for instance [32,33] and papers referenced therein, hydrate can never touch the surface of minerals. From a thermodynamic point of view this can easily be verified in calculated and experimental structures of water. The density of the first adsorbed water layer on minerals such as calcite, hematite, and kaolinite is on the order of three times the density of liquid water. The associated chemical potential is very low and far lower than is ever possible for water in hydrate. Even the average chemical potential of adsorbed water as calculated over the first five layers of water adsorbed on Hematite may be as low as 3.4 kJ/mole lower than liquid water [32]. The reason is the strong Coulombic interactions between water partial charges and charges on atoms in the mineral surfaces. For the same reason it is also easy to see that the regular distribution of partial charges on water hydrogens and oxygens in a hydrate crystal can never favorably match the distribution of positive and negative charges on the mineral surfaces.

What we have presented above provides a tool for analyzing how rapidly the hydrate film will grow, and for specific porosity and minerals, when and at what frequency we need to shock the system in order to break these hydrate films during the formation of hydrate in a specific sediment model (or real sediment from coring) in a laboratory set-up. For this purpose we need a higher resolution time scale of the growth, and in Figure 6 we plot the times needed to reach 100 nm and 10 micron (10,000 nm).

As discussed above it is a straightforward mechanical calculation to determine the magnitude of a mechanical pulse. If 100 nm is an acceptable thickness then the frequency of a mechanical shock pulse needs to be on the average 2 min. Depending on the laboratory set-up this pulse can be generated by sound, mechanical shaking of the sediments section etc. The important thing is that the sediments section is kept in place relative to circulation of fluids and that the mechanical pulse that is used to break the hydrate films is moderate according to the forces needed to break the hydrate films. The necessary force needed is sensitive to many factor like for instance pore size, geochemical issues (mineral types and salinity) as well as temperature and pressure. Calculations of necessary forces to break the hydrate films using mechanical pulses are therefore too case specific to discuss in a more general paper like this.

In the Ignik Sikumi test on hydrate production [8,9,11] large amounts of N₂ were added and the thermodynamic driving force for creating a new CO₂ hydrate was practically lost [9]. Instead of mechanically, breaking the hydrate films with efficient surfactants can keep the water/CO₂ interface free of blocking hydrate and permit for more continuous hydrate formation from the injected CO₂. Even methanol has surfactant properties but is too solvable in water to be interesting. And methanol is of course poisonous and not desirable to use. But some general insight can be found from a recent study [21].

2.3. Residual Thermodynamic Modeling of Hydrate Formation from Water and Dissolved Hydrate Former in Water

Hydrate can nucleate and grow from dissolved hydrate formers as discussed in a number of papers [9,10,14,15,16,18,19,21,29,30]:

$$\Delta G^{(H_2)}=\left[\begin{array}{l}{x}_{H_{2}O}^{H_2}\left({\mu }_{H_{2}O}^{H_2}(T,P,{\overrightarrow{x}}^{H_2})-{\mu }_{H_{2}O}^{water}(T,P,\overrightarrow{x})\right)\\ +{\displaystyle \sum _j{x}_j^{H_2}\left({\mu }_j^{H_2}(T,P,{\overrightarrow{x}}^{H_2})-{\mu }_j^{water}(T,P,\overrightarrow{x})\right)}\end{array}\right]$$

(18)

Chemical potential for water as solvent for CH₄ is not affected significantly by the very low concentrations of CH₄. For CO₂ dissolved in water the approximation in Equation (2) is considered to be good enough for the purpose of this work. Models for chemical potential of the hydrate former, which appear in the last term of Equation (18) are discussed in more detail in [19,20,21] and will not be repeated here. These chemical potentials also enter the canonical partition function for guest molecules in the two types of cavities through Equation (4). The actual film thickness that can be added through this route is limited for CH₄ and significant for CO₂.

There are several implications of the calculations illustrated in Figure 7. Dissociation of hydrate leads to a combination of gas bubbles and supersaturation of hydrate former in the surrounding water. This can then in turn lead to reformation of hydrate from dissolved hydrate former in concentration ranges between the blue and the black curve. If water surrounding the hydrate is diluted below the black curve then hydrate will dissociate. Also keep in mind that every unique concentration of hydrate former in between the blue and black curves will lead to a unique hydrate composition due to the chemical potential dependency of hydrate former as a function of concentration. Mathematically this means that an infinite number of different hydrates are possible. Stability will of course vary substantially.

Figure 8 illustrates a bit more clearly the regions for which it is not possible to grow hydrate from dissolved hydrate formers in water, i.e., regions for which the blue curves are lower in free energy. As expected, the differences in hydrate free energies between CH₄ hydrate formed from solution is slightly less stable than CO₂ hydrate formed from solution. But the differences are hardly significant compared to the corresponding differences in free energies formed heterogeneously from gas (or fluid) hydrate former phase and liquid water.

2.4. Hydrate Nucleation from Primary and Secondary Adsorption of Hydrate Formers on Mineral Surfaces

CH₄ is non-polar and cannot adsorb directly on charged mineral surfaces in competition with water, but CH₄ molecules can get trapped in structured water as illustrated in some recent MSc studies at the University of Bergen [34]. These studies indicates that methane gets dynamically trapped in structured water—typically from the third minimum in the water density as a function of distance from the calcite surface. As long there is a separate methane phase outside of the water layer in contact with the calcite a methane/water region, or phase, exists inside the water. Average water/methane distances in this separate phase are close to the characteristics of a small hydrate cavity. But the dynamics of the structure makes it hard to evaluate associated chemical potentials for methane with sufficient confidence. Work is in progress to overcome these challenges. Some illustrative animations from the MD simulations are available from the first author of the current work.

Direct adsorption is possible for polar and slightly polar hydrate formers like H₂S, which has a significant dipole moment and adsorbs well on the same type of minerals as water due to similarities in molecular structures. Although we do have a variety of chemical potentials for different interaction models for H₂S the limitation is that the results cannot be evaluated. Chemical potentials of adsorbed molecules cannot be measured directly but the density profiles of adsorbed molecules can be derived from measured structure factors. One of our examples for adsorbed CO₂ on Calcite [33] is illustrated in Figure 9 below. The finer details of the molecular simulations are too extensive in the context of this paper so the reader is directed to the original work [33], which is available for free download or can be provided by the first author of this paper. In the context of this paper we only need the calculated value for chemical potential as an example. The average value is of course sensitive to temperature and we can only solve for the hydrate stability pressure for this particular temperature, which is 17.0 bar when liquid water chemical potential is used for the water. This can be compared to a stability limit pressure of 15.7 bars for the same temperature and using CO₂ from gas. This is work in progress and more detailed studies of the chemical potentials of CO₂ and water that actually makes the nuclei is needed.

3. Discussion

Hydrates in sediments can form from many different phases. The most important routes to hydrates are, however, initial nucleation on the interface between water and hydrate former phase and nucleation towards mineral surfaces. These two routes are very different in the sense that hydrate cannot touch mineral surfaces but will be able nucleate towards mineral surfaces and then detach. Nucleation of hydrate on the interface between hydrate former phase and water, on the other hand, can lead to a dramatic slow-down in growth due to mass transport limitations of hydrate formers through the hydrate film. Visual observation or other monitoring techniques will therefore often lead to incorrect interpretations of induction times as nucleation times.

The problem of hydrate films that slow down hydrate growth often leads to a typical misunderstanding that CO₂ is not efficient for producing hydrates. The opposite might actually be the case with the right additives to the CO₂, but this is the subject for a subsequent paper, although the papers [19,20,21] give some insights.

Another frequent misunderstanding is that methane hydrate is more stable than CO₂ hydrate above a certain temperature, as illustrated in the combined Figure 2a,b above and Figure 10b below. Temperature and pressure are independent thermodynamic variables and the free energy is what determines phase stability and relative stability. Free energies of CH₄ and CO₂ hydrates along the equilibrium conditions are plotted in Figure 10b. There are consequences of these differences and this will be discussed further in a subsequent paper.

Even if the diffusivity of hydrate formers through hydrate films is faster than is represented by diffusivity coefficients on the order of 10⁻¹⁵ m²/s or less due to local dissociation mechanisms, a solid state transition [5,6] is highly unlikely for temperatures above freezing, as is also discussed elsewhere [11,16,19,35,36]. The thermodynamic driving force for such a transition would be the difference in free energy derived from Figure 10b. In such a system water is either hydrate or frozen water except in the mineral/water interface and the water hydrate interface, which are both nano scale in thickness. The only other exception is water trapped in small pores and unlikely to form hydrate due to mechanical strain in the hydrate lattice and a relatively large impact of mineral surfaces. The experiments at temperatures far below zero have limited value for relevant hydrate reservoirs since the solid state conversion is likely never going to play any role in reservoirs. Another confusion caused by these experiments is the inclusion of CO₂ in small cavities. It is theoretically possible but not thermodynamically feasible at realistic reservoir temperatures above zero.

The only feasible formation concept is, therefore, that injected CO₂ forms a new CO₂ hydrate and that the released hydrate heat of formation dissociates CH₄ hydrate. There is actually excess heat available for this process since the heat of formation of CO₂ hydrate is higher than is needed to dissociate CH₄ hydrate, as illustrated in Figure 11a above [19,20]. Adding large amounts of nitrogen in the CO₂, or using fuel gas dominated by nitrogen is not a good solution, as discussed earlier in several papers, e.g., [9], and illustrated in Figure 11b for 146 bar. At 276 K more than 30% CO₂ is needed in order to make new hydrate with water. At temperatures above slightly more than 286 K chemical potential for liquid water is also lower than chemical potential of water in hydrate from pure CO₂ and water, which is also illustrated by the equilibrium curve in Figure 2b. A limited amount of nitrogen can be useful for increasing the permeability of the injection gas. And since nitrogen can enter the 25% small cavities of a structure I hydrate primarily stabilized by CO₂, any addition of nitrogen reduces the thermodynamic driving force for making a new CO₂ hydrate. Large amounts of N₂ (or even air) will lead to dissociation of CH₄ hydrate because the N₂ phase (or air) is undersaturated with CH₄ relative to the chemical potential of CH₄ in the hydrate, and if the N₂ (or air) is dry relative to water saturation with reference to water in hydrate that will be an additional driving force for CH₄ hydrate dissociation. Figure 1 in the recent paper by Darnell et al. [37] is therefore wrong since the N₂ hydrate stability limit curve (green curve) is not relevant for the actual processes. The CO₂ hydrate curve, as well as the CO₂ fluid curves is also wrong since CO₂ goes through a phase transition. This has been discussed elsewhere in this work. See also Figure 2. See also Figure 10 and Figure 11 which illustrate the need for a thermodynamic analysis rather than plots in independent thermodynamic variable P and T. A simulation model without real phase transition analysis in terms of free energy based functions are therefore of limited value. Lack of heat transport analysis is also a limitation of this work. The free energy difference between CH₄ in the hydrate and in the outside N₂, as one example, will enter as the thermodynamic benefit but this driving force is implicitly coupled to mass—and heat—transport dynamics. See also Kvamme [19,20] and Kvamme et al. [21] for more details.

The key to efficient additives is illustrated in the paper by Kvamme et al. [21]. Small amounts of methanol in water will have a hydrate-activating effect due to its surfactant properties. First of all the concentration of methanol at the water/CH₄ interface will keep the interface open and free of blocking hydrate films. A second effect is increased transport rates for CH₄ into water. The third effect is a higher solubility of CH₄ in the water below the methanol concentrated layer. The hydrate activator effect of small amounts of methanol is well-known from industrial experience and the Tommeliten study [38]. Methanol is of course not a desirable additive since it will dissolve in water, and will eventually also disappear from the interface, but adding small surfactants to CO₂ is a good approach. Large water/CO₂ surfactants, on the other hand, are not useful additives for enhancing the injection of CO₂. Heavy surfactants on the order of molecular weights of some few thousand is like to cause problems of surfactant clogging and blocking of the flow pathways in between hydrate and mineral surfaces. Thermal breaking of hydrate films are also an alternative. The dynamics of hydrate film formation and growth are still useful tools together with energy balance for calculation of heating/cooling frequencies and necessary heat impulses for breaking desired film thickness.

4. Conclusions

Hydrates can nucleate from many phases. Except for cases of extremely small thermodynamic driving forces it is nano-scale phenomenon in both time and volume, as illustrated in this work using classical nucleation theory. The values for critical hydrate sizes obtained in this work are very much in accordance with earlier results that we have calculated from far more theoretically advanced approaches. Transport through hydrate films established through heterogeneous hydrate formation on the interface between liquid water and hydrate former phase is extremely slow. Simplified mass transport limitation models are able to verify the long induction times that are frequently observed in hydrate experiments. This type of information is critical since it can be used to design experimental equipment for efficient formation of hydrates in the laboratory. In particular it can be used to calculate the necessary frequencies for mechanical shock pulses for breaking mass transport-limiting hydrate films and ensure efficient conversion of water and hydrate-formers in the sediment pores. These calculations will also be useful if cycles of heating/cooling are used to break hydrate films and reform hydrate. The existence of these mass transport limiting films is also a critical issue in the use of CO₂ for combined production of CH₄ from hydrate and safe storage of CO₂ in the form of hydrate. Adding nitrogen in small amounts increases gas permeability but it is absolutely critical that the amount of nitrogen be balanced so that there is still a significant driving force for creating a new CO₂-dominated hydrate from the injection gas. Adding specially designed small surfactants will keep interfaces open and free of blocking hydrate films.