Rainfall-induced seepage damage coupling mechanism and stability evaluation of high and steep slopes in deep open-pit mines

1. Introduction

The stability of open-pit mine slopes is a fundamental scientific and engineering necessity for ensuring safe mining operations. Global field investigations and case studies have demonstrated that rainfall was a primary external trigger for the open-pit slope instability. Particularly, high and steep slopes composed of water-sensitive and easily disintegrated soft rocks (e.g., schist) exhibit pronounced mechanical sensitivity to moisture conditions: rainfall infiltration reduces the effective stress and shear strength of potential sliding surfaces via elevated pore water pressure. It also accelerates microcrack propagation and macrostructural degradation of surface rock masses through induced wetting-drying cycles, forming a positive feedback loop of “seepage intensifying damage → damage promoting seepage”. However, existing studies were predominantly focused on the individual impacts of seepage or damage processes on the slope stability, while the dynamic coupling mechanism between these processes remains almost unclear, hindering accurate characterization of the progressive evolution path of slope instability under rainfall conditions. Therefore, it is crucial to reveal systematically the coupling mechanism between rainfall seepage and rock mass degradation for improving the accuracy of open-pit slope instability early warning systems and enhancing disaster prevention and mitigation capabilities.

Seepage-damage coupling theory serves as the core framework for analyzing the rainfall-induced instability mechanism of high and steep slopes in open-pit mines. Its essence lies in establishing the dynamic correlations among the seepage field, stress field, and damage field, as well as quantifying the mutual feedback effects of multi-field parameters. Research in this field has advanced from single-medium to complex rock masses and from steady-state to unsteady-state seepage, yielding a series of coupling models tailored to open-pit mine scenarios. To address the widespread fractured rock slopes in open-pit mines, Wang et al. [1] developed an unsaturated seepage-stress-damage three-field coupling model. They generated random fracture networks using MATLAB, combined with COMSOL simulations, and derived three-field coupling control equations under unsaturated conditions, and validated the model effectiveness through the PDE module. Reliability analysis clarified that hydrodynamic pressure and shear stress are the dominant factors in the coupling process, with cohesion exerting a more significant influence on slope stability than the internal friction angle. Thereby a methodological framework was provided for precise risk assessment of high and steep fractured rock slopes in open-pit mines. Xu et al. [2] extended the seepage-stress-damage coupling theory to long-term engineering performance prediction. After validating initial settlement data via finite element models, they successfully predicted the damage evolution law for slopes over 50 years of operation and proposed the conclusion that “early settlement is significant, whereas long-term damage is controllable”, offering a conceptual transfer for the long-term stability evaluation of high and steep slopes in open-pit mines. Regarding algorithm optimization, Tang et al. [3] proposed an improved hydraulic coupling algorithm of Discrete Element Method-Finite Difference Method (DEM-FDM), which effectively balances the efficiency and accuracy of unsteady seepage simulations and resolves the technical bottleneck of excessive computational time in transient response simulations during rainfall seepage of open-pit mine slopes. International scholars also contributed valuable extensions to coupling model construction: Sanei et al. [4] integrated the Darcy’s law and the cap elastoplastic model to establish an iteratively coupled reservoir-geomechanical model, revealing the irreversible effect of strain-dependent permeability on seepage. This model can be extended to seepage-damage coupling simulations of high-stress slopes in open-pit mines, improving the prediction accuracy of permeability changes under high-stress conditions. Cleto et al. [5] incorporated the Mesh Fragmentation Technique (MFT), High Aspect Ratio (HAR) elements, and cohesive interface elements into the Finite Element Method (FEM), establishing a simulation model for the interaction between hydraulic fractures and natural fractures that provides key methodological references for seepage evolution analysis of fractured slopes in open-pit mines.

Instability of high and steep slopes in open-pit mines often manifests as “local triggering → global instability”. The seepage-damage response of critical zones (e.g., soil-rock interfaces, weak interlayers, and tailings) and special media constitutes the core entry point for elucidating the instability mechanism. Existing studies clarified the respective coupling pathways and controlling factors based on the material properties of different zones. The Coarse Particle Enriched Zone (CPEZ) at the soil-rock interface is a sensitive area for rainfall seepage. Focusing on this zone, Song et al. [6] integrated a bonding strength degradation model and a viscous force model into the discrete element software MatDEM. Then they successfully simulated the three-stage evolution process of seepage erosion (initial, transitional development, and destruction) under heavy rainfall, which provides a basis for identifying the critical state of seepage damage in interface areas. Marinho FAM Jr [7] quantified the effects of soil-rock interface roughness and matrix suction on shear strength (cohesion, friction angle) through direct shear tests. His results offer direct experimental evidence for analyzing the interface strength weakening of soil-rock mixed high and steep slopes in open-pit mines under rainfall (due to suction reduction). Ma et al. [8] investigated the instability process of slopes with weak interlayers under rapid water level rise, revealing the mechanisms of interlayer strength weakening, slope toe crack propagation, and segmental collapse induced by water level rise. These results allowed making the instability prediction of high and steep slopes with weak interlayers in open-pit mines caused by pit water level rise. For special lithologies such as red-bed mudstone and gypsum rock, Ye et al. [9] used X-ray visualization technology to uncover the dual mechanical behaviors (nonlinear seepage and abnormal diffusion) of red-bed mudstone under hydraulic action and established a quantitative relationship between the hydraulic damage diffusion coefficient and seepage, providing a methodological basis for predicting hydraulic damage evolution in red-bed high and steep slopes in open-pit mines. Experimental studies by Caselle et al. [10] and Noel et al. [11] confirmed that water saturation significantly weakens the strength and creep characteristics of gypsum rock and reduces the compressive strength and fracture toughness of carbonate rocks (positively correlated with porosity). These results offer data support for evaluating rainfall-induced strength attenuation in high and steep slopes containing such lithologies in open-pit mines. The seepage-piping mechanism of tailings pond slopes as a special type of high and steep slopes in open-pit mines garnered considerable attention. Through laboratory tests and discrete element simulations, Wang et al. [12] analyzed the piping-seepage mechanism of tailings with varying fine particle contents, clarifying the regulatory effects of fine particle content on tailings permeability, critical water head, and microstructures (force chains, porosity), which provides a micromechanical basis for piping prevention in high and steep slopes of open-pit tailings ponds. Gupta et al. [13] reviewed numerical simulation studies on waste dump slopes (a typical type of high and steep slope in open-pit mines), highlighting the instability risks posed by increasing mine depth and waste volume and constructing a systematic framework for the stability analysis of waste dump slopes in open-pit mines.

The seepage-damage coupling process of high and steep slopes in open-pit mines is frequently influenced by the superposition of multiple environmental factors. In addition to rainfall, such factors as water level fluctuations and freeze-thaw cycles are also considered for their complex multi-field coupling effects. The instability mechanism of slopes in special regions (e.g., loess and seasonal frozen soil) exhibits unique characteristics. The combined effects of water level changes and rainfall are common contributors for the slope instability in open-pit mines. Based on unsaturated seepage principles, Su et al. [14] conducted finite element simulations to quantify the impacts of rainfall intensity, rainfall type, and water level change rate on the slope pore water pressure and safety factor, clarifying the stability critical thresholds under their combined action. For open-pit mine slopes in loess areas, Tang et al. [15] compared the effects of two infiltration modes (unsaturated infiltration and dominant channel infiltration) on slope instability. He confirmed the dominant role of the channel infiltration in accelerating wetting front migration and inducing surrounding rock instability. The results provide references for a rainfall seepage path analysis and instability early warning of loess high and steep slopes. Xuan et al. [16] established a stress-seepage coupling model based on the Tuanjie Tunnel case, quantifying the differential distribution of pore water pressure and the stabilizing effect of pipe roof support. Their results can be applied to a stress-seepage coupling analysis and supporting optimization of loess high and steep slopes in open-pit mines. Open-pit mine slopes in seasonal frozen soil regions (e.g., the Qinghai-Tibet Plateau) face the challenge of synergistic effects between freeze-thaw cycles and rainfall. Taking the Xiangda Village as a case study, Sandri et al. [17] revealed the synergistic instability mechanism of slopes in this region: “ice lens formation (weak surface development) during the freezing period → pore water pressure increase during the rainfall period”. This study provided theoretical references for preventing and controlling freeze-thaw-rainfall coupled instability of high and steep slopes in open-pit mines in the Qinghai-Tibet region. Additionally, geological structures such as faults alter seepage pathways. Taking the Guanling landslide as a case study, Gui et al. [18] revealed a four-stage mechanism by which buried reverse faults act as dominant seepage zones to induce slope instability (fault development → rainfall infiltration → sliding surface formation → final collapse). This result provides a geological basis for seepage risk identification of high and steep slopes with fault structures in open-pit mines. In terms of multi-scale coupling simulation, Cheng et al. [19] employed a CFD-DEM multi-scale coupling model to explore the seepage erosion mechanism around sheet pile walls in gap-graded soils, offering references for interface seepage control between supporting structures and slope bodies of open-pit mine slopes. An et al. [20] quantified the mechanical response of tunnel linings under heavy rainfall through on-site monitoring and fluid-solid coupling simulations and proposed early warning thresholds and engineering recommendations, which can be further used in rainfall response evaluations of anchoring structures in open-pit mine slopes.

The analysis of seepage-damage coupling mechanisms is inseparable from innovations in testing technologies and numerical methods. In recent years, the application of visualization techniques, novel numerical approaches, and intelligent algorithms has significantly enhanced the accuracy and efficiency of researches on high and steep slope instability mechanisms in open-pit mines. In terms of testing and observation technologies, Sun et al. [21] proposed the Enhanced X-ray Image Digital Radiography (EXIDR) technique. They made intuitive observation of fluid migration behavior in rock fractures and provided a visualization method for revealing the seepage path evolution of fractured rock slopes in open-pit mines. He et al. [22] applied the Material Point Method (MPM) to reproduce the entire process of compound reactivated landslides, demonstrating its advantages over the Limit Equilibrium Method (LEM) and Finite Element Method (FEM) in kinematic simulation. So, an efficient tool was created for dynamical simulation of the instability process of high and steep slopes in open-pit mines. Garsole et al. [23] reviewed the application of Artificial Intelligence (AI) methods in seepage prediction, analyzed key challenges such as data quality and model interpretability, and proposed corresponding solutions. Their AI model optimization strategies offer references for intelligent seepage prediction of open-pit mine slopes. In terms of parameter determination and model verification, Vincent et al. [24] established a porosity-dependent seepage parameter calculation model based on the Darcy’s law, providing a foundational method for rapidly determining seepage parameters of open-pit mine slopes. Al-Janabi et al. [25] compared physical, mathematical, and numerical models (SEEP/W) for seepage in homogeneous embankment dams to verify the model compatibility. This model comparison approach can guide the selection and verification of seepage prediction models for open-pit mine slopes.

In geotechnical engineering measurements, Wang et al. [26] combined field surveys and PFC2D modeling to investigate surface subsidence and sinkhole evolution induced by underground goafs. They revealed a coupled mechanism driven by geological structures and mining activities, with sinkholes propagating in a directional, three-stage rupture-slip pattern along weak planes. This study provides a mechanistic framework linking theory and practice, advancing the understanding and mitigation of mining-related ground disasters. Niu et al. [27] integrated large-scale shaking table tests and root dynamic shear effects, demonstrating that low-frequency seismic inputs (< 18 Hz) dominate slope displacement, with natural frequency/damping shifts indicating damage stages. Failure occurs at the crest and shallow surface via sliding and crushing, while the bamboo-vegetation system withstands up to 0.6 g shaking. So, it requires 2.1-3.6× elevation amplification, because roots mitigate near-surface response and dissipate energy. This work quantifies plant-structure synergy, offering a science-based approach for ecologically resilient design against shallow landslides and advancing eco-geotechnical slope protection. To address sparse in-situ geotechnical data and limited spatial resolution, Chen et al. [28] proposed a hybrid framework combining conditional variational autoencoders (CVAE) with kriging, guided by geological priors from the random field theory. The model accurately reconstructs spatial distributions of soil properties from limited CPT data, with uncertainty increasing with distance from samples and site complexity. It enables fine-scale subsurface modeling, deepens understanding of testing uncertainty, and provides practical insights for early-phase engineering decisions under data-scarce conditions. Geotechnical investigations often lack quantitative methods for optimal test-point layout, limiting efficiency and cost-effectiveness. Crisp et al. [29] applied a genetic algorithm to optimize borehole placement in single/multi-layer soils, achieving up to AUD 4.2 M in savings (4.2 % of project cost). Their results show that optimized layouts can outperform regular grids in specific geologies, highlighting strong practical potential. To overcome labor-intensive and error-prone geotechnical data acquisition, Ma and Guo [30] developed an IoT-based intelligent framework integrating BPNN and SVM for automated, low-human-intervention, multi-task analysis. Case studies in the Ganjiang River Basin demonstrated markedly improved testing efficiency and accuracy, overcoming traditional hardware constraints. This advancement propels geotechnical experimentation toward intelligent information systems and enables more reliable assessments of foundation bearing capacity and mechanical properties.

This study focuses on the high and steep southern slope of a deeply excavated open-pit mine, primarily composed of quartz-mica schist, mica-quartz schist, and other rocks. Field investigations indicate that slope deformation intensifies significantly during rainfalls. To address this issue, the research adopts a multi-step methodology: (1) systematic rock mechanics tests (including slake durability tests and in-situ direct shear tests) are conducted to obtain reliable calculation parameters; (2) Geo-Studio software is employed to establish both a saturated-unsaturated seepage model and a slope stability analysis model. Based on these provisions, the study conducts an in-depth analysis of the dynamic response of the slope seepage field and the variation law of stability under rainfall infiltration conditions. Ultimately, the corresponding slope instability model and recommendations for control measures are proposed. The technical roadmap is illustrated in Fig. 1.

Fig. 1Technical roadmap

2. Project overview

2.2. Slope seepage mechanics model

The micro-seepage pathways within the rock mass are assumed to follow stochastic distributions. Based on the macroscopic schistosity structure of schist, the pore channels can be generalized as equivalent micro-fracture networks. The surrounding rock is idealized as a series of thin plates, with fluid flowing freely between them (Fig. 3).

Under pore pressure, the aperture of micro-fractures changes, leading to variations in fluid flow behavior. Let $b$ denote the width of an equivalent fracture channel. Then the permeability coefficient of a single fracture is expressed as:

1

$k_f=\frac{b^2}{12},$

where $k_f$ is the permeability of a single-fracture medium, and $b$ is the fracture width. For fractured rock mass, the permeability coefficient is expressed by:

2

$K_f=k_f\frac{\rho g}{\mu }.$

Fig. 2Real-world view of southern slope. Photo by Jia Haipeng in 2025 at the Longyu Open-pit Mine

Substituting the expression for $k_f$ yields, the following equation will be derived:

3

$K_f=\frac{b^2\rho g}{12\mu },$

where, $\rho$ is the fluid density, $\mu$ is the dynamic viscosity, and $b$ is the fracture width.

When the fracture aperture increases from $d$ to $D$ under the seepage pressure, the permeability $K_f$ increases by a factor of $D^2/d^2$. The seepage velocity is $v=ki$, where $k$ is the permeability coefficient and $i$ is the hydraulic gradient. With $i=\mathrm{\Delta }h/L$, the velocity becomes:

4

$\nu =\frac{D^2\rho g}{12\mu }\frac{\mathrm{\Delta }h}{L}.$

Fig. 3Schematic diagram of seepage pathways in rock

This indicates that, under constant other parameters, the seepage velocity increases with fracture aperture. Furthermore, the relationship between the rock mass permeability coefficient and effective stress is expressed as:

5

$K_f=k_0{\sigma }^{-\alpha },$

where $k_0$ is the reference permeability coefficient at unit stress, $\sigma$ is the effective stress, and $\alpha$ is the fractal dimension of fracture density. Variations in permeability correspondingly alter the effective stress distribution within the rock mass.

3. Mechanical properties and seepage parameters of rock mass

Three types of tests were conducted on the southern slope containing quartz-mica schist, mica-quartz schist, marble, skarn, and gneissic granite:

(1) Laboratory Geotechnical Tests.

As per the Standard for Engineering Rock Mass Test Methods (GB/T 50266-2013), the following parameters were measured: bulk density, uniaxial compressive strength, shear strength (cohesion and internal friction angle), elastic modulus, and Poisson’s ratio of rock samples under both dry and saturated conditions. Each test group included five parallel runs, and the average value was adopted.

The results of the laboratory geotechnical tests are summarized in Table 1, presenting the measured values of bulk density, uniaxial compressive strength, shear strength parameters, elastic modulus, and Poisson’s ratio for each lithology of the target slope.

Fig. 4Laboratory geotechnical tests. Photo by Jia Haipeng in 2025 at the Rock Mechanics Laboratory, School of Intelligent Construction, Henan University of Urban Construction

a) Compressive strength test

b) Brazilian disc tensile splitting test

c) Shear strength test

As indicated in Table 1, the mechanical properties vary considerably among lithologies. Marble and gneissic granite exhibit the highest uniaxial compressive strength (57-58 MPa), with cohesion values of 6.19-6.75 MPa and internal friction angles of 39°-43°, indicating superior mechanical performance. In contrast, quartz-mica schist demonstrates the lowest performance, with a saturated uniaxial compressive strength of only 9.78 MPa and a cohesion of 1.66 MPa. Furthermore, its strength attenuation rate upon water saturation reaches 35 %, reflecting high water sensitivity and identifying it as a potential weak zone for the slope stability.

(2) Slake Durability Test.

Slake durability tests were conducted in strict accordance with Specifications for Rock Tests in Hydropower and Water Conservancy Engineering (DL/T 5368-2007). Specimens of quartz-mica schist and mica-quartz schist were prepared, with six groups per lithology and masses controlled within the range of 200-300 g, as shown in Fig. 5. The test procedure was as follows: 1) Specimens were dried at 105-110 °C to a constant mass ($m_0$), then cooled, and weighed. 2) They were freely immersed in water for 12 hours. 3) After immersion, specimens were placed in a sieve bucket of an HNB-1 rock slake durability tester and rotated at a specified speed for 10 minutes per cycle. 4) After each cycle, residual specimens were collected, redried, and weighed (mass after the Nth cycle denoted as $m_N$). 6) The above process was repeated for 11 cycles, as shown in Fig. 6 and Fig. 7.

Fig. 5Schematic diagram of slake durability test process. Photo by Jia Haipeng in 2025 at the Longyu Open-pit Mine and the Rock Mechanics Laboratory, School of Intelligent Construction, Henan University of Urban Construction

a) In-situ sampling

b) Water immersion

c) Slake durability testing

d) Oven-drying

e) Residue after slaking

Fig. 6Morphological evolution of residual blocks during slake disintegration of quartz-mica schist

a) As-received specimen

b) The third cycle

c) The seventh cycle

d) The eleventh cycle

Fig. 7Morphological evolution of residual blocks during slake disintegration of mica-quartz schist

a) As-received specimen

b) The third cycle

c) The seventh cycle

d) The eleventh cycle

The slake durability index $I_{dN}$ was calculated using Eq. (6):

where $m_0$ is the initial oven-dried mass of the specimen and $m_N$ is the residual mass after $N$cycles.

Morphological analysis of the residual blocks after cyclic testing revealed three primary types of disintegration products (Fig. 8): 1) Flaky fragments: 5-10 mm thick and 10-30 mm long, resulting mainly from splitting along schistosity planes. 2) Strip-like and medium-thick fragments: 5-20 mm thick with a particle size of 10-30 mm, also caused by failure along schistosity. 3) Irregular blocky fragments: 1-6 mm in size, generated primarily by spalling failure of the rock surface layer.

Fig. 8Three morphological types of slake disintegration products

Based on the morphology of the disintegration products, two distinct disintegration mechanisms were identified for the two schist types (Figs. 9-10). Quartz-mica schist is characterized by single and multiple splitting along schistosity planes, while mica-quartz schist exhibits an additional surface spalling mechanism. For quartz-mica schist, the edges and corners of disintegration products round rapidly with increasing cycle number, eventually forming flaky fragments, indicating a structure highly susceptible to thorough damage.

Fig. 9Schematic diagram of slake disintegration mechanism for quartz-mica schist

Given the significant proportion of mica-quartz schist and quartz-mica schist in the rock mass, slake durability tests were conducted on these lithologies. The results are presented in Table 2.

To quantitatively evaluate the slake durability of the two schists, their slake durability index (SDI) values were classified according to Gamble`s classification criteria and rock slake disintegration degree standards (Table 3). A statistical analysis was also performed to assess data variability and reliability.

The statistical results indicate that both rock types exhibit increasing variability in SDI with cycle number, as reflected by the rising coefficient of variation. Quartz-mica schist shows greater variability in later cycles, supporting its classification as a low-durability rock.

Fig. 10Schematic diagram of slake disintegration mechanism for mica-quartz schist

The relationship between the slake durability index of the schist specimens and the number of cycles is illustrated in Fig. 11.

Fig. 11Graph of relationship between slake durability index and cycle number for two schist types

As shown in Fig. 11, the slake durability index (SDI) of the schists decreases gradually with increasing cycle number, accompanied with intensified rock disintegration. Notable differences are observed between the two schist types: at the same cycle number, quartz-mica schist exhibits the lowest SDI, indicating severe fragmentation. Based on Gamble’s classification: Mica-quartz schist (1000-1200 m elevation) has an average SDI of 82.23 % that classifies it as a moderate-durability rock with a moderate disintegration state. After 11 cycles, the disintegration products are 3-8 mm thick flaky fragments. Quartz-mica schist (1200-1400 m elevation) has an average SDI of 58.96 % that classifies it as a low-durability rock with a severe disintegration state. Disintegration begins as early as the second cycle, yielding in final products as 1-5 mm thin flakes. This suggests that quartz-mica schist is prone to continuous fragmentation under rainfall, increasing the risk of surface instability on slopes.

(3) In-situ Direct Shear Test.

A total of eleven test pits were excavated on platforms at elevations of 1350 m, 1100 m, and 950 m along the southern slope. Each test block had a size of 60 cm×60 cm×60 cm. Strain-controlled loading was applied to obtain shear stress-shear displacement curves under varying normal loads. Based on the reliability theory, shear strength parameters corresponding to an 80 % guarantee rate were derived. This rate was selected because it optimally balances accuracy and conservatism, ensuring engineering safety without overestimating the material strength. Meanwhile the economic feasibility is maintained by avoiding excessively conservative designs. This approach aligns with established geotechnical reliability-based design practices, where an 80 % threshold is widely adopted for slope stability assessments under uncertainty.

Fig. 12In-Situ direct shear apparatus. Photo by Jia Haipeng in 2025 at the Rock Mechanics Laboratory, School of Intelligent Construction, Henan University of Urban Construction

Eleven sets of direct shear tests (Set 1 to Set 11) were performed, with five specimens per set. Each specimen underwent two types of tests: a shear failure test under different normal loads and a post-failure friction test. Peak shear failure strengths were extracted from the shear stress-displacement curves, and the relationships between peak shear strength and normal stress were established.

In accordance with the Mohr-Coulomb criterion, and considering internal pores within the rock mass, the jack pressure gauge readings exhibited a gradual increase during shear fracture. The corresponding “climbing strength” during this “climbing phenomenon” was identified as the most appropriate parameter for calculating shear strength, as it indicates the onset of sliding and failure. The resulting shear stress-deformation and shear strength-normal stress curves are shown in Figs. 13-15.

With the help of a reliability analysis model, shear strength parameters for different lithologies under varying guarantee rates were determined, as summarized in Table 4.

The analysis reveals that shear strength values decrease as the guarantee rate increases, implying that higher precision requirements lead to more conservative parameter estimates. Therefore, an appropriate guarantee rate should be selected based on specific engineering needs. In this study, an 80 % guarantee rate was adopted to balance safety and economic considerations.

Fig. 13Stress-strain-strength relationship curves for test pit at 1350 m platform

a) Shear stress-shear strain relationship curve

b) Shear strength-normal stress relationship curve

Fig. 14Stress-strain-strength relationship curves for test pit at 1100 m platform

a) Shear stress-shear strain relationship curve

b) Shear strength-normal stress relationship curve

Fig. 15Stress-strain-strength relationship curves for test pit at 950 m platform

a) Shear stress-shear strain relationship curve

b) Shear strength-normal stress relationship curve

Given that three test sets were conducted on a schist layer at the 1180-1320 m elevation, three on a schist layer within 940-1180 m, two on a marble layer, and two on a skarn layer. Then the average shear strength parameters were used for each lithology in subsequent numerical simulations. The final rock mass shear parameters are presented in Table 5.

As shown in Table 5, shear strength parameters decrease with increasing guarantee rates. Under the 80 % guarantee rate, the cohesion and internal friction angle of quartz-mica schist at 1180-1320 m elevation are 24.5 kPa and 34.4°, respectively. Only 37 % of those belong to the lower mica-quartz schist ($c$ = 65.7 kPa, $\varphi$ = 30.9°). Moreover, the fault zone exhibits a cohesion of merely 8.9 kPa, indicating a high potential for sliding surface formation. Based on these results, the parameters at an 80 % guarantee rate were selected as inputs for numerical simulations.

To quantitatively address uncertainties in the direct shear tests, a comprehensive analysis was conducted focusing on the following aspects: 1) Spatial variability: The coefficient of variation (COV) of cohesion for schist at different elevations ranged from 15 % to 28 %, while the internal friction angle showed lower variability (COV: 3 %-8 %). 2) Test measurement uncertainty: Normal stress application exhibited ±5 % control precision, while shear stress measurement uncertainty was quantified as ±3 % of reading values. 3) Statistical uncertainty: The 80 % guarantee rate reduced the cohesion values by 12 %-38 % compared to the mean values (50 % guarantee rate), representing a systematic treatment of parameter uncertainty. The specific uncertainty quantification results are summarized in Table 6.

The stress-strain-strength relationship curves (Figs. 12-14) visually demonstrate the scattering of test results, particularly in the pre-peak region, which reflects the inherent uncertainty in rock mass behavior. The reliability-based approach adopted in this study provides a quantitative framework for incorporating these uncertainties into engineering design parameters.

4. Influence of rainfall on slope stability

4.1. Basis for constitutive model

(1) Geological Profile Simplification.

The Mohr-Coulomb model was adopted to describe the elastoplastic mechanical behavior of the soil. Based on the exploratory cross-sectional profile of the slope (Fig. 16), the strata were simplified into four layers: surface highly weathered schist (permeability coefficient: 2.09×10^-5 m/s), moderately weathered schist (3.2×10^-6 m/s), underlying skarn (1.85×10^-6 m/s), and gneissic granite (1.62×10^-6 m/s). The basal boundary and lateral boundaries above the water table were assigned as no-flow boundaries, whereas the slope surface was specified as the unit-flux boundary. The computational domain was discretized using a hybrid triangular-quadrilateral mesh, comprising 2,160 nodes and 1,922 elements. Three characteristic slope positions-the toe (A), mid-slope (B), and crest (C)-were defined as representative locations for analysis.

Fig. 16Analysis section and initial groundwater level diagram

(2) Design of Computational Scenarios.

Six computational scenarios were formulated by combining two rainfall intensities (42.5 mm/day and 81.7 mm/day) with three durations (24 h, 48 h, and 72 h), as summarized in Table 7. The mechanical parameters of the rock mass were adopted from previous studies (Table 4). The corresponding saturated moisture contents were specified as follows: slate, 0.20; strongly weathered slate, 0.10; and gneissic granite, 0.05.

Table 7Model calculation schemes

Parameter	Unit	Scenario 1	Scenario 2	Scenario 3	Scenario 4	Scenario 5	Scenario 6
Intensity	mm/d	42.5	42.5	42.5	81.7	81.7	81.7
Duration	h	24	48	72	24	48	72

(3) Initial Seepage Field Definition.

The initial seepage field was established based on the following considerations: 1) Geological exploration data were used to determine the stratigraphic structure of the slope and the physical properties of the soils. 2) Long-term groundwater level monitoring data indicated an initial phreatic surface located 2.5 m above the slope toe. 3) Based on field sampling and laboratory tests, the hydraulic parameters of the unsaturated soils were calibrated. The inferred distributions of in-situ pore water pressure and volumetric water content (Figs. 17-18) were employed for back-analysis verification, ensuring the consistency between the initial field and actual conditions.

(4) Soil-Water Characteristic Curve.

The relationship between unsaturated hydraulic conductivity and matric suction was described using the Van Genuchten model (Fig. 19):

7

$k\left(\theta \right)=k_S{\left[\frac{1-\alpha {\psi }^n{(1+\alpha {\psi }^n)}^{-m}}{{(1+\alpha {\psi }^n)}^m}\right]}^2,$

where $\alpha$, $n$, and $m$ are all fitting parameters, and $\psi$ denotes the matric suction.

Fig. 17Pore water pressure distribution map under natural conditions

Fig. 18Volumetric water content distribution map under natural conditions

Fig. 19Soil-water characteristic curve (SWCC)

a) Relationship curve between volumetric water content and negative pore pressure

b) Relationship curve between permeability coefficient and negative pore pressure

As illustrated in Fig. 19, both the volumetric water content and the hydraulic conductivity exhibit a nonlinear decreasing trend with increasing matric suction.

(5) Model Validation.

A comparison between the simulated and field-measured pore water pressures under natural conditions revealed a discrepancy of less than 5 %, indicating that the constitutive model and the initial seepage field definition were reasonable and effectively represented the actual hydrogeological conditions of the slope.

4.2. Vertical distribution patterns of pore water pressure

During rainfall infiltration, notable variations in pore water pressure were observed within the near-surface zone of the slope (at 0-8 m depth), whereas regions below the groundwater table (> 8 m) remained largely unaffected (Fig. 20).

Fig. 20Volumetric water content distribution map under natural conditions

a) Scenario 1: rainfall intensity of 42.5 mm/d with a duration of 24 h

b) Scenario 2: rainfall intensity of 42.5 mm/d with a duration of 48 h

c) Scenario 3: rainfall intensity of 42.5 mm/d with a duration of 72 h

d) Scenario 4: rainfall intensity of 81.7 mm/d with a duration of 24 h

e) Scenario 5: rainfall intensity of 81.7 mm/d with a duration of 48 h

f) Scenario 6: rainfall intensity of 81.7 mm/d with a duration of 72 h

For example, under the rainfall intensity of 42.5 mm/d, the pore water pressure at the slope toe increased from –12 kPa to –3 kPa after 24 h, and further rose to –1 kPa at 48 h, and finally approached the saturated condition (0 kPa) after 72 h. In contrast, the pore water pressure at the slope crest exhibited a relatively minor increase over the same period, rising from –15 kPa to only –8 kPa after 72 h. These observations indicate that the lower part of the slope is more susceptible to rainfall-induced changes, thereby exacerbating the potential for slope instability.

4.3. Temporal effect analysis

Time-history curves of pore water pressure at representative slope locations (toe A, mid-slope B, and crest C) are presented in Fig. 21. The results indicate that the pore water pressure at all monitoring points exhibits a trend of “rapid initial increase followed by a gradual deceleration in the rate of increase” as rainfall duration extends. A comparison between Scenario 1 (24 h) and Scenario 3 (72 h) demonstrates an approximately 1.8-fold increase in the affected depth, accompanied with an expansion of the affected zone toward the interior of the slope mass, thereby exacerbating the slope instability.

Fig. 21Variation of pore water pressure with time at characteristic points on slope surface

a) Rainfall intensity of 42.5 mm/day

b) Rainfall intensity of 81.7 mm/day

5. Conclusions

This study focuses on a representative steep, highly fractured rock slope at the Luanchuan Longyu open-pit mine, addressing the complex scenario in which lithological heterogeneity, intense fracturing, and rainfall-driven hydrological perturbations interact in a non-classical manner. By integrating multiphysics monitoring, laboratory testing, and parameter inversion, the following conclusions are drawn:

Coupled Seepage-Fracture Feedback Mechanism: In the steeply inclined and highly fractured rock slope on the southern highwall of the Luanchuan Longyu open-pit mine, a coupled seepage-fracture positive feedback mechanism governs the progressive failure process, even without a macroscopically continuous slip surface. Under rainfall infiltration or elevated groundwater conditions, increased pore water pressure drives fracture propagation, which in turn enhances both permeability ($K_f$↑) and seepage velocity ($v$↑). This leads to a reduction in effective stress and concomitant weakening of the rock mass, establishing a self-accelerating damage chain characterized by “seepage enhancement → fracture extension → rock mass degradation.” The proposed mechanism quantitatively elucidates the intrinsic failure processes of heavily jointed rock masses under atypical sliding modes, providing a critical theoretical foundation for the stability assessment of steep heterogeneous slopes.

Lithological Heterogeneity and Weathering Susceptibility: The control exerted by lithological heterogeneity and weathering susceptibility on the spatial distribution of slope weak zones was identified and quantified. Results indicate that marble and gneissic granite exhibit superior mechanical properties (unconfined compressive strength, UCS ≈ 57-58 MPa; cohesion $c$≈ 6.2-6.8 MPa; friction angle $\varphi$ ≈ 39°-43°). In contrast, quartz-mica schist displays markedly lower saturated strength (UCS = 9.78 MPa, $c$ = 1.66 MPa) and pronounced water-induced degradation (strength reduction rate = 35 %; slake durability index, SDI, declines to 58.96 % after 11 cycles), identifying it as a typical highly weathering-sensitive weak layer. By integrating the parameter back-analysis at the 80 % confidence level ($c$ = 24.5 kPa, $\varphi$ = 34.4°), this lithology is established as the dominant controlling weak stratum within the slope, while fault zones ($c$ = 8.9 kPa) constitute potential slip surfaces. These findings elucidate the essential mechanism by which lithology-structure interactions govern the failure pattern of the slope.

Multi-Field Coupled Evolutionary Pathway: The multi-field coupled evolutionary pathway through which rainfall infiltration induces the instability in steep slopes is elucidated as follows: rainfall modulates volumetric water content and hydraulic conductivity via changes in matric suction, producing a nonlinear hydrological response; preferential flow at the slope toe exacerbates localized hydraulic discontinuities. As the infiltration proceeds, pore water pressure propagates progressively into deeper strata, attenuating the effective stress field and triggering a transition from surficial stratified degradation to overall slope failure. This work constructs a dynamic “surface-to-depth” progressive failure model, providing a theoretical basis and decision-support framework for intelligent slope hazard mitigation strategies based on rainfall thresholds and real-time monitoring.

(4) A systematic comparison between numerical simulation results, in-situ direct shear tests, and disintegration tests was conducted to validate the reliability of the proposed model in capturing two key aspects: (i) unsaturated seepage behavior, and (ii) spatial variability of rock mass strength parameters. For instance, the simulated pore water pressure at the slope toe deviates from the measured values by less than 5 %, demonstrating high reliability in reproducing transient hydrological processes. Furthermore, the strength parameters of quartz-mica schist, evaluated at the 80 % confidence level, exhibit excellent agreement with the statistical outcomes derived from laboratory tests. These findings collectively indicate that the developed “seepage-damage” coupled model is capable of realistically simulating the progressive failure mechanism of natural slopes under coupled hydro-mechanical loading conditions.