Study on the Initiation of Interface Crack in Rock Joints

The interfacial fracture of rock joints is an important although easily ignored issue in jointed rock engineering. To conduct this study, an interface crack model of rock joints was proposed. By analyzing the ratio of stress intensity factor to fracture toughness, the fracture mode of the interface crack was studied. Based on the Mohr-Coulomb criterion, an interface fracture criterion considering T-stress was established. To verify the proposed fracture criterion, laboratory and numerical tests were conducted. Finally, the effect of relative critical size α, internal friction angle φ and cohesion c on the initiation of an interface crack was comprehensively discussed. It is concluded that the proposed fracture criterion can predit the initiation of the interface cracks properly. With an increase in cohesion c, mode II fracture toughness KIIC also clearly increases. When the absolute value of KI is small, the effect of α is much larger than that of φ. In addition, with an increase in the absolute value of the mode I stress intensity factor, the φ of the joint plays a more important role in the initiation of the interface crack.


Introduction
Natural rock masses are composed of a rock matrix and discontinuities. The existence of discontinuities, such as rock joints, has a significant influence on the strength and deformation properties of rock masses. Generally, it can be found that the prime failure modes of rock masses are the shear slips along the rock joints [1,2]. There are also often flaws inside the rock joints, which have a serious influence on the shear mechanical properties of rock joints [3]. To better understand the failure process of rock joints caused by cracks, the mesomechanics method (i.e., fracture mechanics) should be adopted, as it has been used in rock mechanics [4][5][6]. Therefore, the study on the interfacial fracture of rock joints has an important theoretical and practical significance.
Until now, there have been three methods to study the failure of rock joints, namely laboratory experiments, numerical modeling, and theoretical research. Experimental research tends to focus on the influence of the filling material, surface roughness or matedness of the contact surfaces and the coupling of external environments on the failure of rock joints by laboratory tests. For example, She and Sun [7] conducted compressive shear tests on cement-filled joints and investigated the peak shear strength of these joints. Tang and Wong [8] conducted direct shear tests on artificial rock joints with different contact states and analyzed the variation of the peak shear strength. Zhao et al. [9] carried out several groups of direct shear tests on wetting-treated rock joint specimens and studied their shear failure behavior. However, large numbers of tests are time-consuming and it is a difficult to create some specific conditions in laboratory tests, such as setting the interface crack in the rock joint. Moreover, this kind of study can only reveal some common macro-mechanical properties of the rock joints.
Meanwhile, different from experimental technique, numerical simulation is a suitable method to investigate the mesoscopic failure of rock joints. For example, Guo and Qi [10] simulated rock specimens with unfilled joints using the finite difference method and investigated the progressive failure mode of jointed rocks with various undulate joints from a mesoscopic view. Tian et al. [11] proposed an interface model to simulate the shear failure of cemented concrete-rock joints using the finite element method and analyzed the effect of gradual bond failure processes of cohesive interface elements on the postpeak performance of joints. Park and Song [12] simulated the direct shear test to study the influence of meso-parameters on the shear failure of rock joints based on the bonded particle model. Asadi et al. [13] carried out the shear test simulations on rough rock joints using a particle flow code and studied the effects of the meso-properties as well as the geometrical features on the fracture shear behavior. Compared with the continuum models, the discontinuous models, such as the discrete element method, have some unique characteristics, especially the advantage of simulating crack propagation [14]. Although the failure of rock joints have been investigated from a mesoscopic perspective, the effect of an interface crack in the rock joints, which is the main factor in joint failure, has not been researched.
Theoretical research is a key aspect to the study of the failure of rock joints. So far, the effect of surface roughness of joint surfaces on the failure of rock joints [15][16][17] or establishing strength model which consider various physical parameters (e.g., dilation angle and morphology of rock joints) [18][19][20] is hot topic of the theoretical research. However, generally, in these studies, a semiempirical method has been used, which cannot explain the failure of the joint interface from the mesomechanics. In fact, the failure of rock joints always starts at its positions of weakness (i.e., flaws), which can be studied to good effect by fracture mechanics. To study the initiation mechanism of cracks, there are three main theories: maximum tangential stress criterion [21], strain energy density criterion [22], and energy release rate criterion [23]. Among these, the maximum tangential stress criterion (MTS) is the most widely used in rock material for its simple form and precise physical meaning [24][25][26]. In order to predict crack initiation angle more accurately, T-stress near crack tip was introduced to the MTS criterion [27]. The initiation of a crack which is in single medium can be well solved by the theories above-mentioned. When a crack is in the rock joint interface, which is termed as interface crack, the theories are difficult to apply to such case due to the discontinuity of strength parameters along the joint. For interface cracks, the related theoretical research is mainly focused on the interface crack of composite materials (bimaterials) with strong interface, for which the stress field at the crack tip exhibits oscillatory singularity induced by mismatching of the dissimilar materials on both sides of the interface [28,29]. For example, in Cartesian coordinates, Deng [30] first derived analytical solution of the crack-tip stress field for the dissimilar materials. Zhou and Li [31] obtained full elastic fields of the interface crack tip and corresponding energy release rate for orthotropic bi-material. Banks-Sills [32] presented failure criteria for a crack between two dissimilar materials based on energy release rate. Mega et al. [33] proposed several two-and three-dimensional mixed-mode interface failure criteria for predicting delamination failure in multidirectional, laminate composites. Zhang et al. [34] proposed a modified critical energy density criterion for peridynamic interface bond failure analysis. In fact, for the cracks in the rock joint interface, both sides of the interface are the same materials and there is a fracture process zone at the crack tip of rock-like materials [35], for which classical interfacial fracture mechanics is also unsuitable for determining the fracture of the interface crack of rock joints. Moreover, the strength parameters of the rock joint interface are much lower than that of rock matrix, which is termed as weak interface. The weak interface will affect the occurrence of mode I fracture [36] which is always occurred as the cracks is in rock matrix [37]. Therefore, establishing a fracture criterion for the interface crack in the rock joint is very urgent.
In this paper, the interface crack model of rock joints is proposed for the first time. Moreover, by using fracture mechanics, the fracture mode of the interface crack is determined, and a new fracture criterion for interface cracks of rock joints is proposed. Lastly, the new fracture criterion is verified by the experimental and numerical tests using the discrete element method.

Interface Crack Model
Generally, there are many defects (weaknesses) in the rock joints. Hence, a model of a rock joint containing a pre-existing crack is shown as in Figure 1. In this model, for joints with a small thickness, the joint can be simplified to a mathematical interface where the interface crack exists.
termined, and a new fracture criterion for interface cracks of rock joints is proposed. Lastly, the new fracture criterion is verified by the experimental and numerical tests using the discrete element method.

Interface Crack Model
Generally, there are many defects (weaknesses) in the rock joints. Hence, a model of a rock joint containing a pre-existing crack is shown as in Figure 1. In this model, for joints with a small thickness, the joint can be simplified to a mathematical interface where the interface crack exists. For the interface crack, Deng [30] derived general expressions of the crack-tip stress field for dissimilar materials, which is given as: where KΙ and KΙΙ are the stress intensity factors for crack mode Ι and II, respectively. For the interface crack, Deng [30] derived general expressions of the crack-tip stress field for dissimilar materials, which is given as: Re η I n k n r iε n √ 2πσ I n (θ) + Im η II n k n r iε n √ 2πσ where ε is the oscillation index, r and θ are polar coordinates when the ordinate origin is assumed to be at the crack tip, and k n , η I n , η II n ,σ I n (θ),σ II n (θ) andσ III n (θ) are expressions in the derivation. Due to the complexity of the expressions for these parameters, the specific expressions are not listed here. The details can be found in [30].
For bimaterials, the interface crack-tip stress field exhibits oscillatory features. However, for rock joint, both sides of the joint interface are the same materials. Therefore, the bimaterial becomes homogeneous, which is illustrated by ε = 0. Crack-tip stress fields with oscillatory features are simplified to those of traditional fracture mechanics problems. Therefore, the stress field at the interface crack tip of the joint is as follows [38]: where K I and K II are the stress intensity factors for crack mode I and II, respectively. The interface crack model is shown in Figure 1. In Figure 1, the stress on a crack surface whose dip angle is β can be described as: where k is the ratio of minimum principal stress to maximum principal stress, and σ y and τ xy are the normal and shear stresses of the joint plane, respectively.
Then, according to the basic theory of fracture mechanics [38], the stress intensity factors of mode I and II cracks are given as where a is the half-length of the crack. In Equations (3) and (4), σ y is a positive value, which represents the tensile stress. Under compression conditions, the value of the stress intensity factor K I should be negative, which illustrates that K I is able to restrain the initiation of cracks [39].

Determination of Crack Initiation Angle
According to the maximum circumferential stress criterion, the crack will initiate along the direction of maximum tensile stress. However, for interface crack in rock joints, the crack exists in a joint plane with a strength much lower than that of rock block [40], which means that the fracture toughness of the joint plane is much weaker than that of rock block. Therefore, the rock near the interface crack tip cannot be regarded as homogeneous material. In this study, the fracture toughness of mode I and mode II fractures in the direction of joint extension (θ = 0) is defined as K IC_joint and K IIC_joint , respectively, and as K IC_matrix and K IIC_matrix for the other directions (θ = 0), respectively.
When the interface crack is under biaxial compression, the stress field at the crack tip in the polar coordinate can be described as [38]: The fracture modes (here only mode I and mode II fractures are considered) depend on σ θ and τ rθ , regardless of which kind of stress is applied [41]. Accordingly, whether mode I or II fractures occur at the initiation angle θ IC or θ IIC depends on the critical value of the circumferential tensile stress σ θ or the shear stress τ rθ . Here, different joint dip angles β (15 • , 30 • , 45 • , 60 • ) and ratios of minimum principal stress to maximum principal stress k (0, 0.1, 0.3, 0.5, 0.7) were chosen in order to analyze the initiation condition of the interface crack.
In order to better analyze the relationship between stress intensity factors and fracture toughness, the stress intensity factors in different direction angles θ near the crack tip can be defined as [41]: Equation (5) was substituted into Equation (6), resulting in: Whether a mode I or mode II fracture occurs depends on whether K I max or K II max takes precedence to reach its critical value, K IC or K IIC [41]. Therefore, the fracture mode can be determined by: where K I(II)max is the maximum of the mode I or II stress intensity factor around the crack tip. The K IIC /K IC ratio is about 2.6 [41] and the fracture toughness in other directions (θ = 0) is about 4 times that in the direction of joint extension (θ = 0) [42]. In order to determine the initiation angle and the representative conclusions, the ratio of fracture toughness in other directions and that in the direction of joint extension ( ) was selected to be 4. The variation in in various conditions was calculated as in Figure 2.
It is noted that to better facilitate the calculation and illustrate the results, the non-dimension of K I(II) (θ) (i.e., K I(II) (θ)/σ √ πa) was used to analyze the stress intensity factor distribution and mode I fracture toughness of rock matrix K IC_matrix was set set as 1.
As can be seen from Figure 2a-d, when k is set to 0.1 for the mode I stress intensity factor, the maximum value increases first and then decreases with an increase in the crack dip angle β (the angle between the crack and horizontal direction). The angle (θ) corresponding to the maximum value gradually decreases with increase in the crack dip angle β, which is in agreement with the experiments by Rao et al. [41]. Although the mode I stress intensity factor reaches its maximum when the crack dip angle β is about 60 • , the mode II fracture still occurs due to a lower fracture toughness in the joint plane. From Figure 2c,e-h, it is concluded that the mode I stress intensity factor maximum decreases as confining pressure increase and the mode I stress intensity factor is negative in all directions when k is equal to 0.7, which illustrates that the existence of confining pressure inhibits the growth of the mode I stress intensity factor. Considering the above analysis, the conditions such that the crack dip angle β is about 60 • and the confining pressure is 0 are most suitable for a mode I fracture to occur. However, as shown in Figure 2e, the maximum of K I (θ)/K IC and |K II (θ)|/K IIC are 0.3779 and 0.6662, respectively. According to Equation (8), a mode II fracture will still occur and the interface crack will initiate along the joint plane.
Therefore, for interface cracks in rock joints, a mode II fracture often occurs, which is different from traditional rock fracture mechanics, where a mode I fracture more commonly occurs and the initiation angle is about 70.5 • [43].
Equation (5) was substituted into Equation (6), resulting in: Whether a mode Ι or mode ΙΙ fracture occurs depends on whether KΙ max or KΙΙ max takes precedence to reach its critical value, KΙC or KΙΙC [41]. Therefore, the fracture mode can be determined by: where KΙ(ΙΙ)max is the maximum of the mode Ι or ΙΙ stress intensity factor around the crack tip. The ratio is about 2.6 [41] and the fracture toughness in other directions (θ ≠ 0) is about 4 times that in the direction of joint extension (θ = 0) [42]. In order to determine the initiation angle and the representative conclusions, the ratio of fracture toughness in other directions and that in the direction of joint extension was selected to be 4. The variation in in various conditions was calculated as in Figure 2. It is noted that to better facilitate the calculation and illustrate the results, the ) was used to analyze the stress intensity factor distribution and mode Ι fracture toughness of rock matrix _ C m a tr ix K  was set set as 1.

Crack Initiation Criterion
Since a mode II fracture generally occur for interface cracks in rock joints, the crack initiation stress cannot be determined by the maximum circumferential stress criterion. Therefore, in this study, a new fracture criterion is proposed, which is used to predict the initiation stress of the interface crack under compression-shear conditions.
When traditional fracture mechanics theory is used to study the crack initiation mechanism of rock, the Williams expansion is frequently used, which can be written as [44]: where r and θ are polar coordinates when the ordinate origin is assumed to be at the crack tip. In the first term, A 1 can be regarded as the stress intensity factors of K I or K II when the constant terms are absorbed into the trigonometric function f n ij terms. In the second term, A 2 represents the stress that is directly applied to the crack line on the normal plane. The third term is the high order term.
In traditional fracture mechanics theory, only the singular stress term r −1/2 in Williams expansions is considered, which incurs a difference between the theoretical and experimental results [45]. Moreover, for rock material, there is a fracture process zone at the crack tip [46]. For the stress field of the fracture process zone at the crack tip, the proportion of singular stress term decreases and the T-stress term increases [47]. Thus, the effect of T-stress should be considered.
Based on Equation (5), when T-stress is considered, the following is derived: The T-stress can be described as [45]: According to the above analysis, for the interface crack, a "shear fracture" (i.e., mode II fracture) along the joint will generally occur. Moreover, the initiation stress of the interface crack is determined by the stress states of the fracture process zone [48], which illustrates that the initiation stress is a function of the stress in the fracture process zone. Therefore, the Mohr-Coulomb criterion, which is used to describe the shear failure of materials, can be adopted to describe the initiation of the interface crack. The formula can be written as where c is the cohesion, ϕ is the internal friction angle, τ is the shear stress and σ is the normal stress of joints. Equation (10) is substituted into Equation (12) and θ set as equal to 0 (the initiation angle of the interface crack in the rock joint is 0 • ), deriving the following: Both sides of Equation (13) are multiplied by √ 2πr, and the following equation can be obtained: For Equation (14), r = r c and, for simplification, the critical radius r c is presented in the dimensionless form of α = √ 2r c /a. Moreover, the stress field of the fracture process zone is considered. Then, Equation (11) is substituted into Equation (14), deriving: In Equation (15), the term c √ πa is consistent with the definition of stress intensity factor in form, and as such is the mode II fracture toughness produced by the cohesion c. When the average mode II fracture toughness of the fracture process zone is considered, the following is obtained: Substituting α = √ 2r c /a into the right side of Equation (15) and combining with Equation (16) obtains: Equation (17) is substituted into Equation (15), and the fracture criterion of the interface crack can be obtained:

Verification by Experimental Study
Tan and Xu [49] conducted double-edge notched single-edge compression tests (DNSCTs) on the layered concrete, as shown in Figure 3. In their study, the influence of the pouring interval of two layers on the bond characteristics was investigated. For specimens with the same pouring interval, the cohesion c can be determined by direct shear tests and the mode II fracture toughness can be determined by the DNSCT. For the specimens with different pouring intervals, the relationship of K IIC and c is shown in Figure 4. As shown in Figure 4, there is a linear relationship between K IIC and c, where the fitting coefficient R 2 is 0.9980. the pouring interval of two layers on the bond characteristics was investigated. For specimens with the same pouring interval, the cohesion c can be determined by direct shear tests and the mode ΙΙ fracture toughness can be determined by the DNSCT. For the specimens with different pouring intervals, the relationship of KΙΙC and c is shown in Figure 4. As shown in Figure 4, there is a linear relationship between KΙΙC and c, where the fitting coefficient R 2 is 0.9980.  In this test, the interface crack is in a pure shearing state. For this condition, the proposed fracture criterion can be generated as the following expression: From the above equation (Equation (19)), the theoretical relationship between K IIC and c can be obtained. Therefore, the theoretical result is in good agreement with the real test result. Hence, the criterion can be verified by real tests.
Generally, it is difficult to control the test conditions of mode II fractures. Therefore, it is very challenging to obtain the mode II fracture toughness from real tests. However, by using the theoretical study in this paper, the quantitative relation between mode II fracture toughness and the shear strength parameter c is determined, avoiding the difficulty of conducting the real test.  In this test, the interface crack is in a pure shearing state. For this condition, the proposed fracture criterion can be generated as the following expression: From the above equation (Equation (19)), the theoretical relationship between KΙΙC and c can be obtained. Therefore, the theoretical result is in good agreement with the real test result. Hence, the criterion can be verified by real tests.
Generally, it is difficult to control the test conditions of mode ΙΙ fractures. Therefore, it is very challenging to obtain the mode ΙΙ fracture toughness from real tests. However,

Verification by the Numerical Study
Due to the difficulty of conducting the test in the laboratory using the interface crack specimens, to comprehensively and thoroughly verify the proposed fracture criterion of the interface crack, a numerical test is used here.
The numerical simulation material is siltstone [50]. The numerical model is generated by particle flow code (PFC2D), and the contact model between particles is the bondedparticle model (BPM). The size of the two-dimensional model is such that the width is 50 mm and the height is 100 mm. There is a consecutive joint in the model. For a comprehensive study, four models, for which the inclination angles of the joint are 35.8 • , 45 • , 54.2 • , and 63.4 • from the horizontal, have been constructed, as shown in Figure 5. The interface crack is in the middle of the joint, which has a length of 12 mm. The particle diameters of these models satisfy a uniform distribution, with a range between 0.17 mm and 0.3 mm. For the different numerical models, there are different numbers of particles, contacts, etc. For example, for the numerical model with a joint angle of 45 • , there are 28,296 particles, 76,451 contacts, 75,946 parallel bonds and 289 smooth-joint contacts. The particle density of siltstone is 2400 kg/m 3 , and the average porosity is 7%, which is the same as the physical parameters of real specimens [50]. by using the theoretical study in this paper, the quantitative relation between mode ΙΙ fracture toughness and the shear strength parameter c is determined, avoiding the difficulty of conducting the real test.

Verification by the Numerical Study
Due to the difficulty of conducting the test in the laboratory using the interface crack specimens, to comprehensively and thoroughly verify the proposed fracture criterion of the interface crack, a numerical test is used here.
The numerical simulation material is siltstone [50]. The numerical model is generated by particle flow code (PFC2D), and the contact model between particles is the bonded-particle model (BPM). The size of the two-dimensional model is such that the width is 50 mm and the height is 100 mm. There is a consecutive joint in the model. For a comprehensive study, four models, for which the inclination angles of the joint are 35.8°, 45°, 54.2°, and 63.4° from the horizontal, have been constructed, as shown in Figure 5. The interface crack is in the middle of the joint, which has a length of 12 mm. The particle diameters of these models satisfy a uniform distribution, with a range between 0.17 mm and 0.3 mm. For the different numerical models, there are different numbers of particles, contacts, etc. For example, for the numerical model with a joint angle of 45°, there are 28,296 particles, 76,451 contacts, 75,946 parallel bonds and 289 smooth-joint contacts. The particle density of siltstone is 2400 kg/m 3 , and the average porosity is 7%, which is the same as the physical parameters of real specimens [50]. In this study, the smooth-joint contact model (SJM) is used to simulate the behavior of the joints. As shown in Figure 6, the SJM suggests that the particles will overlap with each other along the smooth interface. Hence, this contact model can simulate the behavior of the joints more realistically, which avoids the unrealistic mechanicals properties of joints induced by the bumpiness or roughness of the interface surfaces [51]. Therefore, in In this study, the smooth-joint contact model (SJM) is used to simulate the behavior of the joints. As shown in Figure 6, the SJM suggests that the particles will overlap with each other along the smooth interface. Hence, this contact model can simulate the behavior of the joints more realistically, which avoids the unrealistic mechanicals properties of joints induced by the bumpiness or roughness of the interface surfaces [51]. Therefore, in this way, the initiation of an interface crack under uniaxial compression can be simulated accurately. In this study, the smooth-joint contact model (SJM) is used to simulate the behavior of the joints. As shown in Figure 6, the SJM suggests that the particles will overlap with each other along the smooth interface. Hence, this contact model can simulate the behavior of the joints more realistically, which avoids the unrealistic mechanicals properties of joints induced by the bumpiness or roughness of the interface surfaces [51]. Therefore, in this way, the initiation of an interface crack under uniaxial compression can be simulated accurately.

Calibration of Micro-Parameters
For the particle flow code, micro-parameters are used and their calibration is crucial. In order to simulate the behavior of single-jointed rock, the numerical model of the singlejointed rock specimen is calibrated based on the stress-stain curve obtained by a real experiment [50], as shown in Figure 7. From Figure 7, it is clear that the mechanical behavior of numerical model is in agreement with that of the real specimen. For the singlejointedrock specimen, the joint inclination angle is 68 • . The calibrated micro-parameters of the numerical model are listed in Table 1.

Calibration of Micro-Parameters
For the particle flow code, micro-parameters are used and their calibration is crucial. In order to simulate the behavior of single-jointed rock, the numerical model of the single-jointed rock specimen is calibrated based on the stress-stain curve obtained by a real experiment [50], as shown in Figure 7. From Figure 7, it is clear that the mechanical behavior of numerical model is in agreement with that of the real specimen. For the single-jointedrock specimen, the joint inclination angle is 68°. The calibrated microparameters of the numerical model are listed in Table 1.

Results of Numerical Study
To reproduce the real static response of the specimen, the displacement loading from the upper wall is applied, which simulates the lab uniaxial compression. By fixing the bottom wall and moving the upper wall down, causing an increase in upper loading, the initiation of the interface crack can be clearly simulated. According to [52], the quasistatic loading rate can be set to 0.1 × 10 −3 mm/min. By using numerical studies, the initiation stress of the interface crack specimens can be obtained, as shown in Figure 8. It should be noted that the initiation stress is determined when the first fractures occur around the crack tip. In Figure 8 Table 2.

Results of Numerical Study
To reproduce the real static response of the specimen, the displacement loading from the upper wall is applied, which simulates the lab uniaxial compression. By fixing the bottom wall and moving the upper wall down, causing an increase in upper loading, the initiation of the interface crack can be clearly simulated. According to [52], the quasistatic loading rate can be set to 0.1 × 10 −3 mm/min. By using numerical studies, the initiation stress of the interface crack specimens can be obtained, as shown in Figure 8. It should be noted that the initiation stress is determined when the first fractures occur around the crack tip. In Figure 8 Table 2.    Liu [45] pointed out that for metal and polymethyl methacrylate materials, the test results are ideal when α is equal to 0.1. For rock material, the critical radius r c is larger than that of metal and polymethyl methacrylate material and related research [47] recommended that relative critical size α of less than 1 is appropriate. Therefore, the proposed fracture criterion with a varied relative critical size α, which range from 0.1 to 0.9, are compared with the numerical results. The mode I and II stress intensity factors and the interface crack fracture criterion with a varied relative critical size α are summarized in Figure 9. the interface crack fracture criterion with a varied relative critical size α are summarized in Figure 9. As shown from Figure 9, based on the new fracture criterion, the theoretical curves (μ = tan20° = 0.36) are linear. As the relative critical size α increases, the slope and intercept increase monotonically, which illustrates that the mode ΙΙ fracture toughness and initiation stress of the interface crack increase. When α = 0.6, the numerical test results are in better agreement with the theoretical results. Moreover, according to Equation (14), the critical size rc can be determined as 1.08 mm.

Discussion
In Equation (18), there are three controlling parameters, which are the cohesion c, internal friction angle φ and relative critical size α. The influence of those parameters on the initiation of interface crack is analyzed as follows. Moreover, some comparisons with other studies were also discussed.

Effect of Cohesion c
To study the effect of cohesion c, the parameter c is taken as 7.7 MPa, 17.7 MPa and 27.7 MPa, and the other parameters, such as internal friction angle φ and relative critical size α, were fixed as 20° and 0.6, respectively. According to these conditions, the corresponding results are summarized in Figure 10. It can be seen that with the increase in cohesion c, the mode ΙΙ fracture toughness KΙΙC is 0.87 MPa‧m 0.5 , 1.54 MPa‧m 0.5 and 2.31 MPa‧m 0.5 , respectively, which illustrates that the mode ΙΙ fracture toughness KΙΙC increases with the increase in cohesion. The initiation stress of the interface crack increases with the increase in cohesion. It should be noted that the slope of the relationship curves between mode Ι and ΙΙ stress intensity factors remains unchanged, which indicates that there are no differential effects of varying cohesion c on the initiation of the interface crack. As shown from Figure 9, based on the new fracture criterion, the theoretical curves (µ = tan20 • = 0.36) are linear. As the relative critical size α increases, the slope and intercept increase monotonically, which illustrates that the mode II fracture toughness and initiation stress of the interface crack increase. When α = 0.6, the numerical test results are in better agreement with the theoretical results. Moreover, according to Equation (14), the critical size r c can be determined as 1.08 mm.

Discussion
In Equation (18), there are three controlling parameters, which are the cohesion c, internal friction angle ϕ and relative critical size α. The influence of those parameters on the initiation of interface crack is analyzed as follows. Moreover, some comparisons with other studies were also discussed.

Effect of Cohesion c
To study the effect of cohesion c, the parameter c is taken as 7.7 MPa, 17.7 MPa and 27.7 MPa, and the other parameters, such as internal friction angle ϕ and relative critical size α, were fixed as 20 • and 0.6, respectively. According to these conditions, the corresponding results are summarized in Figure 10. It can be seen that with the increase in cohesion c, the mode II fracture toughness K IIC is 0.87 MPa·m 0.5 , 1.54 MPa·m 0.5 and 2.31 MPa·m 0.5 , respectively, which illustrates that the mode II fracture toughness K IIC increases with the increase in cohesion. The initiation stress of the interface crack increases with the increase in cohesion. It should be noted that the slope of the relationship curves between mode I and II stress intensity factors remains unchanged, which indicates that there are no differential effects of varying cohesion c on the initiation of the interface crack. In essence, the fracture toughness is the characteristic constant of materials, which represents the resistance to the propagation of a crack [38]. When there is no normal stress on the joint plane (pure shearing), the resistance to the propagation of crack (i.e., mode ΙΙ fracture toughness) is provided by the cohesion, which is generated by the mutual attraction between the molecules. Therefore, when the cohesion c increases, the mode ΙΙ fracture toughness increases linearly.

Effect of Internal Friction Angle φ
To study the effect of the friction angle φ, parameter φ was taken as 20°, 25°, 30° and 35°, and the other parameters, such as cohesion c and relative critical size α, were fixed as 7.7 MPa and 0.6, respectively. According to these conditions, the corresponding results are summarized in Figure 11. It can be seen from Figure 11 that the relationship curves between mode Ι and ΙΙ stress intensity factors monotonically increase, and with the increase in internal friction angle φ, the slope of the curves enlarges. However, the intercept term of the curves remains unchanged. This indicates that the internal friction angle affects the interface crack initiation stress and has no influence on the mode ΙΙ fracture toughness KΙΙC. In essence, the fracture toughness is the characteristic constant of materials, which represents the resistance to the propagation of a crack [38]. When there is no normal stress on the joint plane (pure shearing), the resistance to the propagation of crack (i.e., mode II fracture toughness) is provided by the cohesion, which is generated by the mutual attraction between the molecules. Therefore, when the cohesion c increases, the mode II fracture toughness increases linearly.

Effect of Internal Friction Angle ϕ
To study the effect of the friction angle ϕ, parameter ϕ was taken as 20 • , 25 • , 30 • and 35 • , and the other parameters, such as cohesion c and relative critical size α, were fixed as 7.7 MPa and 0.6, respectively. According to these conditions, the corresponding results are summarized in Figure 11. It can be seen from Figure 11 that the relationship curves between mode I and II stress intensity factors monotonically increase, and with the increase in internal friction angle ϕ, the slope of the curves enlarges. However, the intercept term of the curves remains unchanged. This indicates that the internal friction angle affects the interface crack initiation stress and has no influence on the mode II fracture toughness K IIC .
The internal friction angle is another shear strength parameter, which represents frictional properties and the shear strength of the contact area [53]. When there is no normal stress on the joint plane, the mode II fracture toughness is provided by the cohesion. Therefore, the mode II fracture toughness K IIC will not be affected by the internal friction angle. When there is normal stress on the joint plane, the mode I stress intensity factor K I is negative in the fracture process zone. With an increase in the internal friction angle, more driving force is needed to overcome the frictional strength in the fracture process zone. Therefore, the slope of the relationship curves between mode I and II stress intensity factors increases with the increase in the internal friction angle. The internal friction angle is another shear strength parameter, which represents frictional properties and the shear strength of the contact area [53]. When there is no normal stress on the joint plane, the mode ΙΙ fracture toughness is provided by the cohesion. Therefore, the mode ΙΙ fracture toughness KΙΙC will not be affected by the internal friction angle. When there is normal stress on the joint plane, the mode Ι stress intensity factor KΙ is negative in the fracture process zone. With an increase in the internal friction angle, more driving force is needed to overcome the frictional strength in the fracture process zone. Therefore, the slope of the relationship curves between mode Ι and ΙΙ stress intensity factors increases with the increase in the internal friction angle.

Effect of Relative Critical Size α
To study the effect of relative critical size α, parameter α was taken as 0.4, 0.5, 0.6 and 0.7, and the other parameters, such as cohesion c and internal friction angle φ, were fixed as 7.7 MPa and 20°, respectively. According to these conditions, the corresponding results are summarized in Figure 12. It can be seen from Figure 12 that the relationship curves between mode Ι and ΙΙ stress intensity factors monotonically increase, and with the increase in relative critical size α, both the slope and intercept of the curves enlarge. Comparing Figures 11 and 12, it can be concluded that when the absolute value of KΙ is small, the effect of α on the initiation of the interface crack is much greater than that of φ. Then, with the increase in K  , the internal friction angle φ of the joint plays a more important role in the initiation of the interface crack.

Effect of Relative Critical Size α
To study the effect of relative critical size α, parameter α was taken as 0.4, 0.5, 0.6 and 0.7, and the other parameters, such as cohesion c and internal friction angle ϕ, were fixed as 7.7 MPa and 20 • , respectively. According to these conditions, the corresponding results are summarized in Figure 12. It can be seen from Figure 12 that the relationship curves between mode I and II stress intensity factors monotonically increase, and with the increase in relative critical size α, both the slope and intercept of the curves enlarge. Comparing Figures 11 and 12, it can be concluded that when the absolute value of K I is small, the effect of α on the initiation of the interface crack is much greater than that of ϕ. Then, with the increase in |K I |, the internal friction angle ϕ of the joint plays a more important role in the initiation of the interface crack. The relative critical size is also regarded as the characteristic constant of materials in the field of rock fracture mechanics [45]. When the relative critical size increases, more driving force is needed to overcome the cohesion and frictional strength. Therefore, both the slope and intercept of the relationship curves between mode Ι and ΙΙ stress intensity factors enlarge. The relative critical size is also regarded as the characteristic constant of materials in the field of rock fracture mechanics [45]. When the relative critical size increases, more driving force is needed to overcome the cohesion and frictional strength. Therefore, both the slope and intercept of the relationship curves between mode I and II stress intensity factors enlarge.

Comparison with Other Studies
To further study the characteristics of this fracture criterion, it was compared with several traditional crack initiation criterion, which are listed in Table 3. The classic MTS criterion [21] and the MTS criterion considering T-stress [45] were chosen for this comparison. As shown in Table 3, although the application scope and the fracture mode of the proposed criterion are different from that of the two other traditional crack initiation criterions, the initiation conditions are essentially the same, which depend on whether the absolute value of the stress intensity factor reaches its fracture toughness. Although the classical MTS criterion is concise, the MTS criterion considering T-stress is more accurate in predicting the initiation angle. For the proposed criterion, the number of key parameters is six, which is relatively concise and easy to understand. Moreover, based on the fixed initiation angle, the initiation stress can be predicted well.

Conclusions
In order to solve the interfacial fracture of rock joints, the interface crack model of joints was proposed for the first time. Moreover, based on the proposed interface crack model, a new fracture criterion considering T-stress was proposed. Then, the proposed fracture criterion was verified by the laboratory tests and numerical study. Finally, the effect of relative critical size α, internal friction angle ϕ and cohesion c on the initiation of the interface crack was discussed, with the following conclusions drawn: (1) A Mode II fracture generally occurs in the interface crack of rock joints, which is different from the fracture mode that often occurs in the rock matrix. (2) The theoretical results calculated by the proposed fracture criterion are in good agreement with the experimental and numerical results. Hence, the proposed fracture criterion could be verified by the test results. (3) The effect of T-stress was considered in the fracture criterion. It is shown that the cohesion c and internal friction angle ϕ increase with the increase in relative critical size α. (4) The relative critical size α, internal friction angle ϕ and cohesion c all affect the initiation of the interface crack. With an increase in cohesion c, mode II fracture toughness K IIC clearly increases. When the value of K I is small, the effect of α on the initiation of the interface crack is much greater than that of ϕ. Then, with an increase in |K I |, the ϕ of the joint plays a more important role on the initiation of the interface crack.
However, the proposed fracture criterion is primary study of the interface fracture in rock joints. Only the case of pure mode II fracture was verified by experiments, and the other conditions were verified by numerical simulations based on experimental data. Therefore, designing more experiments and simulation conditions for the fracture criterion remain challenges for our future research. Moreover, the proposed fracture criterion combined with cohesive zone model can be written into distinct element code modeling, which can be used to analyze the interface fracture of jointed rock in geotechnical engineering.

Data Availability Statement:
The raw/processed data required to reproduce these findings cannot be shared at this time as the data also forms part of an ongoing study.