An analytical solution for analysis of toppling-slumping failure in rock slopes

An analytical solution for analysis of toppling-slumping failure in rock slopes

Journal Pre-proof An analytical solution for analysis of toppling-slumping failure in rock slopes Hadi Haghgouei, Ali Reza Kargar, Mehdi Amini, Kamra...

5MB Sizes 2 Downloads 259 Views

Journal Pre-proof An analytical solution for analysis of toppling-slumping failure in rock slopes

Hadi Haghgouei, Ali Reza Kargar, Mehdi Amini, Kamran Esmaeili PII:

S0013-7952(19)31542-X

DOI:

https://doi.org/10.1016/j.enggeo.2019.105396

Reference:

ENGEO 105396

To appear in:

Engineering Geology

Received date:

13 August 2019

Revised date:

28 October 2019

Accepted date:

7 November 2019

Please cite this article as: H. Haghgouei, A.R. Kargar, M. Amini, et al., An analytical solution for analysis of toppling-slumping failure in rock slopes, Engineering Geology (2019), https://doi.org/10.1016/j.enggeo.2019.105396

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2019 Published by Elsevier.

Journal Pre-proof An analytical solution for analysis of toppling-slumping failure in rock slopes Hadi Haghgouei1 [email protected], Ali reza Kargar1 [email protected], Mehdi Amini1,* [email protected],

Kamran Esmaeili2 [email protected]

1

School of Mining Engineering, College of Engineering, University of Tehran, Tehran, Iran

2

Department of Civil & Mineral Engineering, University of Toronto, Ontario, Canada *

of

Corresponding author at: School of Mining Engineering, College of Engineering, University of Tehran, Northern Kargar Street, Tehran, Iran

ro

Abstract:

When a slope consists of a horizontal weak rock or soil overlaid by some vertical strong slender rock

-p

columns, the slope is prone to a secondary toppling failure which is known as “toppling-slumping” failure. In the slope, rock columns impose pressure on the lower continuous weak rock or soil mass which leads to

re

a differential settlement at the base of each rock column or a circular shear failure in the continuous mass. If the differential settlements reach a special threshold, rock column overturn and toppling-slumping will

lP

happen. However, the circular shear failure of the soil or weak rock mass which is a well-recognized failure mechanism will occur when shear stress in the slope exceeds its shear strength. In this paper, an attempt has

na

been made to a better understanding of the “toppling -slumping” failure mechanism through developing an analytical solution. For this purpose, first an Airy stress function is acquired for the lower continuous soil mass and its stress distribution and displacement are accordingly determined analytically due to the

ur

overlaid rock columns loads. Then, differential settlements at the base of all rock columns and shear stresses in the soil or weak rock mass are computed and compared with the toppling and shear failure

Jo

criteria to assess slope failure. Finally, several design-charts are developed that determine the stability of a slope against toppling-slumping failure mechanism based on geometric characteristics of the slope and shear strength of its constituent material. Also, a Mathematica package code was developed to evaluate the stress distribution within the slope, and the safety factor. The results show that the vertical stress distribution plays a significant role in the settlement of the rock columns. Results also manifested that at a certain distance from the slope crest, the set of rock columns will act like a single column on a half-plane. The shear or toppling failure of the slope not only depend s on the material’s strength but also on the geometry and the normalized column distance. Keywords: Slope Stability, Secondary toppling failure, Toppling-slumping, Airy stress function, Differential settlement, Sliding.

Journal Pre-proof List of symbols A = Parameters of transformed airy stress function a = Distance between first rock column and slope crest a = Distance between the crest and the last point of loading a/x= Normalized column distance B, C = Parameters of transformed airy stress function c = Cohesion c = Real part of the complex number in linear integration path 𝐶2𝑛, 𝐶2𝑛−1 = Parameters of Fillon-Simpson numerical integration method D = Parameters of transformed airy stress function 𝑒𝑖𝑗 = Strain tensor E = Young’s modulus F(z) = Loading function in transformed space

ro

-p

lP

K, L = Parameters of radial and hoop stress function p = Uniform pressure R= Parameter of radial and hoop stress function r = radial in polar coordinate

  z  = Airy stress function in transformed space  = friction angle  = Density

 ij = stress tensor

 ij = Kronecker delta  = Slope angle  = Poisson’s ratio ,, ,  = Parameters of Filon numerical integration

method  = displacement difference  = Pi number

Jo

ur

na

S 2 n , S 2 n 1 = Parameters of Fillon-Simpson numerical integration method T= Parameter of radial and hoop stress function u = Step function U= Truncated upper bound in numerical integration

1. Introduction

 r ,  ,  r = Stress component in polar coordinate

re

J 2 , J 3 = Second and third deviatoric stress invariant

X  = Distance from the left side of the first column Y= Y-axis of Cartesian coordinate y= Y-axis in complex space Z= Direction of increasing in depth in Goodman and Brown’s model z= Complex number θ = Angle that measured from the bisector of the wedge θ = Lode angle  = half of the wedge angle α = Slope angle in Goodman and Brown’s model  = Angle of polar coordinate in Goodman and Brown’s model  = Airy stress function

of

f  r  = Loading function in real space g = Gravity h = height of the slope h = Subinterval of Fillon-Simpson numerical integration method i = imaginary number i = Counter in FDM I 1 , I 2 , I 3 = Fisrt, second and third stress invariant

V= Vertical displacement X= X-axis of Cartesian coordinate x= X-axis in complex space

In civil and mining projects, engineers encounter serious problems associated with rock and soil slopes. Instability of these slopes can cause significant economic, safety and social impacts as well as physical damage. In order to mitigate these losses, understanding the type and mechanism of slope instability is of great importance. One of the common types of slope instability is toppling of rock blocks. In 1976, Goodman and Bray presented a comprehensive category of all types of toppling failures for the first time, as presented in Fig. 1, (Goodman and Bray 1976). According to this classification, toppling failure is divided into four general

Journal Pre-proof categories, which are flexural, blocky, block-flexure, and secondary modes. Secondary toppling failure mode has more complicated failure mechanisms compared with the three main modes, and it has so far been investigated through numerical, experimental, and analytical methods (Tsesarsky and Hatzor 2009; Frayssines and Hantz 2009; Alejano et al. 2010; Amini et al. 2012; Alejano et al. 2015; Amini et al. 2017; Amini et al. 2018). Goodman and Bray examined the cause of slope instability when several vertical rock columns overlaid a weak horizontal layer of

of

soil or rock (category 4 in Fig. 1), (Goodman and Bray 1976). In this circumstance, either the

ro

shear failure of the weak horizontal soil/rock layer or differential settlement of vertical rock

-p

columns within the continuous horizontal soil/rock layer or a combination of the two may lead to a toppling failure of vertical long rock columns. This mode of failure is called “toppling-

re

slumping” failure (Fig. 2). In 1981, Evans provided a simple model for determining the

lP

settlements on either side of a rock column that overlies a weak material and showed that the

na

rock columns could topple due to the differential settlements (Evans 1981). After this study, toppling in such structures have been reported in several formations such as in combinations of sandstone-shale,

sandstone-colluvium,

ur

shale-basalt,

and

limestone-sandstone

(Evans

1981;

Jo

Pasuto and Soldati 2013; Borgatti et al. 2015; Spreafico et al. 2016; Spreafico et al. 2017; Wyllie 2018). Despite the abundance of these formations in the nature, the mechanism of this instability has not been fully investigated. One of the issues examined in these structures is lateral spreading. In this case, rock columns are subjected to toppling and rock fall, especially when undermining occurs in the weak lower rock/soil layers (Borgatti et al. 2015; Spreafico et al. 2016; Spreafico et al. 2017). Spreafico et al. (2016) utilized a back-analysis method and remote sensing technique to investigate instability of slope in which the loose and weak layer of soil were overlaid by long rock columns. They argued that toppling failure was the main failure

Journal Pre-proof mechanism in the slope’s instability, and undermining at the base of the rock columns is the key factor in this type of failure. To be more precise, toppling-slumping failure is dominant when a weak horizontal layer of rock or soil overlaid by several vertical stronger rock columns. In this case, the weaker layer settles or slides due to the above rock column’s weight. Because of differential settlement or shear failure

of

in the weak layer, the rock column can be toppled. The mechanism of this settlement or shear failure is strongly related to stress distribution in the weak rock/soil layer. Unlike the horizontal

ro

ground surface or half-plane, stress distribution under rock columns near the slope crest is

-p

asymmetrical. Therefore, the stress distributed inside the weak layer under the rock columns and

re

the two sides of the rock columns is not equal. The non-equality of stress on both sides of a

lP

column causes a non-uniform settlement and, consequently, leads to a toppling failure which is known as toppling-slumping. Stress distribution due to external loading on the top of the half-

na

plane is well-known in geotechnical engineering. However, this is not clearly understood when

ur

an external load applied on the above surface of a slope. The present paper aims to introduce an analytic approach for toppling-slumping slope failure. In addition, an attempt will be made to

Jo

achieve a better understanding of the failure mechanism and the factors contributing to the toppling-slumping slope failure. To this aim, the material forming slope was supposed as homogenous material with linear elastic behavior. Furthermore, the rock cliffs are assumed as vertical columns which apply uniform pressure into the slope. The effect of water table and slope saturation are note considered in this work.

Journal Pre-proof To make this investigation, first stress distribution in a slope caused by the upper surface load will be examined by means of the theory of elasticity. Second the corresponding displacement

Jo

ur

na

lP

re

-p

ro

of

and shear failure will be discussed.

Fig. 1 Goodman-Bray’s classification of toppling failures (after Goodman and Bray, 1976).

of

Journal Pre-proof

-p

ro

Fig. 2 Schematic view of toppling-slumping failure in a rock slope.

re

2. Analytical method

lP

2.1. Mellin transformation

na

One of the most widely used methods in the theory of elasticity is the use of the Airy stress function to solve the plane strain/stress problem. In the Airy stress function, the general

ur

formulation is reduced to a governing equation by an unknown variable and the governing

Jo

equation is solved by mathematical methods (Sadd 2009). Equation (1) is called the Bi-harmonic equation in the Polar coordinate system that can be used for the calculation of stress distribution. This equation satisfies the equilibrium and compatibility relations. In Equation (1),  represents the Airy stress function: 2

 2 1  1 2    2   0 r r r 2  2   r

(1)

By solving Equation (1) and achieving the Airy stress function, the stress components can be calculated using Equations (2):

Journal Pre-proof  

 2 r 2

r 

1  1  2  r r r 2  2

 r  

  1     r  r  

(2)

Consider the wedge of Fig. 3, which is loaded at its edges. The values of r and θ can be changed as follows: 0r 

ur

na

lP

re

-p

ro

of

    

Fig. 3 Symmetric pressure applied to the edge of the wedge (after Tranter, 1951).

Jo

Using the Mellin Transform, Tranter proposed an Airy stress function in the transformed space and calculated the stress distribution in a situation where the angle  is 90 degrees, that is the wedge becomes half-plane (Tranter 1951). Equations (3) and (4) represent the transformation and the inversed transform of the Airy stress function: 

Φ  z      r  r z 1dr 0

(3)

Journal Pre-proof 

1

2 i 

(4)

c i 

r  z Φ  z  dz

c i 

In this equation, Φ  z  denotes the transformed Airy stress function of  in the complex space. Equation (5) represents the Tranter proposed stress function for solving the wedge problem:

Φ  A  sin  z    B  cos  z    C  sin   z  2    D  cos   z  2  

(5)

In Equation (5), the parameters A to D depend on the values of z and  and they should be

of

obtained according to the boundary conditions of the problem. It should be noted that z

ro

represents the complex number. Tranter showed that stress can be calculated by finding the stress

2 i

1



re

z  z  1  r  z 2dz

c i 

c i 

[

c i 

2 i 

c i 

d 2  z  ] r  z 2dz 2 d

 z  1 c i 

lP

1

c i 

(6)

d   z 2 r dz d

Jo

 r 

2 i 

na

r 

1

ur

 

-p

function of Equations (5) (Tranter 1951):

Tranter solved the problem for a situation where the angle 2α is 180 degrees (half-plane) and a symmetric pressure is applied to the edge of the wedge, as follows (Tranter 1951):

 

 2   p ,    , 0  r  a r 2

(7)

In Equations (7), P represents the amount of uniform load applied on the edge of the wedge. The wedge in Fig. 3 can be taken as Fig. 4. This figure could represent a slope loaded at its surface. If

Journal Pre-proof this loading is replaced by the weight of several vertical rock columns, located on the surface of

-p

ro

of

the slope, the model has the potential to replicate the toppling-slumping failure mode.

re

Fig. 4 Schematic representation of a slope loaded at its surface.

 

 2  0 ,    r 2

(8)

Jo

 r  

na

 2  f r  ,    r 2

ur

 

lP

The boundary conditions for such a model can be expressed as Equations (8):

  1      0 ,    r  r  

It should be noted that any form of load such as linear, trapezoidal, parabolic etc. can be implemented in the equation (8) by adjusting load function. Tranter pointed out that when r tends to infinity, it is assumed that the stress function is such that all of the values in Equations (9) tend to zero (Tranter, 1951):

Journal Pre-proof r z n

 n , n  0,1, 2,3 r n

rz

 n , n  1, 2  n

r z 1

 3 r  2

(9)

If Equation (2) is multiplied by r z 1 and integrated in terms of r in the range of 0 to infinity, the Equation (10) can be obtained:

0

r z 1

  2 dr   r z 1f  r  dr 2 o r

(10)

of





(11)

-p

    |   z  1 [r z  |  z  r z 1 dr ] 0 0 r 0

re

r z 1

ro

With twice integration by parts on the left side of Equation (10):

As mentioned, when r tends to infinity, it is assumed that the stress function is such that all of the

lP

values in Equation (9) tend to zero. Hence non-integral terms in Equation (11) are equal to zero

na

and, with the Transform definition in Equation (3), Equation (11) is simplified into Equation (12):

ur

    z  1  r z   z  r z 1 dr   z  z  1 Φ  z    0 r

Therefore:

(12)

Jo

r z 1



F  z    r z 1f  r  dr  z  z  1 Φ  z    0

(13)

Where, F  z  is equal to the transformation of the applied load f  r  . By performing similar operations on the other part of Equation (2), we have:

Journal Pre-proof

 z  1

dΦ  z  0, θ α d

z  z  1 Φ  z   0 , θ  α

 z  1

dΦ  z   0 , θ  α d

(14)

By setting the values of Equations (13) and (14) in Equation (5), parameters A to D are defined as presented in the Appendix A. Therefore, the parameters of the Airy stress function are obtained in the transformed state for the boundary conditions shown in Fig. 4 and Equation (8). Now, using Equations (6), the stress can

of

be calculated. By setting the stress function in Equations (6), Equations (15) will be obtained:

a 

z 2

 a   cos  z   cos   z  2     cos   z  2    cos  z    sin  z   sin   z  2     sin   z  2    sin  z    dz      c i  2 i r    z  1 sin  2   sin  2  z  1    z  1 sin  2   sin  2  z  1    (15) 1

c i 

z 2

ur

na

 r 

lP

re

-p

c i 

ro

 z cos  z   sin   z  2     cos   z  2    sin  z   sin   z  2    cos  z    z sin  z   cos   z  2     z 2  z 2         dz c i  4 i r    z  1 sin  2   sin  2  z  1    z  1 sin  2   sin  2  z  1       cos   z  2    sin  z    z  4 cos  z   sin   z  2    z  4 sin  z   cos   z  2     sin   z  2    cos  z    z 2  1 c  i   a   z 2 r   z 2      dz c i  4 i r    z  1 sin  2   sin  2  z  1    z  1 sin  2   sin  2  z  1      1

An important element in complex analysis, such as calculation of integrals of Equation (15), is

Jo

the points at which the function will be singular and are referred to as poles. Therefore, first, the poles of the function must be determined. According to the functions, only the point z = -1 is the pole of these meromorphic functions. However, it must be noted that this point makes the numerators of   and  r zero as well. In other words, this point for these two functions is zeroorder or removable singularity point that the discussed singularity will be solved by eliminating the ambiguity. However, the point z = -1 for  r is a first or simple order pole. Thus, the integrals of Equation (15) can be calculated by linear integration in the z = -1 path from -∞ to + ∞ and

Journal Pre-proof subtracting the residue of the function multiplied by πi. Therefore, the path of integration can be represented as Equations B.1 and B.2 in Appendix B. By performing mathematical operations, the stresses are obtained as Equation (16).

(16-

a)

re

   a   sin  yLog  r    a   dy     r    [  0 Ty  R      r 1 y 2    

-p

ro

of

   sin     cosh     y  sin     cosh     y          y sin 2  sinh  2 y         a    sin  yLog     dy   a  0   sin     cosh     y  sin     cosh     y      r         r         r    y sin 2   sinh 2  y             sin cos   sin cos      sin 2  2  sin 2  2  

   a   cos  yLog  r          Ry T   dy 0   1 y 2    

(16-b)

na



lP



Jo

ur

     a   a  sin  yLog     cos  yLog          sin Cos   sin cos   r   r      Ky  L    dy    Ly  K   dy   ] 2 2 0 0     1 y 1 y sin 2  2 2  sin 2        

 r

 sin     sinh     y  sin     sinh     y  y sin  2   sinh  2 y  a    2r   0  sin     sinh     y  sin     sinh     y   y sin  2   sinh  2 y  

   cos  yLog  a   dy      r     

In Equation (16-b), the values of R, T, K, and L are defined in Appendix C. 2.2. Computing the integrals

(16-c)

Journal Pre-proof The integrals of Equation (16) do not have an exact solution and therefore, they should be solved numerically. There are many numerical methods to solve the complex equations (Rabczuk and Belytschko 2004; Rabczuk et al. 2010; Ren et al. 2016). A suitable numerical solution method for such integrals is the Filon method (Filon 1930). In the Filon method, the integral interval is divided into 2n parts that the length of each is equal to h. In Equation (17), the Filon numerical

f  x  cos  kx  dx  h[Ω f x 2 n sin  kx 2 n   f x 0 sin  kx 0   C 2 n  C 2 n 1 ]

(17)

re

x0

f  x  sin  kx  dx  h[Ω f x 0 cos  kx 0   f x 2 n cos  kx 2 n   S 2 n  S 2 n 1 ]

ro



x 2n

x0

-p



x 2n

of

integration method is shown:

lP

All the parameters are defined in Appendix D.

As is clear, in the numerical solution, the integral boundary should be finite because the

na

replacing series cannot be summed up to a point at infinity. Therefore, the upper boundary of the

ur

integral should be considered as a finite value. Here, the upper boundary is considered as a point after which the integral value is equal to a small value. Since these integrals represent the stress

Jo

associated with the distribution of external loads in the slope, the value of the upper boundary is considered to be such that the integral after that point is maximally equal to 0.01. This is stated in Equations (18): 

U



0

0

U

 f  y  dy   f  y  dy   f  y dy 

 f  y  dy  0.01 U

(18-a)

(18-b)

Journal Pre-proof In the integral shown in Equation (18-b), the asymptotic analysis of function f (y) can be used. In general, the integrant functions of the Equation (16) with asymptotic point of view can be written as a linear combination of me dy in which m and d, are constant values. 2.3. Discrete loading The loading on the slope can be due to rock columns that are spaced apart with certain distance.

na

lP

re

-p

ro

of

To achieve such a case, consider the schematic loading as shown in Fig. 5.

ur

Fig. 5 Schematic representation of discrete loading on a slope.

Jo

To model the loading shown in Fig. 5, it is necessary to introduce a new loading function. For such a loading, the step function can be used. The definition of the step function is given in Equation (19):

0 u  r  ri    1

r  ri r  ri

(19)

Suppose that the spacing between the three rock columns is equal to the width of the columns in Fig. 5 and each column applies a unit load to the slope. In this case, this type of loading could be simulated as several step functions as shown in Fig. E.1 in Appendix E and Equation (20).

Journal Pre-proof f  r   u  r   u  r  1  u  r  2  u  r  3  u  r  4   u  r  5

(20)

Equation (20) can be written as series as in Equation (21): 5

f  r     1 u  r  ri  i

(21)

-p

    i F  z    r z 1f  r  dr     1 u  r  ri   r z 1dr 0 0  i 0 

ro

With respect to Equation (13), F  z  is written as follows:

of

i 0

re

Assuming that the load is applied from the point r=0 to the point r = A at the edge of the wedge

lP

(   ) , and the compressive load is considered as negative:

(22)

ur

na

A  A  1 A i i F  z      1 u  r  ri  r z 1dr   1 ri z 2  0 z  2 i 0  i 0 

intervals.

Jo

Therefore, using this loading function, any kind of uniform load can be applied with arbitrary

2.4. Gravitational loading

In Equation (1), which expresses the Bi-harmonic equation, it is assumed that the body force is zero. This assumption is expected in elasticity (Sadd 2009). However, in slope stability analysis, the stress distribution caused by gravitational loading is of great importance. It is clear that for a shallow slope, gravitational loading has little effect on the stress distribution within the slope. However, for a high slope, stress distribution caused by gravitational loading is significant.

Journal Pre-proof Goodman and Brown presented the Equation (23) based on the theory of elasticity for the gravitational stress distribution in the slope (Goodman and Brown 1963). An important feature in theory of elasticity is the establishment of the principle of superposition. Therefore, due to the importance of gravitational loading in high slopes, the stress distribution caused by loading in the proposed analytical model can be summed up by the stress distribution of the gravitational load proposed by Goodman and Brown (1963). The gravitational stresses provided by Goodman and

 sin α     [Y α  sin α  cos α    2sin α  cos 2 α   X sin α  Y cos α   log    α  tan α   sin α  

g

ro

 xx 

of

Brown are shown in Equations (23).

 sin α     [Y sin α  cos α   2sin 3 α   X sin α  Y cos α   log   α  tan    sin α   (23)

g

re

 yy   gy 

-p

 sin α  Y sin α  1  2 cos 2 α    X cos   1  2sin 2 α  

lP

  sin α  Y sin α  1  2 cos 2 α    X cos α  1  2sin 2 α   

 sin α      g [Y sin 2 α   2 cos α  sin 2 α   X sin α  Y cos α   log   α  tan α   sin α  

na

 xy 

Y X

Jo

  tan 1

ur

  sin α  1  2sin 2 α    X sin α  Y cos α  

3. Evaluation of the analytical model by numerical method Based on the proposed analytical method in section 3, a Mathematica package was developed to evaluate

the

stress

distribution

within

the

slope

(the

package

is

available

on

www.inscribe.ir/slope-stability) and its results were compared with Abaqus finite element software. Fig. 6 shows the mesh and size of the model. It should be noted that 10800 triangular elements were used in this model.

of

Journal Pre-proof

ro

Fig. 6 The size and mesh of the numerical modeling.

-p

3.1. Single rock column loading case

re

To evaluate the validity of the analytical model results, three slopes at 30, 60 and 90 degree slope

lP

angles are considered. Suppose a rock column, as in Fig. 7, located at the top of the slope which

Jo

ur

na

applies a uniform loading, p, to the slope.

Fig. 7 A single rock column located at the top of the slope.

Stress distribution caused by a single rock column in different points of the slope is shown in Figs. 8 to 10. The results of the analytical solution were compared with the finite element numerical model. It should be noted that, in these figures, Y represents the depth of the slope.

of

Journal Pre-proof

(a)

Jo

ur

na

lP

re

-p

ro

(b)

(c)

Fig.8 Stress distribution in a 90 degree slope due to a single rock column loading at Y = 0.5 m: (a) vertical stress; (b) horizontal stress; (c) shear stress.

(a)

Journal Pre-proof

-p

(c)

ro

of

(b)

Jo

ur

na

lP

re

Fig. 9 Stress distribution in a 60 degree slope due to a single rock column loading at Y = 0.5 m: (a) vertical stress; (b) horizontal stress; (c) shear stress.

(a)

(b)

of

Journal Pre-proof

ro

(c)

-p

Fig. 10 Stress distribution in a 30 degree slope due to a single rock column loading at Y = 0.5 m: (a) vertical stress; (b) horizontal stress; (c) shear stress.

re

Accordingly, the analytical model has a high accuracy in predicting the stress distribution in the

lP

slopes. The highest relative error is about 4%. However, by increasing the model dimension in the numerical model and decreasing the mesh size simultaneously, the numerical model solution

na

tends to the analytical model and the relative error will decrease. As stated previously,

ur

gravitational loading has not been considered in this model, so stress distribution is solely due to the presence of external force at the top of the slope. Fig. 11 shows the stress distribution in the

Jo

slope up to a depth of 6 m below the surface. From about a depth of 3 m, the analytic model estimates lower values for vertical stress than the numerical model and the relative error becomes more than 10%. As the depth increases, the error also increases. Fig. 12 presents the vertical stress distribution in the model after summing up the proposed analytical method with the Goodman and Brown method. The same picture compares the analytical method with the numerical model. It is evident that in higher depths, gravitational stress distribution should be taken into consideration using the Goodman and Brown method. As it can be seen, by adding

Journal Pre-proof gravitational loads to the results of proposed method good agreement can be achieved between

-p

ro

of

outcomes of numerical and analytical methods.

(b)

Jo

ur

na

lP

re

(a)

(c) Fig. 11 Stress distribution in a 10 degree slope due to single rock column loading at X/x=0.1 (a) vertical stress (b) horizontal stress (c) shear stress.

of

Journal Pre-proof

(b)

ur

na

lP

re

-p

ro

(a)

(c)

Jo

Fig. 12 Stress distribution in a 10 degree slope by considering the gravitational loading at X/x=0.1 (a) vertical stress (b) horizontal stress (c) shear stress.

3.2. Discrete loading case In Section 2.3, the discrete loading method was discussed. To validate the results of the proposed method, a simple numerical model was developed using Abaqus software based on the slope case presented in Fig. 5 and compared with the analytical model. The results of analytical and numerical models are shown in Fig. 13. As can be seen, the analytical model predicts the stress distribution with a very high accuracy.

of

Journal Pre-proof

(a)

Jo

ur

na

lP

re

-p

ro

(b)

(c)

Fig. 13 Stress distribution in a 60 degree slope due to discrete loading condition at Y = 0.5 m: (a) vertical stress; (b) horizontal stress; (c) shear stress.

4. Evaluation of toppling-slumping and shear failure In the previous sections assessing the stress distribution in the slope was discussed. To evaluate the toppling-slumping and shear failure, two distinct criteria were used. Displacement criterion was used for evaluating the toppling, and stress-based criterion was used for assessing the shear failure. These two criteria will be discussed in the next section.

Journal Pre-proof 4.1 Toppling-slumping failure criteria According to Figs. 8 to 13, in both single and discrete loading cases, the stresses on the two sides of the rock column(s) are not equal. Hence, it can lead to differential settlement on two sides of a rock column and cause toppling. To determine the probability of this type of failure, the displacement must be calculated in the two column ends. Based on the theory of elasticity, when

of

the stress distribution function is known, the displacement field caused by the stress distribution can be determined through the relationship between strain and displacement. Since Equations 16

ro

have no exact solution, they must be solved numerically. Therefore, stress and strain

-p

distributions cannot be cited as a comprehensive function of the slope coordinates so that the

re

displacement distribution could be obtainable by integrating the strain equations. Thus, here the

lP

Finite Difference Method (FDM) is used to determine the displacement at both ends of the rock column. This is shown schematically in Fig. E.1 in Appendix E. In this figure, V represents the

na

vertical displacement. In each of the grid points, the stress is specified (by Equations (16)) and so

1   ij  kk ij  E

(24)

Jo

e ij 

ur

the strain value can be calculated through the Equation 24.

In Equation (24), E represents the Young’s modulus and ν denotes the Poisson’s ratio. By determining the strain at any grid point and using Equation (25) (Jovanović and Süli 2013), it is possible to calculate the displacement in each of the grid points; therefore, the vertical displacement difference ( Δ ) in the both ends of the rock column can be computed by Equation (26):

e yy 

V V  y  ih  V  y 0   O h  y ih

(25)

Journal Pre-proof  V 1 V 2

(26)

Based on Equation (25), the truncation error of the displacement is the order of h. By reducing the value of h, the truncation error will decrease. The value of h depends on the number of grids (n). The higher number of grids lowers the value of h. In order to reduce the computational time and to determine the appropriate solution, the value of h is somehow considered, after which the

of

error of calculating the displacement in two successive values of h would be less than 1%. At this

ro

value, as shown in Fig. 14, the obtained answer by Abaqus becomes equal with that of Equation

Jo

ur

na

lP

re

-p

(25).

Fig. 14 Comparing the displacement calculated from Equation (25) at different grids (n) with FEM (Abaqus) answer.

It is known that toppling occurs when Equation (27) holds (Wyllie 2018):

x  tan p y

(27)

Journal Pre-proof In Equation (27) x and y are width and height of the rock column respectively and  p denotes the dip angle of the slope. This equation can be written in terms of vertical displacement difference () as follow. According to Fig. 15, it is clear that:   x 

 p  sin 1 

of

Therefore, the toppling condition based on the displacement difference will be (Fig. 15):

(28)

Jo

ur

na

lP

re

-p

ro

  x    x sin  tan 1     y  

Fig. 15 Schematic representation of differential settlement and definition of toppling failure based on the displacement.

4.2 Shear failure criteria The probability of toppling failure can be determined by above-mentioned procedure. The question that can be posed here is whether the slope slides before or at the same time as toppling. To investigate the safety of the slope against shear sliding, the Mohr-Coulomb behavioral model was considered first for the slope. The effective parameters of the slope stability against sliding in this behavioral model are cohesion, friction angle, and stress state in the slope. Equation (29)

Journal Pre-proof represents the behavioral model of Mohr-Coulomb based on stress tensor invariants (de Souza Neto et al. 2011). It should be noted that in this equation the compression considered as positive.

I1 sin θ sin    sin   c cos   J 2  cosθ  0 3 3  

(29)

Where,

ro

1  ii  jj   ij  ij  2

-p

I2 

of

I 1   ii

1 2  I 1  3I 2  3

J3 

1  2I 13  9I 1I 2  27I 3  27

(29-b)

(29-c)

(29-d)

(29-e)

Jo

ur

J2 

na

lP

re

I 3  det  ij 

(29-a)

 3 3J 3  1 θ  sin 1  3/2   3  2J 2 

(29-f)

In the above equations I 1 , I 2 , and I 3 are the first, second and third stress invariants, J 2 and J 3 are the second, and third deviatoric stress invariants and θ is the Lode angle. Due to plain strain condition, the stress perpendicular to the model’s surface is considered as   1   3  , where,  1 and  3 are the maximum and minimum principal stresses. By sorting the Equation (29) based on

J 2 , the Equation (30) will be obtained:

Journal Pre-proof I1 sin   c cos  3 J2  sin θ sin  cosθ  3

(30)

By determining the stresses at any points of the slope, and according to Equation (30), the

J2

value will be obtained at each point based on the strength characteristic of the slope. According

J 2 can be also calculated from the values of stress at each point, which

of

to Equation (29-d),

ro

indicates the stress state in the slope. In order to verify the slope stability against shear sliding, the safety factor was calculated along predefined sliding surfaces. The method is that, first a grid

-p

of 20  20 was considered at the side of the slope. From each of the points of this grid, circles

re

were drawn with different radius and the safety factor along these curves was calculated using

 J  J 



2

Strength

Stress

(31)

ur

2

na

FS 

lP

Equation (31):

Jo

In this equation, the strength and stress values were calculated using the Equation (30) and Equation (29-d), respectively. It should be pointed out that in this section we just used the elastic stress which was introduced in section 3. In fact, stress redistribution will occur after one point achieves failure state. Ideally, elastic-plastic analysis should be used in order to consider the stress redistribution effect. However, for safety factor calculation, it is reasonable to use the elastic stress for a wide range of slope conditions encountered in engineering practice (Stianson et al. 2004). It should be noted that the aforementioned Mathematica package can be used to find the safety factor of the slope.

Journal Pre-proof 5. A typical example In order to investigate the toppling-slumping or shear failure in a slope by proposed analytical model, a typical example is studied. The configuration of loading is shown in Fig. 16. In this model, 5 rock columns are spaced at 1% of the width of each column. It should be noted that as shown in Fig. 16, X  was measured from the left side of the first column along the slope

ur

na

lP

re

-p

ro

of

crest in all models.

Fig. 16 Schematic representation of a set of rock column on the slope.

Jo

The load applied by these rock columns to the slope depends on the dimensions of the columns. The height of these columns can be very high (for example, more than 100 m (Evans 1981; Spreafico et al. 2016)). The width of these columns will be determined by the spacing between the vertical joints. Spreafico et al. (2016) have stated that the spacing of the joints in the studied structure is very low and the average spacing between the three measured joints is 0.5, 0.7 and 1.9 m. According to these case studies, the height and width of the rock columns for the typical example was considered as 100 m and 1 m, respectively. The slope is assumed as being completely dry and the effect of water on this model is not considered. These rock columns are

Journal Pre-proof considered as a row and the distance between the first column in this row and the slope crest (a) is variable in different models. In addition, the slope angle (ψ) is also changed in these models. It should be noted that the negative values of a/x represent the undermining at the base of first column. The stress values in the slope are calculated by Equation (16), and then, using Equations (24) and (25), displacements are calculated in the two ends of each rock column. Clearly, the strain and displacement will depend on the Poisson’s ratio and the Young’s modulus.

of

In the examined models, the Poisson’s ratio is considered as 0.3. As mentioned before, toppling

ro

occurs due to differential settlement in the loose layer. This layer can be considered as soil or

-p

weathered rock (Evans 1981). Therefore, the Young’s modulus of the loose layer should have a

re

low value. Accordingly, the value of Young’s modulus varied from 20 MPa to 10 GPa.

lP

In Fig. 17, the displacements in both ends of column () in the slope is plotted for a Young’s

Jo

ur

na

modulus of 500 MPa and with different slope angles and a/x values.

(a)

(b)

(d)

ro

of

Journal Pre-proof

Jo

ur

na

lP

re

-p

(c)

(e)

(f)

of

Journal Pre-proof

(h)

Jo

ur

na

lP

re

-p

ro

(g)

(i)

(j)

Fig. 17 Difference of displacement in both ends of column at (a) negative and zero a/x in 90 degree slope; (b) positive a/x in 90 degree slope; (c) negative and zero a/x in 80 degree slope; (d) positive a/x in 80 degree slope; (e) negative and zero a/x in 75 degree slope; (f) positive a/x and 75 degree slope; (g) negative and zero a/x in 70 degree slope (h)positive a/x in 70 degree slope; (i) negative and zero a/x in 60 degree slope; (j) positive a/x in 60 degree slope.

According to Fig. 17-g to 17-j, for the 70 degree and 60 degree slopes, the displacement difference on both ends of the first column will be negative. This indicates a more vertical movement on the right side of the first column. In Fig. 18-a, the vertical stress distribution under the first column is shown for 80, 70 and 60 degrees slopes with a/x = -0.5. As can be seen, unlike the 80 degree slope, the vertical stress induced on the left side of the rock column will be less

Journal Pre-proof than the right in 70 degree and 60 degree slopes. Hence, at these slope angles, more vertical stress on the right side of the first column caused the more vertical displacement. As can be seen in Figs. 17-g and 17-i by decreasing the amount of undermining (i.e., increasing the value of a/x) at these angles, the absolute value of the displacement difference in the first column will increase. The reason is illustrated in Fig. 18-b. By decreasing the amount of

of

undermining, the vertical stress distribution will be more asymmetric under the first column and the vertical stress on the right side of the rock column will be greater than its left side. Therefore,

ro

vertical stress distribution will play a key role in the displacement of the rock column and thus

-p

increases the possibility of toppling failure.

re

In all models, at a/x with values of -0.5 and -0.4, the displacement difference of the second

lP

column is greater than the third one. By increasing the value of a/x, this difference decreases and at positive values of a/x, the displacement difference of the second column will be less than the

na

third column. This is because at a/x values of -0.5 and -0.4, the induced vertical stress under the

ur

second column is affected by the stress distribution caused by the first column. By increasing the

Jo

a/x, the induced stress under the rock columns is such that the row of columns can be considered as a single column. Fig. 18-c is presented to clarify this issue. In Fig. 18-c, two vertical lines indicate the location of second column. Accordingly, by increasing the value of a/x, the stress distribution becomes more symmetrical under the row of rock columns. At a higher a/x, which represents greater distance from the edge of the slope, the stress distribution is completely symmetric indicating that the rows of columns act like a single column on a horizontal ground surface. In fact, by increasing a/x, the effect of the slope geometry is removed and behavior similar to the stress distribution of a column on a

Journal Pre-proof half-plane will occur. By decreasing the slope angle, this behavior (the change in the stress distribution under the columns, as in the stress distribution under a single column on a half plane) will occur at a lower a/x. Fig. 18-d shows the stress distributions under the rows of rock columns at different angles and with a/x equal to 0. As the slope angle decreases, the stress distribution becomes more symmetric and the column set tends to behave like a single column on

of

a half-plane. In general, the displacement difference in column 1 strongly depends on the slope angle and by

ro

moving to the end columns, the stress distribution under the columns will not be a function of the

-p

slope angle. This is shown in Fig. 18-e. This implies that at farther distance from the slope crest,

Jo

ur

na

lP

re

the effect of geometry will disappear.

(a)

(c)

(b)

Journal Pre-proof

ro

of

(d)

re

-p

(e) Fig. 18 (a)Vertical stress distribution under the first column in the 90, 80 and 70 degree slope with the a/x=-0.5 at the Y=0.5m, (b)Vertical stress distribution under the first column in the 70 degree slope with the a/x=-0.5, -0.3, and -0.1 at the Y = 0.5 m, (c)Vertical stress distribution under the columns in the 90 degree slope with the a/x = -0.5, -0.3, 0 and +3 at the Y = 0.5 m, (d)Vertical stress distribution under the columns at the 90, 80, 75, 70 and 60 degree slope with the a/x=0 at the Y = 0.5 m, (e)Change in displacement difference with slope angle for rock columns.

lP

As shown in Fig. 17-a, all columns will tilt to the left, as the amount of displacement difference on both ends of all columns is positive, and thus all the columns tend to forward toppling at 90

na

degrees for negative and zero a/x values. This is also true for slopes with an 80 degree angle for

ur

a/x = -0.5 to -0.2 and a 75 degree angle for a/x = -0.5. For the rest of the slope cases, the columns near the slope crest tend to backward toppling and the end columns tend to forward toppling. Fig.

Jo

19 represents the definition of forward and backward toppling. As stated above, by increasing the a/x, the vertical stress distribution induced below the row of rock columns tend to be symmetric, and, consequently, at a certain distance from the crest the displacement difference of the columns 1 and 5 become equal in magnitude but in the opposite directions (e.g., curve that represent the a/x = +1 in Fig. 17-j). In the case of Fig. 17-j, column 3 will have a displacement difference equal to zero. Thus, when the vertical stress distribution is symmetric, the set of columns are in an equilibrium state and toppling will not occur.

Journal Pre-proof

(b)

(a)

of

Fig 19. Definition of the: (a) forward toppling (b) backward toppling.

ro

Therefore, the toppling of the row of columns depends on the displacement difference of the first

-p

and last columns. Here are two possible scenarios. In the first case, the first and the last columns are both tilted to the left, if the value of their displacement difference satisfied Equation (28),

re

toppling will occur, otherwise if the displacement difference did not satisfy Equation (28), the

lP

toppling will not happen. In the second case, the first column tilts to the right and the last column to the left. In this case, two different behaviors can be observed. If the displacement of the first

na

and the last columns is equal in size but opposite in direction, the column set is in equilibrium

ur

condition and no toppling occurs (even if the displacement difference in each of the columns

Jo

satisfied Equation (28)). The second behavior occurs when the displacement in the first column according to Equation (28) is insufficient for toppling but the displacement of the last column is sufficient. In this case, the toppling of the last column can cause the instability of the column set. As noted before, the displacement in a slope is related to the Young’s modulus and the Poisson’s ratio of the weak rock/soil layer in the lower part of the slope. In order to reduce the computational time, the Poisson’s ratio is considered as 0.3 for all the models and the Young’s modulus is assumed to be variable. Fig. 20 shows the displacement difference graph for last column, plotted at different a/x values, for a range of Young’s modulus from 20MPa to 10 GPa.

Journal Pre-proof The dotted line on the graph represents the amount of displacement difference needed to topple a column with a height to width of 100:1, calculated by Equation (28). The curves on the right side

Jo

ur

na

lP

re

-p

ro

of

of the dotted line represent situations where toppling occurs.

(a)

Jo

ur

na

lP

re

(b)

-p

ro

of

Journal Pre-proof

(c)

Jo

ur

na

lP

re

(d)

-p

ro

of

Journal Pre-proof

(e) Fig. 20 Displacement difference graph for last column at different a/x values and Young’s modulus (a) 90 degree slope (b) 80 degree slope (c) 75 degree slope (d) 70 degree slope (e) 60 degree slope.

According to Fig. 20, toppling will occur in a slope with the Young’s modulus of less than 500 MPa at all angles and range of the a/x parameter.

Journal Pre-proof Considering the graph of the Young’s modulus of 500 MPa, Fig. 20 indicates that for slope angles of 90 and 80 degree, toppling will not occur at a/x of greater than 1.5 (the curve for Young’s modulus of 500 MPa meets the dotted line in a/x that is below the 1.5). In addition, toppling will not occur for slope angles of 75 and 70 degrees with a/x values of more than 1. This value is equal to 0.5 for a 60 degree slope. Therefore, for the examined slope conditions, the Young’s modulus of 500 MPa is considered critical. Generally, rocks have Young’s modulus

of

higher than 500 MPa, but completely and highly weathered rocks can have a Young’s modulus

ro

of less than 500 MPa (Erguler 2009; Heidari et al. 2013). Soil are also highly vulnerable to this

-p

form of toppling, in terms of differential settlement, due to their low Young’s modulus. It should

re

be noted that in the nature most of the rock columns are sub-vertical instead of vertical. In these cases, based on the inclination of individual rock columns (dipping into or out of the slope face),

lP

the rock columns tend to forward or backward toppling. This inclination can be taken into

na

account by developing the proposed solution through considering the non-uniform loading.

ur

5.2. Shear sliding of the slope

Jo

It was found that with the Young’s modulus, the slope angle and a/x ratio are effective parameters on the toppling of the rock columns, and in a critical composition, the columns located on the top of the slope are toppled. As shown in Section 6.1, the stress distribution in the slope is a function of the slope angle and the a/x ratio. Therefore, in order to investigate the slope stability against sliding, the slopes with different cohesion, friction angle, a/x ratio and angles were considered. The specifications of the developed analytical models are shown in Table 1. In Table 1, the number between the colons is the interval step. It should be noted that the sensitivity and uncertainty analysis were not considered in this research. However, these analyses can be

Journal Pre-proof taken into account based on approaches developed by Vu-Bac et al. using MATLAB toolbox (Vu-Bac et al. 2016). Table 1 The specifications of the developed analytical models

Slope angle ( )

a/x

c (MPa)



90-80-75-70-60

-0.5: 0.5 :2

0.1: 0.1 :1

20: 5 :35

J 2 values based on strength and stress for different slopes. In

of

Fig. 21 shows the distribution of

Jo

ur

na

lP

re

-p

(29-d), respectively at a/x=0, c=0.2 MPa and   30 .

ro

this figure, the strength and stress values are calculated using the Equation (30) and Equation

(a)

(b)

of

Journal Pre-proof

(d)

(e) (f)

Jo

ur

na

lP

re

-p

ro

(c)

(g)

(h)

of

Journal Pre-proof

Fig. 21 Contour of

(j)

ro

(i)

J 2 based on the (a) stress and (b) strength at the 90 degree slope angle; (c) stress and (d) strength at the 80

-p

degree slope angle; (e) stress and (f) strength at the 75 degree slope angle; (g) stress and (h) strength at the 70 degree slope angle (i) stress and (j) strength at the 60 degree slope angle.

re

As can be seen, by decreasing the slope angle the induced stress in the slope decreased. Also, it

lP

is clear that by increasing the slope angle, the warm color contour represents the higher induced stress, moves from the beneath of the rock column to the slope face. Fig. 22 show the safety

na

factor in the slopes with different a/x ratios in which there are three areas of stable, sliding and toppling for the 60 degree slope. The charts for other slopes are presented as Fig. E.3 to E.6 in

ur

Appendix E. In these figures, the safety factor with value of 1 represents the sliding boundary.

Jo

Also, for each slope, according to what was obtained in the previous section, a critical amount of a/x was considered for the toppling boundary. If the slope is not listed in any of the two conditions, it is considered to be stable. It should be noted that, the Young’s modulus was considered to be 500 MPa. In the previous section, it was observed that this value for the case studied is a critical amount which the toppling can occur in the slope. The Young’s modulus does not affect the stability calculations against the slope sliding. This value will only determine the toppling and the stable boundary of the slope (i.e., affecting the a/x determination). As it is observed in Fig. 22 and Figs. E.3 to E.6, with increasing a/x, the safety factor of the slope will

Journal Pre-proof increase against sliding. This is especially noticeable in the variation of a/x from -0.5 to 0. By increasing the slope angle, the stability of the slope will increase. Also, with a decrease in the slope angle, the possibility of toppling will be more dominant than the shear failure. Thus, the

re

-p

ro

of

higher the strength of rock/soil mass, the greater the possibility of toppling.

(b)

Jo

ur

na

lP

(a)

(c)

(d)

Fig. 22 The boundary between shear sliding, toppling and stable for the 60 degree slope with: (a) friction angle of 20; (b) friction angle of 25; (c) friction angle of 30; (d) friction angle of 35.

6. Conclusion This paper addresses the stability analysis of rock slopes prone to toppling-slumping failure mode. This failure mode occurs when a weak horizontal layer of soil or rock is overlaid by

Journal Pre-proof vertical rock columns. The slope failure can occur due to the shear failure of the weak horizontal soil/rock layer or as a result of differential settlement of the vertical rock columns within the continuous horizontal soil/rock layer or as a combination of the two mechanisms. In this study, the stress distribution in a slope was obtained through the theory of elasticity, and subsequently the displacement under the individual vertical rock columns was calculated. The

of

results indicate that the vertical stress distribution plays an important role in the settlement of the rock columns. It was also found that by increasing the normalized column distance, a/x, and

ro

decreasing slope angle, the stress distribution under the rock columns is more symmetrical, and

-p

at a certain distance from the slope crest, the set of rock columns will act like a single column on

re

a horizontal ground surface. Therefore, stress distribution and displacement will be similar to the

lP

stress distribution and displacement on a half plane. Also, it was shown that the differential settlement of the first column strongly depends on the slope angle, while, by moving to the end

na

column, the stress distribution under the column is not depend on the slope angle.

ur

It was shown that differential settlement could lead to the toppling of rock columns in a slope. In

Jo

addition, the a/x and slope angle play a significant role on the toppling of the vertical rock columns. Weathered rock or soil slopes are susceptible to the toppling due to differential settlement. In this research, it was assumed that the slope is homogeneous and isotropic; therefore, Young’s modulus and Poisson’s ratio were remained the same for all points of the slope. As stated by Evans (1981) and Wyllie (2018) studies have shown that weathering degree in rock slopes decreases by increasing the distance from the free face. Therefore, Young's modulus near the slope crest can be lower than other points inside the slope. This can cause more settlement of the column near the slope crest and thus increase the probability of toppling of the

Journal Pre-proof columns. In addition, weathering and erosion can decrease the a/x value and put the stable slope in a critical situation where toppling can occur. The stress based criterion was also proposed and the probability of shear sliding of the slope was investigated. Based on the proposed displacement and stress criterion, several design charts were developed to determine the boundary between the toppling, shearing and stability conditions of

of

the slope according to the geometric characteristics and shear strength parameters of the slope. The results indicate that the shear or toppling failure of the slope not only depends on the

ro

material’s strength but also on the geometry and the normalized column distance (i.e. a/x). It can

-p

be concluded that considering the normalized column distance in slope stability is of great

re

importance.

lP

Appendices

na

Appendix A

ur

Parameters A to D in Airy stress function:

 z  2  cos   z  2    F  z  2z  z  1  z  1 sin  2   sin  2  z  1  

B

 z  2  sin   z  2    F z  2z  z  1  z  1 sin  2   sin  2  z  1  

C 

F z  cos  z   2  z  1  z  1 sin  2   sin  2  z  1  

Eq. (A.3)

D

F  z  sin  z   2  z  1  z  1 sin  2   sin  2  z  1  

Eq. (A.4)

Jo

A

Eq. (A.1)

Eq. (A.2)

Journal Pre-proof Appendix B In the positive part of the imaginary axis:

x  1  c

0 y 

Z  1  iy

Eq. (B.1)

dz  idy

y y

Z  1  iy

ro

  y  0

dz  idy

y  y

Eq. (B.1) Eq. (B.2)

lP

re

Appendix C

-p

x  1  c

of

In the negative part of the imaginary axis:

Values of R, T, K, and L related to Equation (16):

K 

L

na

 y sin 2  sinh  2 y  

ur

T 

sin     cosh     y  sin     cosh     y

cos     sinh     y  cos     sinh     y

 y sin 2  sinh  2 y  

Jo

R

cos     sinh     y  cos     sinh     y

 y sin 2  sinh  2 y  

sin     cosh     y  sin     cosh     y

 y sin 2  sinh  2 y  

Appendix D Filon numerical integration parameters:

Eq. (C.1)

Eq. (C.2)

Eq. (C.3)

Eq. (C.4)

Journal Pre-proof h

sin  2  2sin 2   Ω    2 2 3

  kh

x 2n  x 0 2n

1

 sin  

1  cos 2   sin  2     2  3  

  4

  2

 

n

S 2 n  [f x 2 i sin  kx 2i  0.5 f x 2 n sin  k x 2 n   f x 0 sin  kx 0 ] i 0

n

C 2 n  [f x 2 i cos  kx 2i  0.5 f x 2 n cos  k x 2 n   f x 0 cos  kx 0 ]

ro -p re lP na ur Jo



cos     2 

n

S 2 n 1  f x 2 i 1 sin  kx 2i 1  i 1

n

C 2 n 1  f x 2 i 1 cos  kx 2i 1 

of

i 0

3

i 1

Eq. (D.1)

Journal Pre-proof

Jo

ur

na

lP

re

-p

ro

of

Appendix E

Fig. E.1 Schematic representation of discrete loading by use of the step function.

ro

of

Journal Pre-proof

Jo

ur

na

lP

re

-p

Fig. E.2 Schematic representation of Finite Difference network in a slope.

(a)

(b)

Journal Pre-proof

of

(c)

(d)

ur

na

lP

re

-p

ro

Fig. E.3 The boundary between shear sliding, toppling and stable for the 90 degree slope with: (a) friction angle of 20; (b) friction angle of 25; (c) friction angle of 30; (d) friction angle of 35.

(b)

Jo

(a)

(c) (d)

Journal Pre-proof

ro

of

Fig. E.4 The boundary between shear sliding, toppling and stable for the 80 degree slope with: (a) friction angle of 20; (b) friction angle of 25; (c) friction angle of 30; (d) friction angle of 35.

(b)

ur

na

lP

re

-p

(a)

Jo

(c)

(d)

Fig. E.5 The boundary between shear sliding, toppling and stable for the 75 degree slope with: (a) friction angle of 20; (b) friction angle of 25; (c) friction angle of 30; (d) friction angle of 35.

of

Journal Pre-proof

(a)

ur

(c)

na

lP

re

-p

ro

(b)

(d)

Jo

Fig. E.6 The boundary between shear sliding, toppling and stable for the 70 degree slope with: (a) friction angle of 20; (b) friction angle of 25; (c) friction angle of 30; (d) friction angle of 35.

References

Alejano L, Carranza-Torres C, Giani G, Arzua J (2015) Study of the stability against toppling of rock blocks with rounded edges based on analytical and experimental approaches Eng Geol 195:172184 Alejano LR, Gómez-Márquez I, Martínez-Alegría R (2010) Analysis of a complex toppling-circular slope failure Eng Geol 114:93-104 Amini M, Ardestani A, Khosravi MH (2017) Stability analysis of slide-toe-toppling failure Eng Geol 228:8296 Amini M, Majdi A, Veshadi MA (2012) Stability analysis of rock slopes against block -flexure toppling failure Rock Mech Rock Eng 45:519-532 Amini M, Sarfaraz H, Esmaeili K (2018) Stability analysis of slopes with a potential of slide-head-toppling failure Int J Rock Mech Min Sci 112:108-121 Borgatti L et al. (2015) The 27 February 2014 San Leo landslide (northern Italy) Landslides 12:387-394

Journal Pre-proof

Jo

ur

na

lP

re

-p

ro

of

de Souza Neto EA, Peric D, Owen DR (2011) Computational methods for plasticity: theory and applications. John Wiley & Sons, Duncan CW (2018) Rock slope engineering: civil applications. Fifth edn. CRC Press, New York Erguler ZA (2009) Field-based experimental determination of the weathering rates of the Cappadocian tuffs Eng Geol 105:186-199 Evans R (1981) An analysis of secondary toppling rock failures —the stress redistribution method Quarterly Journal of Engineering Geology and Hydrogeology 14:77-86 Filon LNG (1930) III.—On a Quadrature Formula for Trigonometric Integrals Proceedings of the Royal Society of Edinburgh 49:38-47 Frayssines M, Hantz D (2009) Modelling and back-analysing failures in steep limestone cliffs Int J Rock Mech Min Sci 46:1115-1123 Goodman L, Brown C (1963) Dead load stresses and the instability of slopes Journal of the Soil Mechanics and Foundations Division 89:103-136 Goodman R, Bray J (1976) Toppling of rock slopes In: Rock engineering for foundations and slopes; proceedings of a specialty conference, Vol. 2 p201-233 Am Soc Civ Eng New York Heidari M, Momeni A, Naseri F (2013) New weathering classifications for granitic rocks based on geomechanical parameters Eng Geol 166:65-73 Jovanović BS, Süli E (2013) Analysis of finite difference schemes: for linear partial differential equations with generalized solutions vol 46. Springer Science & Business Media, Pasuto A, Soldati M (2013) Lateral spreading Treatise on Geomorphology 7: 239-248 Rabczuk T, Belytschko T (2004) Cracking particles: a simplified meshfree method for arbitrary evolving cracks International Journal for Numerical Methods in Engineering 61:2316-2343 Rabczuk T, Zi G, Bordas S, Nguyen-Xuan H (2010) A simple and robust three-dimensional crackingparticle method without enrichment Computer Methods in Applied Mechanics and Engineering 199:2437-2455 Ren H, Zhuang X, Cai Y, Rabczuk T (2016) Dual‐horizon peridynamics International Journal for Numerical Methods in Engineering 108:1451-1476 Sadd MH (2009) Elasticity: theory, applications, and numerics. Academic Press, Spreafico MC, Cervi F, Francioni M, Stead D, Borgatti L (2017) An investigation into the development of toppling at the edge of fractured rock plateaux using a numerical modell ing approach Geomorphology 288:83-98 Spreafico MC et al. (2016) Back analysis of the 2014 San Leo landslide using combined terrestrial laser scanning and 3D distinct element modelling Rock Mech Rock Eng 49:2235-2251 Stianson JR, Chan D, Fredlund D Comparing slope stability analysis based on linear elastic or elastoplastic stresses using dynamic programming techniques. In: Proc., 57th Canadian Geotechnical Conf, 2004. pp 23-30 Tranter CJ (1951) Integral transforms in mathematical physics. Tsesarsky M, Hatzor YH (2009) Kinematics of overhanging slopes in discontinuous rock Journal of geotechnical and geoenvironmental engineering 135:1122-1129 Vu-Bac N, Lahmer T, Zhuang X, Nguyen-Thoi T, Rabczuk T (2016) A software framework for probabilistic sensitivity analysis for computationally expensive models Advances in Engineering Software 100:19-31

Journal Pre-proof Conflict of Interest

The authors have no affiliation with any organization with a direct or indirect financial interest in

ro

of

the subject matter discussed in the manuscript.

-p

Authors: Hadi Haghgouei ([email protected])

2.

Ali Reza Kargar ([email protected])

3.

Mehdi Amini ([email protected])

4.

Kamran Esmaeili ([email protected])

Jo

ur

na

lP

re

1.

Journal Pre-proof Highlights of the manuscript can be summarized as follows: 

A fundamental mathematical solution is proposed to determine the stress distribution in the slope due to the surcharge loading.



A new analytical approach is presented for analysis of toppling-slumping failure in rock slopes.



A Mathematica package code was developed to evaluate the stress distribution within the slope and toppling-slumping. Results of the presented analytical method are compared with those of a numerical

of



technique.

ur

na

lP

re

-p

toppling-slumping failure mechanism.

ro

Several design-charts are developed for determine the stability of a slope against

Jo