## Abstract

This paper presents the theory and application to modify the conventional simulator to describe the effects of gas adsorption and gas slippage flow in shale gas. Because of the local desorption of gas and the assumptions of gas desorption instantaneously with the decrease in pore pressure, we define one fictitious immobile “pseudo” oil with dissolved gas. The dissolved gas–oil ratio is calculated from the Langmuir adsorption isotherm constants and shale gas properties. Additional modifications required in the input data are the porosity and relative permeability curves to account for the existence of “pseudo” oil.

#### Authors:

##### Yongmao Hao^{1}, Wendong Wang^{1}, Bin Yuan^{2}, Yuliang Su^{1}, Jie An^{3}, Huai Shu^{4}

^{1}China University of Petroleum (East), No. 66, Changjiang. ^{2}University of Oklahoma, Norman, OK 73019, USA. ^{3}Oil & Gas Technology Research Institute of Changqing Oilfield Company, PetroChina, Xi’an 710018, China. ^{4}Hancheng Branch of PetroChina Coalbed Methane Co., Ltd., PetroChina, Hancheng 715400, China.

Received: 19 January 2017 / Accepted: 17 June 2017

© The Author(s) 2017.

The input rock table considers the changes of rock permeability versus pressure to describe the gas slippage flow effects. In addition, dual-porosity dual-permeability models coupled with local grid refinement method are used to distinguish the impacts of natural fractures and hydraulic fractures on shale gas production with the comparison of vertical well, fractured vertical well, horizontal well, and multistage fractured horizontal well production. This proposed simulation approach shows enough accuracy and outstanding time efficiency. Results show that ignoring gas desorption and slippage flow effects would bring significant error in shale gas simulation The existence of natural fractures also imposes great effects on the productivity of shale gas.

## Introduction

Recently, tight gas and shale gas have achieved great attention because of advancements of horizontal well drilling and large-scale fracturing technology. Recoverable reserves of shale gas in USA are estimated to be 862 TCF (Kuuskraa et al. 2011; Yuan and Wood 2015a, b). One of the issues yet not solved is the lack of predicative tools for shale gas production. The complex shale gas flow mechanisms in nanopores cannot be simply described as Darcy flow equation. In addition, due to large-scale fracturing in shale gas reservoirs with natural fractures, the conventional single porosity model is not enough to simulate the characteristics of shale gas reservoirs. Furthermore, as advanced simulation methods summarized by (Wood et al. 2015), rigorous simulation techniques, such as direct simulation, Monte Carlo method, and Molecular dynamic simulation, are very computationally expensive (Moghanloo et al. 2015; Yuan et al. 2015). Hence, there are two alternative approaches to model shale gas: (1) application of analytical methods to characterize the primary characteristics of shale gas (Ghanbarnezhad Moghanloo and Javadpour 2014), and (2) extending the conventional simulator to effectively model flow from shale gas reservoirs. In this paper, we present the theory and application to modify the conventional simulator to describe the effects of Langmuir gas desorption and gas slippage flow effects (Shabro et al. 2011; Waqas 2012). We also present the dual-porosity dual-permeability model coupled with local grid refinement method to distinguish the impacts of natural fractures and hydraulic fractures on shale gas production.

As for the flow regimes of shale gas, Knudsen number (*K _{n}*), the ratio of molecular mean free path to characteristic length, is commonly used. The Knudsen number can also be used to check the validity of continuum assumption:

- When
*K*< 0.001, the Navier–Stokes equation with no-slip boundary conditions are presented to describe the fluid flow. This regime is called continuum flow regime._{n} - When 0.001 <
*K*< 10, the Navier–Stokes equation will not fully perform well and the flow regimes of rarefied gas, which is neither continuum flow nor free-molecular flow._{n} - When 0.001 <
*K*< 0.1, gas molecules tend to slip at the surface of rock (slippage condition), and as a result, gas flux will increase. However, the slippage is just limited to the vicinity adjacent to pore surfaces. Knudsen layer are defined between the fluid and pore surface with the order of one molecular mean free path. Therefore, this Knudsen layer can be neglected compared with the size of pore throat, flow velocity is linearly extrapolated toward to the pore surface. The finite nonzero flow velocity exists on the wall of pore throats. The shale gas flow will deviate from Darcy flow equation with the modification of Klinkenberg’s modification._{n} - When 0.1 <
*K*< 10, the Knudsen layer will increase. The Navier–Stokes equation will not applicable. Furthermore, Knudsen diffusion increases because the molecular mean free path becomes comparable with the size of pore throat; in other words, the molecule/wall collisions will dominate particle/particle collisions. (Fathi et al. 2012) introduced a correction of double slip to Klinkenberg’s slippage model._{n} - For
*K*> 10, due to the occurrence of free molecule flow, kinetic effects become dominate. The collisions among gas molecules are smaller that of molecule/wall collision. Continuum equations becomes fully ineffective, and the ballistic flow happens (Hadjiconstantinou 2006)._{n}

In this paper, the gas slippage (Klinkenberg effect) will be described with the commonly used Klinkenberg gas slippage factor within the Darcy flow equation. In addition, in tight gas reservoirs, there are significant incidences of natural complex mineralized or plugged fractures. Fracturing of horizontal wells with hybrid fluid or slick water (SWF) (Britt and Smith 2009) not only increases fracture formation contact areas by creating SRVs, but also reactivates natural fractures to some extent, which is vital to tight gas production (Ozkan et al. 2009; Clarkson 2013).

The generation of large-scale and dense SRVs is enhanced by large density of preexisting natural fracture in tight reservoir. Complexity of fracture networks (Fig. 1) can increase the fracture–reservoir contact area due to high fracture density and SRV size.

**Fig. 1** Microseismic fracture map with complex fracture network in shales (Warpinski et al. 2008).

In this paper, we present conventional modified simulation approach to simulate the performance naturally fractured shale gas reservoirs. The main advantage of our work is to provide an equivalent approach to simulate shale gas in conventional black oil model by considering gas desorption, gas flow slippage, all of which are the particular characteristics of shale gas different from conventional gas reservoirs.

The comparison with previous conventional methods indicates the accuracy of our proposed equivalent methods. This approach can also cover the primary geological characteristics of natural fractured shale gas, and the difference of the natural fractures around the well nearby and those hydraulic fractures.

Based on the theory of fluid flow in naturally fractured reservoirs developed by Barenblatt et al. (1960), Warren and Root (1963) presented the concept of dual-porosity models to describe the matrix as the source of interconnected fractures. Kazemi et al. (1976) introduced the dual-porosity model into numerical model to describe fluid flow on large scale. During the large-scale hydraulic fracturing process, there is reactivation of preexisting (cemented) natural fractures and activation of planes of weakness, and even connected with hydraulic primary and secondary fractures, which is so-called fracture networks or stimulated reservoir volume.

Based on the concept of dual-porosity model, it is applicable to introduce the dual-porosity model to describe the characteristics of naturally fractured tight oil reservoirs. The increase in shape factor means the enhancement of fracture density of natural fractures, as shown in Eq. 1.

where, S* _{d}* is shape factor reflecting the geometry of matrix grids, and controlling the flow between matrix and fracture, m

^{−2}. Lx ; Ly ; Lz are the average size of cleaved matrix blocks.

## Mathematical model and simulation approximation

### Gas desorption

There are mainly two steps for shale gas flow, firstly, gas desorption from the pore throat surfaces, and then flow through the fractures and matrix. The slower rate of one of those above two steps determines the rate of shale gas production. If the desorption of shale gas is very fast compared with the rate of flow in shale reservoirs, the shale gas production will be dominantly controlled by the flow in reservoirs but not the desorption process. Based on this assumption, we define one fictitious immobile “pseudo” oil with dissolved gas. The dissolved gas–oil ratio is calculated from the Langmuir adsorption isotherm constants and shale gas properties.

Additional modifications required in the input data are the porosity and relative permeability curves to account for the existence of “pseudo” oil. The amount of gas adsorption per unit mass of rock at the pore throat surface at a given pressure is analogous to the amount of gas dissolved in an “equivalent” unit volume of black “pseudo-” oil at the given pressure. The Langmuir isotherm of adsorption will become the dissolved gas–oil ratio of new defined immobile oil. However, the introduction of “pseudo” oil will require the increase in porosity of the modeled shale gas reservoirs and the changes of phase saturation. The gas–water–oil relative permeability curve will also be modified and the immobile oil property is defined. The governing equations for two phases (oil and gas) in dual-porosity and dual-permeability black oil model have been added as below, which are the basic models for ELCIPSE 100. However, in this paper, we modify the definitions of parameters in conventional black oil model to mimic the functions of composition simulator for shale gas.

1-Matrix System:

Oil Phase:

Gas Phase:

2-Fracture System:

Oil Phase:

Gas Phase:

where the subscript *m* denotes the simulation model.

The amounts of gas in the simulation model and actual shale gas reservoirs must be equal to meet with the material balance (Seidle and Arri 1990), and therefore, the saturation is related to each other by:

where the subscript *m* denotes the simulation model.

The sum of phase saturation must be unity, and therefore,

In addition, in the real shale gas reservoirs, assuming Φ is hydrocarbon filled porosity:

Equations 2a and 2b are simplified to be

This equation relates the shale gas reservoirs porosity with porosity value in the simulation model. The immobile oil saturation can be arbitrary, but once chosen, it must be constant throughout the simulation process. The gas saturation in the model becomes:

Because we have defined oil phase, we have to define the relative permeability curve to describe the flow of gas in the shale gas model. The relative permeability corresponding to the actual gas saturation must be altered to the new gas saturation in model. Actually, in shale gas reservoir, there does not exists mobile oil, and therefore, the immobile oil defined in our model is characterized by very small relative permeability (close to zero) and very large viscosity, in the paper, 10^{7} cp.

In the physical shale gas reservoirs, the amount of gas adsorbed in a unite volume is defined as:

where *V* is the gas content, SCF/ton, *V*_{L} is the Langmuir volume in *SCF/ton*, *b* is the Langmuir pressure constant in psia^{−1}, and *p* is the pressure in psia.

where *P*_{L} is Langmuir pressure in psia.

In the simulation model, V_{m} , the amount of gas dissolved in the same unit volume is:

where *Rs* is the solution gas–oil ratio in MCF/STB; *B*o is the oil formation volume factor in RB/STB. To conserve the mass of phase, the *B*o must be constant; here, we define as 1.0; Combining Eqs. 7 and 9, to give the relation of dissolved gas in an immobile oil.

As noted, the immobile oil saturation can be arbitrary, but the corresponding dissolved gas depends on the choice of defined immobile oil saturation to maintain the material balance of shale gas. In the example of this paper, those above parameters are defined as Tables 1 and 2.

**Table 1** Values of parameters in shale gas simulation.

**Table 2** Equivalent dissolved gas–oil ratio with amounts of gas adsorption (PVTO input values).

The other shale gas parameters are chosen as the data provided in the literature review (Kale et al. 2010).

To demonstrate the accuracy of this equivalent simulation method to modify conventional black oil simulator, the comparison of production rate between the cases using conventional modified model and compositional model is shown as Fig. 2. As shown, the gas production calculated from both models is in good agreement; however, the new proposed method has more advantage of calculation efficiency, which verifies the applicability of this new simulation approach. It is noticed that the shale gas reservoir simulation with new method is processed in short time and need very little data. Furthermore, the modified method from the black oil model is reliable, efficient, and applicable to a wide range of reservoirs, including shale gas and coal-bed methane.

**Fig. 2** Comparison shale gas production rate with respect different simulator.