Fluid-Solid Coupling Numerical Simulation on Natural Gas Production from Hydrate Reservoirs by Depressurization

Natural gas hydrate is considered to be the most prospective energy in the 21st century for its various advantages, such as extensive distributing, vast amount of quantity, high energy density, and environmental-friendly.[1] Since the end of 1980s in last century, numerous countries have been playing highly attention on the exploration of natural gas hydrate and made a long-term plan from the point of energy stratagem.

Nowadays, numerical simulation on the productivity of natural gas hydrates is increasingly becoming a spirited issue for scholars. Various numerous stimulation presented generally are only considering two-phase flow between gas and water and the formation permeability change by hydrate decomposition at the isothermal condition, while the change of mechanical properties and the formation temperature due to hydrates decomposition, the overall realization on the fluid-solid coupling seepage features of hydrates reservoir with the particular phase transition changing are ignored. What's more, the present models place a value on productivity analysis rather than the stability of the decomposing zones nearby wells with weakly cemented, low strength and high porosity and permeability.
The objectives of this study are finally to construct a theoretical analysis system for the whole exploiting process of natural gas hydrates production by depressurization with fluidsolid coupling numerical simulation, which are based on physicochemical theory. Moreover, on the basis of the dynamic change law of stresses status and physical properties for hydrates formation and considering some features of dissociation zones nearby wellbore with weakly cemented, low strength, high porosity and permeability etc, this research does relating research on fluid-solid coupling numerical simulation for natural gas hydrate production, which provides theoretical basis and necessary technological support for the industrial production of natural gas hydrate.
Fluid-Solid Coupling Numerical Simulation on Natural Gas Production from Hydrate Reservoirs by Depressurization 149 morphology change of hydrates and the change of physical and mechanical properties, pore pressure and temperature of gas hydrates formation, which is caused by hydrates dissociation. Additionally, the mingled interaction between seepage field and solid deformation field is considerable. Consequently, the exploitation of natural gas hydrates by depressurization is a complex process influenced by numerous factors.

Research methods of fluid-solid coupling
Oil-bearing formation is a geological system which consists of cemented solid rock skeleton, pore space, and inner fluids (gas and water) in multi-phase porous media. Cemented rock skeleton and pore space make up of rock formation that is a basic unit for appraising the level of mechanical performances and deformation. A distinguished feature of solid-fluid coupling in seepage process, namely solid field and fluid filed mingle and include together while hard to differentiate clearly. Thus, fluid phase and solid phase should be seem as overlapping continuous media and the reaction can happen between different phases [3] . This feature contributes the establish of control equations for fluid-solid coupling model need to aim at special physical phenomena .Meanwhile, the effects of fluid-solid coupling are reflected by control equations, namely some items in control equations of fluid flow can describe solid deformation and vice verse.
According to the multiple definition of field, geological system is classified into interacting multiple solid deformation fields and multiphase fluid seepage fields in math. A solid deformation media field and fluid seepage field are connected by fluid system pressure and multiphase seepage fluid is connected by capillary pressure. The physical quantities of describing solid and fluid parameters are reflected by the representative elementary volume (REV) in continuous medium mathematic model.
Thus, porous media can be regarded as a large number of particles with enough porous and skeleton. Therefore, porous particles can define their materials parameters, such as fluid density, solid density and strength etc. Meanwhile, particles can bear stress and fluid pressure, namely particles can be defined as state variables. When the particle is enough small comparing the whole seepage filed, material parameters and variables are the function of points in space and are changing constan t l y w i t h t h e c o n t i n u o u s c h a n g e o f t h e i r positions. Consequently, porous media can virtually be replaced by a kind of hypothetical continuous media.
The analysis of coupling effect between fluid seepage and rock deformation requires to systematically utilize diverse theories , for instance, rock mechanics, permeation fluid mechanics etc, along with the research of coupling rules for seepage field and deformation and establishing relating control equations. Two meanings are included in this study. For one thing, under the effect of reservoir deformation field, this paper present a research about how to obtain its multiphase seepage rules (oil, gas and water); For another, the research about the stress, strain and strength problems etc. for formation rock under the circumstance of seepage field is discussed.
The solid-fluid coupling numerical simulation models for the exploitation of natural gas hydrate by depressurization mainly include four parts: temperature field, seepage field, deformation field and the relations of seepage and deformation coupling model.

www.intechopen.com
Advances in Natural Gas Technology 150 As previously mentioned the exploitation of natural gas hydrates by depressurization is a complex solid-fluid coupling process including phase transition and is influenced by diverse elements. Some basic assumptions about solid-fluid coupling numerical simulation models in this study are as follows: 1. Formation rock abide generalized Hook law and Drucker-Prager yield criterion, and deformation of rock matrix is generally small displacement deformation; 2. Seepage system is made up of two phases(water and gas) at saturation without concerning the chemical effect between formation rock and fluid; 3. The flow of each phase in fluid satisfied generalized Darcy's law comparing rock skeleton particles; 4. Relative permeability of two phases (gas and water) and capillary pressure are the single function of water saturation; 5. Temperature only impact the decomposition of hydrates but not rock skeleton and fluid properties; 6. Formation rock skeleton density and solid hydrate density are constant during the process of deformation; 7. The relative velocity between solid hydrate and formation is zero and the relative velocity between gas and fluid is zero either.

Kinetic equation of gas hydrates decomposition
Kinetic equation of gas hydrates decomposition is proposed by Kim-Bishnoi [4] , as follows: () mK M A ff g rd g dec e g   (1) or () eg mK M A P P g rd g dec e g The intrinsic dissociation rate constant K rd is a function of the system temperature.
The specific surface area of the hydrate decomposition, A dec , is given by The rate of water generated is calculated according to the rate of the gas generated from the hydrate decomposition (m g ) and Where, N H is the hydrate number, the value of which is 6.0 for methane hydrate.

Energy conservation equations of natural gas hydrate
With regard to heat conduction, convection and supplement of heat from outside and without consideration of kinetic energy and thermal radiation, energy conservation equations would be calculated by the following equation: The left side of Eq. (6) means the increment of internal energy and the first item at the right side of Eq. (6) reflects energy accessing the unit porous media system by heat conduction, the second item at the right side of Eq. (6) illustrates energy taken away from gas and water through hydrate compositing, the third item at the right side of the equation shows energy supplement to hydrate formation from outside.
Without considering throttling effect of influence, energy conservation equation of hydrate can be signified with temperature and illustrated by the following equation after equation groupies derived: Heat of hydrate decomposition causes formation temperature change and limitation for decomposition rate of hydrate. Therefore, the calculation of temperature field is significant to the research of hydrate reservoir simulation.

Rock matrix continuity equation of gas hydrate reservoirs
Rock matrix continuity equation is as follows: If rock particles density is ρ s , then: Assume the density of rock matrix in hydrate reservoir remains unchanged during the deformation, namely ρ s is constant, then: www.intechopen.com Eq. (10) reflects the relationship between the porosity changes and matrix displacement velocity.

Fluid-solid coupling seepage motion equations of gas hydrate reservoirs
Similar with the normal reservoirs, during production from hydrate reservoirs by depressurization, not only fluid particle seepage take place in the porous media, but the reservoirs will undergo deformations and the rock particles will move for the rigid displacement and deformation displacement.
In the deformed porous media, the real velocity of fluid motion in the fluid-solid coupling seepage is [5] : The seepage takes place in the porous media, so the real velocity of fluid relative to the rock skeleton is: The Darcy velocity of fluid seepage is: There is a relationship as follows: So, the relationship between Darcy velocity a u  and the real velocity relative to the rock ra v  can be deduced as follows: Then the relationship among the real velocity, Darcy velocity and solid rock skeleton velocity is: Eq. (16) is the fluid-solid coupling seepage motion equation. While the velocity of rock skeleton particle is:

Fluid-solid coupling seepage continuity equations of gas hydrate reservoirs
Because the object of this study is the law of fluid-solid coupling seepage, we must use the real velocity and the real mass.
To one of the phases in the seepage system  , the mass of per unit volume is: aa a mS   ; While the mass flow rate is: aaa Sv   ; Then we take constituent i for example, and establish continuity equation.
Assume ig X as the mass fraction of constituent i in the gas phase, iw X as the mass fraction of constituent i in the liquid phase.
According to the principle of mass conservation, the continuity equation of constituent i can be obtained as follows: The above equation can be written in abbreviation as follows: ,, For the hydrate reservoirs, there are three phases in the reservoir pore space: liquid phase, gas phase and solid phase. Considering that there will be gas and water products when the solid hydrates decompose, we need to add one output in the continuity equation to reflect the decomposition products of hydrate.
Ignoring the dissolution of gas in the water, we can get the fluid-solid coupling seepage continuity equations of gas phase, fluid phase and solid hydrate in the hydrate reservoirs.
Gas phase: Fluid phase: Solid hydrate: Where, Φ is the absolute porosity of reservoir; ρ g , ρ w , ρ H are the densities of gas, water, hydrate, respectively; S g , S w , S H are the saturations of gas, water, hydrate; g v  , w v  , s v  are the real velocities of gas, water, rock skeleton, respectively; q g , q w are the source or sink terms of gas and water, while set source as positive, sink as negative; m g is the local mass rate of gas generation per unit volume of porous media; m w is the mass rate of water generated from the hydrate decomposition, and m H is the mass rate of methane hydrate decomposition.

Fluid-solid coupling seepage equations of gas hydrate reservoirs
The seepage of gas phase and water phase follow the Darcy law: Gas phase: Water phase: Where, g  is the acceleration of gravity, the value is   Water phase: Compared with normal gas/water seepage equations, the fluid-solid coupling seepage equations, (28) and (29), have two extra parameters, m l and   Sv lls .The first parameter m l reflects that the production from hydrate reservoirs by depressurization is a www.intechopen.com reflects the influence of reservoir skeleton deformation on the seepage field, meanwhile, that the porosity and permeability changed with reservoirs stress state reflects the coupling effect of seepage field and deformation field.
To simplify the questions, we assume that the deformation of reservoir skeleton is a steady deformation with small displacement and ignore the inertia force of rock skeleton. Then the individual derivative is similar to the partial derivative. So: Using Eq. (30), the continuity equation of reservoir rock matrix, namely Eq. (10), can be simplified as follows: According to the solid small deformation theory, we get: The first terms of the right side of Eqs. (34) and (35) are fluid convection flow along with matrix displacement. According to the small displacement theory, these can be ignored. The gas/water two-phase seepage equations (28) and (29) can be simplified further as follows: Gas phase: Water phase: The relationship between gas phase pressure P g and water phase pressure P w is as follows: Where, P c is capillary pressure, which is the function of water saturation.
Eqs. (45) and (47) are the general form of fluid-solid coupling seepage equations of gas hydrate reservoir with gas phase pressure and water saturation as their basic solution variables.

Initial and boundary conditions of seepage field equations
Initial and boundary conditions should be complemented to solve seepage field equations of fluid-solid coupling model. Also, for flow equations, the continuity condition and boundary condition must be satisfied in continuous domain.

Boundary flow continuity
Where, n  is bounding surface unit vector; g q and w q are gas and water phase flow per unit area on bounding surface.

Constant pressure boundary condition
The reservoir boundary pressure or bottom-hole pressure are known as following [6] : Advances in Natural Gas Technology

 
,,, f xyzt P is the given pressure function on boundary G at moment t.

Constant flow boundary condition
This means that the pressure derivative of the reservoir external boundary or the production in internal boundary is known, as following: Where, n is the direction of normal line;

 
,,, f x y zt q is the known function on specify boundary.

Initial condition
The initial conditions of seepage field are mainly the initial pressure and saturation of gas hydrate reservoirs, as following [6] : Where, P l (x, y, z) is the known pressure function, and S l (x, y, z) is the known saturation function.
The water-gas two-phase fluid-solid coupling seepage equations of gas hydrate reservoirs, namely Eqs. (47), can't be solved separately. Several simultaneous equations should be cited, including the kinetic equation of gas hydrates decomposition (1), the continuity equation solid phase hydrate (22), temperature field equations (17) and fluid-solid coupling deformation field equation in the next section. What's more, relevant state equation, assistant equations and the initial and boundary conditions, etc still need to be supplemented.

Deformation field equations of fluid-solid coupling model
The unknown parameters or unknown item in the seepage equations are affected by rock skeleton deformation, so the fluid-solid coupling seepage equations have highly nonlinear characteristics, and need to simultaneity the deformation field equations to solve the question. For porous media deformation field, the basic component equations mainly include equilibrium equation, geometric equation and constitutive equation, and the corresponding definite conditions. In this portion, the porous media deformation field mathematical model is constructed.

Equilibrium equations of reservoir skeleton
On the basis of the conventional force equilibrium principle of cell cube, we can get the deformation field balance equation with tensor form as follows: Fluid-Solid Coupling Numerical Simulation on Natural Gas Production from Hydrate Reservoirs by Depressurization

159
As the total stress supported by formation are composed of pore pressure and skeleton stress which is known as skeleton effective stress. The deformation and strength properties of reservoirs are controlled by the skeleton effective stress instead of the total stress.
We appoint that the tensile stress is positive and the compressive stress is negative. The Terzaghi effective stress formula is combined with: As to water-gas two-phase flow, the calculation of equivalent pore pressure is as follows: Zinkiewicz pointed that, for saturated fluid, the Biot coefficient can be calculated by: Where, K is rock volume modulus, and K s is rock skeleton volume modulus.
Substituting the modified effective stress equation (55) into Eq. (54), we can obtain the fluidsolid coupling deformation equilibrium equations based on the effective stress, which are the basic equations to solve the porous media deformation field.
Where, f i is the gravity term, in which f x is equal to f y , and the value of them is zero, f z is calculated by the formula  

Geometric equations of reservoir skeleton
Based on the solid small deformation theory, the relationship between strain component and displacement component can be described by following geometric equations with tensor form:

Linear elastic constitutive equations
For linear elastic media, stress and strain have a linear corresponding relation. The linear elastic constitutive relation can be expressed by generalized Hook's law with matrix form as follows: In the above equation, D e is the elastic constitutive matrix, which has 36 elements in 3 D coordinate. Elastic mechanics has proved it is symmetric matrices with only 21 independent constant. For the isotropic orthogonal media, there are only two independent elastic material constants, with elastic modulus E and Poisson's ratio signified.
Then the elastic constitutive matrix of rock skeleton can be expressed as:

Elastoplastic constructive equations
For elastic-plastic deformation, the strain is not only related to the current stress, but also to the loading history, loaded/unloaded condition, loading paths and microscopic structure of rock deformation, and so on.
In this study, the elastoplastic constructive relation is adopted to describe the relationship between strain and stress, which is expressed by the incremental form as follows: Eq. (62) can be represented by the matrix form as follows: The elastoplastic matrix can be calculated by the following equation.
The elastic matrix can be got by Eq. (61), and the plastic matrix can be deduced from the elastic-plastic associative flow rule, as follows: Where, A is the hardening index, and F is the yield function.
www.intechopen.com Rock is friction type material, the plastic yield of which need to consider both the shear yield and the volume strain yield. Therefore, the yield function can be expressed using the stress invariant as follows: Currently the Mohr-Coulomb criterion and the Drucker-Prager criterion are used for rock material in the engineering territory.

Modified Mohr-Coulomb criterion
The Mohr-Coulomb criterion can be deduced using the first stress tensor invariant I 1 , the second stress deviator J 2 and Lode stress angle θ σ , as follows: Where: Where, C is the cohesion of rock; φ is internal friction angle; σ 1 ,σ 2 ,σ 3 are three principle stresses.

Drucker-Prager yield criterion
The expression of the Drucker-Prager yield criterion is as follows: Where, α and K are the yield function parameters, the expressions of which are as follows: Comparing the two yield criterions above, the Mohr-Coulomb criterion neglects the effects of the middle principle stress σ 2 , and needs to determine the magnitude of each principle stress when it is used, which is not convenient. So we adopt the Drucker-Prager yield criterion in this research.

Initial and boundary conditions of deformation field equations
Assume that the space region occupied by rock skeleton is d  (the value of d is two or three), and the boundary of the region is  , in which the displacement boundary is u  The known stress boundary: The known displacement boundary: If the effective stress is adopted, so the expressions are: Stress boundary:

Permeability and capillary pressure
The effective permeability of gas hydrate reservoir is the function of hydrate saturation, and it meets the model proposed by Masuda (1997) [7,8] , as follows: Where, K 0 is the absolute permeability of the formation, of which the hydrate saturation is zero; S H is the hydrate saturation; N is the decline exponent of permeability, which is related to the types of hydrates combined in porous space, and the value is 2~15.
The relative permeability and capillary pressure curve can be obtained from the modified models proposed by Van Genuchten (1980) [9] and Parker (1987) [10] .
The relative permeability of the water phase: The relative permeability of the gas phase: The capillary pressure between liquid phase and gas phase: When the saturation of hydrates is 0.4, the relative permeability curve of gas and water phase and the capillary pressure curve are shown as Figs. 2-1 and 2-2.

Effective coefficient of heat conduction
The coefficient of heat conduction for gas hydrate reservoirs is related to the component. Assume that the heat conductions coefficient of rock skeleton, hydrates, water and gas are constant, so the effective coefficient of heat conduction K c can be calculated by the following equation [11][12][13] .

Specific heat capacity
The specific heat capacities of rock skeleton and hydrates are seemed as constants, while the specific heat capacities of water and gas change with temperature. When the pressure is 15MPa, and the temperature range from 273.15K～373.15K, we can get the following relationship (American National Standard) [ Where, C pg and C pw are the specific heat capacities of methane and water, J/Kg/K, which are related to temperature; C ph and C ps are the specific heat capacities of hydrates and rock skeleton, J/Kg/K.

Decomposition heat of methane hydrate
The decomposition heat of methane hydrate can be calculated by Masuda (1999) [15] model as follows: Where, Q H is the decomposition of hydrates; M H is the molar mass of hydrates; c and d is the experiment coefficients, the values are proposed 56599J/mol and -16.744J/mol/K for methane hydrates [16] .

Density of methane
The density of methane is calculated by Peng-Robinson state equation as follows: PR state equation is a mostly used cubic equation of state, which can calculate the gas volume and gas compression factor. PR equation can be expressed with compressibility factor Z, as follows: Where: Where, α(T) and b are the parameters of PR equation; T c is the critical temperature; P c is the critical pressure, and ω is the acentric factor.
Eq. (91) has three real roots, in which the maximum one is the gas compressibility factor, the minimum one is the liquid compressibility factor.
After obtaining the methane compressibility factor under the conditions of temperature T i and pressure P g , the density of methane ρ g can be calculated by the following equation:

Phase equilibrium equation
The phase equilibrium of methane hydrates, water and methane is calculated by Makogon (1997) model [17][18][19] , and the expression is: Where, P e is the equilibrium pressure, Pa; the value of T 0 is 27.15K; A, B, C are experiment coefficients, and the value is 0.0342K -1 , 0.0005K -2 , 6.4804, respectively.
What's more, the viscosities of water and gas are assumed as constant, which don't change with temperature and pressure.

Comprehensive dynamic models of property parameters for gas hydrate reservoirs
To resolve fluid-solid coupling numerical model for gas hydrate reservoirs not only needs equations of seepage fluid and deformation field but also needs relevant supplementary equations. The letter mainly include dynamic models of physical and mechanical parameters such as permeability, porosity and elastic modulus ,which reflect the real time change of the property at the process of exploitation for gas hydrate reservoirs.
Factors affecting physical property parameters of conventional reservoirs include volumetric strain of rock, effective stress, temperature, and so on. Usually one or more factors and property parameters including permeability,porosity and elastic modulus are used to build up a certain relation, that is to build up the dynamic models of physical property parameters.
For gas hydrate reservoirs, hydrate decomposition is the most direct and notable factor to affect physical properties. So in order to build up comprehensive dynamic models of physical properties for gas hydrate reservoirs, two factors must be considered: firstly, it must take example from dynamic models of conventional reservoirs to reflect the relationship between physical property parameters including permeability, porosity and elastic modulus and physical quantity including volumetric strain, effective stress and temperature; secondly, it need to highlight the influence of gas hydrate decomposition effect on related physical parameters.

Comprehensive dynamic model of permeability
Experiments of permeability stress sensitivity for loose sandstone done by Yizhong Zhao [20] indicate that where far away from borehole, the permeability is mainly controlled by principal stress and the fit accuracy of permeability and effective stress is high. The model built by them is indicated as Eq. (95); since sensitivity of nearly well reservoirs is mainly controlled by stress difference, relativity between permeability and volumetric strain is better. So for nearly well reservoir we could adopt the model of permeability and volumetric strain, as the following equation (96) shows.
Where, K σ is permeability of some stress state, 10 -3 μm 2 ; K is the initial permeability, 10 -3 μm 2 ; σ is the effective stress, MPa; c, d and e are regression coefficients of experiment.
Where, ε v is the volumetric strain; φ e is the initial effective porosity; f, j and l are the regression coefficients of experiment.
Combined with the permeability computation model Eq. (77) for gas hydrate reservoirs and considering the difference of permeability effect mechanism for each part of the reservoir, this study builds up two comprehensive dynamic models of permeability for gas hydrate reservoirs. The formula (97) is fit in reservoirs which are far away borehole; however, the formula (98) is fit in the relative analysis of nearly well reservoirs.
The model built up by this study roundly reflects the influence of gas hydrate decomposition effect and change of reservoirs stress state on permeability.

Comprehensive dynamic model of porosity
JianJun Liu, Chinese academy of sciences Seepage institute, and JishunQin, China University of Petroleum, who conduct a lot of experimental research on the relationship between the porosity of low permeable porous medium and the effective stress, proposed that the index relationship could be used to describe the changing regularity between them.
Where, Φ σ is the porosity at some stress state;Φe is the initial effective porosity; m and n are experiment coefficients; σ is the effective stress, MPa.
Palmer suggests that the relationship between porosity of unconsolidated sandstone and the effective stress can be described by index relationship, whose form is similar with the formula (99). When Yizhong Zhao studied the reservoir of unconsolidated sandstone, they optimize the index relationship which is similar with (99) as dynamic models of porosity.
Supposing the distribution of porosity is uniform in gas hydrate reservoirs and ignoring the influence of temperature, so the following relationship can be obtained.
Where, Φ is the absolute porosity when the saturation of hydrates is zero; S H is the saturation of hydrates.
In conclusion, the comprehensive dynamic model of porosity for hydrate reservoirs built for this study is: This model above reflects the effect of hydrate decomposition effect and the influence stress state change on the porosity.

Comprehensive dynamic model of elastic modulus
The research of Guangquan Li, Jianjun Liu etc shows that the relation between elastic modulus and effective stress accords with power function form, as follows: www.intechopen.com Where, E σ is elastic modulus at some stress state, GPa; E is the initial elastic modulus, GPa; a and b are the fitting coefficients of experiment; σ is the effective stress, MPa.
By making sensitivity experiments of elastic modulus for unconsolidated sandstone, Yizhong Zhao points that the regression relation of quadratic polynomial between elastic modulus and effective stress is more accurately, the expression of which is as follows: Where, A, B, C is the fitting coefficients by experiment.
In recent years,through experimental research,the petroleum research institute of Commonwealth Scientific and Industrial Research Organization (CSIRO) and Tohidi, from hydrate research center of Heriot-Watt University point that the decomposition of hydrates will result in the increase of porosity and the decline in the elastic modulus correspondingly, based on which they built the relation between elastic modulus change and the variation of porosity, as follows: Where, is the fitting coefficient by experiment, the value of which is -9.5; ΔE is the decrease of elastic modulus, GPa; ΔΦ is the increment of porosity.
Based on the comprehensive dynamic model of property for gas hydrate reservoirs which we built for this research, the mathematical relation between the hydrate decomposition effect and the increase of porosity resulting from the change of stress condition is as follows: Where, Φ is the absolute porosity when the hydrate saturation is zero; S Hi is the initial saturation of hydrates; S H is the current hydrates saturation; σ is the effective stress, MPa; m and n is the experimental matched coefficients.
Combining Eq. (104) into Eq. (105), the change of elastic modulus is as follows: At the same time, there is the following relationship: Where, E 0 is the elastic modulus before hydrates decomposition, GPa.
Allying Eqs. (106) and (107)  Because the change of Poisson ratio is little, the Poisson ratio is seen as constant during the fluid-solid coupling numerical simulation.

Comprehensive dynamic model of cohesion
During the exploitation of gas hydrate reservoirs, hydrates decomposition can cause the cementation of the reservoirs loosened and the remarkable change of the cohesive force must be taken into consideration.
In 2005, Clennell, from the petroleum research institute of Commonwealth Scientific and Industrial Research Organisation's (CSIRO), and Tohidi, from Heriot-Watt University hydrate research center, proposed that as the decomposition of hydrates, the reservoirs cementation will become weaker, and the cohesion will decrease continuously. So they built the model of cohesion decrease for reservoirs, which was widely adopted by scholars. The mathematical expression is [21][22][23] : Where, C is the cohesive force of the reservoirs after hydrates decomposition, MPa; C 0 is the cohesive force before hydrates decomposition, MPa; ΔΦ is the increment of porosity as a result of hydrates decomposition.
Eq. (109) is also taken as the dynamic model of cohesive force for the production from hydrate reservoirs by depressurization.

Fluid-solid coupling simulation studies on productivity of gas production from hydrate reservoirs by depressurization
Gas production from hydrate reservoirs by depressurization is a complex fluid-solid coupling process with hydrates decomposition, which is affected by the hydrates decomposition effect, fluid-solid coupling and borehole effect. In this process, the reservoir physical parameters shows the complex characteristics, reservoir porosity characteristics change and fluid-solid coupling have the dynamic effect for gas production from hydrate reservoirs by depressurization.
Based on the gas and water fluid-solid coupling mathematical model and the physical and mechanical models build above, this section will solve these models. Taking a gas hydrate reservoir of Gulf of Mexico as an example, we simulate the productivity of gas production from hydrate reservoirs by depressurization. At the same time, the main influence factors to productivity are also researched.

Basic parameters used in the simulation
The parameters used in fluid-solid coupling numerical simulation include basic property parameters, reservoir fluid and the basic properties of solid hydrate. All of the parameters input in the simulation of this study are from a natural gas hydrate reservoir in the Gulf of Mexico, as is shown in Table 3-1 and Table 3

Geometrical model and boundary conditions
To simplify the problem, the plane strain model is adopted considering the symmetry of the hydrate reservoir. Set the overall dimensions as 100m×100m for this finite element model. With respect to the significance of the near borehole formation to the whole capacity simulation, the model meshing is obtained with gradient grid technology, so the mesh density is higher in formation near borehole. The diagram of mechanical model and grid of the finite element model are shown in Figs. 3-3 and 3-4, respectively.
In the process of simulation the boundary conditions for the seepage field are as follows: Line BC and CD keep constant pressure boundary, Line AB and DE keep free pressure boundary, and at the position of the wellbore, namely Point A, the boundary keeps the bottomhole pressure.
The boundary conditions of the deformation field are listed as follows: The maximum effective horizontal principle stress (σ H ) is applied on the boundary line BC, and the minimum effective horizontal principle stress (σ h ) is applied on Line CD. On the Line AB the displacement is constant as zero along X-axis and free along Y-axis; And on the Line DE the displacement is zero along X-axis and free along Y-axis. At the Point A the displacement is fixed to both X-axis and Y-axis.

Fluid-solid coupling numerical simulation on productivity of gas production from hydrate reservoirs by depressurization
For the further studying the influential mechanism of fluid-solid coupling effects on the gas production from hydrate reservoirs by depressurization, this section also analyses the other simulation results of the two simplified models which is based on the gas-water two-phase non-isothermal fluid-solid coupling model established above.
The first model ignores physical parameters change caused by the fluid-solid coupling effect, but it considers the coupling effect between seepage field and deformation field, and it is named the simplified fluid-solid coupling model. Based on the first model, the second model further ignores the coupling effects of seepage field and solid field, so it is named the non-coupling normal model.
For the convenience of analysis, here the gas-water two-phase non-isothermal fluid-solid coupling model established above is renamed as the comprehensive fluid-solid coupling model. The basic characteristics of the three models are shown in Table 3-3.

Numerical simulation on productivity of non-coupling normal model for production by depressurization
Taking into account the non-coupling normal model with phase transition in Table 3-3, which considers the petrophysical properties change as a result of hydrates decomposition and neglects the coupling effect between the seepage field and the deformation field, the gas and water production rate are simulated. The results are shown Figs. 3-5 and 3-6, the unit of which is STCMD.   shows that the gas production rate curve fluctuates greatly because of the instability of gas-liquid two-fluid phase seepage. The variation tendency of gas production rate curve indicates that the gas production rate increases sharply at the early stage, and become stable lately. Only 43 days is set to exploit which is too short to finish hydrates decompose, but it can be indicated that gas production rate will decrease as hydrate saturation reducing, and finally no gas is produced. In general, the gas production rate curve can be divided into three stages, rapidly rising stage, stable stage and dropping stage. We can summarize the characteristic as follows: 1. The rising stage of gas production rate. As the pressure releases the gas production rate rises rapidly. This is because: a. the initial saturation of free gas is 0.1 and when the pressure decreases, the free gas could swell out quickly; b. the pressure of the reservoir prior to production is close to hydrate decomposition pressure at the initial temperature, so the pressure around wellbore will be below the equilibrium pressure of hydrates rapidly and they begin to decompose immediately. 2. The stable stage of gas production rate. Most of the gas produced after free gas output is from the hydrate decomposition, which is controlled by the rate of hydrate dissociation. 3. The dropping stage of gas production rate. The hydrate saturation become lower and lower as the decomposition front going ahead further, and the gas production becomes less correspondingly, and finally no gas could be produced. As only 43 days is set to simulate the gas production in Fig. 3-5, which is not long enough to lead to hydrate decomposition totally, so the obvious dropping stage of gas production is not reached yet in fact.
From Fig. 3-6 it can be seen that the water production rate curve is similar with the one of gas production. It can also be divided into three stages as follows: 1. Instead of gas, water is the only production in the original depressurization stage, and the time of water producing falls behind the gas producing time. This is because the existence of free gas in the reservoir has better flowing ability than the water when the pressure decreasing. 2. After the output of free gas, the water saturation of reservoirs increases gradually as hydrates decompose. The rate of water production decreases rapidly after reaching a peak under the influence of the relative permeability of water and gas. The reason for this phenomenon is as follows, the free water and produced water accumulate towards the wellbore as a result of the propelling of gas, and they come out rapidly when reach to some extent. Consequently, a water production peak takes place. However, due to the fact that both the gross of the free water and the water comes from the hydrate decomposition in the original stage are relatively limited, the rate of water production decreases speedily after reaching the peak. 3. The water producing rate remains steady level after falling. In this period, as the free water output is almost completely, the production is mainly the decomposing water from hydrates, which is controlled by the decomposition rate of hydrates. And the water production exhibits the similar character to that of gas. The decrease of the hydrate saturation and decomposition rate leads to the reduction of the water producing rate, and finally no water is generated.
Shuxia Lee and Yongmao Hao, who are from the Hydrate Research Center of China University of Petroleum, established hydrates synthesizing and exploiting equipment. They obtained the water and gas production rate characteristic curve through small scale hydrate exploiting experiment by depressurization [24,25] . Figs. 3-7 and 3-9 show the results. The size of the pipe adopted in the simulation experiment by ShuXia Lee is Φ38×500mm and Φ38×800mm, which belongs to small one-dimensional model. However, in the numerical simulation of this study a plane strain model is used. The difference of the two models leads to no relative property of the gas and water production rate value. But it is significant to analyze the comparison of the change regulation. Because the model dimension is large and the time is limited in this research, the dropping stage of gas and water production rates have not been happened. Comparing The cumulative gas production is shown as Fig. 3-9. It can be seen that the cumulative gas production curve is well corresponding with the gas production rate curve. At the early stage of depressurization the gas production rate is bigger, and correspondingly the cumulative gas production rises quickly. When the gas production rate reaches the relatively steady stage, the cumulative gas production increases steadily and slowly. It can be inferred that, with gas production rate gradually decreasing at the later stage of depressurization production, correspondingly the increasing range of cumulative gas production will gradually decrease. Until the hydrates decomposition totally, the cumulative gas production finally stays at a certain level.
The cumulative water production is shown as Fig. 3-10. It is easy to see that, the cumulative water production curve is vary similar with the one of gas production, and it is also well corresponding with the water production rate.

Numerical simulation on productivity of simplified fluid-solid coupling model for production by depressurization
During the fluid-solid coupling effect of production by depressurization, the influence of reservoir deformation field to seepage field mainly reflects in two aspects. On the one hand, the deformation of reservoir skeleton will shrink the pore space, which will increase elastic drive energy. On the other hand, the change of effective stress state will cause petrophysical variations of stress sensitive reservoirs, such as porosity, permeability and so on. All of these will affect the production performance of hydrate reservoirs.
In this chapter, the simplified fluid-solid coupling model in Table 3-3 is adopted to simulate gas production from hydrate reservoirs by depressurization, in which the physical parameters alternation caused by fluid-solid coupling is neglected. By comparing with the simulation results of the non-coupling normal model above, we discuss the influence of fluid-solid coupling effects on the production behaviors of hydrate reservoirs. Fig. 3-11 gives the gas production rate curves of both the non-coupling normal model and simplified fluid-solid coupling model. Fig. 3-13 gives the water production rate curves of the two models. The cumulative gas and water production curves of the two models are shown in Figs. 3-12 and 3-14.
From the comparison above, it can be seen that the results of the simplified fluid-solid coupling model and non-coupling normal model are similar in some aspects, but there is also some obvious difference, as follows: 1. Similar to the non-coupling normal model, the gas and water production curves of the simplified fluid-solid coupling model can be divided into three stages generally: rapidly rising stage, stable stage and gradually reduce stage. This is consistent with the results of experiments by Shuxia Lee and Yongmao Hao in China University of Petroleum.
Fluid-Solid Coupling Numerical Simulation on Natural Gas Production from Hydrate Reservoirs by Depressurization 177 Fig. 3-11. Comparison of gas production rate 2. The major difference between two models lies in that both the average production rate and the accumulative production of gas obtained by the simplified fluid-solid coupling model are 3.45% higher than that of the non-coupling normal model. Similarly, both the average water production rate and the water accumulative production are 4.28% higher than that of non-coupling normal model. This is because reservoir skeleton volume strain can lead to reservoir porosity reducing and elastic drive energy increasing to some extent. This is flavor to raise gas and water production rate, accumulative gas and water production of hydrate reservoirs. The increasing of elastic drive energy has more effect to water production than gas production as a result of the poor liquidity of the gas with respect to the liquid.

Numerical simulation on productivity of comprehensive fluid-solid coupling model for production by depressurization
The analysis above shows that shrinkage of reservoir pore space can raise the elastic drive energy because of the fluid-solid coupling effect, which does help to increase the gas production rate and accumulative gas production. This is one side that the fluid-solid coupling affecting production performance. And on the other side, the change of the reservoir effective stress can lead to physical parameters, such as porosity and permeability, which can also influence production behaviors.
In this chapter, we adopt the comprehensive fluid-solid coupling model in Table 3-3 to simulate gas production from hydrate reservoirs by depressurization. Based on the simplified model, the comprehensive fluid-solid coupling model considers the change of physical parameters caused by the stress state variation. Thus, by comparing with the two models, we can study the influence of the physical parameters change of sensitive stress reservoirs, which is caused by the fluid-solid coupling effect, on the production performance of hydrate reservoirs. The comparisons above show that the gas and water production rate curves of the comprehensive and simplified coupling models can also be divided into three stages, namely the rapidly rising stage, stable stage and gradually reduce stage, which is similar to that of the non-coupling normal model, and matches with the experiment results of Shuxia Lee. Further analysis shows that the numerical size of simulative result has significantly difference, as follows: 1. From Figs. 3-15 and 3-16 it can seen that the average gas production rate and accumulative production which are obtained by the simplified fluid-solid coupling model are obviously higher than that of the comprehensive fluid-solid coupling model. The average gas production rate in stable stage of the former model is about 1.144 times greater than that of the later. Accordingly, the gas accumulative production is about 4.4% higher than that of the later. 2. Similarly, from Figs. 3-17 and 3-18 we can also see that the average water production rate and accumulative production simulated by the simplified fluid-solid coupling model is obviously higher than that of the comprehensive fluid-solid coupling model. The average water production rate in the stable stage of the former model is about 1.159 times greater than that of the later. Correspondingly, the accumulative water production is about 15.9% higher than that of the later.  show that the average gas production rate and accumulative gas production of the non-coupling normal model are obviously higher than that of the comprehensive fluid-solid coupling model. The average gas production rate in the stable stage of the former model is about 1.106 times higher than that of the later. And accordingly, the accumulative gas production of the former is about 15.9% higher than that of the later. Similarly, from Fig. 3-21 and Fig. 3-22 it can be seen that the average water production rate and accumulative water production of non-coupling normal model is also obviously greater than that of comprehensive model. The average water production rate during the steady stage of the former model is about 1.111 times higher than that of the later. Correspondingly, the accumulative water production of the former is about 11.1% higher than that of the later.
In conclusion, by analyzing the production performances of the non-coupling normal model, the simplified coupling model and the comprehensive, we can find that the development index, namely the gas production rate, the accumulative gas production and so on, obtained from the simplified model is the highest in the three models. The results of the non-coupling normal model are slightly lower. But the development index simulated by the comprehensive coupling model is significantly lower than the two former models. The reason for this phenomenon is that fluid-solid coupling has different influencing mechanism to production performance.
During the production by depressurization, the fluid-solid coupling effect has two aspects. On one side, the shrinkage of rock matrix caused by fluid-solid coupling effect can raise the elastic drive energy, which is helpful for gas output. So the development index of simplified fluid-solid coupling model which only considers the skeleton volume strain is higher than that of the non-coupling normal model. This is the good aspect that fluid-solid coupling effect contributes to gas production. On the other side, the shrinkage of rock pore space can also cause the reduction of permeability and porosity, and accordingly the seepage captivity also decreases which causes the adverse effect for production. Relatively, the change of porosity and permeability which is caused by fluid-solid coupling plays the dominated role on affecting the production performance. But the increasing of elastic drive energy is limited. Therefore, the predictive value of production index simulated by the comprehensive fluid-solid coupling model which considers the skeleton volume strain and physical parameters change is obviously lower than the results of the simplified fluid-solid coupling model and non-coupling normal model.
The variation of skeleton volume strain and petrophysical properties caused by fluid-solid coupling in the reservoir is the objective existence. So, the overall effect of fluid-solid coupling causes the development index, namely the gas production rate, the accumulative gas production and so on, is obviously lower than that of the non-coupling normal model. Therefore, from the productivity point of view, the coupling effect between the seepage field and the deformation field must not be neglected in the exploitation of hydrate reservoirs. Although there are many numerical models for production from hydrate reservoir in previous research, they failed to totally consider the coupling effect between the seepage filed and rock deformation field.
Compared with the current numerical models, the fluid-solid coupling model, which this study established for the production from hydrate reservoirs production by depressurization, has the superiority as follows: fully consideration of hydrate decomposition effect, the coupling effect of seepage field and deformation field, the change of the temperature, and the dynamic variation of the physical and mechanical properties of formations under the influence of many factors.
In the next portion, based on the gas-water two-phase non-isothermal fluid-solid coupling model established in this research, we will carry out a sensitive analysis on the main factors which affect the productivity from hydrate reservoirs.

Influences analysis on productivity of gas production from hydrate reservoirs by depressurization
The factors that affect the production performance for hydrate reservoirs by depressurization are very more and the exten t o f t h e m i s v a r i o u s . I n o r d e r t o comprehensively evaluate the influence degree of all factors on the development index and provide the theoretical basis for the production optimization of hydrate reservoirs, we primarily study the influence of the reservoir absolute permeability, the bottomhole pressure and the reservoir temperature on hydrate reservoir production index, such as the accumulative gas production.

The effect of absolute permeability of formation
The absolute permeability of formation is a parameter reflecting the seepage ability of porous medium and it is an important factor affecting the production rate after gas hydrates decomposing. In this portion, we use the non-isothermal fluid-solid coupling model developed in this article to simulate the production behaviors with the condition of three reservoir absolute permeability, namely 100mD, 200mD and 300mD. In the simulation the other conditions of hydrate reservoir keep the same. The emphasis is to analyze the influence of formation absolute permeability on production index, such as the hydrate decomposing rate, the gas production rate and accumulative gas production. It is not difficult to find that the influence of absolute formation permeability is obvious to production performance of hydrate reservoir. The greater the formation absolute permeability is, the faster the propagate speed of the pressure drop is in the reservoir, the From the results of the figures above, it is easy to find that despite of the different permeability of the formation, the gas and water production rate curves are similar, which can be divided into three stages on the whole, namely the rapidly rising stage, the stable stage and the gradually reducing stage. The greater the absolute permeability of reservoir is, the faster the propagate speed of pressure drop is, the faster hydrates break down, the greater the gas and water production rates are, and accordingly the bigger the accumulative gas and water productions are. After the well produces for about 200days, the accumulative gas production with the absolute permeability being 300mD is 7.05% higher than that of the formation with permeability being 200mD. The accumulative water production is 7.4% higher than that of the formation of which the permeability is 200mD. The accumulative gas and water production of the formation with the absolute permeability being 200mD is 5.9% and 6.34% respectively higher than that of formation whose absolute permeability is 100mD.

www.intechopen.com
Fluid-Solid Coupling Numerical Simulation on Natural Gas Production from Hydrate Reservoirs by Depressurization 185

The effect of bottomhole pressure
The value of bottomhole pressure is an important parameter to control the gas production from hydrate reservoirs reasonably, and it not only affects the cumulative gas and water production of hydrate reservoirs, but also determines the propagate speed of hydrate dissociation front. Based on the non-isothermal fluid-solid coupling model above, this portion simulate the production performance under the conditions of the initial formation pressure 16.9MPa and three bottomhole pressure, namely 12.9MPa, 13.9MPa and 14.9MPa respectively. Other conditions in the progress keep the same. Emphatically the influence of bottomhole pressure on production index, such as the hydrate decomposing rate, the gas production rate, the accumulative gas production, and so on.   Comparative analysis show that, similarly to the former, the gas and water production rate curves with different bottomhole pressures generally can be divided into three stages: rapidly rising stage, stable stage and gradually reducing stage. And the lower bottomhole pressure, correspondingly, leads to the faster hydrate dissociation, the bigger gas and water production rate and the greater gas and water accumulative production. After 200 days production, when bottomhole pressure is 12.9MPa, gas accumulative production is 7.9% higher than that with the bottomhole pressure being 13.9MPa, and accordingly, the water accumulative production increase 9.6%. When bottomhole pressure is 13.9MPa, the gas and water accumulative production increase 6.9% and 8.7% respectively compared with that being 14.9MPa.
It should be pointed out that low of the bottomhole pressure is contributive to improve the gas accumulation production and other development indexes. But at the same time, the lower bottomhole pressure can lead to the greater producing pressure difference and the worse stability in near well reservoir. Therefore, the stability of the reservoir near well bore should be considered in optimizing design of bottomhole pressure. To improve the production rate of hydrate reservoirs, the bottomhole pressure can not be reduced without restraint.

The influence of formation temperature
The formation temperature is a key parameter to affect the hydrate dissociation rate. Using the non-isothermal fluid-solid coupling model of the hydrate reservoir, this portion will simulate the production performance when formation temperatures are 289.5K, 290K and 290.5K respectively and other conditions are the same. We analyze emphatically the influence of formation temperature to the hydrate dissociation rate, gas production rate and cumulative gas production, etc.
Figs. 3-35 and 3-36 show the location of hydrate dissociation front and the distribution of the saturation of hydrate (SH) after 43 days production under different formation temperatures. Comparisons indicate that the higher formation temperature leads to the faster hydrate dissociation, the greater dissociation front promoting rate and the further forward distance during the same time. As Fig. 3-35 shows, when the formation temperatures are 289.5K, 290K and 290.5K respectively, the hydrate dissociation fronts are 34.105m, 45.08m and 56.08m respectively far away from borehole after 43 days production similarly. Comparisons indicate that, as the same with the previous, the gas and water rate curves can be divided into three stages. They are the rapidly increasing stage, the stable stage and the gradually decreasing stage respectively. The higher formation temperature leads to the faster hydrate dissociation, the larger gas and water production rate. After 200 days production, the gas accumulative production with the formation temperature of 290.5 K is 7.0% higher than that of 290 K, and correspondingly, the water accumulative production increases 8.8%. And, the gas accumulative production with the formation temperature of 290 K is 5.9% higher than that of 289.5 K, correspondingly, the water accumulative production increase 6.9%.
parameters and the stress state, and at the same time it embodies the influence of hydrates decomposition effect. The factors considered in this model are full. And the model can lay a foundation for the fluid-solid coupling numerical simulation of depressurization for gas hydrate reservoirs. 4. The fluid-solid coupling numerical simulation was carried out, and the influence factors were analyzed for the productivity of hydrate reservoirs by depressurization. Results show that the comprehensive effect of the fluid-solid coupling model leads to the productivity from natural gas hydrate lower than that of the non-coupling model. The influence factors analysis reveals that the permeability and bottomhole pressure affect hydrate decomposition rate by changing the propagation rate of pressure drop, but the formation temperature affects hydrate decomposition rate by changing the hydrate phase balance pressure. The production index such as the gas production rate and the cumulative gas etc rises with the increase of reservoir absolute permeability and temperature, and with the decrease of bottomhole pressure.

Nomenclature
m g = local mass rate of gas generation per unit volume of porous media A dec = the specific hydrate decomposition surface area in the porous media M = the molar mass f e and f g = the fugacities of methane gas at equilibrium pressure (P e ) and gas phase pressure (P g ), respectively Φ e and Φ f = the fugacity coefficients of methane gas at P e and P g , respectively