This study employed PFC2D to numerically simulate the microcrack propagation process in red sandstone. The simulated red sandstone exhibited a medium-grained crystalline structure with a uniform particle distribution and an average density of approximately 2,400 kg/m³, with relevant physical properties referenced from the literature (Yang et al., 2013). To analyze the microcrack propagation characteristics of red sandstone with double fissures, a rectangular numerical model measuring 80 mm × 160 mm was constructed, incorporating two pre-fabricated fissures. The inclination angle of the fissures was set at α = 45°, the fissure length at 2a = 16 mm, the fissure width at d = 2.5 mm, the rock bridge length at 2b = 22 mm, and the rock bridge inclination angle at β = 90°, as shown in Figure 1. The model consisted of particles with a minimum radius of 0.3 mm and a particle size ratio (maximum/minimum) of 1.5, resulting in the generation of 26,298 circular particles, 52,116 bonds, and 68,283 contacts. The micromechanical parameters of the red sandstone numerical model were set according to the parameters provided by Yang et al. (2014), as detailed in Table 1.

In particle flow simulations, the interaction between bonded particles is typically modeled using either a contact bond model or a parallel bond model. The parallel bond model is mechanically more complex than the contact bond model is, enabling the simulation of more refined contact behaviors and fissure mechanisms. It is better suited for representing the bonding characteristics of actual rock materials, whereas the contact bond model is more appropriate for simpler contact mechanics analyses (Lee and Jeon, 2011; Zhu and Yu, 2020). Therefore, to investigate the crack propagation characteristics of red sandstone with double fissures, this study adopted the parallel bond model as the interaction model between particles, as illustrated in Figure 2.

For the numerical simulations, a combination of servo-controlled loading and strain control was employed. Specifically, during the entire compression test, axial stress was applied to the upper and lower axial loading plates at a constant loading rate, whereas the lateral loading plates maintained a constant confining pressure through a servo control mechanism. The test was terminated when the axial stress of the specimen decreased to 80% of its peak stress.

Microcrack propagation is driven by changes in internal stress within the specimen. To investigate the microcrack propagation process under different stress levels, five key stress stages were selected for analysis: 60% of pre-peak stress, 80% of pre-peak stress, peak stress, 90% of post-peak stress, and 80% of post-peak stress. At these stress stages, the distribution characteristics of interparticle contact forces and parallel bond forces were dynamically tracked, recorded, and analyzed to explore the micromechanical mechanisms governing crack propagation.

To validate the rationality of the micromechanical parameters listed in Table 1, uniaxial compression tests were conducted on physical specimens of sandstone with double fissures using an MTS815.02 testing machine. Simultaneously, corresponding compression simulations were performed on numerical specimens. The stress-strain curves for both the physical and numerical specimens are presented in Figure 3, whereas their uniaxial compressive strength (UCS) and elastic modulus values are detailed in Table 2.

As shown in Figure 3, the failure patterns of the numerical and physical rock specimens were highly consistent. Both specimens exhibited tensile splitting failure in the upper portion, with the rock bridge showing tensile coalescence. Additionally, the failure surface inclinations in the lower part of the specimens were nearly identical. These observations indicated that the selected micromechanical parameters not only effectively captured the strength and deformation characteristics of the actual rock specimens but also accurately reproduced their fracture patterns, further validating the rationality of the micromechanical parameter selection.

Although the stress-strain curve from the numerical simulation did not perfectly match that of the physical specimen (Figure 3), the compressive strength and elastic modulus of the numerical specimen were closely aligned with those of the physical specimen (Table 2), demonstrating that the numerical model effectively represented the actual mechanical properties of the rock. The numerical specimen exhibited a lower strain than did the physical specimen, primarily because of the compaction stage in the numerical modeling process. In contrast, the physical specimen underwent a pore closure stage during loading, resulting in differences in the initial deformation characteristics.

Figure 4 shows the microcrack distribution characteristics of the double-fissure red sandstone under different the confining pressures. In the figure, red represents tensile microcracks, black represents shear microcracks, and green represents pre-fabricated fissures. As shown in Figure 4, the number of microcracks gradually increased with loading. The dominant crack type was tensile microcracks, which were mainly concentrated near the failure surface, whereas shear microcracks were more sparsely distributed.

Figure 5 presents the number of microcracks at the peak strength point under different confining pressures. The total number of microcracks increased with increasing confining pressure, with tensile microcracks significantly outnumbering shear microcracks. The number of tensile cracks increased stepwise with increasing confining pressure. Within the confining pressure intervals of 5–10 MPa and 15–20 MPa, the rate of change in tensile microcrack number is approximately six cracks/MPa and 27 cracks/MPa, respectively, indicating a relatively moderate increase. In contrast, within the intervals of 0–5 MPa, 10–15 MPa, and 20–25 MPa, the rate of change increases significantly, with similar values across these three ranges. The average rate of change in tensile microcracks within these intervals is approximately 138 cracks/MPa. For shear microcracks, the rate of change gradually increases with confining pressure, rising from 16 cracks/MPa in the 0–5 MPa interval to 58 cracks/MPa in the 20–25 MPa interval. Overall, the total number of microcracks exhibits higher rates of change in the intervals of 0–5 MPa, 10–15 MPa, and 20–25 MPa, showing a clear upward trend with increasing confining pressure, ranging from 156 cracks/MPa in the 0–5 MPa interval to 186 cracks/MPa in the 20–25 MPa interval.

Figure 6 shows the relationship between the number of microcracks and the stress-strain curve under a confining pressure of 15 MPa. In the initial loading stage, almost no microcracks are observed within the specimen. Microcrack initiation begins only around point a (σ _b = 0.6σ _p). As the load continues to increase, the rate of change in the number of microcracks gradually rises until reaching point e (σ _a = 0.9σ _p). Specifically, prior to the peak stress, the rate of change in total microcrack number with respect to strain increases from 4,262 cracks per 1% strain during the linear elastic stage (ab segment) to 14,763 cracks per 1% strain during the yield stage (cd segment). Upon reaching the peak stress, the rate reaches its maximum: during the de segment, the rate of change in total microcrack number is 44,018 cracks per 1% strain. However, after point e, this rate decreases to 30,786 cracks per 1% strain.

When the load reached point b (σ _b = 0.8σ _p), the number of microcracks inside the specimen began to increase rapidly, while the number of shear microcracks also started to grow steadily. As loading progressed to point c (σ _b = 0.9σ _p), the growth rate of microcracks further accelerated, and the specimen transitioned into the nonlinear elastic stage during the b–c phase. At the peak stress point d, the crack propagation rate reached its maximum, accompanied by the most intense crack expansion. During the c–d phase, the specimen entered the yielding stage, and slight fluctuations appeared in the stress-strain curve. After passing peak point d, the crack propagation rate continued to rise. However, after reaching point e, although the number of microcracks continued to increase, its growth rate started to decline.

The relationships between the MSTR and both confining pressure and peak strength are shown in Figure 7. As illustrated, the MSTR increases linearly with both parameters but consistently remains well below 1. Specifically, at a confining pressure of 0 MPa, the MSTR is 0.023, rising to 0.304 at 25 MPa. Fitting results show that the slopes of the MSTR with respect to confining pressure and peak strength are 0.01106/MPa and 0.005/MPa, respectively, indicating that confining pressure has a more pronounced effect on the generation of shear microcracks. Nevertheless, the number of shear microcracks remains significantly lower than that of tensile microcracks throughout.

Figure 8 shows the distribution of interparticle contact forces within double-fissure red sandstone under different confining pressures. In the figure, black represents compressive stress, red represents tensile stress, and green represents pre-fabricated fissures. Throughout the entire experiment, the contact forces between the particles were predominantly compressive.

At a given confining pressure level, during the pre-peak stage when the applied load was relatively low (e.g., σ _b = 0.6σ _p), a contact force concentration zone formed at the fissure tips, whereas the contact forces in the rest of the specimen were relatively uniformly distributed. As the load increased, the contact force concentration zone gradually shifted away from the fissure tips, and the maximum contact force progressively increased. During the post-peak stage, the contact force chains in the rock bridge region became increasingly sparse, and the maximum contact force first increased but then decreased. For example, under σ ₃ = 0 MPa, when σ _b = 0.6σ _p, the maximum contact force was 42.6 kN. At the peak strength, the maximum contact force increased to 98.7 kN. During the post-peak stage, the maximum contact force initially increased to 110 kN (σ _a = 0.9σ _p) and subsequently decreased to 99.7 kN (σ _a = 0.8σ _p).

In the initial loading stage, the contact force chains near the rock bridge are uniformly distributed. As the load increases toward the peak strength, these force chains begin to concentrate along the line connecting the tips of the rock bridge cracks, while becoming sparser in the surrounding areas. After the peak load, particularly when the post-peak strength drops to 0.8 times the peak value, contact force chains in parts of the rock bridge region disappear, resulting in a more scattered and banded distribution. However, with increasing confining pressure, the contact force chains within the specimen become progressively denser under the same pre-peak load ratio, and the direction of contact forces on both sides of the fissures gradually shifts from vertical to align with the inclination of the pre-fabricated fissures.

Figure 9 shows the distribution of parallel bond forces within double-fissure red sandstone under different confining pressures, where black represents compressive stress, red represents tensile stress, and green represents pre-fabricated fissures.

As shown in Figure 9, when σ _b = 0.6σ _p, the parallel bond compressive stress was concentrated at the fissure tips, whereas the parallel bond tensile stress was primarily distributed on both sides of the pre-fabricated fissures. As the load increases to σ _b = 0.8σ _p, the tensile stress regions on both sides of the fissures extended upward and downward, whereas tensile stress concentration began to appear in the rock bridge region.

At peak strength, the compressive stress concentration at the fissure tips decreased, whereas the compressive stress in the rock bridge region intensified. As loading continued in the post-peak stage, the density of both parallel bond compressive and tensile stress distributions gradually decreased, indicating a continuous attenuation of parallel bond forces in the post-peak phase.

With increasing confining pressure, the parallel bond compressive stress remained primarily concentrated at the fissure tips, suggesting that the location of the compressive stress concentration zone was independent of the confining pressure. As the confining pressure increased, the rock bridge region was predominantly subjected to parallel bond compressive stress, whereas the area of parallel bond tensile stress concentration gradually diminished.

The crack propagation characteristics under different confining pressures indicate that the crack propagation process is closely related to the confining pressure. The confining pressure suppresses the propagation of wing cracks at the fissure tips while promoting the initiation and propagation of reverse wing cracks. Under high confining pressure, this process eventually leads to the formation of a reverse shear plane. This finding is consistent with the studies of Yavari et al. (2021), Wang J. et al. (2021), and Wu et al. (2024) on materials such as concrete, further validating the influence of confining pressure on crack propagation.

The underlying mechanism is that, compared with unconfined conditions (uniaxial compression), a low confining pressure imposes an additional lateral constraint, reducing the tensile stress at the fissure tips and thereby inhibiting wing crack propagation. Consequently, under a low confining pressure, the development of wing cracks is significantly suppressed. However, at higher confining pressures, the tensile stress at the fissure tips is further inhibited (see Figure 9), leading to a shift in the crack propagation mode and the formation of reverse wing cracks.

Furthermore, as the confining pressure increases (see Figure 9), the overall failure mode of the specimen gradually transitions from tensile splitting failure to shear failure. The inclination angle of the primary failure surface on the right side of crack ① progressively changes from 90° to align with the pre-fabricated crack inclination angle (45°). This observation aligns with the findings of Li and Ghassemi (2021) and Du et al. (2021), further confirming the effect of confining pressure on failure modes. Specifically, under confining pressure, as the applied load increases, the local tensile stress exceeds the local tensile strength, leading to the formation of tensile cracks. With the accumulation and coalescence of these tensile cracks, a failure surface extending in the shear direction eventually develops. Notably, although the failure mode under confining pressure appears to be shear failure, the formation of the failure surface is primarily governed by the expansion and coalescence of tensile microcracks.

The variation trend in the number of microcracks with confining pressure (Figure 5) shows that tensile microcracks follow a stepwise increase as the confining pressure increases, indicating that 0 MPa, 10 MPa, and 20 MPa are critical confining pressures at which the number of tensile cracks surges. In the 0–5 MPa, 10–15 MPa, and 20–25 MPa ranges, the growth rate of tensile microcracks is relatively high and remains similar across these intervals. In contrast, the growth rate of shear microcracks gradually increases with increasing confining pressure, leading to an overall increase in the growth rate of total microcracks at the critical confining pressure thresholds. This finding aligns with the results of Wang Z. et al. (2024) in their study on the effect of confining pressure on crack quantity in marble specimens, further confirming the role of confining pressure in the crack propagation mechanism.

The relationship between the microcrack quantity and the stress-strain curve (Figure 6) indicates that in the early to mid-stage of the linear elastic phase (with a pre-peak load ratio of approximately 0.6, point a in Figure 6), the number of microcracks inside the specimen remains relatively low. However, in the late linear elastic phase and continuing until the post-peak load ratio reaches 0.9 (point e in Figure 6), the number of tensile microcracks exhibits an approximately exponential increase, whereas the number of shear microcracks also show nonlinear growth, with an accelerating growth rate. This phenomenon can be attributed to the stress concentration effect at the tips of pre-fabricated fissures. Before point a, the parallel bond strength between particles is lower than the interparticle contact force, resulting in limited microcrack initiation. However, as the load increases, once the parallel bond strength reaches or exceeds the interparticle contact force, microcracks begin to form. The increasing stress concentration near newly formed microcracks further accelerates the generation of additional microcracks.

However, when the post-peak load ratio decreases to less than 0.9, the number of both tensile and shear microcracks continues to increase, but the growth rate decreases. We speculate that this may be due to the coalescence of secondary cracks into a fully developed failure surface in the post-peak stage. At this stage, deformation is primarily concentrated along the failure surface, and the reduction in post-peak strength limits the formation of new microcracks, thereby slowing their growth rate.

The analysis of the MSTR as a function of the confining pressure and peak strength (Figure 7) reveals that the MSTR is positively correlated with both parameters, indicating that shear microcracks are more sensitive to the confining pressure than tensile microcracks are. This finding is consistent with the studies of Zhang et al. (2021) and Yu and Duan (2021) on crack propagation in coal specimens. However, although the MSTR increases with increasing confining pressure and peak strength, and shear microcracks exhibit higher sensitivity to the confining pressure, the MSTR value remains below 0.35, indicating that the number of shear microcracks is significantly lower than that of tensile microcracks. While the overall failure mode of the specimen appears to be shear failure, its fundamental cause is still the coalescence of tensile microcracks.

The distribution of contact forces indicates that when the pre-peak load ratio is relatively low (e.g., σ _b/σ _p = 0.6), the contact force is primarily concentrated at the tips of the double fissures, suggesting that microcracks initiate at the fissure tips. As the pre-peak load ratio increases to 0.8, the overall contact force increases, and the concentration zone gradually shifts away from the fissure tips. This shift leads to the formation of wing cracks (Figures 4a–2,b–2) and reverse wing cracks (Figure 4c–2), while also promoting the further expansion of reverse wing cracks (Figure 4d–2). When the pre-peak load ratio reaches 1 (i.e., at peak strength), the contact force chains in the rock bridge region gradually become sparse, indicating a reduction in compression between particles. This reduction facilitates the progressive propagation, coalescence, and eventual linkage of secondary cracks in the rock bridge region.

These findings further confirm that the redistribution of contact forces is a driving mechanism for crack propagation. Specifically, the shift in the concentration zone of contact forces (caused by increased loading), the intensified concentration of compressive stress, and the deflection in the contact force direction (caused by increased confining pressure) are the key factors promoting the initiation, propagation, and coalescence of microcracks.

The distribution of parallel bond forces reveals that as the load increases, the parallel bond tensile stress on both sides of the fissures transfers toward the rock bridge region, where the tensile stress reaches its maximum. This suggests that the rock bridge region is the primary zone for the initiation of tensile microcracks. However, a compressive stress concentration zone remains at the fissure tips (Figure 9), indicating that the fissure tips still experience significant compressive stress.

Figures 8, 9 show that the formation of wing cracks at the fissure tips is primarily due to the concentration of compressive stress, leading to the fracture of parallel bonds. In contrast, both the wing cracks and reverse wing cracks in Figure 4 mainly consist of tensile microcracks, indicating that parallel bond failure is predominantly tensile in nature. This phenomenon can be explained by the compression-induced tensile effect, whose fundamental mechanism is illustrated in Figure 10. Similar conclusions have been reported by Zheng et al. (2020) and Suknev (2011) in their studies on rock and rock-like materials, further supporting this interpretation.

Additionally, the density of parallel bond forces gradually decreases in the post-peak stage as loading continues. This may be attributed to the progressive coalescence of microcracks, forming macro-scale secondary cracks, coupled with a decrease in specimen strength. This redistribution of stress ultimately leads to the overall attenuation of parallel bond forces.