The Effect of Perforation Spacing on the Variation of Stress Shadow

: When the shale gas reservoir is fractured, stress shadows can cause reorientation of hydraulic fractures and affect the complexity. To reveal the variation of stress shadow with perforation spacing, the numerical model between different perforation spacing was simulated by the extended ﬁnite element method (XFEM). The variation of stress shadows was analyzed from the stress of two perforation centers, the fracture path, and the ratio of fracture length to spacing. The simulations showed that the reservoir rock at the two perforation centers is always in a state of compressive stress, and the smaller the perforation spacing, the higher the maximum compressive stress. Moreover, the compressive stress value can directly reﬂect the size of the stress shadow effect, which changes with the fracture propagation. When the fracture length extends to 2.5 times the perforation spacing, the stress shadow effect is the strongest. In addition, small perforation spacing leads to backward-spreading of hydraulic fractures, and the smaller the perforation spacing, the greater the deﬂection degree of hydraulic fractures. Additionally, the deﬂection angle of the fracture decreases with the expansion of the fracture. Furthermore, the perforation spacing has an important inﬂuence on the initiation pressure, and the smaller the perforation spacing, the greater the initiation pressure. At the same time, there is also a perforation spacing which minimizes the initiation pressure. However, when the perforation spacing increases to a certain value (the result of this work is about 14 m), the initiation pressure will not change. This study will be useful in guiding the design of programs in simultaneous fracturing.


Introduction
Shale is the main target rock for global unconventional oil and gas production, and shale gas occupies an indispensable position in the world energy pattern [1]. However, shale reservoirs are dense. To improve the permeability of shale reservoirs and increase shale gas production, horizontal wells and hydraulic fracturing stimulation technology are widely used [2]. When hydraulic fracturing is applied to reservoir reconstruction, multiple perforations can be used for simultaneous fracturing, which can significantly increase fracturing capacity. However, simultaneous fracturing with adjacent perforations can result in a stress shadow phenomenon, which will affect the state of the geostress field, thereby affecting subsequent fractures propagation [3]. Fortunately, reasonable perforation spacing can effectively utilize the positive effect of stress shadow. Therefore, it is important to investigate the influence of perforation spacing on the variation of stress shadow.
In recent years, many scholars have considered the relationship between perforation spacing and stress shadow size through laboratory tests and numerical simulation. Zhou et al. [4] conducted atrue-triaxial hydraulic fracturing experiment to analyze the effect of stress shadow on hydraulic fractures by adjusting the notch distance and found a phenomenon of adjacent hydraulic fractures expanding in the opposite direction. To avoid the stress shadow effect, Morrill and Miskimins [5] obtained the optimum perforation spacing by combining multiple simulation parameters. Liu et al. [6] investigated the hydraulic fracture mutual interference using the finite element method (FEM), and put forward the optimization method of perforation spacing. In addition, the path of hydraulic fracturing is greatly affected by stress shadow [7]. Kresse et al. [8] analyzed the effect of perforation spacing on hydraulic fracture paths and stress shadow using the unconventional fracture model (UFM). Zeng and Yao [9] revealed the intersection pattern of natural and hydraulic fractures, and showed that perforation spacing affects fracture morphology. Weng et al. [10] adopted the UFM to study the interaction between stress shadow and fracture length under different perforation spacing. It was found that the stress shadow is negatively correlated with the perforation spacing. Vahab et al. [11] proposed a multi-state hydraulic fracturing treatment algorithm based on the extended finite element method (XFEM), and used the algorithm to conduct research on multi-perforation fracturing. Khoei et al. [12] studied the interaction between hydraulic fractures and natural fractures using the XFEM method. Furthermore, the change in stress intensity factor (SIF) can also reflect the effect of stress shadow on fracture propagation. Taghichian et al. [13] considered the influence of stress shadow on SIF, which dominates the change in SIF and leads to its substantial reduction.
In addition, perforation density is also a research hotspot. Qi et al. [14] researched the effect of perforation density on the evolution of fracture network by using a 2D coupled flow-stress-damage model. It was obtained that the interaction between micro-fractures and perforations becomes stronger with the increase in perforation density. Luo et al. [15] put forward the control method of variable density perforation and realized the balanced distribution of flow. Ghaderiet al. [16] considered the effect of fracture density on gas saturation distribution. When the fracture density is large, it will only reach a higher oil production rate in a short period. There is an optimum value for the effect of crack density on the oil recovery, and it also depends on reservoir properties. Three-dimensional numerical models were used to analyze the effect of perforation spacing on fracture paths. Additionally, it was found that non-uniform perforation spacing can promote crack propagation [17]. Meanwhile, the two adjacent fractures would deviate from each other under the influence of stress shadow [18]. At present, the same conclusion has been reached in the study of perforation spacing: the smaller the perforation spacing, the greater the size of the stress shadow. Furthermore, the existence of stress shadow will change the fracture paths, leading to the back propagation of adjacent fractures.
Recent studies have shown that perforation spacing, length, and density are the main parameters concerning the size of the stress shadow effect [8,10,16]. However, there are few studies on the change process of stress shadows. Revealing the variation rule of stress shadows is essential for optimizing fracturing design. Thus, two-dimensional numerical models with different perforation spacing were established. Subsequently, the XFEM was used to simulate the variation of stress shadow during the simultaneous fracturing of two adjacent perforations. The change process of stress shadows under different perforation spacing was analyzed.

Principle of the XFEM
By increasing the degrees of freedom and enhancing functions, the XFE Mimplements discontinuity. In addition, the fracture is described by the jump function and the level set function [19]. The displacement vector function can be represented as follows: where the nodal function is denoted by N I (x). a I and b α I represent the joint improvement in the degrees of freedom of the elements penetrated by cracks and the crack tip element, respectively. The Heaviside step function H(x) can be applied to describe fractures [20]. The asymptotic displacement function F α (x) can be represented as follows: The XFEM has been described in detail previously [21].

Fluid-Mechanical Coupling Principle
The mechanical equilibrium equation of rocks can be obtained from the principle of virtual work.
where V and S are the integral space and integral space surface, respectively; δ ε and δ v represent virtual strain field and virtual velocity field, respectively;e is the rockporosity; and T is the external surface force [21].

Fluid Flow Model
Fracture propagation is primarily driven by the fluid pressure generated by fracturing fluid acting on the fracture surface. Supposing the fluid is continuous and incompressible, tangential flow and normal flow can represent fluid flow in fractures.
The tangential flow on the surface of the fracture element can be simulated by the Newtonian model, which can be expressed as follows [21]: The normal flow corresponds to the engineering phenomenon of leak-off, which can be considered as the volume rate at which fluid flows into the simulated area unit. It is represented as follows: where q represents the flow rates; p is the pore pressure of different surfaces; c is theleak-off coefficients; subscripts t and d represent the top and bottom surfaces of a cracked element, respectively [22].

Failure Criterion
This study used the maximum principal stress criterion to simulate fracture.
where σ 0 max is the maximum allowable principal stress. The symbol 〈〉 represents the Macaulay bracket with the usual interpretation. When f is greater than 1 within the tolerance range, the damage is initiated. The damage evolution criterion was represented by damage variable D based on the effective displacement. It can be expressed as follows: where δ 0 m and δ f m are the effective displacement at damage initiation and the effective displacement at complete failure, respectively. δ max m refers to the maximum value of the effective displacement. The author has given a detailed description previously [21].

Model Validation
The author has previously used the Kristonovich-Geertsma-de Klerk (KGD) model to verify the method, so this paper verifies the accuracy of the model by comparing with the results of previous studies. Wu et al. [23] and Zhang et al. [24] established two perforations for fracturing in the horizontal wellbore. Based on the input data of Wu et al.'s [25] model, the simulation results were obtained in this paper ( Figure 1). As shown in Figure 1, our model is relatively accurate. However, when the fracture is about to propagate to the boundary, there is a slight difference between the two results, which may be caused by the boundary effect.
Energies 2021, 14, x FOR PEER REVIEW The damage evolution criterion was represented by damage variable D the effective displacement. It can be expressed as follows: where and are the effective displacement at damage initiation and t displacement at complete failure, respectively. refers to the maximum v effective displacement. The author has given a detailed description previously

Model Validation
The author has previously used the Kristonovich-Geertsma-de Klerk (K to verify the method, so this paper verifies the accuracy of the model by comp the results of previous studies. Wu et al. [23] and Zhang et al. [24] establishe forations for fracturing in the horizontal wellbore. Based on the input data of [25] model, the simulation results were obtained in this paper ( Figure 1). A Figure 1, our model is relatively accurate. However, when the fracture is abo agate to the boundary, there is a slight difference between the two results, wh caused by the boundary effect.

Numerical Model
Considering the variation of stress shadow at different perforation spac square model containing two perforations of 1 m in length was established, a Figure 2. Perforation spacing was defined as variable L. L was chosen as 2, 4, 14, 16, 18, 20, 22, and 24 m, respectively. Assuming that the model is isotro perforation spacing is the only variable, the influence of the perforation spa stress shadow is focused on.

Numerical Model
Considering the variation of stress shadow at different perforation spacing, a 50 m square model containing two perforations of 1 m in length was established, as shown in Figure 2. Perforation spacing was defined as variable L. L was chosen as 2, 4, 6,8,10,12,14,16,18,20,22, and 24 m, respectively. Assuming that the model is isotropic and the The in situ stress balance and fracturing fluid injection were achieved by setting the geostatic analysis step and soils analysis step. The fracturing time was 60 s. Additionally, the pore fluid response was set to transient consolidation. In addition, the model used aCPE4P element with a total of 10,000 in ABAQUS simulation software. The fracturing fluid was injected through the C-flow function. The simulation parameters (Table 1) are based on the data of laboratory experiments and articles [21].

Stress Distribution
The stress shadow has a great influence on the stress distribution. Therefore, the horizontal S11 stress contour was extracted under different perforation spacing, and the stress distribution characteristics of different perforation spacing under the same reservoir parameters and injection rate were analyzed. Figure 3 illustrates the horizontal stress distribution under different perforation spacing when the injection time is 10 s. According to the contour, although the injection rate and reservoir parameters are identical, there are obvious differences in the fracture path and stress distribution with dif- The in situ stress balance and fracturing fluid injection were achieved by setting the geostatic analysis step and soils analysis step. The fracturing time was 60 s. Additionally, the pore fluid response was set to transient consolidation. In addition, the model used aCPE4P element with a total of 10,000 in ABAQUS simulation software. The fracturing fluid was injected through the C-flow function. The simulation parameters (Table 1) are based on the data of laboratory experiments and articles [21].

Stress Distribution
The stress shadow has a great influence on the stress distribution. Therefore, the horizontal S 11 stress contour was extracted under different perforation spacing, and the stress distribution characteristics of different perforation spacing under the same reservoir parameters and injection rate were analyzed. Figure 3 illustrates the horizontal stress distribution under different perforation spacing when the injection time is 10 s. According to the contour, although the injection rate and reservoir parameters are identical, there are Energies 2021, 14, 4040 6 of 16 obvious differences in the fracture path and stress distribution with different perforation spacing at the same fracturing time. Under different perforation spacing, this contour shows that tension stress concentration occurs at the hydraulic fracture tips, and fracture initiate by tension. However, there is compressive stress on the sides of hydraulic fractures, owing to the stress shadow. The stress shadow zone is formed when the additional compressive stresses intersect in the middle of two hydraulic fractures. At the same injection time, with the increase in perforation spacing, the compressive stress concentration area (blue area) keeps separating and gradually gathers towards the perforation, which demonstrates that the size of the stress shadow reduces with the increase in perforation spacing.
Energies 2021, 14, x FOR PEER REVIEW 6 of 16 ferent perforation spacing at the same fracturing time. Under different perforation spacing, this contour shows that tension stress concentration occurs at the hydraulic fracture tips, and fracture initiate by tension. However, there is compressive stress on the sides of hydraulic fractures, owing to the stress shadow. The stress shadow zone is formed when the additional compressive stresses intersect in the middle of two hydraulic fractures. At the same injection time, with the increase in perforation spacing, the compressive stress concentration area (blue area) keeps separating and gradually gathers towards the perforation, which demonstrates that the size of the stress shadow reduces with the increase in perforation spacing. There is a concentrated zone of compressive stress between the two perforations, which decreases with the increase in perforation spacing. (a-l) represents perforation spacing from 2 m-24 m, respectively.

Stress at the Center of Two Perforations
According to Figure 3, it can be found that the middle of the two perforations is a stress shadow zone. Therefore, detecting the stress variation in the middle of the two perforations can intuitively indicate the change process of stress shadows. Hence, the horizontalS11 stress of the two perforation center points was extracted, and the process of stress variation at the point throughout the whole fracturing process was observed, as shown in Figure 4. According to Figure 4, it can be found that the two perforation center points are in the state of compressive stress during the whole fracturing process. Before the beginning of fracturing, the compressive stress at the center point under different perforation spacing is 7 MPa, which is the initial ground stress. As the fracturing progresses, the compressive stress at the center point increases first and then decreases. That

Stress at the Center of Two Perforations
According to Figure 3, it can be found that the middle of the two perforations is a stress shadow zone. Therefore, detecting the stress variation in the middle of the two perforations can intuitively indicate the change process of stress shadows. Hence, the horizontal S 11 stress of the two perforation center points was extracted, and the process of stress variation at the point throughout the whole fracturing process was observed, as shown in Figure 4. According to Figure 4, it can be found that the two perforation center points are in the state of compressive stress during the whole fracturing process. Before the beginning of fracturing, the compressive stress at the center point under different perforation spacing is 7 MPa, which is the initial ground stress. As the fracturing progresses, the compressive stress at the center point increases first and then decreases. That is to say, with the propagation of the fracture, the compressive stress between two perforations increases continuously until the hydraulic fracture extends to a certain length, and then the compressive stress begins to decrease. That is, the stress shadow effect begins to weaken, which indicates that the size of the stress shadow is affected by the fracture length. The size of the stress shadow changes with the propagation of hydraulic fractures. is to say, with the propagation of the fracture, the compressive stress between two perforations increases continuously until the hydraulic fracture extends to a certain length, and then the compressive stress begins to decrease. That is, the stress shadow effect begins to weaken, which indicates that the size of the stress shadow is affected by the fracture length. The size of the stress shadow changes with the propagation of hydraulic fractures. Since the stress at the center point of the two perforations constantly changes, the size of the stress shadow effect can be characterized by the magnitude of the compressive stress at this point. Therefore, the maximum value of the compressive stress at the center point was extracted during the fracturing process under different perforation spacing. Figure 5 shows the maximum value of the compressive stress and the time to reach the maximum value. With the increase in perforation spacing, the maximum compressive stress value decreases, while the time to reach the maximum compressive stress increases. This shows that the size of the stress shadow effect becomes larger with decreasing perforation spacing. When the perforation spacing is large, it takes a longer injection time to produce mutual interference between the two perforations, which is directly related to the fracture length. Since the stress at the center point of the two perforations constantly changes, the size of the stress shadow effect can be characterized by the magnitude of the compressive stress at this point. Therefore, the maximum value of the compressive stress at the center point was extracted during the fracturing process under different perforation spacing. Figure 5 shows the maximum value of the compressive stress and the time to reach the maximum value. With the increase in perforation spacing, the maximum compressive stress value decreases, while the time to reach the maximum compressive stress increases. This shows that the size of the stress shadow effect becomes larger with decreasing perforation spacing. When the perforation spacing is large, it takes a longer injection time to produce mutual interference between the two perforations, which is directly related to the fracture length. is to say, with the propagation of the fracture, the compressive stress between two perforations increases continuously until the hydraulic fracture extends to a certain length, and then the compressive stress begins to decrease. That is, the stress shadow effect begins to weaken, which indicates that the size of the stress shadow is affected by the fracture length. The size of the stress shadow changes with the propagation of hydraulic fractures. Since the stress at the center point of the two perforations constantly changes, the size of the stress shadow effect can be characterized by the magnitude of the compressive stress at this point. Therefore, the maximum value of the compressive stress at the center point was extracted during the fracturing process under different perforation spacing. Figure 5 shows the maximum value of the compressive stress and the time to reach the maximum value. With the increase in perforation spacing, the maximum compressive stress value decreases, while the time to reach the maximum compressive stress increases. This shows that the size of the stress shadow effect becomes larger with decreasing perforation spacing. When the perforation spacing is large, it takes a longer injection time to produce mutual interference between the two perforations, which is directly related to the fracture length.  When the compressive stress at the center reaches the maximum value, the stress shadow effect can be considered to be the strongest, and the fracture length at this time can be considered as the key parameter affecting the stress shadow. The ratio of fracture length to perforation spacing has an important influence on the stress shadow. Figure 6 shows the final fracture propagation path with different perforation spacing. When two perforations are fractured at the same time, the stress shadow alters the in situ stress field distribution, which leads to the two fractures diverging from the maximum horizontal stress direction extension, showing a reverse expansion trend [26]. With the increase in perforation spacing, the deflection degree of the fracture reduces and the stress shadow weakens. According to the time at which the compressive stress reaches the maximum value obtained in Figure 5, fracture paths with different perforation spacing before that time were obtained, which are shown as the grey shadows in Figure 6. When the fracture propagates to the boundary of the gray shaded area in Figure 6, the compressive stress reaches the maximum value, and the stress shadow is the strongest. Thus, it can be found that the fracture length, corresponding to the maximum compressive stress at the two perforation centers, increases as the perforation spacing increases. This also explains the fact that the maximum compressive stress value reduces with the enlargement of perforation spacing in Figure 5, and the time to reach the maximum compressive stress increases. After the fracture propagates to the boundary of the gray shaded area, the fracture spacing increases while the size of the stress shadow is gradually decreasing. Afterward, the in situ stress field again dominates thefracture propagation.

Fracture Propagation Path
When the compressive stress at the center reaches the maximum value, t shadow effect can be considered to be the strongest, and the fracture length at can be considered as the key parameter affecting the stress shadow. The ratio of length to perforation spacing has an important influence on the stress shadow. Figure 6 shows the final fracture propagation path with different perforati ing. When two perforations are fractured at the same time, the stress shadow a in situ stress field distribution, which leads to the two fractures diverging f maximum horizontal stress direction extension, showing a reverse expansio [26].With the increase in perforation spacing, the deflection degree of the fra duces and the stress shadow weakens. According to the time at which the com stress reaches the maximum value obtained in Figure 5, fracture paths with perforation spacing before that time were obtained, which are shown as the gr ows in Figure 6. When the fracture propagates to the boundary of the gray sha in Figure 6, the compressive stress reaches the maximum value, and the stress sh the strongest. Thus, it can be found that the fracture length, corresponding to t mum compressive stress at the two perforation centers, increases as the pe spacing increases. This also explains the fact that the maximum compressive stre reduces with the enlargement of perforation spacing in Figure 5, and the time the maximum compressive stress increases. After the fracture propagates to the ary of the gray shaded area, the fracture spacing increases while the size of t shadow is gradually decreasing. Afterward, the in situ stress field again domin fracture propagation.

Fracture Length
The fracture length was extracted when the fracture propagates to the boundary of the gray shaded area in Figure 6. Since the stress shadow effect is the strongest at this time, it is crucial to study the relationship between the fracture length and the perforation spacing to reveal the stress shadow effect. Figure 7 shows the length of fracture with different perforation spacing at this time. The fracture length is positively related to the perforation spacing when the stress shadow effect is strongest. The ratio of fracture length to perforation spacing was defined as R. The R-value under different perforation spacing was obtained. When the stress shadow effect is the strongest, the R-value under different perforation spacing is distributed between 2 and 4, with an average value of 2.5. Therefore, a conclusion can be drawn that the maximum size of the stress shadow effect is obtained when the fracture length is 2.5 times the length of the perforation spacing.
Energies 2021, 14, x FOR PEER REVIEW

Fracture Length
The fracture length was extracted when the fracture propagates to the bou the gray shaded area in Figure 6. Since the stress shadow effect is the stronge time, it is crucial to study the relationship between the fracture length and the tion spacing to reveal the stress shadow effect. Figure 7 shows the length of fracture with different perforation spacing at t The fracture length is positively related to the perforation spacing when th shadow effect is strongest. The ratio of fracture length to perforation spacing fined as R. The R-value under different perforation spacing was obtained. W stress shadow effect is the strongest, the R-value under different perforation s distributed between 2 and 4, with an average value of 2.5. Therefore, a conclusio drawn that the maximum size of the stress shadow effect is obtained when the length is 2.5 times the length of the perforation spacing.

Discussion
In this work, the variation process of stress shadow under different pe spacing is deeply analyzed and the compressive stress value of the two perfora ter points directly represents the change in the stress shadow effect. It can be mo tive to study the change rule of the stress shadows effect with the fracture pro in the whole fracturing process. At the same time, the ratio of perforation sp fracture length is proposed when the stress shadow effect is strongest, which h significance for guiding fracturing design.
According to Figure 6, the two hydraulic fractures are separated from each the reservoir area between the two perforations is never fractured. Accordin feature, the fracturing design can be optimized for step fracturing or re-fractu reservoir to increase fracture connectivity [3]. According to the perforation dis in Figure 8, the fracturing sequence can be designed, first fracturing ①③⑤, a fracturing ②④. After the completion of fracturing work, it can be re-perforate middle of the two perforations to achieve re-fracturing, to increase the distribu of the fracture network and optimize fracturing. R represents the ratio of fracture length to perforation spacing. R has a strong concentration, mainly between 2 and 4.

Discussion
In this work, the variation process of stress shadow under different perforation spacing is deeply analyzed and the compressive stress value of the two perforation center points directly represents the change in the stress shadow effect. It can be more intuitive to study the change rule of the stress shadows effect with the fracture propagation in the whole fracturing process. At the same time, the ratio of perforation spacing to fracture length is proposed when the stress shadow effect is strongest, which has great significance for guiding fracturing design.
According to Figure 6, the two hydraulic fractures are separated from each other, so the reservoir area between the two perforations is never fractured. According to this feature, the fracturing design can be optimized for step fracturing or re-fracturing the reservoir to increase fracture connectivity [3]. According to the perforation distribution in Figure 8, the fracturing sequence can be designed, first fracturing 1 3 5 , and then fracturing 2 4 . After the completion of fracturing work, it can be re-perforated in the middle of the two perforations to achieve re-fracturing, to increase the distribution area of the fracture network and optimize fracturing. According to the fracturing design method proposed above, XFEM was us verify it. Five initial fractures were arranged as perforations and numbered (Figur First fracturing ①③⑤, the first fracture path was obtained (Figure 9b). Under the of stress shadow, fractures ① and ⑤ expand backward, while fracture ③ exp along the direction of maximum principal stress. The middle of the fracture is alway unfractured zone. Therefore, after the completion of the first fracturing, the perfora ② and ④ can be re-fractured to obtain the fracture path after the second fractu ( Figure 9c). During the second fracturing, the hydraulic fracture path deflects unde influence of the pre-existing hydraulic fractures. Both fractures ② and ④ will d towards fracture ③, and eventually converge to form a fracture network. Meanw the secondary fracturing can activate the primary fracturing silent fracture, so tha fracture of the first fracturing can be reactivated and expanded, further increasin fracture distribution area. In addition, after the second fracturing, there is still an fractured zone in the middle of the two perforations. Therefore, a new perforation c inserted between the two perforations to perform multi-batch fracturing, althoug number of fracturing processes in multi-batch fracturing needs to be studied fu This fracturing design scheme can effectively obtain a complex fracture network improve the production rate. However, a higher pump pressure is required for secondary fracturing. The i tion pressure of the above model is ② ④ > ③ > ① ⑤.Therefore, the subsequent initiation is more difficult due to the effect of stress shadow, especially the superpos of stress shadow of multiple cracks. According to the research results, the optima foration spacing can be determined in fracturing design. Firstly, the two perforation fractured at the same time with a small pumping rate. After the fracture stops expan a large pumping rate is carried out in the middle of the two perforations to generate According to the fracturing design method proposed above, XFEM was used to verify it. Five initial fractures were arranged as perforations and numbered (Figure 9a). First fracturing 1 3 5 , the first fracture path was obtained (Figure 9b). Under the effect of stress shadow, fractures 1 and 5 expand backward, while fracture 3 expands along the direction of maximum principal stress. The middle of the fracture is always the unfractured zone. Therefore, after the completion of the first fracturing, the perforations 2 and 4 can be re-fractured to obtain the fracture path after the second fracturing (Figure 9c). During the second fracturing, the hydraulic fracture path deflects under the influence of the preexisting hydraulic fractures. Both fractures 2 and 4 will deflect towards fracture 3 , and eventually converge to form a fracture network. Meanwhile, the secondary fracturing can activate the primary fracturing silent fracture, so that the fracture of the first fracturing can be reactivated and expanded, further increasing the fracture distribution area. In addition, after the second fracturing, there is still an unfractured zone in the middle of the two perforations. Therefore, a new perforation can be inserted between the two perforations to perform multi-batch fracturing, although the number of fracturing processes in multi-batch fracturing needs to be studied further. This fracturing design scheme can effectively obtain a complex fracture network and improve the production rate. According to the fracturing design method proposed above, XFEM was used to verify it. Five initial fractures were arranged as perforations and numbered (Figure 9a) First fracturing ①③⑤, the first fracture path was obtained (Figure 9b). Under the effec of stress shadow, fractures ① and ⑤ expand backward, while fracture ③ expands along the direction of maximum principal stress. The middle of the fracture is always the unfractured zone. Therefore, after the completion of the first fracturing, the perforations ② and ④ can be re-fractured to obtain the fracture path after the second fracturing ( Figure 9c). During the second fracturing, the hydraulic fracture path deflects under the influence of the pre-existing hydraulic fractures. Both fractures ② and ④ will deflec towards fracture ③, and eventually converge to form a fracture network. Meanwhile the secondary fracturing can activate the primary fracturing silent fracture, so that the fracture of the first fracturing can be reactivated and expanded, further increasing the fracture distribution area. In addition, after the second fracturing, there is still an un fractured zone in the middle of the two perforations. Therefore, a new perforation can be inserted between the two perforations to perform multi-batch fracturing, although the number of fracturing processes in multi-batch fracturing needs to be studied further This fracturing design scheme can effectively obtain a complex fracture network and improve the production rate. However, a higher pump pressure is required for secondary fracturing. The initia tion pressure of the above model is ② ④ > ③ > ① ⑤.Therefore, the subsequent crack initiation is more difficult due to the effect of stress shadow, especially the superposition of stress shadow of multiple cracks. According to the research results, the optimal per foration spacing can be determined in fracturing design. Firstly, the two perforations are fractured at the same time with a small pumping rate. After the fracture stops expanding a large pumping rate is carried out in the middle of the two perforations to generate new fractures and reactivate the existing fractures to increase the fracturing area. Since the proppant is not considered in this model, the competition between perforations can be characterized by the fracture width ( Figure 10). However, a higher pump pressure is required for secondary fracturing. The initiation pressure of the above model is 2 4 > 3 > 1 5 . Therefore, the subsequent crack initiation is more difficult due to the effect of stress shadow, especially the superposition of stress shadow of multiple cracks. According to the research results, the optimal perforation spacing can be determined in fracturing design. Firstly, the two perforations are fractured at the same time with a small pumping rate. After the fracture stops expanding, a large pumping rate is carried out in the middle of the two perforations to generate new fractures and reactivate the existing fractures to increase the fracturing area. Since the proppant is During the first fracture, perforations ① and ⑤ inhibited th tion ③, resulting in a smaller fracture width. With the developm turing, the growth of the second hydraulic fracture has a great im the pre-existing fracture, and the pre-existing fracture tends to clo requirement for proppant. After the second fracturing, it can be div according to the width of the fracture. The first stage is a compet ondary fracture initiation and increasing fracture width lea pre-existing fractures. Perforation ③ has a significant tendency to the fracture decreases rapidly, thus requiring higher proppant st stage, perforation ③ was completely closed. In the third stage, perforations ① and ⑤ is competitive with that of secondary frac pant is subjected to fatigue loading. As a result, higher proppant s the initial fracture to resist the impact of subsequent fractures on th Many scholars have proved that the propagation of hydraulic the distribution of in situ stress fields. Roussel et al. [3] discussed reversal zone and proposed that it was not suitable for re-perfo whereas the stress reversal caused by the stress shadow also ha shown in Figure 11a, the in situ stress field between the two fractu tical to horizontal after fracturing. Therefore, when refracturing tions, the new perforation can be set to be angled with the exist Figure 10. Fracture width. There is a competitive tendency between perforation first and subsequent fracturing.
During the first fracture, perforations 1 and 5 inhibited the opening of perforation 3 , resulting in a smaller fracture width. With the development of secondary fracturing, the growth of the second hydraulic fracture has a great impact on the width of the pre-existing fracture, and the pre-existing fracture tends to close, so there is a higher requirement for proppant. After the second fracturing, it can be divided into three stages according to the width of the fracture. The first stage is a competitive stage, where secondary fracture initiation and increasing fracture width lead to the closure of pre-existing fractures. Perforation 3 has a significant tendency to close and the width of the fracture decreases rapidly, thus requiring higher proppant strength. In the second stage, perforation 3 was completely closed. In the third stage, the fracture width of perforations 1 and 5 is competitive with that of secondary fractures, where the proppant is subjected to fatigue loading. As a result, higher proppant strength is required for the initial fracture to resist the impact of subsequent fractures on the previous fracture.
Many scholars have proved that the propagation of hydraulic fractures will change the distribution of in situ stress fields. Roussel et al. [3] discussed the range of the stress reversal zone and proposed that it was not suitable for re-perforation in this region, whereas the stress reversal caused by the stress shadow also has a positive effect. As shown in Figure 11a, the in situ stress field between the two fractures changed from vertical to horizontal after fracturing. Therefore, when refracturing between two perforations, the new perforation can be set to be angled with the existing perforation. Especially for layered rock masses, the maximum principal stress is perpendicular to the bedding plane (in the Barnett shale). It is generally believed that an angle of 30 degrees between perforation and maximum principal stress can obtain a better fracturing effect. The secondary in situ stress field between the perforations after fracturing is roughly orthogonal to the perforations (Figure 11b), so the parallel of the perforation angle and in situ stress field will have a better Fracture initiation pressure is an essential key parameter in the field fracturing construction [27], which has an important influence on the selection of construction equipment and the control of injection rate in hydraulic fracturing design [28]. Therefore, the study of fracture initiation pressure is crucial for oil and gas production. Existing studies show that fracture trace angle [29], port number [30], the pressurization rate [31] and viscosity of fracturing fluid [32] all affect fracture initiation pressure. However, there are few studies on the influence of perforation spacing on initiation pressure. Therefore, in this work, the initiation pressure was extracted, and the effect of different perforation spacing on reservoir initiation pressure was obtained ( Figure 12). With the decrease in perforation spacing, the initiation pressure increases. Thus is because, when the perforation spacing is small, the fracture initiation must overcome both the in situ stress field and the secondary stress field induced by the stress shadow, which results in larger initiation pressure. Furthermore, when the perforation spacing increases to a particular value, the stress shadow decreases and the initiation mainly overcomes the in situ stress and tensile strength. Additionally, the in situ stress and strength parameters are fixed, so the initiation pressure no longer changes. However, it is interesting that the initiation pressure does not decrease monotonically with the perforation spacing. Before the perforation spacing reaches a certain fixed value, there is already a perforation spacing which can minimize the initiation pressure (the simulation results of this paper show that the minimum initiation pressure is found when the perforation spacing is 10 m). Therefore, if we can find the perforation spacing with the minimum fracturing pressure before the initiation pressure is constant. Then, under this perforation spacing, the expansion of two hydraulic fractures will be affected by the stress shadow, resulting in a backward expansion phenomenon, increasing the expansion area of the hydraulic fracture. Meanwhile, there is a small initiation pressure. This takes advantage of the positive effects of stress shadows. Fracture initiation pressure is an essential key parameter in the field fracturing construction [27], which has an important influence on the selection of construction equipment and the control of injection rate in hydraulic fracturing design [28]. Therefore, the study of fracture initiation pressure is crucial for oil and gas production. Existing studies show that fracture trace angle [29], port number [30], the pressurization rate [31] and viscosity of fracturing fluid [32] all affect fracture initiation pressure. However, there are few studies on the influence of perforation spacing on initiation pressure. Therefore, in this work, the initiation pressure was extracted, and the effect of different perforation spacing on reservoir initiation pressure was obtained ( Figure 12). With the decrease in perforation spacing, the initiation pressure increases. Thus is because, when the perforation spacing is small, the fracture initiation must overcome both the in situ stress field and the secondary stress field induced by the stress shadow, which results in larger initiation pressure. Furthermore, when the perforation spacing increases to a particular value, the stress shadow decreases and the initiation mainly overcomes the in situ stress and tensile strength. Additionally, the in situ stress and strength parameters are fixed, so the initiation pressure no longer changes. However, it is interesting that the initiation pressure does not decrease monotonically with the perforation spacing. Before the perforation spacing reaches a certain fixed value, there is already a perforation spacing which can minimize the initiation pressure (the simulation results of this paper show that the minimum initiation pressure is found when the perforation spacing is 10 m). Therefore, if we can find the perforation spacing with the minimum fracturing pressure before the initiation pressure is constant. Then, under this perforation spacing, the expansion of two hydraulic fractures will be affected by the stress shadow, resulting in a backward expansion phenomenon, increasing the expansion area of the hydraulic fracture. Meanwhile, there is a small initiation pressure. This takes advantage of the positive effects of stress shadows. The heterogeneity of shale is not considered in this study. However, sh sedimentary rock, is not only inhomogeneous in terms of mineral particles but erogeneous in terms of laminar structure [33]. Both affect the fracture propagat [34][35][36]. Han et al. [21] proposed a method to create a mineral heterogeneity mo ing this model, we compared the effects of isotropy and mineral heterogeneity ture propagation paths ( Figure 13). Mineral heterogeneity determines the initi quence of fractures and increases the roughness of fractures. Of course, more im sults need to be studied. Additionally, modeling methods for heterogeneous m with complex microstructure can also be used [37]. Meanwhile, natural fractures control the permeability of rocks [38,39]. T the natural fracture has become an indispensable part of hydraulic fracturing Based on this, a discrete fracture network (DFN) model was proposed to charact The heterogeneity of shale is not considered in this study. However, shale, as a sedimentary rock, is not only inhomogeneous in terms of mineral particles but also heterogeneous in terms of laminar structure [33]. Both affect the fracture propagation path [34][35][36]. Han et al. [21] proposed a method to create a mineral heterogeneity model. Using this model, we compared the effects of isotropy and mineral heterogeneity on fracture propagation paths ( Figure 13). Mineral heterogeneity determines the initiation sequence of fractures and increases the roughness of fractures. Of course, more impact results need to be studied. Additionally, modeling methods for heterogeneous materials with complex microstructure can also be used [37]. The heterogeneity of shale is not considered in this study. However, shale, a sedimentary rock, is not only inhomogeneous in terms of mineral particles but also h erogeneous in terms of laminar structure [33]. Both affect the fracture propagation pa [34][35][36]. Han et al. [21] proposed a method to create a mineral heterogeneity model. U ing this model, we compared the effects of isotropy and mineral heterogeneity on fr ture propagation paths ( Figure 13). Mineral heterogeneity determines the initiation quence of fractures and increases the roughness of fractures. Of course, more impact sults need to be studied. Additionally, modeling methods for heterogeneous materi with complex microstructure can also be used [37]. Meanwhile, natural fractures control the permeability of rocks [38,39]. Therefo the natural fracture has become an indispensable part of hydraulic fracturing resear Based on this, a discrete fracture network (DFN) model was proposed to characterize distribution of natural fractures [40]. A large number of scholars have researched the fect of natural fractures on hydraulic fracture propagation based on the DFN model a put forward three interaction modes [41,42]. The current research only considers the fects of natural fracturing or mineral heterogeneity on hydraulic fracturing, lacking combination of both. Therefore, the next study will consider the hydraulic fractu Meanwhile, natural fractures control the permeability of rocks [38,39]. Therefore, the natural fracture has become an indispensable part of hydraulic fracturing research. Based on this, a discrete fracture network (DFN) model was proposed to characterize the distribution of natural fractures [40]. A large number of scholars have researched the effect of natural fractures on hydraulic fracture propagation based on the DFN model and put forward three interaction modes [41,42]. The current research only considers the effects of natural fracturing or mineral heterogeneity on hydraulic fracturing, lacking a combination of both. Therefore, the next study will consider the hydraulic fracture propagation characteristics under a combination of mineral heterogeneity and natural fracture.

Conclusions
In this work, the variation process of stress shadow under different perforation spacing was analyzed by the XFEM. According to the stress distribution between two perforations, the compressive stress values at the center of two perforations were extracted. The change process of the stress shadow effect was quantitatively analyzed according to the stress value. In addition, the effect of the relationship between fracture length and perforation spacing on stress shadows was examined. The study reached the following conclusions: 1.
The rock reservoir between the two perforations is always in a state of compressive stress. Additionally, the perforation spacing is small while the compressive stress value is large. The stress shadow effect can be quantitatively analyzed by using the compressive stress values at the two perforation centers.

2.
The size of the stress shadow varies with the fracturing process. Moreover, the size is more likely to be maximized at smaller perforation spacing.

3.
As the two fractures deviate from each other, fracture spacing increases continuously.
When it enlarges to a particular value, the stress shadow effect begins to weaken until the in situ stress becomes the dominant factor again. 4.
The stress shadow effect will be strongest if the fracture length is 2.5 times the length of the perforation spacing.

5.
There is a perforation spacing that minimizes the initiation pressure, so this spacing can be used as a fracture design value. Furthermore, the positive effects of the stress shadow can be exploited by changing the fracturing sequence and controlling the perforation angle, thus increasing the fracture area and productivity.