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.
(a)
(b)
(c)
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 [17–21].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 [31–33]:
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.
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.
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.
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.
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.