Abstract
Volume fracturing is a key technology in developing unconventional gas reservoirs that contain nano/micron pores. Different fracture structures exert significantly different effects on shale gas production, and a fracture structure can be learned only in a later part of detection. On the basis of a multiscale gas seepage model considering diffusion, slippage, and desorption effects, a three-dimensional finite element algorithm is developed. Two finite element models for different fracture structures for a shale gas reservoir in the Sichuan Basin are established and studied under the condition of equal fracture volumes. One is a tree-like fracture, and the other is a lattice-like fracture. Their effects on the production of a fracture network structure are studied. Numerical results show that under the same condition of equal volumes, the production of the tree-like fracture is higher than that of the lattice-like fracture in the early development period because the angle between fracture branches and the flow direction plays an important role in the seepage of shale gas. In the middle and later periods, owing to a low flow rate, the production of the two structures is nearly similar. Finally, the lattice-like fracture model is regarded as an example to analyze the factors of shale properties that influence shale gas production. The analysis shows that gas production increases along with the diffusion coefficient and matrix permeability. The increase in permeability leads to a larger increase in production, but the decrease in permeability leads to a smaller decrease in production, indicating that the contribution of shale gas production is mainly fracture. The findings of this study can help better understand the influence of different shapes of fractures on the production in a shale gas reservoir.
1. Introduction
Shale gas as a typical unconventional resource is hosted in organic-rich shale reservoirs [1, 2]. A shale gas reservoir is tight, and its matrix is mainly nano-micro porous [3, 4]. A flow regime in shale gas reservoirs, which includes not only seepage but also diffusion, slippage, desorption, and absorption, is obviously different from that in conventional reservoirs [5–7]. Moreover, the gas flow in the shale matrix is multiscale and nonlinear. Many theoretical and experimental studies on these multiscale gas flow characteristics of shale reservoirs [8, 9] show that the velocity of shale gas in the matrix is so slow that fracturing development methods must be used to increase the production capacity of shale gas reservoirs [10–15].
Given a special shale bedding structure and fragility, volume fracturing is frequently used in the development of shale gas reservoirs, which causes natural cracks to expand into shear slips of brittle rock, forming a complex fracture network with natural and artificial fractures and increasing volume fractures [16–22] to improve the initial production and ultimate recovery. The traditional hydraulic fracturing which is based on the theory of two symmetrical fractures is inapplicable to volume fracturing networks with natural fractures, bedding, and anisotropic prominent shale gas reservoirs [23–25].
Researching complex fracture networks is very important to have accurate production evaluations for shale gas reservoirs [26, 27]. The two shapes of fracture structures, tree-like and lattice-like, have been widely exploited. For shale reservoirs with a tree-like fracture structure, due to difficulty to express the fractures, an equivalent continuous medium model is always used to deal with a fracturing zone through the equivalent permeability according to the fractal theory [28]. For shale reservoirs with a lattice-like fracture network structure, except the treatment of an equivalent continuous medium, a dual-medium model or a discrete fracture network model is also frequently used to simulate seepage characteristics. Eshkalak et al. [29] established a dual-medium model for shale gas with adsorption-desorption and studied the production of single wells through the finite difference method. By using the finite difference method, Swami and Settari [30] determined pressure characteristics and built a production model for gas wells in consideration of the effects of slippage and permeability changes in a fracturing area. Kim et al. [31] considered stress sensitivity and constructed a porous-flow finite element model for shale gas reservoirs with flow-solid coupling. These studies focused only on 2D problems, and no research examined 3D problems. Besides, effects of different fracture network structures on the production of shale gas reservoirs have not been addressed.
The complicated percolation mechanism such as desorption, diffusion, and percolation of shale gas in nanometer porous medium has been studied extensively by scholars, due to the brittleness of shale rock, fracturing is not controllable, and there are just a few laboratory studies on network morphology and seepage law. Zhu et al. [32] established a productivity prediction model of horizontal well multistage fracturing in the shale reservoir considering the interference between the multistage fracturing zones and the pressure drop in the horizontal wellbore. Lu et al. [33] established orthogonal fracture network models with different structures, conducted a series of single-phase and two-phase flow experiments, and studied the influence of seepage characteristics and model structure on permeability in the orthogonal fracture network model. At present, there is a lack of different fracture network morphology description and parameter characterization methods, which cannot better describe or predict the fracture network complexity (density, fracture range, etc.) formed by hydraulic fracturing. Therefore, it is urgent to establish a mathematical model suitable for shale storage characteristics and reveal the gas seepage law under different fracture network conditions.
Numerical methods such as the finite difference method and the finite element method are the powerful methods to study shale gas reservoirs with various fractures [28–36]. First, the finite element algorithm is presented to study multiscale gas flow in a shale reservoir considering desorption, diffusion, and slip effects. Then, the models for the gas reservoir with tree-like and lattice-like fracture structures are established, respectively, under the condition of equal fracture volumes. The flow characteristics of shale gas in a fracturing zone and the effects of fracture network structures on shale gas production are investigated and analyzed. Finally, the factors of shale properties that influence shale gas production in the development process are also researched through the netted fracture models.
2. Finite Element Model for Seepage of Shale Gas
For a shale reservoir with micron pores, the seepage of gas can be expressed by the following formula [37]: where where is the intrinsic permeability, denotes the degree of anisotropy of the matrix when the reservoir is isotropic, , , represents the slippage coefficient, is the Knudsen diffusion coefficient, denotes the density of gas, denotes the effective porosity that is assumed unchangeable in this paper, denotes the isothermal compressibility of gas, denotes the pore pressure, denotes the Langmuir pressure, denotes the Langmuir volume, represents the viscosity coefficient, where the subscript “sc” denotes the variable under standard conditions, denotes the temperature, and denotes the gas compressibility factor, which represents the deviation degree of real and ideal gases under similar conditions.
The finite element method is used to solve Equation (1), and the partial derivative in time in this equation is dealt with an implicit first-order backward difference method. The finite element algorithm for Equation (1) is deduced as where where denotes the numbering of nodes in each element , denotes the shape function at an internal node of an element, denotes the total number of nodes in an element, denotes the Kronecker delta tensor, is the unit outward normal vector at any point of the boundary , denotes the iterative step of time (), and denotes the iterative step for solving the seepage field at the time level ().
The known pressure boundary (Dirichlet boundary) or the known flux boundary (Neumann boundary) is considered for the initial boundary condition for the nonlinear flow of shale gas. (1)The known pressure boundary, in which the pressure value at a node of the boundary is known, does not require an iterative solution. At the boundary, it can be expressed as(2)The known flux boundary can be expressed aswhere denotes the flow velocity on the boundary.
The relationship between pressure and gas viscosity and compression factor can be considered in this finite model. In this model, we use a data table to describe the relationship between pressure and gas viscosity and gas compression factor. But in this study, their relationship is not very important to discuss, so we set up them as constant.
3. Shale Gas Production with Different Fracture Structures
Shale gas is often exploited through horizontal wells and volume fracturing. A schematic of a fracture distribution in shale gas horizontal wells is shown in Figure 1. Many scholars believed that the form of fractures significantly influences the production of shale gas [9]. Lattice-like and tree-like fractures are two typical structures [9]. In this section, of a fracturing area is regarded as the research object. The gas flows in the lattice-like and tree-like fractures are investigated below, respectively.
Figures 2(a) and 2(b) show the geometric models for lattice-like and tree-like fracturing regions, respectively. The boundaries of these two regions give equal fracture volumes. The width of the models is 40 m, their height is 300 m, the top of the fracture network elevation is 100 m, and the fracture width is 0.01 m. Lines A1B1 and A2B2 indicate the outer boundaries, and the pressure values are all 24 MPa; lines C1D1 and C2D2 indicate the inner boundaries, and the pressure values are 6 MPa; lines A1D1 and B1C1 and lines A2D2 and B2C2 are symmetrical boundaries separately. The thickness of the models is set to 1 m. The parameters of the shale matrix are given as follows:
(a) Lattice-like fractures
(b) Tree-like fractures
The parameters of the fractures which are under the assumption that there is no desorption are as follows:
The parameters of shale gas in the matrix are set as
The finite element models with the lattice-like and tree-like fractures are, respectively, established, and the meshes of the local areas of the fractures are shown in Figure 3. The finite element algorithm in Section 2 is used to calculate the pressure and gas production during the production process.
(a) Lattice-like fractures
(b) Tree-like fractures
Figures 4(a) and 4(b) show the pressure distributions in the lattice-like and tree-like fractures after 10 days of reservoir development, respectively. The internal pressure in the fracture networks releases in both the lattice-like and tree-like fracture models. The change of pressure outside the fracture networks is very small. It can be seen that in the initial production time, the most change of pressure occurs in the fracture networks.
(a) Lattice-like fractures
(b) Tree-like fractures
We take line E-E in Figure 4(a) (which is at the plane) as an example to analyze the pressure distribution. The pressure distributions along this line at different times (2, 10, 20, and 30 d) in the lattice-like fracture model are shown in Figure 5. It is clear that the pressure inside the fracture network is higher than that in the fractures. The pressure in the internal region of the network is gradually released with the passage of time. The initial release is rapid, and the later release is gradual. The pressure release rate outside the network is much lower than that in the fracture network.
Figure 6 presents a pressure distribution diagram after 600 days of development. The pressure field distribution changes in the peripheral fracture network are increased, and the distribution of pressure in the lattice-like and tree-like fractures is very similar.
(a) Lattice-like fractures
(b) Tree-like fractures
Figures 7 show the distribution maps of single-day mass fluxes with time for the lattice-like and tree-like fracture models, respectively. The changing trend of the flow rate for the lattice-like fracture structure is consistent with the flow rate for the tree-like fracture structure . A clear difference in specific values occurs from 0 to 100 days. With an increase in time, the difference between the two structures decreases. Figure 8 shows an increasing ratio of the flow rate with time for the tree-like fractures to the lattice-like fractures . During early exploitation time, the tree-like fractures’ daily production is larger than that of the lattice-like fractures. The largest difference between these two structures occurs on the 20th day. After 300 days, the daily productions of the two fracture structures are almost equal. In the case of the same external boundary conditions and the same area, the structure of fracture networks significantly affects the output of the initial production stage but exerts a minimal effect on the later period. Figure 9 presents a cumulative comparison of production for the lattice-like and tree-like fractures. The cumulative production of the tree-like fractures is larger than that of the lattice-like fractures. The difference between the two cumulative productions becomes steady after 300 days, which is in accordance with the results shown in Figures 7 and 8.
With the model of the lattice-like fractures as an example, the effects of the permeability of fractures , the intrinsic permeability of the matrix , the diffusion coefficient of the matrix , and the Langmuir volume on the flow are investigated.
Figure 10 shows the pressure distributions for , , , and with different values on the cross section (line E-E in Figure 4(a)) after 10 days of development for the lattice-like fracture model. , , and exert significant influences on the pressure distribution in the vicinity of a fracture network. exerts a minimal influence on the pressure in the vicinity of the fractures.
(a)
(b)
(c)
(d)
Figures 11 and 12 show the production changes with development time under different values of , , , and . For , the greater the , the higher the daily production in the early period of development, but in the middle and later periods, a production decline is larger. Figures 11(a) and 12(a) also show that the production will not become larger along with the value of . Thus, the determination of a value of in a hydraulic fracturing design needs to be combined with reservoir conditions and to be optimized. For and , the daily and cumulative productions become higher along with the value of and . For , in the early stage of development, due to the existence of fractures, the gas flow in the fracture networks runs very fast, and the influence of in the matrix on the production is very small; at a later time, because the gas in the fracture networks releases out, the production mainly comes from the desorption and diffusion effects in the matrix, and the influence of on the production gradually becomes obvious.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
4. Summary and Conclusions
Considering the multiscale nonlinear gas seepage theory considering desorption, diffusion, and slip effects of shale, a 3D nonlinear finite element algorithm for shale gas seepage is developed. Finite element models for two different fracture network structures are established under the condition of equal fracture volumes, and the effects on the production of fracture network structures are compared and analyzed. (1)The early tree-like fracture production is higher than the lattice-like fracture production because the angle between the direction of fracture branches and the flow direction is less than 90°. This angle plays an important role in the seepage of shale gas. Part of the lattice-like fractures perpendicular to the flow direction exerts a minimal effect on the permeability of shale gas reservoirs(2)In the middle and later periods, owing to a low flow rate in the matrix, the difference between the two fracture network structures can be disregarded, and the productions of these two structures are nearly the same(3)With the lattice-like fracture model as an example, the influencing factors (, , , and ) for the shale gas production in the development process are analyzed. The Langmuir volume, diffusion coefficient, and matrix permeability are positively correlated with shale gas production, and the matrix permeability exerts a significant influence on production, followed by the gas diffusion coefficient. The production does not become larger along with the value of , and the determination of a value of in a hydraulic fracturing design needs to be combined with the reservoir conditions and to be optimized(4)The accuracy and efficiency of our numerical method have been shown. Therefore, this method can be used to solve more complicated boundary value problems of shale gas seepage and can be extended to more sophisticated analysis for shale gas reservoirs exploited by fracturing
Nomenclature
| : | Intrinsic permeability |
| : | Degree of anisotropy of the matrix when the reservoir is isotropic, |
| : | Slippage coefficient |
| : | Knudsen diffusion coefficient |
| : | Density of gas |
| : | Effective porosity that is assumed unchangeable in this paper |
| : | Isothermal compressibility of gas |
| : | Pore pressure |
| : | Langmuir pressure |
| : | Langmuir volume |
| : | Viscosity coefficient, where the subscript “sc” denotes the variable under standard conditions |
| : | Temperature |
| : | Gas compressibility factor |
| : | Numbering of nodes in each element |
| : | Shape function at an internal node of an element |
| : | Total number of nodes in an element |
| : | Kronecker delta tensor |
| : | Unit outward normal vector at any point of the boundary |
| : | Iterative step of time () |
| m: | Iterative step for solving the seepage field at the time level () |
| : | Flow velocity on the boundary. |
Data Availability
The data used to support the study is available within the article.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
The work is supported by the National Natural Science Foundation of China (Nos. 51704015 and 11602221) and China Scholarship Council (No. 201706465064).