Reliability analysis of tunnel roof in layered Hoek-Brown rock masses

Computers and Geotechnics xxx (xxxx) xxx–xxx

Contents lists available at ScienceDirect

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

Research Paper

Reliability analysis of tunnel roof in layered Hoek-Brown rock masses ⁎

X.L. Yanga, , T. Zhoua, W.T. Lib a b

School of Civil Engineering, Central South University, Hunan 410075, China School of Civil and Environmental Engineering, Nanyang Technological University, Singapore 639798, Singapore

A R T I C L E I N F O

A B S T R A C T

Keywords: Reliability analysis Deep tunnels Collapse mechanism Supporting pressure

This paper develops a two-layer model and upper bound solutions for the shape of collapse block in rectangular tunnel subjected to seepage pressure, and then the reliability-based analysis is performed. According to the layered Hoek-Brown rock masses with different characteristics, a continuous failure mechanism of roof collapse in the deep cavities is derived by using upper bound solutions. The performance function of roof stability is proposed and is adopted for reliability analysis and design. Parameters involved in the performance function are random or deterministic variables. For comparison, improved response surface method (RSM) and Monte-Carlo simulation are utilized to calculate the Hasofer-Lind reliability index and the failure probability. It is found that the improved RSM is a reliable method in probabilistic calculation and the supporting pressure has a significant influence on the reliability index. The normal or lognormal distribution of variables has no effect on the failure probability, while the coefficient of variation (COV) of random variables would greatly influence the failure probability. A reliability-based design is performed to determine the probabilistic tunnel pressure for a target reliability index when different COVs of random variables are considered.

1. Introduction

However, only few publications employed the limit analysis and the reliability methods to evaluate the stability of deep tunnels. One of the difficulties lies in the construction of rational failure mechanism for deep tunnels. A curved failure mechanism on the crown of deep tunnel was proposed by Fraldi and Guarracino [19–22]. In their research, the shape and the region of potential collapse of deep tunnel or cavity can be determined by virtue of upper bound theorem. In tunnel engineering, ground water is an adverse factor which has negative influence on the stability of underground structures. With reference to the above works, pore water pressure was incorporated into the stability analysis of tunnel roof [23–25]. Pore water pressure was introduced into the limit analysis by regarding it as an external force acting on soil skeleton, and the results were in good agreement with previous findings. Based on previous research, Li and Yang [26] evaluated the effect of variable detaching velocity along yield surface in the roof failure mechanism, and the existing failure mode was developed to karst area with Hoek-Brown criterion by Yang et al. [27]. Considering these mechanisms, the tunnel behavior in single stratum was characterized. As mentioned above, most of the existing research focused on the homogeneous and isotropic rock masses in the stability analysis of tunnel engineering. Since deep tunnels or cavities are often constructed in layered rock strata, more attention should be paid to the stability problem of such projects as layered rock masses usually are characterized by different mechanical properties in engineering practice.

In the last decades, deterministic methods were commonly used in the stability analysis of tunnels [1]. The mean values of geotechnical parameters are usually considered in these approaches and a global safety factor based on experience is used to account for the uncertainties of geomaterials. However, this factor cannot reflect the inherent uncertainty of each parameter [2]. With further research on statistic properties of geomaterials [3,4], reliability-based analyses and designs were extensively performed and developed in stability assessment of underground structures [5–13]. For the application of reliability methods to a certain case or problem, the probabilistic analysis function should be derived first. Several researchers attempted to combine reliability methods with limit analysis for their advantages [14,15]. The computation efficiency of limit analysis paves the way for building failure mode and characterizes the engineering behavior of a certain structure, while reliability methods provide practical evaluation by considering the inherent uncertainties of geo-parameters. Massih et al. [16] conducted the reliability-based analysis of strip footings. In this research, a failure mode was derived based on upper-bound limit method, and Hasofer-Lind reliability index was used as an indicator for probabilistic assessment. Based on limit analysis theory, Mollon et al. [17,18] proposed a multi-failure mode of tunnel face and utilized the reliability method to gain further insight into the problem.

⁎

Corresponding author. E-mail addresses: [email protected] (X.L. Yang), [email protected] (T. Zhou), [email protected] (W.T. Li).

https://doi.org/10.1016/j.compgeo.2017.12.007 Received 7 June 2017; Received in revised form 29 November 2017; Accepted 15 December 2017 0266-352X/ © 2017 Elsevier Ltd. All rights reserved.

Please cite this article as: Yang, X.L., Computers and Geotechnics (2017), https://doi.org/10.1016/j.compgeo.2017.12.007

Computers and Geotechnics xxx (xxxx) xxx–xxx

X.L. Yang et al.

variable, and ai is the coefficient to be determined. In order to obtain a tentative response surface, a new set of sampling points is selected by using the following vector projection technique. FORM is used to obtain reliability index β and the corresponding design point x* by using the tentative response surface function. The initial values of coefficients a0,a1,…ai,…an are determined by using sampling points centered at the initial central point x i0 and located in each direction at

Seepage force is another key factor which poses a great challenge to tunnel stability, and hence its effects on tunnel stability should not be ignored. Under seepage action, how to develop the potential collapse mode of a deep tunnel in layered rock mass is a new issue, and the determination of required supporting pressure is of great importance. In this paper, reliability analysis is performed to evaluate roof stability of deep tunnels in two-layer rock masses. It is assumed that the potential collapse blocks pass through two layers. The potential failure mechanism subjected to seepage pressure was established and the performance function of roof stability was then derived. Due to the implicit form of performance function, improved response surface method (RSM) and Monte-Carlo simulation (MCS) were used in reliability analysis. The sensitivity of the failure probability to the distribution type and the coefficient of variation (COV) of random variables were studied. Then, based on the reliability–based design (RBD), the required supporting pressures under different situations were determined with respect to a given target reliability index.

x i0 ± hσxi, i = 1,2,…,n

(4)

where h is generally selected from 1 to 2 in the paper. The initial central points can be the mean value points. It may also be shifted towards the failure domain to expedite convergence. 2.3. MCS MCS is a well-accepted reliability method in which samples are generated subjected to the probability density of random variables. The failure probability of structure is estimated as

2. Reliability concepts

Pf = 2.1. Hasofer-Lind index

1 N

n

∑

I (x i )

i=1

(5)

The reliability index has been widely used in the evaluation of structural safety. Hasofer and Lind [28] proposed their reliability index β and the matrix formulation for normal variables is given by

where N is the number of samples. I (x ) = 1 if g (x ) ⩽ 0 and 0 elsewhere. With respect to the law of large numbers, the accuracy of MCS lies in the large number of samples and trails. The convergence of the failure probability is estimated by its COV

β = minx ∈ F (x −μ )T C −1 (x −μ )

COV (Pf ) =

(1)

where x is the vector of n random variables, μ means the vector of their mean values, C is the covariance matrix, and F presents the failure region represented by the performance function g (x ) ⩽ 0 . According to Eq. (1), the Hasofer-Lind reliability index means the minimum distance in units of directional standard deviations from the mean value of random variables to the limit state surface g (x ) = 0 in units of directional standard deviations. Since the reliability index β considers not only the influences of mean values but also the covariance of relevant variables [29], it has been widely adopted in reliability analysis as shown in Refs. [6,8–10,16–18]. Apart from normal variables, non-normal variables are also used to characterize the property of geotechnical materials or loads in engineering practice. To utilize non-normal variables in the calculation of reliability index, Rackwitz–Fiessler transformation [30] can be used to obtain the equivalent normal mean value μiN and the equivalent normal standard deviation σiN . Failure probability is commonly used to evaluate the safety of engineering structures in another aspect. Based on the first-order reliability method (FORM) and Hasofer-Lind reliability index β , the probability of failure Pf can be calculated from

Pf ≅ 1−Φ(β )

(1−Pf )/(Pf N )

3. Limit analysis and Hoek-Brown criterion 3.1. Upper bound theorem According to the upper bound theorem, the actual collapse load is not greater than the load derived by equating the rate of the energy dissipation to the external rate of work in any kinematically admissible velocity field, when the deformation boundary conditions are satisfied [33]. It takes the form as

∫V σij·εij̇ dv ⩾ ∫S Ti·vi ds + ∫V Xi ·vi dv + ∫V −grad u·vi dV

(7)

where σij and εij̇ are the stress tensor and the rate of strain in a kinematically admissible velocity field respectively, Ti and Xi are the surface force and the body force of the studied object respectively, V is the volume of the collapse block, S is the length of velocity discontinuity, −gradu is excess pore pressure, and vi stands for the velocity along the detaching surface. 3.2. Hoek-Brown criterion

(2) Because of the nonlinear characteristics of geomaterials, HoekBrown failure criterion is widely employed to investigate the nonlinear problems in geotechnical engineering. The relationship between normal and shear stresses of Hoek-Brown criterion is

where Φ(·) represents the cumulative distribution function (CDF) of the standard normal variable. 2.2. Improved RSM

τ = Aσc [(σn−σt )/ σc ]B In reliability analysis, an explicit expression of performance function is required. However, explicit performance function is often unlikely to be determined in complex engineering problems. To solve this problem, RSM is proposed to approximate the actual performance function with the use of sampling points. Improved RSM is suggested [31,32] to determine the proper location of sampling points. In their research, the algorithm is presented as

g (X ) ≈ ∼ g (X ) = a0 +

(6)

where A and B are material constants respectively, σn is the normal effective stress, σc and σt mean the uniaxial compressive strength and tensile strength of the rock mass, respectively. The Hoek-Brown criterion is widely used in rock engineering [34–37]. 4. Collapse mechanism and performance function

n

∑ i=1

ai x i

(8)

4.1. Failure mechanism (3)

According to former studies, the velocity discontinuity curves, f (x ) and g (x ) , are introduced to describe the shape of the collapse block in

where x i is the basic random variable, n is the number of the random 2

Computers and Geotechnics xxx (xxxx) xxx–xxx

X.L. Yang et al.

PD = =

∫0

L1

Di̇ 1 1 + f ′ (x )2 dx +

∫L

L2

1

1 1 1 1 − B1

Di̇ 2 1 + g ′ (x )2 dx 1 1 − B1

{−σ + σ [A B ] (1−B ) f ′ (x ) } vdx (1−B )g′(x ) + ∫ {−σ + σ [A B ] } vdx ∫0

L1

t1

c1

L2

t2

L1

2

c2

−1 1

1 2 1 − B2

−1 2

1 1 − B2

(15)

in which L1 and L2 are the widths of two curves illustrated in Fig. 1. The work rate of failure block created by gravity can be expressed as

Pγ =

∫0

L1

ρ1′ [f (x )−f (L1)] vdx +

∫0

L1

ρ2′ f (L1) vdx +

∫L

L2

1

ρ2′ g (x ) vdx (16)

where ρ1′ and ρ2′ are the buoyant weight per unit volume of the upper and lower rock layers calculated by ρ′ = ρ−ρw , in which ρ is the weight per unit volume of the rock and ρw represents the unit weight of water. With reference to the research of Saada et al. [38], the distribution of excess pore pressure in rock mass is given as

u = p−pw = p−ρw h Fig. 1. Impending collapse blocks in rock layers.

arbitrary cross-section tunnel, as shown in Fig. 1. The plastic flow only occurs along the velocity discontinuity surface, and the energy dissipation rate along the surface is calculated with the associated flow rule. It is found that the failure mechanism is symmetrical with the respect of Y-axis. Since the geomaterials are assumed to obey the associated flow rule, the plastic potential functions Ω based on HoekBrown failure criterion of these two stratums are expressed as

Ω1 = τn1−A1 σc1 [(σn1 + σt1)/ σc1 ]B1

(17)

where p means the pore pressure defined by p = ru ρw h , ru means the pore pressure coefficient and h is the vertical distance between the tunnel roof and the top of the collapse block. The gradient of pore pressure is.

−grad u = −

du = ρw −ru ρ dy

(18)

Therefore, the work rate of seepage force along the detaching curve is

Pu =

(9)

L1

∫0 (ρw −ru ρ1)[f (x )−f (L1)] vdx + ∫0 L + ∫ (ρw −ru ρ2 ) g (x ) vdx L

L1

(ρw −ru ρ2 ) f (L1) vdx

2

Ω2 = τn2−A2 σc 2 [(σn2 + σt 2)/ σc 2 ]B2

During the construction of tunnels, the supporting structure is necessary for tunnel safety and stability. The supporting pressure can be considered as an external force in limit analysis, and its work rate is

where A and B are material constants in H-B criterion, and parameters with the subscripts 1 or 2 represent the upper and lower rock layers, respectively. Making reference to Fig. 1, the normal stress σn of any points along the velocity discontinuity are derived as 1

σn1 = −σt1 + σc1 [A1 B1 f ′ (x )]1 − B1

Pq = L2 σq v cosπ

1

(12)

where f ′ (x ) and g ′ (x ) are the first derivative of the collapse curve f (x ) and g (x ) , respectively. Therefore, the energy dissipation rate of a random point on the detaching surface can be written as

ζ [f (x ),f ′ (x ),x ] = PD−Pγ −Pu−Pq = +

{

1

}

L1

ψ1 [f (x ),f ′ (x ),x ] vdx

ψ2 [g (x ),g ′ (x ),x ] vdx + W

(21)

1

1

ψ1 [f (x ),f ′ (x ),x ] = −σt1 + σc1 [A1 B1]1 − B1 (1−B1−1) f ′ (x ) 1 − B1 −(1−ru ) ρ1 f (x ) (13)

(22)

ψ2 [g (x ),g ′ (x ),x ] = −σt2 +

{

}

∫0

where

1 1 v Di̇ 2 = σn2 εṅ 2 + τn2 γṅ 2 = [1 + g ′ (x )2]− 2 −σt2 + σc2 [A2 B2]1 − B2 (1 w 1 −B2−1) g ′ (x ) 1 − B2

∫L

L2

1

1 1 v Di̇ 1 = σn1 εṅ 1 + τn1 γṅ 1 = [1 + f ′ (x )2]− 2 −σt1 + σc1 [A1 B1 ]1 − B1 (1 w

−B1−1) f ′ (x ) 1 − B1

(20)

where σq is the supporting pressure applied on the tunnel boundary. According to the upper bound theorem, the effective shape of collapse block in a limit state can be determined by minimizing the difference of the total energy dissipation rate and the rate of external work. Based on the virtual work equation, an objective function which characterizes the difference of the rate of energy dissipation and the whole external rate of work is expressed as

(11)

σn2 = −σt 2 + σc 2 [A2 B2 g ′ (x )]1 − B2

(19)

1

(10)

1 σc2 [A2 B2]1 − B2

1 (1−B2−1)g′(x ) 1 − B2 −(1−ru ) ρ2 g (x )

(23)

W = L2 σq v−(1−ru )(ρ2 −ρ1) f (L1) v

(14)

(24)

According to the upper bound theorem, the effective shape of collapse block in a limit state can be determined by minimizing the objective function. By transforming the expression of ψ (x ) into Euler’s equation and integrating the results, the variational equations are expressed as

̇ are normal plastic strain rates respectively, γ̇n1 and γ̇n2 ̇ and εn2 where εn1 are shear plastic strain rates respectively, w is the thickness of the plastic detaching zone, and v is the velocity of the collapse block. With respect to the symmetry of the failure mechanism, the whole rate of energy dissipation along the detaching curve is integrated over two different layers which are characterized by diverse parameters. It is given as

∂ψ1 ∂ψ1 ⎤ ∂ − ⎡ =0 ∂f (x ) ∂x ⎢ ∂ f ⎣ ′ (x ) ⎥ ⎦ 3

(25)

Computers and Geotechnics xxx (xxxx) xxx–xxx

X.L. Yang et al. − 1 B1

∂ψ2 ∂ψ2 ⎤ ∂ − ⎡ =0 ∂g(x ) ∂x ⎢ ⎦ ⎣ ∂g ′ (x ) ⎥

[(1−ru ) ρ1 (h + h 0)−σt1] L1−A1

(26)

B −1

B1+ 1 B1

1 1 1 [(1−ru ) ρ1 ] B1 σc1 B1 L1 1 + B1

1

+ (1−ru )(ρ2 −ρ1) hL1 + L2 σq + ⎡(1−ru ) ρ2 (L2 + z ) B2 −σt2⎤ (L2−L1) ⎣ ⎦

By substituting Eqs. (22) and (23) into Eq. (25) and (26) separately, and integrating the result, the analytical expressions of the detaching curves f (x ) and g (x ) are obtained − 1 B1

f (x ) = A1

− 1 B2

g (x ) = A2

− 1 B2

−A2

1

[(1−ru ) ρ1 / σc1 ]

1 − B1 B1 ⎡x

⎢ ⎣

+

n 0 ⎤B1 −n1 (1−ru ) ρ1 ⎥ ⎦

(27)

Therefore, by combining Eq. (30), (31) and (39), the values of L1, L2 and h can be obtained by virtue of numerical tools. Different rock parameters would lead to the different collapse shapes of rock mass. Based on the values of L1, L2 and h , the expressions of detaching curves f (x ) and g (x ) can be determined.

1

[(1−ru ) ρ2 / σc2 ]

1 − B2 B2

B n2 ⎡x + ⎤ 2 −n 3 ⎢ ⎥ (1−ru ) ρ2 ⎦ ⎣

(28)

where n 0 , n1, n2 and n3 are integration constants respectively, which can be determined by boundary condition. Since the detaching curve f (x ) is symmetrical with respect to the y-axis, the integration constant n 0 is equal to zero. As shown in Fig. 1, the following geometry conditions should be satisfied

g (L2) = 0

(29)

f (L1) = −h 0

(30)

g (L1) = −h 0

(31)

4.2. Performance function This paper aims at performing reliability analysis for deep tunnels or cavities constructed in laminar geomaterials subjected to seepage action. Based on the roof failure mechanism presented above, the shape of collapse blocks can be obtained with reference to the upper bound theorem and the variational principle. As discussed above, the dimension of detaching blocks will change when different supporting pressures and geotechnical parameters are considered. Different scopes of collapse blocks would result in different collapse pressures. To maintain the stability of tunnel roof, the supporting pressure should be no less than the weight of collapse rock mass. Therefore, the performance function used in reliability analysis is proposed with respect to the potential roof collapse of deep tunnels.

To maintain the smoothness of detaching curve, an equation between the first derivations of f (x ) and g (x ) at a certain point should be satisfied.

f ′ (x = L1) = g ′ (x = L1)

(32)

G (x ) = σq L2−Wrock

Though computation and simplification, f (x ) and g (x ) are written 1

(40)

in which

as

f (x ) = k1 x B1 −(h + h 0)

B −1

2 B2 + 1 B2 + 1 1 1 [(1−ru ) ρ2 ] B2 σc2 B2 ⎡ (L2 + z ) B2 −(L1 + z ) B2 ⎤ = 0 1 + B2 ⎣ ⎦ (39)

Wrock =

(33)

∫0

L1

ρ1 [f (x )−f (L1)] vdx +

∫0

L1

ρ2 f (L1) vdx +

∫L

L2

1

ρ2 g (x ) vdx (41)

g (x ) = k2 (x +

1 z ) B2 −k2 (L2

+

1 z ) B2

(34)

The expression of the performance function is proposed as B1+ 1 B1

in which h and h 0 are shown in Fig. 1

k1 =

k2 =

1 − B1 − 1 A1 B1 [(1−ru ) ρ1 / σc1 ] B1

− 1 A2 B2

(35)

1 − B2 [(1−ru ) ρ2 / σc2 ] B2 B2 1 − B2

kB z = ⎡ 1 2⎤ ⎢ B2 k1 ⎦ ⎥ ⎣

1 1 − B1 B1 ⎡ − G (x ) = σq L2−ρ1 ⎢A1 B1 [(1−ru ) ρ1 / σc1 ] B1 L1 1 + B1 ⎣ − 1 B2

−ρ2 ⎧A2 ⎨ ⎩

(36)

∫0

L1

ψ1 [f (x ),f ′ (x ),x ] vdx +

∫L

L2

1

ψ2 [g (x ),g ′ (x ),x ] vdx

5. Reliability analysis of deep tunnels

⎧ + W = v [(1−ru ) ρ1 (h + h 0)−σt1 ] L1 ⎨ ⎩ − 1 B1

−A1

B −1

Before the employment of reliability analysis methods, parameters involved in the performance function G (x ) should be identified as deterministic variables or random variables first. In this paper, the uniaxial compressive strength σci , the tensile strength σti , the rock weight ρi , the supporting pressure σq as well as the pore pressure coefficient ru are characterized as random variables with mean values and coefficient of variations. The subscript i in these parameters represents 1 or 2 separately. Their statistic values and distribution types are listed in Table 1. Apart for them, Ai , Bi and h 0 (as defined in Fig. 1) are regarded as deterministic parameters. Their statistic data are shown in Table 2. Generally, the normal distribution of random variables is commonly utilized to reflect the uncertainty of geometry materials. However, the lognormal distribution is always recommended in the reliability analysis to avoid negative values when the coefficient of variables (COV) is

B1+ 1 B1

1 1 1 [(1−ru ) ρ1 ] B1 σc1 B1 L1 1 + B1

+ (1−ru )(ρ2 −ρ1) hL1 + L2 σq 1

+ ⎡(1−ru ) ρ2 (L2 + z ) B2 −σt2⎤ (L2−L1) ⎣ ⎦ − 1 B2

−A2

B −1

2 B2 + 1 1 1 [(1−ru ) ρ2 ] B2 σc2 B2 ⎡ (L2 + z ) B2 1 + B2 ⎣

−(L1 + z )

B2 + 1 ⎫ B2 ⎤

⎦⎬ ⎭

(42)

Since the analytical solutions of L1 L2 and h are difficult to find, the performance function of the rectangular tunnel cannot be expressed explicitly. Reliability analysis of roof stability will be presented in the next section.

(37)

Then the objective function ζ is obtained as

ζ [f (x ),f ′ (x ),x ] =

1 − B2 B2 + 1 B2 + 1 B2 [(1−ru ) ρ2 / σc2 ] B2 ⎡ (L2 + z ) B2 −(L1 + z ) B2 ⎤ 1 + B2 ⎣ ⎦

1 −k2 (L2 + z) B2 (L2−L1)⎫−(ρ2 −ρ1) L1 h 0 ⎬ ⎭

(1 − B1) B2 L1B1 (1 − B2) −L1

⎤ −HL1⎥ ⎦

(38)

By equating the rate of internal energy dissipation to the work rate of external load, the optimized upper bound solution can be achieved. 4

Computers and Geotechnics xxx (xxxx) xxx–xxx

X.L. Yang et al.

Table 1 Statistical values of random variables in analysis.

Table 3 Reliability indexes and design points for normal and lognormal variables.

ρ1 (kN/ m3)

∗ σc2 (MPa)

∗ (kPa) σt2

ρ2 (kN/ m3)

(a) Normal variables 200 0.133 9.96 210 0.512 9.88 220 0.866 9.80 230 1.198 9.72 240 1.509 9.66 250 1.802 9.59 260 2.078 9.45

150.51 151.87 153.05 154.08 154.98 155.76 156.44

21.09 21.35 21.58 21.78 21.97 22.13 22.27

15.07 15.28 15.45 15.61 15.76 15.88 15.99

200.11 200.39 200.63 200.84 201.01 201.17 201.30

23.00 23.01 23.02 23.02 23.03 23.03 23.03

(b) Lognormal variables 200 0.135 9.97 210 0.519 9.88 220 0.884 9.80 230 1.233 9.72 240 1.566 9.64 250 1.884 9.57 260 2.190 9.50

150.51 151.97 153.35 154.64 155.85 156.98 158.03

21.10 21.38 21.66 21.93 22.21 22.49 22.76

15.08 15.30 15.51 15.73 15.94 16.15 16.36

200.11 200.40 200.67 200.91 201.12 201.32 201.49

23.00 23.01 23.02 23.02 23.03 23.04 23.05

Mean value

Coefficient of variation

Distribution type.

σq (kPa)

σc1 (kPa) σt1 (kPa) ρ1 (kN/m3)

10,000 150 2100

0.1 0.1 0.1

Normal/Lognormal Normal/Lognormal Normal/Lognormal

σc2 (kPa)

15,000

0.1

σt2 (kPa)

200 2300

0.1 0.1

Lognormal distribution Normal/Lognormal Normal/Lognormal

0.2

0.1

Normal/Lognormal

ρ2 (kN/m3) ru

Table 2 The deterministic parameters used in analysis. h0 (m)

A1

B1

A2

B2

2

0.4

0.8

0.65

0.7

∗ σc1 (MPa)

∗ (kPa) σt1

Random variable

β

types increases with the increment of supporting pressure. Thus, it can be concluded that the hypothesis of normal parameters is conservative in comparison with the one of lognormal parameters. The values of design points (σc∗i , σti∗, ρi∗, and ru∗) corresponding to different supporting pressures represent the most probable failure point on the limit state surface (LSS). It is the point where the expanding 7dimensional dispersion ellipsoid is tangent to the limit state surface. As ∗ Table 3 shows, the design points σc1 would increase with the increase of supporting pressure for both normal and lognormal variables. Con∗ versely, the design points σc2 , σt1∗ and σt2∗ are slightly higher than their mean values and increase with the increase of supporting pressure. Thus the different parameters of upper or lower rock layers have different influence on the reliability analysis. It is also worth mentioning that the value of design point ru∗ is almost the same as its mean value. Thus ru∗ is not listed in Table 3.

0.25 or higher. The choice of a lognormal distribution is motivated by the fact that the involved random variables are strictly non-negative. Since the expression of performance function G (x ) is explicitly unknown, improved RSM procedures coded in MATLAB are utilized to calculate necessary numerical results.

5.1. Reliability indexes and design points Reliability index is commonly used to estimate the security of engineering structures. In order to reflect the relationship between the supporting pressure and the roof stability of tunnel, incremental supporting pressures are applied to calculate the corresponding reliability indexes and design points. Fig. 2 presents curves of the Hasofer-Lind reliability indexes versus different values of the supporting pressure applied on the tunnel boundary. The numerical results calculated by improved RSM and MCS are listed together for comparison. Besides, the results of normal variables and lognormal variables are also plotted in Fig. 2. As expected, the reliability index increases with the increase of the applied supporting pressure. The values of supporting pressure are rationally selected until an acceptable reliability index is obtained. The comparison of the results of normal variables with those of lognormal variables shows that the reliability index corresponding to normal variables is smaller than the one of lognormal variables. It is also found that the difference of reliability indexes between different distribution

5.2. Failure probability According to Fig. 2, it is found that little differences exist among the reliability indexes obtained by improved RSM and MCS. Based on that, the reliability indexes calculated by RSM are acceptable. Both RSM and MCS will be used in all subsequent research. However, the reliability indexes computed by MCS are obtained at the expense of significantly higher time consumption when compared to the improved RSM based on FORM. Since the accuracy of MCS results lies in the number of samples, the COV is used to estimate the convergence of failure probability. Thus, the minimum sample size to get a reliable result can be determined. As Eq. (6) shows, the COV of failure probability is related not only with the number of samples but also with the failure probability. In this section, the COV of the failure probability versus the different sample sizes of MCS is obtained and discussed. To have a clear visualization of the convergence of MCS results, the COV of failure probability obtained by different sample sizes is plotted in Fig. 3. As expected, the COV of failure probability become steady with the increase of sample size. For the supporting pressure of 240 kPa, a sample size of 150,000 is big enough to obtain a nearly constant failure probability (the COV of failure probability is smaller than 1%). It is also found that a larger sample size is required when higher supporting pressure is applied. In this paper, the sample size of MCS increases from 50,000 to 1,000,000 to obtain a reliable failure probability. Fig. 4 shows the cumulative distribution function (CDF) of supporting pressure obtained by RSM. It is noticed that the distribution of variables (i.e., normal and lognormal) does not significantly influence

Fig. 2. Reliability index versus different σq for normal, lognormal variables.

5

Computers and Geotechnics xxx (xxxx) xxx–xxx

X.L. Yang et al.

Fig. 3. COV of the failure probability with different sample size. Fig. 5. Comparison of the CDFs of the supporting pressure versus different h 0 .

Table 4 Mean values of random variables in different scenarios. Rock parameters

Case 1

Case 2

Case 3

σc1 (kPa) σt1 (kPa)

10,000 150 2100

10,000 150 2100

15,000 200 2300

10,000 150 2100

15,000 200 2300

10,000 150 2100

0.2

0.2

0.2

ρ1 (kN/m3)

σc2 (kPa) σt2 (kPa) ρ2 (kN/m3) ru

Fig. 4. Comparison of the CDFs of the supporting pressure σq .

the failure probability. The effect of the COV of variables is also presented in Fig. 4. One can conclude that even a small fluctuation of the COV of random variables (i.e., 0.05 in this section) would obviously affect the curve of CDF. In these cases, the greater COV of random variables leads to higher failure probability. Therefore, the determination of the COV of variables should be taken seriously to obtain a reliable failure probability. 5.3. Reliability-based design As discussed above, the dimensionless parameters Ai and Bi , the distance from the tunnel roof to the interface of rock layers h 0 are chosen as deterministic value in reliability analyses. Due to the complexity of engineering geology, the influences of different parameters of rock layers on the probabilistic tunnel supporting pressure are investigated in this section. To study the influence of these parameters on the failure probability, a set of tests were carried out. Fig. 5 shows the cumulative distribution function (CDF) of tunnel roof pressure for different values of h 0 . As Fig. 5 presents, the roof pressure increases with the increase of distance h 0 . The variation of rock parameters in the rock strata is also investigated in 3 different scenarios. In Case 1, only 1 layer is considered. In Case 2 and Case 3, 2 layers with weak rock in the upper layer and strong rock in the lower layer and vice versa are studied. The mean value of rock parameters in these cases are shown in Table 4 and these parameters are normally distributed with COV of 0.1. Fig. 6 presents the CDF of tunnel roof pressure in different scenarios. As can be seen in

Fig. 6. Comparison of the CDFs in different cases.

Fig. 6, the roof pressures of layered rock mass are higher than homogenous rock. Besides, rock strata with weak rock in the upper layer have the highest roof pressure among these cases. The finding indicates that layered rock mass has significant influence of the stability of deep tunnels. Generally, an initial target reliability index β of 2.5 or 3 is used in reliability-based design (RBD). The corresponding values for failure probability of them are equal to 0.62％ and 0.13％ respectively. In this section, RBD is achieved to determine the required tunnel roof pressure for a target reliability index of 2.5. This tunnel pressure may also be called as “probabilistic tunnel pressure” [18,39,40]. Fig. 7 presents the probabilistic tunnel pressure for different values of COV of random 6

Computers and Geotechnics xxx (xxxx) xxx–xxx

X.L. Yang et al.

Acknowledgments Financial support was received from the National Basic Research 973 Program of China (2013CB036004), National Natural Science Foundation (51378510) for the preparation of this manuscript. This financial support is greatly appreciated. References [1] Davis EH, Gunn MJ, Mair RJ, Seneviratne HN. The stability of shallow tunnels and underground openings in cohesive material. Géotechnique 1980;30(4):397–416. [2] Antão AN, Santana T, da Silva MV, Guerra NM. Three-dimensional active earth pressure coefficients by upper bound numerical limit analysis. Comput Geotech 2016;79:96–104. [3] Wang F, Li H. Stochastic response surface method for reliability problems involving correlated multivariates with non-Gaussian dependence structure: analysis under incomplete probability information. Comput Geotech 2017;89:22–32. [4] Liu FQ, Wang JH, Zhang LL. Axi-symmetric active earth pressure obtained by the slip line method with a general tangential stress coefficient. Comput Geotech 2009;36:352–8. [5] Low BK, Tang WH. Efficient spreadsheet algorithm for first-order reliability method. J Eng Mech 2007;133(12):1378–87. [6] Luo N, Bathurst RJ. Reliability bearing capacity analysis of footings on cohesive soil slopes using RFEM. Comput Geotech 2017;89:203–12. [7] Su YH, Li X, Xie ZY. Probabilistic evaluation for the implicit limit-state function of stability of a highway tunnel in China. Tunnel Undergr Space Technol 2011;26(2):422–34. [8] Lü Q, Sun HY, Low BK. Reliability analysis of ground–support interaction in circular tunnels using the response surface method. Int J Rock Mech Min Sci 2011;48(8):1329–43. [9] Lü Q, Low BK. Probabilistic analysis of underground rock excavations using response surface method and SORM. Comput Geotech 2011;38(8):1008–21. [10] Goh ATC, Zhang W. Reliability assessment of stability of underground rock caverns. Int J Rock Mech Min Sci 2012;55:157–63. [11] Zhang W, Goh ATC. Reliability assessment on ultimate and serviceability limit states and determination of critical factor of safety for underground rock caverns. Tunn Undergr Space Technol 2012;32:221–30. [12] Goh ATC, Zhang W, Zhang Y, Xiao Y, Xiang Y. Determination of earth pressure balance tunnel-related maximum surface settlement: a multivariate adaptive regression splines approach. Bull Eng Geol Environ 2016:1–12. [13] Goh ATC, Zhang Y, Zhang R, Zhang W, Xiao Y. Evaluating stability of underground entry-type excavations using multivariate adaptive regression splines and logistic regression. Tunn Undergr Space Technol 2017;70:148–54. [14] Xu JS, Pan QJ, Yang XL, Li WT. Stability charts for rock slopes subjected to water drawdown based on the modified nonlinear Hoek-Brown failure criterion. Int J Geomech 2018;18(1):04017133. [15] Simoes JT, Neves LC, Antao AN. Probabilistic analysis of bearing capacity of shallow foundations using three-dimensional limit analyses. Int J Comput Methods 2014;11(02):1342008. [16] Massih DS, Soubra AH, Low BK. Reliability-based analysis and design of strip footings against bearing capacity failure. J Geotech Geoenviron Eng 2008;134(7):917–28. [17] Mollon G, Dias D, Soubra AH. Probabilistic analysis of circular tunnels in homogeneous soil using response surface methodology. J Geotech Geoenviron Eng 2009;135(9):1314–25. [18] Mollon G, Dias D, Soubra AH. Probabilistic analysis and design of circular tunnels against face stability[J]. Int J Geomech 2009;9(6):237–49. [19] Fraldi M, Guarracino F. Limit analysis of collapse mechanisms in cavities and tunnels according to the hoek–brown failure criterion. Int J Rock Mech Min Sci 2009;46(4):665–73. [20] Fraldi M, Guarracino F. Analytical solutions for collapse mechanisms in tunnels with arbitrary cross sections. Int J Solids Struct 2010;47(2):216–23. [21] Fraldi M, Guarracino F. Evaluation of impending collapse in circular tunnels by analytical and numerical approaches. Tunn Undergr Space Technol 2011;26(4):507–16. [22] Fraldi M, Guarracino F. Limit analysis of progressive tunnel failure of tunnels in Hoek-Brown rock masses. Int J Rock Mech Min Sci 2012;50:170–3. [23] Yang XL, Zhang R. Collapse analysis of shallow tunnel subjected to seepage in layered soils considering joined effects of settlement and dilation. Geomech Eng 2017;13(2):217–35. [24] Xu JS, Yang XL. Effects of seismic force and pore water pressure on three dimensional slope stability in nonhomogeneous and anisotropic Soil. KSCE J Civil Eng 2018; 22. http://dx.doi.org/10.1007/s12205-017-1958-y. [25] Li TZ, Yang XL. Risk assessment model for water and mud inrush in deep and long tunnels based on normal grey cloud clustering method. KSCE J Civil Eng 2018; 22. http://dx.doi.org/10.1007/s12205-017-0553-6. [26] Li TZ, Yang XL. Limit analysis of failure mechanism of tunnel roof collapse considering variable detaching velocity along yield surface. Int J Rock Mech Min Sci 2017;100:229–37. [27] Yang XL, Li ZW, Liu ZA, Xiao HB. Collapse analysis of tunnel floor in karst area based on Hoek-Brown rock media. J Central South Univ 2017;24(4):957–66. [28] Hasofer AM, Lind NC. Exact and invariant second-moment code format. J Eng Mech Div 1974;100(1):111–21.

Fig. 7. Design roof pressure with different COV of variables.

variables σci and σti which follow a normal distribution. It is found that the probabilistic tunnel pressure increases with the greater scatter in random variables σci and σti . Furthermore, the influence of σci on probabilistic tunnel pressure is much higher than that of σti . Taken h 0 = 2 m for an example, the probabilistic tunnel pressure increased from 276.52 kPa to 305.38 kPa with the increment of 0.1 in the COV of variables σci and σti . As a conclusion, the determination of probabilistic tunnel pressure depends on the uncertainties of random variables and the site selection of cavities. RBD can play a significant role in the design of the tunnel pressure.

6. Conclusions This study presents a reliability-based analysis for the roof stability of deep tunnels in layered rock strata. Considering the influence of seepage force and supporting pressure, the upper bound solution of the collapse block was developed, and the performance function of roof stability was proposed. By assuming relative geotechnical parameters in the performance function as random variables, the reliability index and the failure probability can be obtained for stability assessment of tunnel roof. The main results are summarized as follows. The two-layer model is developed from the previous research, and the collapse shape of rock mass is majorly influenced by geotechnical parameters. Based on the collapse failure mode, the performance function against roof stability was derived and then used in the reliability analysis. The reliability results obtained from the improved RSM and those calculated from MCS present good agreement. Since the MCS method is much more time consuming, RSM might be a more preferable method. Besides, the hypothesis of normal parameters was found to be conservative than that of lognormal parameters. Therefore, the statistical analysis of geo-parameters should be conducted rigorous to gain accurate probabilistic results. By varying the COV of random variables, CDF curves of the supporting pressure are obtained from the failure probabilities. The greater the COV of random variables, the higher the failure probability. On this basis, the uncertainties of materials should be accurately determined to obtain reliable failure probability. It was also found that the probabilistic tunnel pressure depends on the uncertainties of random variables as well as the site selection and rock parameters of cavities. In addition, the reliability-based design of supporting pressure was performed. When a target reliability index is given, the required supporting pressure of tunnel can be obtained. The results clearly show that the required supporting pressure decreases as the COV declines, and the failure probability is more sensitive to σci than σti . 7

Computers and Geotechnics xxx (xxxx) xxx–xxx

X.L. Yang et al.

[36] Sofianos AI. Tunnelling Mohr-Coulomb strength parameters for rock masses satisfying the generalized Hoek-Brown criterion. Int J Rock Mech Min Sci 2003;40(5):435–40. [37] Yang XL. Effect of pore-water pressure on 3D stability of rock slope. Int J Geomech 2017;17(9):06017015. [38] Saada Z, Maghous S, Garnier D. Stability analysis of rock slopes subjected to seepage forces using the modified hoek–brown criterion. Int J Rock Mech Min Sci 2012;55(10):45–54. [39] Pan QJ, Dias D. Face stability analysis for a shield-driven tunnel in anisotropic and nonhomogeneous soils by the kinematical approach. Int J Geomech 2016;16(3):04015076. [40] Pan QJ, Dias D. The effect of pore water pressure on tunnel face stability. Int J Numer Anal Meth Geomech 2016;40(15):2123–36.

[29] Low BK, Tang WH. Efficient reliability evaluation using spreadsheet[J]. J Eng Mech 1997;123(7):749–52. [30] Rackwitz R, Flessler B. Structural reliability under combined random load sequences[J]. Comput Struct 1978;9(5):489–94. [31] Kim SH, Na SW. Response surface method using vector projected sampling points. Struct Saf 1997;19(1):3–19. [32] Zheng Y, Das PK. Improved response surface method and its application to stiffened plate reliability analysis[J]. Eng Struct 2000;22(5):544–51. [33] Chen WF. Limit analysis and soil plasticity. Amsterdam: Elsevier; 1975. [34] Serrano A, Olalla C. Ultimate bearing capacity of rock masses. Int J Rock Mech Sci Geomech Abst 1997;31:93–106. [35] Sofianos AI, Halakatevakis N. Equivalent tunnelling Mohr-Coulomb strength parameters from given Hoek-Brown ones. Int J Rock Mech Min Sci 2000;39(1):131–7.

8