*3.1. Laboratory Test*

As one of the preferred candidate locations for geological disposal of high-level radioactive waste, Alxa is located in the west of Inner Mongolia Autonomous Region, with the geographic coordinates of 97◦10 ~ 106◦53 E and 37◦24 ~ 42◦47 N. The NRG-1 deep borehole in Alxa revealed that a large number of porphyritic granites are distributed from the surface to 600 m underground. The granitic cores from different depths of the borehole were processed into the standard specimens with a diameter of 50 mm and a height of 100 mm. Before the uniaxial compression tests, the average density of the specimens was 2650.0 kg/m3. During the testing, the MTS 815 servo-controlled hydraulic testing machine at the Institute of Crustal Dynamics, China Earthquake Administration was employed to compress the granites with an axial loading speed of 0.06 m/min.

The failure characteristics of Alxa porphyritic granites after the tests and the stress-strain curves are shown in Figures 3 and 4, respectively. It can be seen that all the specimens were cut by numerous splitting cracks, which conforms to the typical brittle fracturing pattern. The macro mechanical parameters of the specimens were extracted from the stress-strain curves, and the average values of the UCS, the elastic modulus (E) and the Poisson's ratio (ν) were 159.4 MPa, 48.0 GPa and 0.23, respectively. The tensile test is not conducted in this research; however, the average UTS of this kind of granites is about 9.5 MPa according to our previous studies [60,61].

**Figure 3.** The failure characteristics of Alxa porphyritic granites after the uniaxial compressive tests.

#### *3.2. Material Parameters Calibration*

A rectangular numerical model with the same size as the test specimen was generated in PFC2D to emulate the uniaxial loading responses. The model was comprised of 10,412 circular particles with a radius distribution of 0.3–0.45 mm and was assigned the contact bonds of PBM. The effective modulus and stiffness of particles can inherit from the corresponding parameters of PBM. After the initial parameters were given, the simulations of uniaxial compression and direct tension were carried out to calculate the UCS, E, ν and UTS. Then the relevant parameters were continuously adjusted by the comparison between the experimental results and the simulated results. The final calibrated parameters are listed in Table 1, by which the simulated macro mechanical properties are in good agreemen<sup>t</sup> with the actual values, as shown in Table 2. The stress-strain curves and the spatial distribution of cracks are presented in Figure 4, which also proves the accuracy of parametric calibration.

**Figure 4.** The simulated results under the calibrated parametric condition: (**a**) stress-strain curve and microcrack evolution under uniaxial compression; (**b**) microcrack distribution after uniaxial compression; (**c**) stress-strain curve and microcrack evolution under uniaxial tension; and (**d**) microcrack distribution after uniaxial tension.

Most of the microcracks generated in the numerical granite under the compressive loading are the tension-type modes. Fewer microcracks formed in the middle of model under the tensile loading, and aggregated into a macro fracture.


**Table 1.** Calibrated Micro-Parameters in the Numerical Models of Alax Porphyritic Granite.

**Table 2.** Comparison of the Calibrated and Experimental Results of the Macro Mechanical Properties for Alax Porphyritic Granites.


1 The UTS of Alax porphyritic granite was acquired from Zhou et al. [60,61].

#### *3.3. Validation of Hydraulic Fracturing Model*

The fluid parameters for hydraulic fracturing were decided by the actual working conditions. After constructing the fluid networks and adding the fluid parameters in the PBM, it was still vital to verify whether the fracturing simulation could describe the realistic fracturing behaviors, including the breakdown pressure and the propagation of hydraulic fractures. A classical equation for the breakdown pressure *Pwf* when the fluid is injected into a circular elastic borehole was proposed by Haimson and Fairhurst [62] as following expression:

$$P\_{wf} = 3\sigma\_v - \sigma\_\text{h} - P\_0 + \sigma\_t \tag{12}$$

where σv is the confining pressure in vertical direction; σh is the confining pressure in horizontal direction; *P*0 is initial pore pressure; σ*t* is the UTS of the model.

In order to facilitate the comparison of simulated and theoretical breakdown pressures and avoid the stress interference of multiple injection wells, the granite model with a single injection well was established, as shown in Figure 5a. The square fracturing model with a width and height of 1.0 m consisted of 11,778 circular particles with a radius distribution of 4–6 mm. Furthermore, the vertical and horizontal confining pressures were servo-controlled by moving the walls around the model. An injection well with a radius of 50 mm was excavated in the center of the model to inject viscous fluid. It is noted that the particles surrounding the injection well have been replaced by smaller particles to smooth the well surface and prevent stress concentration (Figure 5b). No fluid flow domains were covered on the border of the fracturing model to give the impermeable boundary conditions. The hydraulic apertures *e*0 and *e*inf calculated by Equations (7) and (8) were 2.2 × 10−<sup>6</sup> and 2.2 × <sup>10</sup>−7, respectively, which were input into the model as basic hydraulic property parameters (Table 1). The fracturing fluid was the liquid water with a bulk modulus of 2.0 GPa and a viscosity of 1.0 × 10−3. The numerical model kept the vertical confining pressure σv and the horizontal confining pressure σh equal. Under the conditions of the initial pore pressure of 0 MPa and the fluid injection rate of 1.0 × 10−<sup>5</sup> m<sup>3</sup>/s, the hydraulic fracturing simulations with confining pressures from 5 MPa to 50 MPa (5-MPa intervals) were implemented.

**Figure 5.** Granite model with a single injection well for hydraulic fracturing: (**a**) the model size and applied boundary conditions; and (**b**) particle densification near the injection well.

The simulated values and theoretical values of breakdown fractures are summarized in Figure 6, and the error between them is less than 20%. The main reason for the deviation is that the analytical model assumes that the entire reservoir is impermeable, while the modified hydraulic fracturing model in PFC considers the leakage of fracturing fluid to the surrounding domains.

**Figure 6.** Comparison of breakdown pressures between the simulated values and the analytical solutions.

Figure 7 shows the representative borehole pressure histories, crack spatial distribution and corresponding fluid pressure field in hydraulic fracturing under different confining pressures. It is found that the borehole pressure increases rapidly at the beginning of injection, and as the borehole pressure approaches the breakdown pressure, the borehole pressure increases slowly and non-linearly, which reflects that the leakage of fracturing faster is greater when the borehole pressure is high. The leakage of fluid from the injection well and fractures to the surrounding fluid networks can also be seen from the figures of the fluid pressure distribution.

**Figure 7.** The representative borehole pressure histories, crack spatial distribution and corresponding fluid pressure field in the single-well fracturing under different in-situ stress conditions: (**a**) σv = σh = 20 MPa; (**b**) σv = σh = 30 MPa; and (**c**) σv = σh = 40 MPa.

Before the borehole pressure reaches the breakdown pressure, a few cracks initiated around the wellbore. After increasing to the breakdown pressure, the borehole pressure dropped sharply to the fracture propagation pressure and remained unchanged. As a result, the number of microcracks grew dramatically along certain directions. Hubbert and Willis [63] proposed that the initiation and propagation of deep hydraulic fracturing fractures follow the direction of maximum principal stress. From Figure 7, we can see that the macro hydraulic fractures propagated along the di fferent directions under various confining pressures. The di fficulty in predicting the propagation direction stems from the fact that the micro defects on the surface of the injection well have more significant impact on the initiation and propagation of hydraulic fractures under hydrostatic pressure conditions.

Another noteworthy issue is that almost all the hydraulic fractures acquired from the fluid-mechanical coupling algorithm are made up of tensile cracks, which seems to be distinct from lots of shear-type seismic events observed in hydraulic fracturing tests for granite [64,65]. Al-Busaidi et al. [35] have demonstrated that the shear failures recorded in the laboratory tests are caused by the slippage of the preexisting cracks in the specimens, and the mechanism of injection-induced fracture is predominantly tensile failure. In brief, the simulation results of granite fracturing using the modified fluid-mechanical coupling algorithm are accurate.

#### **4. Hydraulic Fracturing Process in the Specimens with Two Injection Wells**
