Finite element analysis of slope stability by expanding the mobilized principal stress Mohr’s circles - Development, encoding and validation

Djillali Amar Bouzid

Department of Civil Engineering, Faculty of Technology, University of Blida 1, Blida, 09000, Algeria

Keywords:Slope stability, Finite element analysis,Strength reduction method (SRM), Stress point-based factor of safety (FOS), Limit equilibrium method (LEM), Stress deviator,Mohr’s circle, Plastic strain

ABSTRACT In recent years, finite element analysis is increasingly being proposed in slope stability problems as a competitive method to traditional limit equilibrium methods(LEMs)which are known for their inherent deficiencies. However, the application of finite element method (FEM) to slope stability as a strength reduction method (SRM) or as finite element limit analysis (FELA) is not always a success for the drawbacks that characterize both methods.To increase the performance of finite element analysis in this problem, a new approach is proposed in this paper. It consists in gradually expanding the mobilized stress Mohr’s circles until the soil failure occurs according to a prescribed non-convergence criterion.The present approach called stress deviator increasing method (SDIM) is considered rigorous for three main reasons.Firstly,it preserves the definition of the factor of safety(FOS)as the ratio of soil shear strength to the mobilized shear stress. Secondly, it maintains the progressive development of shear stress resulting from the increase in the principal stress deviator on the same plane, on which the shear strength takes place. Thirdly, by introducing the concept of equivalent stress loading, the resulting trial stresses are checked against the violation of the actual yield criterion formed with the real strength parameters rather than those reduced by a trial factor. The new numerical procedure was encoded in a Fortran computer code called S4DINA and verified by several examples. Comparisons with other numerical methods such as the SRM, gravity increasing method (GIM) or even FELA by assessing both the FOS and contours of equivalent plastic strains showed promising results.

1. Introduction

In both geotechnical engineering practice and academia, the problem of slope stability has been a subject of investigation for at least seven decades.Many solutions proposed in accordance to the available computation tools have been enhanced to reach maturity and gained the acceptance of the civil engineering community.

Limit equilibrium methods (LEMs) of slope stability analysis(Fellenius, 1936; Bishop, 1955; Lowe and Karafiath, 1960;Morgenstern and Price, 1965; Spencer, 1967; Janbu, 1968; Sarma,1973), which are based on the method of slices and assume predefined slip surfaces, are widely documented in the literature(Abramson et al., 2001; Cheng and Lau, 2014).

Although LEMs are still considered as design tools for their simplicity and reliable values of factor of safety (FOS), they are subjected to several disadvantages relevant to the assumptions regarding the shape of failure surface and the forces acting at the inter-slice locations.

Over the last two decades, the displacement-based finite element method (FEM) became increasingly popular in analyzing slope stability problems.Its use in this crucial area may have three different forms. The first called the strength reduction method(SRM) (Zienkiewicz et al., 1975; Naylor, 1981; Donald and Giam,1988; Matsui and San,1992; Ugai and Leshchinsky,1995; Dawson et al.,1999; Griffiths and Lane,1999; Cheng et al., 2007; Liu et al.,2015). A brief description of the SRM is given in the next section.The second,less popular than the SRM,is called the gravity increase method (GIM) (Swan and Seo,1999; Li et al., 2009; Sternik, 2013).The GIM consists of gradually increasing the gravitational acceleration(g),while keeping the geomaterial properties unchanged in a finite element process,which continues until the slope reaches the critical failure state.Swan and Seo(1999)compared results of both SRM and GIM and concluded that neither method showed a clear superiority over the other.However,Sternik(2013)found out that,in some circumstances, the GIM may lead to a significant overestimation of the FOS, especially for low slopes. The third form,through which the FEM is employed to analyze the stability of a slope, is termed finite element limit analysis (FELA) (Sloan, 2013;Krabbenhoft et al.,2007;Tschuchnigg et al.,2015a,b).It is currently gaining popularity as it allows rigorous upper and lower bounds of FOS.However,FELA is restricted to problems where only associated plasticity flow rule is allowed.

As far as the slope stability analysis has always been a crucial concern in geotechnical engineering,with no doubt,the FEM is the most appropriate way to deal with this challenging engineering problem. However, the inherent disadvantages regarding the foregoing methods weaken the reliability of the applied methods when it comes to the determination of FOS and slip surfaces.

Starting from the principle stating that the failure in any deformable body is caused by the changes in stresses rather than reducing or increasing soil properties,the new idea in this paper is based on progressively changing the magnitude of stresses in an elaborate way until a situation of failure occurs. The new finite element procedure uses the expansion of the principal stress Mohr’s circle and consequently the increase of the stress deviator for the determination of FOS.The proposed approach takes profit of the linearity of the Mohr-Coulomb criterion which enables a rigorous definition of a stress point-based FOS at any point from the slope domain.Called stress deviator increasing method(SDIM),the new numerical procedure consists of gradually increasing the stress deviator until failure occurs by a non-convergence criterion. The trial factor termed Mohr’s circle expansion factor which controls the Mohr’s circle expansion is considered as the global slope FOS when the iterative process fails to converge. The equations giving the trial values of stresses in terms of the Mohr’s circle expansion factor are described in detail along with the computational methodology which is supported by a flowchart of a Fortran computer program named S4DINA (soil stability study by stress deviator increasing using numerical analysis).

The performance assessment of S4DINA incorporating both slope FOS value and equivalent plastic strains was carried out by considering four examples of homogeneous slopes and an additional example of a slope containing a soft band of frictional material. The comparative study, which involved other methods such as LEM, SRM, GIM and FELA confirmed the reliability and the robustness of the present method.

2. Short description of the SRM

This numerical procedure consists in reducing progressively the strength of soil by a factor termed strength reduction factor FTrial,till an unstable condition results. This unstable situation is evidenced by failure of the solution to converge. The factor FTrial(Fig.1) can geometrically satisfy the following equation:

Fig.1. Definition of in the SRM.

3. New concept: FOS by expanding the mobilized stress Mohr’s circle up to failure

In order to preserve both the mobilized normal stress σmand shear surface orientation identical at equilibrium and at failure and therefore keep the concept of the stress level and the local FOS consistent, the mobilized Mohr’s circle is brought on the verge of failure in such a way that each trial Mohr’s circle crosses the vertical line passing by the point m at a point where the tangent line is parallel to the failure envelope line (Fig. 2).

The innovative concept in this paper is based on two important key conditions which should be satisfied to fulfill the original definition of the FOS. First, the evolution of the mobilized shear stress corresponds to the same mobilized normal stress(σm)at any point in the discretized medium. Second, the shear stress corresponding to the mobilized normal stress should occur on the same plane, which becomes the sliding plane for the ultimate Mohr’s circle.

Starting from the definition of a factor giving the ratio of trial shear stress τTrialto the mobilized shear stress τmSDIM,the purpose of this subsection is to analytically set the necessary equations permitting to evaluate the trial values of the principal stresses.This ratio is expressed as

Fig. 2. Expanding the mobilized stress Mohr’s circle by maintaining the line defining the mobilized shear stress parallel to the failure line.

The main idea is to expand the initial principal stress Mohr’s circle in such a way that the segment O0i remains parallel to the corresponding segment in the subsequent Mohr’s circles (Ott)(Fig.2),ensuring that both the mobilized shear stress τmSDIMand the trial shear stress τTrialoccur on the same sliding plane and holding the definition of the FOS in terms of shear stresses fulfilled.In these conditions, assuming a compression-positive sign convention, the major principal stress should be increased and the minor principal stress should be decreased. The new set of stresses is

Overall, the SDIM proposed in this article can be viewed as a method which preserves the validity of the FOS definition and it is considered robust for at least two main reasons. Firstly, the procedure maintains the progressive development of the shear stress on the same plane, on which the shear strength will occur at the ultimate state. Secondly, the SDIM deals hence with the actual material by employing its real strength parameters (c and φ) and dilation angle ψ rather than those reduced by a factor.

4. The computer code S4DINA

The analytical equations described in the previous section have been implemented in a Fortran computer code called S4DINA.Applied boundary conditions to the soil are pinned supports at the bottom of the slope with no displacements in both horizontal and vertical directions(u = v = 0)and roller supports on both sides of the mesh with no movement in the horizontal direction (u = 0).

Before giving the flowchart of the computer program S4DINA,it is worth providing some valuable information for a better understanding. The slope geomaterial is treated as a visco-plastic material and its behavior is described according to the visco-plastic procedure given by Zienkiewicz and Cormeau (1974). The procedure is explained in detail by Smith et al.(2013).The finite element computational procedure in S4DINA is detailed in the flowchart of Fig. 3.

5. Comparative study involving the FOS and contours of plastic deformation zones

The primary target of this section is to verify the present formulation against the known slope stability methods such as LEM, SRM and any other rigorous method. The values of FOS provided by the SRM were obtained by the use of an existing Fortran computer program code P64 available in Smith et al. (2013), but with slight modifications performed by the author of this paper to create the same conditions and hence a reliable comparison. This program called MSGP64(modified Smith and Griffiths P64)in this paper underwent three principal modifications with respect to the original version of P64. Firstly, the routine for mesh generation in S4DINA was implemented in MSGP64 to exclude any effect of the nature and the number of finite elements on the results.Secondly,the same values of FTrialwere adopted in both S4DINA and MSGP64. Thirdly and similarly to what has been programmed in S4DINA,a supplementary subroutine for computing the equivalent plastic strains was also encoded in MSGP64 for comparing the plastic deformation zones. The LEM FOSs appearing in the comparison examples were obtained using the computer code SLIDE 6.0(Rocscience Inc., 2016). Bishop simplified method (BSM) (Bishop,1955) and Spencer method (SPM) (Spencer, 1967) were chosen for their accuracy in providing perfectly acceptable FOSs for practical purposes (Duncan and Wright,1980).

For the slope examples involving an embankment and a foundation layer, 1602 finite elements were used to model the entire region in both S4DINA and MSGP64(792 for the embankment and 812 for the foundation layer). In the example where only an embankment is modeled, a number of 1392 finite elements were adopted. These meshes were deemed sufficient to give accurate results (Tschuchnigg et al., 2015b). In all the examples considered,Itersmax= 1000 and TOLconv= 0.0001 were utilized to avoid any premature divergence.

Fig. 3. S4DINA flowchart.

Fig. 4. Unsupported vertical cut in a purely cohesive soil.

Table 1 FOS values obtained by different methods for unsupported vertical cut in a soil having cu = 50 kPa.

5.1. Example 1: vertical cut in a purely cohesive soil under undrained conditions

The stability of unsupported vertical cut in a purely cohesive soil under undrained conditions is among the oldest soil mechanics problems.This problem is governed by a dimensionless parameter called stability number Ns= Hγ/cu,where H is the cut height,cuis the deposit undrained cohesion, and γ is the deposit unit weight.Coulomb was the first who developed a solution to determine the critical cut height that an excavation into a purely cohesive deposit may stand without lateral support(Coulomb,1773;Heyman,1973).He proposed a stability number of Ns= 4,on the assumption that the failure of a slope occurs along a plane surface. This stability number correspond exactly to a FOS equal to 1.

In order to see at what extent the proposed method (SDIM) is accurate in dealing with such a problem, a vertical cut of 10 m in height is considered.The soil deposit has an undrained cohesion of cu= 50 kPa and a unit weight of 20 kN/m3to reproduce a FOS of unity.The adopted geometry along with the other soil strength and deformation parameters are showed in Fig. 4.

The values of FOS for both SDIM and SRM along with those of LEM are listed in Table 1,whereas contour lines ofcorresponding to the non-convergence step of computations are illustrated in Fig. 5a and b for SDIM and SRM, respectively.

The examination of the recorded results in Table 1 together with the equivalent plastic zone distribution in Fig. 5 reveals at least three important observations. Firstly, it is quite obvious that both SDIM and SRM yielded the same value of FOS and the same disribution of plastic zones. This is because the final FTrialis around one,and consequently the mobilized stressses in the SDIM remain unchanged (there is no expanding) and the yield criterion in the SRM remains unaltered (there is no reduction in strength).The displayed data are the results of the common code part of both S4DINA and MSGP64. Secondly, all the FOS values are slightly less than 1, proving thus that Coulomb concept of linear failure was faulty and Ns= 4 is not rigorous. Thirdly, the curved shear band is well traced from the slope toe to its top in both SDIM and SRM, and quite different from the slip line assumed by Coulomb.

Taking the average value of the foregoing bounds Ns= 3.7798 and therefore the corresponding cohesion cu= 52.91 kPa, the previous analysis was repeated keeping the other parameters unvaried. The FOS values are displayed in Table 2 whereas the equivalent plastic zone distribution is illustrated in Fig. 6.

The results show clearly that both S4DINA and MSGP64 produced more rigorous solutions.

Fig. 5.contour lines corresponding to the step of failure for a soil having cu = 50 kPa: (a) SDIM; and (b) SRM.

Table 2 FOS values obtained by different methods for unsupported vertical cut in a soil having.cu = 52.91 kPa.

5.2. Example 2: Homogeneous soil slope where the exact FOS is known

Benchmark examples whose exact FOSs are known have been used by many authors to validate their proposed methods(Dawson et al.,1999; Zheng et al., 2005, 2008). The example of this subsection consists of a homogeneous slope,upon which only the gravity was acted, with geometrical configuration and soil properties as illustrated in Fig. 7. With these parameters, the slope has a FOS of exactly 1 according to the limit analysis solution of Chen (1975)who assumed a log-spiral slide curve.

In order to bypass the influence of cohesive-frictional condition(Zheng et al.,2005),a Poisson’s ratio of 0.4 along with an associated flow rule was adopted in both S4DINA and MSGP64. The results drawn from the literature were also based on an associated flow rule.Dawson et al.(1999)employed FLAC to compute the FOS.The results of analysis are showed in Table 3,and the plastic zones are illustrated in Fig. 8.At first glance,the results reported in Table 3 showed that both SDIM and SRM are extremely accurate and even more precise than FLAC simulations performed by Dawson et al.(1999).

Fig. 7. Geometry and material data for a slope having a FOS equal to1.

Table 3 FOS values obtained by different methods for the homogeneous slope shown in Fig.7.

5.3. Example 3: Homogeneous soil slope subjected to an external surcharge

Fig. 6.contour lines corresponding to the step of failure of a soil having cu = 52.91 kPa: (a) SDIM, and (b) SRM.

Fig. 8. contour lines corresponding to the step of failure for an associated flow rule (φ = ψ): (a) SDIM; and (b) SRM.

Fig. 9.contour lines corresponding to the step of failure for an associated flow rule (ψ = φ): (a) SDIM; and (b) SRM.

Fig.10. contour lines corresponding to the step of failure for a non-associated flow rule (ψ = 0°): (a) SDIM; and (b) SRM.

Fig.11. Geometry of the slope and soil parameters considered in the example 3.

In order to assess the performances of the GIM with respect to the SRM,Sternik(2013)studied a homogeneous cohesive-frictional soil slope subjected to an external distributed load whose geometry is illustrated in Fig.11.The author considered two sets of materials according to two slope inclinations, namely β = 30°and β = 45°,as shown in the table of Fig. 11. The author carried out a limited parametric study where four values of the slope height were given,from which heights of H = 6 m and H = 12 m were selected for this study.In all analyses,both associated and non-associated flow rules have been assumed.Furthermore,it is worth to note that only non-associated plasticity calculations were considered by Sternik(2013) in the GIM and the case of slope without surcharge was not examined either.

Before examining the effect of surcharge loading on the slope FOS,the slope was considered first without surcharge.SDIM results along with the other methods are illustrated in Table 4.For lack of space,contour lines(Figs.12-19) are restricted to only the case of non-associated plasticity flow rule.

Table 4 FOS values obtained by different methods applied to the example 3 without surcharge (q = 0).

Fig.12. contour lines corresponding to the step of failure (non-associated plasticity flow rule) for β = 30°, H = 6 m, and q = 0: (a) SDIM; and (b) SRM.

Fig.13. contour lines corresponding to the step of failure (non-associated plasticity flow rule) for β = 30°, H = 12 m, and q = 0: (a) SDIM; and (b) SRM.

For the results of H = 12 m illustrated in Fig.15, the same observations are noticed here.The SDIM results in a clear shear band with an appropriate FOS, whereas SRM exhibits a large plastically damaged zone with an overestimated FOS.

Now it was the custom in those days for princes and princesses to be brought up by fairies, who loved them as their own children, and did not mind what inconvenience they put other people to for their sakes, for all the world as if they had been real mothers

Fig.14. contour lines corresponding to the step of failure (non-associated plasticity flow rule) for β = 45°, H = 6 m, and q = 0: (a) SDIM; and (b) SRM.

Fig.15. contour lines corresponding to the step of failure (non-associated plasticity flow rule) for β = 45°, H = 12 m, and q = 0: (a) SDIM; and (b) SRM.

Fig.16. contour lines corresponding to the step of failure (non-associated plasticity flow rule) for β = 30°, H = 6 m and q = 90 kPa: (a) SDIM; and (b) SRM.

Results of all methods involved in the analysis of example 3 under surcharge loading are listed in Table 5.

Fig.17. contour lines corresponding to the step of failure (non-associated plasticity flow rule) for β = 30°, H = 12 m, and q = 90 kPa : (a) SDIM; and (b) SRM.

Fig.18.contour lines corresponding to the step of failure (non-associated plasticity flow rule) for β = 45°, H = 6 m, and q = 90 kPa: (a) SDIM; and (b) SRM.

Fig.19. contour lines corresponding to the step of failure (non-associated plasticity flow rule) for β = 45°, H = 12 m, and q = 90 kPa: (a) SDIM; and (b) SRM.

From the close examination of Table 5,it is easy to make at least three observations. Firstly, and most importantly, the FOS values provided by the GIM present a significant overestimation regardless the slope inclination or the slope height.This shows why GIM was abandoned as a slope design method in favor of other more rigorous approaches. Secondly, the surcharge loading, although significant,has a little influence on the FOS in comparison with the slope without external loading. Thirdly, the FOSs obtained by the SDIM are close to those of LEM for the four cases. However, SRM FOSs are slightly less than those of SDIM,except the case of β = 45°and H = 12 m for which FOSSDIMis slightly lower than FOSSRM.

The different distributions offor the analyzed cases are depicted in Figs.16-19.On one hand,in SDIM(Figs.16a,17a and 18a and 19a),the damage areas induced by plastic strains form beneath the location of external loading and spread locally to a limited depth forming irregular plasticity bulbs in the case of a steep slope(β = 45°) or reaching the bottom boundary in the case of a flat slope. The other zones remain undisturbed. On the other hand, inthe SRM (Figs. 16b, 17b and 18b and 19b), spurious plastic strain zones appear in irregular shapes beneath the loaded zone and spread towards all the slope boundaries, especially for flat slopes(β = 30°). The entire slope regions are disturbed by the adopted computational process in the SRM.

Table 5 FOS values obtained by different methods applied to the example 3 with surcharge (q = 90 kPa).

Fig. 20. Soil slope considered in the example 4.

Table 6 FOS values obtained by different methods applied to the example 4.

5.4. Example 4: Homogeneous soil slope where FELA FOS is available

Tschuchnigg et al. (2015b) studied a slope stability problem containing only an embankment and considered several geometrical configurations in which only the slope inclination was given a variety of values.Two angles,β = 45°and β = 60°,were retained for the comparative study performed in the present analysis. The geometrical details and the soil properties are presented in Fig.20.This example was chosen as it can serve as a strong assessment tool for the SDIM developed in this paper for the availability of FELA results.

Results of the comparative analysis are reported in Table 6. At first glance, FOSSDIM= 1.53 for an associated flow rule is in good agreement with the upper bound FOS value of FELA for β = 45°,and FOSSDIM= 1.17 is exactly within the range of the two bounds provided by FELA for the case of more steep slope with β = 60°.However,FOSSDIM= 1.42 for β = 45°and FOSSDIM= 1.11 for β =60°,which seem outside the range of the bounded interval defined by FELA, are reasonable values because FELA is not applied for a non-associated plasticity flow rule. Tschuchnigg et al. (2015b)employed Davis approach in this case.

5.5. Example 5: Slope in cohesive-frictional soil with a weak layer

This example described in Fig.23 deals with a slope of 12.25 m in height and an inclination β of 26.6°supported by 8 m thick layer of similar material(c = 28.5 kPa, and φ = 20°).A thin weak layer(c = 0 kPa, and φ = 10°) with 0.5 m thickness is horizontally located 0.75 m below the slope toe surface.The unit weight of both materials was assumed to be γ = 18.84 kN/m3and the Young’s modulus was taken as E = 100 MPa.

This challenging problem has been considered by many researchers in the past,such as Sloan(2013)and Zhou and Qin(2020).In both papers,the FELA combining the conventional FEM with the upper and lower bounds of classical plasticity theorems was adopted. Since FELA is restricted to an associated plasticity flow rule,both S4DINA and MSGP64 were run adopting a dilation angle(ψ)equal to the friction angle(φ).Furthermore,two different values of Poisson’s ratio which are ν = 0.32 and ν = 0.4 were selected in the present study. Results of comparison are gathered in Table 7.

Fig. 21.contour lines corresponding to the step of failure by the present method (SDIM) adopting an associated flow rule: (a) β = 45°; and (b) β = 60°.

Fig. 22.contour lines corresponding to the step of failure by the present method (SDIM) adopting a non-associated flow rule: (a) β = 45°; and (b) β = 60°.

Fig. 23. Soil slope considered in the example 5.

Table 7 FOS values obtained by different methods applied to the slope with homogeneous layer crossed by a weak layer.

Fig. 24. contour lines corresponding to the step of failure adopting an associated flow rule and ν = 0.32: (a) SDIM; and (b) SRM.

Despite the fact that both SDIM and SRM furnished very close FOSs,different performances in terms ofcan be observed when it comes to the distribution of plastic strains over the slope area, as illustrated in Figs.24 and 25 for ν = 0.32 and ν = 0.4,respectively.In these figures, the mode of failure is dominated by intense shear deformations in the weak layer of cohesionless material.The effect of Poisson’s ratio on the outcomes of SRM is clearly visible in Fig.24b by creating a large zone of plastic deformations which extends to the left and bottom boundaries, whereas its influence on the SDIM is much less pronounced. However, in Fig. 25, the erratic zones of plasticity disappear and both methods exhibit a well traced shear bands. Diversely, Sloan (2013) noticed that a secondary failure mechanism occurs along a plane at a right angle to the slope face,the soil beneath the slope face remains intact in both SDIM and SRM.

Fig. 25. contour lines corresponding to the step of failure adopting an associated flow rule and ν = 0.4: (a) SDIM; and (b) SRM.

Although LEMs for the slope stability analysis have been widely and successfully used, this example is a real challenge for LEMs.Two main drawbacks are encountered.Firstly,the slip surface is not detected automatically, and the users have to specify whether the failure surface is circular or not.Secondly,with each chosen surface shape, the method of surface research should be indicated (simulated annealing for example). If a circular surface is chosen, a significant overestimation in the FOS is obtained,and all LEMs yield a value above 1.5 (Table 7). However, if a non-circular shape is selected, even the most rigorous LEMs (Bishop simplified, and Spencer methods) fail to give accurate results. Indeed, significant variations are noticed in the FOS values. SLIDE 6.0 gives the FOS values of 1.184 and 1.221 for Bishop simplified and Spencer methods, respectively. These discrepancies were also discussed in detail by Sloan (2013).

In order to assess the stability of a slope with three different layers, the properties of the weak layer and those of the top layer were kept unchanged whereas the strength properties of the bottom layer were increased. The computation of the FOS indicated that increasing the bottom layer properties has practically no influence on the value of FOS which remains 1.27. This is because in this special problem, the onset of shear band occurs in the weak layer and spreads upward in the top layer. Upon this observation, the strength properties of the top layer were changed to c = 48.5 kPa and φ = 26°; while the characteristics of the weak and the bottom layers were kept constant, as indicated in Fig. 23. Taking a Poisson’s ratio of ν = 0.4 and adopting an associated plasticity flow rule, results of comparison are summarized in Table 8.

It is clear that the FOS values provided by both the SDIM and SRM are accurate although the latter slightly underestimates the FOS. However, results of the two options of LEM suggest two different tendencies: a large overestimation of the FOS when acircular slip surface is chosen,whereas an underestimation when a non-circular slip surface is selected. These results warn the engineers to be cautious when using LEM in this kind of slopes.

Table 8 FOS values obtained by different methods applied to a slope with two different layers sandwiching a weak layer.

6. Conclusions

Among the most known FEMs to deal with the slope stability problems, the SRM has gained a wide acceptance within the geotechnical community for its easiness to implement in a computer code and the stability of its results when dealing with ordinary slope configurations. However, a close inspection to the SRM background formulation reveals that the stress point-based FOS has been built on an inconsistent definition compared to that of LEM,and consequently the global FOS is affected. Furthermore, the numerical instabilities observed in some situations such as those involving steep slopes with relatively higher angles of friction for instance, may diminish the reliability of the SRM.

A new finite element procedure called SDIM for the assessment of FOS of a slope is proposed in this paper. It consists of progressively increasing a factor called Mohr’s circle expansion factor until the slope failure is reached,thus taking the reverse path of the SRM which consists of reducing the soil strength parameters.The slope failure is considered to occur when the iterative process fails to converge within the prescribed range of the maximum number of iterations.

Bearing the acronym of S4DINA, a Fortran computer program has been written for the analysis of slope stability.Performances of S4DINA were compared with those of another computer program performing the SRM computations, taking obviously the conventional LEMs as reference.The comparison analysis was restricted to the evaluation of FOS and to the distribution of the equivalent plastic strains on the discretized domain of the slope.Five examples involving cohesive-frictional soils were considered to assess the capabilities of S4DINA. Results of comparison confirmed the consistency of the present method,especially when compared to FELA methods,and showed that the SRM is sensitive to a non-associated plasticity flow rule and more prone to the erratic zones of plastic flow.

The proposed method (SDIM) is considered robust and reliable for at least five reasons:

(1) SDIM has been rigorously formulated as it preserves the validity of the FOS definition by imposing a correct stress path to reach failure and maintains the gradual development of the shear stress on the same plane on which the shear strength will occur at failure.

(2) SDIM employs the real soil strength parameters(c and φ)and ψ rather than those reduced by a factor. Consequently, the nodal displacements in the SDIM are more accurate and the slope deformation in each FTrialis well defined.

(3) SDIM avoids all the shortcomings of the SRM and it is convenient for slope problems necessitating non-associated flow rule of plasticity.

(4) SDIM has the potential to be extended to three-dimensional(3D) finite element analysis of slope stability problems.

(5) SDIM is convenient for problems where LEM fails to produce appropriate slip surface, and for problems such as slopes with a weak frictional layer.

Declaration of competing interest

The author declares that he has no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.