The Analysis of Hydraulic Fracture Morphology and Connectivity under the Effect of Well Interference and Natural Fracture in Shale Reservoirs

Abstract

Employing multi-stage fracturing technology in horizontally accessed wells is a well-known way to successfully develop shale reservoirs. The interaction between natural fractures and hydraulic fractures has a significant impact on the fracturing effect. In this study, a coupled model of rock deformation and fluid flow was established using the cohesive zone method to simulate the propagation of hydraulic fractures under the synergistic effect of natural fractures and wellbore interference. The influence of in situ stress, fracture spacing, the number of fracture clusters, and the fracturing methods on the formation of fracture networks was analyzed. Studies on the fracture morphology and connectivity of fracture networks show that when the in situ stress difference is small, multiple fractures can easily form, and when the in situ stress difference is large, they can easily gather into a single fracture. An excessive reduction in fracture spacing may impede the optimal propagation and interconnection of hydraulic fractures. The findings reveal that augmenting the fracture spacing from 5 m to 8 m results in a significant 15.59% increase in the overall extent of fracture propagation. Moreover, it also adds to the complexity of the fracture network. Increasing the number of hydraulic fracturing clusters can improve the fracture length and fracture propagation complexity. When the number of fracturing clusters increased from two clusters to five clusters, the maximum fracture propagation width increased by 25.23%. Comparing sequential fracturing and simultaneous fracturing, the results show that simultaneous fracturing can form a more complex fracture network with better connectivity, which is conducive to increasing oil and gas production. The obtained results can provide a reference for horizontal well fracturing designs of shale reservoirs.

1. Introduction

Hydraulic fracturing technology is a key technique for developing unconventional reservoirs, which has been widely applied in low−permeability and low−porosity reservoirs [1,2,3]. Horizontal well fracturing technology creates a complex fracture network by communicating hydraulic fractures with natural fractures in unconventional reservoirs, which forms effective communication channels and increases oil and gas recovery rates [4,5,6]. Hydraulic fracturing technology forms fractures in the reservoir by injecting high−viscosity fracturing fluid into the reservoir through the wellbore, using high−pressure pumps on the surface. When the injection rate of fracturing fluid exceeds the absorption capacity of the reservoir, a high pressure is formed at the bottom of the well. When this pressure exceeds the fracture pressure of the reservoir rock near the bottom of the well, the reservoir will be compressed and fractures will be initiated. At this time, the fracturing fluid is injected into the reservoir continuously, and the fracture will propagate in the reservoir. In order to keep the fracture open, a sand carrier fluid with a proppant is pushed into the reservoir. The carrying sand liquid in the fracture can cause the fracture to continue to propagate while supporting the fracture that has been pressed open so that it does not close. Then, the displacement fluid is injected to displace all of the sand−carrying fluid in the wellbore into the fracture, and the fracture is supported with quartz sand. Finally, the injected high−viscosity fracturing fluid automatically degrades and drains out of the wellbore, leaving one or more fractures in the reservoir and creating a new fluid channel between the reservoir and the wellbore. After fracturing, oil and gas well production generally increases significantly. Fracture morphology represents the fracture’s geometric parameters, including the fracture shape, the fracture length, the fracture width, the fracture complexity, etc. Fracture morphology is an index for the quantitative analysis of the fracturing effect. During hydraulic fracturing, the formation of hydraulic fractures and natural fractures in the reservoir form a complex fracture pattern. The size and complexity of fractures reflect the effect of hydraulic fracturing. Generally, the larger the size of fractures, the more natural fractures communicated by fractures there are, and the more complex the fracture shape formed, the better the effect of fracturing [7,8].

Horizontal wells have the advantages of a large drainage area and the high output of a single well. According to the development experience of the Barnett shale play in the United States, the ultimate recoverable reserves from horizontal wells are more than three times those of vertical wells, with costs only 1.5 times those of vertical wells. Since drilling, completion, and fracturing operations require a lot of equipment, including drilling, fracturing equipment, etc., in order to save fracturing costs and reduce equipment relocation time, well factory multiple−stage hydraulic fracturing in horizontal wells is a commonly used method for economically developing unconventional reservoirs [9,10], such as shale gas in the Longmaxi Formation and shale oil in continental shale. The “well factory” model focuses on drilling, completion, fracturing, and other operations on a platform to form a development factory with dense well locations to improve drilling and fracturing timelines while achieving the centralized management of development. The “well factory” fracturing mode can improve the volumetric reconstruction efficiency of shale reservoirs and reduce construction costs, but the formation of the complex fracture network is affected by the hydraulic fracturing design, such as fracture spacing [11,12]. If the fracturing design is not reasonable, it may have an inhibitory effect on the formation of the fracture network. Therefore, it is of great significance to study the fracture propagation model of horizontal wells in naturally fractured reservoirs for the well factory fracturing design.

The Lower Silurian Longmaxi Formation in Weiyuan, Sichuan Basin is a high−quality shale formation rich in hydrocarbon gas, making it a hot spot for shale gas exploration and development in China. However, the organic carbon content (TOC) is generally low, ranging from 2% to 5%, and the gas content is also low. Additionally, it is located at a deep burial depth of 2800 to 4000 m, with complex surface conditions, presenting significant challenges for development. The natural fractures in the Weifeng Formation−Long Yi 1 sub−section in the Weiyuan area are relatively developed, with microfine fractures dominating. The main types of fractures present are structural fractures, diagenetic fractures, dissolution fractures, and hydrocarbon generation fractures. Based on the core and imaging logging data from well Wei 202, the bedding fractures are widely developed in the Wufeng−Longmaxi Formation, while high−angle fractures are less developed [13,14]. The values of Young’s modulus are relatively high, while the values of Poisson’s ratio are relatively low. The triaxial compressive strength ranges from 97.7 to 281.6 MPa, with an average value of 213.90 MPa. The values of Young’s modulus range from 11 to 33 GPa, with an average of 21 GPa, and the values of Poisson’s ratio range from 0.17 to 0.29, with an average of 0.20. To achieve the large scale and economically efficient development of shale gas in Weiyuan, Sichuan, a “factory-style” well deployment model with a pad and cluster horizontal well system is adopted [15]. Through geological−drilling−fracturing integration design and construction, the daily gas production can reach 1.4 million cubic meters, while the daily water production is 201 cubic meters. The cumulative gas production is 1.252 billion cubic meters, and the cumulative water production is 803,600 cubic meters. This approach has realized the large scale and profitable development of shale gas.

The key problem of well plant fracturing technology is how to increase the reservoir transformation volume and reduce the fracturing cost as much as possible, which is the problem of fracturing mode and fracturing parameter optimization [16]. However, at present, the proposed well plant fracturing mode is more based on practical field experience, and there is less research on the fracturing mechanism. The fracture propagation mechanism of well plant fracturing technology is the result of the interaction mechanism between hydraulic fractures and natural fractures and the synergistic control of multiple hydraulic fracture propagation mechanisms in horizontal wells. During the study of hydraulic fracture mechanisms, the numerical simulation method is an important technical method used to study complex rock fracture problems and dynamic propagation of hydraulic fractures due to the non−linear propagation of fracture propagation and complex mechanical problems [17,18].

Many researchers have investigated the mechanism of interaction between hydraulic fractures and natural fractures. Wu and Olson et al. [19] established a 3D model based on displacement discontinuity and the finite element method, using Newton–Raphson iteration and Picard iteration to solve the coupled problem of fluid flow and rock mechanics. Xi et al. [20] combined the Mohr–Coulomb criterion with the cohesive fracture model to develop a new fracture model that analyzes the near−well interaction between hydraulic fractures and natural fractures. Dehghan et al. [21] established a 3D numerical model based on the extended finite element method, simulating the propagation process of hydraulic fractures intersecting with natural fractures, and fractally analyzed the effects of natural fracture dip angle, strike, and reservoir properties on hydraulic fracture morphology. Rahman et al. [22] simulated coupled hydraulic fractures and arbitrarily oriented natural fractures using the grid method and finite element method, studying the effects of natural fracture orientation and rock mechanics parameters on induced hydraulic fracture propagation. Zhao et al. [23] combined fractal theory and a cohesive fracture model to establish a random fracture network distributed in a reservoir with fractal characteristics, the hydraulic fracture dynamic propagation model, and fractally analyzed the effects of natural fractures’ fractal characteristics on fracture propagation morphology and connectivity. Fatahi et al. [24] simulated the interaction between hydraulic fractures and natural fractures using the discrete element numerical simulation method and validated the results through laboratory experiments. They analyzed the effects of natural fracture angle on hydraulic fracture shear slip. Liu et al. [25] applied a cohesive fracture model based on core test results to investigate the effects of the development of different shale beddings on hydraulic fracture propagation. In terms of multi−fracture propagation in horizontal wells, Tian et al. [26] used the extended finite element method to conduct a numerical simulation study on three different fracturing modes of multi−stage hydraulic fracturing in horizontal wells, analyzing the factors affecting the induced stress and fracture propagation morphology of different fracturing modes. Haddad et al. [27] established a fully coupled poroelasticity displacement discontinuity method (DDM) model of the plane cohesive zone and an extended finite element model to develop a hydraulic fracturing propagation model for three−cluster horizontal wells, revealing the stress shadow effects between hydraulic fractures on fracture propagation, such as merging, parallel growth, or mutual repulsion. Sesetty and Ghassemi [28] studied zipper and MZF using the displacement discontinuity method (DDM). Their work summarized that the spacing of previous fractures and boundary conditions affect the later fracture paths. Zheng et al. [29] established a numerical model of rock deformation and fracturing fluid flow coupling and studied the dynamic propagation of multi−well fracturing fractures under the influence of well interference but did not consider the effects of natural fractures in the reservoir.

In summary, the current research mainly focuses on the intersection and propagation mechanism of single hydraulic fractures and natural fractures and the interference mechanism between horizontal well fractures, while few studies consider the fracture propagation mechanism and propagation pattern analysis under the action of inter−well interference, multiple hydraulic fracture interference, and natural fracture interference. Under the “well factory” fracturing mode, research on the influence of fracturing parameter design on the fracturing effect is insufficient. Therefore, in order to analyze the fracture propagation law of the “well factory”, this article establishes the coupling rock deformation and fracture hydraulic action, combined with the natural fracture network model and finite element numerical method, and based on the ABAQUS finite element analysis platform, establishes the numerical simulation of multi−fracture propagation under the cooperative disturbance of natural fractures and multiple wells through secondary development. The propagation mechanism of multiple hydraulic fractures under multi−well disturbance in fractured reservoirs is analyzed, and the influence of pressure parameters such as ground stress deviation, fracture spacing, number of fractures, and fracturing methods on the fracture propagation morphology is examined.

2. Methodology

The finite element software ABAQUS 2020 provides a cohesive element, which can be used to set element material properties that are different from the solid unit properties. In this paper, a fluid–solid coupling model of industrial hydraulic fracturing of fractured reservoir wells was established by applying the global embedding method. Considering the seepage effect of reservoirs, the model adopted a bilinear constitutive model, the fracture criterion was the maximum damage fracture criterion, the model size was 100 m × 100 m, the horizontal well spacing was 50 m, and the fracture cluster spacing was 5 m. The fracturing method was simultaneous fracturing.

2.1. Constitutive Model

In this study, the method of globally embedded cohesive elements was used to simulate the free propagation process of fractures in the reservoir. The failure mode of rock− is generally brittle, and cohesive elements are used to simulate the characteristics of random fracture in the rock medium. Therefore, the interface−based traction−separation description method [30] was adopted. The elastic behavior was represented by an elastic constitutive matrix, which links the normal and shear stresses of the fracture element.

$$t=\left\{\begin{array}{c}{t}_n\\ t_s\\ t_t\end{array}\right\}=\left[\begin{array}{ccc}{K}_{nn}& 0& 0\\ 0& K_{ss}& 0\\ 0& 0& K_{tt}\end{array}\right]\left\{\begin{array}{c}{\delta }_n\\ {\delta }_s\\ {\delta }_t\end{array}\right\}=K\delta$$

(1)

In this equation, $t$ represents the traction stress, $t_n$ represents the normal stress, and $t_s$ and $t_t$ represent the shear stress in different directions, corresponding to the separation using ${\delta }_n$, ${\delta }_s$, and ${\delta }_t$. $K_{nn}$ is the stiffness matrix of the initial cohesive strength, which represents the stiffness of the cohesive element in the elastic stage.

The model used a bilinear constitutive model [31] to simulate the fracture process of two−dimensional brittle materials, which is linearly elastic until reaching the damage evolution and then linearly softening. The initial damage criterion is the quadratic nominal stress criterion, which means damage starts when the square ratio of the nominal stress in each direction is equal to 1, and can be expressed as:

$${\left\{\frac{〈t_n〉}{t_n^0}\right\}}^2+{\left\{\frac{t_s}{t_s^0}\right\}}^2=1$$

(2)

where $t_n^0$ and $t_s^0$ are nominal stress peaks when the normal and tangential deformation are perpendicular to the interface, respectively. $t_n$ and $t_s$ are normal and tangential stress components, respectively. The symbol indicates a state of pure compressive deformation or stress that does not cause damage.

2.2. Fracture Criterion

The damage evolution criterion is the discriminant for describing the material’s failure after reaching the critical fracture initiation point. For the propagation of cohesive fractures under the compression−shear mode, a rock fracture is caused by tangential stress, and the strength is determined by the cohesive force, friction angle, and normal compressive stress [32,33]. For the tension−shear mode of cohesive fractures, a rock fracture is caused by normal stress and tangential stress, and the strength is determined by tensile strength and cohesive strength [34,35]. The fracture criterion used in this model is as follows.

$$t_n=\left\{\begin{array}{ll}\left(1-D\right)t_n^{-},& t_n^{-}\ge 0\\ t_n^{-},& t_n^{-}<0\end{array}\right.$$

(3)

$$t_s=\left(1-D\right)t_s^{-}$$

(4)

$$t_t=\left(1-D\right)t_t^{-}$$

(5)

where $D$ is the overall damage of the material, capturing the combined effect of all active mechanisms. The initial value of $D$ is 0. If damage evolution is modeled, $D$ monotonically evolves from 0 to 1 after damage initiation and further loading. These parameters are the stress components predicted by the undamaged elastic traction separation behavior of the current strain.

The parameters $t_n^{-}$, $t_s^{-}$, and $t_t^{-}$ are the stress components predicted by the elastic traction separation behavior of the current strain without damage.

$$D=\frac{{\delta }_m^f\left({\delta }_m^{\max }-{\delta }_m^0\right)}{{\delta }_m^{\max }\left({\delta }_m^f-{\delta }_m^0\right)}, \mathrm{where} {\delta }_m=\sqrt{{〈{\delta }_n〉}^2+{\delta }_s^2+{\delta }_t^2}$$

2.3. Fluid–Structure Coupling Model

We assume that the fracturing fluid is an incompressible Newtonian fluid and the fluid flow in the cohesive units is divided into tangential and normal flow. Normal flow indicates the leakage of fracturing fluid into the formation, and tangential flow is the driving force for fracture propagation [36], as shown in Figure 1.

The model considers the process of fracturing fluid leaking from the fracture to the matrix. Within the framework of Biot’s theory, the governing equations for the coupled solid–fluid deformation in a porous medium were established [24,37].

$$\nabla \cdot \sigma -\rho \ddot{u}+\rho b=0$$

(6)

$$\nabla \cdot \left[k_f\left(-\nabla p-{\rho }_f\ddot{u}+{\rho }_{f}b\right)\right]+\alpha \nabla \cdot \dot{u}+\frac{1}{Q}\dot{p}=0$$

(7)

To simulate fluid flow in discontinuous fractures, the continuity equation for fluid flow in the fracture can be written as:

$$\nabla \cdot \dot{w}+\alpha \nabla \cdot \dot{u}+\frac{1}{k_f}\dot{p}=0$$

(8)

Assuming that the fracturing fluid is an incompressible Newtonian fluid, the fluid flow in the cohesive element can be divided into tangential flow and normal flow.

Fluid flow in fractures satisfies the cubic law:

$$qd=-\frac{d^3}{12\mu }\nabla p$$

(9)

where $q$ represents the volume flow rate vector per unit area of cohesive element, $d$ is the aperture of the cohesive element, $\nabla p$ indicates the pressure gradient along the direction of fluid flow, and $\mu$ represents the viscosity of the fluid. Thus, normal flow in fractures is based on the following principles

$$\left\{\begin{array}{l}{q}_t=c_t(p_i-p_t)\\ q_b=c_b(p_i-p_b)\end{array}\right.$$

(10)

where $q_t$ and $q_b$ represent the fluid flow rates at the upper and lower surfaces of the cohesive element, respectively, and $c_t$ and $c_b$ represent the filtration coefficients of the upper and lower surfaces, respectively. $p_i$ represents the pressure interpolated through virtual nodes in the cohesive element, and $p_i$ and $p_b$ represent the fluid pressures on the upper and lower surfaces of the cohesive element, respectively.

According to the mass balance of the fluid, a portion of the injected fluid will fill the fracture within a certain period of time, and the remaining fluid will leak into the rock matrix of the fracture, as shown in the following equation.

$$\frac{\partial d}{\partial t}+\nabla \cdot q+\left(q_t+q_b\right)=q_{inj}\delta \left(x,y\right)$$

(11)

Therefore, the governing equations of the rock matrix consist of the coupled equations for fluid flow and rock deformation.

$${\sigma }_{ij}-{\sigma }_{ij}^0=\frac{E}{1+\nu }\left({\epsilon }_{ij}+\frac{v}{1-2v}{\epsilon }_{kk}{\delta }_{ij}\right)-\alpha \left(p_w-p_w^0\right){\delta }_{ij}$$

(12)

2.4. Boundary Conditions

The boundary conditions of the model include the displacement boundary, stress boundary, fluid pressure boundary on the fracture surface, and cohesive element’s cohesive force. The displacement boundary and stress boundary [30] are shown as

$$u=\overline{u} \mathrm{on} {\Gamma }_u$$

(13)

$$\sigma \cdot n_t=\overline{t} \mathrm{on} {\Gamma }_t$$

(14)

$$\sigma \cdot n_f=t_c-p_f{n}_f \mathrm{on} {\Gamma }_f$$

(15)

where n_t and n_f are the unit outward normal vector of the outer and inner boundaries of stress, respectively. For the fluid flow part, the boundary conditions include the pressure boundary and the fluid boundary.

$$p=\overline{p} \mathrm{on} {\Gamma }_p$$

(16)

$$v_m\cdot n_q=\overline{q} \mathrm{on} {\Gamma }_q$$

(17)

$$Q_f=Q_0 \mathrm{on} {\Gamma }_{f,inlet}$$

(18)

where $n_q$ is the unit outward normal vector of the fluid boundary.

2.5. Numerical Model

For the analysis of the fracture propagation under the influence of the synergistic effect of inter−well interference and the natural fractures in shale reservoirs, the cohesive elements were implemented via secondary development in Python. A well factory model for hydraulic fracture propagation was established, taking into account the coupling effect of natural fracture deformation and fracturing fluid leakage. The model had a size of 100 m × 100 m, the distance between two horizontal wells was 50 m, the distance between fracture clusters was 5 m, and the fracturing method was simultaneous fracturing, as shown in Figure 2. The numbers in the figure represent the fracture clusters. This model was used for the analysis of the results.

A numerical model for hydraulic fracturing pressure calculation was established based on the natural fracture characteristics and hydraulic fracture propagation model of Weiyuan shale gas in Sichuan. The model considered two horizontal wells with staggered hydraulic fractures and two clusters of perforation in the horizontal wellbore. The hydraulic fracturing pressure was calculated in a reservoir with randomly distributed natural fractures. The solid part of the model was meshed with pore fluid elements (CPE4P), while the fracture part was meshed with a cohesive element (COH2D4P) that considered filtration. The mesh was locally refined in the preset hydraulic fracture region. The model parameters are shown in

Table 1. Numerical simulation parameters.

Parameters	Values	Parameters	Values
Elastic modulus (GPa)	21	Porosity (%)	10
Poisson’s ratio	0.2	Permeability (mD)	0.1
Natural fracture angle (°)	0	Single cluster injection rate (m3/s)	0.001
Maximum horizontal effective stress (MPa)	3	Cluster spacing (m)	5
Minimum horizontal effective stress (MPa)	1	Fluid viscosity (mPa·s)	1

.

3. Fracture Initiation Mechanism

3.1. Stress Intensity Factor Theory

Applying elasticity, Westergaard and Muskhelishvili analyzed the calculation of the stress field and displacement field near the tip of a slit plane in an ideal uniform continuous elastic solid, taking a type I crack as an example [38].

$$\left\{\begin{array}{l}{\sigma }_{xx}=\frac{{\mathrm{K}}_{\mathrm{I}}}{\sqrt{2\pi r}}\cos \frac{\mathsf{\theta }}{2}(1-\sin \frac{\mathsf{\theta }}{2}\cdot \sin \frac{3\mathsf{\theta }}{2})+o(r^{-1/2})\\ {\sigma }_{yy}=\frac{{\mathrm{K}}_{\mathrm{I}}}{\sqrt{2\pi r}}\cos \frac{\mathsf{\theta }}{2}(1+\sin \frac{\mathsf{\theta }}{2}\cdot \sin \frac{3\mathsf{\theta }}{2})+o(r^{-1/2})\\ {\tau }_{xy}=\frac{{\mathrm{K}}_{\mathrm{I}}}{\sqrt{2\pi r}}\sin \frac{\mathsf{\theta }}{2}\cos \frac{\mathsf{\theta }}{2}\cos \frac{3\mathsf{\theta }}{2}+o(r^{-1/2})\end{array}\right.$$

(19)

Using the stress function proposed by Westergaard, the full field formula of type I crack stress was calculated as follows:

$${\mathrm{K}}_{\mathrm{I}}={\sigma }_y^{\infty }\sqrt{\pi a}$$

(20)

The stress intensity factor K_I is the only parameter used to characterize the stress field of a crack tip. The crack shape and force in different rocks are different, but the same K_I has the same stress field of the crack tip.

Fracture toughness, also known as the critical stress intensity factor, can characterize the singularity strength of the tip field of linear elastic cracks, independent of the geometry and load of the crack body. It is a basic property of materials and a key parameter for crack analysis. When the crack expands, the stress intensity factor corresponding to the extension pressure is the fracture toughness K_IC, that is:

$$K_{IC}=p_c\sqrt{\pi L}$$

(21)

According to the linear elastic fracture mechanics, when the stress intensity factor is greater than the fracture toughness of the material, the rock will have tensile fractures and the crack will continue to extend. The Irin criterion [39]:

$${\mathrm{K}}_{\mathrm{I}}\ge K_{IC}$$

(22)

3.2. Induced Stress Field Distribution Model

According to the theory of elastoplastic mechanics, under the assumption of homogeneous, isotropic, and plane strain conditions, the schematic diagram of the induced stress field of a single hydraulic fracture is shown in Figure 3. It is assumed that the fracture is vertical and the fracture profile is elliptical. ${\sigma }_H$ and ${\sigma }_h$ represent the original maximum and minimum principal stresses, respectively.

Sneddon studied the change in stress field caused by the propagation of a single hydraulic fracture and derived the formula for calculating the induced stress of a single fracture [40]:

$$\left\{\begin{array}{l}\frac{1}{2}\left(\Delta {\sigma }_y+\Delta {\sigma }_x\right)=-p\left[\frac{r\cos \left(\theta -0.5{\theta }_1-0.5{\theta }_2\right)}{\sqrt{r_1{r}_2}}-1\right]\\ \frac{1}{2}\left(\Delta {\sigma }_y-\Delta {\sigma }_x\right)=p\frac{2r\sin \theta }{L}{\left(\frac{L^2}{4r_1{r}_2}\right)}^{\frac{3}{2}}\sin \left[3\left({\theta }_1+{\theta }_2\right)/2\right]\end{array}\right.$$

(23)

where $\Delta {\sigma }_x$ is the change in the minimum horizontal principal stress, perpendicular to the hydraulic fracture. $\Delta {\sigma }_y$ is the change in the maximum horizontal principal stress, parallel to the hydraulic fracture; $p$ is the water pressure in the joint; $L$ is the fracture height; and $r$, $r_1$, and $r_2$ are the lengths from the target point to the middle point and end of the crack, respectively. $\theta$, ${\theta }_1$, and ${\theta }_2$ are the angles between the corresponding connection and the fracture.

If there are multiple fractures, the induced stress of multiple fractures is superimposed at any point in the formation, and the induced stress field of multiple fractures interferes with the initiation and propagation of adjacent fractures [41].

According to the superposition principle, the expression of the stress field at the n hydraulic fracture is as follows:

$$\left\{\begin{array}{l}{\sigma }_{h\left(n\right)}^{\prime }+{\displaystyle \sum _{i=1}^{n-1}{\sigma }_{x\left(in\right)}}\\ {\sigma }_{H\left(n\right)}^{\prime }={\sigma }_H+{\displaystyle \sum _{i=1}^{n-1}{\sigma }_{y\left(in\right)}}\end{array}\right.$$

(24)

where ${\sigma }_{h\left(n\right)}^{\prime }{\sigma }_{H\left(n\right)}^{\prime }$ are the horizontal maximum principal stress and horizontal minimum principal stress of fracture n, respectively. ${\sigma }_{x\left(in\right)}{\sigma }_{y\left(in\right)}$ are the induced stress components of fracture i in the x and y directions at fracture n, respectively.

The hydraulic fracture propagation direction and morphology will become more complex under the disturbance of fracture−induced stress, and fracture propagation is not only limited by the original in situ stress field. Under the condition of multiple fracture propagation, the induced stress field becomes more prominent. The induced stress interference has both positive and negative effects on hydraulic fracturing engineering. The induced stress forces the fracture to deflect, and it is possible to communicate adjacent cracks and form a complex fracture network. The induced stress affects fracture deflection, which may communicate adjacent fractures and form a complex fracture network. It may also hinder the initiation of adjacent cracks or cause complex propagation paths that are difficult to predict.

According to the analysis of the theoretical model, the influencing factors of the induced stress field are in situ stress, fracture spacing, fracture size, and fracturing sequence. Considering that the hydraulic coupling process and induced stress field are difficult to accurately describe using analytical solutions, numerical simulation provides an effective way to study fracture−to−joint interference under complex conditions.

4. Results and Discussion

Hydraulic fracturing design is the key to determining the fracturing effect. It is necessary to combine the actual geological conditions during well factory fracturing. Research shows that in situ stress difference is the key factor affecting fracture morphology [42,43,44,45]. The hydraulic fractures’ spacing, the quantity of hydraulic fracturing and the method of fracturing are important components of fracturing design, which have a great influence on the fracture propagation. Therefore, the analysis of the hydraulic fracture propagation form and the connectivity of natural fractures with different in situ stress differences and fracturing design parameters can provide some guidance for fracturing construction design.

4.1. Influence of In Situ Stress Difference

In situ stress difference has an important influence on fracture growth morphology. An analytical model of fracture growth under the well factory mode is established when the in situ stress difference is 2 MPa, 3 MPa, and 4 MPa and there is a group of random natural fractures in the reservoir. The fracture growth results are shown in Figure 4.

It can be seen from Figure 4 that different in situ stress differences affect the fracture propagation morphology when the hydraulic fracture meets between wells and the fracture morphology when the hydraulic fracture intersects with the natural fracture. When the stress difference is small, under the effect of stress concentration, the hydraulic fractures between wells do not easily connect when they meet. When the in situ stress difference is large, the hydraulic fractures meet and converge into one fracture. When the stress difference is 2 MPa, the hydraulic fracture is more likely to be captured by the natural fracture in the process of fracture propagation when it meets the natural fracture. When the stress difference is 4 MPa, the hydraulic fracture is more likely to pass through the natural fracture and form a simple vertical fracture. Usually, when the stress difference is low, the natural cracks can become more connected, complex natural cracks can form, and the transformation area is larger, which is preferred. When the stress difference of the reservoir is high, the stimulation effect of fracturing is poor, and it is not easy to exploit oil and gas.

4.2. Influence of Fracture Spacing

In the well factory model, the fracture spacing has a significant impact on the hydraulic fracture propagation pattern. A fracture propagation analysis model is established with fracture spacing of 5 m, 6 m, 7 m, and 8 m. The results of fracture propagation morphology and fracture mode are shown in Figure 5.

It can be seen from Figure 5 that when the fracture spacing is small, the hydraulic fracture propagation is inhibited under the action of the stress shadow. The inter−well fractures can easily converge and connect, forming a single main fracture. When the fracture spacing is 5 m, only one hydraulic fracture in horizontal well 1 propagates, and the other fracture is inhibited. When the fracture spacing is 8 m, the perforation clusters of horizontal well 1 are propagated, and the effect between fractures is small. When the fracture spacing is 6 m, the hydraulic fractures at perforation 1 and perforation 4 communicate and converge into one main fracture, which affects the fracturing reconstruction volume.

We plotted the comparison curves of the number of natural fractures connected by hydraulic fractures and the total length of fractures under different fracture spacings, as shown in Figure 6 and Figure 7, respectively.

As shown in Figure 6 and Figure 7, under the combined effect of hydraulic fractures and natural fractures, the spacing between fractures has a significant impact on the morphology of fractures. As the spacing between fractures increases, the number of natural fractures connected by fractures increases, indicating that a small spacing between fractures can suppress the formation of complex fractures. The total length of fracture propagation increases by 15.59% because when the spacing between fractures is small, the stress shadow between multiple fractures is large, which limits the propagation of fractures. The maximum width of fractures decreases with the increase in the spacing between fractures, but the variation in fracture width fluctuates. This is because the propagation of fractures is influenced by the competition between hydraulic fractures, and on the other hand, when there are more natural fractures connected, leakage can occur and affect the propagation width of fractures. Therefore, for different reservoirs, optimizing the spacing of hydraulic fracturing in the fracturing design can increase the length of fracture propagation and form a complex fracture network.

4.3. Influence of the Number of Fractures

The number of hydraulic fractures has a significant influence on fracture propagation due to the interference between the fractures. The influence of fracture number on the fracture morphology was analyzed by studying the hydraulic fracture propagation models with two, three, four, and five fractures. The fracture propagation results are shown in Figure 8.

It can be seen from Figure 8 that the hydraulic fractures exhibit an asymmetric and complex fracture pattern, intersecting and converging with natural fractures to form fractures with different directions. With the increase in fracture number, the fracture morphology propagation becomes more complex, and the volume of fracturing reconstruction increases. When there are two hydraulic fracturing fractures, the fracture morphology is singular, and the fractures converge into one main fracture. When there are multiple hydraulic fracturing fractures, under the disturbance of multiple fractures, the hydraulic fractures propagate in more directions, forming multiple multi−direction fracture morphologies.

In order to further compare the influence of different fracture numbers on fracture morphology, the comparison curves of the total length of fractures and the number of natural fractures connected by hydraulic fractures are plotted in Figure 9 and Figure 10.

From Figure 9 and Figure 10, it can be seen that under the combined action of hydraulic and natural fractures, increasing the number of hydraulic fracture clusters has an important influence on fracture morphology. With the increase in the number of fracture clusters, the total length of the fractures increases from 149.15 m to 323.16 m, with a growth rate of 116.67%, which is a significant increase. At the same time, the maximum width of the fractures increases by up to 25.23%, but this change is not linear. This is due to the complex competition between hydraulic and natural fractures. The number of fractures communicating with natural fractures increases from two to nine, which significantly increases the complexity of the fracture network, and the direction of fracture propagation also increases, leading to an increase in the reservoir control area. Therefore, the number of hydraulic fracture clusters plays a very important role in the propagation morphology and complexity of the fractures.

4.4. Influence of Fracturing Method

Simultaneous fracturing is the simultaneous fracturing of two wells. Sequential fracturing is the process of fracturing separately and sequentially, with one well completed followed by another. In the well factory mode, different fracturing sequences have a great influence on the fracture propagation pattern, and we compare and analyze the fracture propagation pattern of simultaneous fracturing and sequential fracturing, as shown in Figure 11.

It can be seen from Figure 11 that the hydraulic fractures’ morphology in both wells is different for the two fracturing sequences. In simultaneous fracturing, the hydraulic fractures propagate simultaneously with less inter−well interference and are affected only in the later stage of hydraulic fracture propagation. In sequential fracturing, there is more inter−well interference, the formation of hydraulic fractures in one well causes stress changes near well 2, and the effect of stress concentration near the hydraulic fractures greatly inhibits the propagation of hydraulic fractures in well 2.

In order to further analyze the fracture growth morphology under different fracturing methods, the comparison curves of the number of hydraulic fractures communicating with natural fractures, the total length of fractures, and the maximum width of fractures under different fracture spacings are plotted and shown in Figure 12 and Figure 13.

From Figure 12 and Figure 13, it can be seen that the total length of fractures in sequential fracturing is 49.61% lower than that in simultaneous fracturing, and the fracture width is 17.83% lower than that in simultaneous fracturing. This indicates that the fracturing method has a significant impact on fracture size, and under the reservoir properties and conditions considered, multiple hydraulic fractures in simultaneous fracturing converge into one fracture, and the propagation of hydraulic fractures in the second well of sequential fracturing is inhibited under stress disturbance from the first well. Comparing the number of natural fractures communicated by the two fracturing methods, it can be seen that simultaneous fracturing communicates with more natural fractures, reaching nine, while the hydraulic fractures in sequential fracturing are all captured by a longer natural fracture, communicating with only three natural fractures. Therefore, simultaneous fracturing is more likely to communicate with natural fractures, forming a more complex fracture network.

5. Conclusions

The well factory fracturing mode is an effective development mode for shale gas development at present. The well factory fracturing mode is a very complex project. Due to the limitations of computers, the model in this paper is only a 2D model, while the actual fractured reservoir is a 3D model. This model is a mechanism model, and the model size has a certain gap with the actual size, which reaches several hundred meters. The preliminary understanding of the well factory fracturing mode still has some limitations, and further research is still needed.

By establishing a numerical model of hydraulic fracturing under the joint disturbance of natural fractures and inter−well interactions, this study analyzed the hydraulic fractures’ propagation morphology and fracture connectivity. Simulation results showed that asymmetric fracture morphologies were formed, which were complex and varied. It was found that the in situ stress difference, the hydraulic fracture spacing, the number of fracturing clusters, and the fracturing method had significant impacts on the propagation morphology of hydraulic fractures during the fracturing design process. The following insights were mainly obtained:

(1)

The effect of in situ stress difference on fracture propagation was analyzed. When the local stress difference is small, the hydraulic fracture can easily communicate with the natural fracture. When the local stress difference is large, the fractures can easily converge into a single main fracture.

(2)

In the “well factory” fracturing mode, the fracture spacing will inhibit the propagation of hydraulic fractures. When fracture spacing increases from 5 m to 8 m, the total length of fracture propagation increases by 15.59%, and the complexity of the fracture network is promoted.

(3)

Under the “well factory” fracturing mode, increasing the number of hydraulic fracturing clusters in horizontal wells is beneficial to increasing the length and complexity of fractures. When the number of fracturing clusters increases from two clusters to five clusters, the maximum fracture propagation width increases by 25.23%. Closely spaced fracturing has a significant effect on improving the effectiveness of modified fracturing.

(4)

Different fracturing sequences will produce different stress redistributions, which have an important influence on the propagation of hydraulic fractures. During simultaneous fracturing, the inter−well stress disturbance is relatively small, resulting in more complex fractures. During sequential fracturing, the suppression effect on inter−well fracture propagation is more obvious and affects the propagation of hydraulic fractures.