Abstract
To study the effects of anisotropy and heterogeneity on the shale failure mode and tensile strength, Brazilian splitting tests were performed from both directions of the bedding and layer thickness. Layers containing different bedding and loading angles and layer thicknesses were obtained separately. The results show that, at 0° and 90° angles, the shale cracks grow “linearly”; at 15°, the shale cracks have “arc type” growth; and at 30°–75°, the shale-splitting displays “broken line” crack propagation. The tensile strength from 0° to 90° exhibits an increasing trend. Water has a significant softening effect on the tensile strength of shale—the higher the water content, the lower the tensile strength. In addition, a 3DEC numerical simulation was used to simulate the tests, establishing shale specimen particles with random blocks. In the shale disc, uneven parallel bedding and uniform parallel bedding were set up with different loading angles and layer thicknesses to generate simulated stress-displacement curves, and the effect of layering on shale cleavage was analyzed from a mesoscopic perspective. The tensile strength of shale with uniform parallel bedding was found to be higher under the same conditions, which is consistent with the experimental results. By comparing the experimental and simulation results, from both the macro- and mesoperspectives, the Brazilian splitting crack growth of shale is affected by bedding, displaying a process from disorder to order. This study is of great significance for further exploration of the mechanical properties of shale under loading failure.
1. Introduction
As an unconventional energy source, shale gas has attracted much attention from researchers. The development of shale gas requires large-scale hydraulic fracturing of a shale rock mass to form a fissure grid structure in the rock mass. In addition, shale bedding is well developed [1], and there are many natural fractures hidden in the rock mass with strong heterogeneity and complicated fracture propagation laws. A thorough understanding of the crack propagation law and its influencing factors when shale is stressed is of great significance. Moreover, it is helpful to control the development of the fracture network and improve the production of shale gas.
The Brazilian splitting test indirectly measures the tensile strength of rocks through compression. Owing to the simple operation and easily prepared test pieces, this test has become a universal method for measuring the tensile strength of materials. Many factors affect the experimental results. Huang et al. [2] established a Brazilian splitting mechanical model under a string load and obtained an approximate analytical solution of the internal stress of a disk. It was found that the greater the loading angle of the platform, the lower the stress concentration at the loading place, and the smaller the tensile zone in the disk. Peng et al. [3] carried out a Brazilian splitting test of rocks at alternate loading rates and found that the rock splitting strength increased with increasing loading rate. In addition, through the Brazilian splitting test, crack initiation and propagation laws when the rock splits can be observed. Teng et al. [4] conducted a Brazilian splitting test of Longmaxi shale in combination with acoustic emission and found that the bedding and loading angles affected the shale failure mode. Du et al. [5] combined digital image correlation technology to perform the Brazilian splitting test of carbon shale at different dip angles and found that, with the increase of bedding and loading angle, the split crack initiation time increased, and the extension to penetration time decreased. Gao et al. [6] collected shale from the Yanchang Formation in Ordos and conducted a Brazilian splitting experiment to measure the brittleness index and horizontal differential stress of the shale. It was also found that fractures propagated along the natural fracture during fracturing of the shale. Furthermore, the combination of experimentation and numerical simulation can verify the crack propagation law. Bilgen Carola et al. [7] conducted a Brazilian splitting test and established two- and three-dimensional Brazilian split disc numerical models, obtaining a phase-field method that can predict crack initiation and propagation.
The discrete element method is a dynamic numerical simulation method based on deformable elements, which can simulate the shape and property changes after a block is loaded. Based on the unique separability of its basic unit, the discrete element method is more suitable for simulating plastic mechanical problems such as deformation and failure of massive and fissured media [8]. Wang et al. [9] combined the Neper program with 3DEC software to establish a demonstration model containing rigid, elastic, and fragile particles and found that the assignment of particles and the contact parameters between particles affect the macroscopic mechanical properties of rock compressive strength. Li et al. [10] used the DFN module in 3DEC to simulate the joints in a rock mass for analysis of the discontinuity of the rock material. Meanwhile, Huang et al. [11] used the PFC program to establish a discrete element model of a center-grooved Brazilian disc with different particle sizes. It was found that the larger the particle radius of the specimen, the greater the tensile strength.
Although a lot of experimental and numerical simulation research has been conducted on the extension of shale tensile cracks, there is no systematic study on the effects of the angle and layer thickness of bedding, bond strength, and natural cracks on the extension of tensile cracks. It is necessary to conduct in-depth research from the macro- and microperspectives. In this study, Brazilian splitting tests were conducted on shale specimens from the Longmaxi Formation in the Sichuan Basin, China, to macroscopically study the laws of shale tensile crack propagation. Further, mesomodelling was performed using 3DEC to observe the changes in shale morphology during the splitting process. Stress and displacement fields were generated for a comparative analysis with the experimental results to reveal the propagation rule of shale Brazilian split cracks under the influence of bedding and layer thickness.
2. Numerical Methodology
2.1. Basic Theory of Discrete Element Method
As shown in Figure 1, 3DEC considers that the interaction force between the two blocks is related to the contact between the blocks [12].
(a)
(b)
For elastic deformation, suppose that the increment of the normal force vector between the blocks is directly proportional to the increment of the normal displacement , and then,where is the normal stiffness coefficient of the joint, and is the contact area between the blocks.
Similarly, the increment of the tangential force vector between blocks iswhere is the tangential stiffness coefficient of joints, and is the increment of the tangential displacement. Index i takes values of 1–3 and denotes the components of a vector in the global coordinate system.
The total normal force and total tangential force arewhere and are the initial normal force and initial tangential force, respectively.
Considering the plastic deformation of rock failure, the normal tension is greater than the tensile strength of the contact . If the Mohr–Coulomb criterion is used to calculate the shear strength, when the tangential force exceeds , slip occurs between the blocks, and the tangential force takes the limit value at this time. This indicates thatwhere and are the contact cohesion and friction angle, respectively.
2.2. Building the Model
The numerical simulation uses a combination of 3DEC and Neper. Neper is a free and open source software, which can realise a 3D Voronoi mosaic model and further convert it into 3DEC format code [13]. 3DEC is a three-dimensional computer numerical program based on discrete element simulation, which can calculate the deformation of rock blocks by establishing deformable blocks and different contacts [8].
First, Neper was used to generate a 3D particle mosaic and synthesize a Brazilian disc model. The simulation specification was Φ50 mm × 1 mm. Then, the disc was imported into 3DEC, and the powerful DFN module was used in 3DEC to generate layering.
As shown in Figure 2, the shale disc was set as block Mat 1, and the upper and lower parallel loading plates were set as block Mat 2. The joint material was divided into four parts. The contact surface between the upper and lower loading blocks and the shale disc is Jmat 1. The preset bedding surface is Jmat 2. The virtual joint in the 3DEC simulation is Jmat 3 [14], and the contact surface between shale particles is Jmat 4.
(a)
(b)
(a)
(b)
(c)
(d)
(e)
As shown in Figures 3(a)–3(d), taking the angle between the bedding direction and the vertical loading direction as the bedding angle, a preset bedding of 0°–90° was established. The bedding was a long strip, and its length did not penetrate the disc, whereas its width was consistent with that of the disc. As shown in Figure 3(e), five observation points were set on the centerline of the disc from top to bottom at equal intervals, numbered 1–5, and point 3 is the center of the circle.
Uniformly and nonuniformly parallel bedding were established, and samples were modelled regarding the bedding direction and bedding spacing. When building the bedding models with different layers, the bedding angle was set to 0°. Set the block model to linearly elastic, while set the joints model to area contact elastic/plastic with Coulomb slip failure. Failure of joint in shear or tension results in the use of cohesion, tension, and friction residual values. The default residual values of cohesion and tension are 0. The lower boundary was fixed, and a loading speed boundary of 0.1 mm/min was set at the upper loading block node, which is consistent with the experimental loading speed. Changes in the stress field and displacement field were recorded at the five observation points.
2.3. Calibration of Microparameters
To build the model, values must be assigned to different blocks and contacts in the model. Because the particle and bedding parameter values cannot be measured experimentally, it is necessary to first assume the parameters that need to be assigned and further adjust them according to the results.
Wang et al. [15] proposed that, in 3DEC, the joint normal stiffness can be expressed aswhere is the normal stress.
Considering the operation efficiency, the stiffness value of the adjacent block is greater than 10 times the bedding stiffness value [12], which iswhere and are the normal and shear stiffness of the bedding, respectively; and are the bulk and shear modulus of the adjacent block, respectively; and is the minimum width of the normal direction of the adjacent block in the bedding.
According to classical rock mechanics, the bulk modulus of a block is given byand the shear modulus iswhere is the bulk Young's modulus, and is the bulk Poisson's ratio.
The calculation results were synthesized, with reference to the relevant literature for the assignment of the microparameters of the weak surface of the rock [14, 16–21] and large numbers to the loading block to prevent deformation of the loading block. Finally, the simulation results conforming to the test rules were obtained. The calibration parameters used are listed in Tables 1 and 2.
3. Results and Discussion
3.1. Model Validation
As shown in Figure 4(a), black homogeneous shales from the Longmaxi Formation in the Sichuan Basin were taken and processed in accordance with the International Society for Rock Mechanics [22] standards and prepared into Φ50 mm × 25 mm standard-specification disc test pieces by coring, cutting, and grinding. The parallelism of the upper and lower surfaces of each test piece was controlled within 0.5 mm, and the flatness of the surface was controlled within 0.1 mm. All test pieces were stored in a dry environment at room temperature.
(a)
(b)
The Brazilian splitting test involves platform loading. To ensure that the loading direction line passes through the center of the disk, the specimen was first marked in the preapplication direction, and then the upper loading point was marked. When placing the test piece, a vertical line was used for calibration. After fixing the test piece, the lower platform was fixed, and the upper platform was applied with a loading speed of 0.1 mm/min. The stress-displacement curve of the test was recorded and observed. When a sudden stress drop was measured, the test piece formed a crack and stopped loading. The shape of the specimen after splitting is shown in Figure 4(b). The splitting crack penetrated the upper and lower loading ends and developed vertically through the center point of the disc.
A homogeneous shale disc model without bedding was created in 3DEC software for loading, and the disc displacement field and stress field evolution were recorded.
As shown in Figures 5 and 6, similar to the experimental splitting results, the splitting crack of the homogeneous model extended vertically through the upper and lower loading ends. After splitting, the displacement field on the right side of the crack was larger than that on the left side. The analysis shows that, because of the damage to the particles at the loading end after splitting, the disk homogeneity was destroyed, and the forces on both sides were uneven, resulting in an asymmetric distribution of the stress field.
As shown in Figure 7(a), the tensile stress-main displacement curves of the center point of the test and the simulated disc were compared. During the test loading process, the vertical tensile stress at the center point increased linearly with the loading displacement. Before the splitting failure, the stress increasing trend rarely decreased. The peak value of the visible tensile stress is the “tensile strength” of the shale specimen [17]. More specifically, the tensile strength of the black shale measured by the test specimen was 7.25 MPa. The maximum stress peak at the center of the circle when the homogeneous model was split by simulation was 6.99 MPa, which was 3.59% different from the tensile strength obtained from the experiment. As shown in Figure 7(b), the split peak stresses of observation points 2 and 4, with the center of the circle as the symmetry point, are approximately 70% of the peak stress at the center of the circle, and the same symmetrical stress peaks at points 1 and 5 are the smallest, approximately 15% of those at the center of the circle. During splitting, the displacements of the five observation points were similar, and the maximum displacement was 26.2 µm at the center point.
(a)
(b)
For the bedding shale, which is also from the Longmaxi Formation, two sets of samples were further drilled from 0°, 45° and 90° bedding directions. The samples were processed to create standard test pieces, and the first group of samples was recorded as A-0, A-45, and A-90, according to the angle from small to large; the second group of samples was B-0, B-45, and B-90.
As shown in Table 3, a shale X-ray diffraction analysis was conducted to analyze the mineral composition of the test piece. The clay minerals and nonclay minerals in the two groups of shale were equivalent, each accounting for approximately 50%. Nonclay minerals have the highest quartz content, while clay minerals have the highest illite content.
As shown in Figure 8, for the specimen in the parallel bedding direction, the cracks penetrated the center of the disc along the bedding direction during splitting. Affected by the mineral composition of the specimen, when the main crack develops vertically, the other cracks will split and produce secondary cracks in the vertical direction. Meanwhile, for the test piece in the vertical bedding direction, the crack connects the upper and lower loading ends and is developed into an arc, reaching the maximum arc end at the lateral diameter. At 1/4 of the disk, cracks easily fork to produce lateral secondary cracks. For the 45° bedding direction, the crack connects the two loading points, forming a broken line through the specimen. The main crack first split down as an arc. During the expansion process, it folded in the parallel bedding direction under the influence of bedding. In the vicinity of the lower loading end, the arc crack developed to the lower loading point. In the tests, when the 45° and 90° bedding specimens were split, the crack did not pass through the center of the disc, indicating that the crack was tensile and damaged, and shear slip was also generated. This phenomenon reflects the influence of bedding direction on crack propagation.
(a)
(b)
(c)
(d)
(e)
(f)
As shown in Figure 9, the experimentally obtained shale tensile strength increased with the increase in the bedding and loading angle. At 90°, the tensile strength can be regarded as the tensile strength of the shale blocks. At 0°, the tensile strength can be considered as that of the bedding in the shale. The ratio of the two was 3.25 in group A and 1.5 in group B, indicating that the anisotropy of shale bedding has a significant effect on tensile strength.
3.2. Different Angles
For different angles, as shown in Figure 10, the stress and displacement of the uniform parallel bedding disc model were changed during the loading process to obtain the corresponding stress field and displacement field. The simulation data were recorded during the loading process at the center point of the bedding, and the corresponding stress-displacement curve was obtained.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
As shown in Figures 10 and 11, for uniformly parallel shale models with different bedding angles, at 0° and 90°, the split crack started to form from the upper loading end and developed in a straight line through the centerline of the disc. At 15°–60°, the crack began forming from the upper loading end and first propagated along the vertical center line under the action of vertical tension. During the expansion process, the crack was shifted by the influence of the bedding angle on the lower half of the disc. The bedding angle dominates the direction, and the crack split near the bedding at the lower loading end; the overall crack presents a “Z” broken line. At 75°, the crack began forming from the upper loading end and spread in a small-angle arc to the adjacent bedding of the lower loading section. For nonuniform parallel shale models with different bedding angles, the splitting crack trend is the same as that of the uniform parallel model. Further, it is easier to see the effect of anisotropy on crack development. When splitting at 0° and 45°, secondary cracks evolved at the lower end of the center of the circle, and the development direction was the same as that of the primary crack.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
As shown in Figures 12 and 13, in the splitting stress field cloud of the shale disc model at different angles, the stress value measured at the splitting crack is larger than that measured at any other place, which basically reflects the splitting crack shape of the shale model at different angles. In general, the shale model with uniform parallel bedding exhibited a larger splitting stress value than that with uneven bedding. For bedding shales with different angles, the overall value of the stress field of the 75° model is the smallest.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(a)
(b)
(c)
(d)
(e)
(f)
(g)
As shown in Figures 14 and 15, the stress-displacement curves generated at different monitoring points have the same trend: the stress first increases linearly with the displacement. Before splitting, the stress increased suddenly and then decreased at the instant of splitting. The cleavage stress was the largest at the center of the circle, followed by the stress values at symmetry points 2 and 4 near the center of the circle; and those at symmetry points 1 and 5 near the two loading ends were the smallest. The symmetrical points of the stress-displacement curves can be seen to approach each other.
(a)
(b)
(a)
(b)
As shown in Figure 16, the stress and displacement processes of the center point of the shale disc in different bedding directions were recorded. It was found that the bedding angle in the shale significantly influences its tensile strength. The tensile strength increased with the angle from 0° to 90° in a zigzag manner. For the models with uniform parallel bedding, the splitting strength was the largest at 90°, of 11.98 MPa, and the splitting strength at 0° was smaller, of 6.73 MPa. Their ratio is 1.78. For the models with uneven parallel bedding, the splitting strength was the largest at 90°, of 9.99 MPa, and the splitting strength was smaller at 0°, of 7.72 MPa. The ratio between the two is 1.29. For the 90° model, the shale particles were mainly tensioned during loading, and the splitting strength can be approximated to the tensile strength of the particle blocks, which have a large strength value. Meanwhile, for the 0° model, the weak layer of the middle shale was loaded when the load was applied. The splitting strength can be approximated to the tensile strength of the bedding surface, and the strength value is small.
(a)
(b)
3.3. Different Layer Thickness
As shown in Figures 17 and 18, for parallel bedding shale models with different layer thicknesses, the main crack propagated vertically through the disc. As the layer thickness decreased, the number of bedding increased, and secondary cracks were generated from the middle of the crack in the bedding direction near the loading end. In Figure 18, “2 and 6” indicate the numbers of beds on the left and right semidisk, respectively, of the shale model with uneven bedding.
(a)
(b)
(c)
(a)
(b)
(c)
As shown in Figures 19 and 20, the stresses at the splitting cracks with different layers of shale were more concentrated, which characterizes the vertical cracks. The numerical value of the stress field of a complete unstratified disk is small. The stress field of a model with uneven bedding is generally smaller than that of a homogeneous bedding.
(a)
(b)
(c)
(a)
(b)
(c)
As shown in Figures 21 and 22, for shales with different bedding layer thicknesses, it can be seen that as the layer thickness decreased, the number of bedding increased, and the difference between the stress-displacement curves of recording points 2 and 4 became larger. For the model with a uniform parallel layer thickness of 3.6 mm, the peaks of splitting stress and splitting displacement at point 4 at the lower end are greater. For the 10 and 18 model with the most parallel uneven layering, the peak splitting stress peak at point 2 at the upper end is greater.
(a)
(b)
(a)
(b)
As shown in Figures 23 and 24, for the shale models with uniform and parallel bedding of different thicknesses, the complete disc model without bedding had the smallest center displacement value when splitting. The model with the layer thickness of 5.0 mm exhibited the largest splitting strength of 7.64 MPa, and the model with the layer thickness of 3.6 mm had the smallest splitting strength of 6.79 MPa. The ratio between the two was 1.13. For the bedding models with uneven parallel and different layer thicknesses, the stress-displacement curves of the center of each circle are not very different, and the maximum cleavage strength was 8.02 MPa in the model with the largest bedding number of “10 and 18,” while the model with the bedding numbers of “6 and 12” achieved the minimum splitting strength of 6.37 MPa. The ratio between the two is 1.26. The curve first decreases and then increases. It was found that the shale first undergoes shear failure along the bedding, and the stress decreases; then, the shale particles are compressed, and the tensile stress increases. After the particle tensile strength is reached, interparticle stretching causes damage.
(a)
(b)
(a)
(b)
3.4. Bonding Strength
Studies [23] have shown that shale mainly contains minerals such as clay and quartz feldspar, and different geological shales in different regions have different proportions of mineral components. Different mineral components inevitably lead to different bonding strengths among the shale particles. Studies [24] have shown that the bonding strength between soil particles is mainly affected by cohesion. A 45° uniform parallel bedding model was established. The cohesion of each contact was changed in the same proportions of 0.8, 1.2, and 1.5 times and compared with the original ratio. The simulated bond strength between the bedding and particles differed, and the results are shown in Figures 25–28.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
(a)
(b)
(a)
(b)
As shown in Figures 25 and 26, as the cohesive force between the contacts increased, the split crack decreased, and the displacement and stress at the corresponding split decreased. The smaller the cohesion, the lower the shale bond strength. With the same loading, on the one hand, the 0.8 times cohesive force model formed a through crack, and the overall disk was stressed. On the other hand, the 1.5 times cohesive force model formed dispersed fine cracks, while only the crack area was significantly stressed, and the stress value was small.
As shown in Figures 27 and 28, the stress-displacement change trends of the different contact cohesive shales were the same before splitting failure. With the increase in cohesive force, the splitting strength increased, and the displacement during splitting increased; the measured ratio of the shale tensile strength between 1.5 times and 0.8 times is 2.29. This indicates that an increase in the cohesive force of each contact leads to an increase in the bond strength between the shale particles and the bedding, so that increased energy is required for shale destruction under the same load, and the corresponding shale tensile strength increases. The corresponding splitting strength change is not proportional to the cohesive force change.
3.5. Natural Cracks
As shown in Figure 29, shales with a 45°-direction crack, two 15° and 75° cross cracks, three 0°, 45°, and 90° cross cracks, and horizontal and vertical randomly distributed fracture grids were established. The models were then loaded.
(a)
(b)
(c)
(d)
As shown in Figure 30, for shale models with natural fractures, the split cracks are tortuously affected by the cracks. In the model with only one crack, the split cracks developed directly along the 45° natural crack. For the models with two or three natural cracks, the main crack split through the crack intersection. For models with natural crack grids, the cracks were clearly affected by the natural cracks and became randomly distributed along the crack grid. Random and small secondary cracks formed in the area near the natural cracks. As shown in Figure 31, the stress concentration areas easily formed at the natural cracks.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
4. Conclusions
In this study, 3DEC and Neper software were used to establish the Brazilian splitting models of shale with different angles and layer thicknesses, combined with experimental verification to analyze the influence of bedding geometry and shale cement strength on the tensile strength. The conclusions are as follows: (1)Bedding angles significantly affect the shale-splitting mode. At 0° and 90°, the crack penetrated the upper and lower loading sections through the center of the disk; at 15°, the crack extended through the upper and lower loading ends with an arc-shaped expansion. At 30°–75°, the crack was initiated from the top loading end and shifted from the center line to a polyline-shaped expansion.(2)Bedding significantly affects the tensile strength of shale. For shales with different bedding angles, the ratio of tensile strength measured in the 90° and 0° directions can reach a maximum of 1.79. For shales with different layer thicknesses, the maximum tensile strength ratio can reach 1.26. Thus, shale is clearly anisotropic.(3)Bond strength affects the tensile strength of shale. As the cohesion between the contacts increased, the shale bond strength increased, and the splitting cracks decreased. The ratio of the shale tensile strength measured by the 1.5 times and 0.8 times cohesion models was 2.29.(4)Natural cracks affect the propagation of Brazilian splitting cracks in shale. The cracks are easily distributed along the natural cracks and form small secondary cracks around the natural cracks; stress concentration areas are easily formed at the natural cracks.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that they have no conflicts of interest to report regarding the present study.
Acknowledgments
The research was funded by National Natural Science Foundation of China (No. 11702296), supported by Natural Science Foundation of China (41941018), and supported by the Fundamental Research Funds for the Central Universities (No. 2462019YJRC011).