Critical Filling Height of Embankment over Soft Soil: A Three-Dimensional Upper-Bound Limit Analysis

Abstract

This paper investigates the critical filling height of embankments over soft soil using three-dimensional (3D) upper-bound limit analysis based on a rotational log-spiral failure mechanism. Soft soils are characterized by low shear strength and high compressibility, making the accurate determination of critical filling height essential for evaluating embankment stability. Unlike conventional two-dimensional (2D) analyses, the proposed 3D method captures the true failure mechanism of embankments, providing more realistic and reliable results. The upper-bound analysis equations are derived using the principle of virtual work and solved efficiently through the genetic algorithm (GA), which avoids the limitations of traditional loop and random searching algorithms. The proposed solution is validated by comparing it with existing studies on slope stability and demonstrates higher accuracy and computational efficiency. Parametric studies are conducted to evaluate the influence of the depth–height ratio (the ratio of soft soil depth to embankment height) on the failure width of the embankment, the critical failure surface, and the critical filling height. Results show that the critical failure surface is tangential to the bottom of the soft soil layer and the critical filling height increases as the depth–height ratio decreases. The findings provide a set of critical filling heights calculated under various soft soil depths, strength parameters, and embankment geometries, offering practical guidance for embankment design.

1. Introduction

Soft soils possess the characteristics of low shear strength and strong compressibility [1,2] and are always encountered in engineering constructions [3]. When embankments are built on soft soils and the filling height exceeds the critical filling height, i.e., the highest filling height that a natural foundation can bear, lateral displacement and collapse can be caused [4]. Therefore, the critical filling height is a valuable indicator to evaluate the stability of the embankment [5]. The limit equilibrium method and finite element method are popularly used to explore the critical filling height of embankments [6,7,8,9,10,11]. Although the limit equilibrium method has clear mechanical concepts and simple calculation methods, it only considers stress equilibrium and neglects the constitutive relationship of the soil. Therefore, the solution solved by the limit equilibrium method is neither a strict upper limit value nor a strict lower-bound value. The finite element limit analysis method can simulate the process of embankment failure and obtain the distribution of stress–strain fields during embankment instability [12]. However, the modeling process of the finite element analysis is not only very complicated but also, the constitutive model parameters utilized are difficult to determine. More importantly, the precision of the selected parameters greatly determines the accuracy of the calculation results [11,13].

Recently, the upper-bound limit analysis was used to determine the critical state of embankments. For example, Chen et al. [14] calculated the critical filling height of a pile-supported embankment by using two-dimensional (2D) limit analysis; Gong et al. [15] evaluated the seismic stability of embankments considering the anisotropy and nonhomogeneity of soils; and Li et al. [15] assumed a rigid block sliding along two logarithmic spiral surfaces and two-dimensional (2D) limit analysis to explore the critical filling height of the embankment. However, the embankment failures are 3D in nature. Recent studies have increasingly focused on incorporating three-dimensional features to improve the accuracy of stability predictions. These efforts include the exploration of complex soil behaviors, irregular geometries, and dynamic loading conditions. While promising, many of these approaches remain computationally intensive and are challenging to apply in large-scale engineering projects. Building on this foundation, this study integrates a 3D upper-bound limit analysis with a genetic algorithm (GA) to provide a practical and efficient solution for embankments over soft soils.

However, due to the difficulty of constructing admissible base failure mechanisms of embankments, the 3D analysis method for embankments is not as common as other methods, and a limited number of attempts have been conducted so far. Giger and Krizek [16] firstly utilized the 3D limit analysis method to explore the stability of a soil slope. Qu [17] conducted 3D finite element analyses on three full-scale embankments to assess their failure mechanisms over soft clay deposits. Zhuang and Wan analyzed four cases of embankments rapidly loaded to failure by using the 3D total stress circular arc stability analysis. Berrabah et al. [18] employed the 3D numerical method to analyze the stability of the embankment constructed over locally weak zones. Yang et al. [11] adopted a limit equilibrium-based 3D rotational failure mechanism to explore the stability of convex embankments. Another worthy point to note is that traditional loop algorithm and the random searching algorithm are the more commonly applied algorithms to search the most critical slip surface in current research [10,19]. However, the loop algorithm exhibits not only limited computational efficiency but also the potential to overlook the global minimum value if the incremental step is insufficiently small [19,20]. Moreover, although the random search technique proposed by Chen [14] could provide high computational efficiency, it would get trapped in local minimums. These limitations make it difficult to apply such methods to large-scale or highly complex embankment systems, necessitating the development of more efficient and versatile approaches. Compared to these traditional algorithms, the genetic algorithm (GA) utilized in this study significantly enhances computational efficiency and avoids local minima, leading to more reliable results. Therefore, the genetic algorithm (GA) is adopted to search for the critical safety factor of embankments over soft soil. Due to parameter limitations and the adjustment of the search direction during the search process, the genetic algorithm can effectively avoid missing the global minimum [21,22].

This study explores the critical filling height of embankments over soft soil using 3D upper-bound limit analysis. An admissible rotational log-spiral failure mechanism proposed by Michalowski and Drescher [23] is utilized to replicate the embankment failure. The log-spiral failure surface is composed of two sub-failure surfaces, one passing through the embankment and the other passing through the soft soil. Surcharge loads that may be caused during embankment construction are also considered. The critical filling height is determined by balancing the energy rates of external loads, including the equilibrium between soil weight and surcharge loads, and the frictional force acting along the failure surface. Compared to existing 2D methods, which often oversimplify the complex three-dimensional nature of embankment failure, the proposed approach provides a more accurate representation of the failure mechanisms by directly incorporating 3D effects. Additionally, traditional algorithms such as loop searching and random searching methods are limited in computational efficiency and often fail to locate the global minimum, highlighting the need for a more robust optimization approach. The GA with the advantage of computation efficiency and accuracy is used to solve the formulation and find the most critical filling height and critical failure surface. The solution is validated by degenerating the proposed embankment failure mechanism into the general slope mechanism and then comparing it with existing solutions. Parametric studies are conducted to investigate the influence of the ratio of the soft soil depth to the embankment height on the failure width of the embankment, the critical failure surfaces, and the critical filling heights.

2. Base Failure of Embankment over Soft Soil

2.1. Log-Spiral Failure Mechanism

The upper-bound limit analysis theorem relies on the conservation of energy principle and is designed to assess the maximum limit load or failure boundary of materials exhibiting elastic perfectly plastic behavior [24]. The energy rate equilibrium dictates that the strain rate of the failure mechanism must adhere to both the boundary condition and the corresponding flow rule of the yield criterion simultaneously (Appendix A). Hence, the 3D rotational mechanism proposed by Michalowski and Drescher [23] is the most reasonable failure model for a given frictional-cohesive soil embankment. In this paper, the 3D failure mechanism presented by Michalowski and Drescher [23] for a single-layer slope is further developed and extended for the base failure mechanism of embankments over soft soil. Moreover, according to the characteristics that the embankment is compacted while the foundation is not compacted, the failure mechanism is established by combining the layering characteristics of the embankment and the foundation. The failure mechanism of an embankment over soft soil is shown in Figure 1.

Figure 1 shows the log-spiral failure mechanism of the embankment over soft soil. Based on the layering characteristics of the embankment and the foundation, the failure mechanism is composed of two curvilinear cones with the same symmetry plane, and the surfaces must be tangential to a cone with apex angle $2\phi$ shown as $P_1,P_2$. The schematic illustrations of the embankment failure and soft soil failure are shown as red and blue curves, respectively. The lower outline of two curvilinear cones forming the failure surface can be formulated as two log spirals (i.e., the lines $\mathrm{AC}$ and $\mathrm{AE}$).

$$\begin{array}{c}{r}_1=r_0{e}^{\left(\theta -{\theta }_{\mathrm{A}}\right)tan{\phi }_1}\\ \left({\theta }_{\mathrm{A}}\le \theta \le {\theta }_{\mathrm{C}}\right)\end{array}$$

(1)

$$\begin{array}{c}{r}_2=r_{\mathrm{C}}{e}^{\left(\theta -{\theta }_{\mathrm{C}}\right)tan{\phi }_2}\\ ({\theta }_{\mathrm{C}}<\theta \le {\theta }_{\mathrm{E}})\end{array}$$

(2)

where $r_0$, i.e., $\mathrm{OA}$, is the initial radius of log spirals; ${\theta }_{\mathrm{A}}$ and ${\theta }_{\mathrm{E}}$ are the initial angle and the final angle describing the rotational mechanism, respectively; ${\theta }_{\mathrm{C}}$ is the angle that the failure surface just passes through the interface between the embankment and soft soil foundation; $r_{\mathrm{C}}$ is the radius of the log spiral when the failure surface just crosses the interface between the embankment and the soft soil foundation; and ${\phi }_1$and ${\phi }_2$, separately, represent the internal frictional angles of the embankment and the soft soil. The auxiliary angle, ${\beta }^{\prime }$, is introduced to ensure that the embankment failure penetrates the soft foundation, and this angle is smaller than the slope angle, $\beta$. The shaded area of the circle denotes the location where the failure mechanism just hit the soft soil surface. In the cross-sectional circles, $d_i\left(i=1,2,3\right)$ represents the separation between the center of the cross-sectional circle and the surface of the embankment; $d_4$ is the distance between the center of the cross-sectional circle and the surface of soft soil; $R_i\left(i=1,2\right)$ are the radii of the rotating circle; and $c_1$ and $c_2$ are the cohesions of the embankment and the soft ground, respectively.

2.2. Energy Rate Balance

In upper-bound limit analysis, the energy rate balance is employed to establish formulations, ensuring that the work performed by external loads equals the energy dissipated along the failure surface [24]. Note that the rotational block is assumed to be rigid and without deformation and volume change [23]. Therefore, the internal energy is only dissipated along the sliding failure surfaces. The effect of surcharge loads is taken into account in this study. During embankment construction, surcharge loads are primarily caused by the weight of large rollers or construction debris, which can significantly affect the stability of the embankment by increasing the external work rates [25]. Hence, the forces contributing to the external work rates consist of the soil weight and surcharge loads, which are illustrated in Figure 2.

Then, the energy rate balance equation can be written as

$$W_{\mathsf{\gamma }}^{3\mathrm{D}}+W_{\mathrm{q}}^{3\mathrm{D}}=D_{\mathsf{\gamma }}^{3\mathrm{D}}$$

(3)

where $D_{\mathsf{\gamma }}^{3\mathrm{D}}$ is the rate of internal energy dissipation of the failure mechanism; $W_{\mathsf{\gamma }}^{3\mathrm{D}}$ is the work rate of the soil weight; and $W_{\mathrm{q}}^{3\mathrm{D}}$ is the work rate due to the surcharge loads. The existing literature Zheng et al. [26] has given the details of deriving $D_{\mathsf{\gamma }}^{3\mathrm{D}}$ and $W_{\mathsf{\gamma }}^{3\mathrm{D}}$.

2.3. Work Rate of Surcharge Loads

The uniform load, $q$, is introduced to represent the surcharge loads on the embankment. The magnitude of uniform load is computed based on the entire failure area involved. The rotating mechanism is created by progressively rotating a circle with an expanding diameter around the axis. Consequently, numerous potential failure regions emerge from the intersection of the embankment’s crest with the rotating mechanism during the search for the failure surface [23,27,28]. The detailed expressions of failure area, $S$, can be written as

$$S={\displaystyle \frac{1}{2}}\pi -2\vartheta {[{\displaystyle \frac{{(r_0\cos {\theta }_{\mathrm{A}}-r_0\sin {\theta }_{\mathrm{A}}\cot {\theta }_{\mathrm{B}})}^2+d_{\mathrm{B}}^2}{2\left(r_0\cos {\theta }_{\mathrm{A}}-r_0\sin {\theta }_{\mathrm{A}}\cot {\theta }_{\mathrm{B}}\right)}}]}^2-{\displaystyle \frac{d_{\mathrm{B}}^2-{(r_0\cos {\theta }_{\mathrm{A}}-r_0\sin {\theta }_{\mathrm{A}}\cot {\theta }_{\mathrm{B}})}^2}{2\left(r_0\cos {\theta }_{\mathrm{A}}-r_0\sin {\theta }_{\mathrm{A}}\cot {\theta }_{\mathrm{B}}\right)}}\sqrt{R_1^2-d_1^2}$$

(4)

where $d_{\mathrm{B}}$ is the length of $\mathrm{BS}$ in Figure 2 and $\vartheta$ is the angle corresponding to the trajectory $\mathrm{ASS}\prime$. The expressions of the two parameters can be calculated as

$$d_{\mathrm{B}}=\sqrt{R_1^2-{\left({\displaystyle \frac{r_{0}sin{\theta }_{\mathrm{A}}}{sin{\theta }_{\mathrm{B}}}}-r_{m1}\right)}^2}$$

(5)

$$\vartheta =\arctan \left[{\displaystyle \frac{(R_1^2-d_{\mathrm{B}}^2)-r_0^2{\left(\cos {\theta }_{\mathrm{A}}-\sin {\theta }_{\mathrm{A}}\cot {\theta }_{\mathrm{B}}\right)}^2}{2r_0\sqrt{R_1^2-d_{\mathrm{B}}^2}\left(\cos {\theta }_{\mathrm{A}}-\sin {\theta }_{\mathrm{A}}\cot {\theta }_{\mathrm{B}}\right)}}\right]$$

(6)

The surcharge loads are considered as a body force. The rate of work of the surcharge load is calculated by multiplying the sum of the forces on the failure surface $AS{S}^{\prime }$ by the distance from the centroid of the failure surface to the center of rotation, $O$. Moreover, the work rate of the surcharge load, $W_{\mathrm{q}}^{3\mathrm{D}}$, is dependent on the angular velocity, $\omega$, and the equation can be written as

$$W_{\mathrm{q}}^{3\mathrm{D}}=\omega Sq\left\{\begin{array}{c}{r}_{0}sin{\theta }_{\mathrm{A}}cot{\theta }_{\mathrm{B}}+\\ {\displaystyle \frac{4[\left(r_0\cos {\theta }_{\mathrm{A}}-r_0\sin {\theta }_{\mathrm{A}}\cot {\theta }_{\mathrm{B}}{)}^2+d_{\mathrm{B}}^2\right]\cos \vartheta }{6\left(r_0\cos {\theta }_{\mathrm{A}}-r_0\sin {\theta }_{\mathrm{A}}\cot {\theta }_{\mathrm{B}}\right)\left(\pi -2\vartheta \right)}}\end{array}\right\}$$

(7)

Furthermore, in addition to the special distribution form of surcharge load described above, the rectangular distribution form of surcharge load is also taken into account. Surcharge load is evenly distributed on the embankment crest, and the diagrammatic drawing is illustrated in Figure 3.

In the case of rectangular distribution, the work rate of the surcharge load can be defined as the scalar product of the uniform load and the velocity vector. The rate of external work due to the surcharge load can be written as follows:

$$W_{\mathrm{q}}^{3\mathrm{D}}=2\omega r_0^2{\sin }^2{\theta }_{\mathrm{A}}{{\displaystyle \int }}_{{\theta }_{\mathrm{A}}}^{{\theta }_{\mathrm{B}}}q{\displaystyle \frac{\cos \theta }{{\sin }^3\theta }}\sqrt{R_1^2-d_1^2}\mathrm{d}\theta$$

(8)

2.4. Critical Filling Height of Embankment

The critical filling height of the embankment is the highest filling height that the natural foundation can bear and is an important metric to evaluate the stability of the embankment [25]. By utilizing the upper-bound theorem of limit analysis and balancing the work rate of external forces with the energy dissipation rate along the failure surface, the stability factor is determined as follows:

$${\displaystyle \frac{\gamma H}{c}}=f\left({\theta }_{\mathrm{A}},{\theta }_{\mathrm{E}},{\beta }^{\prime },D_{\mathrm{cr}},{\displaystyle \frac{r_{01}^{\prime }}{r_{01}}},{\displaystyle \frac{r_{02}^{\prime }}{r_{02}}}\right)$$

(9)

where $\gamma$ is the unit weight of the soil. The embankment height, $H$, is not defined, so $H$ needs to be determined from the geometrical and trigonometric relations as

$$H=r_0\left[e^{\left({\theta }_{\mathrm{A}}-{\theta }_{\mathrm{C}}\right)\tan {\phi }_1+\left({\theta }_{\mathrm{E}}-{\theta }_{\mathrm{C}}\right)\tan {\phi }_2}\sin {\theta }_{\mathrm{E}}-\sin {\theta }_{\mathrm{A}}\right]$$

(10)

The function $f\left({\theta }_{\mathrm{A}},{\theta }_{\mathrm{E}},{\beta }^{\prime },D_{\mathrm{cr}},r_{01}^{\prime }/r_{01},r_{02}^{\prime }/r_{02}\right)$ can reach a minimum value, i.e., the least upper-bound value, when ${\theta }_{\mathrm{A}},{\theta }_{\mathrm{E}},{\beta }^{\prime },D_{\mathrm{cr}},r_{01}^{\prime }/r_{01}$ and $r_{02}^{\prime }/r_{02}$ satisfy the following conditions:

$${\displaystyle \frac{\partial F_{\mathrm{s}}}{\partial {\theta }_{\mathrm{A}}}}=0$$

(11)

$${\displaystyle \frac{\partial F_{\mathrm{s}}}{\partial {\theta }_{\mathrm{E}}}}=0$$

(12)

$${\displaystyle \frac{\partial F_{\mathrm{s}}}{\partial {\beta }^{\prime }}}=0$$

(13)

$${\displaystyle \frac{\partial F_{\mathrm{s}}}{\partial D_{\mathrm{cr}}}}=0$$

(14)

$${\displaystyle \frac{\partial F_{\mathrm{s}}}{\frac{\partial r_{01}^{\prime }}{r_{01}}}}=0$$

(15)

$${\displaystyle \frac{\partial F_{\mathrm{s}}}{\frac{\partial r_{02}^{\prime }}{r_{02}}}}=0$$

(16)

where $\partial F_{\mathrm{s}}$ is the stability factor increment; $\partial {\theta }_{\mathrm{A}}$ and $\partial {\theta }_{\mathrm{E}}$ are the rotation angle increments; $\partial D_{\mathrm{cr}}$ is the increment of the ratio of the soft soil depth to embankment height; and $\partial {\beta }^{\prime }$ is the auxiliary angle increment. When the derivatives of the reduction factor $F_s$ with respect to ${\theta }_{\mathrm{A}},{\theta }_{\mathrm{E}},{\beta }^{\prime },D_{\mathrm{cr}},r_{01}^{\prime }/r_{01}$ and $r_{02}^{\prime }/r_{02}$ are equal to zero, the function $f\left({\theta }_{\mathrm{A}},{\theta }_{\mathrm{E}},{\beta }^{\prime },D_{\mathrm{cr}},r_{01}^{\prime }/r_{01},r_{02}^{\prime }/r_{02}\right)$ reaches the minimum value. Consequently, the critical filling height can be expressed as

$$H_{\mathrm{cr}}\le {\displaystyle \frac{c}{\gamma }}\min \left[f\left({\theta }_{\mathrm{A}},{\theta }_{\mathrm{E}},{\beta }^{\prime },D_{\mathrm{cr}},{\displaystyle \frac{r_{01}^{\prime }}{r_{01}}},{\displaystyle \frac{r_{02}^{\prime }}{r_{02}}}\right)\right]$$

(17)

3. Optimization Procedure

The genetic algorithm (GA), as a machine learning method, is utilized to search for the critical value of $\gamma H/c$ and identify all possible critical failure surfaces based on the energy rate balance equation derived for the failure mechanisms. Compared to conventional traversal search algorithms, the GA significantly improves computational efficiency by leveraging evolutionary strategies such as selection, crossover, and mutation [29]. These strategies allow the GA to effectively explore a larger solution space and converge more quickly to the optimal solution, making it more advantageous than loop searching and random searching methods. The calculation flowchart of the GA for the problem considered is illustrated in Figure 4.

Figure 4 shows the flowchart employed to search the least upper bound to $\gamma H/c$. This scheme presents the search process of the genetic algorithm. For the given constraints of values of the embankment over soft soil $c_1$, $c_2$, ${\phi }_1$, ${\phi }_2$, ${\gamma }_1$, ${\gamma }_2$, and $\beta$, new individuals in a population are generated. Each generated individual corresponds to a value with independent variables (${\theta }_{\mathrm{A}},{\theta }_{\mathrm{E}},{\beta }^{\prime },D_{\mathrm{cr}},r_{01}^{\prime }/r_{01},r_{02}^{\prime }/r_{02}$). The independent variables ${\theta }_{\mathrm{A}}, {\theta }_{\mathrm{E}}, {\beta }^{\prime },D_{\mathrm{cr}},r_{01}^{\prime }/r_{01}$ and $r_{02}^{\prime }/r_{02}$ are introduced as the input variables because the geometric shape of the failure embankment can be uniquely defined by these variables. In order to explore the critical filling height of the embankment over different depths of soft soil, the depth restriction, $D_{\mathrm{cr}}$, is introduced. $D_{\mathrm{cr}}$ is defined as the ratio of soft soil depth to embankment height, which constrains the depth of failure and requires that the failure surface is tangential to the bottom of the soft soil. The depth of failure mechanism satisfies the geometric restrictions, which are presented as follows:

$$m_1=\left\{r_0{e}^{[\left({\theta }_{\mathrm{C}}-{\theta }_{\mathrm{A}}\right)\tan {\phi }_1+\left(\theta -{\theta }_{\mathrm{C}}\right)\tan {\phi }_2]}-{\displaystyle \frac{r_{\mathrm{E}}\sin {\theta }_{\mathrm{E}}}{\sin \theta }}\right\}\sin \theta$$

(18)

$$m_2=\left(R_2-{\displaystyle \frac{\sin {\theta }_{\mathrm{E}}}{\sin \theta }}{r}_{\mathrm{E}}+r_{\mathrm{m}2}\right)\sin \theta$$

(19)

$$d=\max \left(m_1,m_2\right)$$

(20)

$$D_{\mathrm{cr}}={\displaystyle \frac{d}{H}}$$

(21)

A set of sliding surfaces can be obtained when depth is limited. By utilizing genetic algorithms, the most dangerous slip surface can be searched. The objective function in the search procedure is developed as follows:

$$f\left(x\right)={\displaystyle \frac{\gamma H}{c}}={\displaystyle \frac{W_{3\mathrm{D}}\left({\theta }_{\mathrm{A}},{\theta }_{\mathrm{E}},{\beta }^{\prime },D_{\mathrm{cr}},\frac{r_{01}^{\prime }}{r_{01}},\frac{r_{02}^{\prime }}{r_{02}}\right)}{D_{3\mathrm{D}}\left({\theta }_{\mathrm{A}},{\theta }_{\mathrm{E}},{\beta }^{\prime },D_{\mathrm{cr}},\frac{r_{01}^{\prime }}{r_{01}},\frac{r_{02}^{\prime }}{r_{02}}\right)}}$$

(22)

The best individual is chosen through a set of selection, crossover, mutation, and elitism. Finally, the optimal results with the highest fitness are generated through constant updates.

4. Validation

To confirm the log-spiral failure mechanism of the embankment over soft soil, the embankment failure mechanism is simplified to the general slope failure, and the outcomes are juxtaposed with those from Michalowski and Drescher [23] and Gao et al. [30]. It is worth noting that Michalowski and Drescher [23] utilized the loop searching method, which systematically evaluates potential solutions within a predefined range, while Gao et al. [30] adopted the random searching technique, which explores potential solutions by generating random values within a specified domain. Both approaches aim to identify the critical height of $\gamma H/c$, but they differ in terms of their searching strategies and computational efficiency. The comparisons in Figure 5 show the variation in the critical height of $\gamma H/c$ with the ratio of width to the height of the slope for different slope angles, $\beta$. The values of γH/c computed by genetic algorithms (GAs) are comparable to or even lower than the outcomes determined by the loop searching method and the random searching technique. This suggests that the γH/c value obtained via the GA is closer to the actual value for the slope. The above validation demonstrates that the optimization method, the GA, is more efficient in searching the least upper-bound value of the critical filling height.

To further validate the feasibility of the present methodology for failure surface projections, the critical slip surface obtained from the limit analysis has been compared with that of the failure surface obtained from [28]. The associated parameters have been defined in Figure 6, which are listed as follows: ${\phi }_1=30°,{\phi }_2=40°,\beta =20°$$, H=10 \mathrm{m}, D_{\mathrm{cr}}=1,\mathrm{and}B=57\mathrm{m},B_0=90\mathrm{m}$.

In Figure 6, it also has been noted that the critical slip surface of the present study and failure surface generated by Mohapatra [29] are well matched with each other. The validation of failure surface projections demonstrates that the GA and the three-dimensional rotational mechanism utilized by the present study are efficient in searching the critical failure surface.

5. Parametric Study

The effects of the cohesion of soft soil, the depth–height ratio, and the surcharge loads on the critical filling height of the embankment are first explored. Then, the effect of the depth–height ratio on the failure width of the embankment is investigated. Finally, the projections of the maximum critical failure surfaces for different conditions are plotted. These analyses are conducted by using nondimensional charts to efficiently cover a wide range of possible parameters.

5.1. Soft Soil Cohesion

Figure 7 shows the critical filling height for the base failure of embankments with different values of $c_2/c_1$ and $\beta .$ The internal frictional angle and cohesion of embankment soil are taken as ${\phi }_1=30°$ and $c_1=5\mathrm{kPa}$, respectively. The slope angles of the embankment are taken as $\beta =45°,\beta =60°$, and $\beta =75°$ and the permitted embankment width is $L=5H$. Note that because the cohesions of embankment soil and soft soil are different, the equivalent cohesion, $\overline{c}=\left(c_{1}H+c_{2}d\right)/\left(H+d\right)$, is employed to describe the soil cohesion.

The figure illustrates that the critical filling height initially shows a gradual increase, followed by a sharp rise as the cohesion of soft soil increases. When the soft soil has higher strength, the critical filling height becomes significantly larger. For instance, in the case of ${\phi }_2=5$°, the critical filling height remains below 18 when the ratio of $c_2/c_1$ ranges between 1 and 5, but it increases sharply when $c_2/c_1$ exceeds six. In contrast, when φ₂ is relatively small, the critical slip surface tends to pass through the intersection of the embankment and soft soil interface. As shown in Figure 7a–c, for all the slope angles of the embankment considered, the critical filling height of the embankment over soft soil decreases as the slope angle increases, particularly for soft soils with smaller friction angles. Additionally, the use of nondimensional charts in this analysis allows for efficient representation and comparison of a wide range of parameters, providing a clearer understanding of the relationships among key factors affecting the embankment stability.

5.2. Depth–Height Ratio

The variation in the critical height of $\gamma H/\overline{c}$ in terms of the depth–height ratio, $D_{\mathrm{cr}}$, is plotted in Figure 8. The internal frictional angle and cohesion of embankment soil are also taken as ${\phi }_1=30°$ and $c_1=5\mathrm{kPa}$, respectively. The normalized embankment width, $L=5H$, is adopted to limit the embankment width.

According to an overall observation of the evolution of the critical filling height of $\gamma H/c$, the process of change can be divided into three stages, namely a rapid change stage when $d/H\le 1$, a moderate change stage when $1<d/H\le 1.5$, and a stable change stage when $d/H>1.5$. This indicates that the internal energy should dissipate rapidly to expand failure to the depth of soft soil at the initial stage, resulting in dissipating much more rapidly in the early stage of failure occurrence than the later stage. For instance, for ${\phi }_2=30$° and $\beta =45°$, the critical filling height of the embankment is reduced by 26.65% when $d/H$ is at rapid change stage. However, the critical filling height of the embankment is reduced by 3.57% when $d/H$ is at stable change stage. This indicates that the critical failure surface tends to develop in the weaker soil of the embankment.

5.3. Surcharge Loads

The geometric restriction and the surcharge load against the stability of the embankment over the soft foundation are depicted in Figure 9 and Figure 10, respectively. The embankment width, $L=5H$, is used. The depth–height ratio, $D_{\mathrm{cr}}$, is equal to one. The slope angle $\beta$ changes from 45° to 80°, and the ratio $c_2/c_1$ varies from 1 to 10.

Figure 9 shows that the critical filling height decreases dramatically with the increase in the slope angle of the embankment, $\beta$. For ${\phi }_1=45°$, $\gamma H/c$ decreases from 43.87 to 14.96 when the slope angle $\beta$ increases from 45° to 80°. On the contrary, the value of $\gamma H/c$ increases with the increase in ${\phi }_1$ and ${\phi }_2$ of the embankment. For the embankment with low soil strength, base failure is more likely to occur as the friction angle of the foundation decreases. Thus, the parameters of slope geometry, $\beta$, show a significant impact on the stability of the embankment over soft soil, which can be interpreted by the fact that the flatter inclination of the embankment can provide higher resistance [31].

Figure 10 shows a comparison between surcharge loads $q=100\mathrm{kPa}$ and $q=0\mathrm{kPa}$ for three different slope angles, i.e., $\beta =45°,\beta =60°,$ and $\beta =75$°. When the embankment crest has surcharge loads, the external work rate, $W_{3\mathrm{D}}$, is not only provided by the soil but also surcharge loads; hence, the filling height of the embankment under surcharge loads is more critical. Taking $\beta =45°$ as an example, the value of $\gamma H/c$ decreases by 37% for the cases of $q=100\mathrm{kPa}$ and $q=0\mathrm{kPa}$ as compared with the case $c_2/c_1=1$. This indicates that the surcharge loads have a negative effect on improving the critical filling height of the embankment over soft soil.

5.4. Failure Width

Figure 11 shows the variations in failure width of the embankment with the depth–height ratio, where $c_2/c_1=1$, $\beta =45$°, ${\phi }_1=30°$, and ${\phi }_2=5°$. Since the projection width of the failure mechanism has a finite width, $L$, the critical base failure of the embankment needs to be checked if the embankment has inadequate stability. Therefore, the distribution of the projection width should be evaluated so that the critical filling height is determined.

In Figure 11, it is shown that the ratio of $L_{\mathrm{cr}}/L$ increases from 0 to 0.68 when $D_{\mathrm{cr}}$ increases from 0 to 2 and $L/H$ is equal to five. When $L_{\mathrm{cr}}/L$ is smaller than one, the failure width of the embankment is included within the embankment width. To study the impact of the depth–height ratio on the failure width of the embankment over soft soil, $L/H=15, L/H=10, L/H=15,$and $L/H=20$ are employed to describe the permitted width of the embankment. For $L/H=5$, $L_{\mathrm{cr}}/L$ increases by 254% when the depth–height ratio increases from 0.5 to 1.5. This is not surprising, since the increase in the failure width is associated with the change in the curvature of the failure surface. Thus, the critical filling height indicates that the constraint on the depth–height ratio has a significant impact on the critical failure width of the mechanism.

5.5. Projection of Critical Failure Surface

Figure 12a plots the 3D failure mechanism of the embankment over the soft ground, and the projections of the maximum cross-section of the critical failure surface are shown in Figure 12b,c. The slope angle of embankment, $\beta$, is taken as $45$°. In Figure 12a, $D_{\mathrm{cr}}=1$, ${\phi }_1=30°$, $c_2/c_1=5$, and $c_2/c_1=1$. In Figure 11, $D_{\mathrm{cr}}=0.5$, ${\phi }_1=40°$, $c_1=10 \mathrm{kPa}$, ${\phi }_2=5°,$ and ${\phi }_2=20°$. Figure 12c is the corresponding 3D failure mechanism.

Figure 12a shows that the critical failure surface is flat when the soft soil has stronger strength, i.e., $c_2/c_1=5$. Moreover, the slip surface expands with a decreasing ratio of $c_2/c_1$, and the critical failure mechanism changes gradually from embankment failure to soft soil failure as the ratio of $c_2/c_1$ decreases. Apart from the cohesion of soil, the internal friction angle of soft soil is also an important factor for the critical slip failure, as shown in Figure 12b. It can be found that critical failure surface tends to develop in weaker soil. By comparing a family of the projections of failure surface curves in Figure 12a,b, the failure width of the embankment is affected by the depths of soft soil and the soil strength. The deeper the depths of soft soil, the greater the failure width. However, when the soil strength of the embankment is large, it will weaken the increase in failure width to a certain extent.

6. Conclusions

This study investigated the critical filling height of embankments over soft soil using a 3D upper-bound limit analysis with a kinematically rotational log-spiral failure mechanism. The critical filling height was determined by balancing energy rates, and the highly efficient genetic algorithm (GA) was employed to optimize the failure mechanism. The main findings are summarized as follows:

(a) The proposed failure mechanism effectively extends to scenarios involving soft soil. Solutions obtained using the GA are more stringent compared to loop searching and random searching techniques, with critical slip surfaces aligning well with finite element limit analysis results.

(b) The most critical failure surfaces are tangential to the base of the soft soil layer. The study also revealed a negative correlation between the critical filling height and the ratio of soft soil depth to embankment height, with the critical filling height increasing as this ratio decreases.

(c) Constraining the failure mechanism to embankment dimensions prevents an indefinite expansion of the mechanism and ensures achievable minimum upper-bound values.

(d) Soil strength parameters (cohesion and internal friction), the depth–height ratio, and surcharge loads significantly influence the critical failure surfaces, with embankment failure dimensions increasing as both the depth–height ratio and surcharge loads increase.

(e) Compared to conventional 2D approaches, the proposed 3D upper-bound limit analysis provides more realistic predictions of critical filling heights, better capturing the three-dimensional nature of failure mechanisms. Additionally, the genetic algorithm (GA) enhances computational efficiency, outperforming traditional optimization techniques in terms of both speed and accuracy.

However, the method relies on idealized soil models and predefined failure mechanisms, which may not fully capture the complexities of real-world embankments. Experimental validation and field applications are needed to further verify the robustness of the proposed approach.

To ensure practical applicability, the proposed GA-based method can be adapted to different soil conditions by calibrating soil parameters such as cohesion, friction angle, and unit weight, providing engineers with a flexible and efficient tool for embankment design in diverse geotechnical settings.