Study of the Critical Pore Radius That Results in Critical Gas Saturation during Methane Hydrate Dissociation at the Single-Pore Scale: Analytical Solutions for Small Pores and Potential Implications to Methane Production from Geological Media

Abstract

We examine the critical pore radius that results in critical gas saturation during pure methane hydrate dissociation within geologic porous media. Critical gas saturation is defined as the fraction of gas volume inside a pore system when the methane gas phase spans the system. Analytical solutions for the critical pore radii are obtained for two, simple pore systems consisting of either a single pore-body or a single pore-body connected with a number of pore-throats. Further, we obtain critical values for pore sizes above which the production of methane gas is possible. Results shown in the current study correspond to the case when the depression of the dissociation temperature (due to the presence of small-sized pores; namely, with a pore radius of less than 100 nm) is considered. The temperature shift due to confinement in porous media is estimated through the well-known Gibbs-Thompson equation. The particular results are of interest to geological media and particularly in the methane production from the dissociation of natural hydrate deposits within off-shore oceanic or on-shore permafrost locations. It is found that the contribution of the depression of the dissociation temperature on the calculated values of the critical pore sizes for gas production is limited to less than 10% when compared to our earlier study where the porous media effects have been ignored.

1. Introduction

Under appropriate conditions of pressure and temperature, water molecules can self-assemble into forming cavities/cages that can be stabilized by the enclathration of guest molecules with size, such that can fit within the cavities (e.g., methane, ethane, propane, hydrogen sulfide, carbon dioxide, nitrogen, hydrogen, tetrahydrofurane, etc.). Depending on the size of the enclathrated guest molecule different type of cavities can form (e.g., 5¹², 4³5⁶6³, 5¹²6², 5¹²6⁴, and 5¹²6⁸, which are formed primarily by pentagonal or hexagonal rings). By combining different number/type of cavities the formation of different hydrate structures is possible, with structures sI, sII, and sH being the more common encountered in nature (Sloan and Koh [1]). The natural hydrate deposits that are found within off-shore oceanic sediments or on-shore permafrost environments can contain significant amounts of methane gas, stored in pure or mixed methane hydrates, are the main focus of the current study. The properties of such systems have been discussed in a number of recent reviews [2,3,4].

Increasing industrial interest in hydrates arises, primarily, from flow-assurance related problems during oil/gas production and transportation (Kelland [5,6]), where hydrates are considered a nuisance. Furthermore, hydrates are encountered in worldwide accumulations (Collet et al. [7]; Chong et al. [8]), containing predominantly methane, both on-shore (under the permafrost), and off-shore (deposited within marine sediments). Hydrate accumulations amount to a significant reservoir of carbon. Note, however, that the size of the carbon reservoir is still under debate among the scientific community (Kvenvolden [9]; Buffet and Archer [10]; Milkov [11]; Klauda and Sandler [12]; Boswell and Collett [13]; Wallman et al. [14]; Pinero et al. [15]; You et al. [16]).

A significant effort has been undertaken during the recent years by the scientific community in order to delineate, with better accuracy, the amount of methane that is stored in hydrate deposits. This effort is important, since hydrates could represent in the future an important potential source of clean gas fuel (Buffet [17]). For the economically viable and environmentally friendly methane production from hydrate deposits, a number of issues need to be resolved adequately (Boswell and Collett [13]). Such issues include among others, the type of hydrate deposit, the amount of hydrate present in the deposit, and the total amount of hydrate that is accessible with the production approach (Boswell and Collett [13]; Walsh et al. [18]). Furthermore, the hydrate dissociation scheme should result in a process where hydrate dissociates under strict control, otherwise the uncontrolled dissociation can be the source of sever accidents and the possible complete destruction of the off-shore facilities. The economics of methane production through hydrate dissociation can be further improved by the simultaneous carbon dioxide sequestration within the hydrate deposits. In that case we can achieve the simultaneous energy production with environmental benefits.

On the other hand, the uncontrolled dissociation of oceanic or permafrost hydrates is also associated with the release of significant amount of methane to the atmosphere (Ruppel and Kessler [19]). Methane is a gas with a stronger green-house gas potential than carbon dioxide (Nisbet et al. [20]). According to Environmental Protection Agency (EPA) [21], methane is estimated to have a Global Warming Potential (GWP) equal to 28–36 over a period of 100 years, while the GWP becomes equal to 84–87 over a period of 20 years. Note that the particular values are according to the Intergovernmental Panel on Climate Change (IPCC) Fifth Assessment Report (IPCC [22]). In addition, the collapse of oceanic floor, due to uncontrolled dissociation of methane hydrate deposits, can create tsunamis that can be disastrous for coastal communities (Vanneste et al. [23]), while drilling through hydrate deposits can pose a serious hazard to offshore oil/gas platforms (McConnell et al. [24]).

The hydrate equilibrium primarily depends on the pressure and temperature conditions and the composition of the system under consideration. The presence of certain components in a hydrate-forming mixture can act as hydrate inhibitors or promoters, which can shift the equilibrium curve. Furthermore, the presence of electrolytes or confinement (i.e., as a result of hydrate formation within porous sediments) can also result in a shift in the equilibrium curve. The effect of confinement on the thermodynamic behavior of dissociation/formation of clathrate hydrates was addressed in a number of studies. Handa and Stupin [25] were the first to report experimental measurements of hydrate inhibition inside mesoporous silica. Handa and Stupin showed that for a constant temperature, the dissociation pressures of methane and propane hydrates in 70-$\stackrel{o}{\mathrm{A}}$-radius pores are higher than the dissociation pressures in the bulk phase. Additional experiments resulting in similar observations were reported among others by Uchida et al. [26,27], Wilder et al. [28], Sheshardi et al. [29], Smith et al. [30], and Anderson et al. [31]. Additionally, by Kang et al. [32,33], Aladko et al. [34], Seo et al. [35], Lee et al. [36], Zarifi et al. [37], Zhang et al. [38], and Liu et al. [39]. Additional experimental references and a detailed discussion on the quality of the available methane hydrate, three-phase equilibrium (i.e., Hydrate–Liquid water–Vapor; H–L_w–V), experimental data under confinement have been discussed in further detail in Tsimpanogiannis [40]. It should be noted, however, that for the case of mixed gas hydrates the inhibition behavior can become more complex, as has been discussed in more detail in the recent experimental study of Gambelli [41].

Theoretical studies by Clennell et al. [42,43], Henry et al. [44], and Anderson et al. [45] were introduced to explain the inhibition effect occurring within the small pores of the porous sediments (i.e., the lowering of the hydrate equilibrium temperature for a constant pressure, or the increasing of the hydrate equilibrium pressure for a constant temperature). The inhibition effect of a porous medium on the hydrate equilibrium curve was modeled by using a modified version of the van der Waals–Platteeuw-based thermodynamic models (van der Waals and Platteeuw [46]; Parrish and Prausnitz [47]; Holder et al. [48,49], where an extra term was added to account for the lowering of the water activity as a result of the presence of the porous medium. Such works were reported, among others by Clarke et al. [50], Klauda and Sandler [12,51], Peddireddy et al. [52], Sun and Duan [53], Song et al. [54], and Duan et al. [55], Li et al. [56], Zarifi et al. [37], Zhou et al. [57], and Azimi et al. [58]. An additional discussion and related references on the effect of the capillarity in narrow pores on the hydrate equilibrium conditions can be found in the recent study by De La Fuente et al. [59].

Tsimpanogiannis [60] examined, in an earlier study (denoted herein as part I), the dissociation of methane hydrate occurring inside a single pore, as well as the effect of pore size on the ability to produce the gas phase that results from methane hydrate dissociation. The particular study focused on obtaining critical pore sizes for which the production of methane gas was possible. Essentially, an effort was undertaken to identify when the resulting gas phase can arrive at the edge of the system (i.e., indicating arrival at the gas production end). In part I, the hydrate inhibition, which results from the presence of small pore-sizes, was ignored. The objective of the current study is to relax that particular assumption and expand the previous work with the incorporation of hydrate inhibition phenomena on the calculation of the critical pore sizes for which methane gas production is possible from geological media. The main purpose is to quantify the effect of hydrate inhibition, occurring due to confinement, on the calculation of the critical pore radius.

The paper is organized as follows: First, we present a brief description of the problem in the framework of simple geological settings, associated with hydrate deposits. Methane gas is the product of hydrate dissociation. Then, the mathematical analysis is introduced and, in a subsequent step, we present the corresponding results for the critical pore-body sizes that are essential to successfully produce methane gas from the porous systems considered. Finally, the conclusions are presented. The conclusions from the current study can provide tentative insights on potential methane production scenarios during hydrate dissociation from geological media.

2. Background

Consider a geological porous medium (i.e., deposited sediment) saturated with water that also contains a distributed solid methane hydrate phase. The geological medium can be located either beneath the ocean floor (off-shore) or under the permafrost (on-shore). Based on earlier studies, Davie and Buffet [61] reported that the hydrate saturation, S_H (defined as: $S_H=V_H/V_P$, where V_P is the total volume of the initial pore space, while V_H is the total volume that is occupied by the hydrate solid; essentially, S_H indicates the percentage of pore space that is occupied by hydrate), in oceanic porous deposits varies between values as low as 4–8% and as high as 20–35% of the available pore space. The particular values have been estimated based on seismic observations, resistivity measurements, and chlorinity profiles [61]. According to Davie and Buffet [61], the lower limit for hydrate saturation has been reported from studies at the Blake Ridge (Ocean Drilling Program (ODP); ODP Leg 164), while the higher limit for saturation has been reported from studies at the Cascadia Margin (ODP Leg 146). The lower hydrate saturations are often associated with hydrate conversion of biogenic methane that is produced locally (at temperatures lower than 50 °C). On the other hand, the higher hydrate saturations are associated with hydrate conversion of thermogenic methane that has migrated from deeper sediments where it has been produced at temperatures in the range 50–150 °C, which can also be combined with locally produced biogenic methane. Xie et al. [62] in a recent study presented a detailed report on hydrate saturation, obtained from natural gas hydrate deposits around the world, that have been studied during the recent years (e.g., Nankai Trough, Japan; Mallik Site Canada; Krishna-Godavari Basin, India; Bay of Bengal, India; Gulf of Mexico, USA; Ulleung Basin, South Korea; Shenhu Area, China; etc.). An extensive list of references regarding these sites can be found in Xie et al. [62]. Additionally, a detailed overview of hydrate occurrences in Europe has been reported by Minshul et al. [63]. Finally, You et al. [16] presented a comprehensive review of the different mechanisms that are associated with methane hydrate formation in geological porous media.

Depending on the fluid saturation in the geological sediments, different hydrate habits can develop (i.e., cementing at grain contact, grain coating, supporting matrix/grain, and pore filling). A detailed discussion has been presented by Waite et al. [64] and You et al. [16].

It is widely accepted that during hydrate formation in porous media with high initial water saturation, all the water inside a pore can solidify into hydrate with the exception of a thin (e.g., 5–50 nm for the case of mineral grains) adsorbed film that remains liquid (Clennell et al. [42]). Ice is strongly nonwetting and forms convex shapes within the pores that minimize the surface area (Clennell et al. [42]). This hypothesis was also confirmed by the two-dimensional visualization experiments of Tohidi et al. [65]. Thus, during the dissociation of hydrate, the liquid phase is the wetting phase while the produced gas is the non-wetting phase that occupies the center of the pore.

Hydrate dissociation occurs at appropriate equilibrium temperature, T_eq, and pressure, P_eq, conditions. For prevailing pressures in the porous medium that are above the dissociation pressure P_eq, the hydrate remains stable within the porous medium. However, when the pressure drops below P_eq, the hydrate dissociates into gas and either water (for temperatures that are above the freezing temperature of water) or ice (for temperatures that are below the freezing temperature of water).

After hydrate dissociation has occurred, we can identify two distinct steps in the process, with respect to the behavior of the produced gas phase. The gas phase partially occupies a given pore volume during the first stage, essentially staying far from the pore walls. For cases where the gas pressure, P_g, becomes higher than the sum of pressures in the surrounding liquid phase, P_l, plus the capillary pressure of the gas bubble, P_c ($P_c=P_g-P_l$), the growth of the gas bubble is possible. Such bubble growth can continue up to the point that the gas phase reaches the pore walls. Upon the arrival of the gas cluster at the pore walls, any further bubble growth is opposed by the capillary forces at the neighboring pore throats. This is considered the end of the first step. A significant consequence results from the first step. In case that within a particular pore the pressure is such that $P_g<P_c+P_l$, then the gas phase cannot grow any further. This implies that the gas phase is not capable of connecting with other neighboring gas clusters that are located in the proximity; therefore, remaining completely isolated in the particular pore. In such a case, the gas phase connectivity, as well as the gas phase production capability, are negatively affected. This particular problem can be remediated through gas mobilization. The gas phase can be mobilized in the presence of strong viscous/buoyancy forces (i.e., encountered close to a production well).

During the second step, if the pressure in the gas phase, P_g, is high enough to overcome the capillary threshold, P_c, of at least one of the neighboring pore throats (i.e., the largest pore throat), the gas phase can expand to occupy additional pores, until capillary pressure equilibrium is achieved within the system.

Hydrates are non-stoichiometric compounds. Namely, the hydrate crystal structure can retain stability, while at the same time certain cavities of the structure can be empty; therefore, the cavity occupancies are variable and depend on the prevailing pressures and temperatures. Figure S1a of the Supplementary Data shows the large and small cavity occupancies of sI methane hydrate along the hydrate equilibrium curve. Shown also (Figure S1b) is the three-phase (H–L_w–V) equilibrium pressure of the system methane + water system. Both the hydrate equilibrium pressure and the cavity occupancy values were calculated using the CSMGem simulator (Sloan and Koh [1]).

Natural hydrate deposits in sediments are often over-pressured. Namely, the prevailing pressure is higher than the three-phase (H–L_w–V) equilibrium pressure required for the hydrate to be stable. This is clearly shown in Figure S2a of the Supplementary Data, where the pressure conditions of several, known oceanic sites, with confirmed hydrate occurrence, are plotted along with the three-phase (H–L_w–V) equilibrium curve (based on data reported by Makogan et al. [66]). Shown also (Figure S2b) are the corresponding overpressures (Tsimpanogiannis and Lichtner [67]). There are cases where the pressure in the hydrate deposit is 450 bars above the hydrate equilibrium pressure. Earlier Monte Carlo (MC) simulations of sI methane hydrates (Papadimitriou et al. [68,69,70]; Lasich et al. [71]; Brumby et al. [72]) or other gas hydrate ([73,74]), have shown that for a given temperature, if the pressure is higher than the hydrate equilibrium pressure, an increase in the cavity occupancies, ${\theta }_{i,C{H}_4}$ (with i = S, L, indicating the small (5¹²) and large (5¹²6²) cavities, respectively, of the sI methane hydrate), can occur. Namely, it was concluded that ${\theta }_{i,C{H}_4}\to 1$, as the pressure increases. In addition, the calculations for methane that were reported in part I (Tsimpanogiannis [60]), indicate that the results when partial cavity occupancy is considered are very close to the results when 100% occupancy is considered. This is true at higher pressures/lower temperatures, where cavities are close to being fully occupied. Minor differences occur at lower pressures/higher temperatures where cavity occupancies deviate from 100%. Therefore, for the remainder of the paper we report results only for the case of 100% cavity occupancy. A detailed discussion of hydrate-related MC simulations has been presented in a recent review by Tsimpanogiannis and Economou [75].

3. Mathematical Modeling

As discussed in part I (Tsimpanogiannis [60]), the following condition needed to be satisfied for a gas bubble that results from the dissociation of hydrate under confinement, in order for the gas phase to grow:

$$P_g\ge P_c+P_l$$

(1)

where the pressure is denoted by P and the subscripts g and l denote the gas and liquid phases, respectively. The capillary pressure is denoted by $P_c$ and results from the presence of curvature in the gas–liquid interface. The capillary pressure can be, (i) either that of the gas bubble, for the case of a bubble located inside the pore-body, or (ii) that of the largest pore-throat, for the case of a gas bubble that has occupied most of the pore-body and resides inside a number of neighboring pore-throats.

In order to calculate the gas-phase pressure, $P_g$, we used the ideal gas law corrected with the compressibility factor, z, in order to account for the high pressures that the gas phase can reach for the problems under consideration. In the current study, the compressibility factor was calculated using the component-specific Equation of State (EoS) for methane introduced by Duan et al. [76].

Regarding the liquid pressure, $P_l$, which also corresponds to the pore pressure (used interchangeably in the current study) it was essential to obtain the pressure field in the domain of interest. Consider the schematic depicted in Figure 1, where a two-dimensional (2D) domain with a randomly distributed hydrate phase is shown. The entire domain is depicted as an effective porous medium, while the geometrical details of the pore network are not shown. Two characteristic zones can be identified; one zone where the pore pressure is $P_l>P_{eq}$, and, thus, no hydrate dissociation is possible, and a second zone where the pore pressure is $P_l\le P_{eq}$, which corresponded to the area where hydrate dissociation is possible. The hydrate-dissociation front is part of the second zone and the following equality applies at the front: $P_l=P_{eq}$. It should be noted that the depiction of the “dissociation front” as a simple straight line in Figure 1 is a simplistic approach (needed for obtaining analytical solutions to the problem), indicating the presence of a zone where the hydrate had dissociated completely and a zone where dissociation had not occurred. However, in a realistic geological system, a “dissociation front zone” exists instead. The width of the particular zone can depend on the width of the pore size distribution that is present in the system. Furthermore, it should be pointed out that, in a real porous medium with pores that follow a pore size distribution, the dissociation front is expected to be a fractal object, rather than a straight line (note that it can reduce to a straight line if all pores are of the same size or, alternatively, in the absence of the porous domain). The particular problem bears similarities with the “liquid-to-gas” phase change problem that is encountered in the drying/evaporation of porous media, which has been discussed in detailed in [77,78,79]. Additionally, it should be noted that P–T diagrams of real geological media, which also indicate the shift of the hydrate equilibrium curves, can be found in the studies of Liu and Flemings [80] and De La Fuente et al. [59].

Condition for Gas Bubble Growth and Production.

3.1. Case A: Ignoring the Effect of Small Pores on the Hydrate Equilibrium Conditions

In the current study, we followed the same notation as the one used in the case of part I (Tsimpanogiannis [60]). In a similar manner, we examined a methane hydrate-saturated system consisting of a single pore-body connected with a number of pore-throats (the coordination number, $\zeta$, indicates the number of pore-throats connected to the single pore-body), where the solid hydrate dissociated into liquid water and methane gas. We were interested in identifying a critical value (denoted with subscript c) for the pore-body radius, $r_{pc,Si}$, such that the gas saturation after dissociation in the porous medium under examination reached a critical gas saturation ($S_{g,Si}\to S_{gc,Si}$). Essentially, the amount of hydrate that was present within the pore-body was adequate to release enough methane gas (upon hydrate dissociation), such that the gas phase could reach the outside boundaries of the porous medium. Therefore, the production of the gas, resulting from hydrate dissociation, could be practically realized. Alternatively, for the case where $r_p\ge r_{pc,Si}$, the production of gas was possible. However, in the case of $r_p<r_{pc,Si}$, not enough gas arrived at the outer boundaries of the porous system and, consequently, an inhibition of gas production was observed. Gas bubble mobilization due to buoyancy or viscous gradients was ignored in the current study. Note that the viscous gradients close to the production well could be important and result in the mobilization of gas bubbles.

In part I, the following general solution was obtained for the critical pore-body radius, $r_{pc,Si}$. The particular solution could be applied to two different porous media systems denoted with Si (with I = 1 or 2). The systems that were examined consisted of a single pore-body and a number of pore-throats. In particular, S1 corresponds to a system consisting of only a single pore-body, while S2 corresponds to a system consisting of a single pore-body connected with a number, $\zeta$ (i.e., coordination number), of pore throats.

$$r_{pc,Si}=\frac{2{\gamma }_{gl}}{z\lambda R{T}_{eq}^{bulk}}\cdot \frac{S_{gc,Si}}{\left(1-b{S}_{gc,Si}\right)}\cdot F_{Si}\left(S_{gc,Si}\right)$$

(2)

where ${\gamma }_{gl}$ is the gas (methane)/liquid (water) interfacial tension, z is the gas compressibility factor, R is the universal gas constant, $\lambda$ is a parameter that is proportional to the amount of gas that is present in the hydrate cages, $b=\frac{P_l}{z\lambda RT}$, and $F_{Si}\left(S_{gc,Si}\right)$ is a function that depends on the geometric characteristics of the porous system and discussed in detail in part I. Figure 2 shows the effect of temperature (and pressure) on the gas–liquid interfacial tension, ${\gamma }_{gl}$ (Figure 2a), and the compressibility factor, z (Figure 2b). The calculations shown in Figure 2 correspond to the bulk phase and are essential for Equation (2). The calculations shown in Figure 2 are performed at P, T conditions that correspond to the three-phase (H–L_w–V) equilibrium conditions for the binary system: methane + water.

The hydrate dissociation temperature in the bulk is a function of pressure, $T_{eq}^{bulk}=f\left(P_{eq}^{bulk}\right)$. For the case of methane hydrate it can be calculated from a simple expression, $P_{eq}^{bulk}=\exp \left(m+\frac{n}{T_{eq}^{bulk}}\right)$ (where m, n are component-specific constants), as the one discussed by Holder et al. [49], or a more refined correlation, such as the one proposed by Moridis [81]. Alternatively, a method that couples the van der Waals–Platteeuw method (van der Waals and Platteeuw [46]) with a cubic equation of state (e.g., the CSMGem simulator discussed by Sloan and Koh [1]) or a Statistical Associating Fluid Theory (SAFT) type equation of state (e.g., see the discussion by Dufal et al. [82]; El Meragawi et al. [83]), can be utilized.

It was shown in part I that meaningful solutions (i.e., $r_{pc,Si}>0$) are obtained as long as:

$$1-b{S}_{gc,Si}>0\Rightarrow S_{gc,Si}<\frac{z\lambda RT}{P_l}$$

(3)

A detailed discussion regarding the behavior of $r_{pc,Si}$, in the absence of any confinement effects, has been presented in part I, and therefore, it is not discussed further in the current study. The growth behavior of the gas bubbles, resulting from hydrate dissociation, has been explored in detail using the hydrate pore-network simulator that has been developed by Tsimpanogiannis and Lichtner [84,85]. The particular problem has also similarities with the earlier pore-network studies on the critical gas saturation during primary oil production, discussed by Du and Yortsos [86], and Tsimpanogiannis and Yortsos [87].

3.2. Case B: Incorporating the Effect of Small Pores on the Hydrate Equilibrium Conditions

When hydrates are located inside small pores (e.g., confined within a porous sediment with a pore radius smaller than 100 nm) a depression of the hydrate dissociation temperature occurs as a result of the confinement (Clennell et al. [42]; Henry et al. [44]). This shift (inhibition) of the equilibrium conditions becomes stronger as the pore size decreases. The temperature depression of the hydrate equilibrium value ($\Delta T_{eq}=T_{eq}^{pore}-T_{eq}^{bulk}$), when compared to the bulk conditions ($T_{eq}^{bulk}$), that results from the presence of a porous medium can be described by the Gibbs-Thompson equation [40]:

$$\frac{T_{eq}^{pore}-T_{eq}^{bulk}}{T_{eq}^{bulk}}=-\kappa \left(\frac{{\gamma }_{sl}\cos {\vartheta }_{sl}}{{\rho }_s\Delta H_s}\right)$$

(4)

where ${\gamma }_{sl}$ is the solid (hydrate)/liquid (water) interfacial tension, ${\vartheta }_{sl}$ is the contact angle between the solid hydrate–liquid water interface and the pore wall (${\vartheta }_{sl}$ is equal to 180° when measured through the solid hydrate and further assuming that the pores are “water-wet”, while the solid hydrate is nonwetting, thus $\cos {\vartheta }_{sl}=1$). The effect of having solid pore walls with different chemistry can be reflected through the change in wettability and in particular through a different contact angle, ${\vartheta }_{sl}$. An extensive discussion and experimental, as well as molecular, simulation data, on the effect of surface chemistry on contact angles, can be found in the recent review by Liang et al. [88]. ${\rho }_s$ is the solid (i.e., hydrate) density, $\kappa$ is the solid–liquid interface curvature, and $\Delta H_s$ is the latent heat of dissociation/melting. Tsimpanogiannis et al. [89] reported a detailed review on experimental and computational methods for the $\Delta H_s$ of sI methane hydrate. The interface curvature between the solid and liquid can be defined as follows:

$$\kappa =\left(\frac{1}{r_1}+\frac{1}{r_2}\right)$$

(5)

where r₁ and r₂ are the two principle radii of curvature with directions that intersect orthogonally at that point.

For the particular case of a cylindrical pore, Anderson et al. [45] argued that one of the principle radii (i.e., the one along the length of the cylindrical pore) is infinite, when dissociation occurs. In that case, hydrate dissociation throughout the entire pore-length can occur instantly. Therefore, for the case of cylindrical pores, the solid–liquid interface curvature is equal to:

$$\kappa =1/r$$

(6)

while for the case of (i) hydrate formation in cylindrical pores and (ii) hydrate dissociation/formation in spherical/cubic pores, it is equal to:

$$\kappa =2/r$$

(7)

It has been shown by Tsimpanogiannis [40] that the hydrate dissociation temperature inside the porous medium is related to the dissociation temperature in the bulk through the following equation:

$$T_{eq}^{pore}=T_{eq}^{bulk}\left(1-\kappa \alpha \right)$$

(8)

where we have defined the parameter $\alpha$, as:

$$\alpha \equiv \frac{{\gamma }_{sl}\cos {\vartheta }_{sl}}{{\rho }_s\Delta H_s}$$

(9)

As it has been discussed in detail in Tsimpanogiannis [40] the temperature depression, $\Delta T_{eq}$, has also been associated with experimental measurements of three-phase equilibrium conditions under confinement in porous media with pore diameter, d, through the following expression:

$$\frac{\Delta T_{eq}}{T_{eq}^{bulk}}=-{\beta }^{*}\cdot \frac{1}{d}$$

(10)

${\beta }^{*}$ the slope of the line that goes through the origin of the axes in the plot $\Delta T_{eq}/T_{eq}^{bulk}$ vs. $1/d$, for the available experimental data of hydrate equilibrium. A value of $0.1695\pm 0.0006$ (nm) has been reported for ${\beta }^{*}$ (Tsimpanogiannis [40]). By combining Equations (4), (8) and (9), we observe that $\alpha =\frac{{\beta }^{*}}{\kappa \cdot d}$, which finally reduces to $\alpha =\frac{{\beta }^{*}}{4}$ for hydrate dissociation in spherical or cubic pores or $\alpha =\frac{{\beta }^{*}}{2}$ for hydrate dissociation in cylindrical pores (see also the detailed discussion in Tsimpanogiannis [40]).

Tsimpanogiannis [40] has examined the case when the three parameters ${\gamma }_{sl}$, $\Delta H_s$ and ${\rho }_s$ are functions of temperature and their dependence is such that the parameter $\alpha$ is constant (i.e., $\alpha \equiv const$) and independent of temperature. The temperature-dependence of ${\rho }_s$ has been obtained from extensive Monte Carlo simulations (Papadimitriou et al. [68,69,70,90]), The temperature-dependence of $\Delta H_s$ has been obtained from molecular dynamics (MD) simulations (Tsimpanogiannis et al. [89]), while the temperature-dependence of ${\gamma }_{sl}$ has been obtained by a hybrid method that combined the aforementioned calculated parameters with experimental measurements of three-phase equilibrium conditions under confinement (Tsimpanogiannis [40]). Figure S3 of the Supplementary Data shows the temperature dependence of the three parameters, ${\rho }_s$, $\Delta H_s$ and ${\gamma }_{sl}$, for the case of sI methane hydrate.

On the other hand, a number of studies have determined the engineering assumption that the values for the particular three parameters (${\rho }_s$, $\Delta H_s$ and ${\gamma }_{sl}$) are independent of temperature in the short range of interest (275–300 K) to gas hydrate equilibria.

Table 1. Comparison of calculated values for ${\gamma }_{sl}$ (methane hydrate–liquid water), using experimental measurements reported by different studies (adopted from Tsimpanogiannis [40]).

Study	${\mathit{\gamma }}_{\mathit{s}\mathit{l}}(\mathbf{mN}/\mathbf{m})$	${\mathit{\rho }}_{\mathit{s}}(\mathbf{kg}/{\mathbf{m}}^3)$	$\mathbf{\Delta }{\mathit{H}}_{\mathit{s}}(\mathbf{kJ}/(\mathbf{mol}{\mathbf{CH}}_4\left)\right)$
Uchida et al. [27]	17 ± 3	914	54.2
Uchida et al. [27] (a)	34 ± 6	914	54.2
Uchida et al. [26]	39	914	45.92
Uchida et al. [26] (b)	78	914	45.92
Uchida et al. [26] (c)	32.72	914	45.92
Anderson et al. [45]	32 ± 1 to 35 ± 1	907–914	53.2–57.7
Anderson et al. [45] (d)	32 ± 3	–	–
Zarifi et al. [37]	31 ± 3	914	57.8
Liu et al. [39]	32 ± 2	914	54.2

(a) Corrected by Anderson et al. [45], accounting for the curvature of dissociation. (b) Corrected by Tsimpanogiannis [40], accounting for the curvature of dissociation. (c) Corrected by Tsimpanogiannis [40], accounting for (i) the curvature of dissociation and (ii) using a refitted value of the slope ${\beta }^{*}$ (${\beta }^{*}=-0.1948 (\mathrm{nm})$) of the line that goes through the origin of the axes. (d) Average value.

shows the range of values that have been used in various studies for the case of methane hydrate. Based on the values that are reported in Table 1 we observe that $\alpha$ can range between $6.7699\times {10}^{-11}$ and $1.1489\times {10}^{-10}$, if we select the min/max values for ${\rho }_s$, $\Delta H_s$ and ${\gamma }_{sl}$ appropriately. The two aforementioned distinct scenarios with respect to the parameters ${\rho }_s$, $\Delta H_s$ and ${\gamma }_{sl}$ (i.e., temperature -dependent or -independent), both result to having a constant value for the parameter $\alpha$.

Figure 3 shows the effect of pore radius on the temperature depression (absolute values) during methane hydrate dissociation in porous media. Results are shown for the temperature of 275 K and the two aforementioned scenarios regarding the parameters ${\rho }_s$, $\Delta H_s$, and ${\gamma }_{sl}$ (i.e., constant and variable). The corresponding case for the temperature of 310 K is shown in Figure S4 of the Supplementary Data. The curves that are shown in Figure S4 of the Supplementary Data shifted by approx. 3K to higher values than those shown in Figure 3. The limited shift observed between 275 K and 310 K is a result of the small variability for the parameters of Equation (4) that are shown in Figure S3 of the Supplementary Data, in the same temperature range. From the particular figure we can observe that for $T_{eq}^{bulk}=275\mathrm{K}$ there is a shift to lower values by, $\Delta T_{eq}=1\mathrm{K}$, for the case of hydrate dissociation in a pore with an approximate radius of $2.5\times {10}^8\mathrm{m}$, while higher $\Delta T_{eq}\prime s$ result for smaller pore-body sizes. The corresponding pore-body radius for $T_{eq}^{bulk}=310\mathrm{K}$ is approximately $2.8\times {10}^8\mathrm{m}$. It is evident; therefore, that for small pore sizes the depression in temperature needs to be taken into account in the calculations. Sahou et al. [91] estimated that the temperature shift (towards lower values) for the case of bulk methane hydrate equilibrium conditions at 8.8 MPa and 9.8 °C, is equal to 4.8 °C. For the particular calculation, Sahou et al. used the constant values for ${\rho }_s$, $\Delta H_s$ and ${\gamma }_{sl}$, reported by Clennell et al. [42] for the case of ice. For the same bulk equilibrium P, T conditions, we calculated in the current study (using variable parameters for ${\rho }_s$, $\Delta H_s$ and ${\gamma }_{sl}$) a temperature shift to lower values that is equal to 2.52 K.

In part I, we reported results obtained from an analytical solution for the case that the depression of the dissociation temperature was ignored and thus P_eq was independent of the pore size. In this paper the aforementioned assumption is relaxed and the depression of the dissociation temperature is taken into account and therefore P_eq depends on the pore size.

By introducing Equation (4) into (2), we can obtain an expression for the critical pore-body radius which consists of two terms:

$$r_{pc,Si}^{total}=r_{pc,Si}^{}+r_{pc,Si}^{shift}$$

(11)

The first term at the right-hand side (RHS) of Equation (11) is the solution that corresponds to the problem in the bulk (discussed extensively in part I), while the second term at the RHS is the correction that corresponds to the effect of the porous medium. The correction term, $r_{pc,Si}^{shift}$, is given as:

$$r_{pc,Si}^{shift}=\frac{z\lambda R{T}_{eq}^{bulk}}{z\lambda R{T}_{eq}^{bulk}-P_l{S}_{gc,Si}}\cdot \left(f\alpha \right)\equiv {\mathsf{\Omega }}_{c,Si}^{}\cdot (f\alpha )$$

(12)

where f is a constant that has a value equal to 1 for the case of hydrate dissociation in cylindrical pores, while a value equal to 2 for the case of hydrate dissociation in spherical or cubic pores. In Equation (12) we have also defined the function ${\mathsf{\Omega }}_{c,Si}^{}$ as follows:

$${\mathsf{\Omega }}_{c,Si}^{}\equiv \frac{z\lambda R{T}_{eq}^{bulk}}{z\lambda R{T}_{eq}^{bulk}-P_l{S}_{gc,Si}}$$

(13)

4. Results and Discussion

4.1. System S1: Single Pore-Body

Figure 4 shows the correction term, $r_{pc,S1}^{shift}$, as a function of temperature, calculated along the three-phase (H–L_w–V) equilibria line for the methane + water binary system. We observe that the correction term has a very week dependence on temperature. Essentially, there is a change of only 0.002% between the minimum (at 275 K) and maximum (at 305 K) values of $r_{pc,S1}^{shift}$. Furthermore, for the case of S1 we observe that the term ${\mathsf{\Omega }}_{c,S1}^{}$ is approximately equal to one (i.e., for all engineering purposes), as is shown in Figure 5, where we plot $({\mathsf{\Omega }}_{c,S1}^{}-1)$ as a function of temperature. We can observe that $({\mathsf{\Omega }}_{c,S1}^{}-1)$ is in the range $1\times {10}^{-6}–3\times {10}^{-5}$, for the temperature range 275 K–305 K. Therefore, the following approximation can provide accurate results for engineering-type calculations:

$$r_{pc,S1}^{shift}\approx f\cdot \alpha$$

(14)

Figure 6a shows for comparison purposes the calculated $r_{pc,S1}^{}$ as a function of temperature (along the three-phase, H–L_w–V, equilibrium line) without accounting the confinement effect (discussed in part I). Figure 6b shows the contribution of the correction term (i.e., $r_{pc,S1}^{shift}$) that is calculated from Equation (12), that originates from confinement. The figure shows the ratio,

$${\mathsf{\Phi }}_{c,S1}^{}\equiv \left(r_{pc,S1}^{shift}/r_{pc,S1}^{total}\right)$$

(15)

plotted as a function of temperature (and the corresponding hydrate equilibrium pressure), for the case of a spherical, cylindrical or cubic pore-body, without any liquid film remaining at the pore-body walls (i.e., ${\omega }_p=1$). Essentially, it shows the value of $\left(100*{\mathsf{\Phi }}_{c,S1}^{}\right)$, namely, the per cent (%) contribution of confinement on the total value of $r_{pc,S1}^{}$.

The Figure 6a,b show two distinct families of curves. In particular, pore-bodies that were located exactly on the hydrate-dissociation front $(P_l=P_{eq}^{bulk})$ belong to the first family of curves and correspond to the upper limit (denoted as solid lines). On the other hand, pore-bodies that were located at points where the pressure was observed to be at a minimum ($P_l=0$) belong to the second family of curves and correspond to the lower limit (denoted as dashed lines). In addition, for pressures such that $P_{eq}^{bulk}>P_l>0$, an infinite number of curves could be identified between the aforementioned upper and lower limits. Essentially, such pressures can be encountered within the porous domain between the location of the hydrate dissociation front and the location of the initial lowering of the pressure (e.g., at the production well for the case, a pressure depletion scheme was followed).

We observed the following regarding the shift in the hydrate equilibrium temperature that resulted from the confinement:

• At pores located on the hydrate-dissociation front (i.e., where $P_l=P_{eq}^{bulk}$), a shift of less than 2% is observed, with the shift magnitude having values in the following order: cubic > cylindrical > spherical.
• At pores located at points with the minimum pore pressure (i.e., where $P_l=0$) for well production via a pressure depletion scheme, a shift of less than 5% was observed, with the shift magnitude having values in the following order: cubic > cylindrical > spherical.

Figure 7b shows the effect of a thin liquid film that remains on the pore-body walls on the value of $\left(100*{\mathsf{\Phi }}_{c,S1}^{}\right)$, as a function of temperature (and the corresponding hydrate equilibrium pressure) for the case of a spherical pore-body. The existence of thin liquid films remaining on the pore-body walls (i.e., ${\omega }_p<1$), results in shifting the correction term to higher values. In particular, for the case of a spherical pore-body located on the hydrate-dissociation front, to less than 0.8% (for ${\omega }_p=1.00$), 1.0% (for ${\omega }_p=0.90$), and 1.4% (for ${\omega }_p=0.80$). The corresponding values for the case of a spherical pore-body, located far from the hydrate-dissociation front (i.e., $P_l=0$), are less than 2.5% (for ${\omega }_p=1.00$), 3.5% (for ${\omega }_p=0.90$), and 4.0% (for ${\omega }_p=0.80$). Similarly, the top panel of Figure 7 shows for comparison purposes the calculated $r_{pc,S1}^{}$, plotted as a function of temperature (along the three-phase, H–L_w–V, equilibrium line), without accounting the confinement effect (discussed in part I).

4.2. System S2: Single Pore-Body Connected with a Number of Pore-Throats

In the current Section, we consider a porous system consisting of a single pore-body (spherical, cylindrical, or cubic) that is connected with four (i.e., $\zeta =4$) cubic pore-throats (i.e., system S2). The remaining geometrical parameters are defined and have values that are similar to cases considered in our previous study (part I). In particular, the parameter, is defined as the ratio of the radius of the pore-body ($r_p$) over the radius of the pore-throat ($r_t$): $\beta =r_p/r_t\mathrm{with}\beta >1$. Additionally, the parameter $\upsilon$ is defined as the ratio of the length of the Unit Porous Block (UPB), $L_{UPB}$ (defined in part I), over the radius of the pore-body as follows: $\upsilon =L_{UPB}/r_p$.

Figure 8 shows the correction term, $r_{pc,S2}^{shift}$, as a function of temperature, calculated along the three-phase (H–L_w–V) equilibria line for the methane + water binary system. We observe that the correction term for pores located far from the hydrate-dissociation front (i.e., where the pore pressure is a minimum; namely, where $P_l=0$) is essentially constant and independent of temperature. On the other hand, for pores located on the hydrate-dissociation front the correction term has strong temperature dependence. In particular, the correction term changes almost four orders of magnitude between 273 K and 295 K. Such behavior is a direct result of the term ${\mathsf{\Omega }}_{c,S2}^{}$, which is shown in Figure 9, plotted as a function of temperature. Therefore, a similar approximation to Equation (14) cannot be used for the case of system S2 and the complete Equation (12) needs to be used for the correction term, $r_{pc,S2}^{shift}$.

Figure 10a shows for comparison purposes the calculated $r_{pc,S2}^{}$ as a function of temperature (along the three-phase, H–L_w–V, equilibrium line) without accounting the confinement effect (part I). On the other hand, Figure 10b shows the contribution of the correction term (i.e., $r_{pc,S2}^{shift}$) that is calculated from Equation (12), that originates from confinement. The figure shows the ratio, $\left[{\mathsf{\Phi }}_{c,S2}^{}\equiv \left(r_{pc,S2}^{shift}/r_{pc,S2}^{total}\right)\right]$ as a function of temperature (and the corresponding hydrate equilibrium pressure), for the case of a spherical, cylindrical or cubic pore-body, without any liquid film remaining at the pore-body walls (i.e., ${\omega }_p=1$).

We observe that, even though the correction term, $r_{pc,S2}^{shift}$, increases by more than three orders of magnitude between 273 K and 305 K, the term $r_{pc,S2}^{}$, increases also by approximately the same orders of magnitude. Therefore, the ratio ${\mathsf{\Phi }}_{c,S2}^{}$, remains essentially constant in the temperature ranges 273 K–293 K, while increases slightly for temperatures higher than 293 K. Notice, however, that in all cases the correction term, $r_{pc,S2}^{shift}$, is less than 0.5% of the total value of $r_{pc,S2}^{total}$. It is also observed that pores that are located at the hydrate dissociation front exhibit the same shift with those pores that are located at points of minimum pressure. Furthermore, it should be pointed out that there is an asymptotic behavior for the case of pores located at the front (see Figure 10a), and it’s effect is clearly shown in Figure 10b, where we observe that the resulting curves are plotted only up to 291 K (spherical pore-body), 293 K (cylindrical pore-body), and 295 K (cubic pore-body), respectively, for the three pore-body geometries that are considered in the current study. For the three particular cases the magnitudes of the shift are having values in the following order: cubic > cylindrical > spherical.

5. Conclusions

In this study we examined the dissociation of sI methane hydrate, occurring inside simple model porous media with small pore sizes. We performed analytical calculations for the critical pore sizes that are required in order for the gas phase to reach the critical gas saturation, and therefore that can be produced. The current study is an extension of our previous work (Tsimpanogiannis [60]), where hydrate dissociation in larger pores was considered. Namely, in our previous work we have ignored the depression of the temperature of methane hydrate dissociation, which is a result of confinement within the porous medium. Note that the particular effect becomes stronger as the pore size decreases. In the current calculations, on the other hand, we considered the effect of the depression of the temperature of hydrate dissociation, and moved forward to quantify its effect on the calculated critical pore sizes.

Two different types of fundamental porous media have been considered. The first porous system (S1), consisted only of a single pore-body, while the second porous system (S2), consisted of a single pore-body connected with a number of pore-throats. Different geometries have been considered for pore-bodies (spherical, cylindrical, and cubic), while only cubic pore-throats have been considered. It is expected, however, that the use of other pore-throat geometries will not alter significantly the main conclusions.

The following observations can be summarized:

• As pore sizes decrease, the shift in the three-phase (H–L_w–V) equilibrium temperature becomes more important.
• In the current study the parameters ${\gamma }_{sl}$, $\Delta H_s$ and ${\rho }_s$ are considered to be functions of temperature and their dependence is such that the parameter, $\alpha \equiv \frac{{\gamma }_{sl}\cos {\vartheta }_{sl}}{{\rho }_s\Delta H_s}$, is constant (i.e., $\alpha \equiv const$) and thus, independent of temperature. The temperature dependence has been obtained from Monte Carlo simulations (${\rho }_s$), molecular dynamics simulations ($\Delta H_s$), and hybrid methods that combine MD/MC simulations with experimental measurements (${\gamma }_{sl}$).
• For system S1, when accounting for the effect of confinement on the H–L_w–V equilibrium $T_{eq}$ shift, results in a contribution up to 8% of the total value of $r_{pc,S1}^{total}$. The shift magnitude was found to have values in the following order: cubic > cylindrical > spherical.
• For system S2, when accounting for the effect of confinement on the H–L_w–V equilibrium $T_{eq}$ shift, results in a contribution up to 1% of the total value of $r_{pc,S2}^{total}$. The shift magnitude was found to have values in the following order: cubic > cylindrical > spherical.
• The current analysis indicated that the contribution of the effect of the hydrate equilibrium $T_{eq}$ shift (due to confinement in small pores) on the critical value of the pores, $r_{pc,Si}^{total}$, in order to be able to produce the gas phase was less than 10%. The contribution to overcome the capillary thresholds of neighboring pore-throats accounted for the remaining 90%.
• Based on the previous observations we can conclude that for hydrate deposits within coarse-grained sediments or fractured systems, there is not going to be a significant negative effect, resulting from accounting for the confinement in small pores.
• Due to the lack of any experimental measurements (to the best of our knowledge) the provided discussion is limited only to the analytical solutions.
• The current analysis does not take into consideration the possible mobilization of the gas phase due to buoyancy or viscous gradients. While the contribution of buoyancy gradients is expected to be insignificant for such small gas clusters, the contribution of viscous gradients can be important, particularly in the vicinity of production wells. When accounting for such gradients, the gas production from hydrate deposits within formations with smaller pore sizes can become possible, due to gas mobilization.