Influences of the intersection angle between interlayer and in situ stresses during hydraulic fracturing process

Influences of the intersection angle between interlayer and in situ stresses during hydraulic fracturing process

Journal of Natural Gas Science and Engineering 36 (2016) 963e985 Contents lists available at ScienceDirect Journal of Natural Gas Science and Engine...

9MB Sizes 0 Downloads 108 Views

Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Contents lists available at ScienceDirect

Journal of Natural Gas Science and Engineering journal homepage: www.elsevier.com/locate/jngse

Influences of the intersection angle between interlayer and in situ stresses during hydraulic fracturing process Quanshu Li*, Huilin Xing Centre for Geoscience Computing, School of Earth Sciences, University of Queensland, Australia

a r t i c l e i n f o

a b s t r a c t

Article history: Received 9 September 2016 Received in revised form 14 October 2016 Accepted 14 November 2016 Available online 16 November 2016

Hydraulic fracturing is commonly used as a main stimulation method in unconventional reservoirs in which interlayers are common and they may change the local stresses and permeability. Interlayers with different intersection angles with the principal stresses may have different initial stress and permeability, and the hydraulic fracture approaches the interlayer from different angles. In this research, finite element based numerical analysis of hydraulic fracturing has been carried out to study the stress and permeability variation within interlayers of four different intersection angles. The stresses and permeability variation shows that: (1) the intersection angle not only influences the stress and permeability of the interlayer, but also determines the timing for the fracture to reach and/or pass the interlayer; (2) the hydraulic fracturing process also influences the stress and permeability within the interlayer, and the higher the intersection angle is, the higher the stress-coupled permeability is. The results can help us to better understand the influence of intersection angle on the mutual effect between interlayer and hydraulic fracturing quantitatively towards an optimal design of hydraulic stimulation process. © 2016 Elsevier B.V. All rights reserved.

Keywords: Hydraulic fracturing Finite element modelling Interlayer Intersection angle Reservoir geology

1. Introduction Hydraulic fracturing is widely used in unconventional reservoir stimulations. Interlayers are common in unconventional reservoirs and represent a key challenge for hydraulic fracturing, owing to their properties, geometry, geometric distribution, and their influence on localized stress and therefore permeability evolution and distribution (Handwerger et al., 2011; Wang et al., 2012; Li et al., 2015; Raji et al., 2015). The stresses and properties of interlayers are normally different from their adjacent rocks (Chudnovsky et al., 2001). In a layered reservoir, the in situ stress difference between the layered formations may act as a barrier for hydraulic fracturing, as well as the difference in elasticity modulus. Other parameters such as fracture toughness difference, interface strength and fracturing fluid pressure distribution and rheology may also influence hydraulic fracturing (Cleary, 1980). Moreover, the difference in the properties could further influence the stress conditions in the reservoir. For example, Teufel and Clark (1981) found the difference in elastic properties between reservoir layers could influence the vertical

* Corresponding author. E-mail address: [email protected] (Q. Li). http://dx.doi.org/10.1016/j.jngse.2016.11.032 1875-5100/© 2016 Elsevier B.V. All rights reserved.

distribution of the minimum horizontal stress state, which in turn would influence the propagation of vertical growth. Nuller et al. (2001) also found that the difference in elastic properties would influence fracturing process by influencing the energy release rate. By investigating the Young's moduli contrast between the layers, Nuller et al. (2001) concluded that both the hoop stress along the crack plane within the more rigid material, and the normal tensile stress at the interface increased when the hydraulic fracture approached from the softer material. Besides the properties of the interlayers, Daneshy (1978) found that the strong bonding between layered formations would be difficult to contain a hydraulic fracture. The bearing stress of interlayers may also create different fracture types within the interlayers (Warpinski et al., 1998). Wu et al. (2004) furtherly pointed out that the hydraulic fracture deflection close to the interlayer may be caused by the different energy release rate in adjacent materials. In addition, the conductivity in adjacent layers can influence the width of a hydraulic fracture if it grows across the interfaces (Smith et al., 2001; Daneshy, 2009). The above research into the influence of interlayers on hydraulic fracturing mainly focused on the initiation and propagation of hydraulic fracture in vertical direction when it encounters different formations. There is a lack of interpretation on the mutual effect between the interlayer and hydraulic fracture. On one hand, the

964

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

interlayers may have different stress and permeability, which may influence hydraulic fracturing process; on the other hand, hydraulic fracturing procedure may influence the stress and permeability within the interlayers, which may reversely influence hydraulic fracturing itself, and the flow leakage from the interlayers. When the intersection angle of the interlayers is different, the intersection angle between the interlayer and principle stresses is different. Thus, the local stress within the interlayers under different intersection angles may be different, which further leads to a difference in the local stress-coupled permeability. Moreover, when the intersection angle is different, the hydraulic fracture may approach the interlayer with different intersection angles, and the propagation distance and influence on the interlayer may be different. In order to figure out the mutual effect between interlayers and hydraulic fracturing, finite element modelling is used in this research to study hydraulic fracture propagation with interlayers under different intersection angles, and the resulting stress and permeability variation within the interlayers. The results can help us to understand the mutual influences between interlayers and hydraulic fracturing process. 2. Methodology In-house code PANDAS (Parallel Adaptive Nonlinear Deformation Analysis Software) was used to simulate the initial and hydraulic fracturing processes. PANDAS is a finite element code developed and verified for simulating coupled fracturing processes in a complex geological setting in 3D, subjected to the affecting hydraulic and thermal factors for simulating coupled crustal dynamics (Xing and Makinouchi, 2002; Xing et al. 2006; Xing, 2014; Li and Xing, 2015), and some benchmarks have been attached to this paper. The maximum stress criterion and Mohr-Coulomb theory are used for failure of poro-elastic-brittle rock materials through the effective stress concept (Li et al., 2015). The stress dependent permeability is used as an exponential law with stress (Dong et al., 2010). Once the rock skeleton is failed, the stresses stored within the element are released, and the new fractured elements are described with the air material, i.e. the Young's modulus is reduced to a very small value, and the permeability stays at the fractured state. By using the continuum damage method, the hydraulic fracturing can be simulated in true three dimension models and the cost of the simulation is comparatively small. Focusing on the near wellbore region, the dimension of the model is 1000 mm  1000 mm (Fig. 1); the width of the thin interlayers is 10 mm. The well is 20 mm to interlayer 1 and 30 mm to interlayer 2. The compressive horizontal stress sh ¼ 5:8 MPa, and sH ¼ 2:3 MPa; and the injection flow rate is 80 ml/s, which is applied where the well is located. The material parameters are shown in Table 1. 3. Results and discussion 3.1. Hydraulic fracturing when q ¼ 0 The numerical simulation of hydraulic fracturing procedure with q ¼ 0 is shown in Figs. 2 and 3. The initial maximum principal stress in interlayer 1 is 2.34 MPa, in interlayer 2 it is 0.5 MPa, and in the reservoir it is 0.9 MPa (Fig. 2a). After the hydraulic fracture is initiated, the maximum principal stress is concentrated on the fracture tips, which is 1.2 MPa, whereas along the hydraulic fracture a compressive stress area has been created (Fig. 2b). Before the hydraulic fracture reaches the interlayers, the induced stress at the

Fig. 1. Top view of finite element model of hydraulic fracturing in a reservoir with interlayers.

Table 1 Material parameters used in numerical analysis. Parameter

Rock

Interlayer 1 (upper)

Interlayer 2 (lower)

Young's modulus (MPa) Poisson's ratio Tensile strength (MPa)

903 0.20 2.9 0.4

517 0.30 25.5 1.e5

1258 0.18 2.3 0.1

Permeability (1015 m2 )

fracture tip has influenced the stress within the interlayers, and it is more obvious in interlayer 2, where the stress has been increased to 0.3 MPa (Fig. 2b). At 86 s, the hydraulic fracture reaches interlayer 1, and due to the longer propagation distance, it reaches interlayer 2 at 100 s. After the hydraulic fracture has reached the interlayers, the induced stress from the hydraulic fracture has further influenced the stress within and beyond the interlayers (Fig. 2c): the maximum principal stress within interlayer 1 has been increased to 0.2 MPa where it is close to the fracture tip; beyond the fracture tip, the stress in interlayer 1 has been increased to 1.7 MPa, and beyond the interlayer, the stress in the reservoir has been increased to 0.7 MPa within a small area; in interlayer 2, the stress is increased to 1.2 MPa, and the stress gradually decreases to 0.9 MPa into the reservoir. At 223 s, the hydraulic fracture passes through interlayer 2, and at 493 s it propagates beyond interlayer 2 while been contained by interlayer 1 with a significant stress concentration (Fig. 2d). The stress concentration spot beyond interlayer 1 has been expanded compared to Fig. 2c, and the compressive stress area along the hydraulic fracture has been shifted and linked to interlayer 1; as in interlayer 2, beyond the fracture, the stress has been gradually increased to 0.4 MPa (Fig. 2d). The initial shear stress in the reservoir is 2.3 MPa, in interlayer 1 it is 1.8 MPa, and in interlayer 2 it is 2.6 MPa (Fig. 3a). After the hydraulic fracture is initiated, shear stress is concentrated at 3 MPa

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

965

Fig. 2. Top view of maximum principal stress (MPa) distribution during the hydraulic fracturing stages for q ¼ 0 : (a) 0 s, initial maximum principal stress; (b) 42 s, the hydraulic fracture approaching the interlayer; (c) 110 s, the hydraulic fracture crossing the interlayer; (d) 493 s, the hydraulic fracture crossed the interlayer.

at the tips of the fracture (Fig. 3b), while along the fracture, the stress is decreased to 1.7 MPa from the initial reservoir stress. The stress concentrations have reached the interlayers before the fracture reaches them (Fig. 3b). When the hydraulic fracture reaches the interlayers, the stress in interlayer 1 has been increased to 3 MPa, and in interlayer 1 it is increased to 2.1 MPa (Fig. 3c). After the hydraulic fracture passes through interlayer 2, the fracture has been contained by interlayer 1 (Fig. 3d). When the fracture is arrested by interlayer 1, the stress in interlayer 1 has been increased to 3 MPa, and the stress concentration has influenced the stress distribution beyond interlayer 1; whilst the fracture passes through interlayer 2, and the stress is decreased to the initial stress in interlayer 2 from the fracture body. However, beyond interlayer 2, the stress decrease is less intense because the reservoir has a higher permeability than interlayer 2. Due to the extremely low permeability, the pore pressure within interlayer 1 is as small as 6  1011 MPa during the hydraulic fracturing process (Fig. 4 iea), but the maximum principal stress

has been increased to 1.0 MPa at 493 s (Fig. 4 ii-a), and the shear stress has been increased to 3.5 MPa at the same time (Fig. 4 iii-a). The stress change leads to a change in permeability, which has been decreased to 1:22  1020 m2 at the central and 1:260  1020 m2 aside, from the initial 1:267  1020 m2 (Fig. 4 iv-a). For interlayer 2, the pore pressure reaches 20 MPa (Fig. 4 ieb), and the maximum principal stress has been generally increased before the fracture reaches the interlayer (Fig. 4 ii-b). Nevertheless, after interlayer 2 has been fractured, the stress is released to zero within the fractured area, while close to the fracture the stress has been increased due to the leakage from the fracture; at further area, the stress has been decreased because of the deformation caused by the fracture (Fig. 4 ii-b). The maximum shear stress has been increased within interlayer 2 before the hydraulic fracture reaches it, and it has been released within the fracture area after the fracture has passed through it, while beyond the fracture, the stress has been decreased from the initial 2.5 MPa (Fig. 4 iii-b). The permeability in interlayer 2 is initially 1:375  1016 m2, which is greatly

966

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 3. Top view of maximum shear stress (MPa) distribution during the different hydraulic fracturing stages for q ¼ 0 : (a) 0 s, initial maximum shear stress; (b) 42 s, the hydraulic fracture approaching the interlayer; (c) 90 s, the hydraulic fracture crossing the interlayer; (d) 493 s, the hydraulic fracture crossed the interlayer.

increased to 1:0  1013 within the fractured area. However, within the area where the stresses have been decreased, the permeability has been reduced to 1:35  1016 m2 (Fig. 4 iv-b). 3.2. Hydraulic fracturing when q ¼ 30 When the intersection angle of the interlayer increases from 0 to 30 , in the minimum horizontal stress direction, the distance from the well to interlayer 1 becomes 23 mm, to interlayer 2 it is 34.6 mm. Since the interlayers are not perpendicular to the propagation direction of the hydraulic fracture, the impact from the fracture is unsymmetrical on the interlayers, and the local stress around the hydraulic fracture is unevenly distributed. (Fig. 5, Fig. 6). The initial maximum principal stresses within interlayer 1 and 2 are close to each other, both are at about 1.3 MPa, and in reservoir rock it is 0.9 MPa (Fig. 5 a). After the hydraulic fracture has been created, the stress concentration from the fracture tip is influencing on the stresses within the interlayers (Fig. 5b). However, the

influence is not symmetric. In interlayer 1, on the left side, the stress is decreased to 1.7 MPa within a small region, while on the right side, there is no such area; in interlayer 2, the stress belt from 1.34 MPa to 1.14 MPa is wider on the right side (Fig. 5b). Because of the longer distance from the well to the interlayers, the hydraulic fracture reaches the interlayers later: at 105 s, the hydraulic fracture reaches interlayer 1, and at 113 s, it reaches interlayer 2, whereas when the intersection angle is 0 , the times are 86 s and 100 s respectively. After the fracture reaches the interlayers, it is more obvious that the stress is decreased in an area on the left wide in interlayer 1, and on the right side, the stress has been increased to 1.0 MPa in a small area; the compressive stress area along the hydraulic fracture has reaches interlayer 2 (Fig. 5c). Under the continuous fluid injection, the hydraulic fracture passes through interlayer 2 at 247 s, while it is contained by interlayer 1. The stress has been decreased further to 1.8 MPa on the left in interlayer 1, and the compressive stress area created along the hydraulic fracture has passed both interlayer 1 and 2, although from the left of interlayer 1 and the right of interlayer 2 (Fig. 5d). It

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

967

Fig. 4. Variations within the interlayers when q ¼ 0 , during hydraulic fracturing process: (i) pore pressure; (ii) maximum principal stress; (iii) maximum shear stress; (iv) permeability.

should be noted that the induced stress area beyond interlayer 1 from the fracture tip has been expanded over time (Fig. 5). As for the maximum shear stress, the initial maximum shear stress in the reservoir and interlayer 2 is 2.3 MPa, and in interlayer 1 is 2.1 MPa (Fig. 6a). The stress is concentrated at 2.9 MPa at the fracture tips, and decreased to 2.0 MPa along the fracture (Fig. 6b). On the right side in interlayer 1, the maximum shear stress is increased to 2.3 MPa before the fracture reaches the interlayer (Fig. 6b), then further increased to 2.9 MPa due to the stress concentration at the fracture tip (Fig. 6c); on the left side, it is further decreased to 2.0 MPa when the fracture is contained by interlayer 1 (Fig. 6d). The induced stress shadow from the fracture is gradually connected to the stress belt of interlayer 2, and passes beyond it (Fig. 6 bed). Fig. 7 further investigates the variations within the interlayers during the hydraulic fracturing process. The pore pressure in interlayer 1 is very small (6:0  1011 MPa) because of the low permeability (Fig. 7 iea), but the stress within interlayer 1 is changed significantly (Fig. 7 ii-a, iii-a). Due to the intersection angle, the maximum principal stress is gradually decreased to 1.8 MPa on the left side of the well while it is asymmetrically

increased to 0.2 MPa on the right in interlayer 1 (Fig. 7 ii-a). The shear stress has been greatly increased to 3.6 MPa on the right in interlayer 1, and reduced to 1.9 MPa on the right (Fig. 7 iii-a). Meanwhile, the permeability in interlayer 1 is unevenly varying due to the state of the stress; on the right side, it is firstly decreased to 1:285  1020 m2, then increased to 1:315  1020 m2, and on the left, it is increased to 1:315  1020 m2 then decreased to 1:28  1020 m2. However, the pore pressure in interlayer 2 is increasing significantly to 15 MPa because the much higher permeability (Fig. 7 ieb), causing the increase in maximum principal stress until it reaches the strength (Fig. 7 ii-b), while the shear stress is increasing until the stresses have been released after fracture has been produced (Fig. 7 iii-b). At the same time, the permeability in interlayer 2 has been increased to 1:0  1013 m2 within the fractured area from the initial 1:35  1016 m2, and slightly decreased to 1:33  1016 m2 on the left side (Fig. 7 iv-b). 3.3. Hydraulic fracturing when q ¼ 60 When the intersection angle increases to 60 , in the minimum

968

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 5. Top view of maximum principal stress (MPa) distribution during the hydraulic fracturing stages for q ¼ 30 : (a) 0 s, initial maximum principal stress; (b) 42 s, the hydraulic fracture approaching the interlayer; (c) 192 s, the hydraulic fracture crossing the interlayer; (d) 492 s, the hydraulic fracture crossed the interlayer.

horizontal stress direction, the distance from the well to interlayer 1 becomes 40 mm, to interlayer 2 it is 60 mm. The initial maximum principal stress in interlayer 1 is 0.6 MPa, in interlayer 2 it is 1.3 MPa, and reservoir 0.9 MPa (Fig. 8a). The induced stress along the fracture is intersecting with the interlayers at an early state (Fig. 8b). Meanwhile, the tensile stresses from the fracture tips is increasing the stresses within the interlayer (Fig. 8b). The hydraulic fracture reaches interlayer 1 at 347 s, interlayer 2434 s, and passes interlayer 2 at 900 s. After the hydraulic fracture has reaches the interlayers, the contact areas between the hydraulic fracture and the interlayers are wide, and the stress has been increased to 0.5 MPa at interlayer 1 which has arrested the hydraulic fracture

(Fig. 8 ced). The tensile stress induced from the fracture tip intersect with the interlayers, and the affected region is overlapped by the compressive stress created along the (Fig. 8 bed). When the intersection angle is 60 , the initial maximum shear stress in interlayer 1 is 2.1 MPa, in interlayer 2 it is 2.9 MPa, in the reservoir 2.3 MPa (Fig. 9a). The maximum shear stress in increased to 3.5 MPa at the fracture tip, and decreased 1.8 MPa along the fracture (Fig. 9b). The decreased stress area is intersecting with the interlayers and decreased the maximum shear stress in interlayer 1e1.8 MPa, in interlayer 2 it is decreased to 2.5 MPa. In the meantime, the maximum stress in the interlayers also have been increased by the stress concentrations at the fracture tip; in

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

969

Fig. 6. Top view of maximum shear stress (MPa) distribution during the different hydraulic fracturing stages for q ¼ 30 : (a) 0 s, initial maximum shear stress; (b) 42 s, the hydraulic fracture approaching the interlayer; (c) 192 s, the hydraulic fracture crossing the interlayer; (d) 492 s, the hydraulic fracture crossed the interlayer.

interlayer 1, the maximum shear stress has been increased to 3.0 MPa at the fracture tip, and in interlayer 2, it is 4 MPa (Fig. 9 d). Further investigations show that the maximum pore pressure in interlayer 1 is only 6:0  1010 MPa (Fig. 10 iea); in interlayer 2 it is increased to 16 MPa (Fig. 10 ieb). The principal stress is not evenly distributed in the interlayers: in interlayer 1, on the right side, the peak stress has been gradually increased to 0.5 MPa overtime, and on the left, it is decreased to 1.3 MPa (Fig. 10 ii-a); in interlayer 2, the maximum principal stress keeps increased until the interlayer has been fractured with the stress been released (Fig. 10 ii-b). However, on the right of the interlayer, the maximum principal stress is decreased to 1.5 MPa (Fig. 10 ii-b). The maximum shear stresses in interlayers have a similar trend

to the maximum principal stresses, which is decreased to 1.9 MPa on the right and increased to 3.0 MPa on the left in interlayer 1 (Fig. 10 iii-a), and increased to 4 MPa on the left, then decreased to 2.8 MPa on the right in interlayer 2 (Fig. 10 iii-b). The initial permeability in interlayer 1 is 1:39  1020 m2, in interlayer 2 it is 1:29  1016 m2. Due to the stress variation within the interlayers, the permeability in interlayer 1 decreases to 1:35  1020 m2 where the stress has been decreased, and increased to 1:42  1020 m2 where the stress has been increased (Fig. 10 iv-a); in interlayer 2 the permeability has been drastically increased to 1:0  1013 m2 within the fracture area, and decreases to 1:26  1016 m2 where the stresses are decreased (Fig. 10 iv-b).

970

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 7. Variations within the interlayers when q ¼ 30 during hydraulic fracturing process: (i) pore pressure; (ii) maximum principal stress; (iii) maximum shear stress; (iv) permeability.

3.4. Hydraulic fracturing when q ¼ 90 When q ¼ 90 , the initial maximum principal stresses in the reservoir and interlayers are homogeneously distributed, which is 0.9 MPa (Fig. 11a). After hydraulic fracture has been formed, the stress is concentrated at the fracture tips at 1 MPa, and the stress along the fracture has been decreased to 1.55 MPa (Fig. 11 bec). The initial maximum shear stress in the reservoir is 2.3 MPa, in interlayer 1 it is 1.4 MPa, and interlayer 2 it is 3.45 MPa (Fig. 12a). After hydraulic fracture has been created, the stress is concentrated at the fracture tips at 3 MPa, and reduced along the fracture to 2.1 MPa (Fig. 12 bec). The area of the stress increase from the fracture tips has expanded through the reservoir, while the decrease from the fracture body has reached beyond the interlayers (Fig. 12 bec). The pore pressure in interlayer 2 is slight changed (Fig. 13 ieb) but in interlayer 1 it stays at zero (Fig. 13 iea). However, the maximum principal stress in interlayer 1 has been decreased to 1.5 MPa at the central (Fig. 13 ii-a). It is resulted from the induced compressive stress along the hydraulic fracture. At the other hand, the maximum principal stress is increased to 0.7 MPa beside the decreased area (Fig. 13 ii-a). It also can be seen that the peaks of the stress are getting closer over time, which means the

influence from the induced stress tends to be stable over time (Fig. 13 ii-a, iii-a). Accordingly, the permeability within the interlayer has been decreased to 1:41  1020 m2 where the stresses have been decreased, and it has been increased to 1:445  1020 m2 where the stresses have been increased. A similar trend can be found in interlayer 2 (Fig. 13 ieb ~ iv-b). 3.5. Hydraulic fracturing of different intersection angles of interlayers under a same time scale By comparing the states at the same time, it is clearer to know the influence of the intersection angle. Fig. 14 illustrates the hydraulic fracturing at 138 s for different intersection angles of the interlayers. At 138 s, in interlayer 1, the pore pressure has reached the highest 2.0E-13 m2 when the angle is 0 , and the peak value of 30 is 1.0E-13 m2, and for the other two angles, the permeability is close to zero (Fig. 14 iea). At the same time, the initial maximum principal stress is 2.34 MPa when the angle is 0 , 1.3 MPa of 30 , 0.6 MPa of 60 , and 0.9 MPa of 90 . This means that the initial principal stress in interlayer 1 is increasing when the intersection angle is increasing; however, after it proceeds 60 , the value of the stress has dropped because its direction has been changed (Fig. 14 ii-a). It is found that the peak values of 0 and 90 are

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

971

Fig. 8. Top view of maximum principal stress (MPa) distribution during the hydraulic fracturing stages for q ¼ 60 : (a) 0 s, initial maximum principal stress; (b) 90 s, the hydraulic fracture approaching the interlayer; (c) 472 s, the hydraulic fracture crossing the interlayer; (d) 1622 s, the hydraulic fracture crossed the interlayer.

pointing to each other. As mentioned above, when the intersection angle is higher, the induced compressive stress from the perpendicular direction of the fracture has a more significant influence. When the intersection angle is 0 , interlayer 1 is perpendicular to the hydraulic fracture, thus the compressive stress from the fracture body has a weakest influence; as when the intersection angle is increasing, the impact from the compressive stress is increasing; when the angle gets to 90 , interlayer 1 is parallel to the hydraulic fracture, the induced compressive from the hydraulic has a

strongest influence while the tensile stress from the fracture tip has a weakest influence. This trend also applies to the shear stress (Fig. 14 iii-a). Due to the different stress condition, the permeability is different. When the intersection angle changes from 0 to 30 , 60 and 90 , the initial permeability is 1:27  1020 m2, 1:31  1020 m2, 1:39  1020 m2 and 1:44  1020 m2 (Fig. 14 iv-a). There is a protrude in each permeability curve, and when the intersection angle is 0 , the protrude is located at the central; when the intersection angle is increasing, the protrude is shifting to the

972

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 9. Top view of maximum shear stress (MPa) distribution during the hydraulic fracturing stages for q ¼ 60 : (a) 0 s, initial maximum principal stress; (b) 90 s, the hydraulic fracture approaching the interlayer; (c) 472 s, the hydraulic fracture crossing the interlayer; (d) 1622 s, the hydraulic fracture crossed the interlayer.

left, and close to the protrude, the permeability on the right is 1:0  1022 m2 larger than the right when the intersection angle is 30 , and 2  1022 m2 larger than the right when the intersection angle is 60 ; when the angle reaches 90 , it is located in the middle again (Fig. 14 iv-a). As for interlayer 2, the pore pressure increases first when the intersection angle is 0 because the distance from the well to the interlayer is the shortest. Because of its high initial permeability, the pore pressure reaches 5 MPa (Fig. 14 ieb). At 138 s, the peak

maximum principal stress has been increased to 1.2 MPa of angle 0 , 0.5 MPa of angle 30 . As for the angle is 60 , the stress partly increased to 1.0 MPa and partly decreased to 1.5 MPa of angle 60 (Fig. 14 ii-b). Considering the pore pressure is almost zero, thus the increase of the stress is from the induced tensile stress from the fracture tip, and the decrease is from the compressive stress created along the fracture. Because of the high intersection angle, the compressive can cause a significant influence on the interlayer. The initial maximum shear stress is 2.5 MPa for 0, 2.3 MPa for 30,

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

973

Fig. 10. Variations within the interlayers when q ¼ 60 during hydraulic fracturing process: (i) pore pressure; (ii) maximum principal stress; (iii) maximum shear stress; (iv) permeability.

2.9 MPa for 60 and 3.45 MPa for 90 (Fig. 14 iii-b). The maximum shear stress has been decreased in different extend for different angles. Due to the different stress states, the initial permeability for the angles from 0 to 90 is 1:38  1016 m2, 1:35  1016 m2, 1:29  1016 m2 and 1:26  1016 m2 respectively. Moreover, the permeability has been significantly increased to 1:67  1016 m2 when the intersection angle is 0 . At 588 s, hydraulic fracture has passed through interlayer 2 for the intersection angles of 0 and 30 , and the stresses are released within the fractured area (Fig. 15 ii-b, iii-b). In the meantime the pore pressure has reached 20 MPa (Fig. 15 i-b), also the permeability has been increased significantly to 1:67  1016 m2 within the fractured area (Fig. 15 iv). With the above numerical simulation, it is noticed that when the intersection angle of the interlayer is different, the following relations have been changed as well: (1) the initial stress within the interlayers; (2) the intersection angle between the hydraulic fracture and the interlayers, which further affects: (3) the propagation distances of hydraulic fracture from the well to the interlayers; (4) the distance between the fracture and interlayer perpendicular to

the propagation direction; (5) and the contact area between the hydraulic fracture and interlayers. Take interlayer 1 for example, when the intersection angle changes from 0 to 30 , 60 and 90 , the propagation distance from the well to the interlayer has been increased from 20 mm to 23 mm, 40 mm and infinite; at the meantime, at the direction perpendicular to the fracture propagation, the distance between the well and the interlayer is decreasing, which is infinite at 0 , 40 mm at 30 , 23 mm at 60 and 20 mm at 90 . These two types of distances are important because the first is related to the tensile stress induced from the hydraulic fracture, and the latter is related to the compressive stress induced along the fracture. 4. Conclusions In this paper, FEM has been used to study hydraulic fracturing initiation and propagation with interlayers of different intersection angles with the interlayer. A coupled stress and fluid constitutive has been applied with a continuum damage mechanism, and the permeability is coupled with stress state. The results lead to the

974

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 11. Top view of maximum principal stress (MPa) distribution during the different hydraulic fracturing stages for q ¼ 90 : (a) initial maximum principal stress; (b) 41 s; (c) 490 s.

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

975

Fig. 12. Top view of maximum shear stress (MPa) distribution during the different hydraulic fracturing stages for q ¼ 90 : (a) initial maximum shear stress; (b) 41 s; (c) 490 s.

976

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 13. Variations within the interlayers when q ¼ 90 during hydraulic fracturing process: (i) pore pressure; (ii) maximum principal stress; (iii) maximum shear stress; (iv) permeability.

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

977

Fig. 14. Variations within the interlayers of different intersection angle at 138 s: (i) pore pressure; (ii) maximum principal stress; (iii) maximum shear stress; (iv) permeability.

978

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 15. Variations within the interlayers of different intersection angle at 588 s: (i) pore pressure; (ii) maximum principal stress; (iii) maximum shear stress; (iv) permeability.

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

following conclusions: 1. With different intersection angles, the initial stress within the interlayer is different. However, the maximum principal stress and shear stress are not proportionally increasing with the intersection angle because the directions of the principal stresses may be changed. With higher intersection angle, the initial permeability is higher due to the stress state, which means the fluid leakage through the interlayer may be higher; 2. The stress of the interlayer has been greatly influenced by the induced stress from the hydraulic fracture before the influence of the pore pressure. When the intersection angle is increasing, the induced compressive stress from the fracture has a wider contact area with the interlayer, also the distance for the compressive stress area to reach the interlayer has been reduced, thus it has a greater impact when the intersection angle is higher; 3. The induced compressive stress from the fracture reduces the permeability of the interlayer, and the induced tensile stress increases the permeability before the hydraulic fracture reach the interlayer, which may further influence the propagation of hydraulic within the interlayers. When the intersection angle is increasing, the induced compressive stress has an increasing impact on the stress within the interlayer, because the intersection angle between the hydraulic fracture and the interlayer is decreasing. Accordingly, the permeability is decreased by the induced compressive stress, while it is increased in the induced tensile stress area; 4. The propagation distance from the well to the interlayer is larger with a higher intersection angle. Thus, the hydraulic fracture needs more time to reach and/or propagate through the interlayers.

Acknowledgement



kAðpb  pa Þ mL

(1A)

where Q is the flow quantity, m3 =s; k, the permeability, m2 ; m is the fluid viscosity, Pa$s. L is the length of the fracture, m; A, the cross sectional area, equals to l$a, in which l is the width of the fracture, and a the height, m; pa and pb are the pore pressure at the two ends of the fracture.

Fig. 1A. Conceptual model of a fracture.

In order to verify the fluid flow simulation, a numerical model of a fracture has been built (Fig. 2A). The length of the fracture is 10 m, the height is 1 m, and the width is 0.1 m Pa is 4:05  107 Pa, and Pb is 4:05  105 Pa. The permeability is 8:88  1011 m2 , and viscosity 8:9  104 Pa$s. The flow quantity is 0.04 m3 =s, thus, the fluid velocity ¼ flow quantity/flow cross section area, which is 0.04/ (1  0.1) ¼ 0.4 m/s. The numerical simulation results have shown that the pore pressure.

Fig. 2A. FEM model of a fracture: pa ¼ 4:5  105 Pa; pb ¼ 4:5  107 Pa.

The authors thank Dr. Stephanie Hamilton for providing help in the language. The authors also thank the anonymous reviewers for improving the paper.

Appendix A. Benchmark of darcy flow PANDAS (Parallel Adaptive Nonlinear Deformation Analysis Software) is a fully 3D finite element code, and capable of solid-fluid coupling simulation. It was developed and verified for simulating coupled fracturing processes in a complex geological setting in 3D, subjected to the affecting hydraulic and thermal factors for simulating coupled crustal dynamics. Although extensive benchmarks of the simulation capabilities of PANDAS have been done (Xing and Makinouchi, 2002; Xing et al. 2006, 2007; Zhang, 2012; Xing, 2014), here are some basic benchmarks and verification related with this research on hydraulic fracturing simulation. Since in this research, Darcy's flow is considered for fluid flow in porous media (i.e. rocks), so the first benchmark is to verify the Darcy's flow simulation. Assume the fluid flow through the fracture (Fig. 1A) follows the Darcy's law, which can be expressed as:

979

Fig. 3A. Pore pressure vs distance from Pa.

980

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 4A. Simulated results of pore pressure (Pa).

Fig. 5A. Simulated results of fluid velocity (m/s).

the stress components at the infinite distance along the axial x, y and z, respectively; Rw the well radius; r the distance from the well axis.

Appendix B. Benchmark of stress near a well The theoretical solutions of the stress around a well can be expressed as (Fjaer et al., 1992):

sr ¼

s0x þ s0y 2

2  t0xy

"

sz ¼

! þ

s0x  s0y !

2

14

sin 2 q þ pw

! R2w R4w þ 3 cos 2 q r2 r4

R2w r2

! ! s0x  s0y R2w R4w  cos 2 q 1 þ 3 2 r2 r4 ! R4 R2 1þ3 w sin 2 q  pw w 4 r r2

s0x þ s0y

s0z

R2w r2

R2 R4 þ3 w 14 w 2 r r4

þ t0xy

sq ¼

1

(2A)



 R2  R4 w  n 2 s0x  s0y cos 2 q þ 4t0xy w sin 2 q 2 r r4

(3A)

# (4A)

where ðx; y; zÞ denote the Cartesian coordinate system, and ðr; q; zÞ denote the cylinder polar coordinate system (Fig. 6A). s0x ; s0y ; s0z are

Fig. 6A. Coordinate systems near a wellbore.

In order to measure the stress near a wellbore, a 3D model of 10 m  10 m  10 m has been built. A well of radius 0.1 m is located at the centre of the model, penetrating through the model (Fig. 7A, Fig. 8A). The far field stresses are 5.8 MPa, 2.3 MPa and 2.6 MPa for s0x ; s0y ; s0z respectively; and three mutually orthogonally intersected faces are fixed to the normal axial movement. In order to inspect the results, several nodes are picked on the path with q ¼ 0 (Fig. 7A).

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 7A. Overview of the model for testing the near wellbore stress state, with the boundary conditions and FEM mesh.

Fig. 8A. Top view of the mesh near the wellbore opening.

The simulated stress results are showing a good match with the theoretical results (Fig. 9Ae12A).

Fig. 9A. Comparison of theoretical and simulated stresses near the wellbore.

981

982

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

Fig. 10A. Simulated stress distribution of sx near the wellbore.

Fig. 11A. Simulated stress distribution of sy near the wellbore.

Fig. 12A. Simulated stress distribution of sz near the wellbore.

Appendix C. Benchmark of hydraulic fracturing propagation In this benchmark the propagation of the hydraulic fracture will be verified by an indoor large scale hydraulic fracturing experiment (Fu et al., 2013). The dimension of the rock sample is 762 mm  762 mm  762 mm, with a wellbore of 27 mm in diameter (Fig. 13A). The propagation of hydraulic fracture can be

tested by the sonar receiver placed on the faces of the rock sample (Fig. 14A). The parameters of the rock sample is shown in Table 1A. It was tested that the fracturing hydraulic pressure was 24 MPa, the initial stresses were 21 MPa, 14 MPa, and 7 MPa respectively (Fig. 16A). The sonar receiver detected sonar signals during the hydraulic fracturing process, and the profile of the hydraulic fracture was decided by the dense sonar signal area (Fig. 15A).

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

983

Fig. 13A. Overview of the rock sample used in the hydraulic fracturing test.

The numerical model is based on the indoor experiment, the boundary conditions follows the experiment (Fig. 16A). In order to save the simulation cost, 1/4 model has been built (Fig. 17A). In the simulation, when an element has been fractured, the stress is released to zero, thus, the profile of the hydraulic fracture can be depicted by the stress distribution. The simulation result of the propagation of the hydraulic fracture showed a good match with the experiment (Fig. 18A).

Fig. 14A. Sonar receiver locations.

Fig. 16A. Boundary condition of the indoor experiment and numerical simulation. Table 1A Parameters of the rock sample. Young's Modulus (MPa)

Poisson's ratio

Porosity (%)

Permeability 1015(m2)

Tensile strength (MPa)

Compressive Strength (MPa)

15500

0.19

14.474

1.036

8

200

Fig. 15A. Sonar detected profile of the hydraulic fracture.

Fig. 17A. Numerical model of the hydraulic fracturing: (a) 1/4 model; (b) finite element mesh.

Fig. 18A. Profile of hydraulic fracture: (a) 60s; (b) 100s; (c) 600 s; (d) 800 s.

Q. Li, H. Xing / Journal of Natural Gas Science and Engineering 36 (2016) 963e985

References Chudnovsky, A., Fan, J., Shulkin, Y., Dudley, J.W., Nichols, W.B., Wong, G.K., 2001. Hydraulic Fracture Containment in Layered Media, Experiment and Computer Simulation. American Rock Mechanics Association. Cleary, M.P., 1980. Analysis of Mechanisms and Procedures for Producing Favourable Shapes of Hydraulic Fractures. Society of Petroleum Engineers. Daneshy, A.A., 1978. Hydraulic Fracture Propagation in Layered Formations. Daneshy, A.A., 2009. Factors controlling the vertical growth of hydraulic fractures. In: SPE Hydraulic Fracturing Technology Conference. Society of Petroleum Engineers, The Woodlands, Texas. Dong, J.J., Hsu, J.Y., Wu, W.J., Shimamoto, T., Hung, J.H., Yeh, E.C., Wu, Y.H., Sone, H., 2010. Stress-dependence of the permeability and porosity of sandstone and shale from TCDP Hole-A. Int. J. Rock Mech. Min. Sci. 47 (7), 1141e1157. Fjaer, E., Horsrud, P., Raaen, A., Risnes, R., Holt, R., 1992. Petroleum related rock mechanics. Elsevier. Fu, H., Cui, M., Zou, J., Peng, Y., Liu, Y, 2013. Experiment of hydraulic fracturing based on sonar detect technics. J. Northeast Pet. Univ. 37 (2), 96e101. Handwerger, D.A., Keller, J., Vaughn, K., 2011. Improved Petrophysical Core Measurements on Tight Shale Reservoirs Using Retort and Crushed Samples. Society of Petroleum Engineers. Li, Q., Xing, H., 2015. Numerical analysis of the material parameter effects on the initiation of hydraulic fracture in a near wellbore region. J. Nat. Gas Sci. Eng. 27 (3), 1597e1608. Li, Q., Xing, H., Liu, J., Liu, X., 2015. A review on hydraulic fracturing of unconventional reservoir. Petroleum 1 (1), 8e15. Nuller, B., Ryvkin, M., Chudnovsky, A., Dudley, J.W., Wong, G.K., 2001. Crack Interaction with an Interface in Laminated Elastic Media. American Rock Mechanics Association. €cke, D.R., Greenwell, C., Gluyas, J., Cornford, C., 2015. The Effect of Raji, M., Gro Interbedding on Shale Reservoirs Properties. Marine and Petroleum Geology.

985

Smith, M.B., Bale, A.B., Britt, L.K., Klein, H.H., Siebrits, E., Dang, X., 2001. Layered modulus effects on fracture propagation, proppant placement, and fracture modelling. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers Inc, New Orleans, Louisiana. Teufel L. W. and Clark J. A. (1981). Hydraulic Fracture Propagation in Layered Rock: Experimental Studies of Fracture Containment. SPE 1981 Low-permeability Gas Symposium. (Denver). Wang, S., Li, Y., Liu, H., Jiang, M., 2012. Hydraulic Fracture Propagation in Unconventional Reservoirs: The Role of Bedding Plane, 8 (3), 173e191. SL: Structural Longevity. Warpinski, N.R., Branagan, P.T., Peterson, R.E., Wolhart, S.L., 1998. An interpretation of m-site hydraulic fracture diagnostic results. In: SPE Rocky Mountain Regional/Low-permeability Reservoirs Symposium. Society of Petroleum, Engineers, Inc, Denver, Colorado, 1998 Copyright 1998. Wu, H., Chudnovsky, A., Dudley, J.W., Wong, G.K., 2004. A Map of Fracture Behavior in the Vicinity of an Interface. American Rock Mechanics Association. Xing, H., 2014. Finite element simulation of transient geothermal flow in extremely heterogeneous fractured porous media. J. Geochem. Explor. 144 (Part A), 168e178. Xing, H., Makinouchi, A., 2002. Three dimensional finite element modelling of thermomechanical frictional contact between finite deformation bodies using R-minimum strategy. Comput. Methods Appl. Mech. Eng. 191 (37), 4193e4214. Xing, H., Mora, P., Makinouchi, A., 2006. A unified friction description and its application to the simulation of frictional instability using the finite element method. Philos. Mag. 86 (21e22), 3453e3475. Xing, H., Zhang, J., Yin, C., 2007. A finite element analysis of tidal deformation of the entire earth with a discontinuous outer layer. Geophysical Journal International 170 (3), 961e970. Zhang, J., 2012. Forward and Inverse Modelling of Enhanced Geothermal Processes. PhD Thesis. University of Queensland.