Analysis of Interaction Between Hydraulic and Natural Fractures

Abstract

The behavior of natural fractures at the hydraulic fracturing (HF) treatment is one of the most important considerations in increasing the production from this kind of reservoirs. Therefore, considering the interaction between the natural fractures and hydraulic fractures can have great impact on the analysis and design of fracturing process. Due to the existence of such natural fractures, the perturbation stress regime around the tip of hydraulic fracture leads to some deviation in the propagation of path of hydraulic fracture. Increasing the ratio of transverse stress to the interaction stress results in a reduction in the deviation of hydraulic fracturing propagation trajectory in the vicinity of natural fracture. In this study, we modeled a hydraulic fracture with the extended finite element method (XFEM) using a cohesive-zone technique. The XFEM is used to discrete the equations, allowing for the simulation of induced fracture propagation; no re-meshing of domain is required to model the interaction between hydraulic and natural fractures. XFEM results reveal that the distance and angle of natural fracture with respect to the hydraulic fracture have a direct impact on the magnitude of tensile and shear debonding. The possibility of intersection of natural fracture by the hydraulic fracture will increase with increasing the deviation angle value. At the approaching stage of hydraulic fracture to the natural fracture, hydraulic fracture tip exerts remote compressional and tensile stress on the interface of the natural fracture, which leads to the activation and separation of natural fracture walls.

1. Introduction

Hydraulic fracture plays different roles in naturally fractured (NF) reservoirs compared with non-fractured reservoirs. Hydraulic fracture propagation in fractured porous media leads to an alternation next to all mechanisms of fracture and propagation characteristics because of the interaction between the weak interface of the natural fracture. The behavior of hydraulic fracture at the collision stage to the natural fracture may lead to the intersection, diversion, or containment. The hydraulic fracture interaction at the interface of the natural fracture is an important factor on the fracture propagation further direction in rock. Natural fracture orientation has to be determined before the operation. Many parameters can influence on the properties of crack at the fracturing process such as pore pressure of the reservoir, rock and fluid properties, state of stress, and many other factors. The coupling of hydraulic fracturing and complex fracture network deformation plays a key role in the naturally fractured porous media.

The creation of complex fracture networks in the naturally fractured media depends on the mechanism of interaction between the induced fracture and the preexisting fracture. Many authors have investigated the effect of NF on the geometry and propagation of induced fracture and some solutions have been provided by most of them for predicting the interaction mechanism [1–4]. Hubbert and Willis specified the minimum fracturing pressure with respect to the state of stress concentration around the borehole [5]. The difference of pore pressure and fracturing fluid pressure effects on the fracturing direction initiation [6]. Other research established correlations between the magnitude of horizontal stress and natural fracture characteristic in the generation of complexity network during fracturing operation [7]. The displacement of the adjacent blocks in edge-to-edge contact in shearing and slipping mode is a significant parameter in hydraulic fracturing propagation in naturally fractured porous media [8]. The extended finite element method (XFEM) is a novel technique for tracking the fracture propagation in naturally fractured reservoirs, which have been extensively used in order to facilitate fracturing interaction mechanism as no–re-meshing domain [9]. Belytschko and Zi used the extended finite element method (XFEM) and cohesive modeling to model crack-propagation paths by the division of the crack tip into the cracked and uncracked regions [10]. Taheri Shakib et al. showed that hydraulic and natural fracture characteristics and situations affect the production rate of fractured reservoirs [11]. They also showed the effect of horizontal stress orientation in stochastic fracture distributed at the hydraulic fracture operation [12]. Also, the interaction scenario of hydraulic fracture propagation in orthogonal and non-orthogonal approaching angle has been investigated [13].

The present paper aims to model the propagation of hydraulic fracture in the naturally fractured reservoirs by the implementation of XFEM. The governing equation of XFEM has been described in this paper. We represent the propagation of hydraulic fracture and the interaction between the hydraulic and natural fractures by XFEM. The research results will provide a theoretical and industrial basis for the application of hydraulic fracturing technology in the effective development of naturally fractured reservoirs.

2. XFEM as a component of interaction between fractures engineering

The methodology of extended finite element method was first proposed by Belytschko and Black to simulate plane-strain fracture propagation problems and modeling discontinuities by using enriched function with the degree of freedom. In extended finite element method, the re-meshing technique in order to track hydraulic fracture propagation and capture the evolved fracture surface is eliminated [14]. The extended finite element form in order to compute the displacement field can be expressed by the following expression:

$u{\left(x\right)}_{xfem}={\displaystyle {\displaystyle \sum }_{i\in \Omega }}{N}_i\left(x\right)u_i+{\displaystyle {\displaystyle \sum }_{j\in {\Omega }_{\tau }}}{{\displaystyle \widehat{N}}}_j\left(x\right)a_{j}H\left(x\right)+{\displaystyle {\displaystyle \sum }_{k\in {\Omega }_t}}{N}_k\left(x\right){\displaystyle {\displaystyle \sum }_{l=1}^4}{F}_l\left(r,\theta \right)b_k^l$E1

where Ω is the set of all nodes, Ω_τ is the nodes intersected by the fracture, Ω_t are the nodes contained and intersected by the tip of the fracture, ${\displaystyle \widehat{N}}\left(x\right)$ and N(x) are the nodal shape function to capture jump across the discontinuities; u_i is the nodal displacement vector of degree of freedom, F_l(r, θ) are the asymptotic function of fracture tip which is used to capture singularity of the strain around the hydraulic fracture tip, a_j and $b_k^l$ are the additional enriched degree of freedom, and H(x) is the Heaviside jump function, which can be expressed as the following form:

$H\left(X\right)=\{\begin{array}{c}1\left(x-x^{}\right).n\ge 0\left(\text{Abovethecrack}\right)\\ -1\left(x-x^{}\right).n<0\text{Otherwise}\left(\text{Belowthecrack}\right)\end{array}$E2

where x is the sample point on the fracture, x is the closest located on the fracture, and n is the unit normal to the fracture. The asymptotic function to model the displacement field at the tip of the fracture can be approximated by

$F_l\left(r,\theta \right)=\left[\sqrt{r}\sin \frac{\theta }{2},\sqrt{r}\cos \frac{\theta }{2},\sqrt{r}\sin \theta \sin \frac{\theta }{2},\sqrt{r}\sin \theta \cos \frac{\theta }{2}\right]$E3

Here, (r, θ) are the polar coordinate locations at the fracture tip.

In this study, fracture modeling is carried out using cohesive behavior at the crack tip. Barenblatt modeled cracking as a cohesive behavior in a model that predicted a nonlinear zone at the crack tip to overcome the limitation within Griffith’s theory. This model can estimate the uncracked structure behavior which is a defect in many other models. Moreover, the cohesive model does not regard singularities in stress behavior as necessary, and removes them from the initial consideration which is a great advantage [15]. In this model, the cohesive crack zone is specified by the relation between the displacement of the fracture face and the cohesion stress applied to the interface. By assuming singular crack propagation within a fractured medium and an advancement of the crack at the tip of the hydraulic fracture, cohesive modeling can be used to calculate nonlinear fracture behavior. The criterion for a fracture to propagate at the cohesive zone is that the energy-release rate must overcome the dissipation-energy rate [16, 17]. Assuming the cohesive zone within the propagated hydraulic fracture, three distinct zones will contribute in the fracturing stage, which are the fully opened zone, partially damaged zone and non-damaged zone.

Figure 1.

The elastic model of mass rock applied by the dynamic modulus.

The fully opened zone is the section that fully separates the upper and lower parts of the crack from the fluid flow. The partially damaged zone or process zone is located around the crack tip where the total stress acts to this zone lower than the critical stress (Figure 1). In addition to the two mentioned zones, non-damaged zone is located beyond the process zone and with no possibility of fracturing fluid.

The tension at the cohesive fracture tip ${\tau }_f$ can be expressed by the cohesion law [18–20]:

${\tau }_f=\frac{\partial \psi }{\partial x^{\prime }}$E4

Here,

is the displacement jump and ψ is the energy density.

The hydraulic fracture propagation in cohesive zone model can be applied by the traction-separation law:

${\left\{\frac{t_n}{t_n^o}\right\}}^2+{\left\{\frac{{\tau }_s}{{\tau }_s^o}\right\}}^2+{\left\{\frac{{\tau }_t}{{\tau }_t^o}\right\}}^2=1$E5

where t_n is the normal traction component, τ_s and τ_t are the shear and tangential traction component, respectively, $t_n^0$ is the peak value of tensile strength of the interface, ${\tau }_s^o$ and ${\tau }_t^o$ are the shear and tangential strength of the rock interface. Macaulay bracket < > represents a pure compressive deformation. Stress component is influenced by the evolution of the damage $(D^{\prime })$ which can be expressed by

$t_I=\{\begin{array}{c}\left(1-D^{\prime }\right){\displaystyle \acute{t_I}},{\displaystyle \acute{t_I}}\ge 0\\ {\displaystyle \acute{t_I}},\text{nodamage}\end{array}I=\text{normal},\text{shear},\text{and tangential}$E6

where

is the stress component without any damage initiation. Damage factor increases from 0 to 1 for non-damage to the fully cracked one. Damage variable developing linearly can be approximated by the following equation:

$D^{\prime }=\frac{{\delta }_r^f\left({\delta }_r^{\max }-{\delta }_r^0\right)}{{\delta }_r^{\max }\left({\delta }_r^f-{\delta }_r^0\right)}$E7

Here, ${\delta }_r^0,{\delta }_r^f$, and ${\delta }_r^{\max }$ represent the initial, final, and the maximum relative displacement, respectively. Total displacement in mixed-mode stage (δ_m) can be obtained by

${\delta }_m=\sqrt{{\delta }_n{}^2+{\delta }_s{}^2+{\delta }_t{}^2}$E8

3. Interaction between induced and natural fractures

Hydraulic fracture in naturally fractured reservoirs is faced with a unique situation which may increase the possibility of deviation from symmetrical propagation. Experimental results reveal that three scenarios may occur at the propagation stage and beyond the collision stage of fluid-driven in hydraulic fracture interaction with the natural fracture, namely diversion, penetration, and containment. Diversion is the situation in which the collided hydraulic fracture has an effective stress too low to initiate new fracture at the front wall of preexisting joint, and as a result the fluid-driven propagates along the natural fracture axis. Many studies have been investigated in order to specify the possibility of occurrence of these scenarios.

Hanson et al. and later Shaffer et al. represented that the magnitude of difference between the young modulus of the two intersected interface has significant influence on increasing the possibility of arresting hydraulic fracture [21, 22]. Based on their experimental reports, as the hydraulic fracture propagates from higher modulus into lower interface, the arresting phenomena increase. In addition to the young modulus, experimental results and numerical analysis reveal the effect of the frictional coefficient on the containment of hydraulic fracture. These results show that if the hydraulic fracture propagates from higher frictional coefficient pathway and collides to lower frictional coefficient interface at the natural fracture, the strain increases parallel to the hydraulic fracture due to the increase in the motion rate at interface region. This increase may result in an abrupt fracture seizing. Daneshy also discussed about the possibility of seizing the growth of hydraulic fracture at the intersection stage based on the opening interface of the natural fracture [23]. Another significant parameter that can influence the crossing criteria of hydraulic fracture is the approaching angle. Blanton using different angle-approaching experiments concluded that the presence of high differential stress and high intersection angle can improve the crossing of hydraulic fracture.

The hydraulic fracture can keep on planar propagation beyond the collision point. However, because of the energy dissipation at the contacting stage, the crossing criteria cannot exactly determine if the hydraulic fracture will penetrate through the other side of weakness plane. The fluid-driven energy must be high enough in order to separate the natural fracture bonding at the intact side of the wall. However, breakage at the other side of the wall might have some offset with the collision point, which originates from the preexisted flaw or mini-cracks along the intact side. Based on Blanton’s results, the reduction of the stress anisotropy and treatment pressure may lead to increase in the possibility of diversion and dissipation of fluid-driven along the natural fracture path and also to complex natural fracture network [24, 25]. Later, Beugelsdijk using laboratory experimental results concluded that at high principle stress difference, the hydraulic fracture may have no interaction with the preexisting discontinuities and may turn around them [26]. In addition to the mentioned scenarios, hydraulic fracture may also cause dilation, long slippage along the natural fracture interface, or may turn around and bypass discontinuities. Inclined weakness plane at the propagation path of induced fracture has high tendency to divert the fluid-driven. However, all of the mentioned scenarios can only be estimated and visually represented using an experimental method. The containment stage is the only stage which can approximate the interaction on the natural fracture and fluid-driven. Beyond this stage, no other method can exactly approve the crossing criteria or diversion.

Hydraulic fracture propagation in the naturally fractured reservoirs plays a different role than the conventional porous media. As the hydraulic fracture passes beyond the induced stress of drilled well, the hydraulic fracture propagation reorientates through the maximum stress principle. The hydraulic fracture propagation in homogeneous porous media is approximately near to the straight path; however, in a real reservoir rock media, because of discontinuities and inhomogeneity, the induced fracture trajectory waver is perpendicular with the minimum compressional stress. The hydraulic fracture tip tends to propagate through the local direction, which has the maximum energy release rate and minimum resistance. Still, there is the possibility of curving and increasing the deviation of hydraulic fracture from straight trajectory by increasing the shearing intensity factor. As long as the induced fracture propagates in opening mode, its fracture trajectory is near to the straight line. When the fracture faced the two materials with different Young’s modulus, the angle of deflection tends to rematch the tip direction in accordance with the lower Young’s modulus material. By increasing the hydraulic fracture length by the propagation of the tip of the hydraulic fracture away from the wellbore, the curvature of hydraulic fracture tends to be decreased. In addition to the rock mechanic properties, the fracturing fluid properties and flow rate injection also have a great impact on the straightness stability. Also, increasing the fracturing fluid viscosity will decrease the leak-off rate and tortuosity of the fracture, but it requires a higher rate of treatment pressure [27]. However, increasing the fluid viscosity in fracturing treatment leads to an abrupt increase in fluid pressure at the fracture path and reduces the flow rate at the fracture tip, because of the uniformity in pressure profile within the hydraulic fracture path. High rate of pressure difference between the fracture tip and the mouth region causes an inhomogeneity in the geometry of the fracture path and lowers the rate of growth [28]. Unlike the high viscosity, lower viscosity will cause a uniform pressure profile within the hydraulic fracture path increasing fluid leakage rate to the adjacent layer. Increasing the fluid leak-off rate will cause a perturbation in the local stress regime and increase the possibility of zigzag fracture pattern. Natural fractures have different response in alteration of the rate of injection and fracturing fluid properties. In the naturally fractured reservoir, increasing the flow rate injection will increase the leak-off rate to the adjacent layer and subsequently cause debonding of the natural fracture in tensile mode [29]. From the studies, reducing the fluid flow injection rate and viscosity of fracturing fluid in fractured media will greatly reduce the possibility of complex fracture network generation [30].

4. Coalescence of hydraulic and natural fractures

After initiation and propagation stage of hydraulic fracture beyond the far-field stress region, the hydraulic fracture tries to rematch its orientation by the maximum stress principle. The hydraulic fracture direction is almost parallel with the orientation of maximum stress principle but not exactly perpendicular to the minimum compressional stress, because it tends to orient its trajectory in porous media along the path of minimum resistance. Despite the stress direction in the local field, the induced fracture trajectory may have a wavy shape because of the inhomogeneity of the porous media along its path. The local stress component at the neighborhood of the fracture tip can be expressed by the following equation:

$\begin{array}{c}{\sigma }_r=\frac{1}{{\left(2r\right)}^{\frac{1}{2}}}\cos \frac{\alpha }{2}\left[K_I\left(1+{\sin }^2\frac{\alpha }{2}\right)+\frac{3}{2}{K}_{II}\sin \alpha -2K_{II}\tan \left(\frac{\alpha }{2}\right)\right]\\ {\sigma }_{\alpha }=\frac{1}{{\left(2r\right)}^{\frac{1}{2}}}\cos \frac{\alpha }{2}\left[K_I\left({\cos }^2\frac{\alpha }{2}\right)-\frac{3}{2}{K}_{II}\sin \alpha \right]\\ {\tau }_{\alpha r}=\frac{1}{2{\left(2r\right)}^{\frac{1}{2}}}\cos \frac{\alpha }{2}\left[K_I\left({\cos }^2\frac{\alpha }{2}\right)+K_{II}\left(3\cos \alpha -1\right)\right]\end{array}$E9

where (r ∧ α) are the polar coordinates of the location from the crack tip, and K_I ∧ K_II are the opening and shearing mode intensity factors, respectively, which is applied to state the stress around the crack tip, for a penny-shaped crack. This crack can be expressed as

$K_I=2\sqrt{\frac{L}{\pi }}\left[P_{ff}-\left({\sigma }_1{\sin }^2\beta +{\sigma }_3{\cos }^2\beta \right)\right]$E10

$K_{II}=\sqrt{\frac{L}{\mu }}\left({\sigma }_1-{\sigma }_3\right)\sin 2\beta$E11

where P_ff is the fluid pressure within the fracture, L is half the length of the fracture, and β is the angle between the fracture and the far-field stress [31]. As the hydraulic fracture propagates all the way into the natural fracture interface, it exerts compressional and tensional stress to the natural fracture, which may lead to reactive the natural fracture in the opening or shearing mode.

In numerical modeling, we can only predict the local displacement within the natural fracture only at the coalescence level. Prior to intersection and activation stage, the natural fracture interface is approximately closed; however, the permeability of the preexisted joint is higher than the matrix element. The stress regime alteration around the induced fractured tip results from the perturbation of location by the natural fracture. Induced fracture tip at the approaching stage, exert compressional and tensional force to the natural fracture interface which lead to debonding in normal and shearing mode. Coalescence of the hydraulic fracture to the natural also may cause other results such as increasing the fluid leak-off rate around this area, lowering the fluid pressure within the induced fracture path, and changing the permeability axis around this area. At the intersection point of, we can observe a higher rate of hydraulic fracture aperture size because of the alteration and increasing the fluid at this region. In tensile-failure mode, for normal displacement to occur at the natural fracture interface, the tensile stress acting on the interval must exceed the tensile strength of the natural fracture. The effective stress (σ_e exerted on the natural fracture interface is given by

${\sigma }_e={\sigma }_n-P_p$E12

where P_p is the pore pressure. In the cohesive model, the normal displacement at the fracture interface results from the cohesive interface forces $\left(F_n^{Cs}\right)$:

$F_n^{Cs}=K_n\text{Δ}{u}_n$E13

where K_n is the normal stiffness and Δu_n is the normal displacement. In addition to tensile stress, shear stress results in the slippage of natural fracture walls, as specified by shear displacement. At the shear-slippage threshold, the shear stress acts at the natural fracture surface to dominate the shear strength of the natural fracture:

${\tau }_s\ge {\tau }_0+\eta \left({\sigma }_n-P\right)$E14

where τ_s is the shear stress, η is the friction coefficient, and τ₀ is the shear strength. The normal displacement at the crack interface caused by the shear forces $\left(F_n^{Cs}\right)$ can be defined as

$F_s^{Cs}=K_s\text{Δ}{u}_s$E15

where K_s is the shear stiffness and Δu_s is the shear displacement along the crack interface in shear slippage [32]. Normal and shear displacement both simultaneously take place at the natural fracture interface, however the shear slippage is a complex function of normal displacement which have nonlinear incremental behavior; this displacement is due to the rule that increasing the normal displacement of the natural fracture interval will decrease the natural fracture wall interaction and increase shear displacement. All of the mentioned phenomena are located at the coalescence stage, which means the fluid front at the time of collision has no invasion through the natural fracture debonding interface. At the touching time of the natural fracture by hydraulic fracture, the rate of compressional or tensioning at the natural fracture interface highly depends on the collision angle. From the result, at the inclined approaching angle, if the acting stress not enough in order causes shearing failure at the natural fracture tip, the lower part of the natural fracture will react in the negative direction way.

5. Interaction between induced fracture and natural fracture with various positions

As mentioned earlier, when the hydraulic fracture propagates through the 90° natural fracture, at the early stage of approaching, the natural fracture is almost closed. By approaching the hydraulic fracture to the natural fracture interface, some activation may occur which may change the local physical properties at that region. In addition to the hydraulic fracture acting stress, the natural fracture also perturbs the stress regime around its area, which is directly proportional to its length. In reality, we cannot represent that if the approaching angle is 90°, then the collision angle is orthogonal too. This is due to the fact that the local perturbation and acting stress in coalescence process are mutual. Natural fracture by acting stress to the tip of the hydraulic fracture will cause deviation on its overall propagation, which may lead to deviation from the 90°. The magnitude of this stress can be expressed by the following equation [33]:

$I=C_s\left({\sigma }_1+P\right)$E16

where σ₁ is the maximum principle stress, C_s is the fracture shape factor, and P is the pressure within the natural fracture.

From Figure 2, assume that the approaching angle is the same as collision angle which is 90°. As seen in Figure 2, the hydraulic fracture approaches the natural fracture in an orthogonal angle. The tensile and shear debonding can be evaluated at the approaching stage of the hydraulic fracture tip to the natural fracture interface in a, b and c areas. a and b areas are located, respectively, at 10- and 5-cm distances from the 50-cm length natural fracture interface, and c area is precisely located at the collision point of the hydraulic fracture to the natural fracture. Stress condition is assumed to be isotopic.

Figure 2.

Evaluated areas for debonding of natural fracture when induction fracture is 90° angle.

The maximum opening and shearing displacement in perpendicular approaching stage approximately occurs at the 20-cm distance from the north tip of the natural fracture. The maximum tensile and debonding size and location in the orthogonal approaching stage are the same. Moreover, debonding evaluation indicates that the minimum debonding size occurs at the 30-cm distance from the north of the natural fracture tip. As already mentioned, in the realistic-induced fracture propagation, debonding displacement alteration in tensile and shearing mode happens because changing the propagation angle at the perturbed stress region is not monotonic.

Perturbation of stress regime around the approaching hydraulic fracture tip will lead to the activation of natural fracture interface prior to the collision stage. In normal opening mode prior to the collision stage, debonding occurs at the time that the pore pressure within the natural fracture dominates the normal closure stress of the natural fracture (P > σ_n). At the “a”-approaching region, the remote stress acting on the natural fracture interface is not sufficient to cause tensile opening of the interface. As the induced fracture propagates at the perpendicular direction to the preexisting fracture interface, the remote normal opening stress increases. Tensile stress at the natural fracture interface is maximum as the tip of hydraulic fracture reaches the “c” region. In shearing slippage, the remote shearing stress, which acts to the natural fracture from the region “a”, is not sufficient to debond in shearing mode. During the induced fracture propagation and approaching the “b” and “c” zones, the tensile and shear stress acting on the natural fracture can overcome the threshold dilation stress, which leads to increasing the aperture and permeability of fractured blocks. Shear permeability of natural fracture interface will increase dramatically by approaching the acting stress to the dilation. If the pore pressure of the natural fractures cannot expose sufficient stress to dominate closure stress, the fracture interface will have displacement in shear rather than debonding in tension. At the 90° angle, the maximum debonding occurs near the collision point. Figure 3 represents the debonding size along the natural fracture and is independent of the rock tightness; the stress dominant of the induced fracture is independent of rock elastic properties. To investigate the phenomenon of debonding, when a natural fracture with a 45° angle is placed in the hydraulic fracture path as same as 90°, we assumed in three regions (Figures 3 and 4).

Figure 3.

Tensile and shear displacements along the deboned zone shown in Figure 2.

Figure 4.

Approaching stage of hydraulic fracture and shear dilation caused by remote stress around induced fracture.

Another main approaching angle, which can be investigated in our study, is an inclined natural fracture with the 45° angle with respect to the propagated hydraulic fracture. In an inclined mode, the lower rate of energy is required in order to reactivate the natural fracture interface at the same distance compared with the perpendicular mode (Figure 5). Unlike many earlier models, the hydraulic fracture is propagated through the interface of the natural fracture, which means that the touching moment of the left side is the same as the right side. Tensile and shear displacements along the debonded crack (45°) are shown in Figure 6. When the induction with 45° angle is close to the natural fracture in the c area, tensile failure phenomenon is such that the natural fracture had an angle of 90°, because the middle area of the natural fracture becomes debonded and the maximum value of debonding occurs at the collision point. But with less distance between the natural and induced fractures, the condition is slightly different. When the hydraulic fracture approaches the 10-cm distance from the natural fracture, the 12-cm distance from the north tip of the natural fracture becomes compressed and the other part becomes debonded. The maximum value of debonding is at the collision point but the symmetry of the debonding zone in the natural fractures with 90° angle does not take place here. After the cutoff point, the natural fracture by hydraulic fracture (c area) of the upper part of the kink point becomes debonded and the lower part becomes compressed (Figure 6). In 45° angle propagation angle, the shear displacement magnitude has a higher value than the tensile opening. In this case, the lower part of the coalescence point has the tendency to bind because of the compression and the upper part in tension turns into debonding (Figure 7).

Figure 5.

Areas of study for debonding investigation when natural fracture with a 45° angle relative to the hydraulic fracture spread.

Figure 6.

Tensile and shear displacements along the debonded zone shown in Figure 5.

In low approaching angle (45°) at the isotropic stress ratio, the shearing displacement is much larger than the tensile mode; however, with an increase in the stress ratio the difference between shearing and tensile opening remains closed to each other [34]. The natural fracture length increases the remote stress caused by the tip of the hydraulic fracture that has a tendency to increase the debonding of the natural fracture [35, 36].

Figure 7.

Debonding induced by the approaching hydraulic fracture to natural fracture.

The approaching stage of the hydraulic fracture was not fully investigated and carried out in a numerical way. As mentioned previously, considering stress regime perturbation around the natural fracture location will cause a deflection on the approaching angle of the hydraulic fracture. As the hydraulic fracture grows toward the natural fracture, influenced by the interaction stress of the natural fracture, the nearest tip edge will be active in a shorter time leading to the propagation of hydraulic fracture in a mixed mode. By increasing the shearing intensity factor, the hydraulic fracture path tends to be more kinked and deviates through the natural fracture interface. The following equation can compute the deflection angle of induced fracture (α) under mixed-mode propagation:

$\tan \frac{\alpha }{2}=\frac{1}{4}\left(\frac{K_I}{K_{II}}\pm \sqrt{\frac{K_I}{K_{II}}+8}\right)$E17

The curvature of the hydraulic fracture by the propagation of the hydraulic fracture will dramatically increase in stress-perturbed zone [33]. If the opening mode dominates in the tip of the hydraulic fracture, the fracture trajectory will tend to be more singular and straight. The rate of the hydraulic fracture deflection highly depends on the treatment pressure, leak-off rate, length of the natural fracture, and stress anisotropy. In this study, we assume that the hydraulic fracture is subjected to an isotropic principle stress. At the early stage of deviation, the natural fracture walls tend to stick together and are almost completely closed. In parallel natural fracture case, in addition to the distance parameter, the alteration of the approaching angle is another factor which was considered. Figure 8 shows the distance from the deviated hydraulic fracture tip on the natural fracture at 10 (a) and 5 m (b) and the exact coalescence (c) of the hydraulic and natural fractures. When the hydraulic fracture reaches the point a, the natural fracture reaches the activation threshold. When the hydraulic fracture approaches the natural fracture (Figure 8b and c), normal displacement occurs, and the natural fracture interface nearly fully separates. As seen in Figure 8, the approaching of the induced fracture will lead to an abrupt increase in the propagation angle and oriented near to perpendicularly. Increases in the values for the deviation angle and interaction stress increase the possibility of natural fracture collision.

Figure 8.

Approaching stage of induced fracture and shear dilation caused by remote stress around hydraulic fracture.

If the collision point in the approaching stage of the hydraulic fracture is assumed to lie at the midpoint of the natural fracture in the isotropic principle stress situation, the tensile displacement is as shown in Figure 9. At the approaching stage, the shear displacement increases nonlinearly because, at a constant shear stress, the shear displacement is also a function of the normal displacement. By increasing the normal displacement of natural fracture interface, the shear displacement has lower resistance to shearing. Moreover, because of continuous changing of approaching angle besides the distance, the shearing, and opening displacement both of them have non-monotonic behavior. As the hydraulic fracture approaches the natural fracture, the approaching angle of the hydraulic fracture increases with respect to the natural fracture location, which leads to a decrease in shearing compression. Surprisingly, the influence of the approaching angle on the shear slippage as the hydraulic fracture approaches the natural fracture is greater than the influence of the distance. As Figure 9 shows, the approaching angle of the hydraulic fracturing tip is 66 (a), 49 (b), and 34° (c). As seen in Figure 10, the deviation of the intersection angle from the perpendicular will result in discrepancies in the natural fracture tip displacement. As the hydraulic fracture interacts with the natural fracture, the pore pressure within the natural fracture changes, which leads to compression and extension within the natural fracture.

Figure 9.

Tensile and shear displacements along the debonded zone shown in Figure 8.

Figure 10.

Deviation of the intersection angle from the perpendicular will result in discrepancies in natural fracture tip displacement.

6. Discussions

Formerly, re-meshing technique has been greatly implemented in order to align the mesh with the tip of the hydraulic fracture for tracking the propagating direction. However, in our study by utilizing the XFEM as no-re-meshing tools can greatly track the hydraulic fracture trajectory to capture the stress and strain field around the tip of the hydraulic fracture. The accuracy of fracture propagation trajectory by refining the mesh around the crack tip can be improved. Stress singularity at the fracture tip is eliminated by the implementation of cohesive zone model in XFEM. Refining the mesh can provide more accurate calculation in the propagation of hydraulic fracture through natural fractures based on shearing or opening mode by computation stress concentration around the fracture tip. The number of iteration to reach convergence in our fracture tip is 5–7. The error between our numerical result and the analytical result is lower than 1%.

7. Conclusion

Natural fractures can have a significant effect on the hydraulic fracture growth and achieve successful treatment. Spacing and trajectory of natural fractures in fractured blocks with respect to the induced fracture propagation has a significant effect on the accuracy of interaction prediction. Numerical analysis of hydraulic fracturing propagation in the naturally fractured reservoir and the interaction between the induced fracture and the natural fracture are the main objectives of this paper. Numerical simulation can be used as a tool to solve this engineering problem.

In this paper, the extended finite element method (XFEM) has been implemented to simulate the coalescence stage of hydraulic fracture and natural fractures. Analysis of interaction between the induced and natural fractures in the fractured reservoirs was discussed in this study. The interaction between the induced and natural fractures depends on the collide angle. Induced fracture causes the opening of the preexisting natural fractures. The tensile and shear debonding of natural fractures in 90 and 45° displayed different behavior caused induced and variations in stresses at the natural fractures. A critical point in interaction between the hydraulic fracture and the natural fractures is the dilation caused by shearing and opening from the northing to the southing along the natural fracture in both degrees which play different scenarios. Decreasing the approaching angle from perpendicular to 45° intensifies the displacement by shearing much more than tensile. In low collision angle, the top stage of the interception point has the maximum debonding in shearing mode and the lower stage has the maximum bonding.