Analysis of Water Migration and Spoil Slope Stability under the Coupled Effects of Rainfall and Root Reinforcement Based on the Unsaturated Soil Theory

Abstract

Root reinforcement is an effective slope protection measure due to root water absorption and soil suction. However, the coupled effect of rainfall and root reinforcement remains unclear, resulting in a challenge to evaluate slope stability in complex environments. This paper regards the root–soil composite as a natural fiber composite and quantifies its reinforcement effect using direct shear tests. The unsaturated soil seepage–stress theory was employed to simulate the effect of rainfall on water migration and the stability of spoil, overburden, and vegetated slopes. Field measurements and pore water pressure tests verified the simulation results. Furthermore, the influences of the slope angle, rainfall parameters, and vegetation cover thickness on slope stability were analyzed. The results showed the following: (1) The root reinforcement enhanced the soil’s ability to resist shear deformation, substantially improving soil shear strength. The cohesion of the root–soil composite (c_rs = 33.25 kPa) was 177% higher than that of the engineering spoil (c_es = 12 kPa) and 32.21% higher than that of the overburden soil (c_os = 25.15 kPa). (2) The overburden and vegetated slopes had lower permeability coefficients and a higher shear strength than the spoil slope, and the effect was more pronounced for the latter, resulting in lower landslide risks. The water migration trend of the vegetated slope was characterized by substantial runoff and a low sediment yield. The safety factors of the spoil slope, overburden slope, and vegetated slope were 1.741, 1.763, and 1.784 before rainfall and 1.687, 1.720, and 1.763 after rainfall, respectively, indicating a significantly higher safety factor of the vegetated slope after rainfall. (3) The slope angle significantly affected slope stability, with lower safety factors observed for higher rainfall intensities and durations. Under these conditions, the slope angle should be less than 30°, and the soil thickness should be 0.5 m for herbaceous vegetation and shrubs and 1.0 m for trees.

1. Introduction

Large spoil slopes are created during construction projects due to extensive engineering construction in China, resulting in increased risks of soil erosion, landslides, collapses, and other disasters. Landslide hazards in abandoned dreg sites pose a high risk. Rainfall and groundwater contribute to slope damage. Rainfall-induced landslides are caused by groundwater recharge and the advance of the wetting front due to rainfall infiltration. Changes in moisture disrupt the mechanical equilibrium [1]. Vegetation influences soil mechanics, soil suction, and hydrology by reinforcing the soil structure and decreasing the permeability coefficient [2]. Numerous studies have demonstrated that a combination of factors causes landslides [3,4]. Although engineering measures can prevent landslides, they are costly and often incompatible with sustainability. Root reinforcement presents a viable alternative for preventing shallow landslides and diminishing slope instability. Numerous researchers have integrated root reinforcement with effective vegetation management and engineering strategies to mitigate shallow landslides and enhance slope stability [5,6,7]. However, the hydrological response mechanism of the root–soil composite to rainfall and the effect of root reinforcement remains unclear, complicating the evaluation of slope stability in complex environments.

Previous research has focused on the influence of single factors on slope stability. Most studies have analyzed some rainfall patterns (uniform, antecedent, and delayed rainfall) and the effects of changes in rainfall parameters, including rainfall threshold, intensity, and duration, on the stability of unsaturated slopes [8,9,10]. However, long-term rainfall intensity can cause significant errors if the ponding depth is considered. Therefore, studies have investigated surface runoff and subsurface seepage. For example, numerous studies have coupled the surface flow equation and the groundwater flow equation to analyze the sensitivity of the influencing factors of infiltration [11]. A multilayer creep–tension mechanism was proposed to analyze slope stability. The combined effect of rainfall and groundwater on multiple weak layers on the slopes was assessed [12]. Furthermore, the slope angle is closely related to the soil’s internal friction angle under rainfall conditions. The larger the slope angle is, the more likely the slope is to collapse rapidly [13].

The mechanical reinforcement effect of plant roots plays an important role in slope stability. The root system can change the soil structure and improve the slope’s mechanical behavior by reinforcing the control of shallow soil movement [14,15]. Researchers usually quantify the root mechanical reinforcement effect through experimental research, theoretical calculation models, and numerical simulations [16]. The Wu model incorporates the root reinforcement effect [17]. The embedded beam element (EBE) model considers root pull-out and fracture damage modes in modeling the root–soil interaction [18]. The fiber bundle and root bundle models consider the asymptotic damage of roots and stress redistribution [2,19]. Researchers have developed a three-dimensional model to assess the effect of root morphology on soil stability [20]. An elastomer was utilized to substitute the root–soil composite. Slopes with plants with different root morphologies showed different mechanical properties. Furthermore, the plant configurations and their spatial variations affect slope stability. Researchers proposed a slope stability model for the regional distribution of shallow landslides, considering the influence of plant roots by analyzing local plant species and their spatial variations [21,22]. In addition to mechanical reinforcement, the root system changes the hydrological response of the soil. Matric suction increases due to the root uptake of water during transpiration and water retention after rainfall events, thereby contributing to slope stability [23]. Various factors affect root water uptake rates and hydraulic properties, including soil properties and hydraulic resistance. Researchers have analyzed rainfall seepage on unsaturated vegetated slopes using Green’s function and water–force coupling [24].

Although considerable work has been conducted to research the effects of a single factor, such as rainfall or vegetation, the coupled influence of rainfall and root reinforcement on slope damage remains unclear. Analyzing a single rainfall factor can improve our understanding of the rainfall-induced landslide mechanism. However, few studies have used the unsaturated seepage–stress theory for practical applications and synthesized hydrogeological, topographic, geomorphological features, and physical properties to study the destabilization of spoil slopes and the potential for landslides. The effect of root reinforcement and water absorption should be considered in slope stability models. The dynamics of the root–soil composite under the coupled effect of rainfall and root reinforcement remain to be addressed. Therefore, it is necessary to evaluate the landslide mechanism under complex multi-physical field coupling accurately.

This study analyzes an abandoned dreg site in Western Sichuan. The effect of root reinforcement on soil shear strength is quantitatively analyzed by direct shear tests. The seepage–stress theory was used to simulate the instability of unsaturated slopes in different restoration stages under the coupled effect of rainfall and root reinforcement. The water migration in the spoil, overburden, and vegetated slopes are comprehensively analyzed to elucidate the influence of root reinforcement on the slope safety factor. Field tests are performed to verify the reliability of the numerical model and method. The influences of the slope angle, rainfall amount, and vegetation cover thickness on slope stability are investigated. The results provide new insights to enable quantitative slope stability evaluations to prevent disasters and manage abandoned dreg sites.

2. Materials and Methods

2.1. Materials and Test Methods

2.1.1. Study Site

The study area is located in Kangding, which is situated in the western part of Sichuan Province, China. The terrain is steep, with rocks primarily consisting of metamorphic and sedimentary types. Field investigations have revealed severe weathering of the rocks in this region, posing a potential risk of landslides. Kangding is located in the sub-humid climate zone, characterized by elevated soil erosion and landslide risks attributable to its prolonged rainy season and high precipitation levels. The average annual temperature is 7.1 °C. Annual precipitation ranges from 800 to 950 mm, with April to September accounting for 90% of the annual precipitation. The spatial distribution of precipitation is uneven, decreasing from southeast to northwest; it exceeds 800 mm in the east.

2.1.2. Soil Samples

Three types of samples were analyzed based on the different stages of restoration at the abandoned dreg site: engineering spoil (Figure 1a,b), overburden soil (Figure 1c), and root–soil composite (Figure 1d,e). The engineering spoil is a fine, granular material similar to soil that is dumped into a slope. Then, overburden soil, which was stripped from the surface soil, was paved on the slope during restoration. Finally, a vegetated slope surface was generated after restoration. Slopes at three stages of ecological restoration encompass spoil slopes, overburden slopes (i.e., spoil slopes covered with overburden soil), and vegetated slopes (i.e., spoil slopes covered with root–soil composites). The natural density of the soil was measured by extracting a sample from the slope with a ring knife. Herbaceous plants adapted to the local geological and hydrological conditions were selected for slope restoration because they have a dense root system that reinforces the soil. Samples of the root–soil composite were obtained, and a direct shear test was conducted to measure the shear strength (Figure 1f). The particle size distribution was determined by sieve analysis. A 4000 g soil sample was dried and weighed, and the particles were sieved through a series of sieves with different pore sizes (Figure 1b). The particle size distribution curve was then determined based on the pore sizes of the sieve mesh and the mass of the particles underneath the sieve mesh (Figure 1g). The soil samples (3000 g) were weighed and dried to a constant weight in an oven at 105 °C. Soil moisture was determined by gravimetric method (weight difference between the wet sample and after oven drying), i.e., ${\theta }_s=m_w/m_s\times 100\%$, where m_w is the weight of the saturated water, and m_s is the weight of the soil after drying [25]. The van Genuchten (VG) model parameters (α and n) of the engineering spoil and overburden soil were obtained from studies on spoils and bare soil’s soil and water characteristics [26,27]. Pore water pressure sensors, Model CYY2, were installed at a depth of 10 cm under the slope surface to continuously monitor the soil’s pore water pressure via Campbell Scientific CR1000X. The soil’s pore water pressure undergoes continuous fluctuations under rainfall. The instantaneous soil water content was determined by the on-site meteorological station, Environmental Monitoring SN-QXZN-M1-DC-12-4G (Jinan, China). The meteorology sensors for soils were buried at a depth of 10 cm under the slope surface. These data are then employed to analyze the relationship between pore water pressure and water content, thereby validating the accuracy of the selected VG model parameters. Assuming that the gas in the soil is uniformly distributed and their influence is relatively small, the pore gas pressure can be neglected; i.e., the matrix suction of the soil is equal to the negative pore water pressure. The soil water retention curve (SWRC) was obtained by fitting the data with RETC 6.02 software [28] (Figure 1h). The soil’s saturated permeability coefficient was confirmed via soil infiltration tests using the ECA-TR09 soil infiltration tester (Beijing, China), conducted at an average depth range of 10–30 cm. The physical and mechanical parameters of the samples are summarized in Table 1.

2.1.3. Direct Shear Tests

The soil shear strength substantially affects slope stability. Significant errors between remodeled and undisturbed samples were prevented by using in situ samples of the engineering spoil, overburden soil, and root–soil composite in the direct shear test. The standard GB-T50123-2019 was used to determine the soil cohesion and internal friction angle, and a fully automatic direct shear instrument (Beijing Huakan Co., Beijing, China) was used [29]. The specimens were placed into the shear box, which contained filter paper and permeable stone. Vertical pressures of 100 kPa, 200 kPa, 300 kPa, and 400 kPa were applied to each group at a 0.8 mm/min shear rate. The maximum shear displacement was set at 6 mm. According to the standard GB-T50123-2019, the shear stress corresponding to a shear displacement of 4 mm is the shear strength when the shear stress–shear displacement curve has no prominent peak.

2.2. Governing Equations

2.2.1. Equations Describing the Seepage in Unsaturated Soils

The Richards equation is based on single-phase flow. The pores between soil particles are not filled with water but contain air. The effect of solid deformation on pore seepage must be considered to assess coupled seepage–stress effects in soils [24]. Therefore, we used the continuity equation for fluid flow in a porous medium under volumetric strain. The governing equation is defined as follows:

(1)$\rho \left(\frac{C_m}{\rho g}+S_{e}S\right)\frac{\partial p}{\partial t}+{\alpha }_o\frac{\partial {\epsilon }_v}{\partial t}\nabla \cdot \rho \left(-\frac{k_s}{\mu }{k}_r\left(\nabla p+\rho g\nabla D\right)\right)=Q_m,$

where ρ is the fluid density; C_m is the specific water capacity; g is the gravitational acceleration; p is the pore water pressure; S_e is the effective saturation; S is the storage coefficient; α_o is the Biot coefficient of the effective saturation; ɛ_v is the volume strain of the medium; k_s is the saturated permeability coefficient; µ is the hydrodynamic viscosity; k_r is the relative permeability coefficient; D is the elevated head; and Q_m is the source term for the seepage. The second term on the left side of Equation (1) describes the effect of soil volume deformation on pore water percolation, which is represented by changing the source term for seepage Q_m in the simulations. S_e, C_m, and k_r are determined by the VG model. Matric suction is a vital parameter affecting the permeability coefficient and saturation in unsaturated soils. The VG model is commonly used to describe the SWRC [30]:

(2)$S_e=\left\{\begin{array}{l}\frac{1}{{\left[1+{\left|\alpha H_p\right|}^n\right]}^m},H_p<0\\ 1,H_p\ge 0\end{array}\right.,$

(3)$C_m=\left\{\begin{array}{l}\frac{\alpha m}{1-m}\left({\theta }_s-{\theta }_r\right){S_e}^{\frac{1}{m}}{\left(1-{S_e}^{\frac{1}{m}}\right)}^m,H_p<0\\ 0,H_p\ge 0\end{array}\right.,$

(4)$k_r={S_e}^l{\left[1-{\left(1-{S_e}^{\frac{1}{m}}\right)}^m\right]}^2,$

where α, m, and n are the curve-fitting parameters; $n={\left(1-m\right)}^{-1}$; H_p is the pressure head; θ_s is the saturated volume water content; and θ_r is the residual water content.

2.2.2. Equations Describing the Stress in Unsaturated Soils

The stress equilibrium equation for soils can be expressed as follows based on the static equilibrium since unsaturated soil is a continuous porous medium:

(5)${\sigma }_{ij,j}+F_i=0,$

where σ_ij is the stress tensor; and F_i is the volume force. When the soil pores are filled with water, the pore pressure reduces the effective stress between the soil particles. In unsaturated soils, an increase in pore water pressure reduces the matric suction, thereby reducing the effective stress. In saturated soils, an increase in pore water pressure directly reduces the effective stress. Therefore, the constitutive equation describing the effect of pore water pressure is rewritten as follows:

(6)${{\sigma }_{ij}}^{\prime }=D_{ijkl}{\epsilon }_{kl}-{\alpha }_{o}p{\delta }_{ij},$

where σ_ij′ is the effective stress tensor; D_ijkl is the elastic parameter tensor; and δ_ij is the Kronecker symbol (i.e., δ_ij = 1 when i = j; otherwise, δ_ij = 0). The second term on the right side of Equation (6) describes the effect of the pore water pressure on the soil’s effective stress. Equation (6) describes the influence of the pore water pressure on the effective stress and the pore seepage in unsaturated soils. The pore water pressure represents the coupled seepage–stress effect as the external stress in the simulation. The geometric equation is expressed as follows:

(7)${\epsilon }_{ij}=\frac{1}{2}\left(u_{i,j}+u_{j,i}\right),$

where ɛ_ij is the strain tensor, and u is the displacement.

When considering transient seepage in unsaturated soils, the unit weight of the unsaturated soil varies rather than remaining fixed (e.g., dry or wet weight) due to the effect of the pore water gravity, which is reflected in the volumetric force F_i. Therefore, the effect of pore water gravity should be considered when analyzing the influence of rainfall infiltration on unsaturated soil slope stability in numerical simulations:

(8)$\gamma ={\gamma }_s+{\gamma }_w\cdot n\cdot s_e,$

where γ_s is the unit weight of dry soil; γ_w is the unit weight of water; n is the porosity of the soil; and γ is the unit weight considering the pore water gravity—this value is not fixed.

2.3. Numerical Solution Method

2.3.1. Geometric Parameters and Model

After conducting a site survey and performing slope measurements, the numerical slope model was established using the 2D plane strain elements in the finite element software COMSOL Multiphysics 6.0 (Figure 2). The slope length and height are 4 H and 2 H, respectively, where H = 5 m. The slope angle is 30°, and the groundwater level height at the left and right boundaries are 3 H and H, respectively. Herbaceous vegetation was used for the restoration of the abandoned dreg site. Therefore, the thickness of the overburden layer in the simulation is 0.5 m, and the overburden of the root–soil composite is assumed to be an isotropic medium. The sample parameters in Table 1 are the shear strengths under total stress conditions. The effective stress accurately reflects the interaction of soil particles without the influence of pore water pressure; therefore, it is closer to the actual soil stress. The effective stress was defined as ${\sigma }^{\prime }=\sigma -p$, where σ′ is the effective stress, and σ is the total stress. This method was used to describe the mechanical properties of the soil in the simulation.

The sample parameters used in the simulation are listed in Table 1. The sample is uniform; the selected spoil slope meets the required isotropy and elastoplastic conditions. The samples’ elastic modulus was determined from the literature since sufficient data were obtained in Western Sichuan, China [31]. The shear strength parameters were determined via direct shear tests. The three slope types included the spoil slope (Case I), the overburden slope (spoil slope covered with overburden soil, Case II), and the vegetated slope (spoil slope covered with a root–soil composite, Case III). Since the accuracy of the numerical simulation depends on the mesh size, the lower layers of the model had a sparse mesh, and the surface soil was densely meshed to reduce computational time while improving calculation accuracy. The damaged surface of layered slopes may not be circular. Therefore, the strength reduction was calculated to determine the slopes’ safety factors [32].

2.3.2. Rainfall Infiltration

The rainfall event in the simulation was selected from the summer 2022 rainstorm event, which was compared to field investigations and data from the weather station. It can represent the typical rainfall situation in Western Sichuan. Moreover, the area often experiences similar rainfall, and the study results apply to similar events. The actual and simulated rainfall amounts are shown in Figure 3. The rainfall simulation was carried out using the software COMSOL for different restoration stages, considering the vegetation’s water requirements during the growth period. The crucial aspect of modeling rainfall-induced slopes using COMSOL is to address the coupling of unsaturated soil seepage stress. As rainwater infiltrates, it alters the pore water pressure within the soil, subsequently affecting the effective stresses between soil particles and leading to slope deformation. In COMSOL, the pore water pressure within the solid mechanics interface is represented as an external stress to realize the influence of the seepage field on the stress field. The equation $S_{\mathrm{e}}\rho \left(\frac{\partial u_x}{\partial t}+\frac{\partial u_y}{\partial t}\right)$ is embedded within the mass source to influence fluid flow within the porous medium via volumetric strain. It realizes the influence of the stress field on the seepage field, thereby achieving fluid–solid coupling. The negative pore water pressure was assumed to be linear in the vertical direction on the slope to simulate the initial unsaturated conditions. A realistic negative pore water pressure was simulated based on the saturation degree of the slope surface to simulate the matric suction on the slope. The pore water pressure was the control condition. Probes in COMSOL were used to measure the instantaneous rainfall and pore water pressure on the slopes. When the pore water pressure on the slope surface was less than 0, the standard flow rate was used as the rainfall intensity, i.e., the flow boundary. Conversely, when runoff and water ponding occurred on the slope surface, a pressure boundary was used. The standard flow rate on the slope changed nonlinearly as k_s·k_r, which was the unsaturated permeability coefficient.

2.4. Experimental Validation Method

The geomorphic features of the landslide were obtained by on-site investigations, including its length, width, and depth. The landslide surfaces of the failures were determined by on-site measurements using a drone and a straightedge. Accurate meteorological data, including rainfall data, were acquired from the weather station established in the researched slope. The failure mode of slopes was simulated via the described numerical simulation method in this paper. Meanwhile, pore water pressure sensors were buried at a depth of 10 cm on the slope surface of the abandoned dreg sites to verify the accuracy of the slope water migration simulation (Figure 4). The simulation results were compared to the on-site data.

3. Results

3.1. Shear Strength for Different Specimens

The shear stress versus shear displacement curves for different specimens and vertical loads are shown in Figure 5. The trends are the same for the three specimens. The constitutive relation is elastoplastic. The shear stresses of the engineering spoil, overburden soil, and root–soil composite increase with increases in the shear displacement and vertical stress. When the shear displacement is small, the shear stress increases linearly, and the constitutive relation is approximately elastic. As the shear displacement increases, the shear stress increases nonlinearly, and the growth rate decreases and stabilizes. Figure 5d shows the fitted curves of the vertical load σ versus the shear strength τ (i.e., the shear stress at a shear displacement equal to 4 mm). The soil cohesion is the y-axis intercept when x = 0. The cohesion of the root–soil composite (c_rs = 33.25 kPa) is 177% higher than that of the engineering spoil (c_es = 12 kPa) and 32.21% higher than that of the overburden soil (c_os = 25.15 kPa).

3.2. Effective Saturation and Pore Water Pressure

As shown in Figure 6, the rainwater infiltration increased the effective saturation of the surface engineering spoil (Case I), overburden (Case II), and the root–soil composite (Case III) from 0.10, 0.55, and 0.80 to 0.46, 0.77, and 0.81, respectively. Although the permeability coefficient is larger than the rainfall intensity in all three cases, the depth of rainwater infiltration and the amount of water transport on the slopes differ. For Case I, the wetting front first reaches the groundwater level at the slope’s foot, and then rises slightly. The wetting area expands continuously until it reaches the groundwater level, and the saturation increases. The water migration trend in Case II was similar to that in Case III. In these cases, rainwater struggled to infiltrate due to the excellent water-holding capacity of the topsoil, combined with the relatively high initial water content and low actual infiltration rate. Specifically, the wetting front remained at the soil stratification interface during rainfall in Case II, causing rainwater to infiltrate slowly after accumulating on the surface. Finally, the saturation of overburden soil increases slightly (Figure 6d). It was observed that infiltration was challenging during rainfall in Case III, with the majority of the rainfall flowing away as runoff along the slope surface (Figure 6f). Compared to the spoil slope, the overburden and root–soil composite demonstrated a superior surface water-holding capacity, resulting in a higher initial moisture content, lower infiltration capacity, and limited rainwater infiltration.

Pore water pressure is a critical indicator of slope stability. Figure 7 depicts the slopes’ pore water pressure before and after rainfall. The rainfall infiltration causes the pore water pressure of the surface layer to increase from −98 kPa to −18.31 kPa in Case I. Rainfall infiltration occurs continuously, and the infiltrated rainwater rapidly changes the soil structure; meanwhile, the pore water pressure increases. The pore water pressure of the surface layer increases from −14.70 kPa to −4.82 kPa in Case II. The infiltration is concentrated in the surface layer of the overburden slope. The pore water pressure of the root–soil composite layer increases from −7.35 kPa to −6.99 kPa in Case III. Runoff occurs in the surface layer.

3.3. Safety Factor

Slope stability analysis is essential for predicting soil settlement, deformation, and slippage due to various loading and environmental conditions. The slope stability analysis is critical for determining the safety of the abandoned dreg site. The slope safety factor can be determined by adding the Mohr–Coulomb criterion to the solid mechanics interface in COMSOL software. The material properties of the Mohr–Coulomb model are parametrically evolved according to a safety factor. The model does not converge at a certain value, which indicates slope instability, where the safety factor is the value when the model does not converge. Figure 8 shows the safety factors before and after rainfall. The safety factors for Cases I, II, and III are 1.741, 1.763, and 1.784 before the rainfall and 1.687, 1.720, and 1.763 after the rainfall, respectively. The greater the rainfall infiltration depth, the larger the decrease in the safety factor. Case I has the largest decrease in the safety factor, and Case III has the smallest one. The reduction in the safety factor is 0.021 for Case III, 0.043 for Case II, and 0.054 for Case I. The actual reduction of the safety factor in all cases in absolute terms is not significant. However, the difference in reduction between the two appears roughly twice and is, therefore, significant. The reduction in the safety factor for vegetated slopes is 61.11% and 51.16% less than for spoil and overburden slopes. The safety factors for the three cases meet the requirements specified in standard GB51018-2014 [33], which mandates that the safety factor of the first level of the abandoned dreg site must not be less than 1.35. Notably, Case III (i.e., vegetated slope) exhibits the highest safety factor and represents the final stage of recovery for the spoil slope. This underscores the effectiveness of vegetation in enhancing the slope’s safety status.

3.4. Factors Affecting the Stability of the Vegetated Slope

This subsection discusses the effect of the slope angle, rainfall parameters, and vegetation parameters on the safety factor of the vegetated slope. The spoil and overburden slopes are not discussed because the vegetated slope is the final slope morphology in engineering.

3.4.1. Slope Angle

The slope angle ranged from 20 to 45° due to the site survey. If the slope height is constant, the soil’s mechanical properties and hydraulic characteristics are the same as in the above model. The influence of the slope angle on the vegetated slope’s stability was investigated for a rainfall intensity of I_r = 72 mm/h and a rainfall duration of T = 9 h. Figure 9 shows the safety factor decreases nonlinearly with an increasing slope angle. The slope angle significantly affects slope stability, and decreasing the slope angle improves slope stability. The fitted equations for the safety factor under natural and rainfall conditions were expressed as follows:

(9)$F_{\mathrm{nc}}=3.06{\mathrm{e}}^{-\frac{\theta }{12.42}}+1.47,$

(10)$F_{\mathrm{ar}}=3.10{\mathrm{e}}^{-\frac{\theta }{13.06}}+1.47$

where F_nc is the safety factor of natural conditions, F_ar is the safety factor after rainfall, and θ is the slope angle.

3.4.2. Rainfall

The slope’s geometric model and hydraulic characteristics were the same as above. The rainfall intensities were 0.5·I_r, I_r, 2·I_r, and 3·I_r. The rainfall intensity was larger than (i.e., 3·I_r) or less than (i.e., 0.5·I_r and I_r) the saturated permeability coefficient of the root–soil composite. The rainfall duration was 15 h. The effects of the rainfall intensity and duration on slope stability were investigated. As shown in Figure 10, the slope safety factor decreases as the rainfall intensity and duration increase. The relationship between rainfall intensity and the saturated permeability coefficient significantly impacts the slope safety factor. When the rainfall intensity is less than the saturated permeability coefficient, i.e., 0.5·I_r, the safety factor of the vegetated slope is reduced by 0.09. However, when the rainfall intensity is larger than the saturated permeability coefficient, i.e., 3·I_r, the safety factor of the vegetated slope is reduced by 0.42.

3.4.3. Vegetation

Herbaceous vegetation, shrubs, and trees are usually planted on slopes during restoration projects. Despite using a fixed soil cover thickness of 0.5 m, as mentioned above, the overburden thickness varied from 0.2 to 2.0 m. Trees have shallow (0.3 m) and deep (1–5 m) roots. Studies have used a range of 0.2–3 m, sufficient for restoration [34]. The slope types’ geometric model and hydraulic characteristics were the same as above. The rainfall intensity was I_r, and the duration was 10 h. The effect of the cover thicknesses on the slope stability was investigated. Figure 11 shows that the safety factor of the slopes differs in the three restoration stages as the cover thickness increases. The presence of roots significantly increases the slope safety factor. When the root depth is small (i.e., less than 1 m), the influence of the cover thicknesses on the safety factor is small. The safety factor increases significantly as the root depth increases from 1 m to 2 m.

3.5. In-Site Experimental Validation

Figure 12 shows two damage types to spoil slopes. Different rainfall intensities and durations cause different damage patterns. The deep landslide in Figure 12a occurred due to long-term, low-intensity rainfall. They were 4.8 m, 5.4 m, and 3.2 m, respectively. Temporary surface runoff and internal seepage after rainwater convergence are the leading causes of shallow slope erosion, as shown in the insert of Figure 12b. A simulation of long-duration, low-intensity rainfall was carried out. The simulated damage mode was consistent with the actual one, as shown in Figure 12b. The sliding surface is deep because the engineering spoil’s large permeability coefficient indicates a large amount of water infiltration. There was no runoff on the slope, and deep sliding occurred.

Three pore water pressure sensors were, respectively, positioned at a buried depth of 10 cm under the slope surface (Figure 4a–c). The measurements under rainfall conditions were juxtaposed with simulated results, as depicted in Figure 13. Notably, the observed pore water pressure data exhibited a trend similar to the simulated results. The pore water pressure on the slope surface increased from 0 to 8 h after the rainfall started and reached the peak. After the rainfall stopped, i.e., from 8 to 24 h, the surface pore water pressure decreased rapidly and stabilized. A pore water pressure test can only provide data at a specific location, whereas the numerical simulation provides continuous data for the infiltration surface and the stratification, as shown in Figure 7. The numerical results have errors because it is difficult to simulate boundary conditions precisely. Moreover, the porosity and infiltration routes on the slope are not uniform.

4. Discussions

4.1. Root Reinforcement Effectiveness

In general, the root–soil composite has the highest cohesion. The root system affects φ, but no trend is evident. A higher c value indicates higher soil strength. The slopes of the stress–shear displacement curves in Figure 5 indicate that the root system improves the soil shear strength. The root system enables the soil to resist shear deformation. The reason is that the shear is generally larger in the root–soil composite than in soils without vegetation [2]. Therefore, the strength and cohesion of the root–soil composite are significantly higher. The reinforcement effect of plant roots significantly increases the soil strength, and the increased cohesion increases the shear strength of the root–soil composite [35]. Root reinforcement increases soil shear strength and slope stability. Therefore, the vegetated slope has the maximum safety factor. Meanwhile, the roots improve slope stability by reducing the soil infiltration rate and improving shear strength.

4.2. Relationship between Landslide and Water in the Soil

As the water infiltrates the slope, the pore water pressure increases continuously, decreasing the matric suction inside the slope. The internal stresses are distributed unevenly, and the shear strength of the medium decreases, resulting in slope damage. Figure 6 shows that the infiltration capacity is lower for the slopes with a cover layer (Case II and Case III). Rainfall only infiltrates the surface layer, or runoff occurs. The water content is higher in the surface layer than in the interior, and the change in soil porosity is small. Thus, the pore water pressure only changes in the shallow area of the slope. Moreover, Figure 7 shows the pore water pressure has a linear distribution, and the infiltration into the deep layer is minimal, indicating that the soil cover of the overburden slope and the vegetation and soil layers in the vegetated slope increase the slope stability by reducing rainwater infiltration, thereby preventing landslides.

The rainfall saturates the pore space by displacing air during infiltration. The infiltration of large amounts of rainwater causes changes in the slope’s seepage field. As the saturation rises, the pore water pressure increases, and the matric suction decreases. Furthermore, the increasing soil saturation causes soil softening, increasing the shear stress and decreasing the shear strength at the slip surface of the slope [12]. Figure 6 indicates that the overburden and vegetated slopes have better water migration and lower rainfall infiltration, resulting in lower landslide risks. The water migration trend of the vegetated slope is characterized by substantial runoff and low sediment yield [36].

4.3. Influence of Slope Angle, Rainfall, and Vegetation on the Vegetated Slope Stability

According to the fitting Equations (9) and (10) for the safety factor to the slope angle, it can be seen that the effect of rainfall is slightly larger on gentle slopes than on steep slopes because gentle slopes allow more rainwater infiltration, and the water remains longer on the slope surface. Therefore, the pore water pressure increases, resulting in a lower safety factor [24]. Considering engineering safety requirements, we recommend a slope angle below 30° after the restoration to achieve a safety factor of at least 1.78. The higher the rainfall intensity and the longer the duration, the larger the slope safety factor reduction. As the duration increases, the reduction in the slope safety factor reaches a threshold value because the soil is saturated, and no more water can infiltrate the soil, as shown in Figure 10. Figure 11 shows that, when the root depth exceeds 2 m, the safety factor reaches the maximum value and stabilizes [37]. The deformation and displacement of the slope toe decreases as the cover thickness increases. Therefore, the thickness of the root–soil composite should be 0.5 m for herbaceous vegetation and shrubs and 1.0 m for trees to increase slope stability and reduce costs.

4.4. Slope Failure Mechanism

The permeability coefficient of the engineering spoil is much larger than the rainfall intensity. The pore water pressure increases when the rainwater infiltrates the deeper part of the slope. Meanwhile, the matric suction decreases, and the slope’s self-weight increases, ultimately leading to the collapse and a severe landslide (Figure 12a). During short-term, high-intensity rainfall, temporary surface runoff and internal seepage occur. Rainwater is the leading cause of shallow slope erosion, converging at the foot of the slope (Figure 12b). Subsequently, cracks form at the foot of the slope and propagate, causing soil erosion [38].

Our study focuses on slope stability. The field measurements verified the accuracy of the numerical simulation. Soil erosion is a complex process. The slope failure generally started with shallow surface erosion, followed by deep sliding and a landslide of the entire slope. However, the surface soil erosion was not simulated. Future research will focus on the development and mechanism of shallow surface soil erosion. A shallow erosion simulation must consider the influence of dry cracks, which was ignored in this study.

5. Conclusions

This paper analyzed abandoned dreg sites in West Sichuan Province, China. The root reinforcement mechanism of the root–soil composite was quantified in terms of mechanical characteristics. According to the seepage–stress theory of unsaturated soil, the coupled effect of rainfall and root reinforcement on slope instability in different restoration stages was numerically investigated and verified by field measurements and tests. The influences of the slope angle, rainfall parameters, and vegetation parameters on the safety factor of the vegetated slopes were determined. The main conclusions are as follows:

(1) Root reinforcement enabled the soil to resist deformation and stress. The root system significantly increased the shear strength of the root–soil composite. Direct shear tests showed that the cohesion of the root–soil composite (c_rs = 33.25 kPa) was 177% higher than that of the engineering spoil (c_es = 12 kPa) and 32.21% higher than that of the overburden soil (c_os = 25.15 kPa).

(2) The plant roots improved slope stability by increasing the soil’s shear strength and decreasing the infiltration rate. The root system improved the slope’s mechanical strength. The safety factors of the spoil, overburden, and vegetated slopes were 1.741, 1.763, and 1.784 before the rainfall and 1.687, 1.720, and 1.763 after the rainfall, respectively. The reduction in the safety factor for vegetated slopes is 61.11% and 51.16% less than for spoil and overburden slopes. The slope safety factor met the code requirements, and the field test verified the model’s accuracy.

(3) The slope angle significantly affected slope stability. It decreased nonlinearly with an increase in the slope angle. The slope safety factor decreased significantly as the rainfall intensity and duration increased. After the restoration of the spoil slope, the slope should have an angle of less than 30°, and the cover thickness should be 0.5 for herbaceous vegetation and shrubs and 1.0 m for trees.

Author Contributions

All authors contributed to the study conception and design. J.H. and H.S. were responsible for the conception and design of the study. H.S. and J.H. designed and performed the numerical simulation; J.H., H.S., Z.Z. and Q.J. analyzed and interpreted the data. J.H., L.L., C.H. and Y.Z. were responsible for the critical revision of the manuscript for important intellectual content and improved the language. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data used are primarily reflected in this article. Other relevant data are available from the authors upon request.

Conflicts of Interest

Authors Lanhua Liu, Caisong He, and Yang Zhou were employed by the company China Academy of Railway Science Co., Ltd. The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Table

Figure 1. Soil samples and properties: (a) sampling engineering spoil; (b) particle size analysis; (c) sampling overburden soil; (d) sampling root–soil composite; (e) cutting redundant stems; (f) direct shear test; (g) particle size distribution curves; and (h) soil water retention curves.

Figure 2. Geometric parameters and mesh of the numerical model.

Figure 3. Actual and simulated rainfall.

Figure 4. Experimental setup for pore water pressure measurements: (a) overview of sensor configuration, (b) sensor layout, (c) sensor detail, (d) weather station, (e) pore water pressure device, (f) data collection and transmission, and (g) data acquisition system.

Figure 5. Shear displacement versus shear stress: (a) engineering spoils, (b) overburden soil, (c) root–soil composite, and (d) vertical stress versus shear stress.

Figure 6. Saturation of the slopes before and after rainfall.

Figure 7. Pore water pressure on the slopes before and after rainfall.

Figure 8. The slope safety factor before and after rainfall.

Figure 9. Relationship between the slope safety factor and the slope angle.

Figure 10. Relationship between slope safety factor and rainfall.

Figure 11. Relationship between slope safety factor and vegetation cover thickness.

Figure 12. Comparison of observed and simulated rainfall-induced damage on spoil slopes: (a) deep landslides, and (b) sliding surface.

Figure 13. Comparison of measured and simulated pore water pressures.

Table 1

Soil material characteristics.

Material	Density ρ (kg·m−3)	Elastic Modulus E (MPa)	Poisson Ratio μ	Cohesion c (kPa)	Internal Friction Angle φ (°)	Saturated Permeability Coefficient ks (m·s−1)	Natural Water Content θo (%)	Saturated Water Content θs (%)
Engineering spoil	1695	20	0.25	12	23.78	6.0 × 10−4	4.05%	36.1%
Overburden soil	1662	10	0.20	25.15	17.33	5.55 × 10−5	32.47%	48.9%
Root–soil composite	1281	13.8	0.24	33.25	12.81	4.17 × 10−5	50.14%	62.3%

References

1. Chen, R.; Ge, Y.H.; Chen, Z.K.; Liu, J.; Zhao, Y.R.; Li, Z.H. Analytical solution for one-dimensional contaminant diffusion through unsaturated soils beneath geomembrane. J. Hydrol.; 2019; 568, pp. 260-274. [DOI: https://dx.doi.org/10.1016/j.jhydrol.2018.10.057]

2. Meijer, G.J. A generic form of fibre bundle models for root reinforcement of soil. Plant Soil; 2021; 468, pp. 45-65. [DOI: https://dx.doi.org/10.1007/s11104-021-05039-z]

3. Ogila, W. Analysis and assessment of slope instability along international mountainous road in North Africa. Nat. Hazards; 2021; 106, pp. 2479-2517. [DOI: https://dx.doi.org/10.1007/s11069-021-04552-9]

4. Li, S.H.; Cui, P.; Cheng, P.; Wu, L.Z. Modified Green–Ampt model considering vegetation root effect and redistribution characteristics for slope stability analysis. Water Resour. Manag.; 2022; 36, pp. 2395-2410. [DOI: https://dx.doi.org/10.1007/s11269-022-03149-6]

5. Li, Y.P.; Wang, Y.Q.; Ma, C.; Zhang, H.L.; Wang, Y.J.; Song, S.S.; Zhu, J.Q. Influence of the spatial layout of plant roots on slope stability. Ecol. Eng.; 2016; 91, pp. 477-486. [DOI: https://dx.doi.org/10.1016/j.ecoleng.2016.02.026]

6. Ji, J.N.; Kokutse, N.; Genet, M.; Fourcaud, T.; Zhang, Z.Q. Effect of spatial variation of tree root characteristics on slope stability. A case study on Black Locust (Robinia pseudoacacia) and Arborvitae (Platycladus orientalis) stands on the Loess Plateau, China. Catena; 2012; 92, pp. 139-154. [DOI: https://dx.doi.org/10.1016/j.catena.2011.12.008]

7. Su, L.J.; Hu, B.L.; Xie, Q.J.; Yu, F.W.; Zhang, C.L. Experimental and theoretical study of mechanical properties of root-soil interface for slope protection. J. Mt. Sci.; 2020; 17, pp. 2784-2795. [DOI: https://dx.doi.org/10.1007/s11629-020-6077-4]

8. Wang, Y.X.; Chai, J.R.; Cao, J.; Qin, Y.; Xu, Z.G.; Zhang, X.X. Effects of seepage on a three-layered slope and its stability analysis under rainfall conditions. Nat. Hazards; 2020; 102, pp. 1269-1278. [DOI: https://dx.doi.org/10.1007/s11069-020-03966-1]

9. Liu, Z.Z.; Yan, Z.X.; Qiu, Z.H.; Wang, X.G.; Li, J.W. Stability analysis of an unsaturated soil slope considering rainfall infiltration based on the Green-Ampt model. J. Mt. Sci.; 2020; 17, pp. 2577-2590. [DOI: https://dx.doi.org/10.1007/s11629-019-5744-9]

10. Palazzolo, N.; Peres, D.J.; Creaco, E.; Cancelliere, A. Using principal component analysis to incorporate multi-layer soil moisture information in hydrometeorological thresholds for landslide prediction: An investigation based on ERA5-Land reanalysis data. Nat. Hazards Earth Syst. Sci.; 2023; 23, pp. 279-291. [DOI: https://dx.doi.org/10.5194/nhess-23-279-2023]

11. Khan, K.U.J.; Wang, C.; Khan, M.W.J.; Liang, Z.; Li, S.; Li, B. Influence of rainfall infiltration on the stability of unsaturated coal gangue accumulated slope. J. Mt. Sci.; 2021; 18, pp. 1696-1709. [DOI: https://dx.doi.org/10.1007/s11629-020-6567-4]

12. Ren, Y.; Li, T.B.; Dong, S.M.; Tang, J.L.; Xue, D.M. Rainfall-induced reactivation mechanism of a landslide with multiple-soft layers. Landslides; 2020; 17, pp. 1269-1281. [DOI: https://dx.doi.org/10.1007/s10346-020-01357-y]

13. Kumar, S.; Roy, L.B. Experimental and Numerical Analysis of Unsaturated Soil Slope Stability with Rainfall and Jute Fibre Reinforcement Condition. Geotech. Eng.; 2023; 54, pp. 15-32.

14. Mao, Z.; Wang, Y.; Mccormack, M.; Rowe, N.; Deng, X.B.; Yang, X.D.; Shangwen, X.; Nespoulous, J.; Sidle, R.; Guo, D. et al. Mechanical traits of fine roots as a function of topology and anatomy. Ann. Bot.; 2018; 122, pp. 1103-1116. [DOI: https://dx.doi.org/10.1093/aob/mcy076] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/29846521]

15. Li, Y.P.; Hu, C.; Jian, L.; Zhao, R.; Li, C. Evaluation of the stability of vegetated slopes according to layout and temporal changes. J. Mt. Sci.; 2021; 18, pp. 275-290. [DOI: https://dx.doi.org/10.1007/s11629-020-6022-6]

16. Hao, G.L.; Wang, L.G.; Liu, X.F. Methods for Studying the Effect of Plant Roots on Soil Mechanical Reinforcement: A Review. J. Soil Sci. Plant Nutr.; 2023; 23, pp. 2893-2912. [DOI: https://dx.doi.org/10.1007/s42729-023-01330-3]

17. Wu, T.H.; Mckinnell, W.P.; Swanston, D.N. Strength of tree roots and landslides on Prince of Wales Island, Alaska. Can. Geotech. J.; 1979; 16, pp. 19-33. [DOI: https://dx.doi.org/10.1139/t79-003]

18. Zhu, J.; Leung, A.; Wang, Y. Modelling root–soil mechanical interaction considering root pull-out and breakage failure modes. Plant Soil; 2022; 480, pp. 675-701. [DOI: https://dx.doi.org/10.1007/s11104-022-05613-z]

19. Thomas, R.; Pollen-Bankhead, N. Modeling root-reinforcement with a fiber-bundle model and Monte Carlo simulation. Ecol. Eng.; 2010; 36, pp. 47-61. [DOI: https://dx.doi.org/10.1016/j.ecoleng.2009.09.008]

20. Ying, Z.X.; Liu, Y.J.; Shi, C.; Lan, Z.C.; Zhang, H.; Ge, G. Simulating root distribution of plant individual with a three-dimensional model. Ecol. Model.; 2021; 455, 109649. [DOI: https://dx.doi.org/10.1016/j.ecolmodel.2021.109649]

21. Mao, Z.; Roumet, C.; Rossi, L.; Merino-Martín, L.; Nespoulous, J.; Taugourdeau, O.; Boukcim, H.; Fourtier, S.; Rey Granado, M.; Ramel, M. et al. Intra- and inter-specific variation in root mechanical traits for twelve herbaceous plants and their link with the root economics space. Oikos; 2022; 2023, e09032. [DOI: https://dx.doi.org/10.1111/oik.09032]

22. Masi, E.B.; Tofani, V.; Rossi, G.; Cuomo, S.; Wu, W.; Salciarini, D.; Caporali, E.; Catani, F. Effects of roots cohesion on regional distributed slope stability modelling. Catena; 2023; 222, 106853. [DOI: https://dx.doi.org/10.1016/j.catena.2022.106853]

23. Qin, C.B.; Wang, R.; Chen, W.K.; Shi, Y.S.; Sun, H.X.; Tang, J.J.; Wang, L.Q. Stability of Ficus virens-Reinforced Slopes Considering Mechanical and/or Hydrological Effects. Forests; 2024; 15, 133. [DOI: https://dx.doi.org/10.3390/f15010133]

24. Wu, L.Z.; Cheng, P.; Zhou, J.T.; Li, S.H. Analytical solution of rainfall infiltration for vegetated slope in unsaturated soils considering hydro-mechanical effects. Catena; 2022; 217, 106472. [DOI: https://dx.doi.org/10.1016/j.catena.2022.106472]

25. Terzaghi, K.T. Theoretical Soil Mechanics; USA Theoretical Soil Mechanics INC: New York, NY, USA, 1943; 49.

26. Liu, Q.; Huang, J.K.; Zhang, Z.W.; Liu, G.M.; Jiang, Q.O.; Liu, L.H.; Khan, I. Effect of freeze-thaw cycles on the shear strength of root-soil composite. Materials; 2024; 17, 285. [DOI: https://dx.doi.org/10.3390/ma17020285] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/38255453]

27. Zhao, M.Y.; Xie, Q.; Ren, X.H.; Zhu, L.; Chen, P.F. Experimental research on the deformation and failure characteristics of coarse-grained soil slopes in seasonal frozen region along Sichuan-Tibet Railway. IOP Conf. Ser. Earth Environ. Sci.; 2020; 455, 12135. [DOI: https://dx.doi.org/10.1088/1755-1315/455/1/012135]

28. Araújo-Neto, R.; Magalhães, I.; Lyra, G.; Maia, S.; Lyra, G. Physical and hydrological characteristics and modelling of the soil water retention curve in the brazilian semi-arid region. Biosci. J.; 2022; 38, e380311-13. [DOI: https://dx.doi.org/10.14393/BJ-v38n0a2022-54195]

29. GB/T 50123-2019 Ministry of Housing and Urban Rural Development of the People’s Republic of China. Standard for Geotechnical Testing Method; China Planning Press: Beijing, China, 2019; pp. 122-125.

30. Hou, J.; Yang, J.Y.; Tan, J. A new method for revealing spatial relationships between shrubs and soil resources in arid regions. Catena; 2019; 183, 104187. [DOI: https://dx.doi.org/10.1016/j.catena.2019.104187]

31. Ran, T.; Huang, Z.Q.; Meng, L.C.; Wang, C.; Li, X.L. Mechanisms and stabilization measures of layered slope: A Study case of a layered slope in Sichuan, China. Geotech. Geol. Eng.; 2022; 40, pp. 785-794. [DOI: https://dx.doi.org/10.1007/s10706-021-01924-7]

32. Bao, Y.D.; Li, Y.C.; Zhang, Y.S.; Yan, J.H.; Zhou, X.; Zhang, X.D. Investigation of the role of crown crack in cohesive soil slope and its effect on slope stability based on the extended finite element method. Nat. Hazards; 2022; 110, pp. 295-314. [DOI: https://dx.doi.org/10.1007/s11069-021-04947-8]

33. GB 51018–2014 Ministry of Housing and Urban Rural Development of the People’s Republic of China. Code for Design of Soil and Water Conservation Engineering; China Planning Press: Beijing, China, 2014; pp. 20-24.

34. Fiolleau, S.; Uhlemann, S.; Falco, N.; Dafflon, B. Assessing probability of failure of urban landslides through rapid characterization of soil properties and vegetation distribution. Geomorphology; 2022; 423, 108560. [DOI: https://dx.doi.org/10.1016/j.geomorph.2022.108560]

35. Feng, S.; Liu, H.W.; Ng, C.W.W. Analytical analysis of the mechanical and hydrological effects of vegetation on shallow slope stability. Comput. Geotech.; 2020; 118, 103335. [DOI: https://dx.doi.org/10.1016/j.compgeo.2019.103335]

36. Gholami, V.; Sahour, H.; Khaleghi, M.; Yousefi, A. Evaluating the effects of vegetation and land management on runoff control using field plots and machine learning models. Environ. Sci. Pollut. Res.; 2022; 30, pp. 31202-31217. [DOI: https://dx.doi.org/10.1007/s11356-022-24347-0] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/36445518]

37. Gao, Q.F.; Zeng, L.; Shi, Z.N. Effects of desiccation cracks and vegetation on the shallow stability of a red clay cut slope under rainfall infiltration. Comput. Geotech.; 2021; 140, 104436. [DOI: https://dx.doi.org/10.1016/j.compgeo.2021.104436]

38. Yang, S.R.; Huang, L.J. Infiltration and failure behavior of an unsaturated soil slope under artificial rainfall model experiments. Water; 2023; 15, 1599. [DOI: https://dx.doi.org/10.3390/w15081599]