2.1 Temperature-dependent functions

The change in temperature can lead to various alterations in the physical properties of unsaturated soil, thereby affecting slope stability. This paper mainly explores the impact of temperature on capillary water. The Young-Laplace equation [25] is employed to represent matrix suction.

(1)

where represents the matrix suction; represents the water-air surface tension; represents the wetting coefficient; represents the average radius of the water-air interface; and respectively represent the pore air pressure and the pore water pressure. The partial derivative of matric suction with respect to temperature is [26]:

(2)

From this partial differential equation, it can be observed that the temperature affects capillary suction primarily in two aspects: 1) the wetting coefficient; 2) the water-air surface tension. A linear non-isothermal is proposed function to describe the effect of temperature on the water-air surface tension :

(3)

where both and are fitting parameters for the linear non-isothermal functions. GRANT and SALEHZADEH [20] performed linear regression on the tension data of the reference interface to fit into the linear non-isothermal function, obtaining estimates for the two fitting parameters:

(4)

(5)

GRANT and SALEHZADEH [20] investigated the influence of changes in the wetting coefficient of porous solids on the temperature sensitivity of the capillary pressure function, proposing a partial differential equation for the wettability coefficient that is temperature-dependent.

(6)

where represents the enthalpy of immersion per unit area. WATSON [27] proposed a method for simulating the relationship between temperature and enthalpy.

(7)

where represents the enthalpy of immersion per unit area at the reference temperature. Substituting Eq. (3) into Eq. (5), the non-isothermal function form of the wetting coefficient can be solved.

(8)

where is a constant determined based on the wetting coefficient at the reference temperature, given by the following expression [20]:

(9)

GRANT and SALEHZADEH [20] substituted the non-isothermal expressions for the wetting coefficient and the water-air surface tension into Eq. (2), resulting in:

(10)

where is a coefficient proposed by GRANT and SALEHZADEH [20], which is equal to the negative ratio of the unit-area immersion enthalpy to the constant . Equation (10) can be solved to obtain a non-isothermal explicit expression for matric suction [28], which reflects the effect of temperature on soil matric suction.

(11)

where represents the magnitude of matric suction in unsaturated soil at the reference temperature.

2.2 Modified shear strength of unsaturated soil

In the past, the Mohr-Coulomb yield criterion was employed to estimate the shear strength of soil. However, it did not account for soil unsaturation. FREDLUND and XING [29] bridged the classic Mohr-Coulomb yield criterion from saturated to unsaturated soil and proposed a method for calculating the shear strength of unsaturated soil.

(12)

where , and respectively represent the shear strength, effective cohesion, and effective internal friction angle of unsaturated soil; represents the normal stress on the failure plane; represents the friction angle corresponding to the variation in matric suction when remains constant.

However, in subsequent experiments, the shear strength of unsaturated soil consistently exhibited a nonlinear behavior with the increase of matric suction, indicating that the strength envelope curve should not be merely a straight line [28, 30, 31]. From Figure 1, it is visually apparent that the strength envelope exhibits linear characteristics only when the matric suction in unsaturated soil is less than the air-entry value (AEV). In the subsequent stage, the strength envelope begins to show a curved trend, indicating the non-linear variation of shear strength with matric suction. Once the matric suction surpasses the residual suction, the strength envelopes curves of various soils start to exhibit distinct non-linear developments. FREDLUND et al [30] proposed a nonlinear shear strength prediction equation:

Figure 1全屏查看

Correlation between the unsaturated shear strength boundary and the soil-water characteristic

(13)

where represents the apparent cohesion of the soil; represents a fitting parameter; and respectively represent the volumetric water content of the soil in the actual state and in the saturated state. The parameter is ascertained through the formula introduced by FREDLUND and XING [29], articulated as follows:

(14)

where is the fitting parameter; is the angle of inclination at the inflection point; indicates the matric suction at the transition point; is residual matric suction; is matric suction. By substituting Eq. (11) into Eq. (13), the shear strength of unsaturated soil at different temperatures can be calculated, requiring only the matrix suction at the reference temperature. It is important to note that when the groundwater level is above the soil surface, the soil is considered saturated. The main focus of this paper is unsaturated soil; therefore, the following discussion pertains only to situations where the groundwater level is below the ground surface.

Table 1 provides parameters for various soil types as summarized by ZHANG et al [17], used to analyze the influence of different unsaturated soil types on shear strength. Here, SWCC1 to SWCC4 respectively represent sand, silt, clay and extremely fine-grained soil.

3.1 3D rotational failure mechanism with a vertical crack

Establishing a kinematically admissible failure mechanism is a prerequisite for using limit analysis [4, 33, 34]. The failure mechanism for the two-step slope is illustrated in Figure 2. The total height of the two-step slope is . The vertical height from the crest to the step plane is , the vertical height from the step plane to the bottom of the slope is . Obviously,

(15)

Figure 2全屏查看

3D rotational wedge failure mechanism with a crack

where and are both depth coefficients. HE et al [24] proposed a failure mechanism applicable to single-step slopes with open cracks at the crest, which incorporates a vertical crack into the existing wedge failure mechanism. Apply this failure mechanism to the two-step slope, considering both open cracks and formation cracks simultaneously. As shown in the Figure 2, the cross-section of the mechanism is composed of slope surfaces FB, BC, CD and DE, vertical crack FG, and sliding surface AD, rotating rigidly about the constant angular velocity at center O. According to the normality rule, it can be inferred that the upper and lower profiles of the failure mechanism and are two logarithmic spiral curves, represented as:

(16)

(17)

where , and are the initial lengths and angles of the logarithmic spirals, and is the effective internal friction angle of the unsaturated soil. The cross-section of the failure mechanism in the rotational direction is a circle with a radius of , and the failure mechanism can be considered as a collection of circles rotating with varying radii between logarithmic spiral curves. The distance between the centers of these circles and the rotation center O is .The expressions for and are as follows:

(18)

(19)

As shown in the Figure 2, there is a vertical crack in the failure mechanism, and the position and length of this crack are unknown. Studying the influence of cracks on slope stability often involves identifying the crack that is most detrimental to stability among all possible cracks and determining its position and length through optimization. By introducing a new variable to represent the position and length of the vertical crack, where is the intersection point of the logarithmic spiral AE and the vertical crack FG. The crack can be expressed as:

(20)

where is the length of the vertical crack; is the length of the slope crest after introducing the crack(see Appendix); and are the two angles of the two-step slope; is the maximum angle of rotation for the mechanism, as shown in the Figure 2. Figure 3 illustrates that by incorporating a planar strain block with a width of into the mechanism, the crack failure mechanism can transition to a 2D failure mode when the slope width approaches infinity [32, 35, 36]. The width of the inserted block can be represented as:

(21)

Figure 3全屏查看

Rotational wedge failure mechanism with a crack incorporating a planar block

where represents the maximum width of the rotating mechanism, and represents the maximum width of the slope.

From the above description, it can be understood that the crack failure mechanism can be controlled by four independent variables:, , and . To accommodate a reasonable and acceptable failure mechanism, the four independent variables should satisfy the following constraints:

(22)

where represents the maximum allowable depth of the vertical crack. The angles , , and represent the polar angles passing through points F, B, C and D, respectively. Their detailed expressions are provided in the Appendix.

3.2 Work rate balance equation

When utilizing the upper bound limit analysis method for solving, it’s essential to first establish a balance between the external power and internal dissipation rates within the failure mechanism. This section will focus on detailing how to construct the power equilibrium equation for the two-step slope with a crack.

In this study, internal energy is dissipated due to the resistance of unsaturated soil, while external power is the work done by the slope soil’s self-weight. Therefore, the power equilibrium equation can be expressed as:

(23)

where and respectively represent the rates of internal energy dissipation generated by the effective cohesion and the apparent cohesion in the soil; represents the rate of external power generated by the self-weight of the soil in the slope failure block.

The formula for calculating the external power generated by the self-weight of the soil in the slope is as follows:

(24)

where represents the weight density of the soil in the slope; represents the linear velocity of each point in the failure block; represents the volume of the slope failure block. For the calculation of for the 3D two-step slope with a crack, the slope failure block should include both the 3D failure part and the planar inserted part:

(25)

where and respectively represent the external power of the 3D failure part and the planar inserted part of the slope. In the calculation, special attention should be paid to and : in slopes containing a vertical crack, the block AFG does not perform work, so the contribution of this part to the external power should be ignored. The calculation formulas for and are in the Appendix.

The internal energy dissipation rate of the slope can be calculated by the following equation [37]:

(26)

where represents the lower logarithmic spiral surface of the failure mechanism.

The type of crack formation determines the dissipation rate along the crack surface. Cracks can be classified into open cracks and formation cracks. Since open cracks form before slope failure occurs, they do not contribute to the total internal energy dissipation, i.e.,

(27)

Formation cracks form during the slope failure process, and their contribution to the internal dissipation rate can be represented as [1]:

(28)

where S_c represents the crack surface.

In the calculation of internal energy dissipation for the unsaturated soil slope with a crack, the effects of effective cohesion and apparent cohesion need to be considered separately. This can be computed as the sum of the dissipation rates along the sliding surface GE of the slope and along the crack surface FG, that is, . The dissipation rate along the sliding surface GE can be regarded as the overall contribution of the dissipation rates along the sliding surface AE, denoted as , and along the sliding surface AG, denoted as , that is, .It should also be noted that each of the internal energy dissipation terms mentioned above consists of two parts: 3D block dissipation term and the planar insert block dissipation term .

Firstly, calculate the internal energy dissipation rate generated by the effective cohesion of the soil:

(29)

where

(30)

(31)

For the formation crack,

(32)

For the open crack,

(33)

where (i=1, 2, 3, 4, 5, 6) represents the angle between the radius of the circle corresponding to the point on the failure surface in different segments and the horizontal radius of the circular section; ,,,, and are dimensionless parameters (see Appendix).

Then calculate the internal energy dissipation rate generated by the apparent cohesion . In this study, the groundwater level is located below the toe of the slope, so it is necessary to calculate the internal energy dissipation in each segment, similar to . The expression for is as follows:

(34)

where

(35)

(36)

For the formation crack,

(37)

For the open crack,

(38)

where ,,,,, are dimensionless parameters, and their detailed expressions are provided in the Appendix.

3.3 Maximum crack depth

Vertical cracks at the crest of a slope cannot be infinitely long; they must not exceed a certain limit. The depth of the crack depends on the stability of the vertical crack, which can be considered as a vertical slope itself. There is a maximum limit to the length of the crack; beyond this constraint, the crack will not remain stable. Under 2D conditions and at a constant temperature, the formula for calculating the maximum depth of cracks in a saturated soil slope is as follows [1]:

(39)

However, this study is based on a 3D failure mechanism and takes into account the temperature effects of unsaturated soil, so the aforementioned formula cannot be used directly. As previously mentioned, methods for solving the stability of 3D slopes can be applied to determine the maximum depth of vertical cracks. The location and depth of existing cracks are assumed to be uncertain. Among all possible cracks, there must be one that approaches the critical condition with the greatest impact on slope stability, and its value can be obtained by minimizing the critical slope height.

From the power balance equation, we can derive:

(40)

In this context, , and represent the soil gravity power, effective cohesion dissipation power, and apparent cohesion dissipation power of the vertical crack block AFG region, respectively. The derivation methods for the expressions of , and are similar to those described earlier, and thus the results are presented directly here:

(41)

(42)

(43)

Thus, the maximum crack depth can be calculated by the following formula:

(44)

Based on the upper bound theorem, the solution derived from Eq. (44) represents an upper limit of the actual solution. Through iterative optimization of the variables , and , one can identify the minimum upper bound value that meets the constraints, which is the maximum depth of the crack.

3.4 Safety factor

The collapse of a slope is a direct consequence of the interplay between inducing forces and resisting forces. Slope stability analysis typically employs two approaches: the first is the strength reduction method [1], which progressively diminishes resistance until the slope fails; the second is the gravity increase method [12], which steadily increases inducing forces until failure occurs. The deficiency of the gravity increase method is that it tends to slightly overestimate the minimum upper bound solution, but it can provide an explicit expression for the safety factor that is readily accepted by geotechnical engineers [12]. In contrast, the strength reduction method can only provide an implicit equation for safety factor that is not easily directly applied to engineering practice. Given that this study incorporates the temperature effects on unsaturated soil and the 3D aspects of the slope, the resulting expressions for internal energy dissipation and external work become more intricate, necessitating a longer computation time. Consequently, the gravity increase method is preferred for calculating the safety factor due to its ability to provide explicit results.

In the gravity increase method, the safety factor is defined as:

(45)

where F_S represent the safety factor; and respectively represent the unit weights of the soil in the failure and initial states. By substituting the calculated external power and internal energy dissipation rate into Eq. (45), the safety factor can be determined. It is important to note that the F_S values obtained using the limit analysis method are upper bounds of the actual solution; hence, one must identify the minimum value among all possible solutions as the final result. In this study, as illustrate in Figure 4, an optimization scheme was developed using the MATLAB software platform, defining four optimization variables ,, and for the crack failure mechanism, taking into account the temperature effect on unsaturated soil. The specific optimization steps are as follows:

Figure 4全屏查看

Flowchart of the optimization program

1) First, input the initial parameters, including 3D two-bench slope parameters, crack types and temperature and soil parameters:

2) Construct the 3D horn mechanism with cracks to calculate the external power, internal energy dissipation rate, and maximum crack depth. The objective function can be expressed as:

(46)

3) Determine the initial values of the four optimization variables ,, and based on the constraints given in Eq. (22).

4) Based on the particle swarm optimization (PSO) algorithm [38], perform a global optimization search for the minimum upper bound solution of the objective function to identify the optimal outcome.

5.1 The impact of 3D effect on two-step slope with different crack conditions

Figure 5 illustrates the trend of stability variation for stepped slopes under different aspect ratios. The results indicate that considering the 3D effects is essential for slope stability: at the same horizontal coordinate, when , the corresponding safety factor is the highest; the larger the value of , the lower the corresponding safety factor; when (corresponding to a slope undergoing 2D failure), the safety factor reaches its minimum value. For ease of understanding, let’s take as an example to discuss how the safety factor of the slope changes with the formation of cracks as changes: when increases from 1 to 2, the safety factor decreases by 21.79%; when increases from 2 to 5, the safety factor decreases by 12.35%; when increases from 5 to 10, the safety factor decreases by 4.58%; and when approaches infinity (indicating a slope with 2D failure), the safety factor decreases by 5.57%. Additionally, it can be observed from the Figure 5 that a slope with no cracks has the best stability among the three, a slope with formation cracks has intermediate stability, and a slope with open cracks has the worst stability. This also demonstrates that considering the presence of cracks significantly reduces slope stability.

5.2 Impact of slope parameters on two-step slope with different crack conditions

Figure 6 presents the impact of the four control parameters of two-step slopes on slope stability: the upper slope angle , the lower slope angle , the depth coefficient , and the step width coefficient . From Figure 6, it can be observed that the stability of the slope deteriorates as the slope angle increases, and the lower slope angle has a greater impact on slope stability compared to the upper slope angle . Considering the instance of formation cracks with the horizontal coordinate , Figure 6(a) illustrates that an increase in the upper slope angle from 45° to 55° results in a 10.27% reduction in the safety factor. Similarly, in Figure 6(b), raising the lower slope angle from 45° to 55° leads to a more significant 24.22% decrease in safety. As illustrated in Figures 6(c) and (d), enhancements in the depth coefficient and the step width coefficient are both beneficial for bolstering slope stability. The influence of cracks on these parameters is generally similar: the presence of open cracks is the most detrimental, while no cracks is optimal. This insight suggests that in the engineering design and reinforcement of stepped slopes, stability can be enhanced through strategies such as sealing existing cracks, expanding step surfaces, and increasing slope height.

5.3 Impact of temperature effects on various soils under different crack conditions and matrix suctions

To explore the pattern of unsaturated soil slope stability changes with temperature under different matrix suctions, this study selected four typical unsaturated soils that have been previously studied as research subjects and used a temperature range close to actual engineering projects (-5°T35°), with results shown in Figures 7-10. The stability patterns of slopes under the three crack conditions align with the conclusions from previous chapters: slopes exhibit maximum stability in the absence of cracks and are at their most precarious when open cracks are present. Furthermore, Figures 7-10 demonstrate that higher matrix suction leads to improved stability for slopes composed of sandy soil.

Figure 7 illustrates the trend of safety factor variation for sand slopes under three types of crack conditions as a function of temperature. It is evident that for saturated sands, where matrix suction is zero, the safety factor of the slope does not change with temperature. For unsaturated sands, as temperature gradually increases, the safety factor of the slope decreases in all crack conditions, and this decrease is approximately linearly related to temperature: taking the case of forming cracks with a matrix suction of 120 kPa as an example: when T rises from -5°C to 5°C, the F_S decreases by 1.11%; to 15°C it decreases by 4.40%; to 25°C by 7.05%; to 35°C by 8.67%. Additionally, Figure 7 indicates that an increase in matrix suction significantly enhances the stability of sand slopes.

Figure 8 depicts the trend of the safety factor for silt slopes as it changes with temperature. For saturated silts, where matrix suction is zero, the slope stability is at a relatively low level and is unaffected by temperature. In the case of unsaturated silts, the presence of matrix suction enhances slope stability, but this enhancement does not significantly continue to improve with further increases in matrix suction: for instance, at 25 ℃ with formation cracks, an increase in matrix suction from 0 to 40 kPa results in a 10.56% increase in F_S; however, an increase from 40 to 160 kPa only leads to a 1.76% increase in F_S, which is significantly less than the initial increase. Moreover, the stability of unsaturated silt slopes does not exhibit significant changes with rising temperature, only minor fluctuations, suggesting that the impact of temperature can be considered negligible.

Figure 9 illustrates the trend of the safety factor for clay slopes as it varies with temperature. For saturated clay, similar to saturated silt, the slope stability is lower compared to the unsaturated condition. However, for unsaturated clay, its slope stability is the most distinctive among the four typical soils, both with changes in matrix suction and with temperature variations. As can be seen from Figure 9, as the temperature rises, matrix suctions of 80 and 120 kPa correspond to curves that first increase and then decrease, exhibiting a peak value. A matrix suction of 40 kPa corresponds to a monotonically decreasing curve, while a matrix suction of 160 kPa corresponds to a monotonically increasing curve. The overall trend of these two curves is also initially increasing and then decreasing, but since the peak occurs outside the horizontal axis range, they appear to be monotonic when viewed on the graph. This indicates that for unsaturated clay slopes, there often exists a maximum safety factor at a certain temperature when changing with temperature. The greater the absolute difference between the actual temperature and this temperature, the lower the stability of the slope.

Figure 10 investigates the stability of extremely fine-grained soil slopes at different temperatures. Compared to the other three types of unsaturated soils, unsaturated extremely fine-grained soils offers the greatest improvement in slope stability: taking the case of formation cracks at T=25 ℃ as an example, when matrix suction increases from 0 to 160 kPa, the F_S for sand slopes increases by 21.79%, for silt slopes by 11.48%, for clay slopes by 53.33%, and remarkably, for extremely fine-grained soil slopes by 279.03%. Under the three types of crack conditions, the stability of extremely fine-grained soil slopes decreases with rising temperature, exhibiting an overall near-linear trend. Similar to sand slopes, enhancing the matrix suction of unsaturated extremely fine-grained soil slopes can significantly improve their stability.

Overall, the finer the soil particles, the greater the F_S; as temperature increases, the F_S of sand and extremely fine-grained soil slopes decreases, with a trend that is approximately linear, the F_S of silt slopes remains almost unchanged, and the F_S of clay slopes first increases and then decreases. The trend of unsaturated soil slope stability with temperature is not uniform; different soils have their unique temperature variation curves. This underscores the importance of understanding soil types in analyzing slope stability.