Analysis of embankment slope stability: the comparison of finite element limit analysis with limit equilibrium methods

This paper presented the analysis of embankment slope stability by considering the problem of embankment slope stability with special effects that it was filled with sand and was placed on purely cohesive clay. The finite element limit analysis of two-dimensional plane strain was employed to analyze the stability of this problem. The embankment slope height (H), the depth factors (d/H) and the embankment slope angle () for the finite element limit analysis of sand was modeled as a volume element with the properties of Mohr-Coulomb material in drained condition. And the clay was modeled as a volume element with the properties of Tresca material in undrained condition where the parameters were soil unit weight (), undrained shear strength (su) and friction angle (). Parametric studies consisted of three dimensionless variables including depth factors (d/H), friction angle () and embankment slope angle (). Results were summarized in the form of the dimensionless stability number (su/H(FS)) and the design chart and application were presented. In addition, the comparison of the solution of stability number with the limit equilibrium methods and the failure mechanisms were also proposed in this paper.


Introduction
Embankment slope stability is commonly found in various road construction projects and at present, the prediction of the embankment slope stability is one of the important roles for civil engineers. In general, there are several effects on the embankment slope stability, such as embankment slope's physical properties, strength parameters, and embankment slope geometry and the construction process is also affected to the embankment slope as well. Moreover, the embankment slope can be categorized into natural embankment slope or construction processed embankment slope. Also, most of the fill slopes are often found in construction projects, which are characteristics of the embankment slope [1,2].
Limit equilibrium method (LEM) is one of the conventional methods which are used to evaluate embankment slope stability because LEM is a simple and inexplicable method. However, there are some limitations in this method, such as the shear assumption that occurs in the division of the components and the assumption of the failure curved lines. In general, LEM bases on the principle of moment and the equilibrium of force in order to evaluate the stability of the embankment slope. There are many important researches such as in [3,4,5,6,7] which analyzed the problems of embankment slope stability by using LEM, but only in 2D analysis system. In the analysis of the embankment slope stability, it is necessary to assume the failure curved line and then divides the components of the failure curved line, which must base on the principle of moment and the equilibrium of force in order to base on the basic principles of LEM. Thereafter, the assumptions present in these researches should check for ensuring that the assumptions are correct and consistent with the principle of LEM.
Limit analysis can also use to analyze embankment slope stability problems and the use of these methods has been incorporated into the research such as in [8,9,10,11]. On the other hand, some researchers have studied in order to find solutions to the case of layered soil slope and cohesive-frictional slope such as in [12,13]. However, there are not many researches on the solution which is obtained from the analysis of embankment slope which is filled of sand and is placed on purely cohesive clay by using finite element limit analysis. Therefore, the research in this topic is still important and should place for the further research in this topic.
For the above reason, this paper proposed the analysis of embankment slope stability by considering the problems of embankment slope stability with special effects that it was filled of sand and was placed on purely cohesive clay by using finite element limit analysis both upper and lower bound on plane strain (2D Plane Strain). Moreover, the comparison of the solutions of stability numbers to the LEM and its applications were presented in this research.

Method of analysis
This research used finite element limit analysis program (OptumG2) [14] to simulate and analyze upper bound (UB) and lower bound (LB) solutions. This was an analytical method based on the principle of plastic limit theorem, which provided accurate analysis results and could apply to common parameters in geotechnical engineer. Finding the answer with plastic bound theorems by using those principles was that if the answer limited to a certain period of time and the computational system was improved until the whole range of answers was closest to each other, the solution would close to the most exact solution. The plastic bound theorems were divided into two theorems: lower bound theorem and upper bound theorem.
In this study, the unit weight of sand was equal to the unit weight of purely cohesive clay. The Eq. (1) was discovered and was highlighted as a major discovery in the study in [15]. However, the friction angle value () was not included in Eq. (1), which illustrated by numerical solutions in the form of stability graph by considering the difference in . In addition, failure mechanisms, which were found, had the appearance of rotational failure (Base slide) and this kind of failure mechanisms occurred in the soft clay layer that was placed on the stiff clay layer. Therefore, it showed that the clay layer was related to failure mechanism of fill embankment slope.
In this study, the results were presented by using the dimensionless stability number, as shown in Eq. (1) Embankment slope height (H), depth factors (d/H) and embankment slope angle () were used for finite element limit analysis. The embankment slope, which was filled of sand, was modeled as a volume element with the properties of Mohr-Coulomb material in drained condition. And the purely cohesive clay was modeled as a volume element with the properties of Tresca material in undrained condition. The parameters were soil unit weight (), undrained shear strength (su) and friction angle () as shown in Fig. 1.
At the boundary condition of this model, the LB was neither horizontal nor vertical (ux = uy = 0). The left and the right bound could not move in a horizontal plane (ux = 0). The ground surface of this model was a free model, without any axial anchorage, as shown in Fig. 2.  The model was a plane strain which the soil element of the finite element mesh of LB was a triangular element with three points of nodal unknown stresses. And the UB was a triangular element with six points of nodal unknown displacement node, as shown in Fig. 3.

Fig. 3. Triangular elements for lower and upper bound analysis
The researcher considered for using OptumG2 program which this program used finite element limit analysis. Moreover, there was a special function, which the model added more elements or adaptive mesh itself and it made the solution of the problems to be closer to the exact solution. The model was designed for five adaptive iterations per the solution analysis of one problem. This function increased the number of elements in areas with high stress and strain changes, also known as effective area, such as the position that occurred to the slip-line field. Five adaptive steps were selected to obtain an accurate solution, where an initial mesh with the number of 5,000 elements, was automatically adapted an increased to a final mesh with the number of 30,000 elements.
The embankment slope angles () were 15, 22.5, 30 and 40 and were considered to use for analyzing and the friction angles () were 25, 30, 35 and 45 and set the value of  = . It was necessary to analyze of embankment slope stability which was filled of sand and was placed on purely cohesive clay that was a rigidperfectly plastic with associated flow rule and undrained young's modulus (Eu/su = 300). For the case of drained young's modulus (E = 2,800N), the undrained poisson's ratio (u) was 0.495 and the drained poisson's ratio of filled of sand () was 0.35, including the value of depth factors (d/H) which was between 1.5, 2, 3, 4 and 5. When considering the depth factors (d/H) was 1, it was equal to the homogeneous embankment slope which the solution of the stability number in this case corresponded to the solution given by [16]. The solution from the finite

Upper Bound Element
Displacement node u y u x element limit analysis was summed in the form of the dimensionless stability number, su/H(FS), as shown in Eq. (1) For LEM, the solutions of this research used AutoSLOPE program by [17] to simulate and to analyze for comparing the solutions in this research. The program automatically found the Self-Adaptive Failure Surface Generation (SAFEGEN), which was a subset of algorithm in the AutoSLOPE program. In addition, the program also analyzed the embankment slope stability by using the LEM to find the factor of safety (FS). Figure 4 showed the results of the last adaptive mesh of embankment slope which was filled of sand and was placed on purely cohesive clay and on two-dimensional plane strain for the case of  = 15,  = 25, d/H = 3 which compared to the solution of the LB and UB. After the 5th adaptive mesh, the number of elements of the model was increased in the effective area or in the shifted section of high shear stress. During the transition period of high shear stress, it showed as a curve through the entire base slide which the position and appearance of the shear stresses were occurred when the soft soil horizon was placed on hard soil horizon. Figure 5 showed the location of the incremental shear strain dissipation of embankment slope which was filled of sand and was placed on purely cohesive clay for the case of  = 15,  = 25, d/H = 3 which compared to the solution of the LB and UB. And when compared to the same mesh, it could see that the displacement as well as the shear strain was in the same position as Fig. 4. The lines in Fig. 10(a) and 10(b) had different values of . The value of the top line was  = 25 and the value of the bottom line was  = 45 which the bottom line in this case had the maximum value of friction angle that was caused the embankment slope to remain stable. When compared with the solution of the LEM, it was found that the solution had a consistent tendency. However, the solution obtained by the LEM was either less than or equal to the UB or LB and it was not in the range of lower and upper areas, but there were some cases as well, as shown in Fig. 10(b) for the case of  = 22.5,  = 25. However, for the case of embankment slope angle values () were 30 and 40 which had the high embankment slope, there was friction angle value () that could maintain the stability of the slope (FS  1) which the values of  were between 35 and 45 for the case of embankment slope angle values were  = 30 and  = 40 as shown in Fig.  10(c) and 10(d), respectively. Figure 10(c) and 10(d) illustrated the solution from finite element limit analysis by using the OptumG2 program. The solution was the relationship between dimensionless variables of the stability number (su/H(FS)) and depth factors (d/H) from 1.5 to 5 which had the embankment slope angle () of 30 and 40. Thus, the graph in this figure was divided into two types: the solution of the LB was a dashed line, while the solution of UB was a solid line. From the figure, it could see that the solution had non-linear relationship which the LB and the UB were nearly the same, almost the same line, indicating that the solution had a close value with the exact solution. The lines in Fig. 10(c) and 10(d) had different friction angle values (). The value of the top line was  = 35 and the value of the bottom line was  = 45 which the bottom line in this case had the maximum value of friction angle that caused the embankment slope to remain stable. When compared with the solution of the LEM, it was found that the solution had a consistent tendency only in the depth factors from d/H = 3 -5. However, the solution obtained by the limit equilibrium method was not in the range of the UB or LB in the depth factors (d/H) from 1.5 to 3, as shown in Fig. 10(c) and 10(d). On the other hand, when the value of embankment slope angle was increased, the embankment slope stability would also decrease which all the solution in these figures were  = 25 and  = 30 respectively. When compared with the solution of the LEM, it was found that the solution was not consistent in the case of embankment slope angle was  = 22.5 while the value of friction angle was  = 25. The solution was not in the range of the LB and UB, in some cases, as shown in Fig. 11(a). Figure 11(d) showed the relationship between dimensionless variables of the stability number (su/H(FS)) and depth factors (d/H) as shown in Fig. 11(c) but the lines in this figure showed the different value of embankment slope angle (), which was as follows 15, 22.5 and 30. The solution was very close to the same line and the embankment slope stability would decrease from the bottom line to the top line which every solution in this figure was  = 35. When compared with the solution of the LEM, it was found that the solution had a consistent tendency. However, the solution obtained by the LEM was either less than or equal to the UB or LB and it was not in the range of lower and upper areas, but there were some cases as well, as shown in Fig. 11(c) for the case of  = 30,  = 35 and in Fig. 11(d) for the case of  = 40,  = 45. From that relationship, Fig. 10 to Fig. 11, it showed that most of the solutions of stability, su/H(FS) and when considering the embankment slope angle () between 15 to 40 and depth factors (d/H) between 1.5 to 3, it showed a decrease in dimensionless stability number (su/H(FS)). This meant that when the embankment slopes which was filled of sand to increase the height of the original soil layer, it would reduce the embankment slope stability. In addition, when considering the depth factors (d/H) which the value was between 4 to 5, it was found that the dimensionless stability number (su/H(FS)) slightly changed and the embankment slope was still stable (FS>1). Moreover, it showed that when it was filled of sand for the construction of the embankment slope, the depth factors would not affect to the embankment slope stability.

Results and discussion
When compared with the solution of the LEM to the finite element limit analysis, it was found that the solution had a consistent tendency in some cases. However, the solution obtained by the LEM was either less than or equal to the UB or LB and it was not in the range of lower and upper areas, but there were some cases as well.

Application examples
A case study used to demonstrate the convenience of using the stability charts to design the embankment slope which was filled of sand. For the first case, the embankment slope angles was  = 15, the friction angles was  = 45, undrained shear strength of purely cohesive clay under the filled sand was su = 20 kPa, soil unit weight was  = 18 kN/m 3 , embankment slope height was H = 5 m, and depth factors was d/H = 4. The parameters were then taken to find the value of dimensionless stability number (su/H(FS)), as shown in Fig. 10 The second case, The embankment slope angles was  = 22.5, the friction angles was  = 35, undrained shear strength of purely cohesive clay under the filled sand was su = 20 kPa, soil unit weight was  = 18 kN/m 3 , embankment slope height was H = 5 m, and depth factors was d/H = 2. The parameters were then taken to find the value of dimensionless stability number (su/H(FS)), as shown in Fig. 11(c), su/H(FS) = 0.171 (LB), su/H(FS) = 0.169 (UB), which was 0.170 (average). After that, the obtained values were calculated for factor of safety (FS) from Eq. (1) which the value was FS = 1.31.

Conclusions
This paper presented the analysis of embankment slope stability by considering the problem of embankment slope stability with special effects that it was filled with sand and was placed on purely cohesive clay. The finite element limit analysis of two-dimensional plane strain was employed to analyze the stability of this problem. The embankment slope height (H), the depth factors (d/H) and the embankment slope angle () for the finite element limit analysis of sand was modeled as a volume element with the properties of Mohr-Coulomb material in drained condition. And the clay was modeled as a volume element with the properties of Tresca material in undrained condition where the parameters were soil unit weight (), undrained shear strength (su) and friction angle (). Parametric studies consisted of three dimensionless variables including depth factors (d/H), friction angle () and embankment slope angle (). Results were summarized in the form of the dimensionless stability number (su/H(FS)) as a function of those input variables.
The solution presented in the form of failure mechanism which consisted of the adaptive mesh and incremental shear strain dissipation. From the analysis result, it was found that the number of elements increased significantly at high shear strain dissipation and it showed as a curve through the entire base slide which the nature of the position was consistent with the research by Lim et al., (2015) [13]. That meant the curvature of failure occurred when the soft soil layer was placed on the hard soil layer. And when the depth factors (d/H) was 1.5 to 3, it showed that the dimensionless stability number (su/H(FS)) tendency decreased (FS  1). This meant that when the embankments slope which was filled of sand to increase the height of the original soil layer, it would reduce the stability of the embankment slope. In addition, when considering the depth factors (d/H) which the value was between 4 to 5, it was found that the dimensionless stability number (su/H(FS)) was slightly changed and the embankment slope was still stable (FS>1). Therefore, when it was filled with the frictional fill materials, such as sand and gravel for the embankment slope construction with the depth factors was d/H  4, it would not affect to the embankment slope stability.
The solution was shown in the graph of the relationship between the dimensionless stability number (su/H(FS)) and the depth factors (d/H). From the comparison, the dimensionless stability number (su/H(FS)) significantly reduced when the depth factors (d/H) was between 1.5 to 3 and friction angle () increased while embankment slope angle () reduced. Moreover, when the depth factors (d/H) was between 4 to 5, it was found that the dimensionless stability number (su/H(FS)) slightly changed and the embankment slope was more stable compared to the range of the depth factors (d/H) was between 1.5 to 3.
When compared with the solution of the LEM to the finite element limit analysis by analyzing the dimensionless stability (su/H(FS)) from the problem of the embankment slope which was filled of sand and was placed on purely cohesive clay, it was found that most of the solution, when embankment slope angle () was 15 to 40 and the depth factors (d/H) was between 1.5 to 5, the dimensionless stability (su/H(FS)) had a consistent tendency. In contrast, when embankment slope angle () was in a very high area and also the value of depth factors (d/H) was very high, it was found that the dimensionless stability (su/H(FS)) had an obvious differences which it showed that the problem analysis of embankment slope stability which was filled of sand and was placed on purely cohesive clay used the LEM only in some cases.