Shock and Vibration

Shock and Vibration / 2021 / Article
Special Issue

Technology, Research, and Outlook of Coal Rock Dynamic Disasters in Deep Mining

View this Special Issue

Research Article | Open Access

Volume 2021 |Article ID 9997597 | https://doi.org/10.1155/2021/9997597

Zhiyong Niu, Shiquan Wang, Hongrui Ma, Songbao Feng, Hengjie Luan, Yibin Yang, "A Numerical Investigation of the Hazardous Injection Area of Induced Earthquake during Hydraulic Fracturing", Shock and Vibration, vol. 2021, Article ID 9997597, 21 pages, 2021. https://doi.org/10.1155/2021/9997597

A Numerical Investigation of the Hazardous Injection Area of Induced Earthquake during Hydraulic Fracturing

Academic Editor: Feng Du
Received04 Apr 2021
Accepted20 May 2021
Published28 May 2021

Abstract

Hot dry rock (HDR) geothermal energy has many advantages, such as being renewable, clean, widely distributed, and without time and weather limitations. Hydraulic fracturing is usually needed for the exploitation of HDR geothermal energy. It has many hidden faults in the reservoir/caprock sequences. Injecting fluid into underground formations during hydraulic fracturing often induces fault slip and leads to earthquakes. Therefore, to well understand the induced fault slip and earthquakes is important for the applications and development of HDR geothermal exploitation. In this study, we investigated the hazardous injection area of the induced earthquakes during hydraulic fracturing. The study was based on a hydraulic fracturing test in Qiabuqia geothermal field in China. According to the field, a fault-surrounding rock-fracturing region system was developed to study the influences of fluid injection on the stability of the specific fault. A total of 60 hydraulic fracturing regions and 180 numerical experiments were designed. The results revealed that the hazardous injection regions that threaten the fault’s stability were near to the fault and concentrated on the following four areas: (a) above the top of the fault in underlying strata; (b) above the top of the hanging wall of the fault in underlying strata; (c) near to the fault planes in both footwall and hanging wall; (d) at the bottom of the footwall of the fault in underlying strata. The hazardous injection area can be controlled effectively by adjusting the injection pressure.

1. Introduction

Geothermal energy has great potential to simultaneously solve the energy crisis and climate problems due to its significant advantages. Firstly, the geothermal resource is rich. According to the estimation of the World Energy Council, it has about of geothermal energy (more than 10,000 times the available energy from all total fossil fuels worldwide), which can meet the energy demand for thousands of years according to the present global energy demand [1]. Secondly, it is renewable. It is mainly generated from the radioactive decay of materials within the Earth [2]. Thirdly, it is clean and green. The generation of geothermal energy produces no pollutants. Fourthly, it has no time and weather limitation. It can supply energy throughout the year due to the constant flow of heat from the Earth [3]. In addition, geothermal is widely distributed. Therefore, full exploitation of geothermal resources plays an important role in energy supply and environmental protection for all countries in the world.

According to the temperature and storage depth, geothermal energy can be divided into three categories: shallow geothermal energy, underground thermal water resources, and HDR resources [4]. Among them, HDR geothermal energy has the greatest potential. It is abundant and available almost in any area up to a certain depth in the world. It has about  JHD energy just within an accessible zone less than 10 km, becoming a strategic choice to solve the energy crisis and reduce environmental pollution [5]. Exploiting geothermal energy usually needs heat carrier fluids to bring the underground heat to the surface. Cold water can be injected into the hot dry rocks and then be heated. Heat energy can be then obtained through the heated water. For developing HDR energy economically on a commercial scale, the adequate fluid always needs to be supplied in the targeted HDR masses. However, the natural HDR mass (as shown in Figure 1(a)) is usually compact and impermeable. There is little or no water or steam in it. In addition, due to its low permeability, it does not allow fluid injecting or flowing in it. In order to solve these problems, an Enhanced Geothermal System (EGS) is usually needed [6]. In such systems, the permeability of rocks of the target reservoir is enhanced by using artificial methods (usually hydraulic fracturing). Cold fluid then can be injected into the reservoir and be heated by the hot rocks (Figure 1(b)).The hot water can be extracted to the surface to release the captured heat for power generation or other uses (Figure 1(c)). In the course of the fluid being injected into underground formations, seismicity is frequently induced. This phenomenon has been reported in many countries and regions, such as Australia [7], Germany [8], America [9], France [10], South Korea [11], Philippines [12], Switzerland [13], New Zealand [14], and Japan [15]. Large-scale, big, destructive earthquakes have exerted a great influence on the geothermal industry worldwide [16]. Induced seismicity by industrial activities often raises public concern and worry. Massive fluid injection during hydraulic fracturing and the subsequent heat mining also often induce seismicity, which may delay or hinder the acceptance of EGS. A better understanding of fluid-injection-induced fault reactivation is necessary for improving the public acceptance of EGS.

The massive fluid injected into underground formations can often induce fault reactivation and slip. Several studies have explored this effect, and several different mechanisms have been proposed to explain it [1721].The pore-pressure increase is considered to be the main mechanism to explain the occurrence of induced seismicity in underground fluid-injection projects. Injecting fluid into reservoir decreases the effective stress of faults which will cause them to to slip. [22]. However, the induced fault slip can also occur when the fault does not exist within the subsurface targeted zone (no longer the problem of effective stress) [23]. When the fluid is injected, the injection zone is expanded which may cause stress perturbation. The stress perturbation can be transmitted to the surrounding faults and leads them to slip. In Oklahoma, USA, the injection-linked seismicity occurs up to 35 km away from the injection wells [24]. It revealed that induced earthquakes may extend far away from the hydraulically connected region. In summary, both pore-pressure increase (direct pressure effect) and stress perturbation (indirect pressure effect) due to fluid injection can cause fault slip and induce earthquakes. Therefore, a full assessment of the stability of the fault beyond the targeted injection zones is also very important.

Numerical simulators have become important tools to solve the problems that arise in fluid-injection engineering in the deep subsurface [22, 25, 26]. In this paper, the influence of fluid injection during hydraulic fracturing on the fault’s stability was fully investigated by ABAQUS. The study took the geothermal potential area, i.e., Qiabuqia Geothermal Field (located in the east of Gonghe Basin in Northwest China) as a case. The numerical model was based on the field measurements (e.g. length, orientation, and inclination of faults) and the experimental data in field test and laboratory test (Density, Young’s modulus, Poisson’s ratio, and Porosity, etc.). In this paper, we mainly investigated the indirect pressure effect on the stability of the fault. Therefore, we only explored the stability of faults that exist not within the targeted injection zone (as shown in Figure 2). Our research focused on the spatial extent of induced seismicity as the following three aspects. (1) The influence of the spatial distribution of injection zone on induced fault slip and magnitude of seismicity. For a practical EGS project, faults naturally exist, but the injection zones are completely determined by human-made factors. For a given fault in the site, it is necessary to find out which injection zone is more likely to lead the fault to slip and induce seismicity. (2) The influence of the change of distance between the injection zone and fault on fault slip. (3) The dangerous injection region significantly affects the fault. (4) How to reduce the EGS-induced fault slip and the magnitude of seismicity.

2. Theoretical Background

2.1. Linear Elastic Material Model

In the Qiabuqia Geothermal Field, the HDR formations mainly consist of hard granite. Granite is an elastic material with high elastic modulus and low permeability. The linear elastic material model can be used to describe its mechanical properties and deformation characteristics. Therefore, the simulation of the granite reservoir and fault was based on the linear theory of elasticity. The total stress is defined from the total elastic strain as [27]:, where is the total stress; is the fourth-order elasticity tensor; and is the total elastic strain.

2.2. Porous Elastic Material Model

Hydraulic fracturing is usually used for establishing EGS in hot dry rocks. Artificial fracture networks can be created in the targeted region with high porosity and permeability. The porous elastic material model can be used to describe the fracturing region as shown in the flowing equations [28, 29]:, where is the shear modulus; is the strain; is the stress; and are the pore-pressure and first stress invariant from the prefaulting state, respectively; and are the undrained and drained Poisson’s ratios, respectively; is the Kronecker delta; is the change in fluid mass per unit volume, kg/m3; is the density of the pore fluid in the reference state, kg/m3; is Skempton’s coefficient; and are drained bulk modulus and solid grain’s bulk modulus; and are the unjacketed pore compressibility and pore fluid bulk modulus, respectively; and is porosity.

The hydromechanical coupling is accomplished through the governing equation of mass conservation for the infiltrating pore fluid:, where is the divergence of the water flux; is the direction of fluid mass flow rate; is fluid mass; and is time.

2.3. Contact Formulations for Fault Modeling

The fault and its surrounding rock are two deformable bodies. Finite-sliding formulation can be used to describe the interaction between them. Master-slave contact method was used to accomplish the interaction analysis of the fault and its surrounding rock. In our study, the surface of the fault-surrounding rock was the master surface and the surface of the fault was the slave surface. The slip of the slave surface along the master surface can be solved by the following equation:, where is the slip; is the coordinate of the point on the contact segment; and are the interpolation functions of the point on the segment; and are the tangent and normal, respectively; and is the curvature of the contact segment. For the detailed derivation of this equation, refer to ABAQUS theory manual [27].

2.4. Moment Magnitude Scaling

The energy released by fault slip determines the magnitude of the earthquake. Seismic moment M0 proposed by Aki was used to quantify the released energy, as follows [30]:, where μ is rigidity of the fault, Pa; A is fault rupture area, m2; and d is the average slip of the rupture area, m.

Then, can be used to calculate the moment magnitude scale of the earthquake [3133]:

3. Model Setup

3.1. Conceptual Model

ABAQUS was used as the numerical tool for simulating our fault-surrounding rock-fracturing region system (Figure 2) [27]. A 2D plane strain model was developed to simulate the system. As shown in Figure 3, the size of the basic model was 6,000 m × 3000 m. The length of the fault was 1000 m. Its dip angle was 60°. The distance between the top of the fault and the top boundary of the model was 1000 m. The hydraulic fracturing region was simulated by an elliptic zone with the major axis of 300 m and minor axis of 100 m. In the basic model, we designed 12 fracturing regions (TML1, TM, TMR1, FM, FT, FB, HM, HT, HB, BML1, BM, and BMR1) as shown in Figure 3. Among them, regions TML1, TM, and TMR1 were in the overlying strata of the fault. The horizontal distance between region TM and TML1 or TMR1 was 100 m. Among them, region TM was just above the fault. The distance between the bottom of region TM and the top of the fault was 300 m. Regions TML1, TM, and TMR1 were in underlying strata of the fault. The horizontal distance between region BM and BML1 or BMR1 was 100 m. Among them, region BM was just under the fault. The distance between the top of region TM and the bottom of the fault was 300 m. Regions FM, FT, and FB were in the footwall of the fault. Among them, region HM was in the middle of the footwall. Regions HM, HT, and HB were in the hanging wall of the fault. Among them, region HM was in the middle of the hanging wall. The horizontal distance between region FM and the fault was 500 m. The aim of the design was to investigate the influence of spatial distribution of hydraulic fracturing region on the stability of the fault. The material properties of the model were mainly obtained from the logs and well data from the field as well as laboratory experiments, as shown in Table 1.


PropertyDensity (kg/m3)Young’s modulus (GPa)Poisson’s ratio (—)Porosity (—)Friction coefficient (—)Permeability (m/s)

Granite2700460.14
Hydraulic fracturing region210060.350.361 × 10−7
Fault2100100.3

The left and right sides of the model were fixed in the horizontal direction, and the bottom was fixed in the vertical and horizontal direction. The top was allowed to deform freely.

The stress within the fault-surrounding rock-fracturing region system should be equal to the lithostatic stress before fluid injection for hydraulic fracturing. Pore pressure within fracturing region should be hydrostatic, with the initial value of 0 MPa. The initial stresses should be geostatic by the following steps: (1) the gravity was applied as an external force to obtain the geostatic stresses under the applied boundary conditions of the model; (2) the geostatic stresses were then inputted to the model as initial internal force before gravity; and (3) repeated step 2 until a mechanical equilibrium was reached between the initial stress fields and the gravity.

3.2. Numerical Testing
3.2.1. Fluid Injecting into the Hydraulic Fracturing Regions in the Basic Model

After the geostatic stresses balance, fluid was injected into the ellipse region to simulate hydraulic fracturing. According to the hydraulic fracturing test in Qiabuqia Geothermal Field, we injected fluid at a constant pressure of 100 MPa to the model plane over 216000 s (60 hours or 2.5days).

3.2.2. Location Change of the Fracturing Regions in the Vertical and Horizontal Directions

In order to investigate the influence of the location of the fracturing region on the fault’s stability, we designed more fracturing regions by changing the regions, as shown in Figure 4.

Firstly, we investigated fracturing regions in the overlying strata. We first changed the location of the regions in a vertical direction. By moving the regions TM, TML1, and TMR1 down 200 m, regions TB, TBL1, and TBR1 were obtained, respectively. By moving the regions TM, TML1, and TMR1 up 200 m, regions TT, TTL1, and TTR1 were obtained, respectively. We then changed the location of the regions in a horizontal direction. By moving region TBL1 400 m to the left each time, regions TBL2 and TBL3 were obtained. By moving TML1 400 m to the left each time, regions TML2 and TML3 were obtained. By moving TTL1 400 m to the left each time, regions TTL2 and TTL3 were obtained. Similarly, regions TBR2, TBR3, TMR2, TMR3, TTR2, and TTR3 can be obtained by moving regions TBR1, TMR1, and TBR1 400 m to the right each time.

Secondly, we investigated the location change of the fracturing regions in underlying strata. Using the similar method above, we obtained the regions BBL3, BBL2, BBL1, BB, BBR3, BBR2, BBR1, BTL3, BTL2, BTL1, BT, BTR3, BTR2, BTR1, BML3, BML2, BMR3, and BMR2.

Thirdly, we investigated the location change of the fracturing regions in the footwall and hanging wall of the fault. By moving regions FT, FM, and FB 400 m to the left each time, the regions FTL, FML, and FBL were obtained, respectively. By moving regions FT, FM, and FB 400 m to the right each time, the regions FTR, FMR, and FBR were obtained, respectively. Among them, regions FTR, FMR, and FBR were closer to the fault. By moving regions HT, HM and HB 400 m to the left each time, the regions HTL, HML, and HBL were obtained, respectively. By moving regions HT, HM, and HB 400 m to the right each time, the regions HTR1, HMR1, and HBR1 were obtained, respectively. Among them, regions HTL, HML, and HBLwere closer to the fault.

Above all, a total of 60 hydraulic fracturing regions were deigned (Figure 4).

3.2.3. Injection Pressure Change of the Fracturing Regions

Injection pressure is an important engineering factor that affects the stress perturbation in the fracturing region and adjacent rocks. It directly affected the induced fault slippage and the seismic moment magnitude . Based on the design of the fracturing regions in Figure 3, we change the injection pressure to 120 MPa and 80 MPa to study its influences on fault stability. The injection time remained constant at 60 hours.

Above all, a total of 180 numerical experiments were designed. After the finish of the injection, the slippage of the fault along the interface can be obtained. Based on the slippage, the seismic moment magnitude can be calculated.

4. Results and Discussion

4.1. Fracturing Region in the Overlying Strata of the Fault
4.1.1. Fracturing regions of TM, TML1, and TMR1

Figure 5 shows the fluid-injection-induced fault sliding distance when the fracturing regions were in overlying strata. As we can see, the induced slips of the three injection schemes were different. When the fracturing region was TM or TMR1, the value of induced sliding distance was positive. It revealed that the fault slid upward along the fault plane. When the fracturing region was TML1, however, the value of induced sliding distance was negative. The fault slid downward along the plane. The detailed results are summarized as follows. (1) When the fracturing region was TM (the black line with the square symbol), the fault slid upward along the whole plane. The upper part slid more significantly. At 150 m of the plane, the sliding distance was the largest, about 23.00 mm. From 50 m to 350 m, the sliding distances were all larger than 20 mm with an average sliding distance of 21.67 mm. Along the whole fault plane (0 m–1000 m), the average sliding distance was 10.14 mm. According to equations (5) and (6), an 2.8047 earthquake can be induced. (2) When the fracturing region was TMR1 (the blue line with the triangle symbol), the induced fault slip displayed a similar tendency (upward along the whole plane). However, the significant slip region moved down to the middle part. At 370 m of the plane, the sliding distance was the largest, about 26.42 mm. From 270 m to 470 m, the sliding distances were all larger than 22.93 mm, with an average sliding distance of 22.02 mm. Along the whole plane, the average sliding distance was 14.53 mm, and the of the induced earthquake was 2.9089. (3) When the fracturing region was TML1 (the red line with the circle symbol), the induced fault sliding distance displayed a contrary tendency (downward along the whole plane). The significant slip region still occurred at the upper part of the fault. At 70 m of the plane, the sliding distance was the largest, about 11.26 mm. From 20 m to 120 m, the sliding distances were all larger than 8.33 mm, with an average sliding distance of 10.60 mm. Along the whole plane, the average sliding distance was 6.46 mm with the of an induced earthquake of 2.6742. Note that region TML1 was on the left of TM (over the footwall of the fault) and region TMR1 was on the right of TM (over the hanging wall of the fault). The results may reveal that injecting fluid above the hanging wall of the fault may induce a larger earthquake than injecting fluid above the footwall of the fault.

4.1.2. Fracturing Regions Changed in the Vertical Direction

Based on the regions of TM, TML1, and TMR1, we studied the influence of location change of the fracturing region in the vertical direction on the induced fault sliding distances and earthquakes.

From Figure 6, we can see the induced fault sliding distance when the fracturing region changed from TM to TB and TT. When the fracturing region was changed to TB, the induced fault sliding distance increased. As was shown in the red line with the hollow circle symbol, the largest sliding distance was up to 78.24 mm. From 50 m to 350 m, the average sliding distance was 62.94 mm, over three times that of region TM. It showed that the upper part of the fault would slide violently. Along the whole plane, the average sliding distance was 30.94 mm, and the of the induced earthquake was up to 3.1277. When the fracturing region was directly over the top of the fault, it was likely to induce a felt earthquake. When the fracturing region was changed to TT, the induced fault sliding distance decreased. As we can see from the blue line with the solid circle symbol, the largest sliding distance was reduced to 15.25 mm. Along the whole plane, the average sliding distance was only 5.80 mm, and the of the induced earthquake was 2.6430. Note that the distance between the bottom of region TB and the top of the fault was 100 m, and the distance between the bottom of region TT and the top of the fault was 500 m. When the fracturing region changed from TB to TT, the average induced sliding distance reduced by more than 80 percent (from 30.94 mm to 5.80 mm). It revealed that the induced fault sliding distance and the moment magnitude of earthquake decrease rapidly with an increase of the distance between the region and fault.

From Figure 7, we can see the induced fault sliding distance when the fracturing region changed from TML1 to TBL1 and TTL1. When the fracturing region changed to TBL1 (the red line with the hollow circle symbol), the induced fault sliding distance also increased. The largest sliding distance was 16.24 mm, and the average sliding distance was 9.74 mm. The of the induced earthquake was 2.6742. When the fracturing region changed to TTL1 (the blue line with the solid circle symbol), the induced fault sliding distance decreased. The largest sliding distance was only 5.52 mm, with an average sliding distance of 3.44 mm. The of the induced earthquake was 2.4917. When the fracturing region changed from TBL1 to TTL1, the average induced sliding distance reduced by about 65 percent (from 9.74 mm to 3.44 mm). The decreased range was less than that when the fracturing region changed from TB to TT.

From Figure 8, we can see the induced fault sliding distance when the fracturing region changed from TMR1 to TBR1 and TTR1. When the fracturing region changed from TMR1 to TBR1, the sliding mode became different. As shown from the red line with the hollow circle symbol in Figure 8, the fault slip along the whole plane was presented as two parts: sliding downward from 0 m to 140 m and upward from 140 m to 1000 m. It caused the fault to be compressed at 140 m. Besides the change in the sliding mode, the amount of sliding distance changed obviously. The largest sliding distance was up to 44.37 mm, occurring at 470 m. The average sliding distance along the whole plane was up to 22.12 and the of the induced earthquake was up to 3.0305. When the fracturing region changed from TMR1 to TBR1, the sliding mode remained unchanged. The fault slid upward along the whole plane. The largest sliding distance was 10.92 mm, and the average sliding distance was 6.85 mm. The of the induced earthquake was 2.6911. When the fracturing region changed from TBR1 to TTR1, the average induced sliding distance reduced by about 70 percent (from 22.12 mm to 6.85 mm). The decreased range was also less than that when the fracturing region changed from TB to TT.

The above results revealed that the induced sliding distance decreased with the vertical distance between the fault and fracturing region. The influence of vertical variation of the fracturing region on the fault slip was largest when the regions were above the top of the fault.

4.1.3. Fracturing Regions Changed in the Horizontal Direction

In this section, we studied the influence of horizontal variation of the fracturing region’s location on induced fault slip when the fracturing regions were in overlying strata of the fault.

Figure 9 shows the induced fault sliding distance of the fracturing region of TBL1-TBL3 and TBR1-TBR3. When the fracturing region shifted away from TBL1 (black dotted line with the square symbol) to TBL2 (red line with the circle symbol) and TBL3 (blue line with the triangle symbol), the tendency of induced fault slip remained similar, but the sliding distance significantly reduced. When the fracturing region was TBL2, the largest sliding distance was reduced to 3.15 mm, and the average sliding distance was only 2.52 mm. The of the induced earthquake was 2.4016. When the fracturing region was TBL3, the largest sliding distance was reduced to 0.64 mm, and the average sliding distance was 0.49 mm. The of the induced earthquake was only 1.9275. When the fracturing region shifted away from TBR1 (green dotted line with the hollow square symbol) to TBR2 (purple line with the hollow circle symbol) and to TBR3 (orange line with the hollow triangle symbol), however, the induced fault slip presented different tendencies. When the fracturing region shifted away from TBR1 to TBR2, the slip along plane was also presented in two sections. However, the zone of sliding downward of the fault became larger. As we can see, the fault slid downward from 0 m to 600 m. It resulted in the fault to be compressed at 600 m. When the fracturing region was shifted to TBR3, however, the fault slid downward along the whole plane.When the fracturing region was TBR2, the largest sliding distance was 5.60 mm, and the average sliding distance was 3.22 mm. The of the induced earthquake was only 2.4726. When the fracturing region was TBR3, the largest sliding distance was 3.50 mm, and the average sliding distance was 2.48 mm. The of the induced earthquake was only 2.3970. We knew that when the fracturing region TBL1 and TBR1 shifted away from the fault in the horizontal direction, the induced fault slip and seismicity both decreased.The horizontal variation of the fracturing region can affect both the slip pattern and the amount of sliding distance. Comparing the results of TBR2 with TBL2 (or TBR3 with TBR3), we knew that injecting fluid over the hanging wall of the fault was more likely to cause a larger fault slip than injecting fluid over the footwall under the same conditions.

Figure 10 shows the induced fault sliding distance of the fracturing region of TML1-TML3 and TMR1-TMR3. When the fracturing region shifted away from TML1 (black dotted line with the square symbol) to TML2 (red line with the circle symbol) and TML3 (blue line with the triangle symbol), the variation characteristics of the induced fault slip were similar to that when the fracturing region shifted away from TBL1 to TBL2 and TBL3. When the fracturing region was TML2, the largest sliding distance was 3.39 mm, and the average sliding distance was 2.51 mm. The of the induced earthquake was only 2.4005. When the fracturing region was TML3, the largest sliding distance was 0.62 mm, and the average sliding distance was 0.50 mm. The of the induced earthquake was 1.9333. When the fracturing region shifted away from TMR1 (green dotted line with the hollow square symbol) to TMR2 (purple line with the hollow circle symbol) and TMR3 (orange line with the hollow triangle symbol), the variation characteristics of the induced fault slip were also similar to that when the fracturing region shifted away from TBR1 to TBR2 and TBR3. When the fracturing region was TMR2, the largest sliding distance was 4.10 mm, and the average sliding distance was 2.65 mm. The of the induced earthquake was 2.4162. When the fracturing region was TMR3, the largest sliding distance was 2.72 mm, and the average sliding distance was 1.88 mm. The of the induced earthquake was 2.3168. We knew that when the fracturing region TML1 and TMR1 shifted away from the fault in the horizontal direction, the induced sliding distance and seismicity decreased.The horizontal variation of the fracturing region would affect both the slip pattern and the amount of sliding distance. Comparing the results of TMR2 with TML2 (or TMR3 with TMR3), we knew that injecting fluid over the hanging wall of the fault was more likely to cause a larger fault slip than injecting fluid over the footwall under the same conditions.

Figure 11 shows the induced fault sliding distance of the fracturing region of TTL1-TTL3 and TTR1-TTR3. When the fracturing region shifted away from TTL1 (black dotted line with the square symbol) to TTL2 (red line with the circle symbol) and TTL3 (blue line with the triangle symbol), the variation characteristics of the induced fault slip were also similar to that when the fracturing region shifted away from TBL1 to TBL2 and TBL3. When the fracturing region was TTL2, the largest sliding distance, average sliding distance, and of the induced earthquake were 2.64 mm, 1.80 mm, and 2.3042, respectively.When the fracturing region was TTL3, the largest sliding distance, average sliding distance, and of the induced earthquake were 0.44 mm, 0.35 mm, and 1.8301, respectively. When the fracturing region was TTR2, the largest sliding distance, average sliding distance, and of the induced earthquake were 2.45 mm, 1.85 mm, and 2.3121, respectively. When the fracturing region was TTR3, the largest sliding distance, average sliding distance, and of the induced earthquake were 1.71 mm, 1.15 mm, and 2.1745, respectively. We knew that when the fracturing region TTL1 and TTR1 shifted away from the fault in the horizontal direction, the induced sliding distance and seismicity decreased. The horizontal variation of the fracturing region would affect both the slip pattern and the amount of sliding distance. Comparing the results of TTR2 with TTL2 (or TTR3 with TTR3), we knew that injecting fluid over the hanging wall of the fault was more likely to cause a larger fault slip than injecting fluid over the footwall under the same conditions.

The above results revealed that, when the fracturing regions were in overlying strata of the fault, the horizontal variation of the fracturing region may affect both the slip pattern and the amount of the sliding distance of the fault. The sliding distance decreased along with the increase of the horizontal distance between the fault and fracturing region. Injecting fluid over the hanging wall of the fault was more likely to cause a larger fault slip than injecting fluid over the footwall under the same conditions. Therefore, for the overlying strata of fault, the area close to the top of the hanging wall of the fault was the more hazardous injection area.

4.2. Fracturing Region in the Footwall and Hanging Wall of the Fault

Figure 12 shows the induced fault sliding distance when the fracturing region was in the footwall and hanging wall of the fault. The horizontal distances between the regions (FT, FM, FB, HT, HM, or HB) and the fault were 500 m. When the fracturing region was FT (black line with the square symbol), the fault slid downward along the whole plane. The largest sliding distance was 5.66 mm, and the average sliding distance was 4.37 mm. The of the induced earthquake was 2.5610. When the fracturing region was HT (green line with the hollow square symbol), the induced fault slip along the whole plane presented in two sections: from 0 m to 810 m, fault slid downward; from 810 m to 100 m, fault slid upward. It caused the fault to be compressed at 810 m. The largest sliding distance, average sliding distance, and of the induced earthquake were 8.62 mm, 5.50 mm, and 2.6276, respectively. When the fracturing region was FM (red line with the circle symbol), the fault also slid downward along the whole plane. The largest sliding distance, average sliding distance, and of the induced earthquake were 7.50 mm, 5.19 mm, and 2.6108, respectively. When the fracturing region was HM (purple line with the hollow circle symbol), the fault slid downward along the whole plane. The largest sliding distance, average sliding distance, and of the induced earthquake were 8.51 mm, 5.91 mm, and 2.6484, respectively. When the fracturing region was FB, the fault slip along the whole plane became complicated. As shown in the blue line with the triangle symbol in Figure 12, the fault slid upward from 0 m to 240 m and downward from 240 m to 1000 m. This caused the fault to be tensioned at the point of 240 m, which was prone to failure. The largest fault sliding distance was 6.04 mm, and the average sliding distance was 3.94 mm. The of the induced earthquake was 2.5310. When the fracturing region was HB, the fault also slid downward along the whole plane (orange line with the hollow triangle symbol in Figure 12). The largest sliding distance, average sliding distance, and of the induced earthquake were 5.38 mm, 4.81 mm, and 2.5888, respectively.

When the horizontal distances between the regions and the fault were reduced to 100 m, the induced fault sliding distances were shown in Figure 13. The sliding pattern and the amount of sliding distance of fault both changed remarkably. When the fracturing region was FTR (located on the right of region FT), the induced fault slip along the whole plane (black line with the square symbol) presented in two sections: from 0 m to 60 m, the fault slid upward; from 60 m to 1000 m, the fault slid downward. It caused the fault to be tensioned at the point of 60 m and prone to failure. The largest fault sliding distance was 33.81 mm, and the average sliding distance was 21.43 mm. The of the induced earthquake was 3.0214. When the fracturing region was FMR (red line with the circle symbol), the fault slip also presented in two sections, but the range of sliding upward increased (0 m–380 m). The largest sliding distance, average sliding distance, and of the induced earthquake were 75.35 mm, 32.38 mm, and 3.1409, respectively. When the fracturing region was FBR (blue line with the triangle symbol), the range of fault that slid upward continuously increased (0 m–660 m). The largest sliding distance, average sliding distance, and of the induced earthquake were 123.32 mm, 60.01 mm, and 3.3195, respectively. When the fracturing region was HTL (located on the left of region FT), the induced fault slip showed an opposite trend (green line with the hollow square symbol). The upper part of the fault slid downward (0 m–370 m) and the lower part of fault slid upward (370 m–1000 m). It caused the fault to be compressed at the point of 370 m. The largest sliding distance, average sliding distance, and of the induced earthquake were 115.25 mm, 45.77 mm, and 3.2411, respectively. When the fracturing region was HML (purple line with the hollow circle symbol), the range of sliding upward of the fault increased to 0 m–370 m. The largest sliding distance, average sliding distance, and of the induced earthquake were 52.67 mm, 22.75 mm, and 3.0387, respectively. When the fracturing region was HBL (the orange line with the hollow triangle symbol), the induced fault slip presented only one section (downward). The largest sliding distance, average sliding distance, and of the induced earthquake were 32.84 mm, 19.92 mm, and 3.0002, respectively.

From the above results, we knew that when the fracturing region was close to the fault, no matter it was in the footwall or hanging wall, injecting fluid into it can induce a felt earthquake. For the fracturing region in the footwall of the fault, the induced slip became larger when the region was closer to the lower part of the fault. For the fracturing region in the hanging wall of the fault, on the contrary, the induced slip became larger when the region was closer to the upper part of the fault. It may reveal that the regions that were close to the lower part of the fault in the footwall and to the upper part of the fault in the hanging wall were more dangerous areas.

When the horizontal distances between the fracturing regions and the fault were increased to 900 m, the induced fault sliding distances were shown in Figure 14. With the fracturing region becoming far away from the fault, the sliding pattern of the fault became simple. When the fluid was injected into the six different regions (FTL, FML, FBL, HTR, HMR, or HBR), the fault all slid downward along the whole plane. When the fracturing region was FTL, the largest sliding distance, average sliding distance, and of the induced earthquake were 1.35 mm, 1.12 mm, and 2.1668, respectively. When the fracturing region was FML, the largest sliding distance, average sliding distance, and of the induced earthquake were 2.56 mm, 11.74 mm, and 2.2944, respectively. When the fracturing region was FBL, the largest sliding distance, average sliding distance, and of the induced earthquake were 3.08 mm, 1.86 mm, and 2.3137, respectively. When the fracturing region was HTR, the largest sliding distance, average sliding distance, and of the induced earthquake were 4.50 mm, 3.21 mm and 2.4717, respectively. When the fracturing region was HMR, the largest sliding distance, average sliding distance, and of the induced earthquake were 3.41 mm, 2.60 mm, and 2.4107, respectively. When the fracturing region was HBR, the largest sliding distance, average sliding distance, and of the induced earthquake were 2.40 mm, 2.16 mm, and 2.3570, respectively. The average sliding distances of six different fracturing regions were all smaller than 3.21 mm. From the above results, we knew that when the fracturing region was far away from the fault, the induced fault slip and seismicity were both small. When the fracturing region was in the footwall of the fault, the induced slip also became larger if the region was closer to the lower part of the fault. When the fracturing region was in the hanging wall of the fault, the induced slip became larger if the region was closer to the upper part of the fault. The same conclusion can be obtained as in the previous paragraph.

4.3. Fracturing Region in Underlying Strata of Fault
4.3.1. Fracturing Regions of BM, BML1, and BMR1

Figure 15 shows the fluid-injection-induced fault sliding distance when the fracturing regions were in underlying strata. When the fracturing region was BM (the black line with the square symbol), the fault slid upward along the whole plane, and the sliding distance along the lower part was large. From 500 m to 1000 m, the average sliding distance was 11.30 mm. From 0 m to 500 m, the average sliding distance was only 5.65 mm. Along the whole fault plane, the average sliding distance was 8.33 mm. The of the induced earthquake was 2.7478. When the fracturing region was BML1 (the red line with the square symbol), the fault also slid upward along the whole plane, and the sliding distance became larger than that of BM. The sliding distance along the middle section of the fault plane was much larger. The largest fault sliding distance was 39.08 mm, occurring at 510 m. The average sliding distance was 24.64 mm, and the of the induced earthquake was 3.0618. When the fracturing region was BMR1 (the blue line with the triangle symbol), the fault slid downward along the whole plane. The largest sliding distance, average sliding distance, and of the induced earthquake were 8.82 mm, 5.75 mm, and 2.6405, respectively. The results revealed that injecting fluid under the footwall of the fault may induce a much larger fault sliding distance than injecting fluid under hanging wall of the fault.

4.3.2. Fracturing Regions Changed in the Vertical Direction

We then studied the influence of location change of the fracturing regions in the vertical direction on the induced fault sliding distances and earthquakes.

Figure 16 shows the induced fault sliding distance when the fracturing region changed from BM to BT and BB. When the fracturing region changed to BT (closer to the fault), the induced fault sliding distance increased. As shown from the blue line with the solid square symbol, the largest sliding distance was 23.68 mm (at 970 m). The average sliding distance was 9.84 mm, and the of the induced earthquake was 2.7960. When the fracturing region changed to BB (farther to the fault), the induced fault sliding distance decreased (red line with the hollow square symbol). The largest sliding distance was reduced to 10.10 mm. The average sliding distance was reduced to 6.27 mm, and the of the induced earthquake was 2.6656. When the fracturing region changed from BT to BB, the average induced sliding distance reduced by 36 percent (from 9.84 mm to 6.27 mm).

Figure 17 shows the induced fault sliding distance when the fracturing region changed from BML1 to BTL1 and BBL1. When the region changed to BTL1 (closer to the fault), the induced fault sliding distance increased (the blue line with the solid circle symbol). As shown, the largest sliding distance was increased to 59.72 mm (at 490 m). The average sliding distance was 35 mm, and the of the induced earthquake was 3.1634. When the fracturing region changed to BBL1 (farther to the fault), the induced fault sliding distance decreased (red line with the hollow circle symbol). The largest sliding distance was reduced to 25.40 mm. The average sliding distance was 15.97 mm, and the of the induced earthquake was reduced to 2.9362. When the fracturing region changed from BTL1 to BBL1, the average induced sliding distance reduced by 54 percent (from 9.84 mm to 6.27 mm).

Figure 18 shows the induced fault sliding distance when the fracturing region changed from BMR1 to BTR1 and to BBR1.When the fracturing region changed to BTR1 (closer to the fault), the induced fault sliding distance increased. As shown from the blue line with the solid triangle symbol, the largest sliding distance was increased to 12.53 mm (at 830 m). The average sliding distance was increased to 8.46 mm, and the of the induced earthquake was 2.7523. When the fracturing region changed from BMR1 to BBR1 (becoming farther to the fault), the induced fault sliding distance decreased (red line with the hollow triangle symbol). The largest sliding distance was reduced to 4.31 mm. The average sliding distance was reduced to 3.03 mm, and the of the induced earthquake was 2.4550. When the fracturing region changed from BTR1 to BBR1, the average induced sliding distance reduced by 52 percent (from 9.84 mm to 6.27 mm).

From the above results, we knew that when region BM, BML1, and BMR1 moved up 200 m to BT, BTL1, and BTL1, respectively, the induced fault slip and seismicity all increased. When region BML1 moved up to BTL1, the increase was the largest. Taking the average induced sliding distance, for example, it increased by 10.36 mm (from 24.64 mm to 35.00 mm). It revealed that the influence of vertical variation of the fracturing region’s location on fault slip was largest when the regions were under the footwall of the fault. For underlying strata of the fault, the area that was close to the bottom of the footwall of the fault was a hazardous injection area.

4.3.3. Fracturing Regions Changed in the Horizontal Direction

In this section, we studied the influence of the fracturing region changing in the horizontal direction on induced fault slip when the fracturing regions were in underlying strata of the fault.

Figure 19 shows the induced fault sliding distance of the fracturing region of BBL1-BBL3 and BBR1- BBR3. As shown, when the fracturing region shifted away from BBL1 (black dotted line with the square symbol) to BBL2 (red line with the circle symbol), the induced fault slip presented a similar tendency (upward). However, the amount of the sliding distance changed significantly. The largest sliding distance was reduced to 10.05 mm, occurring at 370 m. The average sliding distance was 6.00 mm. The of the induced earthquake was 2.6528. When the fracturing region changed to BBL3 (blue line with the triangle symbol), the induced fault slip turned into the opposite tendency (downward). The largest sliding distance was only 3.24 mm, and the average sliding distance was 2.67 mm. The of the induced earthquake was 2.4184. When the fracturing region shifted away from BBR1 (green dotted line with the hollow square symbol) to BBR2 (purple line with the hollow circle symbol), the induced fault slip also presented a similar tendency (downward). The largest sliding distance was 4.61 mm, and the average sliding distance was 3.40 mm. The of the induced earthquake was 2.4883. When the fracturing region changed to BBR3 (orange line with the hollow triangle symbol), the induced fault slip became complicated and presented two sections: from 0 m to 600 m, the fault slid upward; from 600 m to 1000 m, it slid downward. The largest sliding distance was 2.69 mm, and the average sliding distance was 1.14 mm. The of the induced earthquake was 2.1720. We knew that when the fracturing region BBL1 and BBR1 shifted away from the fault in the horizontal direction, the induced sliding distance and seismicity decreased. The horizontal variation of the fracturing region would affect both the slip pattern and the amount of sliding distance. Comparing the results of BBR2 with BBL2 (or BBR3 with BBR3), we knew that injecting fluid under the footwall was more likely to cause a larger fault sliding distance than injecting fluid under the hanging wall of the fault under the same conditions.

Figure 20 shows the induced fault sliding distance of the fracturing region of BML1-BML3 and BMR1-BMR3. When the fracturing region shifted away from BML1 (black dotted line with the square symbol) to BML2 (red line with the circle symbol), the induced fault slip became complicated and presented two sections: from 0 m to 800 m, the fault slid upward; from 800 m to 1000 m, the fault slid downward. The largest sliding distance, average sliding distance, and of the induced earthquake were 8.70 mm, 4.56 mm, and 2.5733, respectively. When the fracturing region was BML3, the induced fault slip (blue line with the triangle symbol) presented an opposite tendency with that when the region was BML1. The fault slid downward along the whole plane and the maximum sliding distance occurred at 610 m, about 2.93 mm. The average sliding distance was 2.45 mm, and the of the induced earthquake was 2.3935. When the fracturing region shifted away from BMR1 (green dotted line with the hollow square symbol) to BMR2 (purple line with the hollow circle symbol), the induced fault slip presented a similar tendency (downward). The largest sliding distance was 4.85 mm, and the average sliding distance was 3.86 mm. The of the induced earthquake was 2.5251. When the fracturing region changed to BMR3 (orange line with the hollow triangle symbol), the induced fault slip also presented two sections: from 0 m to 380 m, the fault slid upward; from 380 m to 1000 m, the fault slid downward. The largest sliding distance was 1.60 mm, and the average sliding distance was 0.90 mm. The of the induced earthquake was 2.1350. We knew that when the fracturing region BML1 and BMR1 shifted away from the fault in the horizontal direction, the induced sliding distance and seismicity decreased.The horizontal variation of the fracturing region would affect both the slip pattern and the amount of sliding distance. Comparing the results of BMR2 with BBL2 (or BMR3 withBML3), We knew that injecting fluid under the footwall was more likely to cause a larger fault sliding distance than injecting fluid under the hanging wall of the fault under the same conditions.

Figure 21 shows the induced fault sliding distance of the fracturing region of BTL1-BTL3 and BTR1-BTR3. When the fracturing region shifted away from BTL1 (black dotted line with the square symbol) to BTL2 (red line with the circle symbol), the induced fault slip also became complicated and presented two sections: from 0 m to 570 m, the fault slid upward; from 570 m to 1000 m, the fault slid downward. The largest sliding distance, average sliding distance, and of the induced earthquake were 9.13 mm, 4.01 mm, and 2.5361, respectively. When the fracturing region was BML3 (blue line with the triangle symbol), the fault slid downward along the whole plane. The maximum sliding distance was 2.47 mm. The average sliding distance was 2.10 mm, and the of the induced earthquake was 2.3488. When the fracturing region shifted away from BTR1 (green dotted line with the hollow square symbol) to BMR2 (purple line with the hollow circle symbol), the induced fault slip presented a similar tendency (downward). The largest sliding distance was 4.41 mm, and the average sliding distance was 3.64 mm. The of the induced earthquake was 2.5081. When the fracturing region changed to BTR3 (orange line with the hollow triangle symbol), the induced fault slip also presented two sections: from 0 m to 140 m, the fault slid upward; from 140 m to 1000 m, the fault slid downward. The largest sliding distance was 1.83 mm, and the average sliding distance was 1.10 mm. The of the induced earthquake was 2.1616. We knew that when the fracturing regions BTL1 and BTR1 shifted away from the fault in the horizontal direction, the induced sliding distance and seismicity decreased. The horizontal variation of the fracturing region would affect both the slip pattern and the amount of sliding distance. Comparing the results of BTR2 with BBL2 (or BTR3 withBML3), we knew that injecting fluid under the footwall was more likely to cause a larger fault sliding distance than injecting fluid under the hanging wall of the fault under the same conditions.

The above results revealed that, when the fracturing regions were in underlying strata of the fault, the horizontal variation of the fracturing region may affect both the slip pattern and the amount of sliding distance. The sliding distance decreased along with the increase of the horizontal distance between the fault and fracturing region. Injecting fluid under the footwall of the fault was more likely to cause a larger fault slip than injecting fluid under the hanging wall of the fault under the same conditions. Therefore, for underlying strata of fault, the area close to the bottom of footwall of the fault was the more hazardous injection area.

4.4. Hazardous Injection Area

Risk assessment of the induced earthquakes in the process of hydraulic fracturing for establishing EGS is essential. Therefore, systematic research on seismic activities and spatial evolution characteristics is needed. Hazardous injection areas of induced felt earthquakes can be identified by numerical simulation. It is helpful for the seismic disaster prevention of such underground fluid-injection engineering. Based on the above analysis, the results of induced fault slip and seismicity were given in a tabular form, as shown in Table 2. The induced earthquakes with larger than 3.0 can be felt by people and may cause public panic. Therefore, we selected the fracturing region that induced earthquake as hazardous injection region. The hazardous injection area of the specific fault can be obtained by connecting all the hazardous regions. When the injection pressure was 100 MPa, the hazardous injection regions were TB, TBR1, FTR1, FMR1, RBR1, HTL1, HML1, HBL1, BTL1, and BML1. According to the ten hazardous injection regions, the hazardous injection area can be drawn, as shown in Figure 22. It revealed that when the injection pressure was 100 MPa, the hazardous injection regions were near the fault and concentrated on four areas: (1) above the top of the fault in underlying strata; (2) above the top of the hanging wall of the fault in underlying strata; (3) near the fault planes in both footwall and hanging wall; and (4) the bottom of the footwall of the fault in underlying strata.


Fracturing regionMaximal slippage (mm)Average slippage (mm)Moment magnitude ()

TTL15.523.442.4917
TTL22.641.802.3042
TTL30.440.351.8301
TT15.255.802.6430
TTR110.926.852.6911
TTR22.451.852.3121
TTR31.711.152.1745
TML111.266.462.6742
TML23.392.512.4005
TML30.620.501.9333
TM22.7210.142.8047
TMR126.4214.532.9089
TMR24.102.652.4162
TMR32.721.882.3168
TBL116.249.742.7931
TBL23.152.522.4016
TBL30.640.491.9275
TB78.2430.943.1277
TBR144.3722.123.0305
TBR25.603.222.4726
TBR33.502.482.3970
FTL11.351.122.1668
FT5.664.372.5610
FTR133.8121.433.0214
FML12.561.742.2944
FM7.505.192.6108
FMR175.3532.383.1409
FBL13.081.742.2944
FB6.043.942.5310
FBR1123.3260.013.3195
HTL1115.2545.773.2411
HT8.625.502.6276
HTR14.503.212.4717
HML152.6722.753.0387
HM8.515.912.6484
HMR13.412.602.4107
HBL132.8419.923.0002
HB5.384.812.5888
HBR12.402.162.3570
BTL159.7235.003.1634
BTL29.134.012.5361
BTL32.472.102.3488
BT23.689.842.7960
BTR112.537.152.7036
BTR24.413.642.5081
BTR31.831.102.1616
BML139.0824.643.0618
BML28.704.562.5733
BML32.932.452.3935
BM13.418.332.7478
BMR18.825.752.6405
BMR24.853.862.5251
BMR31.600.902.1035
BBL125.4015.972.9362
BBL210.056.002.6528
BBL33.242.672.4184
BB10.106.272.6655
BBR14.313.432.4909
BBR24.613.402.4883
BBR32.691.142.1720

When the injection pressure was increased to 120 MPa, the induced maximal slippage, average slippage, and moment magnitude of all fracturing regions were as shown in Table 3. According to the results, the hazardous injection area can be obtained (Figure 23). As shown, it contained 12 hazardous injection regions (TB, TBR1, TMR1, FTR1, FMR1, RBR1, HTL1, HML1, HBL1, BTL1, BML1, and BBL1). The hazardous injection area was enlarged with the increase of injection pressure. The spatial distribution of the regions was almost the same as that when the injection pressure was 100 MPa.


Fracturing regionMaximal slippage (mm)Average slippage (mm)Moment magnitude ()

TTL16.59 (80 m)4.102.5425
TTL23.16 (110 m)2.152.3556
TTL30.52 (420 m)0.421.8829
TT18.30 (120 m)6.962.6958
TTR114.05 (400 m)8.542.7550
TTR22.93 (90 m)1.452.2416
TTR32.05 (320 m)1.372.2252
TML113.56 (70 m)7.772.7276
TML24.06 (160 m)3.002.4521
TML30.74 (490 m)0.591.9813
TM34.05 (140 m)14.752.9132
TMR134.70 (370 m)19.923.0002
TMR24.90 (130 m)2.672.4184
TMR33.25 (370 m)2.252.3688
TBL119.51 (85 m)11.702.8461
TBL23.78 (310 m)3.022.4540
TBL30.76 (630 m)0.591.9813
TB102.62 (70 m)40.973.2090
TBR157.25 (470 m)28.953.1085
TBR26.74 (230 m)3.862.5251
TBR34.20 (460 m)2.982.4502
FTL11.62 (60 m)1.342.2188
FT6.80 (65 m)5.242.6136
FTR140.77 (330 m)25.783.0749
FML13.07 (110 m)2.082.3461
FM8.99 (240 m)6.232.6637
FMR190.63 (220 m)38.943.1943
FBL13.69 (310 m)2.602.4107
FB7.29 (510 m)4.762.5858
FBR1146.67 (410 m)72.063.3725
HTL1138.94 (590 m)55.843.2987
HT10.33 (505 m)6.212.6627
HTR15.40 (685 m)3.852.5243
HML162.65 (800 m)27.323.0917
HM10.21 (760 m)7.092.7011
HMR14.10 (815 m)3.122.4635
HBL139.29 (830 m)23.613.0494
HB6.41 (875 m)5.792.6425
HBR12.89 (500 m)2.602.4107
BTL174.91 (490 m)44.243.2312
BTL211.20 (150 m)4.972.5983
BTL32.20 (620 m)1.322.2144
BT27.94 (970 m)12.372.8623
BTR115.06 (830 m)10.072.8027
BTR25.30 (600 m)4.382.5617
BTR32.96 (535 m)2.522.4016
BML150.91 (540 m)32.393.1410
BML211.27 (310 m)5.862.6460
BML31.48 (770 m)1.082.1563
BM17.13 (630 m)10.282.8087
BMR110.63 (895 m)6.892.6928
BMR25.83 (750 m)4.642.5784
BMR33.51 (615 m)2.942.4463
BBL134.34 (560 m)21.743.0255
BBL212.71 (280 m)7.632.7224
BBL33.23 (150 m)1.372.2252
BB12.16 (570 m)7.472.7162
BBR15.19 (910 m)4.322.5577
BBR25.54 (770 m)3.652.5089
BBR33.89 (730 m)3.212.4717

When the injection pressure was decreased to 80 MPa, the induced maximal slippage, average slippage, and moment magnitude of all fracturing regions were as shown in Table 4. According to the results, the hazardous injection area can be obtained (Figure 24). As shown, it contained only five hazardous injection regions (TB, FMR1, RBR1, HTL1, and BTL1). The area was reduced with the decrease of injection pressure. When the injection pressure was 80 MPa, the hazardous injection regions were all near the fault and concentrated on the following areas:(1) above the top of the fault in underlying strata; (2) near the lower part of the fault in the footwall; (3) near the upper part of the fault in the hanging wall; and (4) under the bottom of the footwall of the fault in underlying strata.


Experiment no.Fault F1
Maximal slippage (mm)Average slippage (mm)Moment magnitude ()

TTL14.46 (75 m)2.782.4300
TTL22.12 (110 m)1.442.2396
TTL30.36 (420 m)0.291.7756
TT12.20 (120 m)4.642.5784
TTR18.77 (400 m)5.502.6276
TTR21.96 (90 m)0.972.1252
TTR32.05 (320 m)0.922.1099
TML18.97 (70 m)5.152.6086
TML22.72 (160 m)2.012.3362
TML30.50 (490 m)0.401.8687
TM2058 (140 m)8.702.7604
TMR119.15 (370 m)10.602.8176
TMR23.29 (130 m)1.802.3042
TMR32.18 (370 m)1.512.2533
TBL112.98 (85 m)7.792.7284
TBL22.52 (310 m)2.012.3362
TBL30.51(630 m)0.391.8614
TB54.66 (70 m)21.323.0199
TBR131.63 (470 m)15.482.9272
TBR24.48 (230 m)2.582.4084
TBR32.80 (460 m)1.992.3333
FTL11.08 (60 m)0.892.1003
FT4.53 (65 m)3.412.4892
FTR126.88 (330 m)17.082.9557
FML12.05 (110 m)1.392.2294
FM6.01 (240 m)4.162.5467
FMR159.91 (220 m)25.823.0753
FBL12.47 (310 m)1.742.2944
FB4.80 (510 m)3.062.4578
FBR197.51 (410 m)47.053.2491
HTL188.70 (590 m)34.973.1632
HT6.91 (505 m)4.162.5467
HTR13.60 (685 m)2.562.4062
HML138.93 (800 m)17.482.9624
HM6.81 (760 m)4.732.5839
HMR12.73 (815 m)2.082.3461
HBL126.36 (830 m)15.712.9315
HB4.30 (875 m)3.842.5236
HBR11.92 (500 m)1.722.2910
BTL144.58 (490 m)25.853.0757
BTL27.34 (150 m)3.202.4708
BTL31.47 (620 m)0.882.0970
BT19.51 (970 m)7.932.7335
BTR110.01 (830 m)6.762.6873
BTR23.52 (600 m)2.912.4433
BTR31.97 (535 m)1.682.2842
BML127.56 (540 m)17.112.9562
BML26.82 (150 m)3.462.4934
BML31.28 (770 m)0.722.0389
BM10.30 (630 m)6.622.6813
BMR17.03 (895 m)4.562.5733
BMR23.88 (750 m)3.092.4607
BMR32.34 (615 m)1.962.3289
BBL116.71 (560 m)10.582.8170
BBL27.35 (280 m)4.562.5733
BBL32.15 (150 m)0.912.1067
BB8.05 (570 m)5.102.6057
BBR13.44 (910 m)2.422.3899
BBR23.68 (770 m)2.872.4393
BBR32.59 (730 m)2.142.3543

5. Conclusions

Injecting fluid into underground formations can exercise a great influence on the stability of the natural faults. When the fault is connected to the injection area, the main cause of fault reactivation and slip is the reduction of effective stress and the weakening of the physical and mechanical properties of the fault induced by fluid injection. When the fault is not connected to the injection area, the fault can also be reactivated by the solid stress perturbation transmitted from the injection area. Therefore, for an EGS project, it is necessary to estimate the stability of the fault that is not connected to the injection area during hydraulic fracturing. In this paper, the influence of fluid injection during hydraulic fracturing on fault’s stability was investigated by a 2D numerical model. The study took Qiabuqia Geothermal Field as a case. The research was focused on the following aspects: (1) the effect of the spatial distribution of the injection zone on induced fault slip and the magnitude of seismicity; (2) the influence of the change of distance between the injection zone and the fault on fault slip; (3) the dangerous injection area that threatens the fault’s stability; and (4) the methods to reduce the fluid-injection-induced fault slip and the seismicity. The results revealed that the spatial distribution of the fracturing region had a significant influence on the induced fault slip and of an induced earthquake. Delimiting hazardous injection area is very helpful for the selection of the location of the hydraulic fracturing region. It is useful for the prevention of fluid-injection-induced earthquakes. This study produced the following conclusions:(1)When the fracturing region was in the overlying strata of fault, the horizontal variation of the region may affect both the slip pattern and the amount of the sliding distance. The induced sliding distance decreased with the increase of the horizontal distance between the fault and fracturing region. It also decreased with the vertical distance. The influence of vertical variation of the fracturing region’s location on fault slip was the largest when the regions were above the top of the fault. Injecting fluid over the hanging wall of the fault was more likely to cause a larger fault slip than injecting fluid over the footwall under the same conditions. For overlying strata of the fault, the area that closes to the top of the hanging wall of the fault was the more hazardous injection area.(2)When the fracturing regions were in underlying strata of the fault, the horizontal variation of the fracturing region may also affect both the slip pattern and the value of the sliding distance of the fault. The sliding distance decreased along with the horizontal distance (or vertical distance) between the fault and fracturing region. Injecting fluid under the footwall of the fault was more likely to cause a larger fault slip than injecting fluid under the hanging wall of the fault under the same conditions. For underlying strata, the area that closes to the bottom of footwall of the fault was the more hazardous injection area.(3)When the fracturing region was in the footwall or hanging wall of the fault, injecting fluid into the region can induce felt earthquake if it was close to the fault. For the fracturing region in the footwall of the fault, the induced slip became larger when the region was close to the lower part of the fault. For the fracturing region in the hanging wall of the fault, on the contrary, the induced slip became larger when the region was close to the upper part of the fault. It revealed that the region that was close to the lower part of the fault in the footwall and the region that was close to the upper part of the fault in the hanging wall were more dangerous areas. The hazardous injection area was significantly influenced by the injection pressure. The area would be increased with the increase of the pressure.(4)The hazardous injection regions that threaten the fault’s stability were near the fault and concentrated on four areas: (a) above the top of the fault in underlying strata; (b) above the top of the hanging wall of the fault in underlying strata; (c) near to the fault planes in both footwall and hanging wall; and (d) the bottom of the footwall of the fault in underlying strata.

Note that this study was only focused on one fault in the site. The above hazardous injection area was only determined by the stability of the specific fault. However, it usually has many faults in the site. It has at least two effects on the hazardous injection area. Firstly, the interactions between the faults can directly influence their stability and slippage during fluid injection. Secondly, for an EGS project, the hazardous injection area should be a superimposed area of the hazardous regions of all the faults. Therefore, it is necessary to study the hazardous injection area determined by all the major faults. This involves simulating several faults in a model. It is beyond the scope of this paper, and further research is needed.

Data Availability

The data used to support the findings of this study are included within the article.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this article.

Acknowledgments

The authors are grateful for the financial support from the start-up fund for doctoral research of Suzhou University (no. 2019jb09 and no. 2019jb18) and Dual-ability teaching team project (Geological Engineering Dual-Energy Teaching Team, no. 2020XJSN06). The authors would like to express their sincere gratitude to Mahmoud from the University of Alberta for his help in proofreading this manuscript.

References

  1. W. G. P. Kumari and P. G. Ranjith, “Sustainable development of enhanced geothermal systems based on geotechnical research—a review,” Earth-Science Reviews, vol. 199, Article ID 102955, 2019. View at: Publisher Site | Google Scholar
  2. P. Bayer, L. Rybach, P. Blum, and R. Brauchler, “Review on life cycle environmental effects of geothermal power generation,” Renewable and Sustainable Energy Reviews, vol. 26, pp. 446–463, 2013. View at: Publisher Site | Google Scholar
  3. A. Aghahosseini and C. Breyer, “From hot rock to useful energy: a global estimate of enhanced geothermal systems potential,” Applied Energy, vol. 279, Article ID 115769, 2020. View at: Publisher Site | Google Scholar
  4. Y. Wang, Y. Liu, J. Duo, M. Liu, and M. Zeng, “Geothermal energy in China: status, challenges, and policy recommendations,” Utilities Policy, vol. 64, Article ID 101020, 2013. View at: Publisher Site | Google Scholar
  5. Z. Feng, Y. Zhao, A. Zhou, and N. Zhang, “Development program of hot dry rock geothermal resource in the Yangbajing Basin of China,” Renewable Energy, vol. 39, no. 1, pp. 490–495, 2012. View at: Publisher Site | Google Scholar
  6. R. M. Potter, E. S. Robinson, and M. C. Smith, Method of Extracting Heat from Dry Geothermal Reservoirs, US, Los Alamos, NM, USA, 1974.
  7. E. L. Majer, R. Baria, M. Stark et al., “Induced seismicity associated with enhanced geothermal systems,” Geothermics, vol. 36, no. 3, pp. 185–222, 2007. View at: Publisher Site | Google Scholar
  8. M. D. Zoback and H.-P. Harjes, “Injection-induced earthquakes and crustal stress at 9 km depth at the KTB deep drilling site, Germany,” Journal of Geophysical Research: Solid Earth, vol. 102, no. B8, pp. 18477–18491, 1997. View at: Publisher Site | Google Scholar
  9. R. P. Denlinger and D. R. H. O’Connell, “Evolution of faulting induced by deep fluid injection, Paradox Valley, Colorado,” Bulletin of the Seismological Society of America, vol. 110, no. 5, pp. 2308–2327, 2020. View at: Publisher Site | Google Scholar
  10. R. Baria, S. Michelet, J. Baumgärtner et al., “Creation and mapping of 5000 m deep HDR/HFR reservoir to produce electricity,” in Proceedings World Geothermal Congress 2005, vol. 8, p. 1627, Antalya, Turkey, April 2005. View at: Google Scholar
  11. K. H. Palgunadi, A. A. Gabriel, T. Ulrich, J.-N. López-Comino, and P. M. Mai, “Dynamic fault interaction during a fluid-injection-induced earthquake: the 2017 5.5 Pohang event,” Bulletin of the Seismological Society of America, vol. 110, no. 5, pp. 2328–2349, 2020. View at: Publisher Site | Google Scholar
  12. C. J. Bromley, C. F. Pearson, and D. M. Rigor, “PNOC-EDC, microearthquakes at the Puhagan geothermal field, Philippines—a case of induced seismicity,” Journal of Volcanology & Geothermal Research, vol. 31, no. 3-4, pp. 293–311, 1987. View at: Publisher Site | Google Scholar
  13. N. Deichmann and J. Ernst, “Earthquake focal mechanisms of the induced seismicity in 2006 and 2007 below Basel (Switzerland),” Swiss Journal of Geosciences, vol. 102, no. 3, pp. 457–466, 2009. View at: Publisher Site | Google Scholar
  14. S. Sherburn, C. Bromley, S. Bannister, S. Sewell, and S. Bourguignon, “New Zealand geothermal induced seismicity: an overview,” in Proceedings of the World Geothermal Congress, Melbourne, Australia, April 2015. View at: Google Scholar
  15. K. y. Nishigami, K. Tadokoro, S. Nagai, T. Mizuno, Y. Kano, and Y. Hiramatsu, “Triggering characteristics of induced earthquakes: brief review of water injection experiments in the Nojima fault and other regions,” Journal of Geography (Chigaku Zasshi), vol. 111, no. 2, pp. 268–276, 2002. View at: Publisher Site | Google Scholar
  16. F. Grigoli, S. Cesca, A. P. Rinaldi et al., “The November 2017 5.5 Pohang earthquake: a possible case of induced seismicity in South Korea,” Science, vol. 360, no. 6392, pp. 1003–1006, 2018. View at: Publisher Site | Google Scholar
  17. M. D. Zoback and S. M. Gorelick, “Earthquake triggering and large-scale geologic storage of carbon dioxide,” Proceedings of the National Academy of Sciences, vol. 109, no. 26, pp. 10164–10168, 2012. View at: Publisher Site | Google Scholar
  18. J. E. Streit and R. R. Hillis, “Estimating fault stability and sustainable fluid pressures for underground storage of CO2 in porous rock,” Energy, vol. 29, no. 9-10, pp. 1445–1456, 2004. View at: Publisher Site | Google Scholar
  19. J. Rutqvist, J. Birkholzer, F. Cappa, and C.-F. Tsang, “Estimating maximum sustainable injection pressure during geological sequestration of CO2 using coupled fluid flow and geomechanical fault-slip analysis,” Energy Conversion and Management, vol. 48, no. 6, pp. 1798–1807, 2007. View at: Publisher Site | Google Scholar
  20. X. C. Wei, Q. Li, X.-Y. Li, Y.-K. Sun, and X. H. Liu, “Uncertainty analysis of impact indicators for the integrity of combined caprock during CO2 geosequestration,” Engineering Geology, vol. 196, pp. 37–46, 2015. View at: Publisher Site | Google Scholar
  21. F. Cappa and J. Rutqvist, “Modeling of coupled deformation and permeability evolution during fault reactivation induced by deep underground injection of CO2,” International Journal of Greenhouse Gas Control, vol. 5, no. 2, pp. 336–346, 2011. View at: Publisher Site | Google Scholar
  22. J. Rutqvist, A. P. Rinaldi, F. Cappa, and G. J. Moridis, “Modeling of fault reactivation and induced seismicity during hydraulic fracturing of shale-gas reservoirs,” Journal of Petroleum Science and Engineering, vol. 107, pp. 31–44, 2013. View at: Publisher Site | Google Scholar
  23. T. H. W. Goebel and E. E. Brodsky, “The spatial footprint of injection wells in a global compilation of induced earthquake sequences,” Science, vol. 361, no. 6405, pp. 899–904, 2018. View at: Publisher Site | Google Scholar
  24. K. M. Keranen, M. Weingarten, G. A. Abers, B. A. Bekins, and S. Ge, “Sharp increase in central Oklahoma seismicity since 2008 induced by massive wastewater injection,” Science, vol. 345, no. 6195, pp. 448–451, 2014. View at: Publisher Site | Google Scholar
  25. H. Liu, H. Wang, H. Lei, L. Zhang, M. Bai, and L. Zhou, “Numerical modeling of thermal breakthrough induced by geothermal production in fractured granite,” Journal of Rock Mechanics and Geotechnical Engineering, vol. 12, no. 4, pp. 900–916, 2020. View at: Publisher Site | Google Scholar
  26. D. Zbinden, A. P. Rinaldi, T. Diehl, and S. Wiemer, “Hydromechanical modeling of fault reactivation in the St.Gallen deep geothermal project (Switzerland): poroelasticity or hydraulic connection?” Geophysical Research Letters, vol. 47, pp. 1–10, 2020. View at: Publisher Site | Google Scholar
  27. Hibbitt, Karlsson, and Sorensen, ABAQUS Theory Manual, Version 6.2, Hibbitt, Karlsson, Sorensen Inc., Providence, RI, USA, 2001.
  28. M. A. Biot, “General theory of three-dimensional consolidation,” Journal of Applied Physics, vol. 12, no. 2, pp. 155–164, 1941. View at: Publisher Site | Google Scholar
  29. J. R. Rice and M. P. Cleary, “Some basic stress diffusion solutions for fluid-saturated elastic porous media with compressible constituents,” Reviews of Geophysics, vol. 14, no. 2, pp. 227–241, 1976. View at: Publisher Site | Google Scholar
  30. K. Aki, “Scaling law of seismic spectrum,” Journal of Geophysical Research, vol. 72, no. 4, pp. 1217–1231, 1967. View at: Publisher Site | Google Scholar
  31. H. Kanamori and D. L. Anderson, “Theoretical basis of some empirical relations in seismology,” Bulletin of the Seismological Society of America, vol. 65, no. 5, pp. 1073–1095, 1975. View at: Google Scholar
  32. H. Kanamori, “The energy release in great earthquakes,” Journal of Geophysical Research, vol. 82, no. 20, 1977. View at: Publisher Site | Google Scholar
  33. T. C. Hanks and H. Kanamori, “A moment magnitude scale,” Journal of Geophysical Research, vol. 84, no. B5, pp. 2348–2350, 1979. View at: Publisher Site | Google Scholar

Copyright © 2021 Zhiyong Niu et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Related articles

No related content is available yet for this article.
 PDF Download Citation Citation
 Download other formatsMore
 Order printed copiesOrder
Views200
Downloads313
Citations

Related articles

No related content is available yet for this article.

Article of the Year Award: Outstanding research contributions of 2021, as selected by our Chief Editors. Read the winning articles.