A new study of multi-phase mass and heat transfer in natural gas hydrate reservoir with an embedded discrete fracture model

Studies of the hydrate cores have shown that natural fractures can be frequently observed in hydrate reservoirs, resulting in a fracture-filled hydrate. Therefore, it is highly necessary for industries to predict the gas well productivity of fracture-filled hydrate reservoirs. In this work, an embedded discrete fracture model is applied to characterize the natural fractures of fracture-filled gas-hydrate reservoirs. The non-linear mass and energy conservation equations which are discretized with the finite-difference method are solved by the fully implicit approach, and the proposed model is justified by a commercial simulator. On the basis of the proposed model, we investigate the influences of natural fractures, fracture conductivity, and hydrate dissociation rate on the gas well productivity and the distributions of pressure, temperature, and hydrate saturation. The simulation results show that hydraulic and natural fractures exert significant impacts on the gas well productivity of the fracture-filled hydrate reservoirs, and the cumulative gas production is increased by 45.6% due to the existence of the connected natural fractures. The connected natural fractures can impose a more important influence on the gas well productivity than the unconnected natural fractures. The cumulative gas production is increased by 6.48% as N_nf is increased from 2 to 50, whereas the increase is 43.38% as N_{f_con} is increased from 0 to 4. In addition, A higher hydraulic fracture conductivity can be more favorable than a higher natural fracture conductivity for improving the gas well productivity, and a higher hydrate dissociation rate can lead to a lower temperature along fractures due to a more noticeable reduction of solid hydrate. This study provides a theoretical basis for developing fracture-filled hydrate reservoirs efficiently in the future.

1 Introduction

Natural gas hydrate (NGH) is formed from gas and water under particular pressure and temperature (Sloan et al., 2007), and 150–180 m³ of gas can be approximately released from one volume of NGH at standard conditions (Carroll, 2020). At present, four different classes of gas hydrate reservoirs have been found in nature (Moridis, 2003; Moridis, 2004; Moridis et al., 2005). Class 1–2 hydrate reservoirs are composed of the overburden layer, gas hydrate layer, gaseous hydrocarbon or water layer, and under burden layer. Class 3 hydrate reservoir distinguishes from Classes 1 and 2 hydrate reservoirs due to the absence of the underlying mobile fluid layer. Class 4 hydrate reservoir is characterized by a low hydrate saturation (<10%) and disperse hydrate distribution (Moridis, 2003; Moridis, 2004; Moridis et al., 2005). Abundant natural gas hydrate resources have been observed in ocean margins and permafrost regions, and the estimated global technically recoverable gas hydrate is approximately 10¹⁵–10¹⁸ m³ (Carroll, 2020; Li et al., 2021b). At present, the depressurization method has been widely applied for developing gas hydrate (Liu et al., 2017). Numerous studies have been conducted to evaluate the potential of developing hydrate reservoirs ascribing to the fact that gas hydrates are known as one of the most promising clean energies (Boswell et al., 2012).

Developing gas hydrate is complicated which involves multiphase flow (gas-water) and hydrate phase transition. Moreover, thermal convection and hydrate decomposition occur during the development of gas-hydrate which can significantly affect temperature distribution and gas well productivity. Therefore, a comprehensive insight into the effect of hydrate phase transition on multi-phase mass and heat transfer of gas-hydrate reservoirs is one of the key issues for optimizing hydrate development. The kinetic model proposed by Kim et al. (1987) was widely used to calculate the hydrate decomposition rate caused by depressurization. In their model, the hydrate decomposition rate is related to hydrate surface area, gas fugacity, and temperature; Pooladi-Darvish and Hong. (2004)compared the effect of thermal conduction and thermal convection on the productivity of gas hydrate wells. Their simulation results showed that thermal conduction has a more significant influence in affecting heat transfer than thermal convection; Sun et al. (2005); Esmaeilzadeh et al. (2008) both established a numerical model considering the one-dimensional fluid flow to evaluate gas well productivity in developing gas-hydrate reservoirs by depressurization method, and the results of their studies showed fluid enthalpy and heat transfer coefficient are the most sensitive parameters in developing gas hydrate. Moridis, (2003), Moridis, (2004), Moridis et al. (2005). developed a new numerical submodule named Tough+Hydtate to study the mass and heat transfer of different types of hydrate reservoirs excluding Class 4. They claimed that Class 1 and Class 2 hydrate reservoirs hold significant promise as potential energy sources. Tonnet and Herri (2009) investigated the gas production and dissociation kinetics of hydrate using experimental and numerical methods. Their studies demonstrated that the mass transfers of gas and water impose a remarkable influence on hydrate decomposition for low permeability gas hydrate reservoirs. Zhao et al. (2014) developed a two-dimensional numerical model to analyze the effects of sensible heat transfer on gas well productivity of hydrate reservoirs. The calculated results of their works showed that the hydrate will dissociate rapidly at a high value of thermal conductivity in the early developing stage.

In all the aforementioned studies, the gas-hydrate reservoirs are assumed as a mono-porous media where the gas hydrate can only be observed in the pore space. However, the most recent studies show that natural fractures commonly exist in gas-hydrate reservoirs (Wang et al., 2011; Wang et al., 2019; Ning et al., 2020). Since natural fractures can significantly influence the gas well productivity of hydrate reservoirs, it is highly necessary for industries to develop a mathematical model to predict the gas well performance of the fracture-filled gas hydrate reservoirs considering the effect of natural fractures.

As one of the most applicable and economical method, numerical simulation has been widely used for evaluating the effect of the complex natural fracture network on gas well productivity. However, in terms of capturing the geometries of natural fractures, the traditional numerical simulator requires applying the local grid refinement (LGR) technique or unstructured grids, which can be inconvenient. In order to overcome the deficiencies of the traditional numerical methods, the embedded discrete fracture model (EDFM) method has been introduced to describe complex fractures (Lee et al., 2001). Compared to the methods of LGR and unstructured grids, the EDFM method characterizes the complex natural fractures with unrefined and structured grids. In recent studies, numerous EDFMs have been proposed to investigate the mechanism of mass and heat transfer in various naturally fractured reservoirs.

AlTwaijri et al., 2018 explicitly modeled complex fracture geometries by use of the EDFM method, and their simulation results confirmed that the complexity of the fracture network imposes a considerable influence in estimating gas and water recoveries. Li et al. (2021) combined an extended finite element method (XFEM) and EDFM to model the geothermal systems. The simulation results justified the application of EDFM for modeling fluid mass and heat transfer. Sangnimnuan et al., 2021 upgraded their numerical model by considering the geomechanics effects to model water-oil two-phase flow based on the EDFM. The calculated results indicated that the stress distribution and mass transfer are mainly affected by hydraulic fracture geometries; Shi, (2021) combined the EDFM method with the boundary element method (BEM) to investigate the temperature change considering water-oil two-phase flow. Their calculated results showed that it is highly necessary to consider the temperature field effect on developing fractured geothermal system.

Although the EDFM was used widely to characterize the natural fracture in conventional oil/gas reservoirs, it has not been hitherto applied in gas-hydrate reservoirs to predict the gas well performance (Yan et al., 2016; Yan et al., 2018; Xu et al., 2017; Xu et al., 2018; Xu et al., 2021; Xu et al., 2022). For the scenarios of gas-hydrate reservoirs, multiphase flow, heat transfer, and hydrate dissociation all raise stringent barriers for one to build the EDFM model to evaluate the gas well productivity of the gas-hydrate reservoirs. Therefore, in this work, taking the effect of the natural fractures into consideration, we develop an EDFM model to characterize the gas well performance of the fracture-filled hydrate reservoirs. This proposed model is validated against a commercial simulator. Furthermore, with the aid of the proposed model, we analyze the influences of the complex fracture network, fracture conductivity, number of connected natural fractures, and hydrate dissociation rate on the gas well productivity and the distributions of pressure, temperature, and hydrate saturation.

2 Methodology

2.1 Model assumptions

This study considers a fractured vertical well in a two-dimension (2D) naturally fractured gas hydrate reservoir. In order to consider heat transfer and to characterize the complex fractures, we make the following assumptions to simulate gas-water flow in hydrate reservoirs.

(1) hydrate is treated as an immobile solid phase;

(2) the compressibility of rock and fluids are considered to characterize porosity change;

(3) the reservoir is homogeneous and the flow of water and gas conforms to Darcy’s law;

(4) the capillary pressure is neglected in the model;

(5) the effects of hydrate saturation changes, hydrate occurrence state, and pore structure on the relative permeability of gas and water are neglected;

(6) only thermal conduction and thermal convection in the reservoir are considered, and the influence of thermal radiation is neglected.

2.2 Mass conservation equations

For the matrix system, the mass conservation equation of the gas and water can be expressed as (Hong and Pooladi-Darvish, 2003):

$\nabla {\left({\rho }_l{\nu }_l\right)}^m+{\rho }_l{q}_{l,mf}+m_l=\frac{\partial }{\partial t}{\left({\rho }_l{V}_v{S}_l\right)}^m,l=g,w(1)$

For the intersecting fracture system, the mass conservation equation of the gas and water can be expressed as:

$\nabla {\left({\rho }_l{v}_l\right)}^f+{\rho }_l{q}_{l,w}+{\rho }_l{q}_{l,mf}+{\rho }_l{q}_{l,ff}+m_l=\frac{\partial }{\partial t}{\left({\rho }_l{V}_v{S}_l\right)}^f,l=g,w(2)$

where ν_l is volumetric flow rates, which can be expressed as:

$v_l=A_c\frac{k{k}_{rl}}{{\mu }_l}\nabla P_l(3)$

For hydrate, the mass conservation equations of the matrix and fracture systems are (Hong and Pooladi-Darvish, 2003):

$-{\left(m_h\right)}^{m,f}=\frac{\partial }{\partial t}{\left({\rho }_h{V}_v{S}_h\right)}^{m,f}(4)$

In Eqs 1–4, subscript l represents gas or water; subscript h represents hydrate; the superscripts m and f represent the matrix and the fracture, respectively; k represents absolute permeability; k_rl represents fluid relative permeability; P represents pressure; μ_l represents viscosity; ρ_l represents density; m_l represents the mass accumulation of gas and water due to hydrate dissociation; m_h represents the mass of hydrate dissociation; V_v represent void pore volume; S_l represents gas and water saturation; S_h represents hydrate saturation; q_l,mf represents volumetric flow rates between matrix and fracture; q_l,ff represents volumetric flow rates between fracture and fracture; q_l,w represents volumetric flow rates between fracture and well.

For the EDFM method, unrefined and structured grids are generally used to characterize complex natural fractures, and the finite difference method can be used as the numerical method to discretize mass conservation equations (Yan et al., 2014). For the right-hand-side of Eqs 1, 2, we can have the following discrete form:

$\frac{\partial }{\partial t}\left({\rho }_l{V}_v{S}_l\right)=\frac{V_v{{\rho }_l}^{n+1}\left(c_s+c_l\right){S_l}^{n+1}\left(P^{n+1}-P^n\right)+V_v{{\rho }_l}^{n+1}\left({S_l}^{n+1}-{S_l}^n\right)}{dt}(5)$

where c_s represents rock compossibility coefficient; c_l represents fluid compossibility coefficient.

The discrete form of fluid flow between matrix grids in Eq. 1 can be written as:

$\nabla {\left({\rho }_l{\nu }_l\right)}^m={{\rho }_l}^{n+1}{T}_{m,ij}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}{\left(P_j^m-P_i^m\right)}^{n+1}(6)$

Likewise, the discrete form of fluid flow between fracture grids in Eq. 2 can be written as:

$\nabla {\left({\rho }_l{\nu }_l\right)}^f={{\rho }_l}^{n+1}{T}_{f,ij}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}{\left(P_j^f-P_i^f\right)}^{n+1}(7)$

The discrete form of mass transfer between matrix grid and fracture grid in Eqs 1, 2 can be written as:

${\rho }_l{q}_{l,mf}={{\rho }_l}^{n+1}{T}_{mf}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}{\left({P_i}^m-{P_j}^f\right)}^{n+1}(8)$

The discrete form of mass transfer between intersecting fracture grids in Eq. 2 can be written as:

${\rho }_l{q}_{l,ff}={{\rho }_l}^{n+1}{T}_{ff,ij}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}{\left({P_j}^f-{P_i}^f\right)}^{n+1}(9)$

The discrete form of mass transfer between the wellbore and fracture grid in Eq. 2 can be written as:

${\rho }_l{q}_{l,w}={{\rho }_l}^{n+1}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}\frac{2\pi k_i\Delta h}{\ln \left(r_e/r_w\right)}{\left({P_i}^f-P_{wf}\right)}^{n+1}(10)$

In Eq. 10, r_w is the wellbore radius, and r_e is the equivalent wellbore radius which can be calculated with the Peaceman formulation (Peaceman, 1978), and the details for calculating m_l in Eqs 1, 2 can be found in Supplementary Appendix SA.

In Eqs 6–9, the terms T_m, T_f, T_mf, and T_ff represent the transmissibility factors between different connections of matrix and fracture grids. Figure 1 presents the schematics of different connections of matrix and fracture grids. Three different types of non-neighbor connections (NNCs) in addition to the connections between neighboring matrix grids can be found in Figure 1, including 1) type #1 NNC, the connections between the embedded fracture segments and matrix grids; 2) type #2 NNC: the connections between neighboring fracture segments of a single fracture; and 3) type #3 NNC: the connections between intersected fractures segments. For the fluid flow between matrix grids, taking x-direction as an example, the transmissibility factor between neighboring matrix grid i and j can be written as:

$T_{m,ij}=\frac{2T_{m,i}{T}_{m,j}}{T_{m,i}+T_{m,j}}(11)$

where

$T_m=\frac{k_m\Delta y\Delta z}{\Delta x}(12)$

In Eq. 12, the k_m represents matrix grid permeability; ∆x represents x-direction matrix grid size; ∆y represents y-direction matrix grid size; ∆z represents z-direction matrix grid size.

FIGURE 1

FIGURE 1. Schematic of the EDFM model for fractured hydrate reservoirs.

For type #1 NNC, the matrix-fracture transmissibility factor can be expressed as (Lee et al., 2001):

$T_{mf}=\frac{2A_f{k}_m}{d_{fm}}(13)$

In Eq. 13, the A_f is the area of the cross section between matrix grid and fracture segment, and d_fm is the average normal distance between the embedded fracture segment and matrix grid. For type #2 NNC, the transmissibility factor of fluid flow between neighboring fracture segments of a single fracture is calculated as (Lee et al., 2001):

$T_{f,ij}=\frac{2T_{f,i}{T}_{f,j}}{T_{f,i}+T_{f,j}},T_{f,i}=\frac{k_{fi}{w}_i\Delta h_i}{\Delta l_{fi}}(14)$

In Eq. 14, the k_fi is permeability of ith fracture segment; ∆l_fi, ∆h_i, and w_i represent the length, height, and width of ith fracture segment, respectively. For type #3 NNC, the transmissibility factor of fluid flow between segment i of fracture 1 and segment j of fracture 2 can be expressed as (Lee et al., 2001):

$T_{ff,ij}=\frac{T_{i,1}{T}_{j,2}}{T_{i,1}+T_{j,2}},T_{i,1}=\frac{k_{fi,1}{w}_i\Delta h}{d_{seg1}},T_{j,2}=\frac{k_{fj,2}{w}_i\Delta h}{d_{seg2}}(15)$

In Eq. 15, d_seg is the distance from the fracture segment to the intersection.

The residual mass conservation equations in terms of the finite difference method can be expressed as (Song and Liang, 2009):

$\left\{\begin{array}{l}\mathrm{M}\mathrm{a}\mathrm{t}\mathrm{r}\mathrm{i}\mathrm{x}:{R_{lm,i}}^{n+1}=\frac{{{\left[V_v{{\rho }_l}^{n+1}\left(c_s+c_l\right){S_l}^{n+1}\left(P^{n+1}-P^n\right)+V_v{{\rho }_l}^{n+1}\left({S_l}^{n+1}-{S_l}^n\right)\right]}_i}^m}{dt}\\ -{\displaystyle \sum }_{j=1}^{nc}{{\rho }_l}^{n+1}{T}_{m,ij}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}{\left({P_j}^m-{P_i}^m\right)}^{n+1}-{\displaystyle \sum }_{j=1}^{nf}{{\rho }_l}^{n+1}{T}_{mf}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}{\left({P_j}^f-{P_i}^m\right)}^{n+1}-m_{l,i}^{n+1}\\ \mathrm{F}\mathrm{r}\mathrm{a}\mathrm{c}\mathrm{t}\mathrm{u}\mathrm{r}\mathrm{e}:{R_{lf,i}}^{n+1}=\frac{{{\left[V_v{{\rho }_l}^{n+1}\left(c_s+c_l\right){S_l}^{n+1}\left(P^{n+1}-P^n\right)+V_v{{\rho }_l}^{n+1}\left({S_l}^{n+1}-{S_l}^n\right)\right]}_i}^f}{dt}\\ -{\displaystyle \sum }_{j=1}^{nc}{{\rho }_l}^{n+1}{T}_{f,ij}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}{\left({P_j}^f-{P_i}^f\right)}^{n+1}-{{\rho }_l}^{n+1}{T}_{mf}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}{\left({P_j}^m-{P_i}^f\right)}^{n+1}\\ -{\displaystyle \sum }_{j=1}^{n_{in}}{{\rho }_l}^{n+1}{T}_{ff,ij}\frac{{k_{rl}}^{n+1}}{{{\mu }_l}^{n+1}}{\left({P_j}^f-{P_i}^f\right)}^{n+1}-m_{l,i}^{n+1}\end{array}\right.(16)$

For the hydrate, we can have:

$\left\{\begin{array}{l}\mathrm{M}\mathrm{a}\mathrm{t}\mathrm{r}\mathrm{i}\mathrm{x}:{R_{hm,i}}^{n+1}=\frac{{{\left[V_v{{\rho }_h}^{n+1}{c}_s{S_h}^{n+1}\left(P^{n+1}-P^n\right)+V_v{{\rho }_h}^{n+1}\left({S_h}^{n+1}-{S_h}^n\right)\right]}_i}^m}{dt}+m_{h,i}^{n+1}\\ \mathrm{F}\mathrm{r}\mathrm{a}\mathrm{c}\mathrm{t}\mathrm{u}\mathrm{r}\mathrm{e}:{R_{hf,i}}^{n+1}=\frac{{{\left[V_v{{\rho }_h}^{n+1}{c}_s{S_h}^{n+1}\left(P^{n+1}-P^n\right)+V_v{{\rho }_h}^{n+1}\left({S_h}^{n+1}-{S_h}^n\right)\right]}_i}^f}{dt}+m_{h,i}^{n+1}\end{array}\right.(17)$

2.3 Energy conservation equations

For the matrix system, the energy conservation equation can be expressed as (Hong and Pooladi-Darvish, 2003):

$\frac{{\partial U}^m}{\partial t}={\left({\rho }_l{\nu }_l{H}_l\right)}^m+\nabla {\left({\lambda }_{eff}\nabla T\right)}^m+U_{mf}+U_{diss}(18)$

$U^m={\left[V_v\left(U_l{S}_l{\rho }_l+U_h{S}_h{\rho }_h\right)+U_r{V}_r\right]}^m(19)$

In addition, the energy conservation equation for the fracture system is:

$\frac{{\partial U}^f}{\partial t}={\left({\rho }_l{\nu }_l{H}_l\right)}^f+\nabla {\left({\lambda }_{eff}\nabla T\right)}^f+U_{mf}+U_{ff}+U_{fw}+U_{diss}(20)$

$U^f={\left[V_v\left(U_l{S}_l{\rho }_l+U_h{S}_h{\rho }_h\right)+U_r{V}_r\right]}^f(21)$

In Eqs 18–21, U^m and U^f represent internal energy of matrix and fracture system; U_mf represents energy exchange between matrix and fracture system; U_ff represents energy exchange between different fractures; U_fw represents energy exchange between fracture and well; U_diss represents heat absorbed by hydrate decomposition; U_l represents internal energy of gas and water; U_h represents internal energy of hydrate; U_r represents internal energy of rock; H_l represents enthalpy of gas and water; T represents temperature; λ_eff represents effective thermal conductivity; V_r represents rock volume.

The effective thermal conductivity λ_eff and enthalpy of the water and gas H_l are given as follows (Hong and Pooladi-Darvish, 2003; Esmaeilzadeh et al., 2008):

${\lambda }_{eff}={\varphi }_v\left(S_w{\lambda }_w+S_g{\lambda }_g+S_h{\lambda }_h\right)+\left(1-{\varphi }_v\right){\lambda }_r(22)$

$H_l={\int }_{T_{ref}}^{T_l}{C}_{pl}dT(23)$

In Eqs 22, 23, ϕ_v represents void porosity; S_w, S_g, and S_h represent water, gas, and hydrate saturation respectively; λ_w, λ_g, and λ_h represent water, gas, and hydrate thermal conductivity respectively; T_ref represents reference temperature; C_pl represents fluid heat capacity.

The internal energy U_l is accumulated energy which neglects the mechanical work of a fluid phase. For gas and water, the relationship between U_l and H_l is (Esmaeilzadeh et al., 2008),

$U_l=H_l-\frac{P_l}{{\rho }_l}(24)$

In Eq. 24, P_l is phase pressure and ρ_l is density. Energy exchange U_mf between jth matrix and the ith fracture can be expressed as (Shi, 2021):

$U_{mf}={\rho }_l{q}_{mf}{H}_l+CI{\lambda }_{mf}\left({T_i}^f-{T_j}^m\right)(25)$

In Eq. 25, CI is the connectivity index between the matrix and fractures (CI = A_f/d_fm); λ_mf is the average thermal conductivity between the matrix and the fractures.

Energy exchange U_ff between ith segments of fracture 1 and the jth segments of fracture 2 can be expressed as (Shi, 2021):

$U_{ff}={\rho }_l{q}_{ff}{H}_l+T{I}_{ff}{\lambda }_{ff}\left({T_i}^{f1}-{T_j}^{f2}\right)(26)$

In Eq. 26, λ_ff is the average thermal conductivity between different fractures; TI_ff is the transmissivity factor between different fractures which can be expressed as (Shi, 2021):

$T{I}_{ff}=\frac{2T{I}_{f1}T{I}_{f2}}{T{I}_{f1}+T{I}_{f2}},T{I}_{f1}=\frac{w_{fi}\Delta h}{\Delta l_{fi}}(27)$

Energy exchange U_fw between well and the ith fracture can be expressed as:

$U_{fw}={\rho }_l{q}_{fw}{H}_l(28)$

In general, the discrete form for the energy conservation equation of matrix and fracture grids can be expressed as (Song and Liang, 2009):

$\left\{\begin{array}{l}\text{Matrix}:R_{tm,i}=\frac{{\left[U^{n+1}-U^n\right]}_i^m}{\Delta t}-{\left[{\displaystyle \sum }_{j=1}^n{\rho }_l^{n+1}{\nu }_l^{n+1}{c}_{pl}\left(T^{n+1}-T_r\right)\right]}^m-{\displaystyle \sum }_{j=1}^n{\lambda }_{m,j}{\left(T_j^{n+1}-T_i^{n+1}\right)}^m\\ -{\displaystyle \sum }_{j=1}^{nf}\left[{\rho }_l^{n+1}{q}_{mf}^{n+1}{c}_{pl}{\left(T_j^f-T_r\right)}^{n+1}+CI{\overline{\lambda }}_{mf}{\left(T_j^f-T_i^m\right)}^{n+1}\right]-m_{h,i}^{n+1}{H}_h\\ \text{Fracture}:R_{tf,i}=\frac{{\left[U^{n+1}-U^n\right]}_i^f}{\Delta t}-{\left[{\displaystyle \sum }_{j=1}^n{\rho }_l^{n+1}{\nu }_l^{n+1}{c}_{pl}\left(T^{n+1}-T_r\right)\right]}^f-{\displaystyle \sum }_{j=1}^n{\lambda }_{m,j}{\left(T_j^{n+1}-T_i^{n+1}\right)}^f\\ -{\rho }_l^{n+1}{q}_{mf}^{n+1}{c}_{pl}{\left(T_j^m-T_r\right)}^{n+1}-CI{\overline{\lambda }}_{mf}{\left(T_j^m-T_i^f\right)}^{n+1}\\ -{\displaystyle \sum }_{j=1}^{n_{in}}\left[{\rho }_l^{n+1}{q}_{ff}^{n+1}{c}_{pl}{\left(T_j^f-T_r\right)}^{n+1}+T{I}_{ff}{\overline{\lambda }}_{ff}{\left(T_j^f-T_i^f\right)}^{n+1}\right]\\ -{\rho }_l^{n+1}{q}_{fw}^{n+1}{c}_{pl}{\left(T_i^f-T_r\right)}^{n+1}-m_{h,i}^{n+1}{H}_h\end{array}\right.(29)$

For Eqs 16, 17, 29, values of k_rl, μ_l, ρ_l, and H_l are taken from the phase upstream region and can be solved with a fully implicit coupling (see Supplementary Appendix SB).

3 Model validation

In order to examine the accuracy of the proposed model, we validate the proposed model against a commercial simulator CMG STARS. As the pressure is decreased in the reservoirs, the gas and water can be gradually released from the hydrate. To illustrate the fluid flow of different phases, we utilize the Brooks-Corey equation (Brooks, 1964) to characterize the relative permeability of water and gas which is expressed as:

$k_{rw}=k_{rw,\max }{\left(\frac{S_w-S_{wc}}{1-S_{wc}-S_{gc}}\right)}^{n_w}(30)$

$k_{rg}=k_{rg,\max }{\left(\frac{S_g-S_{gc}}{1-S_{wc}-S_{gc}}\right)}^{n_g}(31)$

In this work, k_rw,max=0.9, k_rg,max=0.7, S_wc=0.2, S_gc=0.05, n_w=n_g=2. The other benchmark values of the parameters used in the simulation are listed in Table 1.

TABLE 1

TABLE 1. Benchmark values of the parameters used for validation.

3.1 Hydrate reservoir with single fracture

In this case, we established a 155×155×30 m³ hydrate reservoir with a single fracture to conduct the validation. Figure 2 gives a top view of the global grid system used to simulate the production of the established model. As shown in Figure 2, the local grid refinement (LGR) was used to simulate the hydraulic fracture. The reservoir domain is discretized into 31 × 31 = 961 matrix grids in the EDFM model. The bottom-hole pressure (BHP) of the fractured vertical gas well which locates at the center of the hydrate reservoir is set to be a constant value of 5 MPa. The fracture length is 35 m, and the fracture conductivity, which is 1,000 mD∙m, used in CMG STARS is consistent to that used in the EDFM model.

FIGURE 2

FIGURE 2. Top view of the grid system used in CMG STARS.

Figure 3 compares the calculated results from the proposed EDFM model to those from CMG STARS. Figures 3A,B illustrate the gas production rates and water production rates, respectively. Figure 3C shows the temperature distribution along the x-direction at y=77.5 m (see the dotted white line in Figure 3) on the 1800th production day. Figure 3D presents the average temperature of the global reservoir model. Figures 3E,F compare the pressure and temperature distribution between LGR and EDFM at 60-day production respectively. As one can see from Figure 3, the results of the proposed model agree well with the results of CMG STARS.

FIGURE 3

FIGURE 3. Validation of the single fracture model: (A) daily gas rate; (B) daily water rate; (C) matrix average temperature; (D) temperature distribution along the x-direction at y=77.5 m; (E) pressure distribution of LGR and EDFM solutions; (F) temperature distribution of LGR and EDFM solutions.

3.2 Hydrate reservoir with complex fracture network

In addition to the single fracture case, the validation was conducted on a complex fracture network model (see Figure 4). As shown in Figure 4, there are 10 natural fractures as well as a hydraulic fracture in the hydrate reservoir model. The conductivity of the natural fractures used both in CMG and EDFM is 1,000 mD∙m.

FIGURE 4

FIGURE 4. Schematic diagram of complex fracture network model used for validation.

Figure 5 compares the results obtained from EDFM against the results from CMG STARS, including gas production rate (Figure 5A), water production rate (Figure 5B), temperature distribution at y=77.5 m on the 1800th production day (Figure 5C), and the average temperature of the global reservoir model (Figure 5D). In addition, Figure 5e∼5F shows a comparison between LGR and EDFM on pressure and temperature distribution at 60-day production respectively. As shown in Figure 5, the results of the proposed model show slight differences from the results of CMG STARS.

FIGURE 5

FIGURE 5. Validation of complex fracture model: (A) daily gas rate; (B) daily water rate; (C) matrix average temperature; (D) temperature distribution along the x-direction at y=77.5 m; (E) pressure distribution of LGR and EDFM solutions; (F) temperature distribution of LGR and EDFM solutions.

Based on the results shown in Figure 3 and Figure 5, one can find that the proposed EDFM model is reliable in simulating the production of gas-hydrate reservoirs considering complex fracture networks and heat transfer. Thus, in this work, the authors will carry out a comprehensive study of the well production of fracture-filled hydrate reservoirs on the basis of this proposed EDFM model.

4 Results and discussion

With the aid of the proposed EDFM model, the authors studied the effects of natural fractures and heat transfer on the production of gas-hydrate reservoirs. The dimension of the reservoir is 310 × 310 × 30 m³ which is discretized into 31 × 31 × 1 grids. The values of the other reservoir and fluid parameters are the same as those in Table 1.

4.1 Effect of complex fracture network

In order to explore the effect of fracture networks on the performance of gas-hydrate production well, three scenarios were investigated: Scenario #1 which considers 30 natural fractures and 1 hydraulic fracture; Scenario #2 which only considers 1 hydraulic fracture; and Scenario #3 which contains no fracture. Figure 6 shows the top view of three scenarios in this section. It is worth noting that, there are 2 natural fractures (green line in Figure 6A) connected with hydraulic fracture (red line in Figure 6A) for scenario#1. The hydraulic fracture length is 50 m, and the conductivity is 1,000 mD∙m for natural fractures and 2000 mD∙m for hydraulic fracture. The BHP is set to be a constant of 5 MPa.

FIGURE 6

FIGURE 6. Schematics of the different scenarios to investigate the effect of fracture networks: (A) scenario with natural fractures and hydraulic fracture; (B) scenario with hydraulic fracture; (C) scenario without natural fractures and hydraulic.

Figure 7A compares the gas production rates of the three scenarios. In Figure 7A, one can find that the daily gas production rate (q_g) and cumulative gas production (Q_g) with natural and hydraulic fractures (scenario #1) are both higher than those only with one hydraulic fracture (scenario #2) and those without fractures (scenario #3). Comparing the cumulative production of scenario #2 to that of scenario #3, one can find that the hydraulic fracturing treatment can significantly improve gas well productivity. In addition, the existence of the natural fracture increases the cumulative production from 4.39 × 10⁶ m³ (scenario #2) to 6.39 × 10⁶ m³ (scenario #1). This indicates that both hydraulic fracture and natural fractures play important roles in influencing the performance of gas-hydrate wells.

FIGURE 7

FIGURE 7. Simulation outputs of different scenarios: (A) Comparison of the production rates (B) pressure, temperature, and hydrate saturation field maps at the 1800th production day.

Figure 7B shows simulation results of different scenarios at the end of the simulation. As shown in pressure map, a higher pressure drop can be observed along the hydraulic fracture and along the connected natural fractures (see scenario#1). This implies that natural fractures can noticeably increase the expanding speed of the drainage area of the production well due to its connection to hydraulic fractures, thus, leading to the highest well productivity among the three studied scenarios. Besides, a heterogeneous distribution can be observed in temperature field maps affected by natural and hydraulic fractures. The low-temperature zone gradually spreads over time along the fractures for the scenario with natural and hydraulic fractures (scenarios #1 and #2) while it gradually spreads along the radial direction for the scenario without fractures (scenario #3). In addition, the existence of natural and hydraulic fractures increases the hydrate dissociation rate. The maximum hydrate dissociation rate occurs near the hydraulic fracture, connected natural fractures, and the position of the wellbore.

4.2 Effect of number of natural fractures

In this section, different numbers of natural fractures (N_nf =2, 10, 30, and 50) and different numbers of connected natural fractures (N_{f_con}=0, 1, 2, 4) were examined to explore their influence on pressure and temperature distribution. As shown in Figure 8A, 0, 8, 28, and 48 natural fractures (blue lines in Figure 8A) are randomly generated in the hydrate reservoir, respectively. In addition to the randomly generated fractures, in all four cases, two natural fractures are deliberately generated to intersect with the hydraulic fracture. In addition, we further investigated the well productivity with different numbers of connected natural fractures (N_{f_con}) to explore the effect of connected natural fractures on gas well productivity. The connected natural fractures indicate a natural fracture that is connected to the hydraulic fracture. As shown in Figure 8B, in the four studied reservoir models, the distribution and number of the randomly generated unconnected natural fracture (blue lines) remain unchanged through the four studied cases, whereas the number of the connected natural fracture (green lines) is increased from 0 to 4.

FIGURE 8

FIGURE 8. Schematics of the different scenarios to investigate the effect of the number of natural fractures: (A) different values of N_nf with constant value of N_{f_con} (B) constant value of N_nf with different values of N_{f_con}.

Table 2 summarizes the cumulative gas production of the different scenarios. As we can see in Table 2, the cumulative gas production rates are slightly increased as the N_nf is varied from 10 to 50. For example, the cumulative gas production is increased by 6.48% as N_nf is increased from 2 to 50. However, the cumulative gas production rate is significantly increased as the value of N_{f_con} is increased. In addition, the well productivity is more likely to be improved with less connected natural fractures. For example, the cumulative gas production is increased by 95.76 × 10⁴ m³ as N_{f_con} is increased from 0 to 1, whereas the increase is only 16.84 × 10⁴ m³ as N_{f_con} is increased from 2 to 4.

TABLE 2

TABLE 2. Comparison of cumulative gas production of different scenarios with different values of N_nf and N_{f_con}.

A comparison of pressure, temperature, and hydrate saturation after 1800 days of production with different values of N_nf and N_{f_con} is shown in Figure 9. It can be observed in Figure 9A that, although the natural fractures are increased from 2 to 50, the pressure, temperature, and hydrate saturation fields illustrate a negligible difference. This denotes that the natural fractures, which are not connected to the hydraulic fracture, exert a small influence on the productivity of the gas-hydrate well. In Figure 9B, one can find that the connected natural fractures can induce remarkable expansions of the drainage area. The expansion can be more significant as N_{f_con} is increased from 0 to 1 than that as N_{f_con} is increased from 2 to 4. The results in Table 2 and Figure 9 imply that the connected natural fractures improve the well productivity by increasing the drainage areas of the production well. The calculated results represent that the connected natural fractures can exert a more important influence on the well productivity than the unconnected natural fractures.

FIGURE 9

FIGURE 9. Pressure, temperature, and hydrate saturation field maps of different scenarios at the 1800th production day: (A) different values of N_nf with constant value of N_{f_con} (B) constant value of N_nf with different values of N_{f_con}.

4.3 Effect of fracture conductivity

In this section, the authors discussed the impacts of dimensionless hydraulic fracture conductivity (C_fD) on the productivity of a fractured vertical gas well. The dimensionless hydraulic fracture conductivity can be described as (Cinco et al., 1978):

$C_{fD}=\frac{k_f{w}_f}{k_m{x}_f}(32)$

In Eq. 32, k_f and k_m is hydraulic fracture permeability and matrix permeability, respectively; x_f and w_f is hydraulic fracture half length and width, respectively.

Figure 10A displays the gas production rates of a fractured vertical well that is calculated with the proposed numerical model with different dimensionless hydraulic fracture conductivity. As we can see from Figure 10A, A larger dimensionless fracture conductivity leads to higher well productivity, which is attributed to a lower flow resistance in the fracture. One also can find that the q_g and Q_g express a slight difference if the dimensionless fracture conductivity exceeds 100. This is because the fracture can be regarded as infinite-conductive if the dimensionless fracture conductivity exceeds 100, and the flow resistance in the fracture can be ignored.

FIGURE 10

FIGURE 10. Simulation outputs of different scenarios with different values of hydraulic fracture conductivity: (A) Comparison of the production rates (B) pressure, temperature, and hydrate saturation field maps at the 1800th production day.

Figure 10B compares the field maps at the end of the 1800th day with different values of C_fD. As we can see from pressure maps in Figure 10B, the maximum pressure drop occurs at the center of the gas hydrate reservoir. As the fracture conductivity is increased, the expanding speed of the drainage area becomes higher (e.g., the yellow area of C_fD = 100 is larger than that of C_fD = 1 at the end of the 1800th day in pressure field maps). In addition, for a higher conductivity fracture, the temperature drops more rapidly near the fractures (see temperature field maps). This is attributed to the fact that the temperature is mainly influenced by thermal convection and solid-hydrate dissociation. Since the gas and water dissociated from solid hydrate is more readily to be produced through the fracture with a higher conductivity, more solid hydrate will be dissociated near the fracture. The solid-hydrate dissociation is endothermic which can result in a reduction of temperature. In hydrate saturation field maps, one can find that the hydrate saturation near the fracture is lower with higher fracture conductivity, which indicates that more solid-hydrate is dissociated near the fracture.

Moreover, the authors varied the conductivity of natural fractures together with the conductivity of the hydraulic fracture to compare their effects on well productivity. Table 3 summarizes the cumulative gas production for different hydraulic fracture conductivity and different natural fracture conductivity at the end of the simulation. In Table 3, both the natural fracture conductivity and the hydraulic fracture conductivity are varied from 25 mD∙m to 250 mD∙m. The calculated results show that the cumulative gas production is increased both as the hydraulic fracture conductivity and the natural fracture conductivity are increased. Besides, hydraulic fracture conductivity exerts a more significant influence on gas well productivity than natural fracture conductivity. For example, the maximum increase of the cumulative gas production that is induced by varying the natural fracture conductivity occurs with C_hf = 250 mD∙m, and the relative increase is 15.14% which is from 506.07 × 10⁴ m³ with C_nf=25 mD∙m to 582.69 × 10⁴ m³ with C_nf=250 mD∙m. However, with C_nf = 250 mD∙m, the relative increase of the cumulative production by varying the hydraulic fracture conductivity is 58.15% which is from 368.43 × 10⁴ m³ with C_hf=25 mD∙m to 582.69 × 10⁴ m³ with C_hf=250 mD∙m.

TABLE 3

TABLE 3. Comparison of cumulative gas production of a vertical well with different values of hydraulic fracture and natural fracture conductivity.

4.4 Effect of hydrate dissociation rate

Hydrate dissociation rate also imposes a significant influence on the mass and heat transfer in developing gas hydrate. The kinetic model of Kim et al. showed that the hydrate decomposition rate is related to decomposition rate constant (K_D), hydrate surface area (A_dec), and pressure difference which can be expressed as,

$m_h=K_D{A}_{dec}\left(P_e-P_g\right)(33)$

The decomposition rate constant K_D is defined as,

$K_D=K_d\exp \left(-\frac{E}{RT}\right)(34)$

where K_d represents the intrinsic hydrate dissociation rate constant, mole/(day·KPa·m²). As we can see from Eqs 33, 34, the hydrate dissociation rate is a time-varying parameter in developing gas-hydrate. To quantitatively examine the effect of hydrate dissociation rate on the gas well productivity of fracture-filled hydrate reservoir, the different values of intrinsic hydrate dissociation rate constant (K_d) were considered, including 0, 1.07 × 10⁴, 1.07 × 10⁹, and 1.07 × 10¹³ mol/(day·KPa·m²). It is noted that an intrinsic hydrate dissociation rate of 0 indicates that the hydrate dissociation is ignored, and the initial hydrate decomposition is more rapid as the increase of K_d with the same values of hydrate surface area and pressure difference. The gas production rates with different values of K_d are illustrated in Figure 11A. As we can see from Figure 11A, the plots show that the gas well productivity is increased as the hydrate dissociation rate is increased, indicating that hydrate dissociation is favorable for the production of gas-hydrate reservoirs.

FIGURE 11

FIGURE 11. Simulation outputs of different scenarios with different values of K_d: (A) comparison of the production rates (B) pressure, temperature, and hydrate saturation field maps at the 1800th production day.

Figure 11B shows the field maps with different values of K_d. An interesting observation is that the pressure drop is most noticeable with the smallest value of K_d = 0, which implies that the drainage area expands most rapidly without hydrate dissociation rate, whereas the temperature is reduced most significantly with the largest value of K_d = 1.07 × 10¹³ mol/(day·KPa·m²) (see temperature field maps in Figure 11B). It is worth noting that, in Sections 4.1–4.4, the highest well productivity comes together with the largest rate of drainage area expansion and the largest rate of temperature reduction. The calculated results in Figure 11B denote that the temperature change can be more highly related to the hydrate dissociation than the drainage area expansion. In addition, the hydrate dissociation rate can lead to a reduction of solid hydrate, thus, one can see that the hydrate saturation is lower near the fracture with a higher hydrate dissociation rate. The minimum matrix temperature along the natural fracture is around 6.3°C with the value of K_d = 1.07 × 10¹³ mol/(day·KPa·m²), whereas the temperature is 7.1°C with the value of K_d = 1.07 × 10⁴ mol/(day·KPa·m²). Furthermore, the thermotactic habit of hydrate nucleation may play a role in the secondary hydrate formation process which would be important in temperature distribution (Yang et al., 2023).

5 Conclusion

In this work, we developed an embedded discrete fracture model to study the gas well productivity in fracture-filled gas hydrate reservoirs. By using a finite-difference scheme and a fully implicit method to solve the non-linear equations, we validated the proposed model against a commercial simulator. In addition, we investigated the influences of natural fractures, fracture conductivity, and hydrate dissociation rate on mass and heat transfer of fracture-filled gas hydrate reservoir. The calculated results lead us to draw the following conclusions.

(1) This proposed numerical model is reliable for evaluating the gas well productivity of fracture-filled hydrate reservoirs considering the effect of natural fractures and heat transfer, and the research results are of great significance to guide the efficient development of fracture-filled hydrate reservoirs.

(2) Both hydraulic fractures and natural fractures can increase gas well productivity. Moreover, the connected natural fractures can exert a more important influence on the well productivity than the unconnected natural fractures. The cumulative gas production is increased by 6.48% as N_nf is increased from 2 to 50, whereas the increase is 43.38% as N_{f_con} is increased from 0 to 4.

(3) The gas well production rate is increased by increasing fracture conductivity, and the conductivity of hydraulic fracture exerts a more significant influence on the gas well productivity than that of natural fracture. The relative cumulative gas production increment is 15.14% and 58.15% as C_hf and C_nf are increased from 25 mD to 250 mD, respectively.

(4) The well productivity is increased as the hydrate dissociation rate is increased. Besides, a higher hydrate dissociation rate can lead to a lower temperature along fractures due to a larger reduction of solid hydrate. The minimum temperature along natural fractures is around 6.3°C with the value of K_d = 1.07 × 10¹³ mol/(day·KPa·m²).

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

Author contributions

YC: methodology, conceptualization, writing original draft, review, and editing. BT: conceptualization, reviewing, and editing. WL: conceptualization, reviewing, and editing. CL: investigation, reviewing. YZ: reviewing.

Acknowledgments

The authors are grateful to the National Natural Science Foundation of China (51991365), China Geological Survey Project (No. DD20211350), and Guangdong Major Project of Basic and Applied Basic Research (No. 2020B0301030003).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/feart.2023.1132970/full#supplementary-material

References

AlTwaijri, M., Xia, Z., Yu, W., Qu, L., Hu, Y., Xu, Y., et al. (2018). Numerical study of complex fracture geometry effect on two-phase performance of shale-gas wells using the fast EDFM method. J. Pet. Sci. Eng. 164, 603–622. doi:10.1016/j.petrol.2017.12.086

CrossRef Full Text | Google Scholar

Boswell, R., Collett, T. S., Frye, M., Shedd, W., McConnell, D. R., and Shelander, D. (2012). Subsurface gas hydrates in the northern Gulf of Mexico. Mar. Pet. Geol. 34, 4–30. doi:10.1016/j.marpetgeo.2011.10.003

CrossRef Full Text | Google Scholar

Brooks, R. H. (1964). Hydraulic properties of porous media. Fort Collins: Colorado State University.

Google Scholar

Carroll, J. (2020). Natural gas hydrates: A guide for engineers. Burlington: Gulf Professional Publishing.

Google Scholar

Cinco, L., H., Samaniego, V., F., and Dominguez, A., N. (1978). Transient pressure behavior for a well with a finite-conductivity vertical fracture. Soc. Pet. Eng. J. 18, 253–264. doi:10.2118/6014-PA

CrossRef Full Text | Google Scholar

Esmaeilzadeh, F., Zeighami, M. E., and Fathi, J. (2008). 1-D modeling of hydrate decomposition in porous media. Int. J. Chem. Mol. Eng. 2, 45–51.

Google Scholar

Hong, H. N., and Pooladi-Darvish, M. (2003). “A numerical study on gas production from formations containing gas hydrates,” in Canadian international Petroleum conference (Calgary, Alberta: Petroleum Society of Canada). doi:10.2118/2003-060

CrossRef Full Text | Google Scholar

Kamath, V. A., and Holder, G. D. (1987). Dissociation heat transfer characteristics of methane hydrates. AIChE J. 33, 347–350. doi:10.1002/aic.690330228

CrossRef Full Text | Google Scholar

Kim, H. C., Bishnoi, P. R., Heidemann, R. A., and Rizvi, S. S. H. (1987). Kinetics of methane hydrate decomposition. Chem. Eng. Sci. 42, 1645–1653. doi:10.1016/0009-2509(87)80169-0

CrossRef Full Text | Google Scholar

Lee, S. H., Lough, M. F., and Jensen, C. L. (2001). Hierarchical modeling of flow in naturally fractured formations with multiple length scales. Water Resour. Res. 37, 443–455. doi:10.1029/2000WR900340

CrossRef Full Text | Google Scholar

Li, T., Han, D., Yang, F., Li, J., Wang, D., Yu, B., et al. (2021a). Modeling study of the thermal-hydraulic-mechanical coupling process for EGS based on the framework of EDFM and XFEM. Geothermics 89, 101953. doi:10.1016/j.geothermics.2020.101953

CrossRef Full Text | Google Scholar

Li, Y., Liu, L., Jin, Y., and Wu, N. (2021b). Characterization and development of marine natural gas hydrate reservoirs in clayey-silt sediments: A review and discussion. Adv. Geo-Energy Res. 5, 75–86. doi:10.46690/ager.2021.01.08

CrossRef Full Text | Google Scholar

Liu, L., Lu, X., Zhang, X., Liu, C., and Du, B. (2017). Numerical simulations for analyzing deformation characteristics of hydrate-bearing sediments during depressurization. Adv. Geo-Energy Res. 1, 135–147. doi:10.26804/ager.2017.03.01

CrossRef Full Text | Google Scholar

Moridis, G. J., Kowalsky, M. B., and Pruess, K. (2005). “Depressurization-induced gas production from class-1 hydrate deposits,” in All days (Dallas, Texas: SPE). SPE-97266-MS. doi:10.2118/97266-MS

CrossRef Full Text | Google Scholar

Moridis, G. J. (2004). Numerical studies of gas production from Class 2 and Class 3 hydrate accumulations at the mallik site, mackenzie delta, Canada. SPE Reserv. Eval. Eng. 7, 175–183. doi:10.2118/88039-PA

CrossRef Full Text | Google Scholar

Moridis, G. J. (2003). Numerical studies of gas production from methane hydrates. SPE J. 8, 359–370. doi:10.2118/87330-PA

CrossRef Full Text | Google Scholar

Ning, F. L., Liang, J. Q., Wu, N. Y., Zhu, Y. H., Wu, S. G., Liu, C. L., et al. (2020). Reservoir characteristics of natural gas hydrates in China. Nat. Gas. Ind. 40, 1–24.

Google Scholar

Peaceman, D. W. (1978). Interpretation of well-block pressures in numerical reservoir simulation(includes associated paper 6988). Soc. Pet. Eng. J. 18, 183–194. doi:10.2118/6893-PA

CrossRef Full Text | Google Scholar

Pooladi-Darvish, M., and Hong, H. (2004). “Effect of conductive and convective heat flow on gas production from natural hydrates by depressurization,” in Advances in the study of gas hydrates. Editors C. E. Taylor, and J. T. Kwan (Boston: Kluwer Academic Publishers), 43–65. doi:10.1007/0-306-48645-8_4

CrossRef Full Text | Google Scholar

Sangnimnuan, A., Li, J., and Wu, K. (2021). Development of coupled two phase flow and geomechanics model to predict stress evolution in unconventional reservoirs with complex fracture geometry. J. Pet. Sci. Eng. 196, 108072. doi:10.1016/j.petrol.2020.108072

CrossRef Full Text | Google Scholar

Shi, J. (2021). A modified embedded discrete-fracture model to study oil-water two-phase heat and mass transfer in the complex fracture network. Int. J. Heat. Mass Transf. 19.

CrossRef Full Text | Google Scholar

Sloan, E. D., Koh, C. A., and Koh, C. A. (2007). Clathrate hydrates of natural gases. 3rd ed. Boca Raton: CRC Press. doi:10.1201/9781420008494

CrossRef Full Text | Google Scholar

Song, Y., and Liang, H. (2009). 2-D Numerical simulation of natural gas hydrate decomposition through depressurization by fully implicit method. Chin. Ocean. Eng. Soc. 23, 529–542.

Google Scholar

Sun, X., Nanchary, N., and Mohanty, K. K. (2005). 1-D modeling of hydrate depressurization in porous media. Transp. Porous Media 58, 315–338. doi:10.1007/s11242-004-1410-x

CrossRef Full Text | Google Scholar

Tonnet, N., and Herri, J.-M. (2009). Methane hydrates bearing synthetic sediments—experimental and numerical approaches of the dissociation. Chem. Eng. Sci. 12.

Google Scholar

Wang, P., Zhu, Y., Lu, Z., Bai, M., Huang, X., Pang, S., et al. (2019). Research progress of gas hydrates in the qilian mountain permafrost, qinghai, northwest China: Review. Sci. Sin. Phys. Mech. Astron. 49, 034606. doi:10.1360/sspma2018-00133

CrossRef Full Text | Google Scholar

Wang, P., Zhu, Y., Lu, Z., Guo, X., and Huang, X. (2011). Gas hydrate in the Qilian Mountain permafrost and its distribution characteristics. Geol. Bull. China 30, 1839–1850.

Google Scholar

Xu, S., Feng, Q., Wang, S., Javadpour, F., and Li, Y. (2018). Optimization of multistage fractured horizontal well in tight oil based on embedded discrete fracture model. Comput. Chem. Eng. 117, 291–308. doi:10.1016/j.compchemeng.2018.06.015

CrossRef Full Text | Google Scholar

Xu, S., Feng, Q., Wang, S., Zha, M., Cui, R., Gao, F., et al. (2017). “A stochastic permeability model for shale formations based on embedded discrete fracture model,” in SPE/AAPG/SEG Unconventional Resources Technology Conference (OnePetro).

Google Scholar

Xu, S., Guo, J., Feng, Q., Ren, G., Li, Y., and Wang, S. (2022). Optimization of hydraulic fracturing treatment parameters to maximize economic benefit in tight oil. Fuel 329, 125329. doi:10.1016/j.fuel.2022.125329

CrossRef Full Text | Google Scholar

Xu, S., Ren, G., Younis, R. M., and Feng, Q. (2021). Revisiting field estimates for carbon dioxide storage in depleted shale gas reservoirs: The role of geomechanics. Int. J. Greenh. Gas. Control 105, 103222. doi:10.1016/j.ijggc.2020.103222

CrossRef Full Text | Google Scholar

Yan, X., Huang, Z., Jun, Y. A. O., and Huang, T. (2014). The embeded discrete fracture model based on mimetic finite difference method. Sci. Sin. Technol. 44, 1333–1342.

Google Scholar

Yan, X., Huang, Z., Yao, J., Li, Y., and Fan, D. (2016). An efficient embedded discrete fracture model based on mimetic finite difference method. J. Pet. Sci. Eng. 145, 11–21. doi:10.1016/j.petrol.2016.03.013

CrossRef Full Text | Google Scholar

Yan, X., Huang, Z., Yao, J., Li, Y., Fan, D., and Zhang, K. (2018). An efficient hydro-mechanical model for coupled multi-porosity and discrete fracture porous media. Comput. Mech. 62, 943–962. doi:10.1007/s00466-018-1541-5

CrossRef Full Text | Google Scholar

Yang, L., Guan, D., Qu, A., Li, Q., Ge, Y., Liang, H., et al. (2023). Thermotactic habit of gas hydrate growth enables a fast transformation of melting ice. Appl. Energy 331, 120372. doi:10.1016/j.apenergy.2022.120372

CrossRef Full Text | Google Scholar

Zhao, J., Liu, D., Yang, M., and Song, Y. (2014). Analysis of heat transfer effects on gas production from methane hydrate by depressurization. Int. J. Heat. Mass Transf. 13.

CrossRef Full Text | Google Scholar

Nomenclature

A_c Contact area, m²

Adec Hydrate decomposition surface area, m²/m³

Ahs Specific area, m²/m³

Af Contact area of matrix and fracture segments, m²

c_s Rock compossibility coefficient, MP^-1

c_l Gas and water compossibility coefficient, MP^-1

C_p Heat capacities of gas, water and hydrate, J/(kg·°C)

d_fm Average normal distance between matrix and fracture grid, m

d_seg Average normal distance from the fracture to the intersection, m

E Reaction activation energy, J/mole

H Enthalpy of the water and gas, J

H_h Reaction enthalpy, J

k Matrix and fracture permeability, mD

k_m, k_f Matrix and fracture permeability, mD

k_r Relative permeability

K_D Hydrate decomposition rate constant, mole/(m³·KPa·day)

K_d Intrinsic reaction rate constant, mole/(m³·KPa·day)

m_h Mass of hydrate dissociation, kg

m_l Mass accumulation due to hydrate dissociation, kg

M_w, M_g, M_h Molecular weight of water, gas, and hydrate, kg/gmole

N_h Hydrate number

P Pressure, MPa

P_e Equilibrium pressure of hydrate, KPa

P_wf Well bottom-hole pressure, MPa

q_l,ff Volumetric flow rates between intersecting fractures, m³/day

q_l,mf Volumetric flow rates between matrix and fracture, m³/day

q_l,w Volumetric flow rates between fracture and well, m³/day

R Gas constant, J/mole-K

r_e Equivalent wellbore radius, m

r_w Wellbore radius, m

S_w, S_g, S_h Water, gas, and hydrate saturation

T Temperature, °C

T_e Equilibrium temperature of hydrate, °C

T_ref Reference temperature, °C

T_f Transmissibility between fracture segments in same fracture

T_ff Fracture-fracture transmissibility

T_m Matrix-matrix transmissibility

T_mf Matrix-fracture transmissibility

T_ref Reference temperature, °C

U Internal energy, J

U^m Internal energy of matrix system, J

U^f Internal energy of fracture system, J

U_l Internal energy of gas and water, J

U_h Hydrate internal energy, J

U_r Rock internal energy, J

U_diss Heat absorbed by hydrate decomposition

U_ff Energy exchange between different fracture, J

U_fw Energy exchange between fracture and well, J

U_mf Energy exchange between matrix and fracture system, J

v_l Volumetric flow rates, m³/day

V_r Rock volume, m³

V_v Void pore volume, m³

λ Thermal conductivity, J/(m·day·°C)

λ_eff Effective thermal conductivity, J/(m·day·°C)

λ_g, λ_w, λ_h Water, gas, and hydrate thermal conductivity, J/(m·day·°C)

λ_ff Average thermal conductivity between fractures, J/(m·day·°C)

μ Gas and water viscosity, mPa·s

ρ_h Hydrate density, kg/m³

ρ_l Fluid density, kg/m³

ϕ Porosity

ϕ_v Void porosity

∆l, ∆w, ∆h Length, width, and thickness of fracture grids

∆x, ∆y, ∆z Length, width, and thickness of matrix grids

Superscripts

m Matrix

f Fracture

Subscripts

l Gas and water

h Solid hydrate