3D pulsed fracture initiation and propagation around wellbore: A numerical study

Computers and Geotechnics 119 (2020) 103374

Contents lists available at ScienceDirect

Computers and Geotechnics journal homepage: www.elsevier.com/locate/compgeo

Research Paper

3D pulsed fracture initiation and propagation around wellbore: A numerical study

T

Yujie Wanga, Bing Zhaob, Xiaobo Heb, Zhennan Zhanga,

⁎

a b

School of Naval Architecture, Ocean and Civil Engineering, Shanghai Jiao Tong University, Shanghai 200240, China Research Institute of Petroleum Engineering, Sinopec Northwest Oilfield Branch, Urumqi, Xinjiang 830011, China

ARTICLE INFO

ABSTRACT

Keywords: Pulsed fracture Numerical simulation Augmented virtual internal bond Three-dimension Hydraulic fracture Reservoir stimulation

Subjected to the extremely high pressurization rate, the pulsed fracture (PF) can branch to form a complex fracture network around the wellbore, which is the most desired in the reservoir stimulation. To study the PF behaviors further, the 3D numerical simulation is conducted. The rock formation is modeled by the augmented virtual internal bond method, which needs no separate fracture criterion in fracture simulation since the micro fracture mechanism has been embedded into the macro constitutive relation via micro bond potential. The simulation results demonstrate that the 3D PF always presents an ear-like profile. Some branched fractures advance towards the minor in-situ stress direction. The higher the pressurization rate is, the more complex the generated fracture network is. The longer isolation section makes the fracturing process easier. The PF speed is governed by the rock modulus. Besides these qualitative observations, some quantitative conclusions have been drawn. The findings of this work can deepen understanding on the PF and provide valuable references for PF operation in a reservoir.

1. Introduction Hydraulic fracture (HF) technique has been extensively applied to the reservoir stimulation. The HF prefers to grow along the maximum in-situ stress direction. Thus, it is very hard for HF to connect the reservoir in the minor in-situ stress direction. And the generated fracture network is usually not complex as desired. As a new potential stimulation method, the pulsed fracture (PF) technique has received much attention in recent years. Under the extremely high pressurization rate [1–3], multiple fractures may radially initiate on the wellbore wall simultaneously. These fractures grow fast and branch unstably. They coalesce with the natural ones to form a complex fracture network, which greatly improves the permeability of reservoir. In addition, the PF technique is water-saving and environment-protective [4]. So, it should be an inspiring stimulation technique for reservoir. For the conventional hydraulic fracturing process is driven by the quasi-static hydraulic pressure, the rock inherent heterogeneity, grain size and natural fractures have influence on the fracturing behavior as reported

⁎

in Huang et al. [5] and Zhang et al. [6]. While the PF is driven by pulse pressure, the PF behaviors are much more complex than the conventional HF [7–11]. To investigate the regularities of PF, some scholars have conducted the experiments. It is found that the pressurization rate has significant impact on the multiple fracture initiation, e.g., [12–18]. Yang and Rasmus [13] and Zhao et al. [19] obtained the critical pressurization rate for multiple fractures generation through the experiment. Zeigler et al. [20] have ever used three techniques, namely the single burn, the consecutive burn and the pulse tailoring approach, to generate pulse pressure. Compared with the single burn approach, both the consecutive burn and pulse tailoring approaches significantly improve matrix fracturing. The pulse tailoring approach is substantially more energy efficient than the consecutive burn approach. The change of pressurization rate can induce the cavitation effect and local back flow in fracture [21]. These experimental studies reveal many regularities of the PF. For the PF process usually occurs in a very short duration, it is hard to precisely measure some quantities. The size of the test sample is usually

Corresponding author. E-mail address: [email protected] (Z. Zhang).

https://doi.org/10.1016/j.compgeo.2019.103374 Received 18 July 2019; Received in revised form 7 October 2019; Accepted 28 November 2019 0266-352X/ © 2019 Elsevier Ltd. All rights reserved.

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

Fig. 1. Modeling method of augmented virtual internal bond model: (a) the representative element volume; (b) bond deformation.

m

within dozens of inches. Thus, it is inadequate to completely get insight into the PF only through the usual lab test. So, it is quite necessary to conduct the numerical simulation study. Some scholars have carried out the 2D simulation on PF and some meaningful findings have been obtained, e.g. [22–27]. In spite of these advances in 2D study, many features of PF cannot be revealed due to the limitation of 2D method itself. Therefore, the 3D simulation is needed to reproduce the full behaviors of PF. Through 3D simulation, Safari et al. [9] have studied the impact of the ductile-to-brittle transition on fracturing. To investigate the PF behaviors around the wellbore, the 3D numerical simulations will be conducted in this paper. The augmented virtual internal bond (AVIB) [28] will be employed to simulate HF. In AVIB, the micro fracture criteria (initiation, growth, branching and coalescence criterion) have been implicitly embedded into the macro constitutive relation via the micro ‘interatomic’ potential. So, AVIB can simulate the dynamic fracture process efficiently without any separate

P

P

k

P i

P

j

Fig. 2. Pulse pressure exerted on the inner faces of a newly failed element.

Fig. 3. Simulation setup of pulsed fracture: (a) simulation object and boundary condition (isolation section is colored in yellow); (b) the wellbore wall excavated through element removing; (c) linearized pressurization scheme at the pre-peak stage of pulse pressure. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

dp

Wellbore

Sz Sy

Sx

(a)

(b)

P Pk

P0 O

tk

t

(c)

2

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

Fig. 4. Fracturing process with in-situ stresses Sx Sy Sz = 125 120 115 MPa at the time (a) t = 5.48 ms \* MERGEFORMAT; (b) t = 5.56 ms \* MERGEFORMAT; (c) t = 5.64 ms \* MERGEFORMAT; (d) t = 5.72 ms \* MERGEFORMAT.

fracture criterion.

where n \* MERGEFORMAT and t \* MERGEFORMAT are the characteristic lengths for normal and tangential separations, n = t = ~t l 0 \* MERGEFORMAT; ~t \* MERGEFORMAT is the strain of bond at the peak uniaxial tensile force, usually taking ~t = t \* MERGEFORMAT with t \* MERGEFORMAT being the strain value at the peak stress of the uniaxial tension stress-strain curve of material; l 0 \* MERGEFORMAT is the undeformed bond length; n \* MERGEFORMAT and t \* MERGEFORMAT are the normal and tangential deformations of bond; the bond energy n \* MERGEFORMAT and energy ratio q\* MERGEFORMAT are related to the macro elastic constants in AVIB as

2. Brief introduction to AVIB The AVIB [28] stems from the VIB [29]. It is a micro-macro constitutive model. The representative element volume (REV) in view of VIB consists of material particles on the micro scale (Fig. 1a). These particles are connected by virtual internal bond. A bond potential is introduced to describe the interactions between particles. To overcome the fixed Poisson ratio problem induced by the interatomic potential, the Xu-Needleman potential [30,31] is adopted as the bond potential in AVIB, which reads

U( ) =

n

n exp

n n

1+

n n

1

q + q exp

2 t 2 t

n

=

V n2 3E · l02 4 (1 2v )

q= (1)

with

3

2 t (1 4v ) 2· n 2(1 + v)

E \*

MERGEFORMAT

(2) being

the

Young’s

modulus;

v \*

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

Fig. 5. Fracturing process with in-situ stresses Sx Sy Sz = 135 120/115 MPa at the time (a) t = 5.20 ms \* MERGEFORMAT; (b) t = 5.28 ms \* MERGEFORMAT; (c) t = 5.36 ms \* MERGEFORMAT; (d) t = 5.44 ms \* MERGEFORMAT.

MERGEFORMAT the Poisson’s ratio and V \* MERGEFORMAT the volume of REV. According to the Cauchy-Born rule, the bond deformations (shown in Fig. 1b) are associated with the macro strain tensor by n 2 t

=[

=

T

T T

Cijkl =

T

) 2] l02

(3)

where \* MERGEFORMATis the unit orientation vector of bond, = [sin cos , sin sin , cos ]T \* MERGEFORMAT in the spherical coordinates (Fig. 1b). The stress tensor of AVIB is derived as ij

=

1 V

U( ij

n,

t)

2 ij

U( kl

n,

t)

(5)

where ij \* MERGEFORMAT is the strain tensor of REV; the operator \* MERGEFORMAT stands for the integral 2 = 0 ( ) D ( , ) sin d d D ( , ) \* \* MERGEFORMAT with 0 MERGEFORMAT being the bond distribution density in terms of the spherical coordinates \* MERGEFORMAT and \* MERGEFORMAT.

l0 (

1 V

3. Simulation method of PF 3.1. Identification of fracture path In AVIB, the fracturing process is a natural outcome of the constitutive relation computation. With the deformation increasing, the stiffness of an material element decays. As the deformation exceeds a certain value, its stiffness automatically reduces to a negligible value.

(4)

and the tangent modulus tensor is

4

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

Fig. 6. Fracturing process with in-situ stresses Sx Sy Sz = 145 120 115 MPa at the time (a) t = 4.94 ms \* MERGEFORMAT; (b) t = 5.02 ms \* MERGEFORMAT; (c) t = 5.10 ms \* MERGEFORMAT; (d) t = 5.18 ms \* MERGEFORMAT.

dz

dx

dy (a) Fig. 7. Definition of fracture extension: (a) dx; (b) dy and dz.

5

(b)

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

\* MERGEFORMAT is a coefficient greater than unit. Usually, when taking = 5.0 \* MERGEFORMAT, the residual strength of element can be weakened to the negligible level for the given bond potential. 3.2. Application of pulse pressure The internal pressure is the driving force for fracture propagation. The distribution of the internal pressure along fracture is determined by many factors. The pumping rate and the fluid rheology significantly affect the pressure distribution and the fracturing process [32–35]. In the PF process, under the extremely high pressurization rate, the turbulent flow may occur, which makes the fluid flow much more complicated. Nilson et al. [36] found that it takes few milliseconds for the gas pressure completely to enter the fracture and few time later the pressure becomes nearly uniform along the fracture, except a small region just behind the tip. If we precisely calculate the pressure distribution along fracture, the algorithm will become very time-consuming and complicated. So, we could assume that the pressure is uniformly distributed along the fracture as an approximate calculation. Lu and Tao [37] also took such treatment method on fluid pressure distribution along fracture. According to such assumption, once an element is identified as failed in the current time step, the pulse pressure is directly applied to this failed element in the next time step, shown in Fig. 2. The equivalent nodal force induced by internal pressure reads

Fi =

(P S ) n n d

(7)

where Fi \* MERGEFORMAT is the equivalent nodal force vector of node i; P \* MERGEFORMAT is the mean internal pressure in a failed element; S \* MERGEFORMAT is the area of element face; n \* MERGEFORMAT is the normal outwards vector; nd \* MERGEFORMAT is the node number of this element face. Taking the element shown in Fig. 2 as example, the equivalent nodal force induced by internal pressure on the triangular area i-j-m can be calculated as

Fi =

(rij × rim) P 2n d

(8)

where rij, rim \* MERGEFORMAT are the edge vectors of triangular area; nd = 3 \* MERGEFORMAT. 4. Simulation examples The simulation sample for 3D PF is a cube whose dimension is 1.0 m× 1.0 m× 1.0 m \* MERGEFORMAT. A central run-through round hole is excavated in the cube as the wellbore with the radius of 0.08255 m (Fig. 3a). To ensure the perfect randomness of the mesh around a wellbore so as to reduce the mesh-dependence of the simulated fracture path, we firstly mesh the object with the randomized tetrahedral elements and then excavate a wellbore in the cube. In other words, the elements originally in the wellbore are removed. The generated wellbore wall is shown in Fig. 3b. The total element number is 1,113,366 and the total node number is 192,631. The cube is subjected to the in-situ stresses Sx \* MERGEFORMAT, Sy \* MERGEFORMAT and Sz \* MERGEFORMAT. The propellant is controlled to detonate in an isolation section of the wellbore. We take the pre-peak stage of the

Fig. 8. Fracture extension under different in-situ stresses: (a) dx \* MERGEFORMAT, (b) dy \* MERGEFORMAT and (c) dz \* MERGEFORMAT.

The macro fracture path is represented by the failed element array. Thus, to recognize the fracture path, we should firstly identify the failed element through the following criterion 1

>

t

(6)

where 1\* MERGEFORMAT is the first principal value of strain tensor;

6

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

Fig. 9. Fracture profile with different pressurization rates: (a) P = 10.4 MPa/ms, t = 10.94 ms ; (b) P = 20.8 MPa/ms, t = 5.58 ms ; (c) P = 41.6 MPa/ms, t = 2.90 ms .

7

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

4.1. Effects of in-situ stresses To investigate the effect of in-situ stress, we set three cases of different in-situ stresses. In the three cases, we keep Sy and Sz fixed, namely Sy = 120 MPa , Sz = 115 MPa , and change Sx , namely Case-1: Sx = 125 MPa , Case-2: Sx = 135 MPa and Case-3: Sx = 145 MPa . The length of isolation section is dp = 0.35 m \* MERGEFORMAT. The material and other parameters are: Young’s modulus E = 40.0 GPa \* MERGEFORMAT, Poisson ratio = 0.2 \* MERGEFORMAT, rock density = 2400.0 kg m3 \* MERGEFORMAT, t = 0.2 × 10 3 \* MERGEFORMAT which can lead to the macro uniaxial tensile strength t = 4.0 MPa \* MERGEFORMAT according to the relationship E t 2 \* MERGEFORMAT presented in [28]. Take the pressurization t rate P = 20.8 MPa/ms \* MERGEFORMAT. The simulated results are shown in Figs. 4–8. Figs. 4–6 show the fracturing processes of Case-1, Case-2 and Case3, respectively. It is seen that the fractures prefer to initiate and propagate along the maximum in-situ stress direction. They always present an ear-like profile around the wellbore. With the in-situ stress difference decreasing, the branched fractures present a more complex profile, which facilitates the connection of pulse fractures with the reservoir in minor in-situ stress direction. The fracture extensions are defined in Fig. 7. Their variations with time are shown in Fig. 8. It is found that the larger the in-situ stress difference is, the earlier the fracture begins to grow. This is because the large in-situ stress difference can reduce the breakdown pressure. It is consistent with the prediction of the Kirch’s solution [38]. The curves in Fig. 8 are almost parallel to each other, which indicates that the fracture speed dose not vary with the in-situ stresses. 4.2. Effects of pressurization rate To investigate the impact of the pressurization rate on PF, we take different pressurization rates, namely P = 10.4 \* MERGEFORMAT, 20.8\* MERGEFORMAT, and 41.6 MPa/ms \* MERGEFORMAT. The insitu stresses are: Sx = 130.0 MPa , Sy = 120.0 MPa , Sz = 115.0 MPa . Take the isolation section length dp = 0.35 m \* MERGEFORMAT. The parameters of material are: E = 38.0 GPa \* MERGEFORMAT, = 0.2\* = 2400.0 kg m3 \* MERGEFORMAT, MERGEFORMAT, t = 0.21 × 10 3 \* MERGEFORMAT leading to the uniaxial tensile strength t = 4.0 MPa \* MERGEFORMAT. The simulated results are shown in Figs. 9 and 10. Fig. 9 shows the PF at the same fracturing time, namely 0.22 ms, since the fracture initiates. Generally, the PFs initiate and extend in the maximum in-situ stress direction. An ear-like fracture profile is formed around the wellbore finally. Compared with the lower pressurization rate case, the higher one presents more complicated profile. Fig. 10 shows the fracture extensions with time since PF begins to grow. It indicates that the fracture speeds in one direction are almost the same although the pressurization rates are different. The extension speed at different directions may be different, which results from the different fracture branching degrees.

Fig. 10. Fracture extension with different pressurization rates: (a) dx \* MERGEFORMAT, (b) dy \* MERGEFORMAT and (c) dz \* MERGEFORMAT.

pulse as the pressure input, which is linearized as Fig. 3c. The pulse pressure is directly exerted on the rock formation of wellbore wall with no casing and cement. The initial pore pressure of reservoir is P0 = 70.0 MPa \* MERGEFORMAT.

8

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

Fig. 11. Fracture profile with different Young’s moduli at \* MERGEFORMAT: (a) E = 40.0 GPa \* MERGEFORMAT; (b) E = 60.0 GPa \* MERGEFORMAT; (c) E = 80.0 GPa \* MERGEFORMAT.

9

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

t = 4.0 MPa \* MERGEFORMAT, P = 19.2 MPa/ms \* MERGEFORMAT. The simulated results are shown in Figs. 11 and 12. Fig. 11 shows that the PFs generally extend along the maximum insitu stress direction. They eventually present the similar fracture patterns in that PFs grow outwards like an ear. The higher the modulus is, the faster the fracture grows, shown as Fig. 12. It is also suggested that the modulus has little influence on the fractures initiation time in Fig. 12. According to AVIB, the fracture energy is proportional to 2 t E \* MERGEFORMAT. Under the condition that t \* MERGEFORMAT is fixed, the fracture energy is decreased with E \* MERGEFORMAT increasing. Thus, the rock with a higher modulus is easier to fracture as observed in these cases.

4.4. Effects of isolation section length The length of isolation section is an important design parameter of PF treatment. To study its influence on PF, we set different length values, namely dp = 0.20 \* MERGEFORMAT, 0.35\* MERGEFORMAT and 0.50 m \* MERGEFORMAT. The pressurization rate is P = 10.4 MPa/ms \* MERGEFORMAT. The in-situ stresses are Sx = 130.0 MPa , Sy = 120.0 MPa , and Sz = 115.0 MPa . The parameters of material are E = 38.0 GPa \* MERGEFORMAT, = 0.2\* MERGEFORMAT, = 2400.0 kg m3 \* MERGEFORMAT, t = 0.21 × 10 3 \* MERGEFORMAT. The simulated results are shown in Figs. 13 and 14. Fig. 13 suggests that the general fracture patterns are the same in the different isolation section length cases. The isolation section length has a significant impact on the fracture initiation time. The longer the isolation section is, the earlier the fracture initiates, shown in Fig. 14. This is because the longer isolation section can cause a larger deformation at the middle of the isolation section under the uniform internal pressure. The fracture is likely to initiate where the larger deformation occurs. 4.5. Discussion It is seen that in different cases, the complete 3D PF pattern generally looks like an ‘ear’. If we take a certain cross section normal to the wellbore, we can find the fracture in this cross section is also branched, which is consistent with the 2D simulation results in [25,27]. But the 2D simulation cannot present the complete profile of the PF. Compared the fracture speeds in different cases, it is found that among these factors, namely the in-situ stress, the pressurization rate, the rock modulus and the isolation section length, only the rock modulus has influence on the fracture speed. Chandar and Knauss [39] found that the crack velocity do not seem affected by the stress intensity factor through experiment. Kou et al. [40] suggested that the fracture speed is almost constant although there exists slight variation under different dynamic load levels. Glover et al. [41] holds that crack velocity is probably an intrinsic property of material closely related to the elastic modulus. Freund [42] holds that the critical crack velocity with which the crack propagates stably is related to the Rayleigh wave speed while the Rayleigh wave speed is determined by the modulus and Poisson ratio. All these theoretical and experimental researches suggest that the crack speed in PF is related to the modulus, not the loading conditions. This is consistent with the observation in our simulation. From the simulated fracture profiles, it is seen that some of them are

Fig. 12. Fracture extension with different Young’s moduli: (a) dx \* MERGEFORMAT, (b) dy \* MERGEFORMAT and (c) dz \* MERGEFORMAT.

4.3. Effects of rock modulus To investigate the effects of rock modulus on PF, we set different Young’s moduli E = 40.0 \* MERGEFORMAT, 60.0 \* MERGEFORMAT and 80.0 GPa \* MERGEFORMAT. The in-situ stresses are: Sx = 130.0 MPa , Sy = 120.0 MPa and Sz = 115.0 MPa . The isolation section length isdp = 0.35 m \* MERGEFORMAT. The other parameters are = 0.2 \* MERGEFORMAT, = 2400.0 kg m3 \* MERGEFORMAT,

10

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

Fig. 13. Fracture profile with different isolation section lengths: (a) dp = 0.20 m, t = 11.95 ms \* MERGEFORMAT; (b) dp = 0.35 m, t = 10.95 ms \* MERGEFORMAT; (c) dp = 0.50 m, t = 10.40 ms \* MERGEFORMAT.

11

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

that one initiated later in another side. From Fig. 3b it is seen that the generated wellbore wall is not a smooth one. This leads to the heterogeneity of the breakdown pressure on the wellbore wall. Different breakdown pressure causes the non-simultaneous initiation of fracture at the two opposite wellbore wall. We can verify this by an additional simulation. We simulate the case of Fig. 9b again with a smooth wellbore wall (Fig. 15a). It is found that the fracture symmetrically propagates at the two opposite sides, shown in Fig. 15b–d. The asymmetry of the fracture profile in the former simulations does not influence the conclusions drawn in each simulation cases. The asymmetrical profile of fracture should be more consistent with the real cases. Because in the real wellbore wall some fissures are distributed non-uniformly, shown in Fig. 16. These non-uniformly distributed fissures can cause the breakdown pressure heterogeneity. 5. Discuss and conclusions Through the 3D numerical simulation, it is found that the 3D pulsed fracture branches. With respect to the branching phenomenon, the 3D PF is consistent with the 2D one. But the 2D simulation cannot capture the spacial characteristics of the PF. The simulation results suggest that the PF initiation and propagation are subjected to the in-situ stresses. Generally, the PF prefers to extend along the maximum in-situ stress direction. Different from the conventional hydraulic fracture, many PFs obliquely extend towards the minor in-situ stress direction, which facilitates the connection of PFs with the reservoir in that direction. The influences of the pressurization rate, the rock modulus and the isolation section length have been studied. The higher pressurization rate can lead to earlier fracture initiation and more complicated fracture profile. With the rock modulus increasing, the fracture growth is sped up obviously, which is likely to generate a larger fractured zone. The isolation section length has significant impact on the critical pulse pressure for fracture to initiate. The longer the isolation section is, the lower the critical pulse pressure is. This suggests that the PF stimulation becomes easier when increasing the isolation section length. This study reveals the basic PF features. The findings in this study can provide some valuable references for the PF operation. Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 14. Fracture extension with different isolation section lengths: (a) dx \* MERGEFORMAT, (b) dy \* MERGEFORMAT and (c) dz \* MERGEFORMAT.

Acknowledgement The present work is supported by National Science and Technology Major Projects (No. 2016ZX05014-005-003) and the National Natural Science Foundation of China (No. 11772190), which are gratefully acknowledged. Thanks to the anonymous reviewers for their constructive suggestions to improve the manuscript.

not symmetrical along the major horizontal in-situ stress direction, e.g., Fig. 9. Such asymmetry mainly results from the non-simultaneous initiation of the fracture at the two opposite wellbore wall along this direction. The fracture initiated earlier at one side grows longer than

12

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

Fig. 15. Pulsed fracture with a smooth wellbore: (a) the smooth wellbore wall; (b) fracture profile at the time t = 5.38 ms \* MERGEFORMAT; (c) t = 5.47 ms \* MERGEFORMAT; (d) t = 5.56 ms \* MERGEFORMAT.

References [1] Al-Nakhli AR. Chemically-induced pressure pulse: a new fracturing technology for unconventional reservoirs. In: Proceedings of the SPE Middle East Oil & Gas Show and Conference. Manama, Bahrain. SPE-172551-MS; 8-11 March 2015. [2] Salazar A, Almanza E, Chira P, Folse K. Application of propellant high-energy gas fracturing in gas-injector wells at El Furrial field in northern Monagas state Venezuela. In: Proceedings of the SPE International Symposium and Exhibition on Formation Damage Control. Lafayette, Louisiana. SPE 73756; 20-21 February 2002. [3] Jaimes MG, Castillo RD, Mendoza SA. High energy gas fracturing: a technique of hydraulic prefracturing to reduce the pressure losses by friction in the near wellbore - a Colombian field application. In: Proceedings of the SPE Latin American and Caribbean Petroleum Engineering Conference. Mexico City, Mexico. SPE 152886; 16-18 April 2012. [4] Gandossi L. An overview of hydraulic fracturing and other formation stimulation technologies for shale gas production. Joint Research Centre of the European Commission: Luxembourg; 2013. [5] Huang L, Liu J, Zhang F, Dontsov E, Damjanac B. Exploring the influence of rock inherent heterogeneity and grain size on hydraulic fracturing using discrete element modeling. Int J Solids Struct 2019;176–177:207–20. [6] Zhang F, Damjanac B, Maxwell S. Investigating hydraulic fracturing complexity in naturally fractured rock masses using fully coupled multiscale numerical modeling. Rock Mech Rock Eng 2019. https://doi.org/10.1007/s00603- 019-01851-3. [7] Page JC, Miskimins JL. A comparison of hydraulic and propellant fracture propagation in a shale gas reservoir. J Can Petrol Technol 2009;48(5):26–30. [8] Buller D, Hughes S, Market J, Petre E, Spain D, Odumosu T. Petrophysical evaluation for enhancing hydraulic stimulation in horizontal shale gas wells. In: Proceedings of the SPE Annual Technical Conference and Exhibition. Florence, Italy. SPE 132990; 19-22 September 2010. [9] Safari R, Huang J, Mutlu U, Glanville J. 3D Analysis and Engineering Design of Pulsed Fracturing in Shale Gas Reservoirs. American Rock Mechanics Association; 2014. [10] Safari R, Gandikota R, Mutlu O, Ji M, Glanville J, Abass H. Pulse fracturing in shale

Fig. 16. Fissure distribution on the wellbore wall based on the acoustic image (left) and electrical image (right) of real wellbore wall reported in [43].

13

Computers and Geotechnics 119 (2020) 103374

Y. Wang, et al.

[11] [12] [13]

[14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26]

reservoirs: geomechanical aspects, ductile/brittle transition, and field implications. Soc Petrol Eng 2015. SPE 168759. Wu C, Zhang X, Wang M, Zhou L, Jiang W. Physical simulation study on the hydraulic fracture propagation of coalbed methane well. J Appl Geophys 2018;150:244–53. Detournay EM, Cheng A. Influence of pressurization rate on the magnitude of the breakdown pressure. In: Rock Mechanics Proceedings of the 33rd U.S. Symposium. Publ by A.A. Balkema; 1992. p. 325–33. Yang DW, Risnes R. Experimental study on fracture initiation by pressure pulses. In: Proceedings of the SPE Annual Technical Conference and Exhibition Presented at the SPE Annual Technical Conference and Exhibition. Dallas, Texas. SPE63035; 1-4 October 2000. p. 259–72. Cuderman JF. Multiple fracturing experiments: propellant and borehole considerations. In: SPE Unconventional Gas Recovery Symposium. Pittsburg, PA, USA; May 1982. Cuderman JF, Northrop DA. A propellant-based technology for multiple-fracturing wellbores to enhance gas recovery: application and results in Devonian shale. SPE Prod Eng 1986;1(2):97–103. Schatz JF, Hanson JM. Multiple radial fracturing from a wellbore conceptual model and experiments. American Rock Mechanics Association; 1 January 1986. Schatz JF, Zeigler BJ, Hanson J, Christianson M. Multiple radial fracturing from a wellbore - experimental and theoretical results. American Rock Mechanics Association; 1 January 1987. Malhotra S, Rijken P, Sanchez A. Experimental investigation of propellant fracturing in a large sandstone block. SPE Drill Complet 2018;33(2):87–99. Zhao H, Chen H, Liu G, Li Y, Shi J, Ren P. New insight into mechanisms of fracture network generation in shale gas reservoir. J Petrol Sci Eng 2013;110:193–8. Zeigler B, Schatz J, Nudd E. A comparison of multiple radial fracturing enhancement techniques - laboratory investigation. American Rock Mechanics Association; 1 January 1988. Garagash D, Detournay E. An analysis of the influence of the pressurization rate on the borehole breakdown pressure. Int J Solids Struct 1997;34(24):3099–118. Zhu Z, Mohanty B, Xie H. Numerical investigation of blasting-induced crack initiation and propagation in rocks. Int J Rock Mech Min Sci 2007;44(3):412–24. Zhu Z, Xie H, Mohanty B. Numerical investigation of blasting-induced damage in cylindrical rocks. Int J Rock Mech Min Sci 2008;45(2):111–21. Banadaki MMD, Mohanty B. Numerical simulation of stress wave induced fractures in rock. Int J Impact Eng 2012;s40–41:16–25. Zhang Z, Peng S, Ghassemi A, Ge X. Simulation of complex hydraulic fracture generation in reservoir stimulation. J Petrol Sci Eng 2016;146:272–85. Bendezu M, Romanel C, Roehl D. Finite element analysis of blast-induced fracture propagation in hard rocks. Comput Struct 2017;182:1–13.

[27] Liu Z, Peng S, Zhao H, Geng Y, Zhang Z. Numerical simulation of pulsed fracture in reservoir by using discretized virtual internal bond. J Petrol Sci Eng 2019;181:106197. [28] Zhang Z, Gao H. Simulating fracture propagation in rock and concrete by an augmented virtual internal bond method. Int J Numer Anal Meth Geomech 2012;36:459–82. [29] Gao H, Klein P. Numerical simulation of crack growth in an isotropic solid with randomized internal cohesive bond. J Mech Phys Solids 1998;46(2):187–218. [30] Xu XP, Needleman A. Void nucleation by inclusion debonding in a crystal matrix. Modell Simul Mater Sci Eng 1993;1:111–32. [31] Xu XP, Needleman A. Numerical simulations of fast crack growth in brittle solids. J Mech Phys Solids 1994;42:1397–434. [32] Bohloli B, De Pater CJ. Experimental study on hydraulic fracturing of soft rocks: influence of fluid rheology and confining stress. J Pet Sci Eng 2006;53:1–12. [33] Cheng W, Jiang G. Numerical investigations of the fracture geometry and fluid distribution of multistage consecutive and alternative fracturing in a horizontal well. Comput Geotech 2017;92:41–56. [34] Liu Y, Guo J, Chen Z. Fracturing fluid leakoff behavior in natural fractures: effects of fluid rheology, natural fracture deformability, and hydraulic fracture propagation. J Porous Media 2017;20(2):155–73. [35] Dontsov EV, Kresse O. A semi-infinite hydraulic fracture with leak-off driven by a power-law fluid. J Fluid Mech 2018;837:210–29. [36] Nilson RH, Proffer WJ, Duff RE. Modelling of gas-driven fractures induced by propellant combustion within a borehole. Int J Rock Mech Min Sci & Geomech Abstr 1985;22(1):3–19. [37] Lu W, Tao Z. A study on fracture propagation velocity driven by gases of explosion products. Explosion Shock Waves 1994;14(3):264–8. (In Chinese). [38] Kirsch G. Die Theorie Der Elastizität Und Die Bedürfnisse Der Festigkeitslehre. Zeitschrift Des Vereins Deutscher Ingenieure 1898; 29: 797–807. [39] Ravi-Chandar K, Knauss WG. An experimental investigation into dynamic fracture: III. On steady-state crack propagation and crack branching. Int J Fract 1984;26(2):141–54. [40] Kou MM, Lian YJ, Wang YT. Numerical investigations on crack propagation and crack branching in brittle solids under dynamic loading using bond-particle model. Eng Fract Mech 2019;212:41–56. [41] Glover AP, Johnson FA, Radon JC. Crack velocity measurements in the toughness testing of polymers. In: Proceedings of International Conference on Dynamic Crack Propagation. Lehigh University; 1972. [42] Freund LB. Dynamic Fracture Mechanics. Cambridge. UK: Cambridge Univ Press; 1990. [43] Abdideh M, Amanipoor H. Fractures and borehole breakouts analysis of a reservoir using an image log (Case study: SW Iran). Liquid Fuels Technol 2012;30:2360–72.

14