A fully coupled hydraulic-mechanical solution of a circular tunnel in strain-softening rock masses

A fully coupled hydraulic-mechanical solution of a circular tunnel in strain-softening rock masses

Tunnelling and Underground Space Technology 99 (2020) 103375 Contents lists available at ScienceDirect Tunnelling and Underground Space Technology j...

977KB Sizes 2 Downloads 39 Views

Tunnelling and Underground Space Technology 99 (2020) 103375

Contents lists available at ScienceDirect

Tunnelling and Underground Space Technology journal homepage: www.elsevier.com/locate/tust

A fully coupled hydraulic-mechanical solution of a circular tunnel in strainsoftening rock masses Qiang Zhanga, Cong Shaoa, Hong-Ying Wangb,c, Bin-Song Jianga, Yu-Jing Jiangb, Ri-Cheng Liua,

T



a School of Mechanics & Civil Engineering, State Key Laboratory for Geomechanics and Deep Underground Engineering, China University of Mining & Technology, Xuzhou 221116, China b Geoenvironmental Laboratory, Graduate School of Engineering, Nagasaki University, Nagasaki 8528521, Japan c College of Energy and Transportation Engineering, Jiangsu Vocational Institute of Architectural Technology, Xuzhou 221116, China

A R T I C LE I N FO

A B S T R A C T

Keywords: Circular tunnel Hydraulic-mechanical behavior Permeability Biot’s effective stress Strain-softening

This paper presents a fully coupled hydraulic-mechanical strain-softening model considering Biot’s effective stress. Both the hydraulic and mechanical parameters are considered as functions of the confining pressure and plastic shear strain. Using the proposed fully coupled hydraulic-mechanical model, the displacement, stress and pore-water pressure around a circular tunnel in both Mohr-Coulomb and Hoek-Brown strain-softening rock masses are derived using a non-associated flow rule. The proposed solutions are validated by the conventional hydraulic-mechanical elasto-brittle-plastic and strain softening rock mass using analytical and numerical methods, and the examples of the fully coupled hydraulic-mechanical strain-softening rock mass are investigated. The results show that the seepage force enlarges the displacement, residual and plastic radii. With the increase in the initial pore-water pressure, the displacement, water inflow, residual and plastic radii increase in an approximately linear manner. The increasing bulk modulus of the solid constituent increases the above four variables, which reach the approximate limit values, corresponding to the Darcy’s effective stress solutions, when Ks = 10 K0. The increasing permeability of the plastic rock mass reduces the increasing rate of the porewater pressure but increases the water inflow. For the proposed fully coupled examples, both the hydraulic parameters of permeability and Biot’s coefficient show a non-linear decrease with the increasing radius.

1. Introduction Underground engineering is commonly conducted below the groundwater table in a pervious rock mass. In such a saturated stratum, the excavation not only induces a mechanical response but also water seepage in the pores of the surrounding rock. During the deformation and seepage process, the hydraulic behavior affects the mechanical response and vice versa. This is the so-called fully coupled hydraulicmechanical coupling effect. For hydraulic-mechanical coupling engineering, the failure and deformation of the surrounding rock are induced by the variations in the effective stress, which depends on the changes in the total stress and pore pressure. Therefore, the predictions on the water inflow, displacement and plastic region are very important to the stability evaluation, design of the support and drainage systems. In the past decades, many achievements have been made regarding the change in effective stress induced by seepage around a cavity. Early analytical solutions for estimating the effective stress around a circular

tunnel in a pressurized conduit were proposed by Bouvard and Pinto (1969) and by Schleiss (1986). Considering the influence of the seepage force, Fernandez and Alvarez (1994) presented an approach that takes into account the seepage force effect by treating the fractured rock mass as a continuous porous elastic medium by using the image well method proposed by Harr (1962). The distribution of the effective stresses along a radius intersecting the springline of the tunnel was estimated from the closed-form solution derived for the simplified condition that neglects the circumferential seepage force. On the basis of conformal mapping, Kolymbas and Wagner (2007) proposed an analytical expression for the estimation of the steady-state groundwater ingress into a deep and/or shallow drained tunnel with a circular cross-section. Park et al. (2008) derived a simple closed-form analytical solution of the seepage force along the circular tunnel circumference for both zero water pressure and constant total head boundary conditions using the conformal mapping technique. Ming et al. (2010) derived the analytical solutions of pore water pressure for a steady groundwater flow into a horizontal



Corresponding author. E-mail addresses: [email protected] (Q. Zhang), [email protected] (C. Shao), [email protected] (B.-S. Jiang), [email protected] (Y.-J. Jiang), [email protected] (R.-C. Liu). https://doi.org/10.1016/j.tust.2020.103375 Received 30 April 2019; Received in revised form 11 February 2020; Accepted 23 February 2020 0886-7798/ © 2020 Elsevier Ltd. All rights reserved.

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

Nomenclature ψ c φ σc m, s,α

List of symbols σ σr , σ θ εr, εθ εrp ,εθp K, K s

Rp, Rs Ru R1 β u pw p0 G μ σp γp γp*

stress radial and tangential stresses radial and tangential strains radial and tangential plastic strains bulk modulus of the matrix and the solid constituent plastic radius and residual radius radius of unaltered region excavation radius Biot’s coefficient radial displacement pore-water pressure initial ground stress shear modulus Poisson’s ratio critical stress shear plastic strain critical shear plastic strain corresponds to the residual

n γ0p ∗, ζ k Q γw ω p0 ωγ p∗

state Dilatancy angle Cohesive strength friction angle uniaxial compression strength the first, second and third strength parameters of the Hoek-Brown rock mass annulus number of the seepage alternated region fitting parameters permeability water inflow unit weight of water variable value at a confining pressure p0 variable value for a shear plastic strain γp*

Subscripts i ' (∙)

variable value at r = Ri effective stress variable for the ith annulus

the hydraulic-mechanical coupling effect by assuming that the fraction of the permeability ratio for the plastic and elastic zones varies from zero to one. However, the method for defining the proposed fraction was not given. Considering the prior given analytical distribution of the pore-water pressure and strain-dependent permeability proposed by Brown and Bray (1982), Fahimifar and Zareifard (2009, 2014b) and Zou and Li (2016) derived the theoretical hydraulic-mechanical solutions of circular tunnels with the assumption of superposition in displacements induced by hydraulic and mechanical behavior for strainsoftening and elasto-brittle-plastic rock masses, respectively. Unfortunately, a single variable of the absolute volume strain implies an increasing permeability for both the contraction and dilation of the rock mass. This is inconsistent with the decreasing permeability with the increasing confining pressure (Wang et al., 2003). The permeability of fractures is much larger than that of a porous rock matrix (Liu et al., 2016), and the fracture degree can be represented by the plastic strain according to the elasto-plastic theory. Similar to the evolutions of mechanical parameters for strain-softening rock mass, the hydraulic parameters, i.e., Biot’s coefficient and permeability, vary with the confining pressure and plastic strain during the formation and through the cracking process (Hu et al., 2010). In view of the limitations of the previous hydraulic-mechanical solutions of circular tunnels, this paper first analyzed the variations in permeability and Biot’s coefficient during the progressive failure process according to the experimental tests and proposed the fully coupled hydraulic-mechanical strain-softening model (HMSS). Then, the analytical solutions of stress and displacement for the surrounding rock in both the Mohr-Coulomb (MC) and Hoek-Brown (HB) HMSS rock masses were derived. Finally, the proposed solutions were validated using analytical and numerical methods for drained and saturated tunnels, and the influences of pore-water pressure, Biot’s coefficient and permeability of fully coupled HMSS rock mass on the plastic radius, displacement and water inflow were further studied.

tunnel in a fully saturated isotropic semi-infinite aquifer. Considering the Biot’s effective stress and the pore-water pressure distribution proposed by Kolymbas and Wagner (2007), Fahimifar and Zareifard (2013) derived a closed-form solution of an unlined pressure tunnel. Among these solutions, the surrounding rock was considered to be an isotropic and homogeneous elastic rock medium, and the permeability was assumed to be constant by neglecting its changes due to fracture. In this regard, the previous elastic seepage induced stress and deformation solutions were only suitable for shallow cavities, in which the stress was small and the surrounding rock mainly works at the elastic state. For the deep underground engineering, the initial ground stress is very high, and it usually induces rock mass fracture, which thereafter affects the seepage of groundwater. So far, many achievements have been made for the cavities in the drained stratum using the elastobrittle (perfectly)-plastic model (Sharan, 2003; Park and Kim, 2006; Zhang et al., 2012), strain-softening model (Brown et al., 1983; Alejano et al., 2010; Lee and Pietruszczak, 2008; Park et al., 2008; Wang et al., 2010; Zhang et al., 2012b, 2019; Cui et al., 2015) and elasto-plastic coupling model (Zhang et al., 2018). For the cavities in the saturated stratum, the seepage force acted as a body force on the surrounding rock and caused significant difficulties in the analytical derivations of displacement and stress. To overcome this issue, two methods, i.e., total stress and the superposition concept, were commonly used. Brown and Bray (1982) considered the effects of the hydraulic-mechanical coupling on the fracture zone analysis. However, the seepage effect on the elastic region was negligible. Assuming that the failure was induced by Terzaghi’s effective stress, Zareifard and Fahimifar (2014) proposed the elastic-brittle-plastic solution of a circular tunnel by adopting Biot’s effective stress and constant hydraulic parameters. However, both the failure and deformation of the surrounding rock were induced by the change of the effective stress. Moreover, the pore pressure distributions of steady state seepage were obtained with the assumption of constant permeability of the surrounding rock, and this was inconsistent with the experimental results, where the permeability depends on both the confining pressure and fracture distribution (Wang et al., 2003; Chen et al., 2007, 2014; Yang et al., 2015; Xu and Yang 2016). In view of the fully coupled effect between the hydraulic and mechanical parameters, i.e., pore pressure influences the effective stress, and mechanical failure enlarges the permeability of the surrounding rock, many efforts have considered strain and/or stress-dependent permeability (Brown and Bray 1982; Min et al., 2004). Fazio and Ribacchi (1984) and Carosso and Giani (1989) considered the effects of

2. Problem definition 2.1. Proposal of a fully coupled HMSS model Many triaxial compression tests demonstrate the strain-softening behavior in the post-failure region for various rock masses (Martin and Chandler 1994; Alejano and Alonso 2005; Zhao and Cai 2010, Zhang et al., 2013). With the increase in deviatoric stress, first the initial 2

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

Therefore, the surrounding rock in the water-rich stratum manifests a fully coupled HMSS behavior, in which both the hydraulic and mechanical parameters vary with the failure degree and confining pressure and vice versa.

cracks are compacted, then new micro-cracks occur and gradually extend, leading to a decrease first and then an increase in the apparent porosity (Chen et al., 2014; Wang et al., 2003). Fig. 1(a) shows the typical evolutions of permeability and deviatoric stress curves for granite. In the elastic stage, the permeability first decreases then increases due to the compaction effect, and the variation in permeability is very small with a maximum amplitude less than 5%. Therefore, the variation in permeability in the elastic state can be neglected. Once the stress reaches the initial yield value, the permeability significantly increases. With a further increase in the axial strain, the increasing rate of permeability gradually decreases, and remains a relative constant value in the residual state. This is induced by the increase in crack aperture, which can be represented by the plastic strain. For example, the permeability values of the residual specimens at confining pressures of 5 MPa, 15 MPa and 30 MPa are 28.05 × 10−13 m2, 7.78 × 10−13 m2 and 0.93 × 10−13 m2, respectively. Obviously, the permeability of the rock mass at a low confining pressure is much larger than that at a high confining pressure. The main reason is that an obvious tensile-shear failure mode with a main through fracture happens for specimens at low confining pressure, and a ductile failure with many potential blind cracks corresponds to specimens at a high confining pressure. A through crack usually leads to greater permeability. On the other hand, the increasing confining pressure reduces the fracture aperture, which is directly proportional to the permeability. Fig. 1(b) shows the permeability evolution with volume strain. A slight decrease in permeability occurs at the initial state; thereafter, a significant increase in permeability follows during the failure process, and it finally reaches a relatively constant value in the residual state. Thus, the permeability evolution function proposed by Brown and Bray (1982) would be no longer suitable for the permeability evaluation of rock mass in the contraction stage. Two main factors, i.e., confining pressure and failure degree, should be considered to describe the permeability evaluation of the elastic and plastic rock masses. The failure and deformation of the hydraulic-mechanical coupling rock mass are induced by the changes in the effective stress. A general Biot’s effective stress, in which the Darcy effective stress is only a special condition, is considered, and can be expressed as follows:

Fig. 3a shows a circular tunnel excavated in a homogenous and isotropic hydraulic-mechanical coupling rock mass. Before the tunnel is excavated, a hydrostatic ground effective stress p0′, corresponding to a pore-water pressure p0w , is imposed throughout the domain. After the excavation, the internal supporting pressure gradually drops from its ′ , and the pore water is also drained from the initial value p0′ to σr,1 periphery of the cavity. For the plane strain condition, the seepage only occurs in the radial direction, leading to an axisymmetric pore-water pressure distribution. In this study, a steady seepage in the saturated surrounding rock is considered, and the gravity is neglected. To obtain the analytical seepage solutions, a much larger seepage unaltered region with a radius r = Ru is assumed. Compared to the infinite condition, the error induced by the above assumption is negligible when Ru ≫ R1 (Li et al., 2004). During the decrease in the internal effective stress, the convergence displacement u1 at r = R1 gradually increases. A plastic region with ′ is lower than a radius r = Rp develops around the tunnel when σr,1 critical value σp′ . Moreover, a residual region with radius r = Rs occurs when the plastic strain exceeds the critical value for the strain-softening rock mass. To obtain the displacement and stress of the fully coupled

(a) 30

K Ks

200

30 MPa

150

15 15 MPa

100

10 15 MPa

(1) 5 MPa

5

50 30 MPa

0 0.0

0.2

Deviatoric stress (MPa)

2 -13

250

Permeability Stress

20

where σ ′ and σ are the Biot’s effective stress and total stress, respectively, pw is the pore-water pressure, and β is the Biot’s coefficient, which depends on the bulk modulus of solid matrix and solid grain material (Biot, 1941; Nur and Byerlee, 1971), and can be expressed as follows:

β=1−

300 5 MPa

25 Permeability (10 cm )

σ ′ = σ - βpw

2.2. Calculation model

0.4

0.6

0 0.8

ε1 (%)

(b)

-13

10

(2) -14

10 2

Permeability (m )

where K and Ks are the bulk modulus of the matrix material and the solid constituent, respectively. The porous and fractured rock masses show a significant confiningpressure effect. For example, the material property parameters, including the bulk modulus, vary with the confining pressure (Brown et al., 1989). With the increase in deformation, the increasing potential for cracks significantly decreases the apparent bulk modulus of the rock mass. According to Eq. (2), a confining pressure dependent K leads to a confining pressure dependent β. Fig. 2 shows the variation in Biot’s coefficient with axial strain under various confining pressures. Obviously, for a constant confining pressure condition, β shows a nonlinear increase with the increase in axial strain. The Biot’s coefficients of the specimens under the confining pressures of 10 MPa, 15 MPa and 30 MPa increase from 0.55, 0.41 and 0.31 to 0.95, 0.64 and 0.42, with an increasing rate of 72.73%, 56.10% and 35.48%, respectively. The increasing rate of β under a low confining pressure is larger than that under a high confining pressure. On the other hand, β gradually decreases with the increase in the confining pressure.

-15

10

σ3′=0.75 MPa

-16

10

σ3′=1.5 MPa -17

σ3′=2.25 MPa

-18

Initial yield point

10

σ3′=3.0 MPa

10

Contraction

Dilation

-19

10

0.4

0.2

0.0

-0.2

-0.4

-0.6

-0.8

-1.0

Volume strain (%) Fig. 1. Variations of permeability during the loading process under different confining pressures: (a) Stress-strain curve and permeability evolution (Chen et al., 2014), and (b) Permeability evolution under various confining pressures (Wang et al., 2003). 3

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

1.2

10 MPa

1.0 Biot's coefficient

seepage unaltered region is assumed to be one region. Each divided annulus is assumed to be an isotropic medium with constant material property parameters, which depends on the plastic strain and confining pressure (radial stress for the plane strain cylinder) of the selected representative point (Zhang et al., 2012). Fig. 3b shows the typical dividing annulus model of ith annulus with inner and outer radii of Ri-1 and Ri, respectively. The outer radius of each annulus is employed as the representative point. The objective of the hydraulic-mechanical problem is to determine the radial effective stress σr′, tangential stress σθ′ and displacement u with an arbitrary radius r. When the seepage body force is not con′ to p0′ with the sidered, the radial stress monotonously increases from σr,1 increase in radius r. In view of this situation, the radial stress can be discretized into a series of prior known values, and the problem then becomes to determine the corresponding radius of the discretized radial stress (Han et al., 2013). For the perfectly plastic MC rock mass, the closed-form solutions can be explicitly expressed. However, those of the HB rock mass can only be expressed implicitly because of the nonlinear yield criterion expression forms. If the radial stress is a monotonic increasing function, the solutions of plastic HB rock mass around a circular tunnel can be calculated in a difference manner. For a stable hydraulic-mechanical tunnel, the seepage force does not induce tensile failure in rock mass for most realistic values of hydraulic-mechanical parameters. Therefore, the requirement of monotonously increasing radial stress will be satisfied, and detailed derivations will be introduced in Section 4.

5 MPa

Singular value

0.8 20 MPa 0.6 30 MPa 0.4

0.2

0.0

0.2

0.4

0.6

0.8

1.0

1.2

ε1/ε1peak Fig. 2. Variations in Biot’s coefficient during the loading process under different confining pressures (Hu et al., 2010).

2.3. Governing equations For the strain-softening behavior, the softening parameter governs the evolutions of the material properties and should be determined first. The most widely used plastic shear strain γp is employed as the softening parameter, and it can be expressed as follows:

γ p = εθp − εrp

(3)

where εrp and εθp are the radial and tangential strains, respectively. For the axisymmetric problem, the seepage behavior only occurs in the radial direction. The seepage continuity differential equation of the ith annulus can be expressed as follows (Bear, 1972):

d 2p(wi) dr 2

+

w 1 dp(i) =0 r dr

(4)

where p(wi) is the pore-water pressure, and the subscript (i) denotes that the corresponding variable is a function of r in the ith annulus. With the aid of Eq. (1), the equilibrium differential equation of ith annulus can be written as follows:

dσ ′r ,(i) dr



σ ′θ,(i) − σ ′r ,(i) r

+ Fi = 0

(5)

where σr′,(i) and σθ′,(i) are the radial and tangential Biot’s effective dp(wi)

stresses, respectively. Fi = βi dr is the seepage body force, and the subscript i denotes the corresponding variable value of the ith annulus at r = Ri. The plastic behavior is governed by the yield criterion. The most widely used linear MC criterion and nonlinear HB criterion are considered, and they can be uniformly expressed as follows:

Fig. 3. Hydraulic-mechanical calculation model of a circular tunnel, (a) overall model, and (b) dividing annulus model.

σθ′,(i) − σr′,(i) = H (σ ′r ,(i), γi p)

HMSS surrounding rock mass around a circular tunnel, a concerned region R1 ≤ r ≤ 10R1 is divided into n1 annuli, and the least seepage altered region 10R1 ≤ r ≤ Ru is divided into n2 annuli. According to Saint-Venant principle (Timoshenko and Gooodier, 1970), a very large seepage altered region with radius Ru = 200R1 is considered, and the

(6)

where H is a function of radial stress and plastic shear strain. Using the superscripts of “MC” and “HB” to distinguish the MC and HB rock masses, respectively, H can be formulated as follows:

H MC (σ ′r ,(i), γi p) = (Ni − 1) σr′,(i) + Yi 4

(7a)

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

H HB (σ ′r ,(i), γi p) = σci ⎛⎜ ⎝

mi σ ′r ,(i) σci

+ si ⎞⎟ ⎠

αi

where ω is any of the hydraulic and mechanical parameters, i.e., c, φ, σc, m, s, α, ψ, E, μ, G, K, k and β, k is the permeability of the surrounding rock, ω0 is the value corresponding to the uniaxial compression state, and ω p0 and ω γ p ∗ are the values of ω corresponding to a confining pressure p0′ and plastic shear strain γp*, respectively.

(7b)

2c cos φ

1 + sin φ

i where Ni = 1 − sin φi , Yi = 1 − sin φi , φi and ci are the friction angle and i i cohesive strength, respectively, of the MC rock mass; mi, si and αi are the first, second and third strength parameters, respectively, and σci is the uniaxial strength of the HB rock mass. The geometry equations between the radial strain, tangential strain and radial displacement for the ith annulus can be written as follows:

du (i)

εr ,(i) =

dr

,

2.5. Boundary conditions The radial effective stress, pore-water pressure and displacement are continuous at the boundaries of the adjacent annuli. Therefore, the boundary conditions can be expressed as follows:

u (i)

εθ,(i) =

r

(8)

w = p1w p(1) ⎧ σ ′r ,(1) = σ ′r ,1, ⎪ σ ′r (i) = σ ′r (i + 1) = σ ′r , i, p(wi) = p(wi + 1) = piw , ⎪ ⎪ u (i) = u (i + 1) = ui (r = Ri ) ⎨ ⎪ p w = pw 0 ⎪ (n) ⎪ σ ′r (n) = p′0 , ⎩

where εr ,(i) and εθ,(i) are the radial and tangential strains, respectively, and u(i) is the radial displacement. According to the plastic theory, the total strain of the plastic rock mass is composed of elastic strain and plastic strain, and they are related in the following:

εr ,(i) = εre,(i) + εrp,(i) ⎫ εθ,(i) = εθe,(i) + εθp,(i) ⎬ ⎭

εθe,(i) =

(r ⩾ Ru ) (r → ∞) (14)

(9)

where εre and εθe are the radial and tangential elastic strains, respectively. For the HMSS problem, both the elastic and plastic strains are induced by the changes in pore-water pressure and boundary pressure. Hence, the generalized Hooke’s law provides the relation between the elastic strains εre,(i) and εθe,(i) and Biot’s effective stresses σr′,(i) and σθ′,(i) in ith annulus as follows (Timoshenko and Gooodier, 1970):

εre,(i) =

(r = R1)

1 [(1 − μi )(σ ′r ,(i) − p′0 ) − μi (σ ′θ,(i) − p′0 )] ⎫ 2Gi 1 [(1 − μi )(σ ′θ,(i) − p′0 ) − μi (σ ′r ,(i) − p′0 )] ⎬ 2Gi ⎭

3. Hydraulic solutions The general solution of p(wr ) can be obtained by the integrals of Eq. (4). Using the second pore-water pressure boundary condition in Eq. (14), the pore-water pressure of each annulus can by uniformly formulated as follows:

(10)





(15)

where γw is the unit weight of water, and Qi is the water inflow through the periphery of the ith annulus. During the seepage process, the water flow obeys the continuity equation, which requires a constant water inflow for each annulus, i.e., Q1 = ∙∙∙ = Qn = Q. By summing both the left and right hands of Eq. (15) for n annuli, the water-inflow can be determined as follows:

where Gi is the shear modulus and μi is Poisson’s ratio. The plastic theory implies that the flow directions and values of the plastic strain components depend on the selected potential function. For the geotechnical material, the MC plastic potential function is widely used. Then, the radial and tangential plastic strains result in the following:

εrp,(i) + κi εθp,(i) = 0

γw Qi ⎛ r ⎞ ln 2πki ⎝ Ri − 1 ⎠

p(wi) = piw− 1 +

Q=

(11)

2π (p0w − p1w ) n

γw ∑

where κi = (1 + sinψi)/(1 − sinψi) and ψi is the dilation angle of ith annulus.

j=2

where ρj =

1 kj

ln ρj

Rj Rj − 1

(16)

.

2.4. Evolution of material parameters

With the aid of Eq. (16), Eq. (15) can be reformulated as follows:

The geotechnical materials generally show brittle-ductile transition characteristics when the confining pressure is increasing. The brittleductile transition effect results in an increasing critical plastic shear strain, from which the rock mass steps into the residual state, with the increase in confining pressure. According to the previous studies (Zhang et al., 2018), an exponential function is employed to represent the confining-pressure dependent critical plastic shear strain γp* as follows:

p0w ⎡ ⎢∑ ⎣ j=2 =

i

γ(pi)∗

=

γ0p ∗ [1

−e

−ζσ ′r ,(i)

]

p(wi)

1 kj

ln ρj +

1 ki

n

r 1 ln R ⎤ + p1w ⎡ ⎢ ∑ kj ln ρj − i⎥ j i = + 1 ⎦ ⎣ n

∑ j=2

1 kj

ln ρj

1 ki

r ln R ⎤ i⎥ ⎦

(17)

4. Stress and displacement solutions

(12)

4.1. Solution for the elastic rock mass

and ζ are the fitting parameters. where As introduced above, both the hydraulic and mechanical material property parameters depend on the plastic shear strain and confining pressure. The evolutions of the hydraulic and mechanical property parameters can be determined through the triaxial hydraulic-mechanical coupling cyclic loading and unloading tests. For convenience, a plane function is considered to represent the evolution of material property parameters. The evolution function can be uniformly expressed as follows:

γ0p ∗

For the seepage-unaltered elastic region, the pore-water pressure remains unchanged, leading to a zero seepage force. Using the second and the fourth boundary conditions in Eq. (14), the stress and displacement can be obtained according to Lame’s solution as follows: 2 2 ⎧ σ ′r ,(n) = −(p′0 − σ ′r , n ) Rn r + p′0 ⎪ 2 2 σ ′θ,(n) = (p′0 − σ ′r , n ) Rn r + p′0 ⎨ ⎪ u (n) = Rn2 (p′0 − σ ′r , n ) (2Gr ) ⎩

ω (σ ′r ,(i), γi p) = ω0 [1 − (1 − ω p0 ω0 ) σ ′r ,(i) p′0 ][1 − (1 − ω γ p ∗ ω0 ) γi p γ p ∗]

(18)

By combining the radial and tangential kinematic equations of Eq. (8), the deformation compatibility equation can be obtained as follows:

(13) 5

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

dεθ,(i) dr

=

εr ,(i) − εθ,(i)

γ(pi)

(19)

r

u (i) 1 [(Ni − Ni μi − μi ) σ ′r ,(i) + (1 − μi ) Yi = (1 + κi ) ⎧ − ⎨ 2Gi ⎩ r

By substituting the radial and tangential strains of Eq. (10) into Eq. (19) then combining with Eq. (5), a second-order ordinary differential equation of radial effective stress can be derived. Using the second boundary condition in Eq. (14), i.e., σr′,(i) = σr′, i − 1 at r = Ri-1 for the ith annulus, the objective radial stress of the ith annulus can be expressed by integrating the ordinary differential equation as follows:

σr′,(i)

γw Qβi 4πki

, Di =

σ ′r , i − σ ′r , i − 1 + Ci ln ρi (1 − μi ) ρi2 − 1

Theoretically speaking, taking the same derivation method used for the radial stress of the MC plastic rock mass, the hydraulic-mechanical solutions of a circular tunnel in the HB rock mass can be obtained. However, the primitive function cannot be explicitly expressed because of the complex expression forms of the HB yield criterion and the seepage force. Since the seepage body force leads to a non-monotonous distribution of the radial effective stress in the elastic region, the radial stress cannot be prior given as a series of known constant values. In other words, a monotonously increasing distribution of radial stress is essential for the calculation method of determining the radii of a series of prior known radial stresses. By substituting Eqs. (6), (7b) and (15) into Eq. (5), the equilibrium differential equation can be rewritten as follows:

(20)

.

By substituting Eq. (20) into Eq. (5), the tangential stress can be obtained as follows:

σθ′,(i) =

1 − 2μi Ci R R 2 ln i − 1 + Di ⎡ρi2 + ⎛ i ⎞ ⎤ + σr′, i − 1 + Ci ⎢ 1 − μi 1 − μi r ⎝r ⎠⎥ ⎣ ⎦

(21)

Substituting Eqs. (20) and (21) into Eq. (10) and using Eq. (8), the displacement of the ith annulus for elastic region can be obtained as follows:

u (i) =

r ⎧ 1 − 2μi R R 2 Ci ln i − 1 + Di ⎡ (1 − 2μi ) ρi2 + ⎛ i ⎞ ⎤ ⎢ − 2Gi ⎨ 1 μ r r ⎠⎥ ⎝ i ⎣ ⎦ ⎩ + (1 − 2μi )(σ ′r , i − 1 + Ci − p′0 ) ⎫ ⎬ ⎭

dσ ′r ,(i) dr

By substituting Eqs. (6) and (15) into Eq. (5), a one-order differential equation of radial stress can be obtained. Integrating on both sides of the ordinary differential equation and using the second boundary condition in Eq. (14), the radial stress of the MC plastic rock mass can be obtained as follows: ⎜

⎟⎜

Ni − 1

+



Yi − Ci 1 − Ni



⎟⎜



Ni − 1

+

(24)

dr

+ κi

u (i) r

= fi (r ) 1 [A1i (r 2Gi

where

fi (r ) = Yi − Ci B1i 1 − N + B2i , i

Ri − 1 ) Ni − 1 + A2i ],

(

A1i = B1i σ ′r , i − 1 −

where Δσ ′ = , σp′ is the critical supporting pressure corren2 sponding to r = Rp, and n2 is the annulus number of the plastic region. By substituting Eqs. (6), (7b) and (15) into Eq. (5), the equilibrium equation can be approximated by a difference equation as follows:

),

H (σ ′r , i, γi p) − Ci σ ′r , i − 1 − σ ′r , i = 1 − ρi ρi

B1i = Ni (κi − κi μi − μi ) + 1 − μi − κi μi , A2i = B2i = Yi (κi − κi μi − μi ) + (1 − 2μi )(1 + κi ) p0′. Integrating Eq. (25) over the range of [Ri-1, Ri] and using the second displacement continuum condition of Eq. (14), the displacement of the ith annulus for the MC plastic rock mass can be obtained as follows: u (i) =

(30)

σ ′ p − σ ′r,1

(25) Yi − Ci 1 − Ni

(29)

σr′, i = σr′,1 + (i − 1)Δσ ′

Substituting Eqs. (8)–(10) into Eq. (11) results in a differential equation of radial displacement of the ith annulus as follows:

du (i)

(28)

Eq. (29) denotes that the seepage constant Ci should be less than the uniaxial compression strength of the plastic rock mass. The condition of Eq. (29) breaks down only when the residual uniaxial compression strength is very small. Even if this condition is violated, a positive firstorder derivative of radial stress may still be valid throughout the analysis since the radial stress is nonnegative. Thus, the requirement of monotonously increasing radial stress will, therefore, be satisfied for most realistic values of hydraulic-mechanical parameters. Assuming the right hand of Eq. (29) is positive, the radial stress of the plastic rock mass can be discretized as a series of known values as follows:

(23)

Yr − Ni Ci 1 − Ni

αi 1 ⎡ ⎛ mi σ ′r ,(i) ⎞⎟ − C ⎤ s σ + ⎜ ci i i⎥ r⎢ ⎠ ⎣ ⎝ σci ⎦

Ci < σci si αi

Substituting Eq. (23) into Eq. (6) and with the aid of Eq. (7a) yields the tangential stress as follows:

Y − Ci ⎞ ⎛ r ⎞ σθ′,(i) = Ni ⎛σ ′r , i − 1 − i 1 − Ni ⎠ ⎝ Ri − 1 ⎠ ⎝

=

If the right side of Eq. (28) remains positive, the radial stress would follow a monotonously increasing distribution form for the plastic rock mass. For the presented problem, the radial stress is positive. Without loss of generality, the following sufficient relation leads to an absolute monotonously increasing radial stress in the plastic region:

(22)

4.2. Solution for the MC plastic rock mass

Y − Ci ⎞ ⎛ r ⎞ σr′,(i) = ⎛σ ′r , i − 1 − i 1 − Ni ⎠ ⎝ Ri − 1 ⎠ ⎝

(27)

4.3. Solution for the HB plastic rock mass

Ci R R 2 = ln i − 1 + Di ⎡ρi2 − ⎛ i ⎞ ⎤ + σr′, i − 1 ⎢ 1 − μi r r ⎠⎥ ⎝ ⎣ ⎦

where Ci =

− (1 − 2μi ) p′0 ] ⎫ ⎬ ⎭

(31)

Eq. (31) leads to the relationship between the inner and outer radii of the ith annulus as follows:

1 [A1i f1i (r ) + A2i f2i (r ) − A1i f1i (Ri ) − A2i f2i (Ri )] + ui (Ri r ) κi 2Gi r κi

ρi =

H (σ ′r , i, γi p) − Ci H (σ ′r , i, γi p) - Ci + Δσ ′r , i

(32)

(26) where f1i (r ) =

1 r κi + Ni , i − 1 κi + Ni RiN −1

f2i (r ) =

r 1 + κi , 1 + κi

where Δσr′, i = σr′, i − 1 − σr′, i is the increment of radial stress for the ith annulus. With the aid of Eqs. (6), (7b), (9) and (10), the deformation compatibility equation of Eq. (19) can also be approximated by a difference equation for the HB plastic rock mass as follows:

and ui is the displacement at

r = Ri. By substituting Eqs. (8) and (26) into Eq. (3) and with the aid of Eqs. (9)–(11), the plastic shear strain can be obtained as follows: 6

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

εθp, i - εrp, i ⎞ H (σ ′r , i, γi p ) Δεθe, i Δεθp, i = −(1 − ρi ) ⎛⎜ + + ⎟ 2Gi ρi ρi ⎠ ⎝ 1 − ρi

Start

(33)

Input the geometric and material parameters

Δεθe,(i)

where is the tangential elastic strain increment of the ith annulus that can easily be obtained using the incremental form of constitutive equation Eq. (10) as follows:

⎨ Δεθe, i = ⎩

1 [(1 − μi )Δσ ′r , i − μi Δσ ′θ, i] 2Gi 1 [(1 − μi )Δσ ′θ, i − μi Δσ ′r , i] 2Gi

Assume the radial effective stress

Calculate the radial effective stress by Eqs. (39) and (40)

(34) Elastic loop

e ⎧ Δεr , i =

Assume the water inflow Q~

where Δσθ′, i is the tangential stress incremental that can be obtained using the HB yield criterion of Eqs. (6) and (7b) as follows:

Δσθ′, i = Δσr′, i + H (σ ′r , i − 1, γi p) - H (σ ′r , i, γi p)

Calculate the tangent effective stress and displacement by Eqs. (21) and (22)

(35)

No

Judge whether stress state satisfies the yield criterion of Eq. (6)

By combining Eqs. (11) and (33), the radial plastic strain increment Δεrp, i can be obtained as follows:

Yes

The radial strain εr , i − 1, tangential strain εθ, i − 1 and plastic shear strain p γi − 1 can be determined as follows: p e ⎧ εr , i − 1 = εr , i + Δεr , i + Δεr , i ⎪ e εθ, i − 1 = εθ, i + Δεθ, i + Δεθp, i ⎨ p ⎪ γi − 1 = γi p + Δεθp, i − Δεrp, i ⎩

(37)

MC rock mass

Calculate the radius and displacement with Eq. (38)

Calculate the shear plastic strain with Eq. (37)

Calculate the shear plastic strain with Eq. (27)

Update the material parameters using Eq. (13)

n2

ui = εθ, i Rp

Discrete the radial effective stress with Eq. (30)

Calculate the radial and tangential effective stress, displacement with Eqs. (23), (24) and (26)

After n2 iterations, the dimensionless radius ρ of each annulus can be obtained. The radius and displacement corresponding to the prior given radial stress can be obtained as follows:

Ri = Rp

HB rock mass

Plastic loop

(36)

Update the assume radial effective stress

Δεrp, i = −βi Δεθp, i

~ Update the assumed water inflow Q

Calculate the plastic radius with Eq. (41)

∏ ρj ⎫ ⎪ j=i+1 n2 ⎬ ∏ ρj ⎪ j=i+1 ⎭

Judge whether and/or

No

(38)

for MC rock mass for HB rock mass Yes

Calculate the required internal pore-water pressure

5. Calculation procedure

with Eq. (17)

Although the stress and displacement expressions have been derived, the radial stress and displacement of the MC rock mass and the radius of the HB rock mass cannot be obtained immediately due to the variations of the hydraulic and mechanical parameters with the confining pressure and plastic shear strain. The solutions can be theoretically obtained by solving the 2n-order equations; however, the equations are complex and the exact results are difficult to obtain if the number n is very large. In this regard, a simple step-by-step calculation method is presented. The exact solutions should satisfy both the hydraulic and mechanical boundary conditions in Eq. (14). Because the effective stress depends on the gradient of the pore-water pressure, the distribution of the pore-water pressure should be determined first. Eq. (15) denotes that the exact solutions can be obtained if the water-inflow and permeability ∼ are given. Thus, a water-inflow Q is assumed to calculate the seepage force. Then, assuming a radial stress σ͠ r′, n at r = Ru, the solutions for the seepage-unaltered elastic region can be obtained. Then, using the radial stress continuum condition at r = Ru for the nth annulus and the outmost seepage-unaltered region, i.e., the third type of Eq. (18) and Eq. (22), the radial stress σr′, n − 1 can be obtained as follows:

No

Yes Calculate the pore-water pressure of each annulus with Eq. (17)

End

Fig. 4. Flow chart for the calculation sequences of HMSS rock mass around a circular tunnel.

where Tn + 1 = p0′ − σ͠ r′, n . Using the radial continuum conditions at the adjacent boundaries of the annuli in alternating the altered pore-water pressure elastic region, i.e., Eq. (22), the recursive equation of radial stress can be uniformly obtained as follows:

σr′, i − 1 C

=

σr′, n − 1 C

[1 + (1 −

⎫ C − (1 − 2μn ) ⎛ n − p′0 ⎞ ⎝ 2 ⎠⎬ ⎭

2σ ′r , i + 1 −iμ ln ρi 1 − ρi2 ⎧ Gi (1 − 2μi ) Ci i ln ρi − Ti + 1 + 2(1 − μi ) ⎨ Gi + 1 2(1 − μi ) 2(ρi2 − 1) ⎩ [1 + (1 − 2μi ) ρi2 ] − (1 − 2μi )(0.5Ci − p′0 )

2σ͠ ′r , n + 1 −nμ ln ρn 1 − ρn2 ⎧ Gn (1 − 2μn ) Cn n ln ρn − Tn + 1 + = 2(1 − μn ) ⎨ Gn + 1 2(1 − μn ) 2(ρn2 − 1) ⎩ 2μn ) ρn2 ]

Judge whether

2μi + 1 ) ρi2+ 1

⎫ ⎬ ⎭

(40)

+ (1 − 2μi + 1 )(σ ′r , i + 0.5Ci + 1 − p′0 ) . where Ti + 1 = 2Di + 1 (1 − If the plastic behavior occurs in the ith annulus, we have Ri1 < Rp < Ri. The elastic stress at the elasto-plastic boundary should satisfy the yield criterion. By substituting the radial stress of Eq. (20)

(39) 7

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

(a)

Table 1 Geometric and material property parameters of elasto-brittle (perfectly)plastic rock mass.

2.1

β =0.5

Scheme

Value

Radius of tunnel, R1 (m) Initial effective stress, p0′ (MPa)

1.0 1.0

Internal effective pressure, σ1′ (MPa) Initial pore-water pressure, p0w (MPa)

0 0.3

Internal pore-water pressure, p1w (MPa)

0

Young’s modulus, E (MPa) Poisson’s ratio, μ (-) Dilation angle, ψ (Deg) Permeability, k (m/s) c0 and c γ p * , (MPa)

5000 0.2 0 10−7 0.276, 0.055

0.6

0.0

1.5

σ′/p0′

φ0 and φ γ p * , (Deg)

35, 30 0.2, 0.05

s0 and s γ p *

10−4, 10−5

Uniaxial strength σc, (MPa)

50

Proposed method Numerical method

1

2

3

4

5

6

r/R0

(b)

0.6 0.5 0.4

ρp2

w

p /p0

w

1 − 2μi Ci ⎡2(σ ′ − σ ′ + ln ρp ⎤ 2 Ci − H (σ ′r , i − 1, 0) r,i r , i − 1) + ⎥ρ −1 ⎢ 2(1 − μi ) 1 − μi ⎦ p ⎣ (41) =0

0.3 Proposed method Numerical method

0.2

Ri . Rp

By replacing ρi in Eq. (40) with ρp and combining with Eq. (41), the plastic radius Rp and critical pressure σp′ can be obtained using the Newton-Raphson method (Burden and Faires 2010). For the ith annulus, by substituting r = Ri-1 into Eq. (23), the recursive equation for radial stress can be formulated as follows:

Y − Ci ⎞ Ni − 1 Y − Ci + i σr′, i − 1 = ⎛σ ′r , i − i ρ 1 − Ni ⎠ i 1 − Ni ⎝

0.9

0.3

and tangential stress of Eq. (21) into the yield criterion of Eq. (6), a nonlinear equation of plastic radius Rp can be established for both the MC and HB rock masses as follows:



β =0

1.2

m0 and m γ p *

where ρp =

β =1.0

1.8

0.1 0.0

1

3

5

7

9

11

13

15

r/R1



Fig. 5. Distributions of tangential effective stress and pore-water pressure for elastic-brittle-plastic rock mass ( p0w =0.3 MPa): (a) Tangential stress, and (b) Pore-water pressure.

(42)

Different from the MC plastic rock mass, the radial stress in the plastic region of the HB rock mass is prior given. The objective is to determine the radii corresponding to these radial stresses. With the aid of Eq. (32), the dimensionless radius ρi can be obtained. Then, using Eqs. (33)–(37) the plastic shear strain γi p− 1, which is used to update the hydraulic-mechanical properties of the (i-1)th annulus, can be determined. ′ for After several iterations, the required internal effective stress σ͠ r,1 ∼ the MC rock mass and required excavation radius R1 for the HB rock ∼ mass can be obtained. If σ͠ r′,1 = σr′,1 for the MC rock mass and/or R1 = R1 for the HB rock mass, then the assumed radial stress σ͠ r′, n is corrected;

otherwise, σ͠ r′, n should be updated until σ͠ r′,1 = σr′,1 for the MC rock mass ∼ and/or R1 = R1 for the HB rock mass. Additionally, according to the updated permeability of each annulus, the required internal pore-water ∼ p1w corresponding to the assumed Q can be determined using pressure ∼ w ∼ Eq. (16). For the exact solution, p1 should equal to the given internal ∼ pore-water pressure p1w . Otherwise, the assumed water inflow Q should be updated and the above calculation procedure should be repeated ∼ p1w = p1w . The correct values of Q and σ͠ r′, n can be identified using until ∼

Table 2 Hydraulic-mechanical solutions of a circular tunnel in an elastic-brittle (perfectly)-plastic rock mass with an initial pore-water pressure p0w = 0.3 MPa (data in parenthesis are analytical results of Zhang et al. (2012)). Rock type

β

Proposed method EBP

FLAC and/or analytical methods EP

EBP

EP

Rp/R1

u1 E0 p′0 R1

Rp/R1

u1 E0 p′0 R1

Rp/R1

u1 E0 p′0 R1

Rp/R1

u1 E0 p′0 R1

MC

1.0 0.8 0.5 0

2.237 2.123 1.973 1.760

7.758 6.781 5.579 4.038

1.229 1.217 1.197 1.165

1.890 1.780 1.619 1.364

2.250 2.113 2.000 1.761

7.882 6.794 5.594 4.044

1.218 1.223 1.201 1.165

1.943 1.794 1.612 1.364

HB

1.0 0.8 0.5 0

1.926 1.851 1.754 1.615

5.472 4.896 4.165 3.184

1.311 1.295 1.273 1.238

2.191 2.060 1.870 1.574

1.895 1.847 1.750 1.614

5.310 4.902 4.168 3.187

1.334 1.300 1.269 1.237

2.250 2.090 1.862 1.573

Note: EBP denotes the elasto-brittle-plastic behavior, and EP denotes the elasto-perfectly-plastic behavior. 8

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

(a) 0.4

0.5 Proposed solution Solution of Zhang et al. (2012)

0.3

R1=3.0, E=10 GPa, μ=0.25, p′0=20 MPa

0.2

ψn=ψ1=3.75°, γ =0.008

Proposed solution Solution of Zhang et al. (2012)

0.3

σ1′ /p0′

0.4

σ1′ /p′0

(a)

cn=1.0 MPa, c1=0.7 MPa, ϕn=2.0, ϕ1=0.6 p*

′ R1=3.0, E=5.7 GPa, μ=0.25, p0=15 MPa

σcn=30 MPa, σc1=25MPa, mn=2.0, m1=0.6

0.2

-3

p*

0.1

0.1

0.0 0

5

10

15

20

25

30

0.0

35

0

2

4

6

0.5

(b)

0.4

Proposed solution Solution of Zhang et al. (2012)

10

12

14

0.4

0.3

Proposed solution Solution of Zhang et al. (2012)

σ1′ /p′0

σ1′ /p0′

0.3

0.2

8

2Gu1/[R1(p0′-σp′)]

2Gu1/[R1(p0′-σp′ )]

(b)

-3

sn=4×10 , s1=2×10 , ψn=15°,ψ1=5°, γ =0.01

0.2

Rp/R1

0.1

0.1 Rs/R1 0.0

1.0

1.5

2.0

2.5

3.0

3.5

4.0

0.0

4.5

Rs/R1, Rp/R1

1.0

1.3

1.6

1.9

2.2

2.5

2.8

Rs/R1, Rp/R1

Fig. 6. Validation for MC strain-softening rock mass: (a) GRC, and (b) evolution of Rp and Rs.

Fig. 7. Validation for HB strain-softening rock mass: (a) GRC, and (b) evolution of Rp and Rs.

the dichotomy method. The detailed flowchart that summarizes the implementation of the proposed solutions is shown in Fig. 4.

Table 3 Hydraulic-mechanical solutions of a circular tunnel in a strain-softening rock mass with an initial pore-water pressure p0w = 3 MPa (data in parenthesis are solutions by Wang et al. (2010)).

6. Example analysis

Rock type

6.1. Verifications for the elasto-brittle (perfectly)-plastic rock mass To verify the validity of the proposed solutions, two groups of typical input data for both the MC and HB rock masses that were used by many researchers (Park et al., 2008; Wang et al., 2010; Zhang et al., 2018) are considered. The input data are listed in Table 1, and the confining pressure-dependent effect is not considered. To guarantee the calculation precision, the concerned region of [R1, 10R1] is discretized into n1 = 1000 annuli and that of [10R1, Ru] is n2 = 5000 annuli. An acceptable absolute tolerance of 5 Pa for both the pore-water pressure and internal effective stress is used in the following analysis. The elastobrittle (perfectly)-plastic behavior is a special case of strain-softening behavior corresponding to a very small (large) critical plastic shear strain. Therefore, γp* = 10−10 and 10 are used for the example analysis of the elasto-brittle-plastic and elasto-perfectly plastic rock masses, respectively. The numerical procedure of the Fast Lagrangian Analysis of Continua (FLAC) is employed in the hydraulic-mechanical calculation. Because the commercial FLAC procedure cannot consider the vibration of permeability induced by fracture, a steady hydraulic-mechanical behavior with constant permeability is considered. The calculated results are compared with numerical results of FLAC and the analytical solutions, as listed in Table 2. The maximum gaps between

β

Proposed method

Numerical method

Rs/R1

Rp/R1

u1 E0 p′0 R1

Rs/R1

Rp/R1

u1 E0 p′0 R1

MC

1.0 0.8 0.5 0

3.820 3.577 3.261 2.829

5.450 5.118 4.688 4.109

38.623 33.360 27.121 19.721

3.644 3.475 3.184 2.732 (2.802)

5.126 4.907 4.322 3.719 (4.093)

41.055 35.875 29.175 19.025 (19.568)

HB

1.0 0.8 0.5 0

2.861 2.659 2.409 2.0790

3.577 3.334 3.033 2.639

20.835 17.573 13.913 9.762

2.530 2.391 2.164 1.936 (2.070)

2.874 2.718 2.432 2.204 (2.638)

18.377 15.474 11.978 8.132 (9.772)

the proposed results and numerical (analytical) results are 1.75% and 2.96% for the plastic radius and displacement, respectively. Taking the MC elasto-brittle-plastic rock mass as an example, the distribution of the tangential stress with different Biot’s coefficients and pore-water pressures is shown in Fig. 5, which also denotes that the results of the proposed method are almost the same with those of the numerical results by FLAC. Therefore, the proposed solutions are in good agreement with the numerical and analytical results. 9

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

is carried out according to the prior obtained pore-water pressure. The material property parameters are visualized in Figs. 6 and 7 for the MC and HB rock masses, respectively. Figs. 6 and 7 show the comparisons between the proposed solutions and analytical solutions with ground response curves (GRCs), evolutions of plastic and residual radii for the drained MC and HB rock masses, respectively. The proposed solutions are in good agreement with those of the numerical solutions by Wang et al. (2010), with a maximum relative error of 5.72% and 1.70% for plastic radius and displacement, respectively. The HMSS solutions of the above three methods with various Biot’s coefficients are listed in Table 3. In addition to the solutions of Wang et al. (2010), the relative error between the proposed results and numerical results of FLAC are slightly larger with maximum relative errors of 19.82% and 13.91% for the plastic radius and displacement, respectively. This is induced by the localization effect of strain-softening behavior in the numerical simulation method when the element size is very small (Varasa et al., 2005). However, a large element size of FLAC also leads to a large gap with the analytical results. This is the disadvantage of finite element method when it is used for strain-softening behavior analysis. The proposed rigorous HMSS semi-analytical solution can provide a benchmark for these numerical simulations in strain-softening analysis.

Table 4 Geometric and material property parameters of the fully coupled HMSS rock mass. M-C rock mass Value

H-B rock mass

Parameters

Values

Parameters

Values

K0/GPa k0, kσ, k γ p * /m∙s−1

K0/GPa k0, kσ, k γ p * /m∙s−1

E0, Eσ, E γ p * /GPa

6.67 10−7, 0.5 × 10−7, 10−6 10, 12, 7

E0, Eσ, E γ p * /GPa

3.8 10−7, 0.5 × 10−7, 10−6 5.7, 6.84, 3.99

μ0, μσ, μ γ p * (-)

0.25, 0.30, 0.175

μ0, μσ, μ γ p * (-)

0.25, 0.30, 0.175

c0, cσ, c γ p * /MPa

1.0, 0.8, 0.7

σ0, σσ, σ γ p * /MPa

30, 25, 22

φ0, φσ, φ γ p * /Deg

30, 25, 22

m0, mσ, m γ p * (-)

2.0, 1.5, 0.6

ψ0, ψσ, ψ γ p * /Deg

3.75, 3.2

s0, sσ, s γ p * /0.001

4, 3, 2

γ0p ∗ ζ

0.01

ψ0, ψσ, ψ γ p * (Deg)

15, 10, 5

0.375

γ0p * ζ

0.02 0.28

6.2. Verifications for conventional HMSS rock mass In the conventional HMSS analysis, both the Biot’s coefficient and permeability are assumed to be constant during the failure process. The most-used numerical simulation procedure of FLAC and analytical method are employed to validate the proposed HMSS solutions. An initial pore-water pressure p0w = 3.0 MPa and internal pore-water pressure p1w = 0 MPa are considered. Because the failure of the rock mass has no influence on the permeability of the rock mass, the hydraulic behavior is first calculated to obtain the pore-water pressure distribution and water-inflow. Then, the mechanical behavior analysis

(b)

0.6 0.5

0.2

0.1

0.1

0

5

10

15

20 25 u1E/R1 p0′

30

35

0.0

40

1.0

1.5

2.0

2.5

3.0 3.5 Rp/R1

(d) 0.6

0.4

0.5

p0w=0 MPa pw0 =1.0 MPa pw0 =3.0 MPa pw0=5.0 MPa p0w=7.0 MPa p0w=0, Results of Zhang et al.(2018)

0.3

σ1′ /p0′

0.3

0.2

0.0

p0w=0 MPa pw0 =1.0 MPa pw0 =3.0 MPa pw0=5.0 MPa p0w=7.0 MPa p0w=0, Results of Zhang et al.(2018)

0.4 σ1′ /p0′

σ1 ′/p0′

0.3

0.6 0.5

p0w=0 MPa pw0 =1.0 MPa pw0 =3.0 MPa pw0=5.0 MPa p0w=7.0 MPa p0w=0, Results of Zhang et al.(2018)

0.4

(c)

For the proposed fully coupled HMSS rock mass analysis, two data sets from Zhang et al. (2018) in the elasto-plastic coupling strain-softening analysis are considered. The geometric parameters are the same as those presented in Section 6.2, and the mechanical and hydraulic

0.2

pw0 =1.0 MPa

3.0 MPa

4.0

5.0 MPa

4.5

5.0

5.5

7.0 MPa

0.4 σ1′ /p0′

(a)

6.3. Examples of fully coupled HMSS rock mass

0.3 0.2

0.1

0.1 0.0

1.0

1.5

2.0

2.5

3.0

3.5

0.0

4.0

Rs/R1

0.0

0.1

0.2

0.3

0.4

0.5

0.6

Q (10-5m3/m⋅s)

Fig. 8. Hydraulic-mechanical solutions of the MC strain-softening rock mass: (a) GRC, (b) evolution of Rp, (c) evolution of Rs, and (d) evolution of water inflow. 10

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

(b)

0.6 pw0 =0 MPa pw0 =1.0 MPa pw0 =3.0 MPa pw0=5.0 MPa pw0 =7.0 MPa pw0 =0, Results of Zhang et al.(2018)

0.5

σ1′ /p0′

0.4 0.3

0.1

0.1 0.0

4

8

12

16

20

24

(d)

0.30 pw0 =0 MPa pw0 =1.0 MPa pw0 =3.0 MPa pw0=5.0 MPa p0w=7.0 MPa pw0 =0, Results of Zhang et al.(2018)

0.20 0.15

1.9

3.5

4.0

2.2

2.5

7.0 MPa

MPa

3.0 MPa

0.3

0.1

1.6

3.0

0.4

0.05

1.3

2.5

5.0 MPa pw0 =1.0

0.2

1.0

2.0

0.6

0.10

0.00

1.5

0.5

σ1′ /p0′

0.25

1.0

Rp/R1

u1E/R1 p0′

σ1′ /p0′

0.3 0.2

0

p0w=0 MPa pw0 =1.0 MPa pw0 =3.0 MPa pw0=5.0 MPa p0w=7.0 MPa p0w=0, Results of Zhang et al.(2018)

0.4

0.2

0.0

(c)

0.6 0.5

σ1′ /p0′

(a)

2.8

Rs/R1

0.0

0

1

2

3 -5

4

5

6

3

Q (10 m /m)

Fig. 9. Hydraulic-mechanical solutions of the HB strain-softening rock mass: (a) GRC, (b) evolution of Rp, (c) evolution of Rs, and (d) evolution of water inflow.

mass. With the decrease in internal pressure, Rp/R1, Rs/R1 and u1 E0 p′0 R1 gradually increase and reach their maximum values whenσ1′ = 0. With the increase in p0w , the variables of displacement, residual and plastic radii and water-inflow gradually increase, and the increasing rates also increase. Taking the condition of p0w = 5 MPa of the MC rock mass as an example, when the internal pressure decreases from 0.2p0′ to 0.1p0′, Rp/R1, Rs/R1, u1 E0 p′0 R1 and Q increase from 1.19, 1.82, 2.91, 3.35 × 10−5 m3/m∙s to 1.74, 2.65, 7.22, 3.63 m3/m∙s, with increasing rates of 46.21%, 45.70%, 147.85% and 8.29%, and the increasing rates of 93.78%, 94.09%, 332.46% and 17.14% correspond to the internal pressure decreasing from 0.1p0′ to 0. An approximate linear increase of Rp/R1, Rs/R1, u1 E0 p′0 R1 and Q with the increasing p0w can be observed in Fig. 10. Compared to the completely drained conditions, Rp/R1, Rs/R1 and u1 E0 p′0 R1 of p0w = 0.7 MPa increase by 18.86%, 20.40% and 58.44% for the MC rock mass, and those of the HB rock mass are 22.52%, 27.93% and 93.12%. The seepage force significantly enlarges the plastic radius and displacement of the surrounding rock. Additionally, the water-inflow also significantly increases with the increasing p0w for both MC and HB rock masses. Fig. 11 shows the influences of Ks on Rp/R1, Rs/R1, u1 E0 p′0 R1 and Q for the MC and HB rock masses. The increasing Ks enlarges the above four variables in an approximate exponent form with a limit value. The increasing rates of above four variables decrease, and gradually trend to nil. For the MC rock mass, Rp/R1, Rs/R1, u1 E0 p′0 R1 and Q increase from 3.01, 4.56, 22.30, 2.53 × 10−5 m3/m∙s to 3.25, 4.96, 28.43 and 2.54 × 10−5 m3/m∙s with increasing rates of 8.08%, 8.64%, 23.62% and 0.26% when Ks increases from 1.5 K0 to 10 K0. However, those of 0.83%, 0.92%, 2.34% and 0.08% correspond to the increase of Ks from 10 K0 to 105K0. For the HB rock mass, the above variables increase by 9.36%, 11.94%, 36.61% and 0.41% when Ks increases from 1.5 K0 to

Table 5 Calculated results of circular tunnels in the fully coupled HMSS rock mass. Rock type

MC

HB

Scheme

p0w (MPa) 0

1.0

3.0

5.0

7.0

Rs/R1 Rp/R1 u1 E0 p′0 R1

2.97 4.51 22.40

3.05 4.63 23.98

3.21 4.88 27.40

3.37 5.15 31.22

3.53 5.43 35.49

Q (10−5 m3/s)

0

0.84

2.54

4.25

6.00

Rs/R1 Rp/R1 u1 E0 p′0 R1

2.22 2.90 12.21

2.28 3.00 13.45

2.41 3.22 16.25

2.55 3.45 19.58

2.72 3.71 23.58

0

0.77

2.31

3.88

5.48

Q (10−5 m3/m∙s)

parameters are listed in Table 4. Figs. 8 and 9 show the evolutions of GRC, Rp, Rs and Q as the decreasing internal pressure under various initial pore-water pressures for the MC and HB rock masses, respectively. Obviously, similar evolutions of Rp, Rs and GRC of drained conditions can be observed. For the drained cases, the Rp/R1, Rs/R1 and u1 E0 p′0 R1 of the proposed method are 2.97, 4.51, 22.40 and 2.22, 2.90, 12.21 for the MC and HB rock masses, respectively, and those proposed by Zhang et al. (2018) are 2.98, 4.50, 21.29 and 2.14, 2.84, 11.09, respectively, as listed in Table 5. Although the relative error is 9.17% for the displacement for the HB rock mass, the relative errors for Rp, Rs and u1 between the above two methods are all less than 5%. Thus, the proposed solutions are in accordance with the analytical solutions proposed by Zhang et al. (2018). Figs. 8 and 9 imply that the GRC, Rp, and Rs of the HMSS rock mass are similar to those of the completely drained strain-softening rock 11

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

(a)

40

5.5 5.0

Rs

u1

Rp

Q

(a)

12 10

35

0.6

Ks/K0=5, kσ/k0=0.5, kγ /k0=10 p

Ks/K0=5, kσ/k0=1.0, kγ /k0=10 p

0.5

Ks/K0=5, kσ/k0=0.5, kγ /k0=1.0 p

20

3.0

0

1

2

3

4

5

6

7

15

w

0

w

0.2

0

0.1 0.0

6

30 u1

Rp

Q

(b)

3

5

7

11

13

15

0.5

5

Ks/K0=10 , kσ/k0=0.5, kγ /k0=10 p

3

0.4

-5

3 2

Ks/K0=1.5, kσ/k0=0.5, kγ /k0=10

w

0.3 1

0

w

10

9

p /p

2.5

1

p

Q (10 m /m⋅s)

15

u1E0/R1 p0′

Rs/R1, Rp/R1

3.0

p

r/R1

4 20

Ks/K0=1.5, kσ/k0=0.5, kγ /k0=10

5

25

3.5

5

Ks/K0=10 , kσ/k0=0.5, kγ /k0=10 p

p (MPa) Rs

0.3

2

w 0

(b) 4.0

p

3

4

3.5

2.5

6

p /p

25

Ks/K0=5, kσ/k0=1.0, kγ /k0=1.0

0.4

-5

4.0

u1E/R1 p0′

Rs/R1, Rp/R1

30

Q (10 m /m⋅s)

8

4.5

2.0

0

1

2

3

4

5

6

7

5

0.2

0

Ks/K0=5, kσ/k0=0.5, kγ /k0=10 p

p0w (MPa)

Ks/K0=5, kσ/k0=1.0, kγ /k0=10

0.1

p

Ks/K0=5, kσ/k0=0.5, kγ /k0=1.0

Fig. 10. Influence of initial pore-water pressure on the plastic radius, displacement and water inflow: (a) MC rock mass, and (b) HB rock mass.

p

0.0 u1

Rp

Q

34

2.555

4.5

28

2.545

4.0

25

3.5

22

2.535

19

2.530

1

10

2

3

10

4

10

5

10

20

2.330

19

2.325

15

2.315 2.310

3

16 2.6

-5

2.320

17

2.9

Q (10 m /m⋅s)

18

u1E0/R1 p0′

Rs/R1, Rp/R1

3.2

14

2.3

1

5

7

9

11

13

15

10 K0. When Ks increases from 10 K0 to 105K0, the increasing rates are 1.03%, 1.31%, 3.53% and 0.12%. Thus, Ks/K0 mainly influences Rp/R1, Rs/R1 and u1 E0 p′0 R1, and has a little influence on Q, which mainly depends on the permeability and pore-water pressure. Meanwhile, Ks/ K0 mainly depends on the porosity of the rock mass. For the intact rock mass, K0 almost equals to Ks, and β trends to zero. The pore-water pressure can be considered as an initial stress imposed on the outer boundary. However, for the porous rock mass, Ks is much larger than K, leading to β approaching to 1.0. In addition, this is the same with the conventional Darcy effective stress solution. Fig. 12 shows the distribution of pore-water pressure under various Ks and k for the MC and HB rock masses. A large permeability leads to a small change in the pore-water pressure with the increasing r. In the residual and elastic regions, the permeability is governed by the confining pressure. With decreasing kσ/k0 from 1.0 to 0.5, the permeability decreases due to the increasing confining pressure, and the difference in the pore-water pressure between the two conditions gradually increases because of the increasing radial stress as radius r. The significant increase in the permeability of the plastic rock mass increases the porewater pressure in the plastic region around the roadway. However, the increase in Ks only induces a slight increase in pw. Taking the MC rock mass as an example, pw p0w = 0.016 and 0.113 at r = 3R1 and 7R1, respectively, which corresponds to the case of Ks/K0 = 5, kσ/k0 = 0.5 andk γ p ∗ k 0 = 10. pw p0w = 0.028 and 0.152 at r = 3R1 and 7R1 correspond to the case of Ks/K0 = 5, kσ/k0 = 1.0 and k γ p ∗ k 0 = 10. However, when the changes in the permeability of plastic surrounding rock are ignored, pw p0w rapidly increases in the plastic region. For

-5

2.540

10

(b) 3.5

10

3

Fig. 12. Distribution of pore-water pressure ( p0w =3 MPa): (a) MC rock mass, and (b) HB rock mass.

Ks/K0

2.0 0 10

1

3

2.550

3.0 0 10

p

r/R1 31

5.0

Rs/R1, Rp/R1

Rs

Q (10 m /m⋅s)

5.5

u1E0/R1 p0′

(a)

Ks/K0=5, kσ/k0=1.0, kγ /k0=1.0

2

3

10

10

Rs

u1

Rp

Q 4

10

13 12

5

2.305 2.300

10

Ks/K0

Fig. 11. Influence of β on the plastic radius, displacement and water inflow ( p0w =3 MPa): (a) MC rock mass, and (b) HB rock mass.

12

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

(a) 1.6 1.4

Softening region

β /β0

1.2 β /β0, K/K0

12

7. Conclusions

10

Considering the evolutions of permeability and Biot’s coefficient during the progressive failure process, this paper proposed a fully coupled HMSS model of rock mass, in which both the hydraulic and mechanical property parameters were taken into account as variables of plastic shear strain and confining pressure. Assuming that the material parameters are constants in a very small region, a recursive method was employed to derive the hydraulic-mechanical solutions according to the ideal elasto-plastic hydraulic-mechanical solution for both MC and HB rock masses. The proposed solutions were validated by the closed-form solutions of the completely drained surrounding rock and numerical simulation solutions of water-rich surrounding rock around a circular tunnel. With the decrease in internal pressure, the residual and plastic radii, displacement and water-inflow significantly increase, and the increasing rates gradually increase. With the increase in the initial porewater pressure, the residual and plastic radii, displacement and waterinflow gradually increase in approximately linear forms. The increasing bulk modulus Ks of the solid constituent enlarges the residual and plastic radii, displacement and water-inflow, and they gradually trend toward the maximum limit values, which agree with the solutions of Darcy’s effective stress with β = 1.0. The pore-water pressure distribution mainly depends on the permeability of the surrounding rock. A large permeability of plastic rock mass decreases the pore-water pressure, and the increasing Ks only indirectly enlarges the pore-water slightly due to the plastic region. For the proposed fully coupled hydraulic-mechanical examples, the hydraulic parameters are only influenced by the confining pressure in the residual and elastic regions, which leads to slight changes. However, a significant variation in hydraulic parameters occurs in the softening region, and this is in accordance with the experimental results. Compared to the conventional HB rock mass, the increase in the third strength parameter α more significantly enlarges the plastic radius and displacement of the surrounding rock. The proposed fully coupled HMSS solutions are an extension of the conventional hydraulic-mechanical solution, in which the hydraulic parameters are considered to be constants.

Elastic region

8

k/k0 K/K0

1.0

6

0.8

4

0.6

2

0.4

1

2

3

4

5

6

7

8

9

k/k0

Residual region

0

r/R1 1.6

12 Residual region

1.4

Softening region

Elastic region

10 β /β0

β /β0, K/K0

1.2

8

k/k0 K/K0

1.0

6

0.8

4

0.6

2

0.4

1

2

3

4

5

6

k/k0

(b)

0

r/R1 Fig. 13. Evolutions of k and β around the circular tunnel ( p0w =3 MPa): (a) MC rock mass, and (b) HB rock mass.

example, pw p0w = 0.130 and 0.254 at r = 3R1 and 7R1, respectively, for the conditions Ks/K0 = 5, kσ/k0 = 0.5 and k γ p ∗ k 0 = 1 of MC rock mass. The changes in the pore-water pressure with Ks/K0 are induced by the change of the plastic region, in which the permeability increases. Thus, the Ks/K0 only influences the pore-water pressure slightly. Fig. 13 shows the distributions of k, β and K along the radius r for the proposed MC and HB HMSS rock masses with Ks = 5 K0. For the employed data sets, k and β gradually decrease, whereas K gradually increases. In the residual region, the above three variables are only governed by the confining pressure, and only slight decreases in k/k0 and β/β0 of 13.30% and 5.17% occur when r increases from R1 to Rs. When r increases from Rs to Rp, k/k0 and β/β0 significantly decrease by 1080.39% and 17.34%. The values of K/K0 of the surrounding rock are 0.54, 0.59 and 1.24 at r = R1, Rs and Rp, respectively. Similar evolutions of k, β and K can be found for the HB rock mass, as shown in Fig. 13(b). Thus, the fracture of the surrounding rock significantly increases the permeability of the plastic rock mass, and the hydraulic parameters are significantly changed with the increasing plastic strain in the plastic region. The strength parameter α has a significant influence on the residual and plastic radii and displacement. For the proposed HMSS rock mass, Rp/R1, Rs/R1, u1 E0 p′0 R1 and Q of the generalized HB criterion, in which ασ/α0 = 1.1 and α γ p ∗ α 0 = 1.2, are 3.97, 5.34, 49.25 and 2.62 × 10−5 m3/m∙s. With respect to those of the conventional HB criterion of ασ/ α0=α γ p ∗ α 0 = 1.0, the above four variables increase by 64.73%, 65.84%, 203.08% and 13.42%, respectively.

CRediT authorship contribution statement Qiang Zhang: Conceptualization, Methodology, Software, Investigation, Writing - original draft, Data curation, Writing - review & editing. Cong Shao: Validation, Supervision, Writing - review & editing, Visualization. Hong-Ying Wang: Validation, Formal analysis, Visualization. Bin-Song Jiang: Resources, Writing - review & editing, Supervision. Yu-Jing Jiang: Resources, Writing - review & editing, Supervision. Ri-Cheng Liu: Resources, Supervision, Data curation, Writing - review & editing. Declaration of Competing Interest The authors declared that they have no conflicts of interest to this work. Acknowledgements This study has been partially supported by the Fundamental Research Funds for the Central Universities, China (Grant No. 2017XKQY048). References Alejano, L.R., Alonso, E., 2005. Considerations of the dilatancy angle in rocks and rock masses. Int. J. Rock Mech. Min. Sci. 42, 481–507. Alejano, L.R., Alonso, E., Rodriguez-Dono, A., Fernandez-Manin, G., 2010. Application of the convergence-confinement method to tunnels in rock masses exhibiting Hoek-

13

Tunnelling and Underground Space Technology 99 (2020) 103375

Q. Zhang, et al.

fractured rock masses: a numerical study. Int. J. Rock Mech. Min. Sci. 41 (2004), 1191–1210. Ming, H.F., Wang, M.S., Tan, Z.S., Wang, X.Y., 2010. Analytical solution for steady seepage into an underwater circular tunnel. Tunel. Undergr. Space Technol. 25, 391–396. Nur, A., Byerlee, J.D., 1971. An exact effective stress law for elastic deformation of rocks with fluids. J. Geophys. Res. 76, 6414–6419. Park, K.H., Kim, Y.J., 2006. Analytical solution for a circular opening in an elastic- brittleplastic rock. Int. J. Rock Mech. Min. Sci. 43, 616–622. Park, K.H., Lee, J.G., Owatsiriwong, A., 2008a. Seepage force in a drained circular tunnel: An analytical approach. Can. Geotech. J. 45 (3), 432–436. Park, K.H., Tontavanich, B., Lee, J.G., 2008b. A simple procedure for ground response curve of circular tunnel in elastic-strain softening rock masses. Tunn. Undergr. Space Technol. 23, 151–159. Schleiss, A., 1986. Design of previous pressure tunnels. Int. Water Power Dam Construct. 38 (5), 21–26. Sharan, S.K., 2003. Elastic-brittle-plastic analysis of circular openings in Hoek-Brown media. Int. J. Rock Mech. Min. Sci. 40, 817–824. Timoshenko, S.P., Gooodier, J.N., 1970. Theory of Elasticity. MCgraw-Hill Book Company, New York. Varasa, F., Alonso, E., Alejano, L.R., Fdez.-Manína, G., 2005. Study of bifurcation in the problem of unloading a circular excavation in a strain-softening material. Tunn. Undergr. Space Technol. 20, 311–322. Wang, S.G., Elsworth, D., Liu, J.S., 2003. Permeability evolution during progressive deformation of intact coal and implications for instability in underground coal seams. Int. J. Rock Mech. Min. Sci. 58, 34–45. Wang, S.L., Yin, X.T., Tang, H., Ge, X.R., 2010. A new approach for analyzing circular tunnel in strain-softening rock masses. Int. J. Rock Mech. Min. Sci. 47, 170–178. Xu, P., Yang, S.Q., 2016. Permeability evolution of sandstone under short-term and longterm triaxial compression. Int. J. Rock Mech. Min. Sci. 85, 152–164. Yang, S.Q., Huang, Y.H., Jiao, Y.Y., Zeng, W., Yu, Q.L., 2015. An experimental study on seepage behavior of sandstone material with different gas pressures. Acta Mech. Sin. 31 (6), 837–844. Zareifard, M.R., Fahimifar, A., 2014a. Elastic-brittle-plastic analysis of circular deep underwater cavities in a Mohr-Coulomb rock mass considering seepage forces. Int. J. Geomech. 04014077. Zareifard, M.R., Fahimifar, A., 2014b. Effects of seepage forces on circular openings excavatedin Hoek-Brown rock mass based on a generalized effective stress principle. Eur. J. Environ. Civ. Eng. 18, 584–600. Zhang, K., Zhou, H., Shao, J.F., 2013. An experimental investigation and an elastoplastic constitutive model for a porous rock. Rock Mech. Rock Eng. 46 (6), 1499–1511. Zhang, Q., Jiang, B.S., Wang, S.L., Ge, X.R., Zhang, H.Q., 2012a. Elasto-plastic analysis of a circular opening in strain-softening rock mass. Int. J. Rock Mech. Min. Sci. 50, 38–46. Zhang, Q., Jiang, B.S., Wu, X.S., Zhang, H.Q., Han, L.J., 2012b. Elasto-plastic coupling analysis of circular openings in elasto-brittle-plastic rock mass. Theor. Appl. Fract. Mech. 60, 60–67. Zhang, Q., Wang, H.Y., Jiang, Y.J., Lu, M.M., Jiang, B.S., 2019. A numerical large strain solution for circular tunnels excavated in strain-softening rock masses. Comput. Geotech. 114, 103142. Zhang, Q., Zhang, C.H., Jiang, B.S., Li, N., Wang, Y.C., 2018. Elastoplastic coupling solution of circular openings in strain-softening rock mass considering pressure-dependent effect. Int. J. Geomech. 18 (1), 04017132. Zhao, X.G., Cai, M., 2010. A mobilized dilation angle model for rocks. Int. J. Rock Mech. Min. Sci. 47 (3), 368–384. Zou, J.F., Li, S.S., Xu, Y., Dan, H.C., Zhao, L.H., 2016. Theoretical solutions for a circular opening in an elastic–brittle–plastic rock mass incorporating the out-of-plane stress and seepage force. KSCE J. Civ. Eng. 20 (2), 687–701.

Brown strain-softening behaviour. Int. J. Rock Mech. Min. Sci. 47, 150–160. Bear, J., 1972. Dynamics of Fluids in Porous Media. American Elsevier Publishing Company, New York. Biot, M.A., 1941. General theory of three-dimensional consolidation. J. Appl. Phys. 12, 155–164. Bouvard, M., Pinto, N., 1969. Amenagement capivari-cacheoeira. Etude du puits en charge. La Houille Blanche, Paris, France, pp. 747–760 (In French). Brown, E.T., Bray, J.W., Landayi, B., Hoek, E., 1983. Ground response curves for rock tunnels. J. Geotech. Eng. 109 (1), 15–39. Brown, E.T., Bray, J.W., 1982. Rock-support interaction calculations for pressure shafts and tunnels. In: Proceedings of ISRM Symposium, Azachen, pp. 555–565. Brown, E.T., Bray, J.W., Santarelli, F.J., 1989. Influence of stress-dependent elastic moduli on stresses and strains around axisymmetric boreholes. Rock Mech. Rock Eng. 22, 189–203. Burden, R.L., Faires, J.D., 2010. Numerical Analysis, 9th ed. Brooks/Cole, Cengage Learning, Boston. Carosso, G., Giani, G. P. (1989) Analytical solutions for potentials and stresses around a cavity under a water table. In: Serrano (Ed.), Tunnels and Water. Balkema, Rotterdam, pp. 1209-1217. Chen, L., Liu, J.F., Wang, C.P., Liu, J., Su, R., Wang, J., 2014. Characterization of damage evolution in granite under compressive stress condition and its effect on permeability. Int. J. Rock Mech. Min. Sci. 71, 340–349. Chen, Y.F., Zhou, C.B., Sheng, Y.Q., 2007. Formulation of strain-dependent hydraulic conductivity for a fractured rock mass. Int. J. Rock Mech. Min. Sci. 44, 981–996. Cui, L., Zheng, J.J., Zhang, R.J., Dong, Y.K., 2015. Elasto-plastic analysis of a circular opening in rock mass with confining stress-dependent strain-softening behavior. Tunn. Undergr. Space Technol. 50, 94–108. Fahimifar, A., Zareifard, M.R., 2013. A new closed-form solution for analysis of unlined pressure tunnels under seepage forces. Int. J. Numer. Anal. Meth. Geomech. 37 (11), 1591–1613. Fahimifar, A., Zareifard, M.R., 2009. A theoretical solution for analysis of tunnels belowgroundwater considering the hydraulic-mechanical coupling. Tunn. Undergr. Space Technol. 24, 634–646. Fazio, L.A., Ribacchi, R., 1984. Influence of seepage on tunnel stability. In: Proceedings of ISRM Symposium on Design and Performance Underground Excavations, British Geotechnical Society, Cambridge, pp. 173–181. Fernandez, G., Alvarez, T.A., 1994. Seepage-induced effective stresses and water pressures around pressure tunnels. J. Geotech. Eng. 120 (1), 108–128. Han, J.X., Li, S.S., Li, S.C., Yang, W.M., 2013. A procedure of strain-softening model for elasto-plastic analysis of a circular opening considering elasto-plastic coupling. Tunn. Undergr. Space Technol. 37, 128–134. Harr, E.E., 1962. Groundwater and Seepage. McGraw-Hill, New York. Hu, D.W., Zhou, H., Zhang, F., Shao, J.F., 2010. Evolution of poroelastic properties and permeability in damaged sandstone. Int. J. Rock Mech. Min. Sci. 47, 962–973. Kolymbas, D., Wagner, P., 2007. Groundwater ingress to tunnels- The exact analytical solution. Tunn. Undergr. Space Technol. 22 (1), 23–27. Lee, Y.K., Pietruszczak, S., 2008. A new numerical procedure for elasto-plastic analysis of a circular opening excavated in a strain-softening rock mass. Tunn. Undergr. Space Technol. 23 (5), 588–599. Li, Z.L., Ren, Q.W., Wang, Y.H., 2004. Elasto-plastic analytical solution of deep-buried circle tunnel considering fluid flow field. Chin. J. Rock Mech. Eng. 23 (8), 1291–1295. Liu, R., Li, B., Jiang, Y.J., 2016. A fractal model based on a new governing equation of fluid flow in fractures for characterizing hydraulic properties of rock fracture networks. Comput. Geotech. 75, 57–68. Martin, C.D., Chandler, N.A., 1994. The progressive fracture of Lac du Bonnet granite. Int. J. Rock Mech. Min. Sci. Geomech. Abstr. 31, 643–659. Min, K.B., Rutqvistb, J., Tsang, C.F., Jing, L.R., 2004. Stress-dependent permeability of

14