Displacement Analyses for a Natural Slope Considering Post-Peak Strength of Soils

: A natural slope undergoing recurrent movements caused by rainfall-induced groundwater table rises is studied using a novel method. The strength and displacement parameters are back-calculated using a force-equilibrium-based ﬁnite displacement method (FFDM) based on the ﬁrst event of slope movement recorded in the monitoring period. Slope displacements in response to subsequent rainfall-induced groundwater table rises are predicted using FFDM based on the back-calculated material parameters. Important factors that may inﬂuence the accuracy of slope displacement predictions, namely, the curvature of the Mohr-Coulomb (M-C) failure envelope and post-peak strength softening, are investigated. It is found that the accuracy of slope displacement predictions can be improved by taking into account post-peak stress-displacement relationship in the analysis. The accuracy of slope displacement predictions is not inﬂuenced by the curvature of the M-C failure envelope in the displacement analysis.


Introduction
Intensive rainfall causes geohazards in the form of slope failures and debris flows such as those reported by Huang et al. [1][2][3][4][5]. To investigate the failure mechanism of slopes, limit equilibrium methods are widely used in which a single value of input soil strength (called operational strength or weighted strength) [6,7] is used to calculate a safety factor (Fs) for the studied slopes [8,9]. The calculated value of Fs reflects the safety margin of the slope against failure, providing no information on the deformation (or displacements) of the slopes. To remedy this shortcoming, Huang et al. [10,11] proposed a force-equilibrium-based finite displacement method (FFDM) that incorporates a nonlinear (hyperbolic) stress-displacement constitutive law for the geomaterial along the sliding surface, as shown in Figure 1. Hyperbolic curves are widely used for simulating the stressstrain and stress-displacement relationships of soils and the interface between soil and other materials [12][13][14], with the drawback that post-peak strength deterioration for some soils, as shown in Figure 2, cannot be properly accounted for. Tatsuoka et al. [14] showed the effectiveness of utilizing generalized hyperbolic equations and an exponential function for simulating the post-peak strength softening of sandy soils. However, a general description of the post-peak stress-strain (or stress-displacement) relationship is currently unavailable. A critical issue in slope stability analyses is the selection of the operation strength, i.e., the selection of a representative value of the friction angle among the peak, fully softened, and residual friction angles for the geomaterial comprising the sliding surface [15,16]. A conventional limit-equilibrium-based approach produces a back-calculated friction angle based on a known (or assumed) value of the cohesion intercept (c) under the condition of Fs = 1.0 (e.g., Huang) [17]. The result of back-analysis is not affected by the extent of slope deformation (or displacement). To correctly reflect the displacement status of

Methodology
Observing the slope failure in natural slope and a series of rainfall tests, the mechanism of slope displacement of rainfall-induced groundwater table rises will be deserving of studying and discussing.

Constitutive Law for Stress-Displacement Relationship
The following hyperbolic equation is used for the normalized shear stress (τ/τf) vs. shear displacement (Δ) relationship along a potential failure surface ( approach that uses FFDM for slopes undergoing periodic slope movements induced by rainfall (or groundwater table rises). Both strength-and displacement-related parameters back-calculated from the first event of slope displacement were found to be valid for predicting slope displacements that occurred in subsequent rainfall events (or groundwater table rises). In their studies, linear Mohr-Coulomb (M-C) envelopes and hyperbolic stress-displacement relationship were used; however, this may be insufficient to properly simulate the behavior of soils. The present study aims to solve the above problem by using curved M-C envelopes and the hyperbolic stress-displacement relationship with post-peak softening.
Various hillside disasters, not only for social and economic factors, pay a heavy price, but also endanger the ecological environment. This study is proposing, as a reference in practical application and early warning work, and it is also expected that the slope disasters in the future can be significantly reduced. Results of this study preliminarily show that there is a good agreement between the calculated and the measured displacements in the studied slopes. However, obtaining more observed data on wellmonitored slopes is needed so that the technique can be further verified.

Methodology
Observing the slope failure in natural slope and a series of rainfall tests, the mechanism of slope displacement of rainfall-induced groundwater table rises will be deserving of studying and discussing.

Constitutive Law for Stress-Displacement Relationship
The following hyperbolic equation is used for the normalized shear stress (τ/τf) vs. shear displacement (Δ) relationship along a potential failure surface ( Figure 1): Various hillside disasters, not only for social and economic factors, pay a heavy price, but also endanger the ecological environment. This study is proposing, as a reference in practical application and early warning work, and it is also expected that the slope disasters in the future can be significantly reduced. Results of this study preliminarily show that there is a good agreement between the calculated and the measured displacements in the studied slopes. However, obtaining more observed data on well-monitored slopes is needed so that the technique can be further verified.

Methodology
Observing the slope failure in natural slope and a series of rainfall tests, the mechanism of slope displacement of rainfall-induced groundwater table rises will be deserving of studying and discussing.

Constitutive Law for Stress-Displacement Relationship
The following hyperbolic equation is used for the normalized shear stress (τ/τ f ) vs. shear displacement (∆) relationship along a potential failure surface ( Figure 1): where k initial : initial shear stiffness R f : ratio between failure strength and asymptotic shear strength τ ult : asymptotic strength at infinite displacement τ f : failure shear strength expressed using the M-C failure criterion The initial shear stiffness can be expressed as a power function of effective normal pressure (σ' n ) on the failure surface: where K: stiffness number (an experimental constant) P A : atmospheric pressure (=101.3 kPa) G: reference shear stiffness (=101.3 kPa/m) Equations (1)-(6) are incorporated into FFDM with additional force and/or moment equilibrium equations and displacement compatibility functions for a vertically sliced failure mass to calculate the shear displacement along the failure surface (discussed below). Figure 3a shows the failure envelopes for two typical sandy soils reported by Wu et al. [18,19]; these results were obtained using a medium-scale direct shear test apparatus. These failure envelopes are curved rather than straight. The following two methods are equally effective for handling the nonlinearity of the M-C failure envelope: (1) expressing the shear stress (τ) vs. normal stress (σ n ) relationship using a nonlinear equation [20] and (2) using a pressure-level-dependent internal friction angle, ϕ [21]. The latter approach is used in the present study. Figure 3b shows the ϕ vs. log (σ n ') relationship for the data shown in Figure 3a. The straight lines can be expressed using the following equations:

Post-Peak Stress-Displacement Relationship
The post-peak segment of the τ vs. Δ curve is simulated in the following using a curve (also called "Versoria" or "the witch of Agnesi") proposed by Grandi in the 1700s [22], shown in Figure 2. To use this curve efficiently as part of the stress-displacement curve, a normalized local coordinate system (X-Y coordinates) is used. The axis of X = 0 (the ordinate) passes through the point of peak stress (τf). The axis of Y = 0 (the abscissa) is an asymptote of the residual stress τr (= (1-t)·τf). A feature of this curve is that a residual state is attained at an infinite value of Δ. In practice, however, an experimental value of Δ (namely, Δr at X = 1) can be used to mimic the value of τr at a finite displacement Δr with a negligibly small error of τr. In general, the post-peak shear stress (τ) can be expressed as: where t: normalized strength deterioration from peak to residual states Y: normalized post-peak shear stress X: normalized post-peak shear displacement parameter Δf: shear displacement at peak stress states Values of t for various soils derived based on the reported stress-displacement curves [18,19,[23][24][25][26][27][28] are summarized in Figure 4a. The values of t vary over a range of 0.1 to 0.5 for various soils. For a specific soil, values of t can be expressed as a linear function of confining stress, σn. Among the test results shown in Figure 4a, those for two direct shear tests conducted in regions of central Taiwan, reported by Wu [18] and Hsu [19] (see Figure At peak condition for the silty sand reported by Wu [18] and Hsu [19], ϕ b = 49.8 • and 45.1 • and ∆ϕ = 7.9 • and 9.1 • , respectively. At residual conditions for the silty sand reported by Wu [18] and Hsu [19], ϕ b = 46.3 • and 40.2 • and ∆ϕ = 6.9 • and 5.2 • , respectively.

Post-Peak Stress-Displacement Relationship
The post-peak segment of the τ vs. ∆ curve is simulated in the following using a curve (also called "Versoria" or "the witch of Agnesi") proposed by Grandi in the 1700s [22], shown in Figure 2. To use this curve efficiently as part of the stress-displacement curve, a normalized local coordinate system (X-Y coordinates) is used. The axis of X = 0 (the ordinate) passes through the point of peak stress (τ f ). The axis of Y = 0 (the abscissa) is an asymptote of the residual stress τ r (=(1−t)·τ f ). A feature of this curve is that a residual state is attained at an infinite value of ∆. In practice, however, an experimental value of ∆ (namely, ∆r at X = 1) can be used to mimic the value of τ r at a finite displacement ∆r with a negligibly small error of τ r . In general, the post-peak shear stress (τ) can be expressed as: where t: normalized strength deterioration from peak to residual states Y: normalized post-peak shear stress X: normalized post-peak shear displacement parameter ∆ f : shear displacement at peak stress states Values of t for various soils derived based on the reported stress-displacement curves [18,19,[23][24][25][26][27][28] are summarized in Figure 4a. The values of t vary over a range of 0.1 to 0.5 for various soils. For a specific soil, values of t can be expressed as a linear function of confining stress, σ n . Among the test results shown in Figure 4a, those for two direct shear tests conducted in regions of central Taiwan, reported by Wu [18] and Hsu [19] (see  Based on the results of medium-or large-scale direct shear tests on various granular materials [19,[23][24][25][26][27][28], the shear displacement at residual states (Δr) can be expressed as a function of the shear displacement at the peak strength state (Δf) as: where Δratio is the ratio between the shear displacement at the residual state and that at the peak state. Values of Δratio obtained in various studies are shown in Figure 5. As can be seen, Δratio can be assumed as a linear function of σn. For a specific soil, Δratio can be expressed as: where Δratio (100): value of Δratio under reference confining pressure of σn = 100 kPa R: experimental parameter (ranging from 0 to 0.013) Based on the results of medium-or large-scale direct shear tests on various granular materials [19,[23][24][25][26][27][28], the shear displacement at residual states (∆ r ) can be expressed as a function of the shear displacement at the peak strength state (∆ f ) as: where ∆ ratio is the ratio between the shear displacement at the residual state and that at the peak state. Values of ∆ ratio obtained in various studies are shown in Figure 5. As can be seen, ∆ ratio can be assumed as a linear function of σ n . For a specific soil, ∆ ratio can be expressed as: where ∆ ratio(100) : value of ∆ ratio under reference confining pressure of σ n = 100 kPa R: experimental parameter (ranging from 0 to 0.013) Figures 6 and 7 compare experimental and simulated shear stress vs. displacement curves for the tests reported by Wu [18] and Hsu [19], respectively. The input parameters used for simulating the stress-displacement curves in these figures are summarized in Tables 1 and 2, respectively. It can be seen that the simulated pre-peak and post-peak curves agree well with the experimental ones.   [18] and Hsu [19], respectively. The input par used for simulating the stress-displacement curves in these figures are summa Tables 1 and 2, respectively. It can be seen that the simulated pre-peak and po curves agree well with the experimental ones.   [18] and Hsu [19], respectively. The input parameters used for simulating the stress-displacement curves in these figures are summarized in Tables 1 and 2, respectively. It can be seen that the simulated pre-peak and post-peak curves agree well with the experimental ones.

Displacement Compatibility Requirements
Figure 8a schematically shows an example of a potential failure surface with seven vertical slices, for which vectors of shear displacements at the slice base and the slice interface are shown in Figure 8b. Based on the hodograph shown in Figure 8c, a general expression for the shear displacement Δi (i = 1-7) at the base of slice i can be expressed as: where Δ0: vertical displacement at the crest of slice 1; α1, αi: inclination angles for slices 1 and i, respectively

Displacement Compatibility Requirements
Figure 8a schematically shows an example of a potential failure surface with seven vertical slices, for which vectors of shear displacements at the slice base and the slice interface are shown in Figure 8b. Based on the hodograph shown in Figure 8c, a general expression for the shear displacement ∆ i (i = 1-7) at the base of slice i can be expressed as: where ∆ 0 : vertical displacement at the crest of slice 1; α 1 , α i : inclination angles for slices 1 and i, respectively Ψ: angle of dilation.
Ψ: angle of dilation. Figure 8b, c are for the case of Ψ > 0°, representing a volume expansion condition at the pre-peak state for dense soils. It can be seen that Δi increases from the crest of the slope towards the toe. Figure 9a,b show the displacement vectors and hodograph, respectively, for the case of Ψ < 0°, representing a volume contraction condition for loose soils. It can be seen that the shear displacement decreases from the crest of the slope toward the toe. The displacement vectors and hodograph for a special case of Ψ = 0° is shown in Figure 10a,b, representing a non-volume change condition at large shear displacements along the failure surface. A feature for the condition of Ψ = 0° is equal horizontal components for Δi (i = 1-7).   8b,c are for the case of Ψ > 0 • , representing a volume expansion condition at the pre-peak state for dense soils. It can be seen that ∆ i increases from the crest of the slope towards the toe. Figure 9a,b show the displacement vectors and hodograph, respectively, for the case of Ψ < 0 • , representing a volume contraction condition for loose soils. It can be seen that the shear displacement decreases from the crest of the slope toward the toe. The displacement vectors and hodograph for a special case of Ψ = 0 • is shown in Figure 10a,b, representing a non-volume change condition at large shear displacements along the failure surface. A feature for the condition of Ψ = 0 • is equal horizontal components for ∆ i (i = 1-7). GeoHazards 2021, 2, x FOR PEER REVIEW 9 of 22

Calculating Safety Factors Using Method of Slices
A conventional limit equilibrium method of slices, namely Janbu's generalized procedure of slices [29], which satisfies force and moment equilibrium requirements, is used to calculate Fs for the studied slope. For a potential sliding surface, shown in Figure  11, the differential horizontal inter-slice forces (ΔEi) can be formulated based on the vertical and horizontal force equilibrium: where Wi: self-weight of slice i (i = 1−n)

Calculating Safety Factors Using Method of Slices
A conventional limit equilibrium method of slices, namely Janbu's generalized procedure of slices [29], which satisfies force and moment equilibrium requirements, is used to calculate Fs for the studied slope. For a potential sliding surface, shown in Figure  11, the differential horizontal inter-slice forces (ΔEi) can be formulated based on the vertical and horizontal force equilibrium: where Wi: self-weight of slice i (i = 1−n)

Calculating Safety Factors Using Method of Slices
A conventional limit equilibrium method of slices, namely Janbu's generalized procedure of slices [29], which satisfies force and moment equilibrium requirements, is used to calculate Fs for the studied slope. For a potential sliding surface, shown in Figure 11, the differential horizontal inter-slice forces (∆E i ) can be formulated based on the vertical and horizontal force equilibrium: where W i : self-weight of slice i (i = 1−n) α i : base inclination of slice i ∆X i : differential vertical inter-slice force for slice i (=X i −X i−1 ) S fi : ultimate shear strength at the base of slice i based on the M-C failure criterion.
where Ci: cohesive resistance at slice i (= ci·li; li: length of the base of slice i).
A summation of Equation (10) yields the following equation: where Eo, En: horizontal boundary forces at slices 0 and n, respectively. By combining Equations (16) and (19), a constant value of Fs can be derived as: Figure 11. Forces acting on vertical slice of potential failure mass (Janbu's slice method).
Based on the moment equilibrium about the center of the slice base and eliminating the differential terms of Ei and Xi, which become negligibly small when slice widths are sufficiently small, one can obtain: To calculate a constant value of Fs for a potential sliding mass with an arbitrary shape of the shear surface, initial values of Fs = 1.0 and ΔXi = 0 (i = 1−n) were used to derive approximate values of Sfi (Equation (18)) and ΔEi (Equation (16)). Improved values of Fs can be calculated (Equation (20)) using values of ΔXi (based on values of Xi from Equation (21)) until Fs converges with a tolerable convergence error of 0.5%.

Force-Equilibrium-Based Finite Displacement Method
The modified slice method introduced above is incorporated into FFDM by using a local factor of safety (Fsi) in lieu of a constant value of Fs in Equations (16), (18), and (20). Equation (5), expressed by the following equation: where C i : cohesive resistance at slice i (=c i ·l i ; l i : length of the base of slice i).
A summation of Equation (10) yields the following equation: where E o , E n : horizontal boundary forces at slices 0 and n, respectively. By combining Equations (16) and (19), a constant value of F s can be derived as: Based on the moment equilibrium about the center of the slice base and eliminating the differential terms of E i and X i , which become negligibly small when slice widths are sufficiently small, one can obtain: To calculate a constant value of F s for a potential sliding mass with an arbitrary shape of the shear surface, initial values of F s = 1.0 and ∆X i = 0 (i = 1−n) were used to derive approximate values of S fi (Equation (18)) and ∆E i (Equation (16)). Improved values of F s can be calculated (Equation (20)) using values of ∆X i (based on values of X i from Equation (21)) until F s converges with a tolerable convergence error of 0.5%.

Force-Equilibrium-Based Finite Displacement Method
The modified slice method introduced above is incorporated into FFDM by using a local factor of safety (F si ) in lieu of a constant value of F s in Equations (16), (18), and (20). By using local factors of safety (F si ) instead of a constant one (F s ), the shear stress vs. shear displacement relationship shown in Equations (1)- (6), and the displacement compatibility function shown in Equations (14) and (15), the force system for the sliced failure mass is in a static indeterminate condition [10]. Note that Equation (1) is a reversed form of F si : or where τ fi , τ i : ultimate (or peak) shear strength and shear stress, respectively, at the base of slice i (i = 1−n) S fi , S i : ultimate (or peak) shear resistance and shear force, respectively, at the base of slice i l i : length of the base of slice i.
The Equation for S fi (Equation (18)) can be rewritten based on the M-C failure criterion as follows: Rewriting Equation (16) to take into account the local factor of safety (F si ) yields: By substituting Equations (14) and (23) into Equation (25) and summing up ∆E i (i = 1−n; Equation (17)), a closed form solution for ∆ 0 can be obtained: Equation (26) can be rewritten by incorporating the M-C failure criterion, i.e., substituting Equation (24) into Equation (26): Note that Equation (27) is a nonlinear equation, requiring a few iterative computations to derive a converged value of ∆ 0 . To calculate ∆ 0 , a small value of ∆ 0 (=0.001 m) is first used on the right-hand side of Equation (27); improved values of F si , ∆E i , and ∆X i are subsequently calculated using Equations (22), (25), and (21), respectively. This process is iterated until ∆ 0 converges. In the following analyses, an allowable convergence error of 2% is used for deriving a final value of ∆ 0 .
The incremental slope displacement at the base of slice No. i (∆ incr i ) induced by a rainfall-induced groundwater table rise (or a porewater pressure increase) can be expressed as the difference between the post-rainfall slope displacement (∆ b i ) and the pre-rainfall slope displacement (∆ a i ):

Case Study: Lu-Shan Slope
The studied slope is located between km markers 44 and 91 along Provincial Highway 14, which winds through a landslide-prone area ( Figure 12). The sliding mass has maximum width, length and depth of 480 m, 820 m and 100 m, respectively, which spans over a height of about 400 m (between elevations of 1085 and 1500 m above sea level) and covers an area of about 30 ha. In situ borehole logging, in situ testing, laboratory testing, and groundwater table measurements were conducted during the monitoring period. The major rock formation in this area is slate (SL) with intensive joints and is covered by alluvium deposits with depths of a few centimeters to tens of meters. Cross-section C-C' is approximately in the direction of landslides, along which boreholes B01, B03, B06 and B09 were drilled. Groundwater tables were measured at boreholes B03, B06 and B09 using electronic groundwater pressure sensors. Slope displacements were measured at borehole B09 using an in-hole inclinometer. Slope displacements at boreholes B03 and B06 were based on Global Positioning System (GPS) data obtained at stations B21 and A25, respectively, which are located in the vicinity of the boreholes.
testing, and groundwater table measurements were conducted during the monitoring period. The major rock formation in this area is slate (SL) with intensive joints and is covered by alluvium deposits with depths of a few centimeters to tens of meters. Crosssection C-C' is approximately in the direction of landslides, along which boreholes B01, B03, B06 and B09 were drilled. Groundwater tables were measured at boreholes B03, B06 and B09 using electronic groundwater pressure sensors. Slope displacements were measured at borehole B09 using an in-hole inclinometer. Slope displacements at boreholes B03 and B06 were based on Global Positioning System (GPS) data obtained at stations B21 and A25, respectively, which are located in the vicinity of the boreholes. Figure 13 shows cross-section C-C' with a detected failure surface (bold solid line). The failure surface reported by Huang et al. [11] is also plotted (dotted line). In the present study, a new failure surface based on time-domain reflectometry (TDR, as shown in Figure 13) and GPS data is drawn; this surface was not considered in the study of Huang et al. [11]. The slip surface is located approximately at the interface of colluviums and highly fractured (or weathered) slate strata with a maximum depth of about 100 m, associated with rock quality designation (RQD) values distributed over a wide range. Deriving the material properties using an undisturbed sampling technique for such a depth is costly and time-consuming. Figure 13 shows the cross-section C-C' of the Lu-Shan slope with a deep-seated slip surface approximately along the interface of colluvium and heavily weathered slate (SL) interlayered with sandy slate (SSL), as shown in the cross-section. An example of the rock quality designation (RQD) profile obtained from borehole B09 is also shown in Figure 13. The slip surface shown in Figure 13 is identical with that reported in a monitoring and remediation project supervised by the governmental authority [17].   Figure 13 shows cross-section C-C' with a detected failure surface (bold solid line). The failure surface reported by Huang et al. [11] is also plotted (dotted line). In the present study, a new failure surface based on time-domain reflectometry (TDR, as shown in Figure 13) and GPS data is drawn; this surface was not considered in the study of Huang et al. [11]. The slip surface is located approximately at the interface of colluviums and highly fractured (or weathered) slate strata with a maximum depth of about 100 m, associated with rock quality designation (RQD) values distributed over a wide range. Deriving the material properties using an undisturbed sampling technique for such a depth is costly and timeconsuming. Figure 13 shows the cross-section C-C' of the Lu-Shan slope with a deep-seated slip surface approximately along the interface of colluvium and heavily weathered slate (SL) inter-layered with sandy slate (SSL), as shown in the cross-section. An example of the rock quality designation (RQD) profile obtained from borehole B09 is also shown in Figure 13. The slip surface shown in Figure 13 is identical with that reported in a monitoring and remediation project supervised by the governmental authority [17].    Table 3.   Table 3.    Table 3.   Figure 15 shows a comparison of the measured groundwater table heights at various locations of the failure mass. It can be seen that all groundwater tables respond to the rainfall event with abrupt rises at the central (borehole B06) and near-top (borehole B03) parts of the slope. The groundwater table heights near the toe of the slope (borehole B09) were a higher level than those measured at the upper portions of the slope (boreholes B03 and B06), suggesting that a perched water table exists around the toe of the slope, undermining the stability of the slope.    For most of the rainfall events, measured slope displacements at various locations are similar, except those for the last rainfall event. For the last event of rainfall (or slope movements), the displacements at the upper portions (boreholes B03 and B06) were greater than those around the toe (borehole B09), suggesting that a displacement field with a decreasing shear displacement toward the toe of the slope (see Figure 9a,b) is applicable here. This pattern of displacements also suggests that a volumetric contraction occurred within the sliding mass, which acted as a buffer in transmitting lateral thrust from the upper portion of the sliding mass.  Figure 16 compares the measured displacements at various locations of the sliding mass. For most of the rainfall events, measured slope displacements at various locations are similar, except those for the last rainfall event. For the last event of rainfall (or slope movements), the displacements at the upper portions (boreholes B03 and B06) were greater than those around the toe (borehole B09), suggesting that a displacement field with a decreasing shear displacement toward the toe of the slope (see Figure 9a,b) is applicable here. This pattern of displacements also suggests that a volumetric contraction occurred within the sliding mass, which acted as a buffer in transmitting lateral thrust from the upper portion of the sliding mass. are similar, except those for the last rainfall event. For the last event of rainfall (or slope movements), the displacements at the upper portions (boreholes B03 and B06) were greater than those around the toe (borehole B09), suggesting that a displacement field with a decreasing shear displacement toward the toe of the slope (see Figure 9a,b) is applicable here. This pattern of displacements also suggests that a volumetric contraction occurred within the sliding mass, which acted as a buffer in transmitting lateral thrust from the upper portion of the sliding mass.   (20)) calculated based on the groundwater tables measured after rainfall events nos. 1-8. In this analysis, cohesion intercept (c) = 65 kPa and ϕ = 32 • were used based on the hypothesis that the slope attained the critical condition (F s = 1.0) during rainfall event no. 1 with a high groundwater table. The value of c = 65 kPa was determined based on the following empirical equation [30,31]: where c: cohesion intercept (kPa) z: depth (m).
GeoHazards 2021, 2, x FOR PEER REVIEW 15 of 22 Figure 16. Comparison of cumulative slope displacements measured at various boreholes. Figure 17 shows variations of Fs (Equation (20)) calculated based on the groundwater tables measured after rainfall events nos. 1-8. In this analysis, cohesion intercept (c) = 65 kPa and φ = 32° were used based on the hypothesis that the slope attained the critical condition (Fs = 1.0) during rainfall event no. 1 with a high groundwater table. The value of c = 65 kPa was determined based on the following empirical equation [30,31]: where c: cohesion intercept (kPa) z: depth (m). The average depth of the slip surface z ≒ 65 m (Figure 13), resulting in c = 65 kPa in the analysis. Note that an arbitrary value of c can be used to derive a calibrated value of φ based on the measured groundwater table height in event no. 1. The essence of the proposed procedure is that accurate predicted slope displacements for subsequent rainfall events (event nos. 2-8) can be obtained regardless of the c value used (to be explained in Figure 18). It can be seen in Figure 17 that values of Fs decrease rapidly during each rainfall event. An increasing trend of Fs in response to decreasing heights of the groundwater table during the period August 2005 to April 2006 can also be seen. In general, the increasing Fs corresponded to the trend of decreasing groundwater table heights during event nos. 1 to 4. However, this increasing trend of Fs (from event nos. 1 to 4) contradicts to the trend of slope displacements (see Figure 13). Therefore, an immediate rainfall-induced groundwater table rise may have a dominant influence on slope movement (and cumulative slope displacement). The influence of a long-term descending groundwater table on (cumulative) slope movement is insignificant. A misleading conclusion may be obtained if only values of Fs are used for slope stability evaluations.     (20)) calculated based on the groundwater tables measured after rainfall events nos. 1-8. In this analysis, cohesion intercept (c) = 65 kPa and φ = 32° were used based on the hypothesis that the slope attained the critical condition (Fs = 1.0) during rainfall event no. 1 with a high groundwater table. The value of c = 65 kPa was determined based on the following empirical equation [30,31]: where c: cohesion intercept (kPa) z: depth (m). The average depth of the slip surface z ≒ 65 m (Figure 13), resulting in c = 65 kPa in the analysis. Note that an arbitrary value of c can be used to derive a calibrated value of φ based on the measured groundwater table height in event no. 1. The essence of the proposed procedure is that accurate predicted slope displacements for subsequent rainfall events (event nos. 2-8) can be obtained regardless of the c value used (to be explained in Figure 18). It can be seen in Figure 17 that values of Fs decrease rapidly during each rainfall event. An increasing trend of Fs in response to decreasing heights of the groundwater table during the period August 2005 to April 2006 can also be seen. In general, the increasing Fs corresponded to the trend of decreasing groundwater table 65 m (Figure 13), resulting in c = 65 kPa in the analysis. Note that an arbitrary value of c can be used to derive a calibrated value of ϕ based on the measured groundwater table height in event no. 1. The essence of the proposed procedure is that accurate predicted slope displacements for subsequent rainfall events (event nos. 2-8) can be obtained regardless of the c value used (to be explained in Figure 18). It can be seen in Figure 17 that values of Fs decrease rapidly during each rainfall event. An increasing trend of F s in response to decreasing heights of the groundwater table during the period August 2005 to April 2006 can also be seen. In general, the increasing F s corresponded to the trend of decreasing groundwater table heights during event nos. 1 to 4. However, this increasing trend of F s (from event nos. 1 to 4) contradicts to the trend of slope displacements (see Figure 13). Therefore, an immediate rainfallinduced groundwater table rise may have a dominant influence on slope movement (and cumulative slope displacement). The influence of a long-term descending groundwater table on (cumulative) slope movement is insignificant. A misleading conclusion may be obtained if only values of F s are used for slope stability evaluations.
assumptions are investigated later. In practice, an infinite number of curves o m) vs. K relationship, as shown in Figure 19, can be used. However, only three are shown in Figure 19 for simplicity. It can be seen that these curves, which use of Rf and n, yield different back-calculated values of K that meet the requirem 633 m) = 0.064 m. The strength parameters, c = 65 kPa, φ = 32° and Ψ = 0°, obta 18 and the displacement parameters, Rf = 0.75, n = 0.1 and K = 200, obtained i used to calculate the incremental displacements in rainfall events 2-8. cumulative displacements for all rainfall events are shown in Figure 20a,c. Th curves of slope displacements agree well with the measured ones in the follow The calculated cumulative displacements are comparable to the observ displacements; (2) The use of Ψ = −12° for constructing the displacement calculated displacements that well agree with the slope displacements obser locations; (3) Three sets of K, n, and Rf generate similar cumulative slope displa   Figure 18 shows the results of a parametric study of the F s vs. ϕ relationship for various values of c (=0 and 65 kPa) and ∆ϕ (=0 • , 1 • , 3 • , and 5 • ) using Janbu's method of slices. The value of ∆ϕ is changed to investigate its effect (or that of the curvature of the M-C failure envelope) on the value of F s . The condition of ∆ϕ = 0 • is a special case of a linear M-C failure envelope. Figure 18 also shows that different values of ∆ϕ result in different values of ϕ that generate F s = 1.0. These back-calculated values of ϕ are used in the following analyses. Figure 19 shows the result of a back-analysis for displacement-related parameters (K, R f , and n in Equations (1)-(4) and (6)), obtained using c = 65 kPa, ϕ = 32 • and Ψ = 0 • . The use of Ψ = 0 • is based on the assumption that the slope is in a large displacement state before the occurrence of rainfall event no. 1 on 4 August 2005. This assumption is not necessarily true; other assumptions are investigated later. In practice, an infinite number of curves of the ∆ i (x = 633 m) vs. K relationship, as shown in Figure 19, can be used. However, only three typical curves are shown in Figure 19 for simplicity. It can be seen that these curves, which use various values of R f and n, yield different back-calculated values of K that meet the requirement of ∆ i (at x = 633 m) = 0.064 m. The strength parameters, c = 65 kPa, ϕ = 32 • and Ψ = 0 • , obtained in Figure 18

Influence of Post-Peak Softening on Slope Displacements
To consider post-peak softening behavior in the displacement analysis using FFDM, a back-analysis for the displacement parameter K based on rainfall event no. 1 similar to that performed in Figure 19 is shown in Figure 21 Figure 22 shows the results of displacement calculations that take into account post-peak softening of the slope material. A comparison of Figure 22 with Figures 20a-c reveals that post-peak softening has a remarkable influence on the calculated slope displacements for event no. 8, for which relatively large displacements were measured. The mechanism of post-peak softening was further investigated using the parameter set shown in Figure 22. The stress-displacement behavior at the base of slice no. 25 (central portion of the slope) and slice no. 49 (near the slope toe) was considered, as shown in Figure 23a. Figure 23b shows the stressdisplacement curves corresponding to pre-rainfall and post-rainfall stress states at the base of slice nos. 25 and 49 obtained without consideration of post-peak softening. The arrows in this figure show the transition of stress (and displacement) states from prerainfall to post-rainfall states. It can be seen that the soil at slice no. 25 (upper portion of the slope) is under the state of stress hardening due to the limitation of using a hyperbolic stress-displacement relationship. To take into account the influence of post-peak softening behavior of geomaterials on the rainfall-induced slope displacements, typical values of post-peak softening parameters obtained in Figures 4b and 5, namely, t1 = 0.00006, t2 = 0.2 and Δratio = 2, are used in the FFDM analyses. The results shown in Figure 23c suggest that the soil at the base of slice no. 25 undergoes a transition from a pre-peak state to a postpeak state during rainfall event no. 1, which was not found in previous analyses.

Influence of Post-Peak Softening on Slope Displacements
To consider post-peak softening behavior in the displacement analysis using FFDM, a back-analysis for the displacement parameter K based on rainfall event no. 1 similar to that performed in Figure 19 is shown in Figure 21. Note that infinite sets of parameters (c, ϕ), Ψ, (t 1 , t 2 , ∆ ratio ) and (K, n, R f ) can be used to the requirement of calculated ∆ i (at  Figure 22 shows the results of displacement calculations that take into account post-peak softening of the slope material. A comparison of Figure 22 with Figure 20a-c reveals that post-peak softening has a remarkable influence on the calculated slope displacements for event no. 8, for which relatively large displacements were measured. The mechanism of post-peak softening was further investigated using the parameter set shown in Figure 22. The stressdisplacement behavior at the base of slice no. 25 (central portion of the slope) and slice no. 49 (near the slope toe) was considered, as shown in Figure 23a. Figure 23b shows the stress-displacement curves corresponding to pre-rainfall and post-rainfall stress states at the base of slice nos. 25 and 49 obtained without consideration of post-peak softening. The arrows in this figure show the transition of stress (and displacement) states from pre-rainfall to post-rainfall states. It can be seen that the soil at slice no. 25 (upper portion of the slope) is under the state of stress hardening due to the limitation of using a hyperbolic stress-displacement relationship. To take into account the influence of post-peak softening behavior of geomaterials on the rainfall-induced slope displacements, typical values of post-peak softening parameters obtained in Figures 4b and 5, namely, t 1 = 0.00006, t 2 = 0.2 and ∆ ratio = 2, are used in the FFDM analyses. The results shown in Figure 23c suggest that the soil at the base of slice no. 25 undergoes a transition from a pre-peak state to a post-peak state during rainfall event no. 1, which was not found in previous analyses.

Conclusions
A natural slope subjected to recurrent rainfall-induced movements was studied using FFDM. In addition to the use of a hyperbolic model of the shear stress vs. shear displacement relationship and an improved displacement compatibility function, curved M-C failure envelopes and post-peak strength softening were taken into account in the slope displacement analyses. The results of the case study showed that: (1) The accuracy of slope displacement predictions was not influenced by the curvature of the M-C failure envelope used in the back-analyses and displacement analyses.
(2) Post-peak soil strength softening was successfully described by the proposed stressdisplacement relationship model. New parameters, including normalized post-peak strength deterioration (t) and residual displacement ratio (Δratio), were studied using direct shear tests on various soil types. Good agreements between the experimental and simulated post-peak stress vs. displacement relationships were obtained.
(3) The accuracy of the slope displacement computations can be improved by incorporating the post-peak stress-displacement relationship in the analysis. The results of a case study on the Lu-Shan slope revealed that a transition from pre-peak to post-peak states for the slip surface near the slope toe may have taken place during the first rainfall event in the monitoring period. This critical feature may be

Conclusions
A natural slope subjected to recurrent rainfall-induced movements was studied using FFDM. In addition to the use of a hyperbolic model of the shear stress vs. shear displacement relationship and an improved displacement compatibility function, curved M-C failure envelopes and post-peak strength softening were taken into account in the slope displacement analyses. The results of the case study showed that: (1) The accuracy of slope displacement predictions was not influenced by the curvature of the M-C failure envelope used in the back-analyses and displacement analyses. (2) Post-peak soil strength softening was successfully described by the proposed stressdisplacement relationship model. New parameters, including normalized post-peak strength deterioration (t) and residual displacement ratio (∆ ratio ), were studied using direct shear tests on various soil types. Good agreements between the experimental and simulated post-peak stress vs. displacement relationships were obtained. (3) The accuracy of the slope displacement computations can be improved by incorporating the post-peak stress-displacement relationship in the analysis. The results of a case study on the Lu-Shan slope revealed that a transition from pre-peak to postpeak states for the slip surface near the slope toe may have taken place during the first rainfall event in the monitoring period. This critical feature may be overlooked when post-peak strength deterioration is not taken into account in the displacement analysis.
Author Contributions: C.-L.L. contributed to the conception, formal analysis, investigation, data curation, writing-review and editing of the study; C.-C.H. contributed to the methodology, validation, project administration, funding acquisition. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.