Abstract

To investigate the crack behaviour of rock or rock-like material in uniaxial loading, a series of numerical simulations were conducted on gypsum specimens containing a single flaw with different inclination angle (0°–90°) and length (10 mm–30 mm). Based on the numerical simulations results, the effect of flaw length on peak strength, the CI stress, and the CD stress were investigated with different inclination angles. The results show that the peak strength decreased initially and then increased with increasing of the flaw angle. Meanwhile, the peak strength decreased gradually when the length of the preexisting flaw increased. When the inclination angle was 30°, 45°, and 60°, the reduction degree of peak strength increased with increasing of the flaw length. The CI stress and CD stress not only depend on the inclination angle but also depend on flaw length. Four types of crack were observed in numerical simulations. The present research facilitates increased understanding of crack behaviour of rock under different conditions.

1. Introduction

The cracking behaviour of rock and rock-like specimens containing artificial flaws under compressive loading has been experimentally investigated for decades. For example, gypsum [15], glass [6], polymethyl methacrylate (PMMA) [79], diastone [13, 7, 1013], sandstone [8, 14], granite [7], and limestone [9] were tested. Similar observations have been made in their experiments using different terms—two types of cracks are observed in compression: tensile or wing cracks and shear or secondary cracks (Figure 1). Wing or tensile cracks were usually found to be the first cracks to appear from the preexisting flaws. Secondary cracks are shear cracks and always initiate from the tips of the flaws after the first cracks. They are initially stable, but they may become unstable during the loading process near crack coalescence or near specimen failure. Some studies revealed that crack inclination angles have significant effect on the crack propagation and coalescence process. Bobet [15] analysed the initiation of secondary cracks in gypsum with two flaws (angle from 30 to 60) under compression. It observed all flaw inclination angles, wing cracks, and secondary cracks initiate at the same stress level. Park and Bobet [16] investigated crack coalescence in gypsum specimens with open and closed flaws under compression and one of their results showed that the mean wing and secondary crack initiation stresses increase with flaw inclination angle and coalescence stress slightly increases as the flaw angle increases. Yang and Jing [14] analysed strength failure and crack coalescence behaviour of brittle sandstone samples containing a single fissure (angle from 0° to 75°) under uniaxial compression. They found that the fissure angle has a key effect on the strength and deformation behaviour of sandstone samples under uniaxial compression, and the mechanical parameters (uniaxial compressive strength Young’s modulus and peak axial strain) of flawed sandstone firstly decrease and then increase with increasing fissure angle. Zhou et al. [17, 18] conducted a series of uniaxial compression experiments and numerical simulation on rock-like materials (formed from a mixture of sand, plaster, limestone, and water at mass ratio of 126 : 9 : 9 : 16) containing multiple flaws (angle from 30° to 60°) and concluded that the crack initiation modes depend on the flaw angle and the nonoverlapping length.

A number of experimental and numerical investigations have been conducted for fractured rock under compression since Nemat-Nasser and Horii [20] investigated the fracturing mechanisms of flaw plates (model material) and obtained that flaw length is one of the parameters controlling the failure mode of sample. Wong et al. [21], Tang and Kou [22], and Tang et al. [23] performed uniaxial and biaxial compression on numerical model samples containing a number of large, preexisting flaws and a row of suitably oriented smaller flaws by RFPA2D software, which shows that the crack growth is stable and stops at some finite crack length with the increase of lateral pressure, while a lateral tensile stress with even a small value will result in an unstable crack growth after a critical crack length is attained. Yang and Jing [14] performed a series of uniaxial compressive tests on brittle sandstone which contained a single fissure (angle from 30° to 60° and length from 0 to 25 mm). Zhang et al. (2017) proposed a modified phase-field model for mixed mode fracture. Zhou et al. (2016) used the nonordinary state-based peridynamics to investigate crack propagation and coalescence in precracked rock-like Brazilian disks. Zhou et al. (2017) used the improved peridynamic method to simulate the failure of rock-like material under the compressive loading. Wang et al. (2016, 2017) also analysed crack propagation and coalescence in rocks under uniaxial compression using peridynamic method.

It was found that the uniaxial compressive strength, Young’s modulus, and peak axial strain of sandstone samples with preexisting single fissure are all lower than those of intact sandstone sample, which is closely related to the fissure length and fissure angle.

The above experimental and numerical results for precracked rock or rock-like material have distinctly shown that inclination angle and length of cracks in rock or rock-like material have an important influence on the crack behaviours. However, it should be noted that some previous investigations focused mainly on a fixed length of precrack with varying inclination angle for rock and the length of crack was not taken in account. Some previous literatures analysed the peak strength and deformation modulus with different crack lengths but lower inclination angle (30°–60°), due to hard fabrication of crack. Less experimental studies were performed on rock with different crack lengths and different inclination angles (0–90°) to investigate the initiation and propagation process of precrack rocks. Therefore, in this paper, to better understand the effects of inclination angle and length of crack on the peak strength, crack initial stress, crack secondary stress, and failure behaviour of precrack rock, a series of numerical simulations with distinct element code (UDEC) on the uniaxial compression tests were conducted for precrack gypsum specimen with varying inclination angle (0–90°) and length (10–30 mm) of cracks. Moreover, the Voronoi block in UDEC is adopted to better describe crack initiation, propagation, and coalescence process.

2. The 2D Distinct Element Method

2.1. Introduction about Distinct Element Method

UDEC is a two-dimensional numerical program based on the distinct element method for discontinuity modelling. It has been widely used in the literature to explore damage of rock [2429], as well as Brazilian test analysis [30, 31]. UDEC especially joint model used in UDEC captures several features that are representative of the physical response of joints. For joint model in UDEC, all blocks have connection with contacts (corner-to-corner contacts, edge-to-corner contacts, and edge-to-edge contact). However, edge-to-edge contact is important, because it corresponds to the case of a rock joint closed along its entire length (Figure 2). If the maximum tensile stress exceeds the tensile strength of the contact or the maximum shear stress exceeds the shear strength of the contact, contact will break. And a microtensile crack or a microshear crack, respectively, will form. The contact behaviour is described by where and are the normal and shear stresses, and are the normal and shear displacements, and are the tensile and residual tensile strengths, is the shear strength, and are the cohesive and the residual cohesive strengths, and are the friction and the residual friction angles, and is the incremental contact shear displacement.

2.2. Procedure of the Simulation

In this study, the UDEC model is the gypsum specimen which was selected from previous literature [1]. The reason is that the experiment tests about gypsum specimen containing single crack under different conditions were conducted by Wong [1], and the physical experiments were also repeated numerically by PFC in previous literatures [32, 33]. The UDEC model was set up with elastic but indivisible Voronoi blocks. Model size was chosen to replicate those of the diametrical cross section of specimens used for laboratory testing (i.e., 50 mm width and 100 mm height) (see Figure 2). Average edge length of Voronoi blocks was selected as 0.002 m which yielded about 1320 Voronoi blocks in total within the model. Every Voronoi block is divided into a plurality of triangular-shaped block. The maximum edge length of the triangular zones is 0.00133 m. Figure 2 shows a typical UDEC model and how individual Voronoi blocks are numerically interconnected.

For a distinct element model such as the UDEC that synthesizes macroscale material behaviour from the interactions of microscale components, the input properties of the microscopic constituents are usually not known. Such intrinsic characteristics imply that it is nearly impossible to calibrate a set of macro output properties from the assembly such that they are exactly equal to that of the physical tests. Only a set of comparable properties can be obtained from the calibrated model. The parameters of gypsum in experimental tests and numerical simulations are shown in Table 1. Figure 3 shows the stress-strain curves of the intact gypsum specimen of UDEC and experiment test under uniaxial compression.

The second stage of simulation involved a series of numerical simulation tests with single crack under different conditions: the crack length of gypsum specimen is fixed at 10 mm, 15 mm, 20 mm, 25 mm, and 30 mm, respectively, while the crack inclination angle β measured horizontally is varied from 0° to 90° (see Figure 2).

In the present numerical study, specimens are loaded under a uniaxial vertical compression in a velocity controlled manner. A sufficiently low loading rate of 0.005 m/s is applied to ensure that the specimen remains in a quasi-static equilibrium throughout the test. It takes approximately 1000000 steps to load one specimen to complete failure.

3. Numerical Simulation Results and Analysis

The results of the simulation work performed as described in the previous section were analysed in different ways to meaningfully describe crack behaviour of gypsum in different conditions under uniaxial compression.

3.1. Stress-Strain Characteristics

Figure 4 shows the axial stress-axial strain curves and failed contact number-axial strain curves of UDEC model containing a single flaw with different inclination angles and lengths under uniaxial loading, respectively. As shown in Figure 4, all axial stress-axial strain behaviour can be divided into four stages. To show this part clearly, an axial stress-axial strain curve for model containing a single flaw with inclination angle 30° and length 10 mm under uniaxial compression is depicted, and four pictures of UDEC model at different stages were included in Figure 5.

Before the initiation of macrocrack from the preexisting crack, only a few microcracks (white colour) initiated at the axial stress of 9 MPa. In this stage, the curve of failed contact number is almost unchanged. As the axial stress reached 17.6 MPa, macrotensile crack (white colour) initiated from both ends of preexisting crack. The number of failed contacts increased steadily with increasing axial strain in this stage. When the axial stress increased to about 25.5 MPa, the secondary crack initiated from left top of preexisting crack. The number of failed contacts increased rapidly when the axial strain is changed a little. The peak value of stress-strain curve is 26 MPa. After the peak point, a large macrocracks’ coalescence has happened and the model failed.

The failed contact numbers under different conditions are shown in Figure 6. It can be observed from Figure 6 that the failed contact number at first increased and then decreased with increasing of inclination angle β under the same length of preexisting crack. When the inclination angle β was 0° and 90°, respectively, the failed contact number is almost unchanged with length of preexisting crack from 10 mm to 30 mm. However, when the inclination angle β was 15° and 60°, respectively, the failed contact number increased gradually when length of preexisting crack increased. But when the inclination angle β was 30° and 45°, respectively, the failed contact number increased gradually when length of preexisting crack decreased. And this decrease trend is most obvious when the inclination angle β was 45°, as shown in Figure 6.

The peak strengths of model containing single flaw are shown in Figure 7. It is shown from Figure 7 that the peak strength at first decreased and then increased with increasing of the flaw angle. The behaviours are akin to those observed by Zhou et al. [17], Yang and Jing [14], Xu et al. [34], and Wasantha et al. [35]. Meanwhile, the peak strength decreased gradually when length of preexisting crack increased. To show the extent of the reduction in peak strength of the model in relation to an intact specimen, the reduction degree (Yang et al., 2012) [17] is defined aswhere is the peak strength of intact model (the value is 33.5 MPa). is the peak strength of model containing single flaw.

The reduction degrees for peak strengths of different models are listed in Table 2 and Figure 8.

For the model with flaw length of 10 mm and 15 mm, increased from 2.6% to 24.2% when the inclination angle β was increased from 0° to 90°. For the model with flaw length of 20 mm, 25 mm, and 30 mm, increased from 2.9% to 50.5% when the inclination angle β was increased from 0° to 90°. For the model with flaw length  mm, increased when the inclination angle β was 30°, 45°, and 60°, respectively. For the model with flaw length  mm, almost unchanged when the inclination angle β was 0° to 90°, respectively. When the inclination angle β was 90°, 75°, and 0°, the flaw length has no influence on . When the inclination angle β was 30°, 45°, and 60°, increased with increasing of the flaw length .

3.2. Crack Initiation Stress and Crack Damage Stress

In most of the literatures studied crack behaviour, crack initiation (CI) stress, and crack damage (CD) stress were the most frequently discussed. Nicksiar and Martin (2012) [19] reviewed the various methods about crack initiation stress in previous literatures. In this study, Poisson’s ratio method was used to obtain the crack initiation stress (see Figure 9(b)). The ratio remains constant until cracking occurs. The stress associated with the change in Poisson’s ratio is referred to as the CI stress [19, 36]. And the volumetric strain method was used to obtain the damage stress (see Figure 9(a)). When the unstable crack begins to grow, the volumetric strain reaches maximum. The stress associated with maximum volumetric strain is referred to as the CD stress [19, 37, 38]. In general, Poisson’s ratio and volumetric strain are calculated based on the lateral and axial strains. In UDEC model, the lateral and axial strains are calculated using the average displacement in the and directions at dozens of locations in some sampling area to monitor the specimen response (Figure 10).

The CI stresses and are shown in Table 3 and Figure 11 by a series of calculations. As shown from Table 3, a general increase of the CI stress with the increase of flaw inclination angle β from 30° to 90° is observed. The CI stress decreased mainly when flaw inclination angle β < 30°. The CI stress decreased gradually with increasing of flaw length . The effect is relatively more prominent for 30°, 45°, and 60°.

displayed the similar law in Figure 11. decreased mainly when flaw inclination angle β < 45°. increased with the increase of flaw inclination angle β from 45° to 90°. increased with the decrease of flaw length . In addition, for the model with flaw inclination angle β of 0°, 15°, 75°, and 90°, changed from about 0.6 to 0.8. For the model with flaw inclination angle β of 30°, 45°, and 60°, changed from about 0.1 to 0.3.

The CD stresses and are also shown in Table 4 and Figure 12 by a series of calculations. The CD stresses and decreased with the increase of flaw length . And the CD stress at first decreased and then increased with the increase of flaw inclination angle β from 0° to 90°. decreased mainly when flaw inclination angle β < 30°. increased with the increase of flaw inclination angle β from 30° to 90°. In addition, when flaw inclination angle β < 15° or β > 75°, the CD stress is substantially equal to the peak strength. With the decrease of flaw length , the effect is relatively more obvious.

3.3. Description of Crack Behaviour

In a flaw-containing specimen, the block sizes have effect on the crack processes [33]. Wong results showed seven crack types with different trajectories and initiation mechanism (Figure 13). Based on Wong results, this section will analyse the crack types of model containing single flaw with different inclination angle and length. The crack behaviours of model containing single flaw with different inclination angle and length were analysed in Figures 14, 15, and 16. Four pictures of each model show the crack initiation, propagation, crack type at peak stress, and sketch of cracks.

For length of 10 mm and flaw angle β of 30°, 45°, and 60°, the crack initiation and crack propagation processes are shown in Figure 14. As can be seen from Figure 14, the first macroscopic cracks which are wing cracks (T1) are initiated from the tips of flaws. Cracks initiating later than the first cracks were usually referred to as secondary cracks. For the flaw angle β of 30°, a portion of the secondary cracks initiated upwards near the left tip, which are antitensile cracks (T3). A portion of the secondary cracks initiated downwards from the left tip, which are mixed tensile-shear cracks (S-T) (Figure 14(30-4)). For the flaw angle β of 45°, the secondary cracks consisted of some microcracks clusters initiated at a distance away from the preexisting flaw, which are tensile cracks (T) (Figure 14(45-4)). For the flaw angle β of 60°, three secondary cracks initiated from the left tip of flaw, two of them initiated downwards are mixed tensile-shear cracks (S-T), the third crack initiated upwards is antitensile crack (T3). One secondary crack initiated downwards from the left tip of flaw is antitensile crack (T3) (Figure 14(60-4)).

Figure 15 shows the crack initiation and crack propagation processes with length of 20 mm and flaw angle β of 30°, 45°, and 60°. For the flaw angle β of 30°, two macroscopic wing cracks initiated from the flaw tip as the first cracks (T1). Three secondary cracks initiated near the left flaw tip, which are substantially parallel to the first crack. The one initiated upwards is antitensile crack (T3). And other cracks are mixed tensile-shear cracks (S-T). Three secondary cracks initiated from right flaw tip, which are type 2 tensile cracks (T2) (Figure 15(30-4)). For the flaw angle β of 45°, two first cracks initiated from the flaw tip, which are wing cracks (T1). Two secondary cracks initiated upwards from right flaw tip are mixed tensile-shear cracks (S-T). Two secondary cracks initiated downwards from left flaw tip are type 2 tensile cracks (T2) (Figure 15(45-4)). For the flaw angle β of 60°, two macroscopic wing cracks (T1) initiated from the flaw tip as the first cracks and then four secondary cracks initiated from the right flaw tip, one of which is type 2 tensile crack (T2). Other secondary cracks are tensile cracks (T) (Figure 15(60-4)).

Figure 16 shows the crack initiation and crack propagation processes with length of 30 mm and flaw angle β of 30°, 45°, and 60°. For the flaw angle β of 30°, two macroscopic wing cracks (T1) initiated from the flaw tip as the first cracks. Four secondary cracks consisted of some microtensile cracks initiated near the left flaw tip region. Two of them extended upwards are type 3 tensile cracks (T3). One secondary crack initiated from right flaw tip is type 2 crack (T2) (Figure 16(30-4)). For the flaw angle β of 45°, two first cracks initiated from the flaw tip, which are wing cracks (T1). Two types of macroscopic cracks initiated from the flaw tip. The one initiated from the left flaw tip is mixed tensile-shear crack (S-T). While the one initiated from the right flaw tip is type 2 tensile crack (T2) (Figure 16(45-4)). For the flaw angle β of 60°, two macroscopic tensile wing cracks (T1) initiated from the flaw tip as the first cracks. One secondary crack as shear crack (S) initiated at distance from the right flaw tip. Two secondary cracks initiated near the left flaw tip are mixed tensile-shear cracks (S-T) (Figure 16(60-4)).

Five crack patterns are observed in single flaw with different geometries shown in Table 5: tensile crack (T), type 1 tensile crack (T1), type 2 tensile crack (T2), type 3 tensile crack (T3), mixed tensile-shear crack (S-T), and shear crack (S). Meanwhile, the crack types changed with different geometries. This indicated that the geometries have great effect on the crack type.

4. Discussion

From the loading of UDEC model, the peak strength, the CI stress, and the CD stress were obtained. The CI stress ratios ( at first decreased and then increased as the flaw inclination angle β increased, which was also observed in some references [7, 32, 39]. With the increasing of flaw length, the variation shape of the CI stress ratios changed from “V” into “U”; note that some points are not within the scope of the law. One plausible explanation was due to definition of the CI stress. In this paper, Poisson’s ratio method was used to obtain the crack initiation stress [19]. But in actual operation, the CI stress is defined with 10% change of Poisson’s ratio in some models. Another plausible explanation was due to random distribution of block sizes.

In the present numerical study, tensile crack (T), type 1 tensile crack (T1), type 2 tensile crack (T2), type 3 tensile crack (T3), mixed tensile-shear crack (S-T), and shear crack (S) were obtained. Compared with Wong literature [2, 3], some of the paths of cracking were different in the present study. One plausible explanation was due to influence of grain size. The grain size in UDEC influences the onset of cracking and possible initiating mechanisms [19]. Table 3 shows the crack process with length of 10 mm and flaw angle β of 30°, 45°, and 60° under different grain sizes. The number of microcracks increased with the increase of grain size. The first crack which is wing crack (T1) has the same trend. With the increase of grain size, the trend becomes better. Table 6 also shows difference about the type of secondary crack. For example, for the flaw angle β of 45°, the secondary cracks for appear to be tensile crack (T). But the secondary cracks for appear to be antitensile crack (T3), mixed tensile-shear crack (S-T), and shear crack (S) together. So, the grain size has effect on type of crack in numerical simulation. In particular, the smaller the grain sizes, the more the time required for calculations. In the future research, to obtain good results, a more reasonable block size should be adopted.

5. Conclusion

A series of numerical simulations on rock-like specimens containing single flaw with different inclination angle and flaw length were subjected to uniaxial compressive loads. The main conclusions can be summarized as follows:(1)The peak strength at first decreased and then increased with increasing of the flaw angle. Moreover, the peak strength decreased gradually when length of preexisting crack increased. When the inclination angle β was 90°and 0°, the flaw length has no influence on . When the inclination angle β was 30°, 45°, and 60°, increased with increasing of the flaw length .(2)The CI stress and the CI stress ratios ( depend on the inclination angle and flaw length. The CI stress, the CI stress ratios (, and the CD stress at first decreased and then increased as the flaw inclination angle β increased. With the increasing of flaw length, the variation shape of the CI stress and the CI stress ratios changed from “V” into “U”; the inclination angle and flaw length almost have no effect on .(3)The crack behaviour also depends on the inclination angle and flaw length. Four types of crack were observed in numerical simulations, involving tensile crack (T), type 1 tensile crack (T1), type 2 tensile crack (T2), type 3 tensile crack (T3), mixed tensile-shear crack (S-T), and shear crack (S).

Conflicts of Interest

The authors declare no conflicts of interest related to the publication of this paper.

Authors’ Contributions

Zhang Zhijun and Dai Bing contributed equally to this work.

Acknowledgments

This study was financially supported by the National Natural Science Foundation of China (51324744, 51774187, and 51374129), the Innovation Foundation of Doctoral Student in Hunan Province (71380100006), Key Research and Development Plan Projects in Hunan (2017SK2280), the Key Research Foundation of Education Bureau of Hunan Province (17A184), and the Innovation Team Project Plan for Prediction and Control of Uranium Mine Geotechnical Engineering Disaster of USC (NHCXTD04).