Boundary effects of rainfall-induced landslides

Computers and Geotechnics 61 (2014) 341–354

Contents lists available at ScienceDirect

Computers and Geotechnics journal homepage: www.elsevier.com/locate/compgeo

Boundary effects of rainfall-induced landslides Abid Ali a,⇑, Jinsong Huang a, A.V. Lyamin a, S.W. Sloan a, M.J. Cassidy a,b a b

ARC Centre of Excellence for Geotechnical Science and Engineering, The University of Newcastle, Callaghan, NSW 2308, Australia Centre for Offshore Foundation Systems and ARC Centre of Excellence for Geotechnical Science and Engineering, The University of Western Australia, Crawley, WA 6009, Australia

a r t i c l e

i n f o

Article history: Received 15 January 2014 Received in revised form 14 May 2014 Accepted 30 May 2014

Keywords: Boundary conditions Landslides Rainfall Drainage Hydraulic conductivity Infinite slope

a b s t r a c t In the study of landslides, it is generally assumed that an impermeable boundary exists at a certain depth and failure occurs at this boundary. In reality this is not always the case and failures can occur at any depth. This paper aims to study the effect of boundary conditions on landslides, using a series of seepage and stability analyses performed over a range of rainfall intensities, and for different failure mechanisms, by studying the failure time and depths corresponding to fully drained, partially drained, and impermeable boundaries. It is shown that these conditions can significantly affect the occurrence and depth of rainfall-induced landslides. Ó 2014 Elsevier Ltd. All rights reserved.

1. Introduction Landslides are responsible for causing loss of life and extensive damage to property in many regions throughout the world. Even though an otherwise stable slope may fail due to human-induced factors, such as excavation at the toe or loading due to construction, many slopes simply fail due to rainfall infiltration [2,6]. Rainfall induced landslides are common in tropical and subtropical regions where residual soils exist in slopes and there are negative pore water pressures in the unsaturated zone above the water table [14]. In an unsaturated soil, these negative pore water pressures contribute toward its shear strength and thus help to maintain stability [7]. The infiltration of rainwater causes a reduction in this negative pore water pressure and an increase in the soil unit weight (due to an increased saturation), both of which have a destabilizing influence [25]. Rainfall-induced landslides are characterized as failures occurring along a plane parallel to the ground surface. Given the physics involved in landslides, various infinite slope models have been used to assess their stability after heavy rainfall [4,8,11,17,21,25,26]. These models are based on the premise that each slice of an infinitely long slope receives the same amount and intensity of rainfall [5], the time required for infiltration

⇑ Corresponding author. Tel.: +61 2 4985 4974. E-mail address: [email protected] (A. Ali). http://dx.doi.org/10.1016/j.compgeo.2014.05.019 0266-352X/Ó 2014 Elsevier Ltd. All rights reserved.

normal to the slope is far less than that required for flow parallel to the slope, the wetting front propagates in a direction normal to the slope [25], and the failure depth is small compared to the length of the failing soil mass. The validity of these assumptions has been checked against the predictions of two-dimensional numerical models, with the conclusion that an infinite slope approximation may be adopted as a simplified framework to assess failures due to the infiltration of rainfall [11,26]. The present study is based on an infinite slope model. Research indicates that several factors affect the stability of a slope subjected to rainfall infiltration. The characteristics of the rainfall (duration, intensity and pattern), the saturated hydraulic conductivity of soil, the slope geometry, the initial conditions and the boundary conditions have all been identified as influential factors. Studies on rainfall-induced landslides have considered different types of boundary conditions: an impermeable boundary condition (e.g. [11,13,25]), fixed water table (e.g. [12,26]) or a drained boundary condition (e.g. [5,18]). Few studies [20,26] have also reported that when a soil slope is underlain by a less permeable layer, the interface (between the soil and the less permeable layer) acts as a pseudo-boundary, (generating positive pore water pressures, thus) causing the slope to fail at the interface. In such a case, the location of the interface could significantly affect the stability of a layered slope [26]. In addition, rainfall-induced landslides could be shallow (i.e. at an intermediate depth in the soil profile) or deep

342

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

Nomenclature a c0 Df dt FN FT FS FS1 I K Kr Ks m n N q S Se Sr t tf

scaling suction effective cohesion depth of failure incremental time step normal force tangential force factor of safety factor of safety at the depth of failure rainfall intensity hydraulic conductivity relative hydraulic conductivity saturated hydraulic conductivity van Genuchten model parameter soil porosity van Genuchten model parameter flux infiltrating the slope degree of saturation effective degree of saturation residual degree of saturation rainfall duration time of failure

(i.e. at the bottom of the soil profile). Some researchers (e.g. [5,13,22]) have reported the generation of positive pore water pressures (at shallow depths) during rainfall infiltration and considered it to be the cause of shallow landslides. However Cheng [3] showed that the positive pore water pressures that are generated during infiltration in a homogeneous soil profile, are in fact due to the high nonlinearities associated with the solution (of seepage flow problem (both spatial and temporal) and if a sufficiently fine mesh was considered, no positive pore pressures would be generated). Therefore, shallow failures during infiltration in homogeneous soils occur due to a reduction in the negative pore water pressures [11,12]. Li et al. [11] reported that if failure does not occur during propagation of wetting front then it will occur at the boundary (as they considered an impermeable boundary in their study) due to a rise in the water table. A study which systematically evaluates the response of a slope having different failure mechanisms, to various boundary conditions has still not been performed. Therefore, this paper addresses the effects of the boundary conditions on rainfall-induced landslides. Three different types of boundaries are considered: fully drained, partially drained and impermeable. Although the factor of safety is an important aspect, it is also imperative to know what would be the location and the timing of the collapse mechanism if a slope was to fail. Thus, in the present study, the time of failure tf and the depth of failure Df are studied for different boundary conditions. The depth of failure can be used to estimate the volume of the failing soil mass. Four case studies are performed to assess the effect of the boundary conditions on the depth of failure and time of failure for different rainfall intensities, slope inclinations and soil properties. The different boundary conditions are approximated by adopting different ratios of the saturated hydraulic conductivity of the overlying soil layer relative to the underlying rock. It is shown that the boundary conditions significantly affect the failure process by influencing either the occurrence of failure or the failure depth for a given rainfall intensity. To obtain the pore water distributions for a slope subjected to rainfall, numerical simulations of 1D seepage have been performed using the HYDRUS 1D [16] software.

u uw ua W D z z0

pore water pressure head pore water pressure pore air pressure weight of failing soil mass slope depth slope normal direction vertical direction slope angle unit weight of soil unit weight of soil solids unit weight of water soil volumetric water content residual water content saturated water content total normal stress effective normal stress soil shear stress soil shear strength effective friction angle coefficient of effective stress

a c cs cw h hr hs

r r0 s sf /0

v

2. Seepage analysis Assuming that the effect of pore-air pressure is insignificant and that water flow due to thermal gradients is negligible, one-dimensional uniform flow in a variably saturated soil can be described by a modified form of Richards equation [15]. Therefore, the flow in an unsaturated infinite soil slope can be described by the 1D equation (e.g. [26]):

   dh d du K ¼ þ cos a dt dz dz

ð1Þ

where h is the volumetric water content, t is time, u is the pore water pressure head, a is the inclination of the slope to the horizontal, K is the hydraulic conductivity and z is the spatial coordinate as shown in Fig. 1. To solve the above equation numerically, the water content h is assumed to vary with the pore water pressure head u according to the van Genuchten [23] model as:

" #m h  hr 1 Se ¼ ¼ hs  hr 1 þ ðauÞN

ð2Þ

where Se is the effective degree of saturation, hs and hr are the saturated and residual water content respectively, a is the suction scaling parameter and N, m are the parameters of the van Genuchten model. Noting that the volumetric water content is related to the degree of saturation S and the porosity n (h = nS), the effective degree of saturation can also be expressed in terms of the degree of saturation S in the following form:

Se ¼

S  Sr 1  Sr

ð3Þ

where Sr is the residual degree of saturation. To complete the description, the hydraulic conductivity K can be estimated as:

K ¼ K sK r

ð4Þ

where Ks is the saturated hydraulic conductivity and Kr is the relative hydraulic conductivity given by van Genuchten [23]: m 2

K r ¼ Se1=2 ½1  ð1  S1=m Þ  e

ð5Þ

343

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

illustrated in Fig. 1. The corresponding shear stress on the failure surface is given by:

α

1 unit

z

s ¼ W cos a sin a

z /cosα

In this study, the shear strength of soil sf (in terms of effective stress) is assumed to be described by the Mohr–Coulomb model:

I O

q=

ð8Þ

sf ¼ c0 þ r0 tanð/0 Þ

α I cos

Df

Soil,

0

where c is the effective cohesion, r is the effective normal stress, and /0 is the effective friction angle. To consider the influence of pore water pressure on the shear strength of a variably saturated soil, Terzaghi’s effective stress principle is modified according to the formulation of Bishop [1]:

D

K s1 Df /cos α

e Plan ure

r0 ¼ ðr  ua Þ þ vðua  uw Þ

Fail

W FT

Roc

ð9Þ 0

ð10Þ

where ua is the pore-air pressure, uw = cwu is the pore water pressure, (ua  uw) is known as matric suction and v is called the coefficient of effective stress and is a constitutive property of the soil that depends on the degree of saturation. For a variably saturated soil, v denotes the proportion of matric suction that contributes to the effective stress and generally varies between 0 (for a perfectly dry soil) and 1 (for a completely saturated soil). Though many mathematical forms of v have been proposed in the past, in the present study v is considered equal to the effective degree of saturation, Se [24]. With ua = 0 and v = Se, Eq. (10) can be written as:

k, K s 2

FN

α

r0 ¼ r  Se uw

ð11Þ

Slope failure occurs when the applied shear stress s exceeds the mobilized soil shear strength sf. The factor of safety FS can then be computed as:

Fig. 1. Limit equilibrium set-up.

sf tan /0 c0  Se uw tan /0 ¼ þ s tan a W cos a sin a

2.1. Slope stability assessment

FS ¼

Once the pore water pressure distribution is obtained through seepage analysis, the factor of safety FS at any given time t can then be determined by limit-equilibrium techniques. The stability of an infinite slope is estimated by using a closed form solution similar to that proposed by White and Singham [25], where the failure is considered to occur along a plane parallel to the ground surface. A soil column of a unit width is considered, where the self-weight W is used to obtain the normal force FN and tangential force FT at any depth. The expression for the factor of safety is derived along the same lines as White and Singham [25]. Referring to Fig. 1, the total normal stress r on the failure surface can be computed as:

where FS = 1 corresponds to a limiting condition for equilibrium, and failure occurs when FS is less than 1.

r ¼ W cos2 a

ð6Þ

where W is the self-weight of the soil column above the failure plane. The unit weight of soil during rainfall infiltration will increase due to an increase in water content; neglecting this effect (i.e. assuming a constant unit weight) could lead to a higher (and thus false) estimate of stability [19,26]. Therefore, to account for the variation in unit weight due to the variation in the water content with depth, W is determined as:

W¼

Z

Df = cos a

0

¼

1 cos a

Z 0

c

dz 1 ¼ cos a cos a

Df = cos a

Z 0

Df = cos a

½ð1  nÞcs þ nScw dz

½ð1  nÞcs þ hcw dz

ð7Þ

where Df is the failure depth, c is the bulk unit weight of soil, cs is the unit weight of soil solids, and cw is the unit weight of water and the water content h is defined in Eq. (2). Note that the integration takes place along the vertical direction as

ð12Þ

3. First failure determination To estimate the time of failure tf, the time discretization must be sufficiently accurate. In the present study, an iterative search method was developed to accurately determine the time of failure and is described below. Step 1. Assume the time duration t for the seepage analysis. The value of t is chosen such that steady state conditions are achieved for a fully drained boundary. (Note: For other boundary conditions the same value of t can be used). Larger values of t are required for low rainfall intensities and as the rainfall intensity increases, t is reduced. Step 2. Discretize the slope depth D into j nodes and the time duration t into 250 time steps, dt (the maximum value of 250 time steps is a limitation of HYDRUS 1D). If t > 0.025 days, dt = t/250, else dt = 104 days. Step 3. Obtain the pore water pressure head u and the water content h at each time step for the entire profile using seepage analysis (from HYDRUS 1D). Step 4. From Eqs. (2) and (7) calculate the effective degree of saturation Se, and the self-weight W respectively, for the entire profile at each time step (using the water content h obtained from seepage analysis). From Eq. (12), calculate FS for the entire profile at each time step. Step 5. Find FS1, Df and tf. FS1 is the factor of safety of the slope when it fails at a time tf and depth Df.

344

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

Discretize slope depth D and time t

t

t - t1

Assign t tf

Perform seepage analysis for t days

Obtain u, θ for the entire profile at each time step dt

t1

10-4 days

No Apply Eq. (12) to calculate FS

Find FS1 , Df and t f. FS1 min (FS(i,j) 1); tf t(i); Df D(j) FS1 cannot be found

FS1 < 0.999

t

tf?

Yes

Calculate t 1

FS1 > 0.999

No failure Accept FS1, tf and Df

t1

10-4 days

Find FSmin and Dmin at time t

end

Fig. 2. Flow chart illustrating the search strategy to find FS1, Df and tf.

i. For the first time step find all depths where FS < 1. If the condition FS < 1 is not met in the first time step then move to the next time step until FS < 1 is true. The ith time step where FS is less than 1 is tf. From all depths where FS < 1, the minimum value (of FS) is FS1 and the corresponding depth is Df. (If failure does not occur, then we report the minimum factor of safety FSmin and the corresponding depth Dmin at time t. The analysis is then terminated).1 ii. Check the value of FS1. If FS1 > 0.999, accept the results (FS1, Df, tf) as satisfactory and exit the analysis. iii. If FS1 < 0.999 and tf < t, reduce t to tf and then rerun the analysis (i.e. go to Step 1.). Step 6. If FS1 < 0.999 and t = tf, check if t can be further reduced to bring FS1 closer to 0.999. This is done by finding t1.

1 At t = tf, the FS may become less than unity at several depths in the slope (or beyond a certain depth as fixed time intervals are considered). The minimum FS (i.e. FS1) is reported as the factor of safety of the slope. The location where FS is minimum (at t = tf) is the point where the failure initiates in the slope. Therefore, we refer to the minimum value of FS as FS1, and the depth corresponding to FS1 as the depth of failure, Df.

t1 ¼

0:999  FS1  dt ðFSði  1; Df Þ  FS1 Þ

ð13Þ

where FS((i  1), Df) is the factor of safety at depth Df for the (i  1)th time step. Step 7. If FS1 < 0.999, tf = t and t1 < 104 days, accept the results (FS1, Df, tf) as satisfactory and exit the analysis. The steps described above are illustrated in Fig. 2. The critical value of factor of safety to find tf is taken as 0.999 so that tf is accurately determined. If FS is less than 0.999, it is possible that there may be a large difference between the actual time of failure and that detected as tf. However, for the Step 7 the minimal time step of 104 supersedes the preference given to FS i.e. if FS < 0.999 and t1 < 104 the results are accepted as satisfactory as the time difference between actual failure and tf is less than 9 s in real time, which is quite small. 4. Case studies A hypothetical slope, at an inclination of a to the horizontal, is considered as shown in Fig. 1. A 1 m thick homogeneous soil layer

345

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

intensity is I, and the flux q infiltrating the soil at the top is given by Icos a. No ponding of water is allowed at the ground surface at any time. The mechanical and hydraulic properties of the soil are given in Table 1 [25]. The soil–water characteristic curve (SWCC) and the permeability curve for the soil are shown in Fig. 3. The numerical model described above was validated by comparing its results with those published by White and Singham [25], when using similar initial conditions and soil properties. The failure time obtained from the analysis, tf = 81.72 h, is very close to the value of tf = 88.08 h reported by White and Singham [25]. The difference in the results (approximately 7%) can be explained due to the different approximations of v; in the present formulation v = Se while White and Singham [25] considered it equal to S. If v = S is considered, tf = 87.48 h which agrees well with the reported time of tf = 88.08 h thereby confirming the validity of the new implementation. Fig. 4 shows the results from the new algorithm (v = Se) for the slope considered by White and Singham [25]. Three different boundary conditions are implemented at the bottom of soil layer namely, fully drained, partially drained and

Table 1 Material properties of the soil [25] considered in the present study. Parameter

Symbol

Value

Porosity Unit weight of solids Unit weight of water Saturated hydraulic conductivity Residual water content Saturated water content Scaling suction Shape parameter Shape parameter Effective cohesion

n

Effective friction angle

/0

0.4 20 10 8.64 0.128 0.4 5 1.5 (N  1)/N 0 (Case-1, 3) 1 (Case-2, 4) 35

cs cw Ks hr hs a N m c0

Units kN/m3 kN/m3 m/day

1/m

kPa deg.

is underlain by rock (i.e. D = 1 m). The water table is 2 m (measured normal to the slope) below the ground surface and is assumed to be constant throughout the analysis. The initial pore water pressure profile is hydrostatic. It is assumed that as rainfall occurs, water infiltrates in a wetting front normal to the slope. The rainfall

SWCC

Permeability curve

Hydraulic conductivity (m/day)

0.45

Volumetric water content [-]

0.4 0.35 0.3 0.25 0.2 0.15 0.1 10

-4

10

-2

10

0

10

0

10

-5

10

2

10

-4

10

Suction head [m]

-2

10

0

10

2

Suction head [m]

Fig. 3. Soil-water characteristic curve (left) and the permeability curve (right) of the soil.

0

-0.2

-0.2

-0.4

-0.4

Depth (m)

Depth (m)

FS1 = 0.9996, Df = -1.00 m, tf = 81.72 hrs 0

-0.6

-0.8

-0.8

-1 -20

-0.6

-1 -15

-10

-5

0

← Failure

0

1

0

10

20

30

2

3

4

Factor of Safety

Pore water pressure (kPa) 40

50

Time (hrs) Fig. 4. Slope results from new procedure.

60

70

80

346

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

Table 2 Boundary conditions investigated in the present study.

Table 3 Summary of parametric study.

Sr. no.

Name

Ks1/Ks2

1 2 3 4 5 6 7 8

Fully drained

Very small (<0.001) 0.01 0.1 1 10 100 1000 Very large (>10,000)

Partially drained

Impermeable

impermeable. The type of boundary condition depends on the ratio Ks1/Ks2, where Ks1 and Ks2 are the respective hydraulic conductivities of the soil and the underlying rock. The ratio Ks1/Ks2 is simply used as an index to represent the relative degree of permeability at the bottom of the soil profile. For a fully drained boundary, there is no accumulation of infiltrated flux as it drains completely. This implies that the initial negative pore water pressure at the bottom of soil remains constant, and it happens when the ratio Ks1/Ks2 is very small (Ks1/Ks2 < 0.001). A fully drained condition is simulated by assigning a constant head at the bottom of the soil layer (the value of head being equal to the initial pressure head at the bottom of soil layer). For a partially drained boundary the pore water pressure at the boundary increases during flow. This is simulated by considering the range of Ks1/Ks2 values described in Table 2. For an impermeable boundary, no flow can occur across the boundary and infiltrated flux accumulates at the boundary. This condition is implemented by assigning a zero flux condition at the bottom of the soil layer. For the cases of a fully drained boundary or an impermeable boundary, only the top 1 m of soil is considered in the flow and stability analysis. For the case of a partially drained boundary, the top 2 m of the soil profile is considered for the flow simulations while the top 1 m is considered for the stability assessment. The mechanical deformations in the soil skeleton due to the flow of water are assumed to be negligible and are not considered in the analysis. It is also assumed that no failure can occur in the rock. The various boundary conditions considered in the analysis are summarized in Table 2. The slope described above was subjected to 24 different rainfall intensities I, defined by the ratio I/Ks1, for each of the three boundary conditions. In general, the influence of the boundary conditions will depend on the inclination of the slope and the material properties. In this section, we summarize the results of the stability assessment of the infinite slope model just described. The section is divided into four parts. Case 1 studies the impact of the boundary conditions on the failure depths for a cohesionless slope, while Case 2 investigates the same effect for a cohesive slope. Cases 3 and 4, respectively, consider the influence of the boundary conditions on the failure depths for a steep cohesionless slope and a steep cohesive slope. Table 3 summarizes the results for all the cases. Except for the cohesion c0 2 and slope angle a, all the parameters are held constant, including the friction angle /0 (=35°).

4.1. Case 1: slope with a = 36° and c0 = 0 kPa To investigate the effect of the boundary conditions on a slope in cohesionless soil, the slope inclination a was set to 36° and the drained cohesion c0 was set to 0 kPa. The slope was subjected 2 In a practical sense, the value of c0 = 1 kPa (Case 2,4) is very low. However, this value was considered to illustrate the different failure mechanisms that will be discussed later.

Case 1 2 3 4

Boundary conditions Table 2 (1–8)

a

c0 (kPa)

Rainfall intensity I/Ks1

36°

0 1 0 1

0.01, 0.02, 0.03, 0.04, 0.05, 0.0625, 0.75, 0.875, 0.1, 0.125, 0.150, 0.175, 0.20, 0.25, 0.30, 0.35, 0.40, 0.45, 0.50, 0.60, 0.70, 0.80, 0.90, 1.00

40°

to 24 rainfall intensities (see Table 3), and the response of the pore water pressure and the factor of safety were recorded. The pore water pressure and factor of safety profiles for selected I/Ks1 values of 0.01, 0.25, 0.30 and 0.40 are shown in Figs. 5–8 respectively. For the boundary conditions where failure does not occur, the results are plotted for steady-state flow conditions. Referring to Fig. 5, I/Ks1 = 0.01 represents the situation where the rainfall intensity is small compared to the saturated hydraulic conductivity of the soil. Previous research [9–11,13,27] indicates that for the same rainfall duration, the reduction in the negative pore water pressure is directly proportional to I/Ks1. In the present study, I/Ks1 = 0.01 is the lowest rainfall intensity considered and thus gives the least reduction in the negative pore water pressures. At this rainfall intensity, the pore water pressure and factor of safety profiles are very sensitive to the type of boundary condition, as its effect can be observed at a shallow depth (0.40 m). Beyond this depth, the pore water pressure and factor of safety profiles differ depending upon Ks1/Ks2. For more permeable boundaries (i.e. the fully drained boundary and those partially drained boundaries where Ks1/Ks2 < 1), the pore water pressure and factor of safety at the boundary tend toward their initial values. Indeed, at a fully drained boundary, the pore water pressure and factor of safety are exactly the same as their initial values. As the boundary becomes less permeable (i.e. Ks1/Ks2 increases), there is a greater reduction in the negative pore water pressure at the boundary and, consequently, a greater reduction in the factor of safety. However when I/Ks1 = 0.01, failure (denoted by cross) occurs only when the boundary is impermeable. The pore water pressure and factor of safety profiles for I/Ks1 = 0.25 are shown in Fig. 6. In this situation, a greater reduction in negative pore water pressure can be observed in the soil profile. Compared to I/Ks1 = 0.01, boundary effects can be noticed only beyond a larger depth (0.90 m). This situation also represents the maximum rainfall intensity where the failure is still deep and no failure is observed for the examples with more permeable boundaries. The pore water pressure and the factor of safety profiles for I/Ks1 = 0.30 are shown in Fig. 7. This value of I/Ks1 represents the onset of shallow failure (i.e. failures which occur above the boundary). Shallow failure occurs only for Ks1/Ks2 = 1, with the deep failures occurring for the less permeable boundaries. In this case, the boundary conditions clearly affect the depth of failure. The pore water pressure and factor of safety profiles for I/Ks1 = 0.40, shown in Fig. 8, indicate that all the failures are shallow and occur independently of the boundary conditions. Comparing these results to those of I/Ks1 = 0.30, it can be concluded that as the rainfall intensity increases, the failure occurs at a smaller depth. In this case, the boundary conditions do not have any effect on the failure depths and times. From Figs. 5–8 it can be observed that a residual negative pore water pressure exists at the depth of failure, implying that failure occurs due to a reduction in the negative pore water pressure. Fig. 9 shows the variation in the depth of failure Df and the time of failure tf with I/Ks1 for various boundary conditions. In this figure

347

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354 Ks1/Ks2 FD

0.01

0.1

1

10

100

1000

IMP

0

0

Initial Profile

-0.2

-0.4

-0.6

-0.4 Failure

Initial Profile

Depth (m)

Depth (m)

-0.2

-0.6

-0.8

-0.8

-1

-1 -15

-10

-5

1

0

1.5

Pore water pressure (kPa)

2

2.5

3

Factor of Safety

Fig. 5. Pore water pressure and factor of safety profiles at I/Ks1 = 0.01 for Case 1. (For interpretation of the figure legends, the reader is referred to the web version of this article).

Ks1/Ks2 FD

0.01

0.1

1

10

0

100

1000

IMP

0

-0.2 -0.9

Initial Profile

-0.95

-0.4

-1 -3

-2 -1 Pore water pressure (kPa)

0

-0.6

-0.4

-0.6

Failure

Depth (m)

-0.2

Depth (m)

Depth (m)

-0.8 -0.85

-0.8

-0.8

Depth(m)

Initial Profile

-0.8

-0.9

-1 1

-1

1.05 1.1 1.15 Factor of Safety

-1 -15

-10

-5

0

1

1.5

Pore water pressure (kPa)

2

1.2

2.5

3

Factor of Safety

Fig. 6. Pore water pressure and factor of safety profiles at I/Ks1 = 0.25 for Case 1.

Ks1/Ks2 FD

0.01

0.1

1

10

100

1000

IMP

0

0

-0.4

-0.2

-0.9

Initial Profile

-0.95 -1 -3

-2 -1 Pore water pressure (kPa)

0

-0.6

-0.4

-0.6

Failure

-0.8

Initial Profile

Depth(m)

Depth (m)

-0.2

Depth (m)

Depth (m)

-0.8 -0.85

-0.8

-0.8

-0.9

-1 1

1.1 1.2 Factor of Safety

2

2.5

-1

-1 -15

-10

-5

Pore water pressure (kPa)

0

1

1.5

Factor of Safety

Fig. 7. Pore water pressure and factor of safety profiles at I/Ks1 = 0.30 for Case 1.

3

348

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354 Ks1/Ks2 FD

0.01

0.1

1

10

100

0

0

-0.2

-0.2

1000

IMP

-0.6

-0.4 Failure

-0.4

Depth (m)

Depth (m)

Initial Profile

-0.6

Initial Profile

-0.8

-0.8

-1

-1 -15

-10

-5

0

1

1.5

Pore water pressure (kPa)

2

2.5

3

Factor of Safety

Fig. 8. Pore water pressure and factor of safety profiles at I/Ks1 = 0.40 for Case 1.

Ks1/Ks2 FD

0

0.01

Boundary effect

0.1

1

10

100

IMP

No boundary effect

0

10

tf (days)

Critical I/Ks1

-0.2

1

tf (days)

-0.4

Df (m)

1000

1.5

-0.6

-1

10

-2

10

0.5

0

0.2

0.4

0.6

0.8

1

I/Ks1

-0.8

-1 0

0.2

0.4

0.6

0.8

1

I/Ks1

0

0

0.2

0.4

0.6

0.8

1

I/Ks1

Fig. 9. Variation of Df with I/Ks1 (left) and tf with I/Ks1 (right) for Case 1.

‘‘FD’’ represents a fully drained boundary, ‘‘IMP’’ represents an impermeable boundary, and different Ks1/Ks2 values represent the different partially drained boundaries. It can be observed that for values of I/Ks1 ranging from 0.01 to 0.25, the failures are always deep and occur only for the less permeable boundaries. A value of I/Ks1 = 0.30 gives the rainfall intensity at which shallow failures start to occur, while for I/Ks1 = 0.40 the failures are always shallow irrespective of the boundary conditions. Thus, for the given slope, there exists a critical rainfall intensity (I/Ks1 = 0.40 in this case) beyond which failures are always shallow and completely independent of the nature of the boundary condition. The depth of failure Df decreases non-linearly with an increase in I/Ks1. Recently, it has been shown [11] that failure may occur due to the propagation of a wetting front alone, without any generation of positive pore water pressure, for cases where the slope angle (a) is greater than the friction angle (/0 ) of the soil. The failure depth Df then depends on the residual negative pore water pressure. As the rainfall intensity I/Ks1 increases, a greater reduction in negative pore water pressure occurs which explains why shallow failures occur at higher rainfall intensities. It can also be seen that the failure time tf decreases exponentially with an increase in I/Ks1, but only up to

a critical value of I/Ks1 after which the failure time is independent of the boundary conditions. For rainfall intensities less than the critical I/Ks1, the failure time decreases as the boundary becomes less permeable. Fig. 10 shows the variation in the failure depth Df and failure time tf with different boundary conditions. For the purposes of plotting, a fully drained boundary and an impermeable boundary are represented by Ks1/Ks2 values of 0.001 and 10,000 respectively. It can be observed from the failure time tf that, for a rainfall intensity of I/Ks1 < 0.25, failure occurs only when the boundary is less permeable. For I/Ks1 > 0.4, the failure depths and times are the same for all the boundary conditions considered. As I/Ks1 increases from 0.25 to 0.4, the failures move gradually from being deep to being shallow. Here the effect of the boundary conditions changes from controlling the occurrence of failure to controlling the depth of failure. From the above discussions it can be deduced that, for more permeable boundary conditions, only shallow failures can occur. As the rainfall intensity increases beyond a critical value, the failures become shallower and completely independent of the boundary conditions.

349

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

0

1.00

0.01

0

I/Ks1 increasing

10

Df (m)

-0.4

-0.6

tf days

I/Ks1 increasing

-0.2

Critical I/Ks1 0.40

0.35

10

-1

0.25

0.30 0.35 0.40

10

-0.8

0.04

0.05

-2

0.30

1.00

-1 -3

0

10

4

10

10

10

-3

10

Ks1/Ks2

0

10

4

Ks1/Ks2

Fig. 10. Variation of Df with boundary conditions (left) and tf with boundary conditions (right) for Case 1.

4.2. Case 2: slope with a = 36° and c0 = 1 kPa

4.3. Case 3: slope with a = 40° and c0 = 0 kPa

To investigate the effect of the boundary conditions on a slope having an identical inclination as before (a = 36°), but considering a cohesive soil (c0 = 1 kPa), we again predicted the behavior for the same 24 rainfall intensities as in Case 1. Fig. 11 shows the variation of the failure depths Df and the failure times tf with I/Ks1. It can be observed that the failures are always deep and there are no shallow failures. The failure time tf still decreases exponentially with I/Ks1, just like in Case 1. Fig. 12 shows the variation of the failure depths Df and the failure times tf with different boundary conditions for Case 2. In these examples, failure occurs only for the less permeable boundaries. Thus, the specific type of boundary condition determines whether failure occurs or not. Observation of pore water pressures at failure leads to the conclusion that failure occurs due to generation of a positive pore water pressure only. This case shows that a change in the failure mechanism (from a reduction in negative pore water pressure to generation of a positive pore water pressure) can significantly affect the depth of failure.

This case is considered to check the effect of boundary conditions on the stability of a steeper slope (a = 40°) than considered in Case 1 for cohesionless soil (c0 = 0 kPa). The variation of the failure depths Df and failure times tf with I/Ks1 are shown in Fig. 13. Compared to Case 1, the onset of shallow failures is much earlier in terms of I/Ks1 i.e. the critical I/Ks1 beyond which the shallow failure occurs is much lower (=0.0625) than the critical value of I/Ks1 for Case 1 (=0.40). The variation of the failure depths Df and failure times tf with the boundary conditions for Case 3 are shown in Fig. 14. As long as I/Ks1 < 0.03, the failures are always deep. Here the boundary conditions control the occurrence of failure, with a value of I/Ks1 = 0.04 signifying the onset of shallow failures. Thus, in this situation, the boundary conditions clearly control the depth of failure, except when I/Ks1 > 0.0625 where the failures are always shallow. This case shows that increasing the slope angle results in the failures becoming shallower, thereby reducing the effect of the boundary condition on the failure depths.

Ks1/Ks2 FD

0.01

0.1

1

10

100

1000

IMP

1.5

0

0

10

1

tf (days)

-0.4

Df (m)

tf (days)

-0.2

Failures due to boundary effect only

-0.6

-1

10

0.5

-0.8

0

0.2

0.4

0.6

0.8

1

I/Ks1

-1

0

0.2

0.4

0.6

I/Ks1

0.8

1

0 0

0.2

0.4

0.6

I/Ks1

Fig. 11. Variation of Df with I/Ks1 (left) and tf with I/Ks1 (right) for Case 2.

0.8

1

350

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

0

0.01 s1

increasing

-0.2

t (days)

-0.4 Failures due to boundary effect only

-0.6

0.125

0.10

-1

10

0.30 0.40

f

Df (m)

0.02

I/K

0

10

1.00

-2

10

-0.8

-1 -3 10

10

0

10

4

-3

0

10

4

10

Ks1/Ks2

10

Ks1/Ks2

Fig. 12. Variation of Df with boundary conditions (left) and tf with boundary conditions (right) for Case 2.

Ks1/Ks2 FD

0.01

0.1

1

10

100

1000

IMP

1.5

0

0

10

tf (days)

-0.2

-0.6

Df (m)

Df (m)

-0.4

-0.6

tf (days)

1 Boundary effect No boundary effect

-0.8

-1

10

-2

10

0.5

-1

0

Critical I/Ks1

-0.8

0.2

0.4

0.6

0.8

1

I/Ks1

0.03 0.04 0.05 0.06 0.07 0.08 I/Ks1

-1 0

0.2

0.4

0.6

0.8

0

1

0

0.2

0.4

I/Ks1

0.6

0.8

1

I/Ks1

Fig. 13. Variation of Df with I/Ks1 (left) and tf with I/Ks1 (right) for Case 3.

0

1.00

10

0

0.01 0.02

0.40

-0.2

tf (days)

I/Ks1 increasing

Df (m)

-0.6

Critical I/Ks1

10

I/Ks1 increasing

-0.4

0.03 0.04 0.05 0.0625

0.30

-1

0.05

-0.8

0.30 0.40

0.625

10

0.04

-2

1.00

-1 -3 10

0.03

0

10

Ks1/Ks2

4

10

-3

10

0

10

Ks1/Ks2

Fig. 14. Variation of Df with boundary conditions (left) and tf with boundary conditions (right) for Case 3.

4

10

351

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354 Ks1/Ks2 FD

0.01

0.1

1

10

-0.8

1000

IMP

0

10

Boundary effect

-0.85

-0.95

1

Critical I/Ks1

-1 0.2

0.25

0.3

0.35

0.4

I/Ks1

-0.6

tf (days)

-0.4

No boundary effect

-0.9

tf (days)

Df (m)

-0.2

Df (m)

100

1.5

0

-0.8

-1

10

0.5

0

0.2

0.4

0.6

0.8

1

I/Ks1

-1

0 0

0.2

0.4

0.6

0.8

1

0

0.2

0.4

I/Ks1

0.6

0.8

1

I/Ks1

Fig. 15. Variation of Df with I/Ks1 (left) and tf with I/Ks1 (right) for Case 4.

4.4. Case 4: slope with a = 40° and c0 = 1 kPa

uw ¼ Fig. 15 shows that in this case the variation of Df and tf is quite different compared to any of the cases considered previously. Shallow failures start to occur after I/Ks1 exceeds 0.30. However, unlike Case 1 and Case 3, where failure occur at the surface (Df = 0.01 m) when I/Ks1 = 1, in this case Df reduces only up to 0.61 m for I/Ks1 = 1. Fig. 16 shows the variation of Df and tf with different boundary conditions. A value of I/Ks1 = 0.20 signifies the onset of shallow failures, which always occur once the rainfall intensity exceeds I/Ks1 = 0.30. The critical value of I/Ks1 is higher in this case than for soil without cohesion (Case 3). Thus it can be concluded that an increase in cohesion increases the critical value of I/Ks1 for shallow failures. It can be also said that an increase in cohesion results in deeper failures, as is evident from Case 2 and Case 4. The observations made in the four cases can be explained with the analytically estimated pore water pressure at failure. Assuming that the soil is fully saturated (i.e. Se = 1) at failure, the pore water pressure corresponding to different failure depths can be estimated analytically. Substituting FS = 1 and Se = 1 into Eq. (12) gives:

  /0 c0  1  tan W cos a sin a tan a

ð14Þ

tan /0

The above expression can be used to estimate the pore water pressure at failure, for all the four cases and thus understand the failure mechanism. The initial and failure profiles of the pore water pressure for all the four cases are shown in Fig. 17. The pore water pressure at failure is shown on a magnified scale in the inset. Since hydrostatic initial conditions were considered, the initial profile of pore water pressure (for the four cases) depends only on the inclination of the slope and is more negative for slopes having gentler inclinations (Case 1 and Case 2). The residual negative pore water pressure (at failure) increases with an increase in depth. This implies that the failure depth will depend on the reduction in the negative pore water pressure and will be the deepest when the reduction in the negative pore water pressure is the least. This is because with an increase in depth, the destabilizing force due to the weight component increases and, therefore in case of a deep failure, the negative pore water pressure need not be reduced to the same extent as that required for a shallow failure. This explains

0 0.01

0

10

0.02

-0.8

-1

10

Critical I/Ks1

0.10

0.125 0.25 0.30 0.40

0.200

1.00

-2

0.40

10

0.30 0.25 0.200 0.175

-1 10

tf (days)

1.00

-0.6 I/Ks1 increasing

Df (m)

-0.4

I/Ks1 increasing

-0.2

-3

10

0

Ks1/Ks2

10

4

10

-3

10

0

Ks1/Ks2

Fig. 16. Variation of Df with boundary conditions (left) and tf with boundary conditions (right) for Case 4.

10

4

352

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

0

-0.1 Depth (m)

-0.2

-0.2 -0.3

-0.4 -0.6 -0.8

-0.4

-1 -2

-0.5

-1 0 1 Pore water pressure (kPa)

-0.6 Initial condition

Initial condition

α = 36

-0.8 -0.9 -1

Case 1 Case 2 Case 3 Case 4 -15

-10

-5

0

Pore water pressure (kPa) Fig. 17. Profiles of pore water pressure – initially and at failure for all the cases; inset, profile of pore water pressure at failure on a magnified scale.

0

0 Initial condition for α = 40

-0.4

-0.2

Depth (m)

Depth (m)

-0.2

Initial condition for α = 36

-0.6

-0.8

-1 -20

Case 1 Case 2 Case 3 Case 4 -15

-0.4

-0.6

-0.8

-10

-5

0

-1 0.5

5

Failure

-0.7

α = 40

From Fig. 17 it can be observed for Case 1 and Case 3 that failures at all depths are due to a reduction in the negative pore water pressure. However, due to a greater inclination, a larger residual pore water pressure can be observed for Case 3. This implies that as the slope becomes steeper, failures will be shallower provided that the reduction in the negative pore water pressure is the same compared to a gentler slope. This explains why failures are shallower in Case 3 at the same values of I/Ks1 compared with Case 1, and the onset of shallow failures in Case 3 at lower values of I/Ks1. Comparing the failure profiles of Case 1 with Case 2 (or Case 3 with Case 4) shows that an increase in cohesion not only increases the pore water pressure required for failure but could also change the failure mechanism. Considering the failure profile of Case 2, it can be seen that due to an increase in cohesion, there is a change in the failure mechanism. Failure will occur (at any depth) only when positive pore water pressures is generated in contrast to Case 1

1

Pore water pressure (kPa)

Case 1 Case 2 Case 3 Case 4 1.5

2

2.5

3

Factor of Safety

Fig. 18. Pore water pressure and factor of safety profiles of a fully drained boundary at I/Ks1 = 0.2 for all cases.

0

0 Initial condition for α = 40

-0.4

-0.2

Depth (m)

Depth (m)

-0.2

Initial condition for α = 36

-0.6

-0.8

-1 -20

Case 1 Case 2 Case 3 Case 4 -15

-0.4

-0.6

-0.8

-10

-5

0

Pore water pressure (kPa)

5

-1 0.5

Failure

Depth (m)

i. Why failures are deep, as the reduction in the negative pore water pressure is not sufficient enough to cause a shallow failure. ii. Why with an increase in I/Ks1 ratio, the failure becomes shallow as now a greater reduction in the negative pore water pressure can be achieved with a higher I/Ks1.

0

1

Case 1 Case 2 Case 3 Case 4 1.5

2

2.5

Factor of Safety

Fig. 19. Pore water pressure and factor of safety profiles of an impermeable boundary at I/Ks1 = 0.2 for all cases.

3

353

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

0

2

-0.2

1.2

Critical I/Ks1 (Case 1)

-0.6

tf (days)

-0.4

1.4

1.6 1.4

tf (days)

Df (m)

Critical I/Ks1 (Case 3)

1.6

1.8

1.2 1

-0.8

1 0.8 0.6 0.4

0.8

0.2

0.6 Critical I/Ks1 (Case 4)

Case 1 Case 2 Case 3 Case 4

0 -2 10

0.4

-1

0

10 I/Ks1

10

0.2

-1 0

0.2

0.4

0.6

0.8

1

0

0.2

0.4

0.6

0.8

1

I/Ks1

I/Ks1

Fig. 20. Variation in failure depths Df and failure time tf of all cases with I/Ks1 for an impermeable boundary only.

where failures occurred due to a reduction in the negative pore water pressure. Since the profile is homogeneous, positive pore water pressure is not generated during infiltration, until the wetting front has reached the boundary. This explains why in Case 2 the failures were always deep even at very high rainfall intensities (I/Ks1 = 1) and why failures occurred only when the boundary was less permeable. The failure profile of Case 4 shows that positive pore water pressure is required for failure to occur at depths less than approximately 0.75 m. This explains why in Case 4, the failure depths reduced only up to 0.61 m even at very high rainfall intensities (e.g. I/Ks1 = 1). However, it must be stated that although it is possible for a slope to fail at depths very close to the boundary (say Df > 0.90 m), the failure depth would depend on the nature of the boundary condition as the failure could be shallow for a particular boundary but deep for another boundary type (see Fig. 7). Therefore, the I/Ks1 ratio should not be confused as the rainfall intensity beyond which the failure is shallow for a given boundary; rather, it is the rainfall intensity beyond which the failure will be shallow for all boundary types. Figs. 18 and 19 show the pore water pressure and factor of safety profiles for a fully-drained boundary and an impermeable boundary, respectively, when I/Ks1 = 0.2. Fig. 20 compares the failure depths Df and failure times tf for all the cases with an impermeable boundary. It can be seen that Case 3 is the least affected by the boundary conditions as far as the failure depths Df are concerned. Case 2, on the other hand, is most affected by the boundary conditions since failure always occurs along the boundary. For Case 1 and Case 4, the boundary conditions affect the behavior of the slope as long as the critical rainfall intensity needed to cause shallow failure is not exceeded. Comparing the failure times tf reveals that, as the slope angle increases, the failure time decreases (as the failure depths are now smaller for the same rainfall intensity). It is also evident that as the cohesion increases, the failure time increases (as the failure depths are now greater) provided the rainfall intensity I/Ks1 stays the same.

and purely frictional soils have been considered. From the results obtained, the following conclusions can be drawn:

Overall, it can be concluded that the nature of boundary plays an important role in the triggering of landslides, by either controlling the depth of failure or the overall occurrence of failure. It should be noted that the present study did not consider scale effects and effect of different initial conditions. These factors may have a significant influence for some cases.

5. Conclusion

Acknowledgements

The effects of various boundary conditions on the failure depths and failure times of a soil slope have been studied for a variety of rainfall intensities I/Ks1 and slope angles. Slopes in both cohesive

The authors wish to acknowledge the support of the Australian Research Council in funding the Centre of Excellence for Geotechnical Science and Engineering.

i. For more permeable boundaries (i.e. fully drained or partially drained boundaries where Ks1/Ks2 < 1), the failures are always shallow. Deep failures only occur for less permeable boundaries (i.e. impermeable or partially drained boundaries where Ks1/Ks2 > 1). This is true for all the four cases that were investigated. ii. When failure occurs due to a reduction in negative pore water pressure there exists a critical rainfall intensity (I/Ks1 value), beyond which the failure is always shallow and independent of the boundary conditions. In this case, the failure depths decrease non-linearly as I/Ks1 increases and the shallowest failure occurs for the highest rainfall intensity. The critical value of I/Ks1 depends on the slope angle and material properties. As the slope inclination increases or the cohesion decreases, the critical value of I/Ks1 decreases. Lastly, once the rainfall intensity falls below the critical rainfall intensity, the boundary conditions effects become increasingly noticeable as the failures occur only at the boundary. iii. The effect of the boundary conditions is more pronounced in gentler slopes than steeper slopes. In the former, the failures occur at a greater depth. Also, the boundary conditions have a more significant effect on the failure depths when the failures occur because of the generation of positive pore water pressure rather than when they occur as a result of a reduction in negative pore water pressure.

354

A. Ali et al. / Computers and Geotechnics 61 (2014) 341–354

References [1] Bishop AW. The effective stress principle. Teknisk Ukeblad 1959;39:859–63. [2] Brand EW. Landslides in South Asia: a state-of-art report. In: 4th International symposium on landslides, Toronto; 1984. [3] Cheng YG. Numerical simulation of unsaturated flow using modified transformation methods. In: Civil Engineering. National University of Singapore; 2008. p. 259. [4] Cho SE, Lee SR. Evaluation of surficial stability for homogeneous slopes considering rainfall characteristics. J Geotechn Geoenviron Eng 2002;128(9):756–63. [5] Collins BD, Znidarcic D. Stability analyses of rainfall induced landslides. J Geotechn Geoenviron Eng 2004;130(4):362–72. [6] Fourie AB. Predicting rainfall-induced slope instability. Proc ICE-Geotechn Eng 1996;119(4):211–8. [7] Fredlund DG, Rahardjo H. Soil mechanics for unsaturated soils. New York: John Wiley & Sons Inc.; 1993. [8] Iverson RM. Landslide triggering by rain infiltration. Water Resour Res 2000;36(7):1897–910. [9] Kasim FB, Fredlund DG, Gan JKM. Effect of steady state rainfall on long term matric suction conditions on slopes. In: Proceedings of the second international conference on unsaturated soils. Beijing, China; 1998. [10] Lee LM, Kassim A, Gofar N. Performances of two instrumented laboratory models for the study of rainfall infiltration into unsaturated soils. Eng Geol 2011;117(1–2):78–89. [11] Li WC, Lee LM, Cai H, Li HJ, Dai FC, Wang ML. Combined roles of saturated permeability and rainfall characteristics on surficial failure of homogeneous soil slope. Eng Geol 2013;153:105–13. [12] Lu N, Godt JW. Infinite slope stability under steady unsaturated seepage conditions. Water Resour Res 2008;44(11):W11404. [13] Ng CWW, Shi Q. A numerical investigation of the stability of unsaturated soil slopes subjected to transient seepage. Comput Geotech 1998;22(1):1–28. [14] Rahardjo H, Lim TT, Chang MF, Fredlund DG. Shear-strength characteristics of a residual soil. Can Geotechn J 1995;32(1):60–77.

[15] Richards LA. Capillary conduction of liquids through porous mediums. Physics 1931;1(5):318–33. [16] Simunek J, van Genuchten MT, Sejna M. The Hydrus-1D software package for simulating the movement of water, heat, and multiple solutes in variably saturated media, Version 4.16, HYDRUS Software Series 3. Riverside, California, USA: Department of Environmental Sciences, University of California Riverside; 2013. p. 340. [17] Tsai TL. The influence of rainstorm pattern on shallow landslide. Environ Geol 2008;53(7):1563–9. [18] Tsai TL. Influences of soil water characteristic curve on rainfall-induced shallow landslides. Environ Earth Sci 2011;64(2):449–59. [19] Tsai TL, Chen HF. Effects of degree of saturation on shallow landslides triggered by rainfall. Environ Earth Sci 2010;59(6):1285–95. [20] Tsai TL, Chiang SJ. Modeling of layered infinite slope failure triggered by rainfall. Environ Earth Sci 2013;68(5):1429–34. [21] Tsai TL, Yang JC. Modeling of rainfall-triggered shallow landslide. Environ Geol 2006;50(4):525–34. [22] Tsaparas I, Rahardjo H, Toll DG, Leong EC. Controlling parameters for rainfallinduced landslides. Comput Geotechn 2002;29(1):1–27. [23] van Genuchten MT. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci Soc Am J 1980;44(5):892–8. [24] Vanapalli SK, Fredlund DG, Pufahl DE, Clifton AW. Model for the prediction of shear strength with respect to soil suction. Can Geotechn J 1996;33(3):379–92. [25] White JA, Singham DI. Slope stability assessment using stochastic rainfall simulation. Procedia Comput Sci 2012;9:699–706. [26] Zhan TLT, Jia GW, Chen YM, Fredlund DG, Li H. An analytical solution for rainfall infiltration into an unsaturated infinite slope and its application to slope stability analysis. Int J Numer Anal Methods Geomech 2012;37(12):1737–60. [27] Zhan TLT, Ng CWW. Analytical analysis of rainfall infiltration mechanism in unsaturated soils. Int J Geomech 2004;4(4):273–84.