Abstract

Currently, there is no proper method to predict the pore pressure disturbance caused by multistage fracturing in shale gas, which has challenged drilling engineering in practice, especially for the infilling well drilling within/near the fractured zones. A numerical modelling method of pore pressure redistribution around the multistage fractured horizontal wellbore was put forward based on the theory of fluid transportation in porous media. The fracture network of each stage was represented by an elliptical zone with high permeability. Five stages of fracturing were modelled simultaneously to consider the interactions among fractures. The effects of formation permeability, fracturing fluid viscosity, and pressure within the fractures on the pore pressure disturbance were numerically investigated. Modelling results indicated that the pore pressure disturbance zone expands as the permeability and/or the differential pressure increases, while it decreases when the viscosity of the fracturing fluid increases. The pore pressure disturbance level becomes weaker from the fracture tip to the far field along the main-fracture propagation direction. The pore pressure disturbance contours obviously have larger slopes with the variation of permeability than those of the differential pressure. The distances between the pore pressure disturbance contours are smaller at lower permeability and higher viscosity. The modelling results of the updated pore pressure distribution are of great importance for safe drilling. A case study of three wells within one platform showed that the modelling method could provide a reliable estimation of the pore pressure disturbance area caused by multistage fracturing.

1. Introduction

Shale gas becomes more and more important worldwide. The shale gas production reached  m3 in America in 2016, which made up more than 40% of the total natural gas production of America [1]. The horizontal well factory and multistage hydraulic fracturing are the two most important technologies for shale gas production commercially [24].

There is over sixty years of history in the study of fracture propagation in porous media, and many models have been developed, such as the PKN model [5, 6], KGD model [7, 8], and some three-dimensional models [9, 10]. Recently, more complex models for multiple-fracture propagation in horizontal wells have been built, such as the Unconventional Fracture Model (UFM) [11, 12]. In order to obtain complex fracture networks in low permeable shale reservoirs through multistage hydraulic fracturing, many studies have been done to analyze the stress distribution around the fractured horizontal wellbore [13], which is mainly focused on the extent of stress reversal [14] and the effect of stress shadow [15]. Besides the analysis of fracture initiation and propagation, the stress redistribution results are further used to optimize the stage spacing [14, 16, 17] and well space [18]. Such kinds of works have greatly helped to improve the efficiency of stimulation treatments in shale gas formation.

In fact, hydraulic fracturing in the pay zone also has a pressure-elevating effect. The pore pressure will increase during and after the fracturing work. The natural gas within the pay zone of the fracturing well would be driven to the nearby well. When the gas invades the nearby wellbore under drilling, it would lead to a gas kick and overflow. It seems easy to understand the pore pressure elevation caused by hydraulic fracturing, but it is not easy to quantitatively describe these updated pore pressure distributions. Currently, there is no proper method to predict or detect such pore pressure disturbances from the industry, and it has not yet gained enough concerns to warrant a scientific study. [19]. Unfortunately, it has challenged the drilling engineering in practice, especially for the infilling well drilling within/near the fractured zones. The infilling well drilling practice in Fuling, the first and largest commercially developed shale gas field in China, indicated that some well drilling works were obviously influenced by troubles related to the large-scale hydraulic fracturing work [20], such as gas kick, overflow, and drilling fluid being polluted by fracturing fluid from neighbouring wells. The resultant percentage of nonproductive time could reach as high as 21.04~52.82% according to a simple summary of six infilling wells.

Pore pressure distribution around a single fracture could be approximately described as an ellipse (Koning, 1985; [21, 22]). But for the multistaged fracturing in shale gas, pore pressure distribution around several fractures should be evaluated simultaneously considering the interactions among fractures. In this study, firstly, a numerical modelling method of pore pressure distribution around a multistage fractured horizontal wellbore was put forward based on the theory of fluid transportation in porous media. Then, a series of numerical modelling experiments were carried out to investigate the effects of formation permeability, fracturing fluid viscosity, and pressure within the fractures on the pore pressure distribution. Finally, the implication of the safe drilling of an infilling well was further discussed. And a case study of infilling well drilling with the aid of numerical modelling showed that the average ROP reached 11.02 m/h and no drilling troubles or accidents occurred.

2. Methods

2.1. Model Description

The horizontal section for shale gas production is usually around 1000 meters. Multistage fracturing can produce a fracture network including main fractures and different scales of subfractures (Figure 1(a)). The length of the fracture network of each stage can be 200~400 meters according to the interpretation of microseismic monitoring. In order to investigate the pore pressure disturbance at the scale of the horizontal well, the modelled size is  m in length and  m in width, representing the half horizontal plane along the horizontal wellbore section (Figure 1(b)). According to the wire-mesh model by Xu et al. [23], the fracture network of each stage can be represented by an elliptical zone. Here, we use an ellipse (semimajor axis  m, short axis  m) with very high permeability of  mD to represent the fracture network (as shown in Figure 1(c)) compared to that of the low permeability of the shale formation ( mD). In this study, we simulate five stages of fracturing. The center points of these five half-ellipses from left to right are (440, 0), (470, 0), (500, 0), (530, 0), and (560, 0). The space between two nearby ellipses is 30 meters equal to the fracturing stage space.

For the transportation of any single phase (fluid or gas) within the porous formation, it is controlled by the continuity equation [24] according to the law of conservation of mass: where is the density of the fluid, is the porosity, and is the percolation velocity. If the fluid percolation follows Darcy’s law, the motion equation should be [25] where is the viscosity of the fluid, is the permeability of the formation, and is the gradient of pressure. Combining equations (1) and (2) gets

For a horizontal well, the influence of temperature on the density of fluid can be neglected. So the density variation with pressure can be expressed by where is the density at the reference pressure () and is the compressibility of the fluid. Because the compressibility of the shale rock is much smaller than that of fluid, the decrease of porosity with pressure is very small, so the right side of equation (3) can be further derived as

Therefore, equation (3) can be written as

For a certain investigated domain () as shown in Figure 1(b), the initial pressure () is equal to the shale gas reservoir pressure which is taken as the reference pressure:

For the pressure boundary of in Figure 1(c), the hydraulic fracturing pressure () is applied:

So the differential pressure () is approximately the power on the fracture surface driving the seepage of fracturing fluid into the shale formation, which basically results in the increase of pore pressure () near the fracture.

For other outer boundaries of the model (Figure 1(b)), they are defined to be zero flow because these boundaries are at the far-field condition:

The numerical calculation for the definite problem described by equations (6)–(9) within the domain shown in Figure 1(b) is carried out by FEM software. The multifrontal massively parallel sparse (MUMPS) method is used to solve the equations with speed-up convergence. Triangular grids are used to mesh the model. In order to ensure modelling convergence and acceptable accuracy, the maximum grid size is 1/(2-3) of the elliptical minor axis. Generally, the number of grids in a model is about 104~105 (as shown in Figure 2). After the pressure field at every time step is calculated, the pressure variation on the output lines (see definitions in Data Processing) can be extracted for detailed analysis.

2.2. Data Processing

The initial reservoir pressure is taken as a reference pressure to define the pore pressure disturbance as . Therefore, the pressure disturbance () is dimensionless and represents the pore pressure elevation due to hydraulic fracturing. In order to quantitatively analyze the pore pressure disturbance, three data output lines are set as shown in Figure 1(b).

(1) Data Output Line 1: y = 0. This line starts from the origin point of the model and extends to the right boundary of the model along the positive x-axis direction. The pore pressure distribution along this line can reflect the variation of pore pressure near the horizontal wellbore

(2) Data Output Line 2: y = 200. This line starts from the point (0, 200) on the left boundary of the model and extends to the right boundary along the positive -axis direction. The pore pressure distribution along this line can reflect the variation of pore pressure near the fracture tip in the -axis direction

(3) Data Output Line 3: x = 440 and y ≥ 200. This line starts from the tip (440, 200) of the left fracture and extends to the top boundary along the positive -axis direction. The pore pressure distribution along this line can reflect the variation of pore pressure in the -axis direction

For each modelling, the pore pressure distributions along these three lines are extracted, and then the corresponding pore pressure disturbances are calculated for further analysis.

2.3. Design of the Numerical Modelling Experiments

To investigate the effects of different factors of permeability (), viscosity of the fracturing fluid (), and differential pressure (, via changes in the hydraulic fracturing pressure at a constant initial reservoir pressure) on the pore pressure disturbance induced by fracturing, three series of modelling experiments are designed as detailed in Table 1. For each series, only one of the controlling parameters of permeability, viscosity of the fracturing fluid, or differential pressure is varied while the others are held constant. Therefore, any resultant change in overall pore pressure disturbance should be directly caused by the change of the varied parameter. Other modelling parameters for all these modellings are shown in Table 2 unless otherwise stated.

3. Results and Discussions

With the range of parameters shown in Tables 1 and 2, the pressure distribution patterns of the modellings are similar. Figure 3 shows the modelling results for one of the realizations of Table 1. The affected zone of pore pressure disturbance increases with the duration of differential pressure. In the following, the pore pressure distributions on the three output lines at different times are analyzed quantitatively for each series of modelling (Table 1).

3.1. The Influence of Permeability

To investigate the influence of permeability on the pore pressure disturbance, a series of modellings (series 1 in Table 1) with different permeabilities varying from 0.2 mD to 2.0 mD are carried out for the constant fracturing pressure and viscosity ( MPa and  mPa·s, respectively). The results are shown in Figure 4.

The results shown in Figure 4 indicate that the pore pressure disturbance zone () expands obviously from the fractures in both the and directions if the permeability of the shale formation increases. Although the permeability of a shale gas reservoir is usually very low, the pore pressure disturbance due to such large-scale, multistage hydraulic fracturing is significant. According to our modelling at  mD, the pore pressure disturbance zone with 5% increase () could be reached as far as 40.3 m in the direction and 48.4 m in the direction from the fracture.

3.2. The Influence of Viscosity

To understand the influence of viscosity on the pore pressure disturbance, a series of modellings (series 1 in Table 1) with varying viscosities (1 mPa·s, 2 mPa·s, and 3 mPa·s) are carried out when the permeability and the differential pressure are held constant ( mD,  MPa). The results are shown in Figure 5.

The pore pressure disturbance zone decreases when the viscosity of the fracturing fluid increases (Figure 5). This is because the increase of viscosity makes the fluid transportation in the shale much more difficult. Slick water with relatively low viscosity (usually <10 mPa·s) is widely used for the fracturing in shale gas reservoirs. According to our modelling at  mPa·s, the pore pressure disturbance zone with 5% increase () could be reached at 24.6 m in the direction and 27.3 m in the direction from the fracture. These results indicate that the influence of viscosity on the pore pressure disturbance is much less than that of formation permeability.

3.3. The Influence of Differential Pressure

In order to study the influence of differential pressure on the pore pressure disturbance, the formation pressure ( MPa), permeability ( mD), and viscosity of the fracturing fluid ( mPa·s) are kept constant and the hydraulic fracturing pressure is set to be 90 MPa, 100 MPa, and 110 MPa. So the differential pressure varies at 60 MPa, 70 MPa, and 80 MPa, respectively. The modelling results are shown in Figure 6.

As shown in Figure 6, the pore pressure disturbance zone expands as the differential pressure increases. But the increasing rate is relatively low. According to our modelling at  mD and  mPa·s, when the differential pressure increases from 70 MPa to 80 MPa, the pore pressure disturbance zone with 5% increase () only increases 1.19 m and 1.41 m in the and directions, respectively. This means, for the pressure level usually applied in the fracturing of shale gas formation, the influence of differential pressure on further expansion of the pore pressure disturbance zone is very small, which can almost be ignored compared to the larger scale of the fracture network zone (semimajor axis  m).

3.4. Implication of the Safe Drilling of the Infilling Well

For the implementation of the “well factory,” it is important to design a safe well spacing to avoid well interference during drilling and fracturing, especially for infilling well drilling. Based on the geometric positions of wells and fractures, the pore pressure disturbance along the main-fracture propagation direction should be carefully considered. Generally, the effective fracture zones are expected to connect to one another between the two closest nearby stimulation stages from the two sides of the two parallel wells [18], based on the view of gas production. However, such a connection between the effective fracture zones should be controlled to a certain extent to avoid troubles and events during infilling well drilling, such as gas kick and lost circulation, which requires a safe well spacing.

For example, there are two parallel horizontal wells (Well A and Well B) within the same plane and both of them have been hydraulically fractured (as shown in Figure 7) at the same condition that the pore pressure disturbance zones are the same. Now an infilling well (Well C) is planned to drill to accelerate shale gas recovery. If the density of the drilling fluid is designed according to the initial reservoir pressure and the ability of well control could deal with pressure disturbances lower than , then Well C would probably undergo gas kick in sections where the pore pressure disturbance zones of Well A and Well B overlapped (). To drill Well C safely, the density of the drilling fluid should be optimized considering the pore pressure disturbance. Therefore, it is important to understand the pore pressure disturbance along the main-fracture propagation direction.

In order to quantitatively investigate the pore pressure disturbance along the main-fracture propagation direction ( direction in this paper), the modelling results on data output line 3 are further analyzed (Figure 8). As shown in Figure 8, the pore pressure disturbance level becomes weaker from the fracture tip to the far field along the main-fracture propagation direction. The pore pressure disturbance contours obviously have larger slopes with the variation of permeability than those of the differential pressure (Figures 8(a) and 8(b)), which indicates that the distance of the pore pressure disturbance in the direction is more sensitive to the change of permeability. The pore pressure disturbance contours first decrease quickly at relatively low viscosity and then go down slowly to a plateau with further increase of the viscosity (Figure 8(c)). The distances between the pore pressure disturbance contours are smaller at lower permeability and higher viscosity (Figures 8(a) and 8(c)), which means that the gradient of pore pressure along the direction is large.

To analyze the sensitivity of pore pressure to different parameters in a universal way, we define two dimensionless parameters along the main-fracture propagation direction: (1)Dimensionless pressure (): (2)Dimensionless distance (): where is the half-length of the fracture, is the seepage distance in the direction from the fracture tip, and is the seepage velocity.

Our modelling results of series 3 (see detailed parameters in Table 1) show that the ~ curves for different differential pressures completely overlap with each other. It means the relationship between and at a certain time is independent to the hydraulic fracturing pressure () and is controlled by the fluid mobility (). This is because the seepage velocity is proportional to the pressure gradient for Darcy’s flow which is consistent with equations (10) and (11). So here we discuss the dimensionless pressure distribution along the data output line 3 at a constant differential pressure ( MPa,  MPa), and the results are shown in Figure 9. On the whole, the dimensionless pressure decreases quickly as the dimensionless distance increases. The lower the fluid mobility is within the formation, the faster the dimensionless pressure drops, which means a smaller pore pressure disturbance area (Figure 9). For different types of mobility, there is a sharp decrease of the dimensionless pressure near the fracture tip (). And following the sharp drop, there is a small section with a slower dimensionless pressure drop or even a slightly dimensionless pressure increase due to the fracture interaction. The nearby fracture can contribute to the pore pressure and the seepage velocity, which leads to the abnormal high value of the dimensionless pressure according to equation (10). This effect becomes more obvious as the mobility increases. The pore pressure disturbance contours are also added to Figure 9 for comparison.

According to the design method of casing program [26], the safety margin of equivalent density for gas kick is 0.05-0.10 g/cm3. This indicates that it would be safe for the infilling well drilling with similar designs (including drilling fluid density and casing program) of the nearby wells in the areas where the pore pressure disturbance is smaller than 1.05 (i.e., ). If the infilling well needs to penetrate zones with a pore pressure disturbance of to meet reservoir engineering or geological design requirements, the casing program or/and the drilling fluid density for the infilling well should be carefully optimized based on the updated pore pressure distribution. At this condition, the numerical modelling of the pore pressure disturbance as mainly discussed in this study could provide helpful and practical results.

4. Field Applications

The spacing of two nearby paralleled horizontal wells is relatively large at the early development stage of the Fuling shale gas field, mainly 600-1300 meters [27, 28]. By contrast, the well spacing within the same drilling unit for shale gas production in America is usually 200-300 m [29]. For such a kind of large well spacing in Fuling, infilling well drilling is a promising way to accelerate recovery. However, the pore pressure of the interested block has been greatly disturbed by numerous multistage hydraulic fracturing, which push the drilling work of infilling wells into great risk.

Well 81-2HF, Well 81-3HF, and Well 81-4HF are horizontal wells within the same well platform no. 81 (Figure 10). These three wells have similar casing program: (3rd spud). Their targets are within the same formation.

Well 81-2HF was drilled to 5005 m ( m) at the third spud using oil-based mud (OBM) with a density of 1.47 g/cm3 (Figure 10). And it was completed on June 12, 2015. Then, the 1000 m horizontal section was hydraulically fractured with fifteen stages. The fracturing work continued for about ten days and finished on Aug. 9, 2015. Nearly 31200 m3 of fluid was injected to the formation. The pump pressure was about 50~95 MPa with a pump displacement of 2~15 m3/min. The tested transient gas production was  m3/d. Then, Well 81-2HF was closed due to an unfinished gas pipeline network.

The pore pressure disturbance caused by the hydraulic fracturing in Well 81-2HF was numerically modelled. Although there was no microseismic detection data to quantitatively describe the fracture half-length of Well 81-2HF, the microseismic data of 104 stages of similar fracturing stimulations from 4 wells in the Fuling Gas Field showed that the average fracture half-length was 189.4 m. Our modelling results based on the fracturing parameters of Well 81-2HF showed that the pore pressure disturbance contour of was about 48.73 m from the facture tip along the fracture propagation direction. However, the nearest distance from the target point of Well 81-2HF to that of Well 81-3HF was only 217.4 m at  m. This means Well 81-3HF penetrated the zones with the pore pressure disturbance of , which has a high probability of undergoing gas kicks.

In fact, the drilling practice of Well 81-3HF verified the above prediction. Well 81-3HF and Well 81-4HF were drilled by the “well factory” mode after Well 81-2HF was fractured. The drilling sequence is shown in Figure 11. The third spud drilling of Well 81-3HF was seriously influenced by the multistage fracturing of Well 81-2HF. When Well 81-3HF was at the third spud, Well 81-2HF was closed with the wellhead pressure as high as 40 MPa. There were gas kicks fifteen times in Well 81-3HF within the depth of 2994~4505 m. Among these fifteen gas kicks, thirteen kicks had to close the well and exhaust the gas by burning. The flames of these burnings reached 2~6 meters, and the average burning time was nearly 70 min. The drilling fluid was weighted seven times, and its density increased from 1.36 g/cm3 to 1.85 g/cm3 (Figure 11). The nonproductive time due to these kicks reached 288 hours, and the average ROP of the third spud decreased to 5.17 m/h. As a comparison, there was no gas kick during the third spud drilling of Well 81-2HF with the drilling fluid density of 1.47 g/cm3 and the average ROP was 6.77 m/h.

Direct evidence of the drilling troubles in Well 81-3HF related to the fracturing of Well 81-2HF was a string sticking at 4129 m on Apr. 21, 2016. It failed to release the sticking by conventional methods. But after the wellhead pressure of Well 81-2HF dropped to 30 MPa at 8:00 am on Apr. 24, the drill string was successfully released (Figure 11). This indicated that the geostress and pore pressure disturbance caused by the hydraulic fracturing were responsible for the sticking.

On the other hand, the nearest distance from the target point of Well 81-2HF to that of Well 81-4HF was about 528.0 m at  m. So according to the modelling results, the trajectory of Well 81-4HF is out of the pore pressure disturbance area of Well 81-2HF. The drilling practice of Well 81-4HF was consistent with this result. There was only one gas kick for the third spud drilling of Well 81-4HF with the drilling fluid density of 1.54 g/cm3 (Figure 10), and the average ROP was 8.14 m/h.

5. Conclusions

A numerical modelling method of pore pressure redistribution around a multistage fractured horizontal wellbore was put forward based on the theory of fluid transportation in porous media. The fracture network of each stage was represented by an elliptical zone with high permeability. Five stages of fracturing were modelled simultaneously to consider the interactions among fractures.

The effects of formation permeability, fracturing fluid viscosity, and pressure within the fractures on the pore pressure disturbance were numerically investigated. Modelling results indicated that the pore pressure disturbance zone expands as the permeability and/or the differential pressure increases, while it decreases when the viscosity of the fracturing fluid increases. The pore pressure disturbance level becomes weaker from the fracture tip to the far field along the main-fracture propagation direction. The pore pressure disturbance contours obviously have larger slopes with the variation of permeability than those of the differential pressure. The distances between the pore pressure disturbance contours are smaller at lower permeability and higher viscosity.

The modelling results of the updated pore pressure distribution is of great importance for safe drilling. A case study of three wells within one platform showed that the modelling method could provide a reliable prediction of the pore pressure disturbance area caused by multistage fracturing.

Data Availability

The data in the tables used to support the findings of this study are included within the article. The data in the figures used to support the findings of this study are available from the corresponding author ([email protected]) upon request.

Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.

Acknowledgments

This study is supported by the State Key Laboratory of Shale Oil and Gas Enrichment Mechanisms and Effective Development, China (No. ZC0607-0016); the National Natural Science Foundation of China (No. 51704309); and the Fundamental Research Funds for the Central Universities, China (18CX07008A).