Numerical Simulation of Shale Gas Multiscale Seepage Mechanism-Coupled Stress Sensitivity

,e complexity of the gas transport mechanism in microfractures and nanopores is caused by the feature of multiscale and multiphysics. Figuring out the flow mechanism is of great significance for the efficient development of shale gas. In this paper, an apparent permeability model which covers continue, slip, transition, and molecular flow and geomechanical effect was presented. Additionally, a mathematical model comprising multiscale, geomechanics, and adsorption phenomenon was proposed to characterize gas flow in the shale reservoir. ,e aim of this paper is to investigate some important impacts in the process of gas transportation, which includes the shale stress sensitivity, adsorption phenomenon, and reservoir porosity.,e results reveal that the performance of the multistage fractured horizontal well is strongly influenced by stress sensitivity coefficient. ,e cumulative gas production will decrease sharply when the shale gas reservoir stress sensitivity coefficient increases. In addition, the adsorption phenomenon has an influence on shale gas seepage and sorption capacity; however, the effect of adsorption is very weak in the early gas transport period, and the impact of later will increase. Moreover, shale porosity also greatly affects the shale gas transportation.


Introduction
Shale gas is an unconventional energy with extensive distribution and large reserves, which can alleviate the pressure of global natural gas supply and demand gap.Furthermore, it can also provide reliable energy supply for the sustainable development of the world economy.However, the understanding of the shale gas flow mechanism is far behind production practice.Many relevant studies show that the gas transportation exhibits the characteristics of multiscale.In addition, the influence of stress sensitivity on seepage flow cannot be neglected.erefore, it is very important to understand the multiscale seepage mechanism of shale gas under stress sensitivity.
Many previous works have been performed to understand the characteristics of gas flow and stress sensitivity in shale gas [1][2][3].Huang et al. [1] shows that gas transportation in shale cannot be simply characterized by the Darcy law, and a method by calculating the Knudsen number was proposed to analyze the shale gas flow regime.Additionally, a unified model taking into account slip and real gas effect was proposed to model gas transport behavior for different gas transport mechanisms in nanopores [4].Zhang et al. [5] studied the stress sensitivity characteristics of shale with experiments, and the research results showed that the stress sensitivity of the shale core used in this experiment was generally higher than that of tight sandstone.Meanwhile, the study also pointed out that the compressibility of shale is much higher than that of sandstone due to the strong compressibility of shale nanopores, microfractures, and clay minerals.Besides, the gas flow mechanism of fractured vertical well with stress sensitivity is also studied [6].However, little work has been done to investigate multiscale gas flow coupling stress sensitivity.
In this paper, a mathematical model comprising multiscale, geomechanics, and adsorption phenomenon was proposed to characterize gas flow in the shale reservoir.Firstly, the characteristics of shale gas reservoirs are analyzed, and the physical parameters of gas at high pressure and high temperature are calculated by using the Peng-Robinson state equation.Secondly, the seepage mathematical model is solved by using nite di erence, and the results are compared with commercial software.Furthermore, the correction coe cient of shale permeability with stress sensitivity and multiscale seepage is discussed under di erent pressures and pore sizes.Finally, to evaluate the impact of stress sensitivity, adsorption, and porosity on the shale gas mechanism, several cases concerning production dynamics of multistage fractured horizontal wells have been carried out.

Shale Reservoir Characterization
With the application of high-resolution electron microscopy in investigation of shale pore structures, some studies have illustrated that the pore structure of shale is complex and behaviors show strong heterogeneity, and the pore types in shale are classi ed as intergranular pores, intragranular pores, some dissolved pores, and microfractures [7,8].Due to the di erence and diversity of pore development, the distribution of the shale pore diameter varies from 1 nm to 1000 nm [9][10][11].Moreover, the results of core analysis show that the permeability of shale varies from 10 −6 mD to 1 mD [11].In a word, the physical property and structural characteristics of the shale reservoir directly a ect the transportation mechanism of shale gas.e physical model of shale gas production and transport is shown in Figure 1.Matrix pores, natural microfractures, and hydraulic fractures consist of the network channel for gas transfer [11].Among them, the matrix accounts for the main reservoir space, and the fracture plays a key role in the transportation channel.Free gas and adsorbed gas transport from natural fractures to hydraulic fractures under pressure di erence.In addition, when gas migrates in the nanoscale pore, the gas molecules collide seriously with the pore roaring wall of the shale, slippage, transition, and free-molecular ow will be included in microscale and nanosacle pore space besides Darcy ow.
erefore, the transportation of shale gas involves di erent spatial scales.
e occurrence state of shale gas can be divided into the form of adsorption gas, free gas, and dissolution gas.And several investigations have illustrated that the adsorption gas may account for 20%∼80% of the total reserves in shale gas [12].us, the adsorption capacity of shale is very important for shale gas reservoir evaluation, and the physical changes of adsorption gas may have an important impact on shale production behaviors.To characterize the relationship between adsorbed gas and pressure in shale reservoirs, a mass of experiments have been done to understand the characteristics of shale gas absorption in recent years.And several models are proposed to characterize the shale gas adsorption phenomenon.Specially, the Langmuir adsorption equation has been widely used because of its simple form and slight error [12,13].Furthermore, the parameters of the Langmuir equation can be easily obtained by tting the experimental data.e mathematical form of rock adsorption volume V E is represented as follows: where V L and P L is the Langmuir pressure and Langmuir volume, respectively, which can be obtained by tting the isothermal adsorption experimental data, and the isothermal adsorption curve of shale measured in the laboratory is shown in Figure 2.
In fact, there is no di erence in the composition between the adsorbed gas and the free gas.According to the component analysis of the shale gas produced, the study shows that the shale gas is a natural gas mainly composed of methane, and the methane content is up to 98% [14].It is a typical dry gas reservoir.With the change of temperature and pressure, the methane state is divided into the solid state, liquid state, gaseous state, and supercritical state [15].As shown in Figure 3, it can be seen from that the actual shale reservoir temperature is much higher than the critical temperature of methane, which means the shale gas is in a supercritical state.
In addition, we calculate the high-pressure physical parameters of methane by the Peng-Robinson equation of state [16], which is expressed as follows: where (3) Besides, the equation can also be written as follows:  Journal of Chemistry where e gas density can be calculated: And gas viscosity is calculated by Lee's empirical model [17]: where K 9.379 + 0.016M g (1.8T) e calculation results are shown in Figures 4-6, including deviation factor, density, and viscosity.is part of work provides gas physical parameters for establishing the transport model of the shale gas reservoir.

Apparent Permeability Model with Stress Sensitivity.
To describe the variation of gas transport in porous media, researchers have proposed to use the Knudsen number to classify the shale gas ow regime [18].Knudsen number was de ned as the ratio of mean free path of gas molecules to pore size of porous media.Flow regimes are divided into four categories, as shown in Table 1. Figure 7 is the Knudsen number corresponding to di erent pore sizes under di erent pressure conditions.It can be seen that the ow regime includes continuous ow, slippage ow, and transition ow in the shale gas reservoir.erefore, the transport mechanism in shale could not be characterized by the classical Darcy law.
So, an apparent permeability model based on Beskok-Karniadakis (BK model) theory is widely used to characterize the multiscale transport of shale gas [19][20][21][22][23][24], and the model can be expressed as follows: Furthermore, when the pressure changes during the transportation of shale gas in porous media, the stress sensitivity and adsorption will cause the variation of shale pore radius, which lead to the change of shale apparent permeability.In view of the above phenomenon, Guo et al. [25] modi ed the Beskok-Karniadakis permeability model (BKG model) based on the capillary bundle model.And the modi ed Beskok-Karniadakis permeability model can be illustrated as follows: where where σ eff is the e ective stress rock, which can be calculated as follows: To clearly understand the e ect of stress sensitivity and adsorption on apparent permeability, the correction coe cient of permeability ξ is de ned in the following form: ξ

Multiscale Transport Model of Shale Gas.
In consideration of the actual reservoir formation containing nature fracture networks and matrix porosity, the matrix and the fracture system are equivalent to a porous medium.Additionally, to make this mathematical model more tractable and easier to understand, the following assumptions are made: (1) e gas component is only methane, and it ows at a constant temperature.(2) e gas exists in the form of adsorption and free.
And the desorption phenomenon can be characterized by Langmuir isothermal adsorption theory; moreover, the process of desorption and adsorption can reach an instant balance.
(3) Gas reservoir model is homogeneous and isotropic, ignoring the compressibility of rock.( 4) e e ect of reservoir stress sensitivity and desorption on shale permeability is taken into account.
Based on the above assumptions and according to the law of mass conservation, the continuity equation of shale gas can be derived as follows: Substituting equation (15) into equation ( 14), we can obtain the governing equation of multiscale migration of shale gas considering stress sensitivity: Initial conditions and boundary conditions should be added to solve the partial differential equation.e initial condition of the gas reservoir is the distribution of pressure at each point.Since the gas is in a static equilibrium state before transport, the pressure at any position is the original formation pressure: When there is no edge and bottom water around reservoir, it can be considered that the gas reservoir is a closed boundary condition, which is mathematically known as the second boundary condition or Neumann boundary condition: According to the working system of the actual shale gas production wells, the inner boundary conditions usually include the conditions of fixed-output gas quantity and fixed-bottom flow pressure production.As we all know, the well radius is far less than the well distance, and thus, the production well is treated as the source in the numerical simulation.When the bottom hole pressure is determined, the Peaceman well model can be adopted to correspond the bottom hole flow pressure and gas production volume, and then, the source term or sink term can be added to the differential equation [26], which can be describes as follows: As for isotropic formations, the equivalent supply radius is r e � 0.14

Numerical Solutions of Mathematical Model.
In this part, the finite difference method is used to discretize the above differential equation of shale gas seepage flow.And local grid refinement is used to describe the hydraulic fractures.Figure 8 represents a simplified numerical simulation flowchart of shale gas seepage in the porous media.To approve the accuracy of the developed numerical model presented in this paper, production behavior for a vertical well in a convention gas reservoir is calculated with a simple case of the presented model, and the results are validated with commercial software, which is an adaptive implicit black oil simulator.e basic parameters for the gas reservoir model are shown in Table 2. Gas adsorption, multiscale flow mechanism, and stress sensitivity are not considered in the model.Figure 9 represents the match result of gas production rate, and Figure 10 shows the distribution of pressure on model profile.Simulation results of this work match well with the commercial simulation software.

Results and Discussion
4.1.Apparent Permeability Analysis in Shale Gas.In this section, the effect of stress sensitivity and adsorption on apparent permeability in the shale matrix is discussed.e basic parameters for the seepage model are shown in Table 3. e relationship between permeability correction coefficient and pressure under a different pore size without stress sensitivities is shown in Figure 11(a).It is can be clearly found that the correction coefficient of permeability is greater than 1 without considering the stress sensitivity and adsorption, which means the apparent permeability of shale is greater than the intrinsic permeability.Additionally, the correction coefficient of permeability decreases with the increase of pore size, and the reason is that the effect of stress and the increase of pore size reduce the Knudsen number; thus, the gas flow in porous media is in continuous flow, and apparent permeability and absolute permeability reach a consensus.
Besides, the relationship between permeability correction coefficient and pressure under a different pore size with stress sensitivity equal to 0.06 MPa −1 is shown in Figure 11(b).It reveals that the apparent permeability is affected by the stress sensitivity; specifically, when the pore diameter is greater than 5 nm, the permeability correction coefficient decreases first and then rises with the decrease of pressure.When the pore diameter is less than 5 nm, the permeability correction coefficient increases with the decrease of pressure.
ese results are in accordance with findings of earlier investigations [25].In a word, compared with the BK model, BKG model can characterize the multiscale seepage mechanism in shale gas when considering the effective stress variation [27,28].

Effect of Stress sensitivity.
In this section, the influence of stress sensitivity for gas transport in shale matrix and fracture is investigated thoroughly, and the basic parameters for the multistage fractured horizontal well numerical simulation are listed in Table 4. e parameters in Table 4 are derived from a horizontal well in Fuling Shale Gas Field, China.In order to increase the calculation efficiency, only three hydraulic fractures are considered.And the 3D shale gas reserve model with a multistage horizontal well is represented in Figure 12.In the basic case, the shale stress Journal of Chemistry sensitivity coe cient is 0, which means there is no stress sensitivity.In addition, we keep the parameters of the reservoir and multistage as a constant while change the stress sensitivity coe cient.
rough these scenarios, we can obtain the impacts of the stress sensitivity on the multistage horizontal well interstitial ow.
e simulation results of the behavior of multistage fractured horizontal well in a 500 days period for the four cases are shown in Figures 13(a) and 13(b).It is noticeable that the gas production rate and cumulative gas production will decrease when the shale gas reservoir stress sensitivity coe cient increases.
e reason is that increasing stress sensitivity coe cient means there are more deformation for the shale rock skeleton, and it can reduce shale gas apparent permeability.e simulation case without stress sensitivity has the biggest cumulative gas production of about 3011 × 10 4 m 3 after 500 days, while for the case of stress sensitivity coe cient 0.12 MPa −1 , the cumulative gas production is about 0.41 times of that value.In conclusion, the stress sensitivity has an important e ect on shale gas transport, and the deformation of the rock skeleton will a ect the behaviors of multistage horizontal well in shale gas.
Specially, we analyzed the multiscale seepage mechanism in shale gas with the case of stress sensitivity, namely, 0.03 MPa −1 , 0.05 MPa −1 , and 0.12 MPa −1 .Figure 14 represents the pressure distribution at di erent times with stress sensitivity coe cient equal to 0.05 MPa −1 .It is clearly seen that the shale matrix pressure evolution of the case of coe cient is 0.05 MPa −1 .
e process of shale gas transportation can be divided into three stages.Firstly, the gas in the fractured shale area ows into the wellbore.Secondly, the matrix gas ows into the hydraulic fracture.en, the gas within the kerogen will undergo desorption and di usion to the matrix pore.Figures 15(a ese results demonstrate that stress sensitivity coe cient can re ect the impact degree from e ective stress, and the shale gas matrix  6 Journal of Chemistry permeability will decrease after gas production, which leads to increase the seepage resistance and decrease the gas owing into wellbore.So, in order to avoid the e ect of stress sensitivity on shale gas production, reasonable production pressure di erence should be investigated considered particularly.

E ect of Adsorption Phenomenon.
In this study, we perform sensitivity analysis of the absorption volume.Four cases are proposed to investigate the impact of the absorption on shale gas multiscale transport.e Langmuir pressure is treated as a constant, while the Langmuir volume is changed.To be speci c, four cases of Langmuir volume, namely, 0 m 3 /kg, 0.004 m 3 /kg, 0.008 m 3 /kg, and 0.012 m 3 /kg, are analyzed.Other reservoir parameters are same, and the simulation results are shown in Figures 16(a) and 16(b); it can clearly be found that gas production rate and cumulative gas production increase with an increase in Langmuir volume.And the cumulative gas produced increases by 1.6% and 4.3% for the Langmuir pressure from 0.004 to 0.008 and 0.012, respectively.e reason is that an increase in the Langmuir volume means an increase in adsorbed gas reserves; moreover, when the pressure drops, the adsorbed gas is desorbed to free gas, which supplements the formation energy.In summary, the adsorption phenomenon has an in uence on shale gas seepage and sorption capacity, which is an important e ect to analyze the multistage fractured horizontal well production behavior; however, the e ect of adsorption is very weak in the early gas transport period, and the impact will increase later.

Impact of Shale Porosity.
ree separate simulation cases were validated with shale porosity of 0.01, 0.05, and 0.10.e of multistage fractured horizontal well behavior are presented in Figures 17(a) and 17(b).It can clearly be found that the signi cance di erence in gas production rate decline is apparent between the models.Furthermore, the cumulative gas production has been increased to 1.6 times with the porosity increased from 0.05 to 0.1.And the di erence of shale matrix porosity can a ect the Knudsen number of grids, as shown in Figures 18(a) and 18(b).Knudsen gradually increases with the decrease of pressure.With the continuous decrease of pressure, the ow pattern of gas transits from continuous ow to slippage and transition ow.Because the Knudsen number is the function of pressure, large pressure drop on the fracture leads to large number of Knudsen around the fracture.ese results highlight the importance of matrix porosity to shale gas productivity.e reason is that the larger the storage space the shale matrix has, the more the free gas it contains, which can delay the decline of formation energy.

Impact of Fracture Length and Conductivity.
In this section, we perform sensitivity analysis of the fracture length   and conductivity.e e ect of fracture length on gas production rate and cumulative gas production is shown in Figures 19(a

Conclusions
In this paper, a model of gas multiscale transportation-coupled stress sensitivity is proposed to investigate the shale gas seepage mechanism.e mathematical model is solved using the nite di erence method.According to the analysis of key parameters, the following conclusions can be deduced: (1) A new apparent permeability model considering multiscale ow and stress sensitivity was proposed and analyzed.Pore size and stress sensitivity have an important in uence on the permeability correction factor.e stress sensitivity can lead to the decrease of shale permeability.(2) e behavior of the multistage fractured horizontal well is simulated with a di erent stress sensitivity coe cient.e stress sensitivity actually has a signi cant impact on shale gas multiscale seepage.e gas production rate and cumulative gas production will decrease when the shale gas reservoir stress sensitivity coe cient increases, and the stronger the stress sensitivity, the lower the recovery factor of shale gas.(3) e Langmuir pressure and shale porosity signicantly a ect the behavior of the multistage fractured horizontal well.e cumulative shale gas increases as the Langmuir volume and porosity increase.However, the variation of porosity has more obvious e ect on gas production than Langmuir volume in the early period.(4) e Langmuir pressure and shale porosity signicantly a ect the behavior of the multistage fractured horizontal well.e cumulative shale gas increases as the Langmuir volume and porosity increase.And cumulative gas production increases with the increase of fracture length and conductivity.And the variation of porosity has more obvious e ect on gas production than Langmuir volume in the early period.

Figure 1 :
Figure 1: Physical model of gas transport.
) and15(b)  show the distribution of the permeability correction factor under di erent stress sensitivity coe cients after 500 days.

Figure 11 :
Figure 11: e relationship between permeability correction coe cient and pressure: (a) without stress sensitivity; (b) with stress sensitivity.

Figure 13 :Figure 14 Figure 12 :
Figure 13: e behavior of the multistage fractured horizontal well: (a) gas production rate; (b) cumulative gas production.
Figures20(a)  and20(b).Cumulative shale gas production increases with the length of fracture.It can be clearly found that cumulative gas production increases with the increase of fracture length and conductivity.

Figure 16 :
Figure 16: E ect of Langmuir volume on behavior of the multistage fractured horizontal well: (a) gas production rate; (b) cumulative gas production.

Table 1 :
Division of ow regimes according to the Knudsen number.

Table 2 :
Gas reservoir parameters for validation.

Table 3 :
Basic parameters for the seepage model.

Table 4 :
Main parameters for gas reservoir simulation.