Mesoscopic instability mechanism and grouting reinforcement of deep roadway surrounding rock with weak interlayers

Abstract

Constructing deep underground tunnels faces severe challenges when the surrounding rock contains weak, unstable layers. This study investigates how the angle and spacing of these weak layers lead to the deformation and collapse of deep tunnel structures. A particle flow simulation, based on a 1,000–m-deep tunnel, was developed to recreate the physical behavior of the surrounding rock and observe how it cracks, shifts, and loses its internal support under various geological conditions. The findings indicate that steeper weak layers cause the rock’s internal support network to break down more severely, resulting in a peak crack density of approximately 1,250 per square meter and a fractal dimension of 1.5848 in the highly fractured zones. Conversely, a larger distance between these weak layers preserves more solid rock, significantly reducing overall damage, restricting the maximum displacement of the roadway roof to between 0.3 m and 0.4 m and maintaining a lower fractal dimension of approximately 1.33. The primary cause of tunnel failure is the widespread internal fracturing which in field conditions led to cumulative floor heaving exceeding 6 m and rib shrinkage exceeding 4 m prior to reinforcement. To control this instability, the injection of a grouting slurry into surrounding rock. This cementing slurry intervention glues the broken rock fragments back together, successfully rebuilding the internal support network and optimizing the fractal dimension of crack distribution. Ultimately, using this grouting reinforcement strategy fundamentally repairs the rock’s internal connections, transforming the broken material back into a stable structure.

1 Introduction

As shallow resources become largely depleted, deep underground excavation under extreme high-stress conditions has become increasingly critical (Kang et al., 2023). However, deep mine roadways are prone to nonlinear large deformations under high stress conditions (Yang et al., 2017; Ma et al., 2025). Compounded by complex geological features such as weak interlayers, frequent disasters like collapses, roof falls, rib failures, and floor heaving occur during roadway excavation and maintenance (Zuo et al., 2013; Ren et al., 2024). These incidents result in devastating economic losses and casualties, presenting a significant challenge to the safe and efficient utilization of China’s coal resources (Zhang et al., 2025).

To address the challenge of maintaining the stability of the roadway surrounding rock in deep underground environments, a great number of researchers have made notable achievements in both theoretical studies and technical applications. Through elastoplastic theory, the stress distribution and zoning characteristics of the surrounding rock along with their evolution patterns were predicted (Zhang et al., 2012). The importance of synergistic load-bearing between anchor cables and the surrounding rock was recognized. A core concept matching active support with stress regulation was proposed, leading to the gradual refinement of anchoring theory (Wu et al., 2024). A series of advanced support materials and technologies have been developed, including high-strength anchor bolt/cable, energy-absorbing steel supports, and negative Poisson’s ratio anchor cables. These high-performance support components have expanded beyond the core active reinforcement methods of anchors, mesh, cables, and shotcrete (Jiao et al., 2013; Kang et al., 2015; Saxena et al., 2016; Jing et al., 2020; Khaleghparast et al., 2023; He et al., 2023). A synergistic control system for the bolting-modification-destressing for roadway surrounding rock has been developed (Kang et al., 2021). Beyond these established active support systems, the development of advanced grouting materials has recently emerged as a critical research direction for deep underground engineering. Current literature highlights the implementation of high performance chemical grouts and nanomaterial reinforced slurries that offer superior penetration and bonding capabilities (Liu et al., 2021; Onaizi et al., 2021). These advanced materials enhance the mechanical properties of deep rock masses through micro scale reinforcement mechanisms such as crack filling and the reconstruction of stable force chain networks (Sharkawi et al., 2018). By effectively bonding fractured rock fragments at a mesoscopic level, these technologies provide a fundamental basis for maintaining structural integrity under extreme geological conditions (Liang et al., 2026). Through numerical simulation, the deformation and failure mechanisms of the roadway surrounding rock corresponding to aforementioned theories and techniques were revealed, thereby providing optimized solutions and countermeasures (Meng et al., 2016; Du et al., 2016; Xin et al., 2024). Obviously, the continuous medium method effectively elucidates the redistribution of rock stress and the evolution patterns of plastic zones (Alber et al., 2009; Guo et al., 2021). However, sedimentary rocks are inherently a nonhomogeneous body, particularly stratified geomaterials containing widespread weak interlayers like mudstone, which further exacerbates the spatial heterogeneity of coal-rock mechanical behavior (Zhu et al., 2018; Kang et al., 2023). Such weak-surface structures alter stress transmission pathways and failure modes, causing support systems to gradually lose anchoring capacity under interfacial shear and interlayer slippage (Wang et al., 2018; Wu et al., 2025b; Yu et al., 2022; Zhu et al., 2025). Continuum theory, constrained by the continuum assumption and the principle of scale separation, struggles to characterize microprocesses such as interface sliding in weak interlayers, interlayer debonding, and shear failure at anchor cables and anchoring interfaces (Bobet and Einstein, 2011; Kang, 2014; Małkowski, 2015). While physical experiments can simulate these processes through similar materials, they are limited by scale effects and monitoring capabilities, which makes it difficult to capture dynamic evolution patterns such as crack initiation and force chain reconstruction within rock masses (Ding et al., 2017; Wang et al., 2023).

To examine the micro-mechanism of deformation and failure in anchored rock masses under shear slip at interlayer interfaces, this research thus takes the stability control of roadway surrounding rock at Kouzidong Coal Mine—where the mining depth reaches 1,000 m—as its engineering practice foundation. To analyze how the dip angle and spacing of weak interlayers affect the displacement field, crack field, force chain network, and stress evolution of surrounding rock, a discrete element model of deep roadway surrounding rock under weak interlayer impacts was established (Zhang et al., 2024; Chen et al., 2025; Peng et al., 2025). A discrete element model was established to analyze how the dip angle and spacing of weak interlayers govern the displacement field, crack field, and force chain network (Tang et al., 2021; Man et al., 2022). By developing a “three media-four interfaces” mesoscale contact model, this study successfully simulated the shear slip processes at weak interlayers and anchor interfaces, revealing the unstable patterns of roadway surrounding rock resulting from multi-stage shear failure of anchor cables triggered by weak interlayers (Wang et al., 2021; Wang et al., 2025). Given these limitations of continuum theory in capturing interface sliding, interlayer debonding, and localised shear failure, a discrete element method approach is therefore adopted. This method explicitly models particle scale interactions and fracture propagation to simulate the deformation and failure processes of the deep roadway surrounding rock containing weak interlayers.

2 Field engineering overview and discrete element model

2.1 Field engineering overview

The Kouzidong Coal Mine, a typical deep soft rock mine, serves as the physical prototype for this study. Its coal-bearing strata belong to the Carboniferous-Permian System of the North China Basin. As a typical kilometer-deep soft rock mine, it features high ground stress, high ground temperature, loose rock properties, and complex geological conditions (Kang et al., 2021). This study takes the mine’s transportation roadway as the research carrier, with its specific engineering conditions serving as the foundation. The roof and floor of the roadway consist of interbedded mudstone and sandy mudstone, with clay minerals such as illite and kaolinite exceeding 50% in content, while quartz accounts for only 39% of the composition. The mudstone exhibits distinct layering and softens readily upon contact with water. The compressive strength of the alternating layers of mud and sand in the roadway roof and floor rock ranges from 25 to 37 MPa, while the cohesive strength ranges from 8 to 12 MPa. In-situ stress testing demonstrates that the horizontal stress in the roof rock reaches 21.8 MPa, and the vertical stress reaches 25.1 MPa.

The roadway adopts a cross-sectional design with straight walls and an arch shape, measuring 5,800 mm in width and 4,100 mm in height, while its support scheme is presented in Figure 1a. The anchor bolts adopted in the scheme was ϕ22 mm with a length of 2.5 m, spaced at 0.8 × 0.8 m intervals, with a preload of 260 Nm; and ϕ21.8 mm anchor cables with lengths of 9.2 m (for the roadway roof) and 6.2 m (for the roadway ribs), spaced at 1.2 × 1.4 m intervals, with a preload of 160 kN. In situ borehole endoscopic observation provides a visual record of the internal roof rock strata as shown in Figure 1b. A digital borehole imaging system was utilized to detect the microstructural integrity and the distribution of fractures within the mudstone layers. The imaging results reveal clear evidence of structural damage within the rock mass. Specifically, the images captured in Figure 1b show the presence of shear fractures and composite fractures which confirm the intense fragmentation of the weak interlayers prior to the application of reinforcement.

FIGURE 1

Table 1 provides a detailed illustration of the lithological distribution of the coal seam’s roof and floor strata. The roof strata consist primarily of interbedded weak mudstones. The upper 37.1 m interval is dominated by mudstone or sandy mudstone, with an 8.2 m thick layer of fine sandstone occurring between 37.1 m and 45.3 m. The floor strata also consist mainly of mudstone and sandy mudstone, with cumulative mudstone thicknesses totaling 27.1 m. Test results for the physical and mechanical properties of the coal body are detailed in Table 2.

The presence of weak mudstone interbedding in the roof and floor strata led to severe deformation and failure of the roadway surrounding rock, despite the implementation of the aforementioned support measures. This was particularly evident during the mining phase, where significant nonlinear deformation characteristics were observed in the roadway walls and floor. The deformation and failure of the on-site roadway’s surrounding rock are illustrated in Figure 2. The coal pillar rib exhibits severe fracturing, with numerous roadway bolts and cables broken, steel bands torn, and support plates overturned. During the whole process of roadway excavation and longwall mining operations, the roadway experienced nearly 10 cumulative bottom lifts, with a rib clearance depth exceeding 3 m, a cumulative bottom bulge exceeding 6 m, and a cumulative rib shrinkage exceeding 4 m.

FIGURE 2

2.2 DEM methodology and mesoscale contact models

With the large-deformation haulage roadway in this mine as the engineering background, this paper investigates the influence of sand-mud interbedding on the stability of the roadway surrounding rock. We set weak interbedding inclinations (15°, 30°, 45°, 60°) and spacings (0.4, 0.8, 1.2, 1.6 m) to establish a two-dimensional discrete element model with a 30 m × 30 m grid. The rock mass within this numerical domain is represented by a distribution of particles with radii ranging from 0.03 m to 0.05 m. This specific size range was chosen to ensure that the thin weak interlayers, which are critical to the stability of the 1,000–m-deep roadway, are discretized by an adequate number of particles to accurately characterize their mechanical behavior and interfacial shear slip. Furthermore, maintaining the particle radii between 0.03 m and 0.05 m provides the high mesoscopic resolution required to capture the evolution of force chain networks and displacement fields while ensuring an optimal balance with overall computational efficiency. Regarding the boundary conditions during the excavation phase, a specific control strategy was implemented to ensure the validity of the unloading simulation. After the model reached the initial geostress equilibrium, constant stress boundaries were maintained on the top and right edges of the domain using a servo mechanism. This ensures that the far-field geostress of the 1,000–m-deep environment remains constant as the roadway is excavated. At the same time, the bottom and left edges were configured as displacement restricting boundaries to simulate the structural constraints of the deep rock mass. This combined boundary configuration allows for an accurate representation of the stress redistribution and the subsequent deformation of the surrounding rock during the excavation and unloading process. The roadway dimensions and anchor bolts/cables layout were configured according to actual dimensions, as illustrated in Figure 3, and the medium of anchor bolts/cables, sandstone, mudstone, and their corresponding anchor interfaces were considered. Calibration was performed according to the micromechanical parameters illustrated in Table 3. The anchor bolts/cables were calibrated to the 335 MPa yield strength steel, while the sandstone and mudstone mechanical parameters were calibrated to field-measured strengths as shown in Table 2. The anchor interface was calibrated to the 60 MPa field strength of resin grout.

FIGURE 3

The detailed procedures for numerical simulation are as follows: First, we assigned values to the surrounding rock of the roadway based on the field-measured rock stress (horizontal stress of 21.8 MPa and vertical stress of 25.1 MPa). Subsequently, the excavation process of the model roadway was simulated, and the excavated surrounding rock was supported in accordance with the on-site support scheme. This approach simulates the evolution of the microstructural behavior of the roadway surrounding rock under excavation unloading effects. In the discrete element model, the positioning generation method is employed to construct the anchor bolts/cables medium. This ensures full contact between the anchor bolt/cable and the coal-rock mass while preventing bending, fracture, or even particle redistribution of the anchor bolt/cable during the model equilibrium process (Wu et al., 2025a).

3 Results and analysis

3.1 Characteristics of displacement deformation

Figure 4 illustrates how the dip angle and spacing of weak interlayers affect the evolution of the displacement field in the roadway’s surrounding rock, providing a clear visual reference for subsequent analysis. As shown in the figure, deformation occurs due to rock mass unloading, with the roadway exhibiting distinct displacement characteristics in both the roof and floor as well as the ribs. If active support can effectively resist rock mass deformation during the geostress equilibrium process, the roadway will not become unstable. It is evident that under this support plan, the surrounding rock of the roadway is unstable, displaying characteristics such as roof falls, floor heaving, and wall convergence. To rigorously validate the accuracy of the model, we performed a quantitative comparison between the numerical results and the field monitoring data recorded at the 1,000–m-deep site. The simulation results are almost identical to the field detection records, where the cumulative floor heaving was measured at 6.2 m and the rib shrinkage was measured at 4.1 m prior to implementation of any reinforcement. However, the dip angle and spacing of the weak interlayers significantly influence these displacement deformation patterns. When the dip angle of the weak interlayer grows, the displacement deformation of the roof and floor increases, while there is a relative decrease in the displacement deformation of the ribs. As the spacing of the weak interlayer increases, the displacement deformation of both the roof and floor and the ribs decreases remarkably.

FIGURE 4

3.2 Distribution evolution and numerical characteristics of cracks

How the dip angle and spacing of weak interlayers influence the evolution of crack fields in the surrounding rock of the roadway is illustrated in Figure 5, effectively reflecting key findings on crack propagation characteristics. As illustrated, cracks invariably initiate at the surface layer of the surrounding rock and the unsupported floor, subsequently propagating deeper. With increasing dip angle of the weak interlayer, the height at which roof cracks develop significantly increases, becoming more densely distributed, while corresponding crack development in the ribs diminishes. As the spacing of weak interlayers increases, the roof, floor and ribs decrease markedly. Notably, the weak interlayers, which are more prone to producing shear cracks, will break the anchor bolts and cables at a certain angle due to the shear slip along the weak interlayers during rock mass unloading. Increased dip angles of weak interlayers lead to higher anchor cable failure rates. Particularly in surrounding rock of the roadway with 60° dip angles, multiple-segment shear failures occur in roof anchors and cables, causing the roof to collapse as a unit. As the spacing between weak interlayers increases, the number of the interlayers decreases, the broken anchor rods and cables under the shear action of the weak interlayers also decrease accordingly. Typically, when the spacing between weak interlayers is 0.4 m, broken anchor bolts and cables appear in both the roadway roof and ribs. When the spacing between weak interlayers is 1.6 m, broken anchor bolts and cables only occur in the roadway ribs.

FIGURE 5

For the purpose of further investigating the evolution patterns of crack fields within the roadway’s surrounding rock, the influence of weak interlayer dip angle and spacing on crack density is illustrated in Figure 6. The radial direction represents crack density, while the circumferential direction indicates crack orientation. The Fig. reveals that as the weak interlayer dip angle increases, crack density in all directions within the surrounding rock of roadway significantly rises. With an increase in the spacing between weak interlayers, the number of cracks in all directions inside the roadway surrounding rock drops significantly, which aligns with the observed crack evolution trends in the study. Regarding the orientation of cracks, they are mainly in the vertical direction—this makes the roadway roof more likely to spall and collapse, and this phenomenon matches the evolution pattern of the displacement field.

FIGURE 6

Quantifying the damage level of surrounding rock is a prerequisite for predicting and assessing roadway instability. This paper constructs a fractal model to characterize structural damage based on crack distribution patterns throughout the surrounding rock of roadway loading process. The fractal dimension that reflects crack distribution characteristics is calculated by means of the box-counting method. The specific steps are as follows: First, the original crack distribution illustrated in Figure 7a is converted into a binary graph shown in Figure 7b using an adaptive thresholding process. To ensure methodological transparency and eliminate subjective bias, the adaptive gray level threshold is determined using Otsu’s method. This algorithm automatically identifies the optimal separation point by maximizing the inter class variance between the pixels representing cracks and those representing the rock matrix. By employing this adaptive approach, we ensured that the binary representation of the fracture network accurately captures the structural damage regardless of variations in local contrast. Within this binary graph, the unit box assignment rule is defined as follows: if a unit box contains at least one crack pixel, it is assigned a value of 1; otherwise, it is assigned a value of 0. Subsequently, the fractal dimension for crack distribution is determined by calculating the slope of the linear fit (as shown in Equation 1) between the total number of boxes and the number of non-empty boxes in a double-logarithmic coordinate system (Figure 7c).

FIGURE 7

In the formula: is the fractal dimension for crack distribution characterizing the degree of structural damage; represents the total count of boxes; represents the count of non-empty boxes.

The above calculation results indicate that, as the primary load-bearing structure of the surrounding rock of roadway, the deep rock mass consistently exhibits a high degree of crack development. To further quantify the evolutionary patterns of fractal characteristics in the fractal dimension of the crack distribution of surrounding rock of roadway, the fractal dimension of crack distribution in surrounding rock of roadway under different geological conditions was calculated using the box-counting method. This yielded evolutionary curves of fractal dimension for crack distribution in anchored roadway rock mass under varying geostress, as illustrated in Figure 8.

FIGURE 8

As shown in Figure 8, the fractal dimension that characterizes rock mass damage gradually increases with the progression of unloading deformation failure in the surrounding rock of roadway. Under the same stress conditions, the fractal dimension also increases as the dip angle of the weak interlayer increases. The fractal dimension decreases as the spacing between weak interlayers increases. However, as illustrated in Figure 8b, the fractal dimension that corresponds to the spacing of 0.4 m is smaller than the fractal dimension that corresponds to the spacing of 1.2 m. It is possible that the stress in the hard rock transmits along relatively continuous paths, resulting in the formation of a large-scale interconnected fracture network in the stress concentration zones, which increases the fractal dimension.

3.3 Force chain distribution and stress evolution

The force chains in the discrete element model are used to characterize the mesoscopic interactions and load-bearing pathways between particles. The force chains which distribute densely indicate that the structure is compact and undamaged, while the force chains with large values characterize the main bearing zone of the structure. The force chain field of the surrounding rock mass can represent the internal stability of the surrounding rock and the anchoring effect of the anchor bolts and cables, however, it is necessary to note that when the values are overly high, it is easy to exceed its ultimate force chain strength, leading to the redistribution of the force chain network. Figure 9 illustrates the influence of weak interlayer dip angle and spacing on the evolution of the force chain field of the surrounding rock of the roadway. The Figure 9 reveals that the surface rock mass exhibits a discontinuous force chain distribution due to unloading deformation failure, while the deep rock mass force chain network remains relatively dense with minimal numerical variation. As the dip angle of the weak interlayer increases, the fracture force chains in the roof and floor expand, while those in the two ribs relatively contract, consistent with the displacement field and crack field. With increasing spacing of the weak interlayer, the fracture force chains in the roof, floor and two ribs all present a shrinking trend. Notably, the high strength chain network often manifests within the anchor bolt anchorage zone. However, driven by the shear slip of the weak interlayer and the shear failure of the anchor bolts, the force chain network subsequently fractures and weakens. It is particularly evident in the surrounding rock of roadway with a weak interlayer dip angle of 60° illustrated in Figure 9, where a distinct force chain fault formed in the roof, leading to the overall collapse of the roof rock mass.

FIGURE 9

To further quantify the stress evolution patterns within the surrounding rock of roadway, Figures 10, 11 respectively illustrate the horizontal and vertical stress evolution patterns at the ribs, shoulders, and roof regions of roadways with weak interlayers at inclinations of 15° and 60°, and spacings of 0.4 m and 1.6 m. As illustrated in Figures 10, 11, both horizontal and vertical stresses in the surrounding rock of roadway drop sharply due to unloading deformation. The timing of stress drop-off is delayed as one moves deeper into the rock mass. For deeper rock masses with minimal deformation, stresses remain largely unchanged. Notably, some stress curves exhibit a subsequent drop after exceeding deeper monitoring points. This phenomenon is attributed to the presence of anchor bolts and cables within these monitoring points, which enhance their load-bearing capacity. This subsequent stress drop is fundamentally driven by localized shear interactions at the anchorage interface. As the weak interlayer undergoes interfacial shear slip, it exerts a concentrated lateral force on the anchor cable crossing the shear plane, inducing high localized shear stress. When this localized stress exceeds the ultimate shear strength of the anchoring system, the anchor cable undergoes sudden shear failure which disrupts the mechanical coupling between the support structure and the surrounding rock. This failure results in a subsequent stress drop as documented in Figures 10aIII,bIII, 11aIII,bIII. Therefore, the stress fluctuations in measurement points containing anchor cable media are significantly greater than those in pure rock mass. This can be interpreted as the effect of high-tensile anchor cable media and their anchoring action on the rock mass (e.g., Figures 10aVI,bVI, 11aVI,bVI). Previously, it was believed that surface rock stress in roadways typically falls below that of deeper rock masses. However, both Figures 10, 11 reveal the stress in the surface rock exceeding that of the deeper rock mass. This phenomenon is attributed to the anchoring effect of the anchor cables resisting the deformation and compression of the rock mass into the roadway. As the deformation and failure of the rock mass in the anchored zone progresses, its stress gradually declines, while the stress in the deeper rock mass remains largely unchanged. From the comparison between Figures 10a,b, it is evident that an increase in the dip angle of the weak interlayer leads to a greater magnitude of rock stress decline and an earlier onset time, corresponding to enhanced stress fluctuations at the anchor rod-equipped measurement points. From the comparison of Figures 11a,b, it is evident that an increase in the spacing of weak interlayers leads to a decrease in the magnitude of rock stress decline and a delay in the timing of decline. This corresponds to a significant weakening of stress fluctuations at monitoring points containing anchor cables.

FIGURE 10

FIGURE 11

Figure 12 illustrates the influence of weak interlayer dip angle and spacing on the number of internal force chains in the surrounding rock of roadway. As the dip angle of weak interlayers increases, the number of force chains decreases significantly, indicating that a large number of force chains broke. As the spacing between weak interlayers increases, the number of force chains also increases. Notably, the surrounding rock of roadway exhibits a marked reduction in the number of force chains due to unloading deformation failure, which is consistent with the evolution pattern of the force chain field in roadway surrounding rock.

FIGURE 12

3.4 Particle failure patterns

Figure 13 illustrates the influence of weak interlayer dip angle and spacing on the failure patterns of the surrounding rock of roadway. As illustrated, the failure patterns of the surrounding rock are consistent with the fracture patterns of the force chain network. As the dip angle of the weak interlayer increases, the particle failure in the roof and floor rock mass intensifies, while the particle failure in the rib rock mass relatively decreases. Meanwhile, particle shear failure traces along the weak interlayer can be observed. The increasing spacing of the weak interlayers leads to a thicker layer of hard rock, resulting in a certain degree of reduction in particle failure in both the roof and floor rock masses as well as the rib rock masses. Since there are no surface support components implemented in the simulation, such as metal mesh or steel bands, the particle failure directly results in rock fall. Both the increasing dip angle of the weak interlayers and the decreasing spacing lead to intensified rock fall.

FIGURE 13

4 Discussion

Figure 14 further provides a microscopic analysis of surrounding rock deformation and failure under the influence of weak interlayers, comparing cases where interlayer inclinations are 15° and 60°, and spacing is 0.4 m and 1.6 m. The force chain vectors of the anchor cables indicate that the surrounding rock, which has been deformed by compression toward the roadway, exerts tensile stress on the cables. Due to the absence of protective support components such as metal mesh and steel bands, the force chain network of the roadway’s surrounding rock has fractured over a large area, failing to provide surface protection. This failure triggers roof falls, rib falls, and floor heaving disasters. As the dip angle of the weak interlayer increases, the shear slip angle at the interface becomes larger; this increment further leads to higher shear stress and deformation at the interface of the roadway roof, which is consistent with the study’s mechanical inference. It is easier to cause shear failure of the anchor cables at the roadway roof (Figure 14a). The decrease in the spacing of the weak interlayer results in a growth in the quantity of shear interfaces, resulting in more shear failures of the anchor cables (Figure 14b). In conclusion, there are two primary factors that cause the instability of the surrounding rock of the roadway under the influence of weak interlayers: First, the large-scale failure of the surface rock force chain network prevents effective surface support. Second, the insufficient shear resistance of the anchor cables leads to anchorage failure. Therefore, for roadway surrounding rock stability control under such geological conditions, the focus should be on enhancing the shear resistance of the support components and strengthening the force chain connection between the surface rock and the deep rock mass.

FIGURE 14

4.1 Implementation of grouting reinforcement

Based on the analysis of the instability mechanisms driven by deep weak interlayers, two core shortcomings in the original support scheme are identified: First, the absence of surface reinforcement components fails to constrain the fractured rock, thereby inducing severe macroscopic deformation, such as roof and rib falls. Second, the support system lacks the capacity to control the shear slip at interfaces induced by weak interlayers. The insufficient shear resistance of anchor cables makes it difficult to withstand these interlayer shear forces, leading to multi-segment shear failure.

To address this instability, an optimized control scheme combining boundary confinement and grouting reinforcement is proposed: (1) The specific scheme is shown in Figure 15. Install protective surface components, such as metal mesh and steel bands, to distribute support stress and enhance structural integrity. (2) Introduce grouting slurry within a 6.5–m radius centered on the roadway arch to permeate the dynamic fracture network. This leverages the bonding effect of the grouting slurry to reconstruct the force chain network of the rock mass. The injected fluid permeates the anisotropic fracture channels governed by pressure-driven seepage. We assume that the primary role of the injected slurry is to provide mechanical bonding between the fractured rock fragments to facilitate the reconstruction of the internal support network. The numerical implementation of this process within the discrete element framework was achieved by adjusting the mechanical properties of the existing particle contacts rather than generating new particles to represent the slurry. Specifically, we modified the normal and tangential bond strengths and the stiffness of the particle interfaces within the 6.5 m grouting radius to simulate the reinforcement effect of the solidified cementing fluid. By increasing these contact parameters, the solidified slurry effectively glues the highly fractured rock matrix back together. This grouting intervention halts the progressive macroscopic deformation and fundamentally forms a continuous and high strength force chain network. To ensure the analysis remains centered on the mechanical reinforcement and structural integrity of the rock mass, complex fluid interactions and specialized rheological behaviors are excluded from the current study. As the fluid solidifies during the curing process, it effectively bonds the highly fractured rock matrix. This grouting intervention halts the progressive macroscopic deformation, significantly optimizes the fractal dimension of crack distribution within the surrounding rock, and fundamentally forms a continuous, high-strength force chain network.

FIGURE 15

4.2 Validation of force chain reconstruction via cementing slurry injection

To simulate the displacement field, crack field, force chain field, and particle failure patterns of the roadway surrounding rock in the excavation and unloading process, and to confirm the effectiveness of the support, the grouting intervention was introduced after the model roadway excavation and the subsequent unloading process reached a state of preliminary equilibrium. By applying the cementing slurry at this stage, the simulation effectively accounts for the reinforcement of a rock mass that has already undergone initial deformation and microstructural damage. This ensures that the slurry permeates the existing fracture network to perform its bonding and reconstruction functions in a realistic geological environment. A two-dimensional discrete element model was constructed according to the optimized support scheme (Figure 16a).

FIGURE 16

Numerical simulation results indicate that after support optimization, deformation and failure in the roadway ribs were significantly reduced. The maximum displacement of the roadway roof ranged between 0.3 and 0.4 m, resulting in effective resolution of roof falls. Since the roadway floor was also reinforced with grouting, it remains in a stable state (Figure 16b).

Due to the significant reduction in roadway deformation, the propagation depth of shallow tensile cracks and deep shear cracks in the roadway roof and ribs decreased markedly (Figure 16c), and the overall number of cracks also decreased significantly (Figure 16d). For the weak interlayers, which are more prone to developing shear cracks, anchor cables oriented at certain angles to these interlayers exhibited shear failure. The number of cracks around the failed anchors was markedly higher than in other locations. This indicates that grouting effectively filled the fissures at the interfaces of the weak interlayers, successfully preventing crack propagation into the deeper rock mass.

A continuous high-strength force chain network formed within the anchor cable anchorage zone (Figure 16e), with a clear force chain transmission path and no significant force chain discontinuities. Noticeable stress concentration phenomena occurred around the steel bands, particularly at the roof, consistent with the displacement field simulation results.

The failure patterns of the surrounding rock largely correspond with the fracture patterns in the force chain network (Figure 16f). There is only some obvious granular failure in localized areas of the roof. This indicates that when the rock bolts are performing their anchoring function, they transfer stress to the surrounding particles through their own deformation under load. Consequently, some particles fail due to bearing excessive stress. However, the extent of this failure is limited and remains under control. This demonstrates that the support system has effectively held back the propagation of granular failure to some degree, which aligns with the expected control effect of the optimized scheme.

4.3 Comparison with existing literature and research contributions

Stability control of deep roadways with weak interlayers remains a critical underground engineering challenge. While extensive research exists, key gaps remain in mesoscopic instability mechanisms and fundamental grouting reinforcement principles.

Most existing studies rely on continuum mechanics frameworks, which effectively predict macroscopic stress redistribution but cannot explicitly characterize discontinuous processes like weak interlayer sliding and anchor-rock interface shear failure (Moussaei et al., 2019; Jia et al., 2025). While some researchers have adopted the Discrete Element Method (DEM), most studies are limited to single weak interlayer scenarios or simplified support systems, failing to quantify the coupled effects of interlayer dip angle and spacing on full-field evolution (Gao et al., 2025). Additionally, current grouting research focuses primarily on material development, with insufficient attention to the mesoscopic mechanism of force chain network reconstruction (Zhijun et al., 2025).

Building on these prior findings, our comprehensive discrete element analysis addresses these gaps by quantifying the coupled influence of weak interlayer dip angle and spacing on roadway stability, showing that steeper interlayers exacerbate force chain fracture while larger spacing limits maximum roof displacement to 0.3–0.4 m. Complementing this, we establish a fractal dimension-based damage assessment framework that more accurately captures fracture network complexity than traditional methods, and reveal that grouting’s core role is reconstructing continuous high-strength force chain networks rather than merely increasing macroscopic rock strength, deepening understanding of grouting principles.

Methodologically, our DEM model explicitly simulates crack initiation, propagation, and interface shear slip, more accurately reproducing progressive failure processes than continuum-based approaches. Furthermore, our “three media-four interfaces” mesoscale contact model precisely characterizes multi-medium interactions, improving simulation reliability compared with existing DEM studies.

Despite these advances, our two-dimensional model cannot fully capture three-dimensional interlayer distribution and anchor cable stress states. Additionally, our grouting simulation only considers solidified slurry mechanics without flow processes, and we do not address time-dependent weak interlayer properties or water-rock interactions, all of which will be incorporated in future three-dimensional coupled CFD-DEM studies.

5 Conclusion

References

• 1

  AlberM.FritschenR.BischoffM.MeierT. (2009). Rock mechanical investigations of seismic events in a deep longwall coal mine. Int. J. Rock Mech. Min. Sci.46, 408–420. 10.1016/j.ijrmms.2008.07.014

  • CrossRef
  • Google Scholar
  • View reference in article

• 2

  BobetA.EinsteinH. H. (2011). Tunnel reinforcement with rockbolts. Tunn. Undergr. Space Technol.26, 100–123. 10.1016/j.tust.2010.06.006

  • CrossRef
  • Google Scholar
  • View reference in article

• 3

  ChenQ.QiaoH.YangT.SongM.RenT.WangG.et al (2025). Study on shear performance between layers of rock-filled concrete based on morphology analysis. Constr. Build. Mater.493, 143151. 10.1016/j.conbuildmat.2025.143151

  • CrossRef
  • Google Scholar
  • View reference in article

• 4

  DingS.JingH.ChenK.XuG.MengB. (2017). Stress evolution and support mechanism of a bolt anchored in a rock mass with a weak interlayer. Int. J. Min. Sci. Technol.27, 573–580. 10.1016/j.ijmst.2017.03.024

  • CrossRef
  • Google Scholar
  • View reference in article

• 5

  DuZ.QinB.TianF. (2016). Numerical analysis of the effects of rock bolts on stress redistribution around a roadway. Int. J. Min. Sci. Technol.26, 975–980. 10.1016/j.ijmst.2016.09.003

  • CrossRef
  • Google Scholar
  • View reference in article

• 6

  GaoJ.GuoZ.ZhuJ. (2025). Model test on the ultimate failure mechanism of automatically formed roadway in weakly cemented soft rock. Tunn. Undergr. Space Technol.166, 106968. 10.1016/j.tust.2025.106968

  • CrossRef
  • Google Scholar
  • View reference in article

• 7

  GuoX.LiC.HuoT. (2021). Shapes and formation mechanism of the plastic zone surrounding circular roadway under partial confining stress in deep mining. Geomechanics Eng.25, 509–520. 10.12989/GAE.2021.25.6.509

  • CrossRef
  • Google Scholar

• 8

  HeM.GuoA.MengZ.PanY.TaoZ. (2023). Impact and explosion resistance of NPR anchor cable: field test and numerical simulation. Undergr. Space10, 76–90. 10.1016/j.undsp.2022.10.001

  • CrossRef
  • Google Scholar
  • View reference in article

• 9

  JiaC.LiangG.ShiC.YuJ. (2025). Investigation on three-dimensional large deformation characteristics of tunnel in stratified rock mass with coupled discrete–continuous method. Tunn. Undergr. Space Technol.164, 106777. 10.1016/j.tust.2025.106777

  • CrossRef
  • Google Scholar
  • View reference in article

• 10

  JiaoY.-Y.SongL.WangX.-Z.Coffi AdokoA. (2013). Improvement of the U-shaped steel sets for supporting the roadways in loose thick coal seam. Int. J. Rock Mech. Min. Sci.60, 19–25. 10.1016/j.ijrmms.2012.12.038

  • CrossRef
  • Google Scholar
  • View reference in article

• 11

  JingH.WuJ.YinQ.WangK. (2020). Deformation and failure characteristics of Anchorage structure of surrounding rock in deep roadway. Int. J. Min. Sci. Technol.30, 593–604. 10.1016/j.ijmst.2020.06.003

  • CrossRef
  • Google Scholar
  • View reference in article

• 12

  KangH. (2014). Support technologies for deep and complex roadways in underground coal mines: a review. Int. J. Coal Sci. Technol.1, 261–277. 10.1007/s40789-014-0043-0

  • CrossRef
  • Google Scholar
  • View reference in article

• 13

  KangH.YangJ.MengX. (2015). Tests and analysis of mechanical behaviours of rock bolt components for China’s coal mine roadways. J. Rock Mech. Geotechnical Eng.7, 14–26. 10.1016/j.jrmge.2014.12.002

  • CrossRef
  • Google Scholar
  • View reference in article

• 14

  KangH.JiangP.WuY.GaoF. (2021). A combined “ground support-rock modification-destressing” strategy for 1000-m deep roadways in extreme squeezing ground condition. Int. J. Rock Mech. Min. Sci.142, 104746. 10.1016/j.ijrmms.2021.104746

  • CrossRef
  • Google Scholar

• 15

  KangH.GaoF.XuG.RenH. (2023). Mechanical behaviors of coal measures and ground control technologies for China’s deep coal mines – a review. J. Rock Mech. Geotechnical Eng.15, 37–65. 10.1016/j.jrmge.2022.11.004

  • CrossRef
  • Google Scholar

• 16

  KhaleghparastS.AzizN.RemennikovA.AnzanpourS. (2023). An experimental study on shear behaviour of fully grouted rock bolt under static and dynamic loading conditions. Tunn. Undergr. Space Technol.132, 104915. 10.1016/j.tust.2022.104915

  • CrossRef
  • Google Scholar

• 17

  LiangX.HeC.FengK.WangC.HuZ.ZhangJ.et al (2026). Grouting effectiveness in shield tunnels through water-rich sandy pebble strata and its impact on structure. Tunn. Undergr. Space Technol.175, 107785. 10.1016/j.tust.2026.107785

  • CrossRef
  • Google Scholar

• 18

  LiuC.SuX.WuY.ZhengZ.YangB.LuoY.et al (2021). Effect of nano-silica as cementitious materials-reducing admixtures on the workability, mechanical properties and durability of concrete. Nanotechnol. Rev.10, 1395–1409. 10.1515/ntrev-2021-0097

  • CrossRef
  • Google Scholar

• 19

  MaB.XieH.ZhangX.ZhouH.ZhouC.SunW.et al (2025). Experimental investigation of rock mineralogical effect on energy transfer and rockbursts induced by tensile fracturing of roof strata. Int. J. Rock Mech. Min. Sci.189, 106087. 10.1016/j.ijrmms.2025.106087

  • CrossRef
  • Google Scholar

• 20

  MałkowskiP. (2015). The impact of the physical model selection and rock mass stratification on the results of numerical calculations of the state of rock mass deformation around the roadways. Tunn. Undergr. Space Technol.50, 365–375. 10.1016/j.tust.2015.08.004

  • CrossRef
  • Google Scholar

• 21

  ManJ.HuangH.AiZ.ChenJ. (2022). Analytical model for tunnel face stability in longitudinally inclined layered rock masses with weak interlayer. Comput. Geotechnics143, 104608. 10.1016/j.compgeo.2021.104608

  • CrossRef
  • Google Scholar

• 22

  MengQ.HanL.XiaoY.LiH.WenS.ZhangJ. (2016). Numerical simulation study of the failure evolution process and failure mode of surrounding rock in deep soft rock roadways. Int. J. Min. Sci. Technol.26, 209–221. 10.1016/j.ijmst.2015.12.006

  • CrossRef
  • Google Scholar

• 23

  MoussaeiN.SharifzadehM.SahriarK.KhosraviM. H. (2019). A new classification of failure mechanisms at tunnels in stratified rock masses through physical and numerical modeling. Tunn. Undergr. Space Technol.91, 103017. 10.1016/j.tust.2019.103017

  • CrossRef
  • Google Scholar

• 24

  OnaiziA. M.HuseienG. F.LimNHASAmranM.SamadiM. (2021). Effect of nanomaterials inclusion on sustainability of cement-based concretes: a comprehensive review. Constr. Build. Mater.306, 124850. 10.1016/j.conbuildmat.2021.124850

  • CrossRef
  • Google Scholar

• 25

  PengM.GuoX.XiaC.ShiZ. (2025). Stability analysis of a rock slope-pile-anchor coupled reinforcement system with bedding weak zones using an improved SPH method. Eng. Geol.357, 108305. 10.1016/j.enggeo.2025.108305

  • CrossRef
  • Google Scholar

• 26

  RenF.SongT.MaK.KarakusM. (2024). Experimental investigation on the influence of weak interlayers on sandstone rockburst and associated microcracking mechanism. Int. J. Rock Mech. Min. Sci.182, 105890. 10.1016/j.ijrmms.2024.105890

  • CrossRef
  • Google Scholar

• 27

  SaxenaK. K.DasR.CaliusE. P. (2016). Three decades of auxetics research − materials with negative poisson’s ratio: a review. Adv. Eng. Mater18, 1847–1870. 10.1002/adem.201600053

  • CrossRef
  • Google Scholar

• 28

  SharkawiA. M.Abd-ElatyM. A.KhalifaO. H. (2018). Synergistic influence of micro-nano silica mixture on durability performance of cementious materials. Constr. Build. Mater.164, 579–588. 10.1016/j.conbuildmat.2018.01.013

  • CrossRef
  • Google Scholar

• 29

  TangB.YeboahM.ChengH.TangY.YaoZ.WangC.et al (2021). Numerical study and field performance of rockbolt support schemes in TBM-excavated coal mine roadways: a case study. Tunn. Undergr. Space Technol.115, 104053. 10.1016/j.tust.2021.104053

  • CrossRef
  • Google Scholar

• 30

  WangQ.JiangB.PanR.LiS. C.HeM. C.SunH. B.et al (2018). Failure mechanism of surrounding rock with high stress and confined concrete support system. Int. J. Rock Mech. Min. Sci.102, 89–100. 10.1016/j.ijrmms.2018.01.020

  • CrossRef
  • Google Scholar

• 31

  WangF.XiaK.YaoW.WangS.WangC.XiuZ. (2021). Slip behavior of rough rock discontinuity under high velocity impact: experiments and models. Int. J. Rock Mech. Min. Sci.144, 104831. 10.1016/j.ijrmms.2021.104831

  • CrossRef
  • Google Scholar

• 32

  WangJ.LiuP.WuC.HeM.GongW. (2023). Mechanical behavior of soft rock roadway reinforced with NPR cables: a physical model test and case study. Tunn. Undergr. Space Technol.138, 105203. 10.1016/j.tust.2023.105203

  • CrossRef
  • Google Scholar

• 33

  WangX.LiY.ZhaoG.ZhangQ.FuQ.FanC.et al (2025). Study on mechanical behavior and reinforcement mechanism of fissured soft rock reinforced by anchor-grouting integration. Constr. Build. Mater.493, 143142. 10.1016/j.conbuildmat.2025.143142

  • CrossRef
  • Google Scholar

• 34

  WuJ.WongH.ZhangH.YinQ.JingH.MaD. (2024). Improvement of cemented rockfill by premixing low-alkalinity activator and fly ash for recycling gangue and partially replacing cement. Cem. Concr. Compos.145, 105345. 10.1016/j.cemconcomp.2023.105345

  • CrossRef
  • Google Scholar
  • View reference in article

• 35

  WuJ.YangS.WilliamsonM.WongH. S.BhudiaT.PuH. (2022a). Microscopic mechanism of cellulose nanofibers modified cemented gangue backfill materials. Adv Compos Hybrid Mater8, 177. 10.1007/s42114-025-01270-9

  • CrossRef
  • Google Scholar

• 36

  WuJ.ZhangW.WangY.JuF.PuH.RiabokonR.et al (2022b). Effect of composite alkali activator proportion on macroscopic and microscopic properties of gangue cemented rockfill: Experiments and molecular dynamic modelling. Int. J. Miner Metall. Mater32, 1813–1825. 10.1007/s12613-025-3140-8

  • CrossRef
  • Google Scholar

• 37

  XieH.LiC.HeZ.LuY.ZhangR.et al (2021). Experimental study on rock mechanical behavior retaining the in situ geological conditions at different depths. Int. J. Rock Mech. Min. Sci.138, 104548. 10.1016/j.ijrmms.2020.104548

  • CrossRef
  • Google Scholar

• 38

  XinX.MengQ.PuH.WuJ. (2024). Theoretical analysis and numerical simulation analysis of energy distribution characteristics of surrounding rocks of roadways. Tunn. Undergr. Space Technol.147, 105747. 10.1016/j.tust.2024.105747

  • CrossRef
  • Google Scholar

• 39

  YangR.LiY.GuoD.YaoL.YangT.LiT. (2017). Failure mechanism and control technology of water-immersed roadway in high-stress and soft rock in a deep mine. Int. J. Min. Sci. Technol.27, 245–252. 10.1016/j.ijmst.2017.01.010

  • CrossRef
  • Google Scholar

• 40

  YuM.ZuoJ.SunY.MiC.LiZ. (2022). Investigation on fracture models and ground pressure distribution of thick hard rock strata including weak interlayer. Int. J. Min. Sci. Technol.32, 137–153. 10.1016/j.ijmst.2021.10.009

  • CrossRef
  • Google Scholar

• 41

  ZhangQ.JiangB.-S.WangS.GeX. r.ZhangH. q. (2012). Elasto-plastic analysis of a circular opening in strain-softening rock mass. Int. J. Rock Mech. Min. Sci.50, 38–46. 10.1016/j.ijrmms.2011.11.011

  • CrossRef
  • Google Scholar

• 42

  ZhangW.HuL.YaoZ.-B.XiongY. R.ZhaoJ.MaT.et al (2024). In-situ and experimental investigations of the failure characteristics of surrounding rock through granites with biotite interlayers in a tunnel. Eng. Geol.343, 107816. 10.1016/j.enggeo.2024.107816

  • CrossRef
  • Google Scholar

• 43

  ZhangB.ZhaoH.JiD.GaoL.XuX.WangH. (2025). Study on the rational layout and control of gas drainage roadways in deep inclined coal seams. Phys. Fluids37, 067143. 10.1063/5.0273236

  • CrossRef
  • Google Scholar

• 44

  ZhijunZ.BangyaoL.HuaimiaoZ.YakunT.LinglingW. (2025). Study on the microstructure and mechanical properties of uranium tailings reinforced by MICP under the action of activated carbon. Sep. Purif. Technol.367, 132947. 10.1016/j.seppur.2025.132947

  • CrossRef
  • Google Scholar

• 45

  ZhuW.JingH.YangL.PanB.SuH. (2018). Strength and deformation behaviors of bedded rock mass under bolt reinforcement. Int. J. Min. Sci. Technol.28, 593–599. 10.1016/j.ijmst.2018.03.006

  • CrossRef
  • Google Scholar

• 46

  ZhuQ.YinQ.TaoZ.HeM. (2025). Cyclic shear responses of rough-walled rock joints subjected to dynamic normal loads. J. Rock Mech. Geotechnical Eng.17, 3289–3297. 10.1016/j.jrmge.2024.11.028

  • CrossRef
  • Google Scholar

• 47

  ZuoJ.WangZ.ZhouH.PeiJ.LiuJ. (2013). Failure behavior of a rock-coal-rock combined body with a weak coal interlayer. Int. J. Min. Sci. Technol.23, 907–912. 10.1016/j.ijmst.2013.11.005

  • CrossRef
  • Google Scholar