Abstract
The slope stability problem is an important issue for the safety of human beings and structures. The stability analysis of the three-dimensional (3D) slope is essential to prevent landslides, but the most important and difficult problem is how to determine the 3D critical slip surface with the minimum factor of safety in earth slopes. Basing on the slope stress field with the finite element method, a stability analysis method is proposed to determine the critical slip surface and the corresponding safety factor of 3D soil slopes. Spherical and ellipsoidal slip surfaces are considered through the analysis. The moment equilibrium is used to compute the safety factor combined with the Mohr-Coulomb criteria and the limit equilibrium principle. Some assumptions are introduced to reduce the search range of center points and the radius of spheres or ellipsoids. The proposed method is validated by a classical 3D benchmark soil slope. Simulated results indicate that the safety factor of the benchmark slope is 2.14 using the spherical slip surface and 2.19 using the ellipsoidal slip surface, which is close to the results of previous methods. The simulated results indicate that the proposed method can be used for the stability analysis of a 3D soil slope.
1. Introduction
Landsides are a common worldwide geological disaster that can cause heavy casualties and huge economic losses [1–3]. The slope is not only an important environment of human existence but also an important part of engineering construction. Slope instability has been called one of the three major geological disasters, along with earthquakes and volcanoes [4]. The slope instability problem is closely associated with safety and economic benefits; as a result, slope stability analysis has a very important practical significance and economic value. A great deal of accumulated experience has been obtained on slope stability analysis, which is based on limit equilibrium theory [5–7]. Limit equilibrium method is the most popular method in assessing slope stability [8, 9]; however, this method cannot provide the relationship of internal stress and strain in rock or soil slopes, and they do not always provide unique factors of safety owing to the inherent assumptions of limit equilibrium analyses [10]. Currently, finite element method is widely used in slope stability analysis because the method can calculate the stress and strain based on the nonlinear constitutive relation of rock and soil mass [2, 11, 12]. However, this method cannot consider both the stress and moment equilibrium [13]. Although the strength reduction method can determine the failure zone and the safety factor [2, 14, 15], it is difficult to link the finite element calculation results with the traditional safety factor of the slope.
Currently, many engineers still adopt two-dimensional (2D) analysis methods due to the various limitations of the 3D slope stability analysis [6]. The three-dimensional (3D) slope is often simplified as a 2D slope problem because the above methods are not convenient for use in the direct 3D slope stability analysis [16, 17]. The most commonly used method in practical engineering to analyze the 3D slope stability is to consider each profile as a 2D slope [18]. Finally, the safety factor of the slope can be calculated by the weighted average of each profile; however, the method often introduces some errors. At present, various existing computational methods have different accuracies in determining the critical slip surface and the corresponding safety factor [2, 19]. Cheng et al. introduced a practical way in using NURBS surface and ellipsoidal surface to simulate a three-dimensional sliding surface [20]. Ahangar-Asr et al. proposed a three-dimensional approach in conjunction with genetic algorithm (GA) to investigate the effect of earthquake force inclination on the minimum safety factor of slope and the shape and direction of the corresponding failure surface [21]. Hajiazizi and Tavana determined three-dimensional nonspherical critical slip surface in earth slopes using an optimization method to obtain the nonspherical critical slip surface, which is more consistent with the actual slip surface in nature by using the three-dimensional alternating variable local gradient [15]. According to the above references, the safety factor using 3D analysis methods is higher than that from the corresponding 2D analysis; that is, the 3D analysis provides a more economic slope design method [20].
In this paper, a 3D slope stability analysis method, combined with the finite element method and the limit equilibrium principle, is proposed to determine the critical slip surface and the corresponding safety factor of 3D soil slopes. Some assumptions are introduced to reduce the search range of the center points and the radius of spheres or ellipsoids. And the proposed method is validated by a classical benchmark slope and compared with other conventional 3D slope stability analysis methods. Some useful conclusions are presented in this paper.
2. Method
In this section, according to the rotational failure in 3D soil slopes, the stress field of the 3D soil slope is calculated using the finite element method. The safety factor of a certain slip surface can be determined by using the elastic-perfectly plastic soil model, the Mohr-Coulomb criteria, and the limit equilibrium principle. And the spherical and ellipsoidal slip surfaces are used to determine the critical slip surface.
2.1. 3D Slope Stability
Landslides (or slope failure) often pose a great threat to human beings and/or manmade structures. Figure 1(a) shows a landslide that occurred in the soil slope at the Zhaotong City, Yunnan, Southwest China, which is a typical form of rotational failure in soils. Rational assessment of the slope stability using analytical methods is very important for the hazard prevention and mitigation of landslides. According to the stability analysis of 3D soil slopes, two key issues should be determined: (a) the location of the critical slip surface and (b) the corresponding safety factor of this critical slip surface.
(a)
(b)
As shown in Figure 1(b), for the rotational sliding in soil slopes, the form of the slip surface is complicated and is influenced by the geological conditions; however, the critical slip surface can be assumed to be a spherical surface or ellipsoidal surface, which can explain most conditions of the actual landslide characteristics in soil slopes. Furthermore, this assumption can ensure the feasibility of determining the safety factor of the critical slip surface.
2.2. Stress Field of a Slope
The stress field of a slope under natural conditions can be simulated by different numerical methods. A commonly used method for determining the stress field of a slope is the finite element method [13]. Figure 2(a) shows the vertical stress distribution of a 3D soil slope solved using the finite element method. A potential slip surface is located behind the stress field of the 3D slope; this surface is the critical slip surface of the slope and must be determined. As shown in Figure 2(a), the slope elements are cut by a spherical surface or an ellipsoidal surface, with several different elements interacting with the cutting surface, and the slip surface is composed of a group of elements. Figure 2(b) shows the stress distribution of a slip surface element in one representative element.
(a)
(b)
As shown in Figure 2(b), the principal stresses , , and of the representative element can be determined via numerical simulation, assuming that the area of the slip surface is and the total stress on the slip surface is : where and are the normal stress and shear stress on the slip surface, respectively.
According to force equilibrium analysis, the total stress, , can be divided into three parts along the axis directions of , , and :where , , and are the angles between the normal direction of the slip surface and the axis directions of , , and , respectively. The relationship between these there angles is given by
The normal stress and total stress on the slip surface can be determined as follows:
The shear stress on the slip surface can be computed as follows:
The shear resistance on the slip surface can be determined according to the Mohr-Coulomb criterion: where and are the cohesion and friction angle of the soil, respectively.
If the geometric information and the stress field of one representative element are known, then the shear stress and shear resistance on this slip surface element can be determined.
2.3. Safety Factor of a Certain Slip Surface
For the stability analysis of 3D soil slopes, the determination of the 3D critical slip surface from a series of potential slip surfaces in global optimization is a difficult task, thus making the solution time tremendously long [19]. Currently, spherical shapes are commonly used for the slip surfaces in the 3D slope stability analysis because of the relative lack of difficulty in performing global optimization. However, the nonspherical slip surface is the most general shape that can be considered in the 3D slope stability analysis process [24]. In this paper, both spherical and ellipsoidal (nonspherical) shapes of slip surfaces are considered for the stability analysis of 3D soil slopes. Some assumptions are introduced during the determination of the safety factor for a certain slip surface, and details are shown in Figure 3.
(a)
(b)
As shown in Figure 3, the sliding directions of the slip surface elements are all parallel to the plane. If the determination of the safety factor of a slip surface is regarded as a two-dimensional problem and the total number of slip surface elements is , then the safety factor () of the slip surface can be computed as follows:where and are the total shear resistance and shear stress on the slip surface, respectively, and and are the cohesion and friction angle of the slope soil, respectively.
However, the safety factor under 3D conditions is difficult to calculate using (7) because the directions of the shear stress and the normal stress are not consistent for each slip surface element. The geometric features of a sphere or an ellipsoid have the same characteristics: a fixed center, where the normal stress for each slip element is pointing to this fixed center and perpendicular to the slip surface. Furthermore, it is assumed that shear stresses are directed consistently with the overall sliding direction, which is all parallel to the plane. Next, the safety factor is defined in terms of moment rather than in the terms of shear stress [10], which can be calculated by integrating all of the slip elements; subsequently, the safety factor of a determined slip surface can be calculated as follows: where is the resisting moment, which is related with the shear resistance on the slip surface; is the slip moment, which is related with the shear stress on the slip surface; and is the vertical distance from the center of a spherical or ellipsoid to the surface of each element. If the slip surface is a spherical style, then is equal to the radius of the sphere; however, if the slip surface is an ellipsoidal style, then is related to the size parameters of the ellipsoid.
2.4. Determination of the Critical Slip Surface
During the process of determining the critical slip surface, the number of calculation steps is high because of the numerous potential slip surfaces; as a result, some assumptions and optimization algorithms should be introduced. The spherical slip surface can be described with the following equation:where are the coordinate values of center point of a sphere and is the radius of the sphere.
However, for the ellipsoidal slip surface, the description equation can expressed as follows: where are the coordinate values of the center point for an ellipsoidal shape and , , and are the radii in , , and directions, respectively.
According to (9), there are four uncertain parameters to determine for the spherical slip surfaces. The soil slope is cut by several different spheres and forms several slip surfaces, which have the corresponding safety factor. The slip surface with the minimum safety factor is the critical one. Figure 4 shows a 3D slope cut via a sphere from a different perspective.
(a)
(b)
As shown in Figure 4, for the spherical slip surface, fortunately, the sphere rotated in any direction has no effect on the shape of slip surface. However, for the ellipsoidal slip surface, nine uncertain parameters should be determined: , , , , , , , , and , where , , and are the rotation angles of the ellipsoid surface around the axis directions , , and , respectively. The increasing number of uncertain parameters will lead to an exponential increase in the computational burden. Figure 5 shows a 3D slope cut by an ellipsoid from different perspectives.
(a)
(b)
To reduce the number of calculations, the following assumptions are introduced. First, the semimajor axis of the ellipsoidal is parallel to the slope (as shown in Figure 5(a), where OP is assumed parallel to MN) and there are no rotations around the directions and . As a result, only one rotation angle exists, where and are equal to 0 and is equal to the inclination of slope. With this assumption, the number of uncertain parameters is reduced to six. Furthermore, it is assumed that a relationship exists among the radius parameters of the ellipsoid; here, the assumption for the relationship among , , and presented by Zhao et al. is used that parameter is the maximum value and the value of is less than [25].
Furthermore, determining the appropriate center coordinate and search range of the sphere or ellipsoid embedded in the three-dimensional slope is a problem [15]. To reduce the number of calculations, the center points are limited in one plane [13, 26], which is the same as the two-dimensional problem. According to the symmetry soil slopes, the center points are located in the vertical plane across midpoint of the slope (details shown in Figure 6).
As shown in Figure 6, the center points are located in the rectangle , and the coordinate is a constant. In this manner, the appropriate search domain can effectively reduce the number of calculations. For a determined geometry of 3D soil slopes, the search range of the center points can be further reduced. The optimization for the search range of the center points can simplify the calculation and improve the computing speed [15]. Figure 7 shows the further optimization assumption for the search range of the center points.
As shown in Figure 7, point is the midpoint of the straight line , and two auxiliary lines are added, which are perpendicular to , respectively. Regarding to the slope geometry and the potential slip surface, the search range of the center points can be greatly reduced into the range of . For the soil slopes, when the sliding surface is tangent to the plane MEF (as shown in Figure 6), the slip surface cuts a small part of the 3D slope, which cannot be the most dangerous slip surface. While the depth of the tangent plane in the slope is greater than MR with the slip surface, the resistance moment of the slip surface will be relatively large [8]. The radius of the spheres or ellipsoids can be limited to a small range [4]. Furthermore, to address the very large number of calculations for the 3D slope, some optimization algorithms can be introduced to improve the computational efficiency [15].
2.5. Calculation Process
The calculation process for determining the critical slip surface of 3D soil slopes from stress field using the finite element method is shown as follows:(1)According to a special 3D soil slope, first, the numerical model of this slope should be established, and, then, the finite element method is used to calculate the stress field of this slope by means of gravity loading; here, the Mohr-Coulomb constitutive model for the soils is used during the numerical simulation process.(2)After the completion of numerical simulation, the computing nodes, elements, mechanical parameters, and stress field of the slope are exported to four data files (nodes.txt, elements.txt, parameters.txt, and stress.txt); the four files can be loaded by the VBA (Visual Basic for Applications) programming language.(3)After the main program has loaded the four exported files, that is, nodes.txt, elements.txt, parameters.txt, and stress.txt, to set the search range of the center points, the rotation angle of the ellipsoid, the range of values of the radius, and other initial conditions, including the boundary condition of the slope, can be automatically detected and the shear strength parameters of each element can be found in the file of parameters.txt.(4)After several random slip surfaces are generated under some certain optimization principles, the safety factor of each slip surface can be determined using (8), where the minimum value is the safety factor of the soil slope and the corresponding slip surface is the critical one.(5)After exporting the geometry information of the critical slip surface, the stress field on the critical slip surface can be visualized using various drawing software packages, such as Surfer or Tecplot.
3. Results
In this section, a classical benchmark soil slope is used verify the proposed stability analysis method for 3D soil slopes [22, 23, 26]. Comparison analysis between the other methods is performed for the location of the critical slip surface and the corresponding safety factor.
3.1. Benchmark Slope Conditions
The 3D benchmark soil slope consists of three different types of soil layers: layers 1, 2, and 3, from top to bottom, where layer 2 is a weak interlayer with a height of 0.3 m [17]. Figure 8 shows the geometric condition and finite element mesh of this 3D benchmark slope, which is commonly used for the verification of the 3D slope stability analysis method.
As shown in Figure 8, the benchmark slope has the following characteristics: a height of 22.2 m, length of 80 m, width of 50 m, and slope inclination of 27°. Here, the commercial finite element software ABAQUA is used to calculate the stress field within the soil slope. The slope mesh consists of solid-brick elements with eight nodes and elastic-perfectly plastic behavior (Figure 8). Table 1 presents the mechanical parameters of the benchmark soil slope.
As shown in Figure 8, the bottom surface is fixed of displacement constraints in , , and directions; the two side surfaces in direction are fixed of displacement constraints in and directions; and the two side surfaces in direction are also fixed of displacement constraints in and directions. After all the initial conditions are defined for the numerical mesh, the elastic stress of the slope is firstly calculated under the only action of gravity loading. When the elastic process of the gravity loading is finished, the associated displacements and strains of the numerical nodes and elements are all reset to zero. The initial stress field for the soil slope is determined, and then all the numerical meshes are used the elastic-perfectly plastic properties. The secondary calculation process of the slope stress field is carried out by using the Mohr-Coulomb model; this stress field is the input data for slope stability analysis. After the completion of numerical simulation, the computing nodes, elements, mechanical parameters, and stress field of the slope are exported to four data files by using the self-designed secondary development codes for the ABAQUS. Figure 9 shows the stress distribution characteristics of the benchmark slope.
(a)
(b)
From Figure 9, the benchmark slope is mainly under compressive stress with a low level; because the slope height is at a small scale, the maximum compressive stress of the benchmark slope is in the range of 0 to 0.45 MPa. Furthermore, small tensile stress exists at the slope surface, and stress dislocation appears around weak layer 2. As the potential slip surface cannot be directly determined from the stress field of the 3D slope, some indirect methods should be introduced to determine the critical slip surface and its corresponding safety factor. The following stability analysis results are performed based on these stress field results for the benchmark slope.
3.2. Spherical Slip Surface
First, a potential spherical slip surface is assumed to exist in the 3D benchmark soil slope, and the center of these potential spherical slip surfaces are located in the middle of the width direction, as shown in Figure 10(a). The benchmark slope is cut by these potential spherical slip surfaces.
(a)
(b)
According to the geometric information of the benchmark slope, the precise search range for the center point of potential spheres can be determined, as shown in Figure 10(b). The details of the search range of the center points can be gained from Figure 10(b); that is, the center point must be located in the range . As a result, the corresponding radius range of the spheres is from 18.3 m to 29.3 m.
Once the search range of the center points and radius for the spherical slip surfaces are determined, several random spheres are generated. The slope is cut by these spheres; for each spherical slip surface, the geometric information and stress distribution of the surface can be computed and the corresponding safety factor can be determined, where the minimum safety factor is the critical slip surface. Figure 11(a) shows the normal stress distribution on the critical slip surface, which is in the range of 6.0 kPa to 145 kPa. The normal stress on the slip surface is related with the stress field of the soil slope and can be used to compute the shear resistance of the slip surface. Figure 11(b) shows the safety factor results of the critical slip surface.
(a)
(b)
From Figure 11(b), the center point of the critical slip surface is “6.1 m, 19.3 m, and 25.0 m” and the radius is 24.5 m. The slip moment and resisting moment are 187.9 kN·m and 401.2 kN·m, respectively, the corresponding safety factor is 2.135. A comparison of the proposed method with other previous methods is presented in Table 2 (previous methods are from [8, 17, 22]).
As shown in Table 2, although there is little difference in the shape of spherical slip surface, the safety factor of benchmark soil slope using the presented method is close to those of other methods. The simulated results indicate that the proposed method can be used for the spherical stability analysis of a 3D soil slope.
3.3. Ellipsoidal Slip Surface
The search of the ellipsoidal slip surface is more complicated than the search of the spherical slip surface. Some assumptions are introduced to reduce the number of calculations; Figure 12(a) shows the 3D soil slope cut by a ellipsoid, where the center of these potential ellipsoidal slip surfaces are located at the middle in the width direction.
(a)
(b)
As shown in Figure 12(a), no rotation angle exists around the axis directions and , but the rotation angle around the axis direction is equal to the inclination of the soil slope and the rotation angles are “27°, 0°, and 0°,” respectively. The initial value of are set as (20 m, 15 m, and 15 m), respectively; these values are no larger than 25 m. The search range of the center points is the same as the spherical slip surface.
Once the search range of the center points and radius for the ellipsoidal slip surfaces are determined, several random spheres are generated. The slope is cut by these ellipsoids, and the safety factor of these ellipsoidal slip surfaces can be computed, where the minimum one is the critical slip surface. Figure 12(b) shows the safety factor results of the critical slip surface. As shown in Figure 12(b), the center point of the critical slip surface is “8.0 m, 15.1 m, and 25.0 m” and the radius values are “24.4 m, 11.7 m, and 20.7 m,” respectively. The slip moment and resisting moment are 189.1 kN·m and 414.6 kN·m, respectively; the corresponding safety factor is 2.193. Comparison of the presented method with previous methods is presented in Table 3 (improved simplified Janbu’s method [23]).
As presented in Table 3, the safety factor determined using the proposed method is 2.193, which is larger than that of Zhang’s result and is also larger than that of the spherical slip surface because the proposed method meets all forces and moment equilibrium conditions base on slope stress field. However, from all of other results and this paper’s results, we found that the safety factor of this benchmark slope is in the range of 2.0 to 2.2. The simulated results indicate that the presented method can suit not only the stability analysis of spherical slip surface but also the ellipsoidal slip surface.
4. Discussion and Conclusion
In geotechnical engineering, the slope stability problem is an important issue for construction safety and long-term operation. Because many theoretical issues remain not fully resolved, the mechanics of the slope instability should continue to be studied [27–29]. The finite element method is the most widely used method for slope stability analysis [4]. However, the critical slip surface and the corresponding safety factor cannot be determined directly based on the stress field and displacement field of the slope; thus, an indirect method should be introduced to solve this problem, such as the strength reduction method [13, 28]. However, this method requires iterative calculations and limits the node unbalanced force, which inevitably result in the calculating efficiency not high. In this paper, a stability analysis method was proposed to determine the critical slip surface of 3D soils slopes from the stress field with the finite element method. Spherical and ellipsoidal slip surfaces are considered for the 3D soil slopes. And the safety factor of the slip surface is computed based on Mohr-Coulomb criteria and limit equilibrium principle.
A classical benchmark soil slope is used to verify the proposed stability analysis method for 3D soil slopes. The simulated results indicate that, for the spherical slip surface, the slip moment and resisting moment are 187.9 kN·m and 401.2 kN·m, respectively, and the corresponding safety factor is 2.135, which is close to the value of the previous three different methods considered [8, 17, 22]. For the ellipsoidal slip surface, the slip moment and resisting moment are 189.1 kN·m and 414.6 kN·m, respectively; the corresponding safety factor is 2.193, which is close to Zhang’s result [23]. From all of the previous results and this paper’s results, the safety factor of this benchmark slope is in the range of 2.0 to 2.2.
The proposed method can resolve some limitations of the existing 3D methods that are based on the limit equilibrium approach. For the limit equilibrium approach, the method cannot consider the internal stress field of the slope, and many assumptions are necessary considered for the slices of slope soils. Besides, the traditional finite element method cannot directly determine the critical slip surface of slope and calculate the corresponding safety factor. The organic combination of the FEM with the limit equilibrium principle can fully account for soil deformation and the influence of elastoplastic stress adjustment on the slope stability. That means the proposed method can determine the critical slip surface via mathematical programming using the stress field of the slope as the input. Based on the FEM and limit equilibrium principle, the slope transfers the element interaction force by itself. When a critical state is reached, the element nodes appear failure. The proposed method can not only determine the location of critical slip surface but also can calculate the corresponding safety factor of the 3D soil slopes. Furthermore, it makes good use of geometric combination to reduce the search range during the calculation process.
Regarding the proposed method for the stability analysis of a 3D soil slope, some reasonable assumptions were introduced to reduce the search range of the center points and radius. In addition, the definition of the safety factor was relatively unique by using the resisting moment and the slip moment. Combined with the assumption of the sliding direction, the safety factor is very convenient for the integral solution in the sliding surface. Therefore the proposed method can provide some realistic guidance for the 3D symmetric soil slope stability.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
The authors gratefully acknowledge the support of the National Natural Science Foundation of China (51209156), Youth Science and Technology Fund of Sichuan Province (2016JQ0011), and the Science Foundation for Excellent Youth Scholars of Sichuan University (2013SCU04A07).