Numerical investigation on cracking behavior of granite with intersecting two-flaws: A flat-joint modeling method

In practical rock engineering, it is of considerable significance for the study of mechanical behavior of rock containing flaws. The uniaxial compression failure mechanism of granite containing intersecting two-flaws is investigated through the flat-joint modeling approach. The flat-joint model is introduced briefly, and numerical models of intact granite specimens are established. After that, a series of micro-parameters are confirmed by comparison with the laboratory tests of intact granite specimen. Then uniaxial compression simulations are conducted on granites with intersecting two-flaws. The numerically simulated results show good consistency with laboratory tests, including the strength properties, deformation properties, and failure modes. The strength and deformation characteristics of the granites are dependent on the relative position length. Furthermore, the relationship between the crack evolution process and the stress-strain curves for the specimens is studied concerning the relative position length. Finally, the crack propagation mechanism and cracking type of granite specimens are revealed by analyzing the micro-force field


INTRODUCTION
Rock, as a typical heterogeneous material, usually contains a large number of flaws (weak surfaces, faults, fissures, and voids, etc.) from the microscopic to macroscopic scale (Lee et al., 2017;Li et al., 2020;Wang et al., 2020;Wu et al., 2020;Yang et al., 2016;Zhou et al., 2014).The flaws can control the mechanical behavior of rock, and have a significant impact on the safety and stability of rock engineering (Liu et al., 2020).Lots of rock engineering practices have shown that rock engineering failures are usually caused by the initiation and propagation of cracks in the rock mass.Hence, it is of considerable significance to investigate the crack evolution behavior of flawed rock in various rock engineerings.
To investigate the failure mechanism of flawed rock comprehensively, many scholars have conducted experimental research on rock and rock-like materials with different fracture geometries (Chiu et al., 2016;Hadjigeorgiou et al., 2009;Qin et al., 2019;Wang and Tian, 2018;Wong and Einstein, 2009;Zhou et al., 2014).Shen et al. (1995) conducted several laboratory tests on gypsum specimens containing pre-existing fractures.They reported that different fracture geometries produced significantly different stress fields in the rock bridge area, and hence, resulted in various failure modes.Mughieda and Alzo'ubi (2004) investigated the influence of bridge inclination on the failure modes of rock-like materials containing flaws.They observed that wing cracks started to propagate before failure because of the tensile stress concentration.Cheng et al. (2015) studied the coalesce behaviors of rock with multiple flaws in laboratory tests and observed that the failure process was similar to that of a flawed pair.Chiu et al. (2016) proposed a rough-joint model to simulate the deformability and strength of rock joint and studied the characteristics of shear displacements, normal closure, and shear dilations of the joints.Zhang et al. (2019) studied the compressive failure process of rock-like materials containing X-type flaws aligned in the loading direction.They revealed that the crack evolution laws and fracture mechanisms of rock mass with X-type flaws.Wang et al. (2020) investigated the failure behavior of jointed rock mass with various joint density and joint distribution.Chai et al. (2020) conducted a series of static and dynamic compression tests on jointed rock samples filled with different filling materials, and revealed the strength characteristics and failure modes of the filled rock joints.
Along with experimental studies, the numerical analysis methods have been extensively used to investigate the failure behaviors of the rock mass (Chen et al., 1998;Feng et al., 2006;Ning et al., 2011;Potyondy and Cundall, 2004;Wu and Wong, 2013).Among these numerical methods, the particle flow code (PFC), which can directly reproduce cracks through bonding fractures between particles without using the complex theory of fracture mechanics, is widely developed to reveal the failure mechanism for the rock (Zhang and Wong, 2012).By adopting PFC2D, Fakhimi et al. (2002) numerically studied the sandstone failure with a circular flaw under biaxial compression.They indicated that the numerical model of PFC2D was able to reproduce the damaged areas observed in laboratory tests.Lee and Jeon (2011) analyzed the fracture evolution of granite containing a combination of horizontal and inclined flaws by PFC2D.Their results revealed that there were reasonable agreements between the results of the laboratory tests and the PFC2D modeling.Khazaei et al. (2015) investigated the damage process of intact rocks based on the acoustic emission energies through the PFC3D simulation and experimental test.Jin et al. (2017) investigated the effects of flaws on the crack initiation, failure pattern, and energy mechanism of the rock-like materials under uniaxial compression using PFC2D.
The literatures aforementioned pay more attention to the mechanical behavior of rock material containing parallel flaws and non-intersecting flaws (Figure 1(a) and (b)).However, there are few works on intersecting flaws (Lajtai et al., 1994;Esterhuizen et al., 2011;Yang et al., 2014).As is well known all, stress redistribution during the construction of rock engineering can generate numerous new flaws, which contains a multitude of intersecting flaws, as shown in Figure 1(c).Therefore, in this research, to extend the previous experimental study, numerical investigation for rock with intersecting flaws is conducted using particle flow code in two-dimension (PFC2D).Firstly, the flat-joint contact model in PFC2D is briefly introduced.After that, the granite micromechanical parameters are calibrated based on the previous laboratory tests (Yin et al., 2016).And then, the PFC2D simulations are conducted for intersecting two-flaws-contained specimens under uniaxial compression.Following that has been simulated by PFC2D, the influence of relative position length on the strength properties, deformation behaviors and ultimate failure modes are evaluated thoroughly.Also, the cracking mechanism of granite specimens is analyzed.Finally, the evolution laws of the micro-force field within the granites containing intersecting two-flaws are revealed.The PFC2D results can provide insightful information for crack propagation of rock containing intersecting flaws in rock engineering.method Peng Liang et al.

Flat-joint model
Over the past decades, the most widely used in the modeling of rock by PFC is the standard bonded-particle model (BPM) (Cho et al., 2007;Potyondy and Cundall, 2004).However, an intrinsic limitation of the standard BPM has been an inability to reproduce the ratio of the unconfined compressive strength (UCS) to the tensile strength (TS), the appropriate internal friction angle, and the nonlinear strength envelope.To overcome the limitation, the flat-joint contact model (FJM) is proposed (Potyondy, 2012).Figure 2 illustrates a schematic diagram of the FJM.The FJM is composed of particles connected through flat-joint contacts.In FJM, the interface can be discretized into several equal-length segments, with each segment initially bonded.As the bonded segments break, the interface behavior changes from the completely bonded states to the wholly de-bonded and frictional states.However, the fully de-bonded interface is not removed, which will continue to resist relative rotation while it cannot work in the standard BPM.

Numerical specimen and simulation procedure
Granite specimens are used in this research.The size of the numerical specimen (70 mm×140 mm) by PFC2D is the same as that of specimens used in laboratory tests (Yin et al., 2016).Each numerical specimen consists of 39803 particles and 109546 contacts.In PFC2D, the determination of micro-parameters for specified material is complicated.For a better understanding of the influence of micro-parameters on the numerical results of PFC2D, an in-depth study of the relationship between the micro-parameters and the macroscopic properties of rock materials is conducted.In the present study, the micro-parameters are confirmed by the trial and error methods, which are used by many researchers (Manouchehrian et al., 2014;Yang et al., 2016;Zhang and Wong, 2012).After multiple calibrations, the micro-parameters for the granite specimens from the numerical tests in this study are proposed, which are listed in Table 1.method Peng Liang et al.Solids and Structures, 2020, 17(4), e280 4/14 The comparison of laboratory tests and numerical results is listed in Table 2.As shown in Table 2, the numerical results of the intact granite specimen achieve a good agreement with laboratory tests.Moreover, it can be seen that the failure pattern of granite simulated is conformed to the failure pattern of the specimen obtained in laboratory tests (Figure 3 (a) and (b)).

Mechanical properties
Experimental results

Numerical results
Uniaxial compression strength, UCS (MPa) 120.66 120.39 Young's modulus, E (GPa) 10.49 10.52 Axial strain (10 -2 ) 1.45 1.21 To generate intersecting flaws, some ball particles are deleted after creating the intact numerical specimen.The flaw geometry is shown in Figure 3(c).The primary flaw ① length 2m is 30 mm, and the secondary flaw ② length n is 15 mm.The inclination angle β is fixed to 45°.And p is the relative position length (RPL, the distance between the intersection of two flaws and the top of the flaw ①), which varies from 0 mm to 30 mm at an interval of 7.5 mm.The opening thickness of pre-existing flaws in numerical models is about 1.0 mm.For the pre-flawed granite specimen, the loads are continuously applied until the failure point is reached.The failure point is defined as the position of the stress-strain curves at which the axial stress of post-peak falls to 60% of the peak strength.

NUMERICAL RESULTS ANALYSIS
Based on the laboratory tests of granite with intersecting two-flaws of relative position length (Yin et al., 2016), we will make a systematic analysis of the granite specimens and make a quantitative comparison between the numerical results and the laboratory tests.

Analysis of strength and deformation characteristics
Figure 4 shows the comparison results of axial stress-strain curves.By comparing with the experimental results, it can be observed that a trend from the numerical simulations has a reasonable agreement with the results from the laboratory tests before peak strength.Furthermore, the post-peak part of stress-strain curves obtained numerically is more ductile than that obtained from laboratory tests, which results mainly from the different fracturing mechanisms between numerical simulation and laboratory tests.It is seen that there are two noticeable stress drops after the peak-strength for p = 30 mm, which is different from the results of laboratory tests.Notably, for all the stress-strain curves obtained by PFC2D, there is stress when axial strain is zero.It is mainly because of the generation process of the specimens, during which locked-in stress must be inputted to make each particle contact strictly with other particles (Cao et al., 2018;Yang et al., 2016).Moreover, there are other slight discrepancies for the shapes of stress-strain curves between the experimental and numerical results in Figure 4.These discrepancies are attributed to the fact that there are some natural micro-cracks in the experimental rock specimens.These natural micro-cracks can affect the mechanical response of the rocks.Besides, the PFC2D results may be influenced by the number of particles.Moreover, the results show differences in values from the laboratory tests and PFC2D simulation.In laboratory tests, the UCS takes an obvious nonlinear relation as the RPL increases.The average UCS firstly decreases from 77.93 to 50.52 MPa when p increases from 0 to 15 mm, whereas UCS increases severely from 50.52 to 68.05 MPa when p increases from 15 to 22.5 mm.Besides, the average UCS decreases from 68.05 to 34.36 MPa when p rises from 22.5 to 30 mm, which is not observed in numerical simulation.In the numerical simulation, the average UCS decreases when p increases from 0 to 30 mm.The maximum average UCS is 72.66 MPa when p is 0 mm, whereas the minimum average UCS is 47.86 MPa at the RPL of p is 30 mm.From Figure 5, it is seen that when p equals to 0, 22.5 and 30 mm, Young's modulus of PFC2D simulation is similarly equal to that from laboratory tests.Furthermore, when p equals to 7.5 and 15 mm, we can conclude that Young's modulus of PFC2D simulation is higher, by 1.65 GPa and 2.24 GPa, than that obtained by laboratory tests, respectively.These differences from the above analysis may be due to the following reasons.The PFC2D model is two-dimensional, which cannot accurately reflect the three-dimensional physical phenomena of rock by laboratory tests.Besides, it is possibly the effect of the heterogeneity of the materials or inherent randomness in the PFC2D models.

Analysis of cracking characteristics
The cracks can cause discontinuities and nonlinearities in the rock.Compared with laboratory research, the PFC2D simulation is more convenient and effective for obtaining mechanical information of rock cracks.However, the overall pattern of the exact micro-crack propagation path in various locations could not be obtained in the laboratory tests.Therefore, the cracking process of the specimens is investigated for the entire deformation process from a microscopic point of view by using PFC2D.Figure 6 illustrates the comparison results of failure modes.In Figure 6, the red and blue short lines stand for the tensile and shear cracks, respectively.Following the compared results, we can see that the results obtained from PFC2D modeling show similar patterns with the ultimate failure modes of the experimental specimens.In other words, compared with laboratory data, the PFC2D model can capture the development of micro-fractures during a compression test.The specimen failure is mainly caused by tensile cracks.Moreover, the cracks in granite specimens emanate from the flaw tips and propagate along the direction of axial stress in most cases.Likewise, we have observed some densely overlapped tensile and shear cracks in a direction between the flaw plane and the axial stress.Also, numerous crack clusters are observed nearby the flaws of the numerical specimens.From the PFC2D simulation, the stress-strain curves and the number of cracks are calculated in the loading process of granite specimens (as shown in Figure 7).In these Figures, the 'Crack-num' represents the total number of failure bonds, 'Tensile crack' represents the number of failure normal bonds, and 'Shear crack' represents the number of failure shear bonds.The number of tensile cracks in the models is larger than that of the shear cracks, and the difference becomes more pronounced as the stress increases.To make a systematic evaluation of the cracking behavior in granite specimens, the stress-strain curves are separated into three stages as follows.Stage I from 0 to point i is the crack initiation stage; stage II from point i to peak strength point v is the stable crack growth stage; stage III from point v to point vii is the post-peak stage.The axial stress at the point i corresponds to the stress of crack initiation.The stress point i, as an important stress threshold in the progressive failure process of rock, is widely studied by many scholars.Martin (1993) suggested using the calculated crack volumetric strain and acoustic emission method to identify crack initiation.However, it is difficult to confirm the stress of crack initiation due to the early appearance and propagation of micro-cracks in the PFC2D models.Especially for flaw-contained rock, the first crack usually initiates near the flaws.Therefore, in the present study, as the crack number of the PFC2D models reaches about 10% of the total crack number at point v, the corresponding stress is identified as the stress of crack initiation (Yang et al., 2016).In this paper, typical specimens are selected for crack evolution analysis.Figure 8 and Figure 9 illustrate the evolution process of cracks for the specimen of p = 0 mm and p = 30 mm, respectively.And some features of the crack propagation can be summarized.It should be noted that the failure points correspond to those shown in Figure 7 For p = 0 mm, as the axial stress applied on the specimen reaches the stress point i, anti-tensile crack 1 a and tensile crack 1 b develop from the primary flaw ① and secondary flaw ②, respectively (Figure 8(a)).With the increase of the axial load, the granite will undergo the stable deformation, i.e., the stable crack growth stage.During stage II, the micro-cracks increase stably, which link to form macro-cracks gradually.However, the granite specimen can still maintain its stability.Furthermore, the number of all cracks is rare when the axial stress levels are lower.Whereas, as the stress levels are higher, the number of shear and tensile cracks increases steadily.With the increase of the stress, crack 1 a and crack 1 b propagate along the load direction (Figure 8(b)).And after that, the shear crack 2 originates from the lower tip of flaw ①, as shown in Figure 8(c), then propagates along the direction of flaw ①.Afterwards, as the stress (σ1) is loaded to 71.04 MPa, the anti-tensile secondary crack 3 initiates from the lower tip of flaw ①, and the initiation direction is almost vertical to flaw ① (Figure 8(d)).The tensile crack 1 b propagates when stress point reaches point v, but only propagates a little.Also, the major anti-tensile secondary crack 4 and some crack clusters are observed around the upper tip of flaw ① (Figure 8(e)).The propagating speed of crack 3 and crack 4 is higher than that of crack 1 a and 1 b .In correspondence Figure 8(d) and Figure 8(e), although the axial load increases slightly, the evolution of cracks can still be seen as the deformation of the granite specimen increases.
After reaching peak strength, the stress drops dramatically.The specimen enters the post-peak stage.At this stage, the cracks grow rapidly, and the crack growth rate is higher than that of the first two stages.As the stress decreases to the point vi, the crack 1 a , 3 and 4 propagate along the stress direction to the bottom edge of the specimen.For the primary flaw ①, two tensile cracks 5 a and 5 b develop simultaneously from the lower and upper tips of the flaw, respectively.Moreover, the starting direction is approximately vertical to the horizontal direction (Figure 8  For p = 30 mm in stage I, the tensile crack 1 b and shear crack 1 a develop from the primary flaw ① and secondary flaw ②, respectively (Figure 9(a)).As the specimen is loaded from point ii to point iii before the peak strength, the tensile crack 1 b appears to propagate in the stress direction.Also, the new anti-tensile crack 2 originates from the lower tip of flaw ① and propagates in the lower direction of the specimen (Figure 9 After the point v, the crack 4 b propagates to the left bottom edge of the specimen.Also, the crack 4 a , crack 1 b , and crack 3 propagate along the load direction.Moreover, several far field cracks 5 are observed near the upper tip of flaw ① (Figure 9(f)).As the stress decreases to point vii, the anti-tensile crack 6 develops from the upper tip of flaw ① and propagates to the downward direction.Also, the crack 1 b and 2 propagate to the boundary of the specimen along the stress direction (Figure 9(g)).
Based on the above analysis of the failure process for granite with intersecting two-flaws, some characteristics of the crack evolution are summarized.The crack propagation modes could change in the loading process, and the specimens show different failure modes due to various RPLs.Besides, the tensile crack number is markedly larger than the number of shear crack, and differences are more obvious as the deformation increases.For the specimens with various RPLs, the cracks of tension and anti-tension mostly appear from the upper tip of flaw ①, and propagate along Latin American Journal of Solids and Structures, 2020, 17(4), e280 11/14 the direction of major principal stress.However, for the flaw ②, cracks of tension usually emanate from the upper tip.
Besides, the anti-tensile cracks are often in the opposite direction of the tensile cracks or secondary tensile cracks.The far field cracks and crack clusters always originate after the appearance of the first crack for all the simulated specimens.

Analysis of micro-mechanical behavior
To reveal the crack evolution mechanism comprehensively, the analysis of micro-force field is critical.However, it is nearly impossible to be obtained in laboratory tests.The distribution of contact force and flat-joint model interface force during crack propagation of specimens are tracked and recorded.Here, we take the specimen for p = 0 mm as an example, and Figure 10 shows the cracking process and the micro-force field of the granite with intersecting two-flaws (p = 0 mm).The point i-vii in Figure 10 correspond to those shown in Figure 7(a).In Figure 10, the blue lines represent the interface force of flat-jointed contact, while the red lines represent the cracks in the blue images.The black images of Figure 10 are the distribution of flat-jointed force chain, in which black color represents compressive force, while the red color represents tensile force.Notably, in the black images, the line orientation and thickness stand for the force orientation and magnitude, respectively.Generally, the magnitude of the force near the flaws is greater than that at the regions far from the flaws, which is similar to the results of the specimen with one flaw or two un-parallel flaws (Zhang and Wong 2012).As the specimen is loaded to point i (σ1 = 22.01 MPa, ε1 = 2.316×10 -3 ), the maximum interface force is 61 kN.With the initiation of anti-tensile crack 1 a and tensile crack 1 b , the concentration of tensile stress at the tip regions of flaw ① and flaw ② is released, which is shown by more red segments at the tips of crack 1 a and crack 1 b than that at the other positions in the flat-jointed force chain field.Interface force mainly concentrates near the two tips of flaw ①.As the stress increases to point ii (σ1 = 45.20 MPa, ε1 = 4.481×10 -3 ), the crack number has a small rise, crack 1 a and crack 1 b propagate continually.Red colors in the flat-jointed force chain field appear at the tips of crack 1 a and crack 1 b , which means that the new crack tips will become a new area of tensile stress concentration.Interface force concentration is also observed near the two tips of flaw ①.At point ii, the maximum interface force is 119 kN.At point iii (σ1 = 57.09MPa, ε1 = 5.582×10 -3 ), the maximum interface force is 155 kN.From point iii, the red color crowds nearby the lower tips of flaw ①, which indicates the initiation of shear crack 2. The interface force concentrates at tips of flaws more obviously.At point iv (σ1 = 71.04MPa, ε1 = 6.873×10 -3 ), the maximum interface force increases from 155 to 194 kN.Nearby the region of flaw tips, the red color increases greatly, which indicates that the new crack tips will serve as the concentration area of tensile stress.At the peak strength point v (σ1 = 72.66MPa, ε1 = 7.129×10 -3 ), the maximum interface force decreases greatly from 194 to 179 kN.At point v, more tensile force concentrations are observed than at point iv due to the propagation of crack 3 and crack 4.Moreover, there are some micro-crack clusters and recognizable macro-cracks.At point vi (σ1 =63.20 MPa, ε1 =7.398×10 -3 ), the maximum interface force is 151 kN.It should be noted that the flat-jointed interface force varies distinctively, and more blank regions appear compared than point i-v.The concentration of tensile stress near the tips of flaw ① is released, resulting from the appearance of crack 5 a and crack 5 b , and the corresponding red line segments also reduce near the tips.At point vii (σ1 = 44.24MPa, ε1 = 7.672×10 -3 ), the maximum interface force increases slightly from 151 to 158 kN.At this time, large amounts of blank areas appear, which indicates the release of concentration of tensile stress on account of the propagation of macro-crack clusters 6.
Based on the above analysis, it can be seen that the concentration region of tensile stress will transfer to the new crack tips as the cracking propagation.Besides, once the tensile or anti-tensile cracks emanates, the concentration of tensile stress will be released obviously.The occurrence of macro-cracks can lead to a reduction in the red and black parts near the flaw tips, indicating that stress concentration has been released.

Figure 1 :
Figure 1: Rock mass containing parallel, non-intersecting and intersecting flaws (arrow and yellow circles refer to the flaws).

Figure 2 :
Figure 2: The flat-joint contact model between particles

Figure 3 :
Figure 3: Results of failure modes and specimen geometry with intersecting two-flaws: (a) Experimental failure pattern; (b) Numerical failure pattern; (c) Flaw geometry.

Figure 4 : 14 Figure 5 :
Figure 4: Comparison results of stress-strain curvesFigure5illustrates the results of strength and deformation properties.The uniaxial compressive strength (UCS) and Young's modulus of granite specimens with intersecting two-flaws are dependent on the RPL.In laboratory tests, Young's modulus exhibits the same trend as UCS with the increase of the RPL, i.e., firstly decreases, and then increases, finally decreases.By contrast, both the UCS and Young's modulus of the specimens in PFC2D models decreases with increasing RPL, and have a minimum value at p = 30 mm.

Figure 6 :
Figure 6: Comparison results of failure modes
(f)).As the method Peng Liang et al.Latin American Journal of Solids and Structures, 2020, 17(4), e280 10/14 stress decreases, cracks from the two tips of flaw ① and the upper tip of flaw ② propagate until the boundary of the granite model.Also, numerous macro-crack clusters 6 are observed at the bottom right of the model (Figure 8(g)).
(b) and (c)).As the deformation increases, the stress reaches the point iv, and several crack clusters 3 are observed near the upper tip of flaw ②.Besides, the crack 1 b and crack 2 propagate steadily along the stress direction (Figure 9(d)).Afterwards, as the deformation increases to peak strength (point v), the secondary tensile crack 4 a and anti-tensile crack 4 b start to propagate in the stress direction (Figure 9(e)).

Table 1
Input micro-parameters for the granite specimen