- Email: [email protected]
Elasto-plastic coupling analysis of circular openings in elasto-brittle-plastic rock mass
Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
Contents lists available at SciVerse ScienceDirect
Theoretical and Applied Fracture Mechanics journal homepage: www.elsevier.com/locate/tafmec
Elasto-plastic coupling analysis of circular openings in elasto-brittle-plastic rock mass Q. Zhang ⇑, B.S. Jiang, X.S. Wu, H.Q. Zhang, L.J. Han State Key Laboratory for Geomechanics and Deep Underground Engineering, School of Mechanics and Civil Engineering, China University of Mining and Technology, Xuzhou 221008, PR China
a r t i c l e
i n f o
Article history: Available online 28 June 2012 Keywords: Circular opening Cyclic loading and unloading Elasto-plastic coupling Analytical solution FLAC3D
a b s t r a c t This paper deals with elasto-plastic coupling solutions for the prediction of displacements around circular openings in elasto-brittle-plastic rock mass subject to hydrostatic stress. Both linear Mohr–Coulomb criterion and nonlinear Hoek–Brown criterion on the basis of a non-associated flow rule are considered. The restrictive condition between strength parameters and Young’s modulus in post-failure region for rock materials is theoretically established. Meanwhile, the proposed restrictive condition is validated by FLAC3D code that the advanced elasto-plastic coupling model (EPCM) based on Mohr–Coulomb criterion is merged in. Finally, the elasto-plastic coupling dimensionless displacements of circular opening in elasto-brittle-plastic rock mass are analyzed with different deterioration degree of Young’s modulus using two kinds of rocks. The results show that the deformation of rock mass increases obviously with the decrease of Young’s modulus in the plastic region, but Young’s modulus has little influence on the plastic region radius and stresses distribution form. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction Although plane strain problems of circular openings in rock mass are such simple cases, they can provide theoretical foundation for geotechnical and mining engineering problems, such as the optimization design of tunnels, boreholes and mine shafts. Lots of researchers show special interest on this axisymmetric problem and have done a great deal of contributions on the closed-form solutions for this problem [1–10]. As for most of the analyses reported in the past, the solutions are obtained considering different models of material behavior, such as ideal elasto-plastic, elasticbrittle-plastic and strain-softening models, with different yield criteria, like the linear Mohr–Coulomb (M–C) and nonlinear Hoek–Brown (H–B) criteria [1]. The stress state and deformation of M–C rock mass with a circular hole were earlier studied [2]. Until to 1983, the closed-form solutions for the stress and radial displacement of circular openings was obtained on the assumption of constant elastic strain in the plastic region, which is equal to that at the elasto-plastic interface [3]. Later, a numerical method for elasto-plastic rock mass was proposed avoiding the defects of assumption of constant elastic strain for plastic rock mass [4]. In order to avoid the difficulty, the so-called self-similarity theory combined with analytical method was adopted as an effective way to predicted rock mass displacements and stresses in ⇑ Corresponding author. Tel.: +86 13641535671. E-mail address: [email protected] (Q. Zhang). 0167-8442/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.tafmec.2012.06.008
elasto-plastic H–B rock mass [5]. Obviously, the above achievements are all based on the ideal elasto-plastic model, however, most of the brittle rock mass shows sudden loss of strength when the stress arrived at its peak strength. Considering the brittle failure characteristics, lots of scholars have also done further investigation. By assuming four types of strain distribution forms in the plastic region, the closed-form solutions of circular openings were analyzed with both the M–C and H–B yield criteria. As we know, the elastic strain should obey the Hooke’s law, even in the plastic region of surrounding rock. So only the third kind of solution is correct [8]. In order to obtain the explicit expressions for displacements, it was assumed that the elastic strain in the plastic region obeys the same law in thick wall cylinder condition. So, the solutions for generalized H–B rock mass are also not exactly coherent with the real displacements [9–11]. Meanwhile, the strain-softening calculation procedure using numerical method [12–14] and semi-analytical and semi-numerical method [15,16] were proposed, and the examples are well consistent with each other. The material properties depend on many factors, such as the stress state and extent of damage. The deterioration of strength parameters are commonly considered in the former studies, but the Young’s modulus was rarely considered. After the excavation of tunnels, the stress state of surrounding rock turns from three dimensional stress state to quasi-two dimensional state. The deterioration of confining condition and extent of broken would lead to the decrease of Young’s modulus. Additionally, numbers of traditional compression tests show that the Young’s modulus of rock
61
Maximum deviatoric stress (MPa)
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
P0
P0
Pi
R0
Axial strain (-)
Rp
Fig. 1. Cyclic loading and unloading test under lower confining pressure.
P0 Plastic region
2. Definition of the problem Fig. 2 shows a circular opening excavated in isotropic, continuous rock mass. The initial stress is subject to a hydrostatic pressure p0. The inner surface of the opening is subject to a pressure pi. The inner pressure pi is equal to the initial stress state p0 before the excavation. After the tunnel is caved, the inner pressure pi decreases gradually and the surrounding rock deformation converges. When pi is lower than the critical value pc, i.e. the elastoplastic support pressure, the plastic zone develops and the circumferential stress shows a sudden drop. Then, the rock mass material properties, such as the strength parameters and deformation parameters, reach the residual value. The coupled elasto-plastic model is shown in Fig. 3 which is required to solve the stress and displacement of rock mass. It should be noted that the residual Young’s modulus of rock mass is lower than that of the intact rock mass for the elasto-brittle-plastic rock mass. In this paper, two of the most commonly-used yield criteria are considered: the linear M–C criterion
Elastic region
P0
Fig. 2. A circular opening in an infinite medium.
Ideal plastic Maximum deviatoric stress (MPa)
sample is deteriorated in the post-failure region, namely the elasto-plastic coupling effect, as shown in Fig. 1. Generally speaking, those researches can be concluded to two kinds. One supposed the Young’s modulus depends on the confining pressure or minor principal stress, i.e. the pressure-dependent Young’s modulus (PDM) model [17–20], and the other supposed the Young’s modulus is the function of radius, i.e. the radius dependent Young’s modulus (RDM) model [21,22]. To some extent, the RDM model can be regarded as an extension of the PDM model when the initial material parameters of surrounding rock at different radius is function of r and radial stress. As we know, the progressive failure of rock mass can be regarded as an accumulation process of irrecoverable deformation [23]. So, the RDM model would be more reasonable, and was adopted to obtain a more rigorous elasto-plastic coupling analytical solutions for stresses and displacements of circular openings in elasto-brittle-plastic rock mass. Additionally, almost all of the analytical solutions are based on the continuum mechanics, and the deterioration degree of material properties are taken objectively, especially of the Young’s modulus [24]. The obtained results would be contrary with the common sense and theoretical hypothesis, when the Young’s modulus of post-failure rock mass is less than a certain value. In another word, the Young’s modulus can’t deteriorate unboundedly with the elasto-plastic coupling model and there is a certain restrictive condition between rock mass strength parameters and Young’s modulus for elasto-plastic coupling rock mass. Finally, the proposed model was merged in FLAC3D code by compiling the elasto-plastic coupling model (EPCM) for a more widely use.
Brittle plastic
Unloading path
Axial strain (-) Fig. 3. Material behavior models.
r1 ¼ Nr3 þ Y
ð1Þ
and the nonlinear H–B criterion
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
r1 ¼ r3 þ mrc r3 þ sr2c
ð2Þ
where r1 is the major principal stress, r3 is the minor principal stress, N = (1 + sin u)/(1 sin u), Y = 2ccos u/(1 sin u), u and c are the friction angle and cohesion of rock mass respectively. rc is the unaxial compressive strength of intact rock mass, m and s are material property parameters. As for the plane strain problem of circular openings, tangential stress (rh) and radial stress (rr) are often the major and minor stresses, respectively. In this way, it is defined that r1 = rh and r3 = rr, then the yield criteria of Eqs. (1) and (2) for both intact and broken rock mass can be expressed as
rh ¼ Nrr þ Y for intact rock mass rh ¼ Nr rr þ Y r for residual rock mass
ð3Þ ð4Þ
62
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
or
rh rh
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ rr þ mrc rr þ sr2c for intact rock mass qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ rr þ mr rcr rr þ sr r2c for residual rock mass
ð5Þ ð6Þ
where the Nr, Yr, rcr, mr and sr are residual values of rock properties compared to the intact rock mass, and Nr = (1 + sin ur)/(1 sin ur), Yr = 2crcos ur/(1 sin ur), in which ur and cr is the residual friction angle and cohesion, respectively. According to the plastic theory, the geotechnical materials generally obey the Mohr–Coulomb plastic potential function, which can be expressed as
g ¼ rh arr
ð7Þ
where a = (1 + sin w)/(1 sin w), and w is the dilation angle of rock mass. The plastic strain can be obtained with the potential theory, as shown in the following equation:
epr ¼ k
@g ; @ rr
eph ¼ k
@g @ rh
ð8Þ
where k is the plastic strain multiplier. The strain of rock mass in plastic region is composed of elastic and plastic parts. Combined the integration results of Eq. (8) with the elastic strain, the total strain of rock mass induced by excavation can be expressed as 1 er ¼ 2fG ½ð1 lÞðrr p0 Þ lðrh p0 Þ ka 1 eh ¼ 2fG ½ð1 lÞðrh p0 Þ lðrr p0 Þ þ k
) ð9Þ
where er and eh are the radial and tangential strain, respectively, f is the deterioration coefficient of Young’s modulus, G is the shear modulus of intact rock mass. The geometric equation of this axisymmetric problem can be denoted as
er
@u ¼ ; @r
eh
u ¼ r
ð10Þ
The boundary conditions of this problem can be described as
r ¼ R0 ;
ðrr ¼ pi Þ
r ¼ Rp ;
ðrrp ¼ rre ; up ¼ ue Þ
r ! 1; ðrr ¼ p0 Þ
9 > = > ;
ð11Þ
where AMC = pi + BMC, BMC = crcot ur. Using the geometric equation of Eq. (10), the compatibility equation for deformation can be expressed as
@ eh er eh ¼0 @r r
ð14Þ
By substituting Eqs. (9) and (13) into Eq. (14), the differential equation for plastic multiplier k can be obtained. And its general solution can be expressed as
k¼
ð1þaÞ Nr 1 T r ð1 N2r Þð1 lÞ r þ AMC 2fG Rp 2fGðN r þ aÞ R0
ð15Þ
where T is the integration constant that should be determined. Then the strain and deformation of surrounding rock mass can be obtained by substituting Eq. (15) into Eq. (9), namely
aþ1 R þ T a rp ð1 2lÞp0r Nr 1 aþ1 R 1 ¼ 2fG gAMC Rr0 T a rp ð1 2lÞp0r
1 er ¼ 2fG Nr gAMC
eh
Nr 1 r R0
ð16Þ
u ¼ r eh where g = (1 l)(1 + Nra)/(Nr + a) l, p0r ¼ p0 þ BMC . Using the continuity of the radial deformation at the elastoplastic interface, the integration constant T can be determined, as shown in the following equation:
T ¼ ðg þ fÞAMC ðRp =R0 ÞNr 1 ð1 þ f 2lÞp0r
ð17Þ
As we know, the brittle failure occurs only when the elastic stresses reach its peak strength. After the sudden drop of stress, the new stress comes to the state of residual strength, so the stresses at the elasto-plastic interface follow the following equations
rhe rre ¼ 2ðp0 pcMC Þ rhp rrp ¼ ðN 1ÞpcMC þ Y
ð18Þ
where rhe and rre are the tangential and radial stress on the elastic side of elasto-plastic boundary, rhp and rrp are the corresponding stress on the plastic side, pcMC is support pressure of M–C rock mass, namely pcMC ¼ ð2p0 þ YÞ=ð1 þ NÞ. The stress and deformation in the elastic region can be easily obtained. Then the radius of plastic region can be obtained by substituting the corresponding stress expression into Eq. (18).
where rrp and rre are the radial stress of rock mass in the plastic region and in elastic region, respectively, and the up and ue are the rock mass deformation in the plastic and elastic region, respectively.
Rp ¼ R0
3. Analytical solutions
3.2. Elasto-plastic coupling solutions of H–B rock mass
3.1. Elasto-plastic coupling solutions of M–C rock mass
In the same way, the stress, deformation and radius in plastic region can be obtained for H–B rock mass.
The differential equation of equilibrium for axisymmetric problems in residual region can be expressed as (not considering the body force of rock mass)
drr rr rh þ ¼0 dr r
ð12Þ
By substituting Eq. (4) into Eq. (12) and using the boundary condition at r = R0, namely, the first equation of Eq. (11), the stresses in plastic region for the M–C rock mass can be denoted as
rr ¼ AMC
Nr 1 r R0
rh ¼ Nr AMC
BMC
Nr 1 r R0
9 > =
> BMC ;
ð13Þ
c 1 pMC þ BMC Nr 1 AMC
ð19Þ
9 =
rr ¼ AHB ln2 Rr0 þ BHB ln Rr0 þ pi
ð20Þ
rh ¼ AHB ln2 Rr0 þ ð2AHB þ BHB Þ ln Rr0 þ BHB þ pi ; where AHB = mrrcr/4, BHB ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi mr rcr pi þ sr r2cr .
h
i9
2 1 er ¼ 2fG d1 ln Rr0 þ ð2d1 þ d2 Þ ln Rr0 þ d2 þ d3 C ar r 1a > > =
h
2 1 eh ¼ 2fG d1 ln Rr0 þ d2 ln Rr0 þ d3 þ Cr 1a u ¼ r eh
i
> > ;
ð21Þ
where d1 = (1 2l)AHB, d2 = 2dAHB + (1 2l)BHB, d3 = (1 2l) (pi p0) + dBHB 2dAHB/(1 + a), d = (1 l)(a 1)/(1 + a) and C is integration constant.
63
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
In the same way, using the continuity of the radial deformation and the stress state at the elasto-plastic interface, the integration constant C and the radius of residual region can be denoted as
Rp 2dAHB C ð1 2l þ fÞ p0 pcHB dBHB 2dAHB ln þ R0 1 þ a qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2
ffi3 BHB þ B2HB 4AHB pi pcHB 5 Rp ¼ R0 exp 4 2AHB a ¼ R1þ p
ð22Þ
ð23Þ
where
pcHB
Based on the generalized Hooke’s law, the stress at the peak point A (in Fig. 4) satisfies the following equation
EeA1 ¼ rA1 2lrA3
ð25Þ
where eA1 is the elastic strain at point A. After the brittle failure occurs, stress state located at point B. As we know, the re-loading path obeys the elastic path before it reaches the unloading point. So, the stress and strain at point B should satisfy
Er eB1 ¼ rB1 2lrB3
0 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 m2 mp 1 m 0 ¼ p0 rc @ þ þs A 2 4 8 rc
ð24Þ
which is the critical support pressure of H–B rock mass.
where Er is the residual Young’s modulus and e is elastic strain at point B. Before and after brittle failure, axial strain e1 and confining pressure r3 for rock sample is uniform. So the minimum value of deterioration coefficient of Young’s modulus (f = Er/E) can be obtained by submitting Eq. (25) into Eq. (26), namely,
4. The restrictive condition of residual material properties Here, the loading-unloading process in the post-failure region of conventionally triaxial compression test is taken to study the restrictive condition for ideal elasto-brittle-plastic rock mass properties, as shown in Fig. 4, and the axial pressure and lateral pressure are r1 and r3 respectively. When the stresses satisfy Eq. (2) or Eq. (4), brittle failure occurs and the strength properties of rock mass and deformation parameter (E) both deteriorate. A series of unloading paths at the point B on the stress-strain curve with different Young’s modulus deterioration coefficient are shown in Fig. 4. The path of L1 shows that the post-failure rock mass has the same Young’s modulus as intact rock mass. However, when the Young’s modulus decreases to a certain value, namely, the unloading path of L2, the accumulated plastic strain turns to be zero after the load is completely unloaded. Moreover, when the Young’s modulus decreases further, such as path of L3, after the load is completely unloaded, the point would locate at point G with a negative plastic strain both in the axial and circumferential direction. Obviously, the unloading path of L3 violate the Drucker assumption, which denotes that the length of rock sample increases and the radius of rock sample decrease during the loading and unloading process. For stable materials, irrecoverable deformation should be consistent with the direction of loads, but the condition of L3 shows that the rock sample becomes thinner and longer than that at the initial stress state after the loading-unloading process, which is against with the common sense. So, the path of L2 should be the limit condition of Young’s modulus deterioration degree.
ð26Þ B 1
fcr ¼
rA1 2lrA3 rA1 2lr3 ¼ rB1 2lrB3 rB1 2lr3
ð27Þ
The rock mass should satisfy the peak and residual yield criteria at point A and B respectively. As for M–C rock mass, the limitation for fcr can be simplified as
fcr MC ¼
ðNr 2lÞr3 þ Y r ðN 2lÞr3 þ Y
ð28Þ
In the same way, the limit of fcr for H–B rock mass can be expressed as
fcr HB ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi mr rcr r3 þ sr r2 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffifficr ð1 2lÞr3 þ mrc r3 þ sr2c
ð1 2lÞr3 þ
ð29Þ
At the boundary of elasto-plastic interface of circular openings in elastic-brittle-plastic rock mass, the radial stress is continuous and the tangential stress drops. So, the rupture process is the same as the traditional compression test of brittle rock mass. By substituting the radial stress at the elasto-plastic interface into Eqs. (28) and (29), respectively, the minimum value of fcr for M–C and H–B rock mass can be expressed as
ðNr 2lÞpcMC þ Y r ðN 2lÞpcMC þ Y pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1 2lÞpcHB þ mr rcr r3 þ sr r2cr pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ ð1 2lÞpcHB þ mrc r3 þ sr2c
fcr MC ¼
ð30Þ
fcr HB
ð31Þ
Maximum deviatoric stress (MPa)
5. Verification examples
B
C
L3
L2 G
5.1. Common compression test with numerical method
A
O
L1 F
E
Axial strain (-) Fig. 4. The unloading path in the post-failure region.
The common compression test to obtain rock mass material properties is generally used in geotechnical engineering. So in this part, the elasto-plastic coupling constitutive model is developed based on FLAC3D code. Meanwhile, the linear Mohr–Coulomb criterion with ultimate tensile strength is considered, and the domains used in the definition of the yield criterion in (r1, r3)-plane is shown in Fig. 5. When the stresses of rock mass are located at the outside of yield surface, the plastic revision would carry out, and the implementation of the elasto-plastic coupling calculation is shown as the following. When a new loading step is carried out, it is supposed that all of the components of the stress component induced by the totalstrain increment is Si (the subscript i = 1, 2, 3 represents the three principal stress)
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
h=0
-
+
Domain 1
f s =0
-
Maximum principal stress (MPa)
64
9
rN1 ¼ rI1 kt a2 > = rN2 ¼ rI2 kt a2 > rN3 ¼ rI3 kt a1 ;
ð37Þ
where kt is the plastic multiplier for the ultimate tensile strength criterion, namely, kt ¼ rI3 rt =a1 . In addition, the intersection of Mohr–Coulomb criterion and maxmum tensile criterion, namely the combined yield criterion can be described as
Domain 2 f t =0
h ¼ r3 rt þ aP ðr1 rP Þ
+
p
ð38Þ
p
where a and r are two constants defined as
+ Minimum principal stress (MPa)
9 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi = 1 þ N2uðcp Þ þ Nuðcp Þ ffiffiffiffiffiffiffiffiffiffiffiffi p ¼ rt Nuðcp Þ 2cðcp Þ Nuðcp Þ ;
aP ¼
rP
ð39Þ
Then the plastic flow direction should be the diagonal between the fs = 0 and ft = 0 in the r1 r3 plane.
@h ¼ aP ; @ r1
Fig. 5. Domains used in the definition of the yield criterion.
9 S1 ¼ a1 Dee1 þ a2 Dee2 þ Dee3 >
= S2 ¼ a1 Dee2 þ a2 Dee2 þ Dee3
> ; S3 ¼ a1 Dee3 þ a2 Dee1 þ Dee2
ð32Þ
where the variations of a1 and a2 are a function of plastic shear strain cp, namely, a1 ¼ Kðcp Þ þ 4Gðcp Þ=3; a2 ¼ Kðcp Þ 2Gðcp Þ=3. Then the tentative elastic stress components, can be calculated as
rIi ¼ rOi þ Si
ð33Þ
where rIi and rOi are respectively the components of the tentative elastic stress and old stress. If the tentative elastic stress does not satisfy the yield criterion, the tentative elastic stress rIi would be at a real stress state, otherwise, the plastic revision implementing. If the tentative elastic stress satisfies the shear yield criterion (stress state in domain 1), the stress flow direction are calculated based on the plastic potential gs.
@g s ¼ 1; @ r1
@g s ¼ 0; @ r2
@g s ¼ Nw @ r3
ð34Þ
where gs is the potential flow rule, namely, gs = r1 Nwr3, and Nw is the dilatancy parameter, which is as a function of w(cp), namely, Nw = [1 + sin w(cp)]/[1 sin w(cp)]. Then the modified stress state can be obtained by the subtraction of the plastic stress increment, namely
9
rN1 ¼ rI1 ks ða1 a2 Nw Þ > = rN2 ¼ rI2 ks a2 ð1 Nw Þ > rN3 ¼ rI3 ks ða1 Nw þ a2 Þ ;
ð35Þ
where ks is the plastic multiplier for Mohr–Coulomb yield criterion, f s ðrI ;rI Þ
and it can be determined by ks ¼ ða1 a2 Nw Þð1 a13Nw þa2 ÞNu . In the other hand, when the minor principal stress reaches the tensile strength (in domain 2), that is to say, ft = r3 rt, the tensile flow rule should be considered. And the flow direction can be calculated with the flow equations corresponding to the tensile plastic potential function gt = r3.
@g s ¼ 0; @ r1
@g s ¼ 0; @ r2
@g s ¼1 @ r3
ð36Þ
Then the modified stress can be calculated by the Eq. (37), namely,
@h ¼ 0; @ r2
@h ¼1 @ r3
ð40Þ
Then the modified stress can be calculated by the Eq. (41), namely,
9
rN1 ¼ rI1 kh ða1 aP þ a2 Þ > = rN2 ¼ rI2 kh a2 ðaP þ 1Þ > ; rN3 ¼ rI3 kh ða1 þ a2 aP Þ
ð41Þ
where kh is the plastic multiplier for the combined strength hðrI ;rI Þ criterion, which can be determined by kh ¼ ða þa aP Þþ1aP ð3a aP þa Þ. 1
2
1
2
When a new loading step is carried out, the calculation steps can be followed with the following steps: Step 1: new loading step; Step 2: calculate the tentative elastic stress, and update the stress state, then judge whether the new stress satisfies the peak yield criterion. If not, go to step 5, or go to step 3; Step 3: modify the stress state, plastic strain and softening index (namely, the plastic shear strain in this paper) using the return mapping method in the principal stress space; Step 4: update the material parameters and deformation parameters. Step 5: next loading step. The elasto-plastic coupling model (EPCM) is developed based on the initial strain-softening model, and the material parameters of rock mass are input according to the table method. So, the EPCM can implemented by setting the equivalent plastic strain corresponding to the residual material parameters sufficiently small. In the past, most of the scholars determines the residual Young’s modulus subjectively based on strain-softening model with FLAC3D code, without any theoretical foundation. Even though we can obtain their value by traditional laboratory tests, the discreteness of geotechnical materials would also induce error results, especially after brittle failure happens. So, the correctness should be validated by the threshold condition, namely, f P fcr. This paper presents the restrictive condition for the elasto-brittle-plastic rock mass, which also exists for the strain-softening rock mass. The similar results
Table 1 Rock sample material properties. E (GPa)
v
c (MPa)
u (°)
cr (MPa)
ur (°)
W (°)
300
0.25
0.3
35
0.1
15
10
65
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
5.2. Circular openings problems In order to study the effect of different Young’s modulus deterioration degree on the surrounding rock mass deformation, the typical examples for kinds of M–C and H–B rock mass are studied [9,14,15,25]. The rock mass properties are list in Table 2.
Maximum principal stress (MPa)
Fig. 7. Numerical calculated model.
=0.15
=0.5 =1.0
Axial strain (-) Fig. 8. Loading and unloading curves with r3 = 0.
Maximum principal stress (MPa)
can easily be drawn for the strain-soften rock mass considering the additional axial strain in the softening part, which isn’t further studied in this paper. Then, we study the elasto-plastic coupling effect through triaxial compression numerical method based on the developed EPCM with FLAC3Dcode. In this example, the M–C criterion is chosen. The rock sample properties are shown in Table 1. The minimum fcr = 0.226 and fcr = 0.336 can be obtained with confining pressure r3 = 0 MPa and r3 = 1 MPa respectively using Eq. (27). The restrictive condition of Eq. (28) shows that the deterioration coefficient for residual rock mass should lie in the region between (Nr 2l)/ (N 2l) and Yr/Y. And in this example, the minimum fcr are 0.226 and 0.376 for the confining pressure of r3 = 0 and r3 ? 1 respectively. Fig. 6 shows the minimum deterioration coefficient for Young’s modulus with different confining pressure. In this example, the residual Young’s modulus increases with confining pressure. The loading-unloading process is studied with different confining pressure using the elasto-plastic coupling model. Fig. 7 shows the numerical model, which is discretized into 15,320 hexahedral elements with 17,056 nodes. The standard rock sample presented by ISRM, 50 mm diameter, 100 mm height is adopted. Meanwhile, the loading rate of 1.0 108 is taken to simulate a static loading process. In order to verify the correction of the proposed restrictive condition, another two sets of deterioration coefficient less and bigger than calculated fcr are considered. Meanwhile, the ideal elastobrittle-plastic model is also taken to compare with the proposed model. Fig. 8 shows the loading-unloading path of unaxial compression sample. We can see that the unloading point with the value f of 0.15 (which is less than fcr) located at the negative part of axial strain. The results illustrate that the length of rock sample increases in the loading and unloading process. The same conclusion can be drawn from triaxial compression tests, as shown in Fig. 9. The confining pressure r3 of 1 MPa is loaded initially after the axial displacement is fixed. Using the generalized Hooke’s law, the axial pressure caused by confining pressure should be r1 = 0.5 MPa. Similarly, the case with f of 0.30 also induces negative plastic strain in axial direction after the loading-unloading process. The plastic strain is a non-decreasing variable for the stable material, so the former two cases that induce negative plastic strain should be wrong. In the other hand, it indicates the correctness of the proposed restrictive condition of elasto-plastic coupling rock mass.
=0.5 =0.15 =1.0
Axial strain (-)
Critical Young’s modulus deterioration coefficient (-)
Fig. 9. Loading and unloading curves with r3 = 1.0 MPa.
Minimum principal stress (MPa) Fig. 6. Restrictive condition between Young’s modulus and confining pressure.
Eqs. (13) and (20) show that the deterioration of Young’s modulus has no influence on the stress of rock mass, so, only the deformation of rock mass is studied here. The limitations of f for different cases are calculated using Eqs. (28) and (29), as shown in Table 3. In this problem, the residual region is determined by the stress state on the elasto-plastic boundary. And the stresses are only related with the strength properties of rock mass both in the elastic and plastic region. This research result is the same as the former study [9,15], as shown in Figs. 10 and 11. Another two sets of Young’s modulus deterioration coefficient are considered, which satisfies the restrictive condition of Eqs. (28) and (29). Figs. 12 and 13 show the influence of Young’s
66
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
Table 2 Rock mass properties [25]. Types of rock
Hard rock
Soften rock
Radius of opening, R0 (m) Initial stress, p0 (MPa) Supporting pressure, r0 (MPa) Young’s modulus, E(GPa) Poisson ratio, v dilation angle, w (°)
1 1 0 50 0.2 0, 30
1 1 0 5 0.2 0, 30
M–C rock mass c (MPa) u (°) cr (MPa) ur (°)
0.173 55 0.061 52
0.276 35 0.055 30
H–B rock mass m s mr sr rc, rcr (MPa)
0.5 1.0E4 0.3 1.0E5 75
0.2 1.0E4 0.05 1.0E5 50
Fig. 12. Dimensionless displacement in M–C hard rock mass.
Table 3 Calculated limitation of f. Types of rock
Limitation of f
Hard rock
Soften rock
H–B
M–C
H–B
M–C
0.726
0.536
0.524
0.414
Fig. 13. Dimensionless displacement in M–C soften rock mass.
Radial and tangential Stresses (MPa)
Fig. 10. Distribution of stresses in M–C hard rock mass.
This paper K.H. Park [13]
Fig. 14. Dimensionless displacement in H–B hard rock mass.
Dimensionless radius r/R 0 (m/m) Fig. 11. Distribution of stresses in H–B soft rock mass.
modulus deterioration with M–C rock mass. When the dilation angle of hard rock mass is 30°, the inner surface deformation displacement increases by 6.09% with the decrease of f from 1.0 to 0.6, and the corresponding value is 6.24% when the rock mass
Fig. 15. Dimensionless displacement in H–B soften rock mass.
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
dilation angle is 0°. For the M–C soften rock mass, when f decreases from 1.0 to 0.5, the air surface displacement increases by 24.90% and 26.37% with dilation angle of 30° and 0°, respectively. Figs. 14 and 15 show the influence of Young’s modulus deterioration with H–B rock mass. Similarly, the deformation on inner surface increases with the deterioration of Young’s modulus. For the H–B hard rock, the displacement on air surface increase by 1.09% and 1.53% with the decrease of f from 1.0 to 0.8, corresponding to the dilation angle of 30° and 0° respectively. And for H–B soft rock, when f decrease from 1.0 to 0.6, the growth of air surface displacement is 13.24% and 15.59% with the dilation angle of 30° and 0° respectively. Additionally, both of M–C and H–B rock mass show that the rock mass displacement growth slightly with the increasing dilation angle and the Young’s modulus sharply influences the deformation of rock mass than that of the dilation. Another thing should be noted is that the restrictive condition for Young’s modulus deterioration degree also depends on the difference between peak strength and residual strength of rock mass.
[2]
[3] [4] [5]
[6]
[7]
[8] [9] [10]
[11]
6. Conclusion [12]
The analytical solutions considering elasto-plastic coupling effect was deduced for circular openings in elasto-brittle-plastic rock mass governed by Mohr–Coulomb and Hoek–Brown criteria subject to a hydrostatic in suit stress. Meanwhile, the restrictive condition between strength parameters and Young’s modulus was presented, which indicates there is a minimum Young’s modulus deterioration degree. The results were further validated by the numerical test in traditional triaxial compression from loading to unloading using FLAC3D code with the developed EPCM. Finally, the displacements around circular openings were studied with different deterioration degree of Young’s modulus using two kinds of rock mass. The results show that the deterioration of Young’s modulus has little influence on stresses, but increases rock mass deformation sharply in the plastic region. The solutions in this paper are an expansion for the traditional elasto-brittle-plastic and ideal elasto-plastic models. As we know, the calculated deformation of rock mass is usually less than the in situ value, however, the deformation of rock mass increases sharply when the elasto-plastic coupling character is considered, but none do the stresses in the plastic region. Acknowledgements This project is supported by the Fundamental Research Funds for the Central Universities (2012QNB23) and also supported by the National Natural Science Foundation of China (Project Nos. 51174196, 51179185 and 51104151.
[13]
[14]
[15]
[16]
[17] [18]
[19] [20]
[21]
[22]
[23]
[24]
[25]
References [1] E. Hoek, C.T. Carranza-Torres, B. Corkum, Hoek–Brown failure criterion-2002 edition, in: Proceedings of the 5th North American Rock Mechanics
67
Symposium and 17th Tunnelling Association of Canada Conference, Toronto 2002, pp. 267–273. A.L. Florence, L.E. Schwer, Axisymmetric compression of a Mohr–Coulomb medium around a circular hole, Int. J. Numer. Anal. Meth. Geomech. 2 (1978) 367–379. E.T. Brown, J.W. Bray, B. Ladanyi, E. Hoek, Ground response curves for rock tunnels, J. Geotech. Eng., ASCE 109 (1) (1983) 15–39. Y. Wang, Ground response of circular tunnel in poorly consolidated rock, J. Geotech. Eng. 122 (1996) 703–708. C. Carranza-Torres, Self Similar Analysis of the Elastoplastic Response of Underground Openings in Rock and Effects of Practical Variables, Ph.D. Thesis, University of Minnesota, 1998. C. Carranza-Torres, C. Fairhurst, The elasto-plastic response of underground excavations in rock masses that satisfy the Hoek–Brown failure criterion, Int. J. Rock Mech. Min. Sci. 36 (1999) 777–809. C. Carranza-Torres, Elasto-plastic solution of tunnel problems using the generalized form of the Hoek–Brown failure criterion, Int. J. Rock Mech. Min. Sci. 41 (2004) 480–481. K.H. Park, Y.J. Kim, Analytical solution for a circular opening in an elasticbrittle-plastic rock, Int. J. Rock Mech. Min. Sci. 43 (2006) 616–622. S.K. Sharan, Elastic-brittle-plastic analysis of circular openings in Hoek–Brown media, Int. J. Rock Mech. Min. Sci. 40 (2003) 817–824. S.K. Sharan, Exact and approximate solutions of displacements around circular openings in elastic-brittle-plastic Hoek–Brown rock, Int. J. Rock Mech. Min. Sci. 42 (4) (2005) 524–532. S.K. Sharan, Analytical solutions for stresses and displacements around a circular opening in a generalized Hoek–Brown rock, Int. J. Rock Mech. Min. Sci. 40 (2008) 78–85. E. Alonso, L.R. Alejano, F. Varas, G. Fdez-Manin, C. Carranza-Torres, Ground response curves for rock masses exhibiting strain-softening behavior, Int. J. Numer. Anal. Meth. Geomech. 27 (2003) 1153–1185. K.H. Park, B. Tontavanich, J.G. Lee, A simple procedure for ground response curve of circular tunnel in elastic-strain softening rock masses, Tunnel. Undergr. Space Technol. 23 (2008) 151–159. Y.K. Lee, S. Pietruszczak, A new numerical procedure for elasto-plastic analysis of a circular opening excavated in a strain-softening rock mass, Tunnel. Undergr. Space Technol. 23 (2008) 588–599. S.L. Wang, X.T. Yin, H. Tang, X.R. Ge, A new approach for analyzing circular tunnel in strain-softening rock masses, Int. J. Rock Mech. Min. Sci. 1 (2010) 170–178. Q. Zhang, B.S. Jiang, S.L. Wang, X.R. Ge, H.Q. Zhang, Elasto-plastic analysis of a circular opening in strain-softening rock mass, Int. J. Rock Mech. Min. Sci. 50 (2012) 38–46. F.H. Kulhawy, Stress deformation properties of rock and rock discontinuities, Eng. Geol. 9 (4) (1975) 327–350. E.T. Brown, J.W. Bray, F.J. Santarelli, Influence of stress-dependent elastic moduli on stresses and strains around axisymmetric boreholes, Rock Mech. Rock Eng. 22 (3) (1989) 189–203. M.R. Asef, D.J. Reddish, The impact of confining stress on the rock mass deformation modulus, Geotechnique 52 (4) (2002) 235–241. M. Verman, B. Singh, M.N. Viladkar, J.L. Jethwa, Effect of tunnel depth on modulus of deformation of rock mass, Rock Mech. Rock Eng. 30 (3) (1997) 121–127. R.T. Ewy, N.G.W. Cook, Deformation and fracture around cylindrical openings in rock – I. Observations and analysis of deformations, Int. J. Rock Mech. Min. Sci. Abstr. 27 (5) (1990) 387–407. P.A. Nawrocki, M.B. Dusseault, Modelling of damaged zones around openings using radius-dependent Young’s modulus, Rock Mech. Rock Eng. 28 (4) (1995) 227–239. P.K. Kaiser, A new concept to evaluate tunnel performance-influence of excavation procedure, in: Preceding of the 22nd US Symposium on Rock Mechanics, Cambridge, MA, 1981, pp. 264–271. C.G. Zhang, J.H. Zhao, Q.H. Zhang, X.D. Hu, A new closed-form solution for circular openings modeled by the unified strength theory and radiusdependent Young’s modulus, Comput. Geotech. 42 (2012) 118–128. T. Ogawa, K.Y. Lo, Effects of dilatancy and yield criteria on displacements around tunnels, Can. Geotech. J. 24 (1987) 100–113.
Contents lists available at SciVerse ScienceDirect
Theoretical and Applied Fracture Mechanics journal homepage: www.elsevier.com/locate/tafmec
Elasto-plastic coupling analysis of circular openings in elasto-brittle-plastic rock mass Q. Zhang ⇑, B.S. Jiang, X.S. Wu, H.Q. Zhang, L.J. Han State Key Laboratory for Geomechanics and Deep Underground Engineering, School of Mechanics and Civil Engineering, China University of Mining and Technology, Xuzhou 221008, PR China
a r t i c l e
i n f o
Article history: Available online 28 June 2012 Keywords: Circular opening Cyclic loading and unloading Elasto-plastic coupling Analytical solution FLAC3D
a b s t r a c t This paper deals with elasto-plastic coupling solutions for the prediction of displacements around circular openings in elasto-brittle-plastic rock mass subject to hydrostatic stress. Both linear Mohr–Coulomb criterion and nonlinear Hoek–Brown criterion on the basis of a non-associated flow rule are considered. The restrictive condition between strength parameters and Young’s modulus in post-failure region for rock materials is theoretically established. Meanwhile, the proposed restrictive condition is validated by FLAC3D code that the advanced elasto-plastic coupling model (EPCM) based on Mohr–Coulomb criterion is merged in. Finally, the elasto-plastic coupling dimensionless displacements of circular opening in elasto-brittle-plastic rock mass are analyzed with different deterioration degree of Young’s modulus using two kinds of rocks. The results show that the deformation of rock mass increases obviously with the decrease of Young’s modulus in the plastic region, but Young’s modulus has little influence on the plastic region radius and stresses distribution form. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction Although plane strain problems of circular openings in rock mass are such simple cases, they can provide theoretical foundation for geotechnical and mining engineering problems, such as the optimization design of tunnels, boreholes and mine shafts. Lots of researchers show special interest on this axisymmetric problem and have done a great deal of contributions on the closed-form solutions for this problem [1–10]. As for most of the analyses reported in the past, the solutions are obtained considering different models of material behavior, such as ideal elasto-plastic, elasticbrittle-plastic and strain-softening models, with different yield criteria, like the linear Mohr–Coulomb (M–C) and nonlinear Hoek–Brown (H–B) criteria [1]. The stress state and deformation of M–C rock mass with a circular hole were earlier studied [2]. Until to 1983, the closed-form solutions for the stress and radial displacement of circular openings was obtained on the assumption of constant elastic strain in the plastic region, which is equal to that at the elasto-plastic interface [3]. Later, a numerical method for elasto-plastic rock mass was proposed avoiding the defects of assumption of constant elastic strain for plastic rock mass [4]. In order to avoid the difficulty, the so-called self-similarity theory combined with analytical method was adopted as an effective way to predicted rock mass displacements and stresses in ⇑ Corresponding author. Tel.: +86 13641535671. E-mail address: [email protected] (Q. Zhang). 0167-8442/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.tafmec.2012.06.008
elasto-plastic H–B rock mass [5]. Obviously, the above achievements are all based on the ideal elasto-plastic model, however, most of the brittle rock mass shows sudden loss of strength when the stress arrived at its peak strength. Considering the brittle failure characteristics, lots of scholars have also done further investigation. By assuming four types of strain distribution forms in the plastic region, the closed-form solutions of circular openings were analyzed with both the M–C and H–B yield criteria. As we know, the elastic strain should obey the Hooke’s law, even in the plastic region of surrounding rock. So only the third kind of solution is correct [8]. In order to obtain the explicit expressions for displacements, it was assumed that the elastic strain in the plastic region obeys the same law in thick wall cylinder condition. So, the solutions for generalized H–B rock mass are also not exactly coherent with the real displacements [9–11]. Meanwhile, the strain-softening calculation procedure using numerical method [12–14] and semi-analytical and semi-numerical method [15,16] were proposed, and the examples are well consistent with each other. The material properties depend on many factors, such as the stress state and extent of damage. The deterioration of strength parameters are commonly considered in the former studies, but the Young’s modulus was rarely considered. After the excavation of tunnels, the stress state of surrounding rock turns from three dimensional stress state to quasi-two dimensional state. The deterioration of confining condition and extent of broken would lead to the decrease of Young’s modulus. Additionally, numbers of traditional compression tests show that the Young’s modulus of rock
61
Maximum deviatoric stress (MPa)
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
P0
P0
Pi
R0
Axial strain (-)
Rp
Fig. 1. Cyclic loading and unloading test under lower confining pressure.
P0 Plastic region
2. Definition of the problem Fig. 2 shows a circular opening excavated in isotropic, continuous rock mass. The initial stress is subject to a hydrostatic pressure p0. The inner surface of the opening is subject to a pressure pi. The inner pressure pi is equal to the initial stress state p0 before the excavation. After the tunnel is caved, the inner pressure pi decreases gradually and the surrounding rock deformation converges. When pi is lower than the critical value pc, i.e. the elastoplastic support pressure, the plastic zone develops and the circumferential stress shows a sudden drop. Then, the rock mass material properties, such as the strength parameters and deformation parameters, reach the residual value. The coupled elasto-plastic model is shown in Fig. 3 which is required to solve the stress and displacement of rock mass. It should be noted that the residual Young’s modulus of rock mass is lower than that of the intact rock mass for the elasto-brittle-plastic rock mass. In this paper, two of the most commonly-used yield criteria are considered: the linear M–C criterion
Elastic region
P0
Fig. 2. A circular opening in an infinite medium.
Ideal plastic Maximum deviatoric stress (MPa)
sample is deteriorated in the post-failure region, namely the elasto-plastic coupling effect, as shown in Fig. 1. Generally speaking, those researches can be concluded to two kinds. One supposed the Young’s modulus depends on the confining pressure or minor principal stress, i.e. the pressure-dependent Young’s modulus (PDM) model [17–20], and the other supposed the Young’s modulus is the function of radius, i.e. the radius dependent Young’s modulus (RDM) model [21,22]. To some extent, the RDM model can be regarded as an extension of the PDM model when the initial material parameters of surrounding rock at different radius is function of r and radial stress. As we know, the progressive failure of rock mass can be regarded as an accumulation process of irrecoverable deformation [23]. So, the RDM model would be more reasonable, and was adopted to obtain a more rigorous elasto-plastic coupling analytical solutions for stresses and displacements of circular openings in elasto-brittle-plastic rock mass. Additionally, almost all of the analytical solutions are based on the continuum mechanics, and the deterioration degree of material properties are taken objectively, especially of the Young’s modulus [24]. The obtained results would be contrary with the common sense and theoretical hypothesis, when the Young’s modulus of post-failure rock mass is less than a certain value. In another word, the Young’s modulus can’t deteriorate unboundedly with the elasto-plastic coupling model and there is a certain restrictive condition between rock mass strength parameters and Young’s modulus for elasto-plastic coupling rock mass. Finally, the proposed model was merged in FLAC3D code by compiling the elasto-plastic coupling model (EPCM) for a more widely use.
Brittle plastic
Unloading path
Axial strain (-) Fig. 3. Material behavior models.
r1 ¼ Nr3 þ Y
ð1Þ
and the nonlinear H–B criterion
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
r1 ¼ r3 þ mrc r3 þ sr2c
ð2Þ
where r1 is the major principal stress, r3 is the minor principal stress, N = (1 + sin u)/(1 sin u), Y = 2ccos u/(1 sin u), u and c are the friction angle and cohesion of rock mass respectively. rc is the unaxial compressive strength of intact rock mass, m and s are material property parameters. As for the plane strain problem of circular openings, tangential stress (rh) and radial stress (rr) are often the major and minor stresses, respectively. In this way, it is defined that r1 = rh and r3 = rr, then the yield criteria of Eqs. (1) and (2) for both intact and broken rock mass can be expressed as
rh ¼ Nrr þ Y for intact rock mass rh ¼ Nr rr þ Y r for residual rock mass
ð3Þ ð4Þ
62
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
or
rh rh
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ rr þ mrc rr þ sr2c for intact rock mass qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ rr þ mr rcr rr þ sr r2c for residual rock mass
ð5Þ ð6Þ
where the Nr, Yr, rcr, mr and sr are residual values of rock properties compared to the intact rock mass, and Nr = (1 + sin ur)/(1 sin ur), Yr = 2crcos ur/(1 sin ur), in which ur and cr is the residual friction angle and cohesion, respectively. According to the plastic theory, the geotechnical materials generally obey the Mohr–Coulomb plastic potential function, which can be expressed as
g ¼ rh arr
ð7Þ
where a = (1 + sin w)/(1 sin w), and w is the dilation angle of rock mass. The plastic strain can be obtained with the potential theory, as shown in the following equation:
epr ¼ k
@g ; @ rr
eph ¼ k
@g @ rh
ð8Þ
where k is the plastic strain multiplier. The strain of rock mass in plastic region is composed of elastic and plastic parts. Combined the integration results of Eq. (8) with the elastic strain, the total strain of rock mass induced by excavation can be expressed as 1 er ¼ 2fG ½ð1 lÞðrr p0 Þ lðrh p0 Þ ka 1 eh ¼ 2fG ½ð1 lÞðrh p0 Þ lðrr p0 Þ þ k
) ð9Þ
where er and eh are the radial and tangential strain, respectively, f is the deterioration coefficient of Young’s modulus, G is the shear modulus of intact rock mass. The geometric equation of this axisymmetric problem can be denoted as
er
@u ¼ ; @r
eh
u ¼ r
ð10Þ
The boundary conditions of this problem can be described as
r ¼ R0 ;
ðrr ¼ pi Þ
r ¼ Rp ;
ðrrp ¼ rre ; up ¼ ue Þ
r ! 1; ðrr ¼ p0 Þ
9 > = > ;
ð11Þ
where AMC = pi + BMC, BMC = crcot ur. Using the geometric equation of Eq. (10), the compatibility equation for deformation can be expressed as
@ eh er eh ¼0 @r r
ð14Þ
By substituting Eqs. (9) and (13) into Eq. (14), the differential equation for plastic multiplier k can be obtained. And its general solution can be expressed as
k¼
ð1þaÞ Nr 1 T r ð1 N2r Þð1 lÞ r þ AMC 2fG Rp 2fGðN r þ aÞ R0
ð15Þ
where T is the integration constant that should be determined. Then the strain and deformation of surrounding rock mass can be obtained by substituting Eq. (15) into Eq. (9), namely
aþ1 R þ T a rp ð1 2lÞp0r Nr 1 aþ1 R 1 ¼ 2fG gAMC Rr0 T a rp ð1 2lÞp0r
1 er ¼ 2fG Nr gAMC
eh
Nr 1 r R0
ð16Þ
u ¼ r eh where g = (1 l)(1 + Nra)/(Nr + a) l, p0r ¼ p0 þ BMC . Using the continuity of the radial deformation at the elastoplastic interface, the integration constant T can be determined, as shown in the following equation:
T ¼ ðg þ fÞAMC ðRp =R0 ÞNr 1 ð1 þ f 2lÞp0r
ð17Þ
As we know, the brittle failure occurs only when the elastic stresses reach its peak strength. After the sudden drop of stress, the new stress comes to the state of residual strength, so the stresses at the elasto-plastic interface follow the following equations
rhe rre ¼ 2ðp0 pcMC Þ rhp rrp ¼ ðN 1ÞpcMC þ Y
ð18Þ
where rhe and rre are the tangential and radial stress on the elastic side of elasto-plastic boundary, rhp and rrp are the corresponding stress on the plastic side, pcMC is support pressure of M–C rock mass, namely pcMC ¼ ð2p0 þ YÞ=ð1 þ NÞ. The stress and deformation in the elastic region can be easily obtained. Then the radius of plastic region can be obtained by substituting the corresponding stress expression into Eq. (18).
where rrp and rre are the radial stress of rock mass in the plastic region and in elastic region, respectively, and the up and ue are the rock mass deformation in the plastic and elastic region, respectively.
Rp ¼ R0
3. Analytical solutions
3.2. Elasto-plastic coupling solutions of H–B rock mass
3.1. Elasto-plastic coupling solutions of M–C rock mass
In the same way, the stress, deformation and radius in plastic region can be obtained for H–B rock mass.
The differential equation of equilibrium for axisymmetric problems in residual region can be expressed as (not considering the body force of rock mass)
drr rr rh þ ¼0 dr r
ð12Þ
By substituting Eq. (4) into Eq. (12) and using the boundary condition at r = R0, namely, the first equation of Eq. (11), the stresses in plastic region for the M–C rock mass can be denoted as
rr ¼ AMC
Nr 1 r R0
rh ¼ Nr AMC
BMC
Nr 1 r R0
9 > =
> BMC ;
ð13Þ
c 1 pMC þ BMC Nr 1 AMC
ð19Þ
9 =
rr ¼ AHB ln2 Rr0 þ BHB ln Rr0 þ pi
ð20Þ
rh ¼ AHB ln2 Rr0 þ ð2AHB þ BHB Þ ln Rr0 þ BHB þ pi ; where AHB = mrrcr/4, BHB ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi mr rcr pi þ sr r2cr .
h
i9
2 1 er ¼ 2fG d1 ln Rr0 þ ð2d1 þ d2 Þ ln Rr0 þ d2 þ d3 C ar r 1a > > =
h
2 1 eh ¼ 2fG d1 ln Rr0 þ d2 ln Rr0 þ d3 þ Cr 1a u ¼ r eh
i
> > ;
ð21Þ
where d1 = (1 2l)AHB, d2 = 2dAHB + (1 2l)BHB, d3 = (1 2l) (pi p0) + dBHB 2dAHB/(1 + a), d = (1 l)(a 1)/(1 + a) and C is integration constant.
63
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
In the same way, using the continuity of the radial deformation and the stress state at the elasto-plastic interface, the integration constant C and the radius of residual region can be denoted as
Rp 2dAHB C ð1 2l þ fÞ p0 pcHB dBHB 2dAHB ln þ R0 1 þ a qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2
ffi3 BHB þ B2HB 4AHB pi pcHB 5 Rp ¼ R0 exp 4 2AHB a ¼ R1þ p
ð22Þ
ð23Þ
where
pcHB
Based on the generalized Hooke’s law, the stress at the peak point A (in Fig. 4) satisfies the following equation
EeA1 ¼ rA1 2lrA3
ð25Þ
where eA1 is the elastic strain at point A. After the brittle failure occurs, stress state located at point B. As we know, the re-loading path obeys the elastic path before it reaches the unloading point. So, the stress and strain at point B should satisfy
Er eB1 ¼ rB1 2lrB3
0 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 m2 mp 1 m 0 ¼ p0 rc @ þ þs A 2 4 8 rc
ð24Þ
which is the critical support pressure of H–B rock mass.
where Er is the residual Young’s modulus and e is elastic strain at point B. Before and after brittle failure, axial strain e1 and confining pressure r3 for rock sample is uniform. So the minimum value of deterioration coefficient of Young’s modulus (f = Er/E) can be obtained by submitting Eq. (25) into Eq. (26), namely,
4. The restrictive condition of residual material properties Here, the loading-unloading process in the post-failure region of conventionally triaxial compression test is taken to study the restrictive condition for ideal elasto-brittle-plastic rock mass properties, as shown in Fig. 4, and the axial pressure and lateral pressure are r1 and r3 respectively. When the stresses satisfy Eq. (2) or Eq. (4), brittle failure occurs and the strength properties of rock mass and deformation parameter (E) both deteriorate. A series of unloading paths at the point B on the stress-strain curve with different Young’s modulus deterioration coefficient are shown in Fig. 4. The path of L1 shows that the post-failure rock mass has the same Young’s modulus as intact rock mass. However, when the Young’s modulus decreases to a certain value, namely, the unloading path of L2, the accumulated plastic strain turns to be zero after the load is completely unloaded. Moreover, when the Young’s modulus decreases further, such as path of L3, after the load is completely unloaded, the point would locate at point G with a negative plastic strain both in the axial and circumferential direction. Obviously, the unloading path of L3 violate the Drucker assumption, which denotes that the length of rock sample increases and the radius of rock sample decrease during the loading and unloading process. For stable materials, irrecoverable deformation should be consistent with the direction of loads, but the condition of L3 shows that the rock sample becomes thinner and longer than that at the initial stress state after the loading-unloading process, which is against with the common sense. So, the path of L2 should be the limit condition of Young’s modulus deterioration degree.
ð26Þ B 1
fcr ¼
rA1 2lrA3 rA1 2lr3 ¼ rB1 2lrB3 rB1 2lr3
ð27Þ
The rock mass should satisfy the peak and residual yield criteria at point A and B respectively. As for M–C rock mass, the limitation for fcr can be simplified as
fcr MC ¼
ðNr 2lÞr3 þ Y r ðN 2lÞr3 þ Y
ð28Þ
In the same way, the limit of fcr for H–B rock mass can be expressed as
fcr HB ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi mr rcr r3 þ sr r2 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffifficr ð1 2lÞr3 þ mrc r3 þ sr2c
ð1 2lÞr3 þ
ð29Þ
At the boundary of elasto-plastic interface of circular openings in elastic-brittle-plastic rock mass, the radial stress is continuous and the tangential stress drops. So, the rupture process is the same as the traditional compression test of brittle rock mass. By substituting the radial stress at the elasto-plastic interface into Eqs. (28) and (29), respectively, the minimum value of fcr for M–C and H–B rock mass can be expressed as
ðNr 2lÞpcMC þ Y r ðN 2lÞpcMC þ Y pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1 2lÞpcHB þ mr rcr r3 þ sr r2cr pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ ð1 2lÞpcHB þ mrc r3 þ sr2c
fcr MC ¼
ð30Þ
fcr HB
ð31Þ
Maximum deviatoric stress (MPa)
5. Verification examples
B
C
L3
L2 G
5.1. Common compression test with numerical method
A
O
L1 F
E
Axial strain (-) Fig. 4. The unloading path in the post-failure region.
The common compression test to obtain rock mass material properties is generally used in geotechnical engineering. So in this part, the elasto-plastic coupling constitutive model is developed based on FLAC3D code. Meanwhile, the linear Mohr–Coulomb criterion with ultimate tensile strength is considered, and the domains used in the definition of the yield criterion in (r1, r3)-plane is shown in Fig. 5. When the stresses of rock mass are located at the outside of yield surface, the plastic revision would carry out, and the implementation of the elasto-plastic coupling calculation is shown as the following. When a new loading step is carried out, it is supposed that all of the components of the stress component induced by the totalstrain increment is Si (the subscript i = 1, 2, 3 represents the three principal stress)
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
h=0
-
+
Domain 1
f s =0
-
Maximum principal stress (MPa)
64
9
rN1 ¼ rI1 kt a2 > = rN2 ¼ rI2 kt a2 > rN3 ¼ rI3 kt a1 ;
ð37Þ
where kt is the plastic multiplier for the ultimate tensile strength criterion, namely, kt ¼ rI3 rt =a1 . In addition, the intersection of Mohr–Coulomb criterion and maxmum tensile criterion, namely the combined yield criterion can be described as
Domain 2 f t =0
h ¼ r3 rt þ aP ðr1 rP Þ
+
p
ð38Þ
p
where a and r are two constants defined as
+ Minimum principal stress (MPa)
9 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi = 1 þ N2uðcp Þ þ Nuðcp Þ ffiffiffiffiffiffiffiffiffiffiffiffi p ¼ rt Nuðcp Þ 2cðcp Þ Nuðcp Þ ;
aP ¼
rP
ð39Þ
Then the plastic flow direction should be the diagonal between the fs = 0 and ft = 0 in the r1 r3 plane.
@h ¼ aP ; @ r1
Fig. 5. Domains used in the definition of the yield criterion.
9 S1 ¼ a1 Dee1 þ a2 Dee2 þ Dee3 >
= S2 ¼ a1 Dee2 þ a2 Dee2 þ Dee3
> ; S3 ¼ a1 Dee3 þ a2 Dee1 þ Dee2
ð32Þ
where the variations of a1 and a2 are a function of plastic shear strain cp, namely, a1 ¼ Kðcp Þ þ 4Gðcp Þ=3; a2 ¼ Kðcp Þ 2Gðcp Þ=3. Then the tentative elastic stress components, can be calculated as
rIi ¼ rOi þ Si
ð33Þ
where rIi and rOi are respectively the components of the tentative elastic stress and old stress. If the tentative elastic stress does not satisfy the yield criterion, the tentative elastic stress rIi would be at a real stress state, otherwise, the plastic revision implementing. If the tentative elastic stress satisfies the shear yield criterion (stress state in domain 1), the stress flow direction are calculated based on the plastic potential gs.
@g s ¼ 1; @ r1
@g s ¼ 0; @ r2
@g s ¼ Nw @ r3
ð34Þ
where gs is the potential flow rule, namely, gs = r1 Nwr3, and Nw is the dilatancy parameter, which is as a function of w(cp), namely, Nw = [1 + sin w(cp)]/[1 sin w(cp)]. Then the modified stress state can be obtained by the subtraction of the plastic stress increment, namely
9
rN1 ¼ rI1 ks ða1 a2 Nw Þ > = rN2 ¼ rI2 ks a2 ð1 Nw Þ > rN3 ¼ rI3 ks ða1 Nw þ a2 Þ ;
ð35Þ
where ks is the plastic multiplier for Mohr–Coulomb yield criterion, f s ðrI ;rI Þ
and it can be determined by ks ¼ ða1 a2 Nw Þð1 a13Nw þa2 ÞNu . In the other hand, when the minor principal stress reaches the tensile strength (in domain 2), that is to say, ft = r3 rt, the tensile flow rule should be considered. And the flow direction can be calculated with the flow equations corresponding to the tensile plastic potential function gt = r3.
@g s ¼ 0; @ r1
@g s ¼ 0; @ r2
@g s ¼1 @ r3
ð36Þ
Then the modified stress can be calculated by the Eq. (37), namely,
@h ¼ 0; @ r2
@h ¼1 @ r3
ð40Þ
Then the modified stress can be calculated by the Eq. (41), namely,
9
rN1 ¼ rI1 kh ða1 aP þ a2 Þ > = rN2 ¼ rI2 kh a2 ðaP þ 1Þ > ; rN3 ¼ rI3 kh ða1 þ a2 aP Þ
ð41Þ
where kh is the plastic multiplier for the combined strength hðrI ;rI Þ criterion, which can be determined by kh ¼ ða þa aP Þþ1aP ð3a aP þa Þ. 1
2
1
2
When a new loading step is carried out, the calculation steps can be followed with the following steps: Step 1: new loading step; Step 2: calculate the tentative elastic stress, and update the stress state, then judge whether the new stress satisfies the peak yield criterion. If not, go to step 5, or go to step 3; Step 3: modify the stress state, plastic strain and softening index (namely, the plastic shear strain in this paper) using the return mapping method in the principal stress space; Step 4: update the material parameters and deformation parameters. Step 5: next loading step. The elasto-plastic coupling model (EPCM) is developed based on the initial strain-softening model, and the material parameters of rock mass are input according to the table method. So, the EPCM can implemented by setting the equivalent plastic strain corresponding to the residual material parameters sufficiently small. In the past, most of the scholars determines the residual Young’s modulus subjectively based on strain-softening model with FLAC3D code, without any theoretical foundation. Even though we can obtain their value by traditional laboratory tests, the discreteness of geotechnical materials would also induce error results, especially after brittle failure happens. So, the correctness should be validated by the threshold condition, namely, f P fcr. This paper presents the restrictive condition for the elasto-brittle-plastic rock mass, which also exists for the strain-softening rock mass. The similar results
Table 1 Rock sample material properties. E (GPa)
v
c (MPa)
u (°)
cr (MPa)
ur (°)
W (°)
300
0.25
0.3
35
0.1
15
10
65
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
5.2. Circular openings problems In order to study the effect of different Young’s modulus deterioration degree on the surrounding rock mass deformation, the typical examples for kinds of M–C and H–B rock mass are studied [9,14,15,25]. The rock mass properties are list in Table 2.
Maximum principal stress (MPa)
Fig. 7. Numerical calculated model.
=0.15
=0.5 =1.0
Axial strain (-) Fig. 8. Loading and unloading curves with r3 = 0.
Maximum principal stress (MPa)
can easily be drawn for the strain-soften rock mass considering the additional axial strain in the softening part, which isn’t further studied in this paper. Then, we study the elasto-plastic coupling effect through triaxial compression numerical method based on the developed EPCM with FLAC3Dcode. In this example, the M–C criterion is chosen. The rock sample properties are shown in Table 1. The minimum fcr = 0.226 and fcr = 0.336 can be obtained with confining pressure r3 = 0 MPa and r3 = 1 MPa respectively using Eq. (27). The restrictive condition of Eq. (28) shows that the deterioration coefficient for residual rock mass should lie in the region between (Nr 2l)/ (N 2l) and Yr/Y. And in this example, the minimum fcr are 0.226 and 0.376 for the confining pressure of r3 = 0 and r3 ? 1 respectively. Fig. 6 shows the minimum deterioration coefficient for Young’s modulus with different confining pressure. In this example, the residual Young’s modulus increases with confining pressure. The loading-unloading process is studied with different confining pressure using the elasto-plastic coupling model. Fig. 7 shows the numerical model, which is discretized into 15,320 hexahedral elements with 17,056 nodes. The standard rock sample presented by ISRM, 50 mm diameter, 100 mm height is adopted. Meanwhile, the loading rate of 1.0 108 is taken to simulate a static loading process. In order to verify the correction of the proposed restrictive condition, another two sets of deterioration coefficient less and bigger than calculated fcr are considered. Meanwhile, the ideal elastobrittle-plastic model is also taken to compare with the proposed model. Fig. 8 shows the loading-unloading path of unaxial compression sample. We can see that the unloading point with the value f of 0.15 (which is less than fcr) located at the negative part of axial strain. The results illustrate that the length of rock sample increases in the loading and unloading process. The same conclusion can be drawn from triaxial compression tests, as shown in Fig. 9. The confining pressure r3 of 1 MPa is loaded initially after the axial displacement is fixed. Using the generalized Hooke’s law, the axial pressure caused by confining pressure should be r1 = 0.5 MPa. Similarly, the case with f of 0.30 also induces negative plastic strain in axial direction after the loading-unloading process. The plastic strain is a non-decreasing variable for the stable material, so the former two cases that induce negative plastic strain should be wrong. In the other hand, it indicates the correctness of the proposed restrictive condition of elasto-plastic coupling rock mass.
=0.5 =0.15 =1.0
Axial strain (-)
Critical Young’s modulus deterioration coefficient (-)
Fig. 9. Loading and unloading curves with r3 = 1.0 MPa.
Minimum principal stress (MPa) Fig. 6. Restrictive condition between Young’s modulus and confining pressure.
Eqs. (13) and (20) show that the deterioration of Young’s modulus has no influence on the stress of rock mass, so, only the deformation of rock mass is studied here. The limitations of f for different cases are calculated using Eqs. (28) and (29), as shown in Table 3. In this problem, the residual region is determined by the stress state on the elasto-plastic boundary. And the stresses are only related with the strength properties of rock mass both in the elastic and plastic region. This research result is the same as the former study [9,15], as shown in Figs. 10 and 11. Another two sets of Young’s modulus deterioration coefficient are considered, which satisfies the restrictive condition of Eqs. (28) and (29). Figs. 12 and 13 show the influence of Young’s
66
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
Table 2 Rock mass properties [25]. Types of rock
Hard rock
Soften rock
Radius of opening, R0 (m) Initial stress, p0 (MPa) Supporting pressure, r0 (MPa) Young’s modulus, E(GPa) Poisson ratio, v dilation angle, w (°)
1 1 0 50 0.2 0, 30
1 1 0 5 0.2 0, 30
M–C rock mass c (MPa) u (°) cr (MPa) ur (°)
0.173 55 0.061 52
0.276 35 0.055 30
H–B rock mass m s mr sr rc, rcr (MPa)
0.5 1.0E4 0.3 1.0E5 75
0.2 1.0E4 0.05 1.0E5 50
Fig. 12. Dimensionless displacement in M–C hard rock mass.
Table 3 Calculated limitation of f. Types of rock
Limitation of f
Hard rock
Soften rock
H–B
M–C
H–B
M–C
0.726
0.536
0.524
0.414
Fig. 13. Dimensionless displacement in M–C soften rock mass.
Radial and tangential Stresses (MPa)
Fig. 10. Distribution of stresses in M–C hard rock mass.
This paper K.H. Park [13]
Fig. 14. Dimensionless displacement in H–B hard rock mass.
Dimensionless radius r/R 0 (m/m) Fig. 11. Distribution of stresses in H–B soft rock mass.
modulus deterioration with M–C rock mass. When the dilation angle of hard rock mass is 30°, the inner surface deformation displacement increases by 6.09% with the decrease of f from 1.0 to 0.6, and the corresponding value is 6.24% when the rock mass
Fig. 15. Dimensionless displacement in H–B soften rock mass.
Q. Zhang et al. / Theoretical and Applied Fracture Mechanics 60 (2012) 60–67
dilation angle is 0°. For the M–C soften rock mass, when f decreases from 1.0 to 0.5, the air surface displacement increases by 24.90% and 26.37% with dilation angle of 30° and 0°, respectively. Figs. 14 and 15 show the influence of Young’s modulus deterioration with H–B rock mass. Similarly, the deformation on inner surface increases with the deterioration of Young’s modulus. For the H–B hard rock, the displacement on air surface increase by 1.09% and 1.53% with the decrease of f from 1.0 to 0.8, corresponding to the dilation angle of 30° and 0° respectively. And for H–B soft rock, when f decrease from 1.0 to 0.6, the growth of air surface displacement is 13.24% and 15.59% with the dilation angle of 30° and 0° respectively. Additionally, both of M–C and H–B rock mass show that the rock mass displacement growth slightly with the increasing dilation angle and the Young’s modulus sharply influences the deformation of rock mass than that of the dilation. Another thing should be noted is that the restrictive condition for Young’s modulus deterioration degree also depends on the difference between peak strength and residual strength of rock mass.
[2]
[3] [4] [5]
[6]
[7]
[8] [9] [10]
[11]
6. Conclusion [12]
The analytical solutions considering elasto-plastic coupling effect was deduced for circular openings in elasto-brittle-plastic rock mass governed by Mohr–Coulomb and Hoek–Brown criteria subject to a hydrostatic in suit stress. Meanwhile, the restrictive condition between strength parameters and Young’s modulus was presented, which indicates there is a minimum Young’s modulus deterioration degree. The results were further validated by the numerical test in traditional triaxial compression from loading to unloading using FLAC3D code with the developed EPCM. Finally, the displacements around circular openings were studied with different deterioration degree of Young’s modulus using two kinds of rock mass. The results show that the deterioration of Young’s modulus has little influence on stresses, but increases rock mass deformation sharply in the plastic region. The solutions in this paper are an expansion for the traditional elasto-brittle-plastic and ideal elasto-plastic models. As we know, the calculated deformation of rock mass is usually less than the in situ value, however, the deformation of rock mass increases sharply when the elasto-plastic coupling character is considered, but none do the stresses in the plastic region. Acknowledgements This project is supported by the Fundamental Research Funds for the Central Universities (2012QNB23) and also supported by the National Natural Science Foundation of China (Project Nos. 51174196, 51179185 and 51104151.
[13]
[14]
[15]
[16]
[17] [18]
[19] [20]
[21]
[22]
[23]
[24]
[25]
References [1] E. Hoek, C.T. Carranza-Torres, B. Corkum, Hoek–Brown failure criterion-2002 edition, in: Proceedings of the 5th North American Rock Mechanics
67
Symposium and 17th Tunnelling Association of Canada Conference, Toronto 2002, pp. 267–273. A.L. Florence, L.E. Schwer, Axisymmetric compression of a Mohr–Coulomb medium around a circular hole, Int. J. Numer. Anal. Meth. Geomech. 2 (1978) 367–379. E.T. Brown, J.W. Bray, B. Ladanyi, E. Hoek, Ground response curves for rock tunnels, J. Geotech. Eng., ASCE 109 (1) (1983) 15–39. Y. Wang, Ground response of circular tunnel in poorly consolidated rock, J. Geotech. Eng. 122 (1996) 703–708. C. Carranza-Torres, Self Similar Analysis of the Elastoplastic Response of Underground Openings in Rock and Effects of Practical Variables, Ph.D. Thesis, University of Minnesota, 1998. C. Carranza-Torres, C. Fairhurst, The elasto-plastic response of underground excavations in rock masses that satisfy the Hoek–Brown failure criterion, Int. J. Rock Mech. Min. Sci. 36 (1999) 777–809. C. Carranza-Torres, Elasto-plastic solution of tunnel problems using the generalized form of the Hoek–Brown failure criterion, Int. J. Rock Mech. Min. Sci. 41 (2004) 480–481. K.H. Park, Y.J. Kim, Analytical solution for a circular opening in an elasticbrittle-plastic rock, Int. J. Rock Mech. Min. Sci. 43 (2006) 616–622. S.K. Sharan, Elastic-brittle-plastic analysis of circular openings in Hoek–Brown media, Int. J. Rock Mech. Min. Sci. 40 (2003) 817–824. S.K. Sharan, Exact and approximate solutions of displacements around circular openings in elastic-brittle-plastic Hoek–Brown rock, Int. J. Rock Mech. Min. Sci. 42 (4) (2005) 524–532. S.K. Sharan, Analytical solutions for stresses and displacements around a circular opening in a generalized Hoek–Brown rock, Int. J. Rock Mech. Min. Sci. 40 (2008) 78–85. E. Alonso, L.R. Alejano, F. Varas, G. Fdez-Manin, C. Carranza-Torres, Ground response curves for rock masses exhibiting strain-softening behavior, Int. J. Numer. Anal. Meth. Geomech. 27 (2003) 1153–1185. K.H. Park, B. Tontavanich, J.G. Lee, A simple procedure for ground response curve of circular tunnel in elastic-strain softening rock masses, Tunnel. Undergr. Space Technol. 23 (2008) 151–159. Y.K. Lee, S. Pietruszczak, A new numerical procedure for elasto-plastic analysis of a circular opening excavated in a strain-softening rock mass, Tunnel. Undergr. Space Technol. 23 (2008) 588–599. S.L. Wang, X.T. Yin, H. Tang, X.R. Ge, A new approach for analyzing circular tunnel in strain-softening rock masses, Int. J. Rock Mech. Min. Sci. 1 (2010) 170–178. Q. Zhang, B.S. Jiang, S.L. Wang, X.R. Ge, H.Q. Zhang, Elasto-plastic analysis of a circular opening in strain-softening rock mass, Int. J. Rock Mech. Min. Sci. 50 (2012) 38–46. F.H. Kulhawy, Stress deformation properties of rock and rock discontinuities, Eng. Geol. 9 (4) (1975) 327–350. E.T. Brown, J.W. Bray, F.J. Santarelli, Influence of stress-dependent elastic moduli on stresses and strains around axisymmetric boreholes, Rock Mech. Rock Eng. 22 (3) (1989) 189–203. M.R. Asef, D.J. Reddish, The impact of confining stress on the rock mass deformation modulus, Geotechnique 52 (4) (2002) 235–241. M. Verman, B. Singh, M.N. Viladkar, J.L. Jethwa, Effect of tunnel depth on modulus of deformation of rock mass, Rock Mech. Rock Eng. 30 (3) (1997) 121–127. R.T. Ewy, N.G.W. Cook, Deformation and fracture around cylindrical openings in rock – I. Observations and analysis of deformations, Int. J. Rock Mech. Min. Sci. Abstr. 27 (5) (1990) 387–407. P.A. Nawrocki, M.B. Dusseault, Modelling of damaged zones around openings using radius-dependent Young’s modulus, Rock Mech. Rock Eng. 28 (4) (1995) 227–239. P.K. Kaiser, A new concept to evaluate tunnel performance-influence of excavation procedure, in: Preceding of the 22nd US Symposium on Rock Mechanics, Cambridge, MA, 1981, pp. 264–271. C.G. Zhang, J.H. Zhao, Q.H. Zhang, X.D. Hu, A new closed-form solution for circular openings modeled by the unified strength theory and radiusdependent Young’s modulus, Comput. Geotech. 42 (2012) 118–128. T. Ogawa, K.Y. Lo, Effects of dilatancy and yield criteria on displacements around tunnels, Can. Geotech. J. 24 (1987) 100–113.