Structural Analysis of Backfill Highway Subgrade on the Lower Bearing Capacity Foundation Using the Finite Element Method

)e present study is to investigate the stability of the backfill subgrade on the lower bearing capacity foundation. A finite element (FE) model was constructed to simulate the high-filled road subjected to the actual self-weight load. )e strength reduction method was adopted to establish an analysis model of slope stability. At the same time, the sensitivity analysis of the factors affecting the slope stability was carried out through parametric studies, including the elastic modulus, cohesion, internal friction angle, and slope rate. )e results showed that the slope stability analysis model established by the strength reduction method can characterize the stability of the slope by calculating the slope safety coefficient. )e mutation point of the relationship curve between the displacement generated in the slope and the reduction coefficient can be used as the criterion. Under the condition of a given strength reduction coefficient, the calculated results obtained through FE modeling can show the development of the equivalent plastic zone in the form of cloud diagrams.


Introduction
Road slope refers to the infrastructure constructed or modified to meet the needs of road construction [1][2][3][4][5]. e slope formed under certain conditions (topography and geology) will break the previous mechanical balance into an unstable body due to changes in external factors. Under the action of its own weight or other loads, it will follow a certain relatively weak surface [4,6,7]. A landslide is a kind of bad geological phenomenon such as sliding downwards as a whole, intermittently slowly, and sometimes even suddenly. e factors that cause landslides include excavation of side slopes, rapid embankment filling, earthquakes, river erosion, sudden drops in reservoir water levels, and heavy rains [6,[8][9][10]. e slope is not only a special geological environment but also an important part of engineering construction. Due to the construction demand for a large number of highway projects in recent years, the stability of the slope during the construction process is a link that cannot be ignored [11][12][13][14].
In the past, there were various approximate methods and specialized methods, which can be used to calculate the global safety factor that depends on the slope height, the steepness, and the constitutive properties of the soil constituting the slope. Nash [15] conducted a comprehensive review of many classic slope stability analysis methods. e limitation of many of these classical methods is that Mohr-Coulomb shear strength behavior with a fixed friction angle is usually assumed because of the inherent application of simple statistical methods to the soil. Approximately, the calculation of stress and related shear strength is a continuous mechanical problem with static uncertainties [16][17][18]. To reduce manual trials and errors associated with these technologies, such as slicing methods, many engineering software packages have been developed, but most software packages still use simple static methods to approximate soil stress and strength [19][20][21][22].
In the FE analysis, the method of increasing the external load or reducing the material strength of the rock-soil mass is the earlier method adopted [23,24]. e stability analysis of these slopes mainly includes the judgment basis of slope instability, the stability of the excavated slope, the comparative analysis of the selection of different software, and the analysis of the internal force of the slope with the supporting structure [25,26]. In recent years, the FE method has been increasingly used to predict the displacement and stress in statically indeterminate slopes, dams, and embankments [5,13,14,[27][28][29][30][31][32][33][34][35][36][37][38][39][40][41][42][43][44][45].
e potential attractiveness of calculating slope stability in a continuum/finite element framework lies in many previous studies [7]. e equilibrium stress, strain, and related shear strength in the soil have been calculated accurately in the studies conducted by Aryal and Duncan [41,44]. e general material models of soil (including Mohr-Coulomb and many other models) have been employed in the studies conducted by Ardah et al. [46]. e FEM method has been applied in 2D or 3D with complex slope combinations and soil sediments to model almost all types of mechanisms [47]. e FEM method has also been extended to address the damage caused by leakage, brittle soil behavior, random "yield soil" properties, and interventions such as engineering geotextiles, soil nails, drainage ditches, and retaining walls [5,38,40,[42][43][44][45]. To fully understand the stability condition of the backfill soil on the high-filled road with lower foundation bearing capacity, this article adopts the classic slope stability analysis method. In this method, the inclined soil is subjected to the actual self-weight load, and the purpose is to find a continuous surface through the soil that has the smallest safety coefficient against the slip or shear failure. e safety coefficient is defined as the degree to which the shear strength of the soil is reduced so that the slope has just reached the critical failure state.

Research Objective and Overview
e present study is to investigate the stability condition of the backfill soil on the high-filled road with lower foundation bearing capacity, as shown in Figure 1. A FE model was constructed to simulate the high-filled road subjected to the actual self-weight load. e strength reduction method was adopted to establish an analysis model of slope stability. At the same time, the sensitivity analysis of the factors affecting the slope stability was carried out through parametric studies, including the elastic modulus, cohesion, internal friction angle, and slope rate. Under the condition of a given strength reduction coefficient, the simulation results obtained by ABAQUS software show the development of the equivalent plastic zone in the way of cloud maps. e slope is stable under the reduction factor when the failure of transfixion or local yield failure occurs. If the reduction coefficient continues to increase, the slope is in a critical failure state when the plastic strain between adjacent iterative steps increases too much or the equivalent plastic strain reaches the yield limit at some nodes, and the strength reduction coefficient is positioned as the minimum safety factor of the overall stability of the slope.

Model Construction
In the present study, the FE method was adopted. In this method, a structure is regarded as a whole composed of a finite number of elements through nodes [5,14,30]. Except for the nodes which are fixed on the boundary, the displacement of each node that can produce displacement can be calculated by using equilibrium conditions, and then the internal forces of each element can be calculated by the node displacement [33][34][35].
In the finite element analysis of the strength reduction of the slope, the stability of the slope usually adopts the nonconvergence of the solution as the failure criterion [48]. Within the maximum number of iterations, if the calculation fails to converge, it means that no stress distribution that can satisfy both the failure criterion and the overall balance is found, which means that the soil has been damaged. Specifically, the actual strength parameters c ′ and φ ′ of the backfilled roadbed are divided by a reduction factor F trial at the same time to obtain a set of new values of c ″ and φ ″ after reduction, i.e., en, the reduced values of c ″ and φ ″ are taken into the trial calculation of finite elements as new material parameters [49,50]. When the finite element calculation converges, the value of F trial is slightly larger, and then the trial calculation is carried out until the finite element no longer converges, which indicates that the soil reaches the critical limit state and the slope shear failure occurs. At this time, the critical slip surface and safety factors are obtained. e gravity load is determined by the following method [51,52]. On each one of the finite elements, the gravitational load generated by the weight of the soil can be obtained by the following equation: in which S is the area of the element; e is the element number. e result of this integration is to take the product of the area of each element and the weight of the soil as the element's gravity load and then distribute it to each node. Various complex constitutive models can be considered in the FE analysis, but the most common in engineering analysis is the ideal elastoplastic model because the results of the ideal elastoplastic model are the most comparable to the results of the limit equilibrium method. For general road slopes, the shear failure of the soil is mainly controlled, and the calculation accuracy of the distribution and size of the plastic zone is relatively high, and the requirements for the displacement are relatively low.
In the present study, the Mohr-Coulomb plastic model was employed to characterize the constitutive property of the soil. e yield criterion of the Mohr-Coulomb model is assumed as follows. When the shear stress acting on a certain point is equal to the shear strength of that point, the point will be broken, and the shear strength has a linear relationship with the normal stress acting on the surface. e yield surface equation of the Mohr-Coulomb model [53,54] can be presented as follows: in which φ(θ, f α ) is the friction angle on the meridian plane; c(ε − pl , θ, f α ) represents the changing process of the cohesive force of the material according to the isotropic hardening model. p is the equivalent compressive stress. q is the Mises equivalent stress. R mc is the deviatoric stress coefficient of the Mohr-Coulomb model [55,56], which can be defined as in which φ is the bevel angle of the yield surface of the Mohr-Coulomb model on the P-RMCQ plane, and it generally refers to the internal friction angle of the material. θ is the direction angle of the generalized shear stress.

Material Properties.
Four types of civil engineering materials were adopted in the present simulation [48]. Table 1 gives the material properties.
e bounding and unbound pavement materials (including the surface layer, base layer, and subbase layer) were simplified to one type of material property [57]. Elastic behavior without plastic deformation characterized by Young's modulus and Poisson's ratio was adopted for the pavement materials. e density of the backfill soil was determined as 1800 kg/m 3 . Regarding the Mohr-Coulomb parameters of the backfill soil, the cohesion was determined as 15, 20, and 25 kPa, respectively. e internal friction angle was selected as 25°and 30°, and two elastic moduli of 50 MPa and 100 MPa were determined. e foundation with the lower bearing capacity was characterized by silty clay and mud clay [48]. e corresponding Mohr-Coulomb parameters are shown in Table 1 as well. Figure 2 gives the physical size of the whole model, including the pavement, subgraded by the backfill soil, and the foundation with lower bearing capacity [48]. ese physical sizes were determined according to earlier studies [36][37][38][39][40][41][42]. e slope angle α is an important parameter considered in the present study, and it varied for different values (tanα  Figure 3 gives the FE meshing for different slopes used in the present study. Element selection is the foundation for the FE analysis. For the same FE model, different calculation results can be obtained by using different elements. Regarding the FE model in the present study, the solid element was the optimized one to simulate one part of the whole model. Since a solid element can be connected to other elements through any of its surfaces, it can construct almost any shape and bear any load in the FE model.

Meshing and Element Selection.
To reduce the computations and improve the calculating accuracy, the 8-node plane-strain-reduced (CPE8R) element was adopted in the present FE simulation. It is a general plane strain element, in fact, a special plane stress element, which can be used to simulate a section of a very long structure, such as a dam. It is understood that the plate or shell element may be more efficient to solve the plane strain problem considering the high efficiency during the application. However, it should be noted that the research objective in the present study is about the subgrade and road pavement, which are not satisfying the premise of plate or shell theory. erefore, the solid element, CPE8R, was employed in the present study. Also, it should be noted that the structure containing the plane strain element should be defined in the global X-Y plane, that is, for all nodes z � 0.

Effects of Slope
Angle. e soil was defined by the Mohr-Coulomb model, which is the so-called rigid plastic model or Saint-Venant model. Regarding the Mohr-Coulomb model, when the stress of the soil is less than the yield stress, the soil will not deform, just like a rigid body. When the stress reaches the soil yield stress, the plastic deformation will increase until the soil failure occurs. erefore, the equivalent plastic strain (PEEQ) should be evaluated first, considering the proposed models. Figure 4 gives the results of PEEQ when the slopes are 1 : 1, 1 : 1.25, 1 : 1.5, 1 : 1.75, and 1 : 2, respectively.
As shown in Figure 4, under the conditions of a given strength reduction coefficient, the calculation results obtained through the FE simulation can show the development of the equivalent plastic zone in the form of cloud diagrams. When the slopes were 1 : 5, 1 : 1.75, and 1 : 2, respectively, failure to achieve penetration or local yield failure indicates that the slope is stable under this reduction factor, while when the slopes equal 1 : 1.25 and 1 : 1, the plastic deformation between adjacent iteration steps at some nodes increases too much or the equivalent plastic strain reaches the yield limit. From the displayed results (Figure 4), the plastic zone has penetrated the top of the slope, and the equivalent plastic strain and displacement have an infinite development trend, with obvious sudden changes, indicating that the slope is already in a critical failure state at this time.  Also, it can be observed that when the slopes are 1 : 1, 1 : 1.25, 1 : 1.5, 1 : 1.75, and 1 : 2, the maximum PEEQ are 1.058e-2, 9.026e-3, 8.999e-3, 6.556e-3, and 6.807e-3. Basically, as the slope decreases, the maximum PEEQ decreases, as shown in Figure 4. It can be indicated that reducing the slope can improve the stability of the road slope. Also, it can be observed that when the slope is lower than 1 : 1.75, the maximum PEEQ tends to converge.

Effects of Cohesion and Internal Friction Angle.
It can be observed that, as the cohesion equals 15 kPa, the maximum PEEQ for internal friction angles 25°and 30°are 1.058 × 10 −2 and 1.022 × 10 −2 , respectively. A higher internal friction angle can increase the value of maximum PEEQ. As the cohesion equals 20 kPa, the maximum PEEQ when the internal friction angles equal 25°and 30°are 1.019 × 10 −2 and 1.023 × 10 −2 , respectively. However, as the cohesion increases to 30 kPa, the maximum PEEQ when the internal friction angles equal 25°and 30°are still 1.019 × 10 −2 and 1.023 × 10 −2 , respectively. is means that, as the internal friction angles equal 25°and 30°, even the cohesion increases, and the maximum PEEQ will remain the same.
is result can also prove that 20 kPa is the split point of the maximum PEEQ. It can further be observed that the increase rate of the maximum PEEQ between 15 kPa and 20 kPa is very small (about 0.1%). erefore, it can conclude that higher cohesion and internal friction angle can decrease the maximum PEEQ, but when the cohesion or internal friction angle is higher than one constant value, the maximum PEEQ will remain the same. e cohesive force of the backfill has a great influence on the overall stability of the slope. With the increase of the cohesion of the backfill, the safety and stability factor also greatly increases, indicating that when the cohesive force of the subgrade soil is high, the stability of the roadbed slope is also high. It can be seen from Figure 5 that the cohesion is increased by 5 kPa, and the stability coefficient can be Advances in Civil Engineering increased by about 10%. When the cohesion is large, the contribution to the safety and stability coefficient gradually decreases. It can also be seen that if the backfill material with greater cohesion is selected, the stability of the roadbed slope can be significantly improved. e safety and stability factor remains unchanged as the cohesive force of the weak foundation changes. e main reason may be related to the high fill height of the roadbed. Because when the height of the roadbed is large, the additional stress generated by the external load has a limited range, that is, the working area of the roadbed does not penetrate deep into the soil foundation, the influence of the change of the soil foundation material parameters on the results can be almost ignored. rough the above analysis, we can know the importance of subgrade soil cohesion in highway engineering, and foundation cohesion has little effect on stability.

Effects of Elastic Modulus.
According to the strength reduction theory, when the finite element calculation does not converge due to the reduction of strength parameters, the slope shear failure occurs. erefore, the reduction coefficient of the strength parameter corresponding to the last convergence calculation can be defined as the safety coefficient of the slope. Based on this determination criterion, Figure 6 gives the relationship between the elastic modulus and safety coefficient when the cohesion and internal friction angle equal 15 kPa and 25°, respectively.
It can be found from Figure 6 that the elastic modulus of the backfill soil has little effect on the safety and stability coefficient of the slope. In the case of a lower modulus of elasticity, the safety and stability coefficient of the slope is higher than that of the high modulus of elasticity. e possible reason is that the plastic expansion zone is generated, and when the elastic modulus of the subgrade that occupies most of the entire embankment is small, the overall settlement displacement of the structure increases, which causes the center of gravity of the plastic zone to decrease, which will reduce the sliding force. On the contrary, the safety and stability coefficient of the slope is increased, and this kind of situation cannot be attributed to the improvement of the safety and stability of the slope. Figure 7 gives the safety coefficients of the road slope for varying cohesion and internal friction angle. e slight effect of the elastic modulus on the safety coefficient of the road slope can also be found since the slope safety coefficients were almost the same for different elastic moduli. However, it should be noted that the slope safety coefficient was obviously improved when the internal friction angle increased from 25°to 30°. e internal friction angle is the representative of the internal friction of the soil,   Advances in Civil Engineering which specifically includes the surface friction between the soil particles and the bite force generated by the interlocking action between the soil particles. e larger the internal friction angle, the greater the internal friction force, so as the internal friction angle increases, the maximum equivalent plastic strain value becomes smaller, resulting in a higher value of the slope safety factor. Figure 8 gives the safety coefficients of the road slope when the cohesion was 15 kPa and 20 kPa, respectively.
A slight effect of the elastic modulus on the safety coefficient of the road slope can also be found because of the similar values of slope coefficients. Besides, the cohesion increased the values of the slope safety coefficient. Various physical and chemical forces between soil particles (6)     Advances in Civil Engineering determine the value of soil cohesion, including Coulomb force, static force, van der Waals force, and cementation. e greater the cohesive force, the greater the suction force between the particles so that only when the external force is given greater pulling force can plastic strain occur, and the corresponding antisliding force is greater.

Conclusions
According to the current construction requirements of the process of backfilling the subgrade on the weak foundation, this paper conducted a FE analysis on the slope of this special engineering condition based on the commonly used strength reduction method in the stability analysis of the road slope.
To accurately investigate the influence of slope and the mechanical properties of the soil on the stability of the backfill subgrade during the construction process, the slope stability model based on the strength reduction method is established in this study. At the same time, the influence analysis of the slope, cohesion, and internal friction angle on the slope stability was carried out. e research process can highlight the following conclusions: (1) e slope stability analysis model established by the strength reduction method can characterize the stability of the slope by calculating the slope safety coefficient. In the analysis process, the mutation point of the relationship curve between the displacement generated in the slope and the reduction coefficient is used as the criterion. (2) As the slopes were lower than 1 : 1.5, failure to achieve penetration or local yield failure indicates that the slope is stable under this reduction factor, while as the slopes were higher than 1 : 1.5, the plastic 40 (c=15kpa; φ=25°)    zone has penetrated the top of the slope, and the equivalent plastic strain and displacement have an infinite development trend, indicating a critical failure state. As the slope decreases, the maximum PEEQ decreases; reducing the slope can improve the stability of the road slope. When the slope is lower than 1 : 1.75, the maximum PEEQ tends to converge. Higher values of internal friction angle and cohesion can increase the value of maximum PEEQ, but when the cohesion or internal friction angle is higher than one constant value, the maximum PEEQ will remain the same. e effects of the elastic modulus of the backfill soil on the safety and stability coefficient of the slope are slight.
Regarding future development, a more precise 3D FEM model should be proposed by considering the dimensional sizes in the longitudinal direction. e reasonable physical sizes of the FEM model should also be optimized and determined in the future. Comparison results between the 2D and 3D models should also be included to evaluate the accuracy and efficiency.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that they have no conflicts of interest to report regarding the present study.