Re-profiling of a squeezing tunnel considering the post-peak behavior of rock mass

International Journal of Rock Mechanics & Mining Sciences 125 (2020) 104153

Contents lists available at ScienceDirect

International Journal of Rock Mechanics and Mining Sciences journal homepage: http://www.elsevier.com/locate/ijrmms

Re-profiling of a squeezing tunnel considering the post-peak behavior of rock mass Kai Guan, Wancheng Zhu *, Xige Liu, Jiong Wei, Leilei Niu Center for Rock Instability and Seismicity Research, School of Resource and Civil Engineering, Northeastern University, Shenyang, 110819, China

A R T I C L E I N F O

A B S T R A C T

Keywords: Squeezing ground Strain-softening rock mass Re-excavation Material point History-dependent behavior

In problems involving large deformation, one of main concerns in rock engineering is the occurrence of a heavily reduction in tunnel section and even shield jamming. In this paper, a rigorous theoretical model is established to characterize the process of tunnel re-profiling works in heavily squeezing ground to re-construct the desired profile, based on the finite strain approach and Lagrangian material description. The post-peak behavior of rock mass and elasto-plastic strains in plastic zone during each excavation process are taken into account. Therefore, the proposed method enables to quantify the repetitive disturbances to a plastified surrounding rock and accumulated degradation of material properties. A numerical procedure is presented, in which the historydependent behavior and path of material point in strain-softening rock mass during the entire excavation pro­ cess are determined via our previous work, in combination with incremental algorithm and spatial interpolation. The discrepancy from the accurate solution induced by the negligence of elastic strains in all plastic zones is investigated by comparing the present results with the existing study. An extensive computations are then carried out to clarify some practical questions, including the effects of post-peak behavior, ground condition and old temporary support on the repetitive re-profiling, the variation of rock pressure developing on a rigid support after excavation, and the required minimum over-excavation size to achieve the desired profile. Also the effectiveness of advancing a pilot tunnel is investigated under different post-peak and geometrical conditions, in order to illustrate the applicability of the construction method in squeezing ground.

1. Introduction In regard to rock mass with large deformability and low strength in deep mining and tunnelling engineering, the hazard of large deforma­ tion can occur easily when the construction techniques and support measures cannot to accommodate the geological environment. The primary manifestations of squeezing ground include the reduction in cross-section of opening and the damage to the installed support, thus re-profiling works and support replacement are necessary to reconstruct the clearance profile after excavation. In particular, if the deformation of rock mass remains unexpected high, shield jamming and cutter-head blocking which seriously delays the construction schedule may take place, and these problems can only be solved by re-excavation of the squeezing rock mass to release tunnel boring machine (TBM).1 Therefore, an accurate prediction of squeezing behavior after excavation and re-profiling works is necessary to the tunnel stability analysis and support design. A wide range of deep rock caverns involving large deformation

underwent re-profiling operations. Yang et al.2 investigated the large deformation failure mechanism of a main ventilation roadway in Xin’An coal mine, and found that the supporting materials were still broken even after several roadway repairs and re-profiling works (see Fig. 1). The large wall convergence of a 10 m diameter hydropower headrace tunnel in India led to the severe violation of the original profile, Hoek3 proposed to re-excavation of the invaded rock mass under the protection of a forepole umbrella. Due to the difficult ground conditions, a tunnel boring machine was trapped in Yacambú-Quibor squeezing tunnel with the maximum cover reaching up to 1270 m, and it was eventually removed by re-excavating the deformed structure manually several years later.4 Other similar large deformation problems necessitating carven re-excavation were also observed in Zhangcun coal mine,5 Uluabat tunnel,6 Saint Martin access adit,7 and Laodongshan tunnel.8 From the perspective of mechanics, the re-profiling work in heavily squeezing rock can be treated as an unloading problem from the deformed state of previous excavation stage in order to obtain the desired size of cross-section. In practice, the re-profiling work is similar to the yielding construction procedures (such as over-excavation and

* Corresponding author. E-mail address: [email protected] (W. Zhu). https://doi.org/10.1016/j.ijrmms.2019.104153 Received 3 April 2019; Received in revised form 30 October 2019; Accepted 1 November 2019 Available online 3 December 2019 1365-1609/© 2019 Elsevier Ltd. All rights reserved.

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

List of symbols

δ, T εr, εt η, ηc

Terms defined in Eq. (27) Total tangential and radial strain Plastic internal variable and corresponding critical value μ Poisson’s ratio σ0 Hydrostatic in-situ stress σa Temporary support pressure exerted on tunnel σD Rock strength σ r, σt Radial and tangential stress σ Rp Radial stress at elasto-plastic interface σs Initial rock pressure acted on a rigid new support ϕ, ϕp, ϕr Friction angle and corresponding peak and residual values ψ , ψ p, ψ r Dilation angle and corresponding peak and residual values ω, ωp, ωr Strength parameter and its peak and residual values ω1, ω2, ω3 Coefficients in strain solution Ω1, Ω2 Coefficients in displacement

a0, a ~0 a c, cp, cr C12 E f(x,y) m ne, np r r(0) Re Rp, Rp0

Initial and current tunnel radius Approximate value of initial tunnel radius Cohesion and corresponding peak and residual values Function of material locations Young’s modulus Function in power series defined in Eq. (28) Function of friction angle Number of small enough annulus in elastic and plastic zone Location of material point Initial undisturbed location Radius of research domain Elasto-plastic interface and corresponding previous location Rp ini , Rpa, Rpb Guess values of plastic radius Tol Calculation tolerance Total tunnel convergence ua uRp0 Total displacement at the previous location of elastoplastic interface β Function of dilation angle Δ Increment of parameter Radial stress increment in plastic zone Δσ r

Superscripts e Elastic component p Plastic component (n) nth excavation stage value Subscripts (j) jth material point value

tunnelling, and obtained stress expression around a circular tunnel considering the rheology of the rock mass. Later, this method was extended to the cases of non-circular profile,12,13 and twin tunnels excavation.14 Peng et al.15 investigated the effects of total stage and excavation thickness on rock behavior using the Schwarz alternative method. However, the above solutions only acquired the elastic response of rock mass, hence disregarding the irreversible energy dissipation involving in the plastic deformation since the rock mass yields in elasto-plastic manner. Moreover, rock mass in squeezing con­ ditions is deforming after excavation and presents the significant char­ acteristic of geometric non-linearity, so the above works built in the context of undeformed configuration instead of the actual deformed one are not suitable to deal with re-profiling operation of tunnel involving large deformation in this paper. The heavily squeezing phenomena in geotechnical engineering can be characterized accurately by finite strain elasto-plasticity.16 Recently, Vrakas and Anagnostou17 presented a finite strain theoretical analysis of the squeezing response to tunnel repetitive re-profiling using the ma­ terial coordinates, and clarified some practical problems related to the ground response curve. However, the elastic strains in their model were all ignored at both initial excavation and subsequent re-profiling stages, and this extreme simplification may create larger discrepancy from the accurate prediction along with the further re-profiling works.18 In addition, the rock behavior in Vraka’s solution was treated as brittle–­ plastic with only cohesion reduced directly from peak value to residual one, thus leading to neglect the progressive strength degradation41 and strain-softening behavior at each material point due to the accumulated plastic strains during repetitive re-profiling works. In general, the rock mass exhibits the strain-softening behavior in post-peak stage, and all of the strength parameters (i.e., cohesion, friction angle and dilatancy) vary with the deformation history.19–23 Kova’ri24 and Vrakas and Ana­ gnostou17 studied the effectiveness of an advancing pilot tunnel in squeezing ground, and concluded that this construction method was not beneficial and economic with respect to the future enlarged full cross-section, which was different from engineering experience.25,26 This inconsistency also motivate the present authors to revisit this problem considering different conditions. This paper extends the tunnel re-profiling model of Vrakas and Anagnostou17 to take into account both of the elastic strains involved in

Fig. 1. Re-profiling of a roadway in Xin’an coal mine.2

advancing pilot tunnel), accounting for that all of them are closely related to the process of repetitive excavation.9 Several theoretical an­ alyses had been presented to investigate the tunnel behavior following repetitive excavation disturbances. Wang et al.10,11 assumed the tunnel radius was time-dependent and grew with the excavation stage to characterize the enlargement process of sequential excavation in 2

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

Fig. 2. Model setup for re-profiling of a squeezing tunnel (Left picture from Vrakas and Anagnostou30).

Fig. 3. Schematic illustration of the repetitive re-profiling process in squeezing rock mass (Note: Dash line denotes the configuration before excavation).

the whole process of (repetitive) re-profiling works and the accumulated degradation of material properties induced by the repetitive excavation disturbances by considering the strain-softening behavior of rock mass. Emphasis of the present study is put on the development of a more rigorous history-dependent theoretical analysis of material responses under the re-profiling works in highly deformed ground conditions. This paper is organized as follows. Section 2 gives the problem description and theoretical model; Section 3 provides the formulation depicting the rock behavior after initial excavation and subsequent re-profiling works

based on the finite strain approach, and the numerical implementation is also presented; Section 4 makes comparison with the existing study and deals with some practical problems considering different post-peak behavior, including the rock pressure developing on a rigid support, the amount of over-excavation to acquire the desired profile, and the motion and stress paths of material points; Section 5 investigates the effectiveness and applicability of driving a pilot tunnel first under different conditions, and Section 6 gives some conclusions.

3

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

where εpt and εpr denotes the tangential and radial plastic strain, respectively. It is worth pointing out, however, that for the material points un­ dergoing the repetitive excavation disturbances, the plastic internal variable η in Eq. (1) should be calculated based on the total plastic strains, because the irreversible plastic strains are accumulated gradu­ ally along with the increasing re-profiling step. The degradation law of strength parameters should be determined by laboratory test, however, a bilinear one is defined as follows for simplicity33–35 (

ωðηÞ ¼

ωp

ωp

ωr

�η

ηc

ωr

0 � η < ηc

(2)

η � ηc

where ω represents one of cohesion c, friction angle ϕ and dilation angle ψ , and the value of ω depends on η. ωp and ωr are the initial (peak) and residual value of a strength parameter. ηc is the critical plastic internal variable from which the residual behavior is first observed, and it de­ termines different constitutive relations which include perfectly-plastic model with ηc being 0, strain-softening model and brittle-plastic model with ηc being infinity.35 Alonso et al.36 pointed out that the magnitude of ηc generally vary between 0.001 and 1. The stresses σt and σ r in plastic zone follow the Mohr-Coulomb failure criterion given by

σt ¼ mðηÞσr þ σD ðηÞ (3)

1 þ sin ϕðηÞ 2cðηÞcos ϕðηÞ mðηÞ ¼ ; σ D ðηÞ ¼ 1 sin ϕðηÞ 1 sin ϕðηÞ Fig. 4. Movement of material points in the first excavation stage: (a) the cur­ rent location after excavation; (b) inversion of the initial location before excavation.

where m(η) and σD(η) are related to the strength parameters and plastic internal variable, and mp, σ Dp and mr, σ Dr denote the peak and residual values, respectively. According to non-associated plastic flow rule, the dilatancy coeffi­ cient β(η) defines the ratio of principal plastic strain rate, that is

2. Problem description and mechanical model Fig. 2 shows the re-profiling model of a circular tunnel in squeezing ground under plane strain, where the shaded area denotes the plastified rock mass needed to be further excavated to acquire the desired profile. The rock mass subjected to an initial hydrostatic stress σ 0 is considered as homogeneous, and the uniform support pressure σ a is exerted on tunnel wall, with Lagrangian material coordinate (r, t) used in this problem. In the following analysis, sign convention is defined as positive for compressive. In order to simplify the theoretical analysis, the time-dependent behavior of squeezing rock is ignored in this study, which can be applied to some cases with respect to the short-term large deformation.27–29 Moreover, the assumption is made that the rotational principal directions keep fixed due to the axisymmetric problem here, which leads to a simplified elasto-plastic constitutive law relating Cau­ chy stresses (acting on the deformed configuration) to deformation tensor within the framework of finite strain using a hypoelastic-plastic model,30 and also provides convenience for theoretical derivation due to the kinematically determinate problem under this condition. Ac­ cording to the rigorous numerical simulation by Vrakas31 and Schuerch and Anagnostou,32 the rock response based on this assumption can still provide good borderlines and average approximation for the real solu­ tion. It is noted that despite of some errors induced by the above as­ sumptions, the approach described here attempts to deal with the effect of the most relevant factor33 (i.e., post-peak behavior) on squeezing responses of rock mass subjected to excavation disturbances. The strain-softening behavior of material points is considered by adopting the plastic internal variable η, which controls the evolution process of the strength parameters and post-peak yield surface according to plastic deformation. η can be written as34

η ¼ εpt

εpr

ε_ pr 1 þ sin ψ ðηÞ ¼ ε_ pt 1 sin ψ ðηÞ

βðηÞ ¼

(4)

The repetitive re-profiling works are illustrated schematically in Fig. 3. The excavation process can be modelled by gradually decreasing the rock pressure σs on the desired excavation boundary a0 to a given support pressure σ a. During the initial excavation stage, the rock mass is

unloaded from the in-situ stress σ 0 (i.e., σs ¼ σ 0 ) until the failure cri­ terion Eq. (3) is attained, and then a increasing thick plastic ring is ð1Þ

developed as the support pressure equals σa . The strength parameters of material points in rock mass are updated after excavation, and an equilibrium state of stresses and strains developing on the current deformed configuration is achieved at the end of the first excavation, ð1Þ

with the tunnel wall squeezed from radius a0 to að1Þ and the plastic radius being Rð1Þ p . It is noted that this equilibrium state constitutes the initial condition of rock mass in the next excavation stage if the desired profile is not obtained and the re-profiling works become necessary. In this case, the deformed tunnel should be enlarged and the rock mass ð1Þ

within the (partially) plastified range ½að1Þ ; a0 � will be excavated in the ð2Þ

next stage, which is characterized by unloading of the rock pressure σs

ð2Þ

corresponding to the radial stress at r ¼

ð2Þ a0

at the end of the first stage.

In general, a0 equals a0 except for the special construction procedures of over-excavation of deformed tunnel as well as advancing a pilot tunnel with respect to the subsequent extension to the full cross-section ð2Þ

ð1Þ

(i.e., a0 > a0 ). A new equilibrium state and updated material prop­ erties are obtained after the redistribution of stress and strain field under ð2Þ

ð1Þ

the support pressure σ a . The current radius of tunnel increases to að2Þ while the plastic zone is extended to the radius Rð2Þ due to the accu­ p mulated disturbances to the surrounding rock. The above re-profiling of ð2Þ

(1) 4

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

in Eq. (6) are always fixed.

Eq. (6) yields the relation between a0 and a(1) when r(1) ¼ a(1): ð1Þ

(7)

ð1Þ

a0 ¼ að1Þ þ Δuð1Þ a

where Δua represents the tunnel convergence increment after the first excavation. ð1Þ

3.2. The re-profiling works Accounting for the generality of re-profiling works, the following analysis is based on the nth excavation with n � 2. 3.2.1. Derivation of the history-dependent rock behavior in plastic zone 3.2.1.1. Stress and deformation of material points. The total displace­ ment u(n) of material point induced by the total of n excavation stages, can be calculated by the difference between the initial locations r(0) and the final location r(n) after the nth excavation, that is

Fig. 5. Concentric annuli in plastic zone after the nth excavation stage.

uðnÞ ¼ rð0Þ

a deformed plastified rock mass is repetitive until the desired profile is maintained.

The location of material point before the nth excavation is donated

by r0 , which is also the final location r(n stage, then we have ðnÞ

3. Finite strain elasto-plastic analysis of repetitive excavation disturbances

ðnÞ

r0 ¼ rðn

The radius of research domain is set as Re (see Fig. 4(a)), which

(npþ1) material points which the radial stress increment Δσ r constant between adjacent points, that is35

ð1Þ

np

1Þ

rðnÞ ¼ rðn

rðnÞ ¼ uðnÞ

uðn

(10)

1Þ

where u(n 1) is the total displacement after the n-1 excavation stages. The numerical method solving the strain-softening behavior of rock mass after re-profiling work is similar to that adopted in the first exca­ vation. The plastic zone is divided into np concentric annuli containing npþ1 material points as shown in Fig. 5, and the calculation process starts from the outmost one to the innermost one.

keeps

For the jth annulus, the outer and inner radius is rj

ðnÞ

(5)

and rðjþ1Þ , ðnÞ

respectively. It is noted that the notation of rj denotes the final location ðnÞ

of the jth material point after the nth excavation. The constant radial stress increment is also applied for adjacent material points, i.e.

where σ Rp is the radial stress at the elasto-plastic interface with σRp ¼ ð2σ 0 σDp Þ =ðmp þ 1Þ, which is independent of the excavation stages. The rock response to the first excavation in a squeezing strainsoftening rock mass can be determined conveniently by our previous work.23 More specifically, the calculation starts from the outermost annulus to the innermost one in plastic zone incorporating the constant radial stress increment Eq. (5) between adjacent points (see Fig. 4(a)), then the stresses, displacement and strains of material point in each small enough plastic annulus can be given in closed forms, which are similar to the existing expressions for the rock mass with brittle-plastic behavior,23 while the strength parameters are updated once the plastic internal variable η in Eq. (1) is known. It is noted, however, that ac­ counting for the axial in-situ stress q ignored in this study, the value of q is set to be 20 when our previous work23 is adopted here. Thus, the current location r(1), updated strength parameters (cohesion c(1), friction angle ϕ(1) and dilation angle ψ (1)), displacement increment Δu(1), stress σ (1) and strain ε(1) of each material points are all obtained after the first excavation via the method presented in Guan et al..23 As shown in Fig. 4(b), the initial location of material point in un­ disturbed rock mass can be obtained as rð0Þ ¼ rð1Þ þ Δuð1Þ

(9)

ðnÞ

should be much larger than accounting for the expanding plastic ring in the process of subsequent re-profiling works. At the end of the first excavation, it is assumed that the elastic zone consists of (neþ1) equally spaced material points, while the plastic zone is composed of

σ Rp

at the end of the (n-1)th

1Þ

ΔuðnÞ ¼ r0

ð1Þ Rp

σ ð1Þ a

1)

According Eq. (8) and Eq. (9), the induced displacement increment during the nth excavation can be calculated by

3.1. The first excavation stage

Δσð1Þ r ¼

(8)

rðnÞ

ΔσðnÞ r ¼

σðnÞ a

σ Rp

(11)

np

thus leading to the following expression to estimate the radial stress: (12)

ðnÞ ðnÞ σ ðnÞ rðjþ1Þ ¼ σ rðjÞ þ Δσ r

where σ rðjþ1Þ and σ rðjÞ corresponds to the radial stress of material points ðnÞ

at rj

ðnÞ

ðnÞ

and rðjþ1Þ , respectively. ðnÞ

It can be observed from Eq. (11) that the radial stress increment Δσr

ðnÞ

may be different for each excavation stage if the support pressure σ a exerted on the tunnel wall is changed. According to the Mohr-Coulomb failure criterion Eq. (3), the ðnÞ

tangential stress at rðnÞ ¼ rðjþ1Þ can be estimated by ðnÞ

(13)

ðnÞ ðnÞ ðnÞ σ ðnÞ tðjþ1Þ ¼ mj σ rðjþ1Þ þ σ DðjÞ

(6)

where mj

ðnÞ

Although material points are moving continuously towards the tunnel wall during the repetitive excavations, their initial locations r(0)

1þsinϕj

ðnÞ

¼

1 sinϕj

ðnÞ

, σ DðjÞ ¼ ðnÞ

2cj cosϕj ðnÞ

ðnÞ

1 sinϕj

ðnÞ

, and cj

friction angle of material point at rðnÞ ¼ The final location 5

ðnÞ rðjþ1Þ

ðnÞ

ðnÞ rj ,

and ϕj

ðnÞ

is cohesion and

respectively.

can be expressed as23

K. Guan et al.

ðnÞ

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

� ! .� ðnÞ ðnÞ σðnÞ 1 mj rðjÞ þ σ DðjÞ � ðnÞ � ðnÞ ðnÞ mj σðnÞ 1 rðjÞ þ Δσ r þ σ DðjÞ

ðnÞ

rðjþ1Þ ¼ rj

In combination with Eqs.23 and 24, Eq. (19) leads to the following expression with respect to displacement increment and current location

(14) pðnÞ

The relation between the plastic radial strain increment Δεr pðnÞ t

of material point at rðnÞ ¼ rðjþ1Þ : ðnÞ

and the

ðnÞ

where

ðnÞ βj

1þsinψ j

ðnÞ

¼

1 sinψ j

ðnÞ

, with

ψ ðnÞ j

being the dilation angle at r

ðnÞ

¼

Since total strain increment Δε can be decomposed plastic Δε and elastic ΔεeðnÞ ones, Eq. (15) yields ðnÞ ðnÞ ðnÞ r þ βj Δ t

ε ¼ ΔεeðnÞ þ r

ðnÞ

ðnÞ C12ðjÞ

ε

ΔεeðnÞ r ΔεeðnÞ t

1 þ μ� ¼ ð1 E

μÞΔσðnÞ r μÞΔσðnÞ t

σ ðn r σ ðn t

Δσ ðnÞ ¼ σðnÞ t t

ðnÞ

ωðnÞ 2ðjÞ

1Þ

ðnÞ

!mðnÞ j

ðnÞ

Ω2ðjÞ N! N þ

ðnÞ � δj

ωðnÞ 3ðjÞ þ ð1 þ μÞð1

1

(27)

h � �NþδðnÞ j ðnÞ T rj

i

(28)

1

� .� � ðnÞ ðnÞ ðnÞ mj 2μÞ 1 þ βj σ DðjÞ

�

ðnÞ

ðnÞ

ðnÞ

�h

.� ðnÞ ðnÞ σðnÞ mj rðjþ1Þ þ σ DðjÞ

�i . o 1 E

�i . 1 E

ðnÞ

(30)

ðnÞ

ðnÞ

ðn 1Þ

previous location rðjþ1Þ of the material point with the current location

1Þ

ðnÞ rðjþ1Þ

ðnÞ

(19)

ðnÞ

ðnÞ

� � ��i ðnÞ1 β þ1 ðnÞ j 1; T rj

(31)

ðn 1Þ

where C12ðjÞ is also related to rðjþ1Þ (see Eq. (26)). Therefore, Eq. (31)

should be solved numerically using the least square method. ðn 1Þ

Once rðjþ1Þ is obtained, the initial undisturbed location rðjþ1Þ of the

(20)

ðnÞ

(21)

material information after the (n-1)th excavation using interpolation technology. The total displacement after n excavation stages is given by Eq. (8) as

and

uðjþ1Þ ¼ rðjþ1Þ

ðnÞ

ðn 1Þ

ð0Þ

εðnÞ t ¼

n X

εðnÞ r ¼

ðnÞ

rð0Þ rð1Þ rðn 1Þ rð0Þ þ ln ð2Þ þ ⋅ ⋅ ⋅ þ ln ðnÞ ¼ ln ðnÞ ð1Þ r r r r

(33)

n X

ΔεðiÞ r ¼ ln

drð0Þ drð1Þ drðn 1Þ drð0Þ þ ln ð2Þ þ ⋅ ⋅ ⋅ þ ln ðnÞ ¼ ln ðnÞ ð1Þ dr dr dr dr

(34)

i¼1

It can be seen from Eqs.33 and 34 that, the total strains of material point depend on the initial location r(0) and current location r(n), which is different from the definition of the incremental strains in Eqs.23 and 24.

can be estimated by spatial interpolation

since the parameters of all material points at the end of the previous

Eq. (33) yields the total tangential stain at rðnÞ ¼ rðjþ1Þ : ðnÞ

excavation are known, and then ω3ðjÞ in Eq. (22) is obtained. ðnÞ

Logarithmic strain is adopted here, and the incremental total strains

ð0Þ

εðnÞ tðjþ1Þ ¼ ln

ðnÞ

� � dð ​ ΔuðnÞ Þ drðn 1Þ ΔεðnÞ ¼ ln ðnÞ ¼ ln 1 þ r drðnÞ dr

ΔεðiÞ t ¼ ln

i¼1

of material point at rj is determined If the previous location rj (which will be dealt with in the following analysis), the corresponding

and Δεr during the nth stage are given as17 � � ΔuðnÞ rðn 1Þ ¼ ln ðnÞ ΔεðnÞ ¼ ln 1 þ t ðnÞ r r

(32)

ðnÞ

rðjþ1Þ

The total strains are represented as the sum of incremental total strains during each stage, i.e.

(22)

ðnÞ ðn 1Þ ðnÞ ðn 1Þ ωðnÞ þ ω2ðjÞ σtðjÞ 3ðjÞ ¼ ω1ðjÞ σ rðjÞ

ð0Þ

material point at the current location rðjþ1Þ can be predicted based on the

ðnÞ

ðn 1Þ

ðnÞ

δj Ω1ðjÞ f j

ðnÞ

1Þ

and σtðjÞ

is determined according to Eq. (25), that is

h ðn 1Þ ðnÞ ðnÞ rðjþ1Þ ¼ rðjþ1Þ C12ðjÞ

before the nth excavation with its current location rj . If the number of material points in plastic zone npþ1 in Fig. 5 is sufficiently large, the error by this manipulation can be ignored. Therefore, Eq. (21) becomes

ðnÞ Δε t

h

(26)

ðnÞ

rðjþ1Þ

Taking into account Eq. (10) at rðjþ1Þ (i.e., Δuðjþ1Þ þ rðjþ1Þ ¼ rðjþ1Þ ), the

(18)

þ ω2ðjÞ σðn t

ðn 1Þ

∞ � �� X ðnÞ 1; T rj ¼

Ω2ðjÞ ¼ ω1ðjÞ þ mj ω2ðjÞ

in Eq. (21) are approximated by the ones at the previous location

stress state σ rðjÞ

�

!βðnÞ þ1 j

ðn 1Þ

1Þ

ðn 1Þ

σ

(25)

(29)

In order to simplify the calculation process, the stresses σr

ðn 1Þ t ðn 1Þ rj

� ��i ðnÞ1 β þ1 ðnÞ j 1; T rj

μΔ σ r

� � ðnÞ ¼ ð1 þ μÞ 1 μ βj μ � � ðnÞ ðnÞ ¼ ð1 þ μÞ βj βj μ μ

ðnÞ ðn ωðnÞ 3ðjÞ ¼ ω1ðjÞ σ r

ðn 1Þ

¼

ðnÞ

rðjþ1Þ

n ðnÞ Ω1ðjÞ ¼ exp

where

ω

�

ðnÞ

rj

ðnÞ � � rj ðnÞ ; T rj ¼ ðnÞ ¼ ðnÞ mj þ 1 rðjþ1Þ

fj

(17)

Substituting Eqs.17 and 18 into Eq. (16) gives � i 1 h� ðnÞ ðnÞ ðnÞ ðnÞ ðnÞ ΔεðnÞ ω1ðjÞ σ r þ ωðnÞ ωðnÞ r þ βj Δεt ¼ 2ðjÞ σ t 3ðjÞ E

ðnÞ 1ðjÞ

!βðnÞ þ1 j

N¼0

� ðnÞ

ðnÞ

ðnÞ

where the stress increment ΔσðnÞ during the nth excavation is calculated as the difference between the current stress σðnÞ of material point and its corresponding stress σ ðn 1Þ after the previous excavation stage, that is ðnÞ Δσ ðnÞ r ¼ σr

¼

� ðnÞ

μΔ σ t

ðnÞ

δj Ω1ðjÞ f j

βj þ 1

ðnÞ δj

For each single stage of excavation, Hooke’s law can be written in the following incremental form37: 1 þ μ� ¼ ð1 E

ðnÞ

Δuj þ rj

pðnÞ

(16)

ðnÞ βj Δ eðnÞ t

h ðnÞ ¼ C12ðjÞ

where

ðnÞ rj .

ðnÞ

Δε

ðnÞ

rðjþ1Þ

(15)

ðnÞ

ΔεpðnÞ þ βj ΔεpðnÞ ¼0 r t

ðnÞ

Δuðjþ1Þ þ rðjþ1Þ

plastic tangential strain increment Δε induced by the current exca­ vation stage for the jth annulus is given by Eq. (4):

(23)

rðjþ1Þ ðnÞ

rðjþ1Þ

(35)

Considering Eq. (18) and the Hooke’s law Eq. (17), the elastic ones of total stain are deduced as

(24)

6

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

ðn 1Þ

Fig. 6. Inversion of the material point locations before the nth excavation stage based on the current material state. ① The previous location rðjþ1Þ before the nth ðn 2Þ

ðn 1Þ

stage is calculated via Eq. (31) and the current location rðjþ1Þ . ② The material point location rðjþ1Þ after the (n-2)th stage is interpolated based on rðjþ1Þ and other ðnÞ

existing material points. ③ Computation is repetitive using spatial interpolation until the initial location rðjþ1Þ is obtained. ð0Þ

n X

εeðnÞ ¼ t

ΔεeðiÞ ¼ t

i¼1

" 1þμ ð1 E ¼

eðnÞ r

ε

n X

¼

Δε

eðiÞ r

i¼1

μÞ

ΔσðiÞ t

μ

i¼1

1þμ� ð1 E

" 1þμ ¼ ð1 E

n X

n X

#

ðn 1Þ

Fig. 6). After the material location rðjþ1Þ at the end of the (n-1)th exca­

Δσ ðiÞ r

σ0

Δσ

ðiÞ r

�

μ

i¼1

μ σ ðnÞ r n X

ðn 1Þ

vation is calculated via Eq. (31), its previous location r0ðjþ1Þ before this

(36)

i¼1

μÞ σ ðnÞ t μÞ

¼

eðnÞ

n X

σ0

excavation (red line in the (n-2)th stage of Fig. 6), which also equals

��

ðn 2Þ

rðjþ1Þ according to Eq. (9), can be obtained by spatial interpolation since

the locations of the given material points during the (n-2)th excavation are all known (black line in the (n-2)th stage of Fig. 6). The computation

# Δσ

ðiÞ r

procedure is repetitive until the undisturbed initial location rðjþ1Þ of ð0Þ

i¼1

1 þ μ� ð1 E

μÞ σðnÞ r

σ0

�

μ σ ðnÞ t

σ0

material point is obtained. It should be noted that, for the material point

��

at rðnÞ ¼ rðjþ1Þ , its stresses, strains, displacement and strength parameters ðnÞ

(37)

before the current excavation stage can be also estimated in the similar way. According to Eq. (2), the accumulated degradation law of strength

eðnÞ

where εt and εr are the total elastic tangential stain and total elastic radial stain induced by n excavation disturbances, respectively.

parameters ω (i.e., cohesion cðjþ1Þ , friction angle ϕðjþ1Þ and dilation angle ðnÞ

ðnÞ ðnÞ ψ ðnÞ ¼ rðjþ1Þ after the total of n excavation ðjþ1Þ ) at the material point r

3.2.1.2. Accumulated degradation law of material strength. For a single excavation stage, the conventional definition in Eq. (1) can be used to

disturbances is given as follows

calculate the plastic internal variable increment Δηðjþ1Þ of the material ðnÞ

point at Δη

ðnÞ ðjþ1Þ

ðnÞ rðjþ1Þ ,

¼ Δε

�

that is

pðnÞ tðjþ1Þ

Δε

(38)

Substituting Eq. (15) into Eq. (38), we now have � � � �� � ðnÞ ðnÞ pðnÞ ðnÞ ðnÞ eðnÞ Δηðjþ1Þ ¼ 1 þ βj Δεtðjþ1Þ ¼ 1 þ βj Δεtðjþ1Þ Δεtðjþ1Þ

i¼1

�

εeðnÞ tðjþ1Þ þ

n h X

ðiÞ

βj

�

ðiÞ

Δεtðjþ1Þ

ωr

� ηðjþ1Þ

> :

ηðnÞ c

ðnÞ

0 � ηðjþ1Þ < ηðnÞ c ðnÞ ðjþ1Þ

η

eðiÞ

�i

Δεtðjþ1Þ

ðnÞ

denotes the critical plastic internal variable for the nth excavation. The critical plastic internal variable ηc in Eq. (41) may be different for each stage of excavation, and tends to be larger along with the increasing excavation stages (see Fig. 7), because more released energy by the further excavation may be used to rock dilatancy instead of rock failure. This prediction will be investigated in our future work.

(39)

(40)

i¼1

ðn 3.2.1.3. Relation among Rp0 , RðnÞ p and Rp ðnÞ

eðnÞ

where εtðjþ1Þ and εtðjþ1Þ are obtained by Eq. (33) and Eq. (36), ðnÞ

respectively.

It can be observed that, in order to calculate

ηðnÞ ðjþ1Þ

the elasto-plastic interface (i.e.,

¼

ðnÞ Rp

. For the material point at

as shown in Fig. 5), the total

ðnÞ

in Eq. (40), both of

uRp ¼ RðnÞ p

point with rðnÞ ¼ rðjþ1Þ must be estimated in advance. Namely, the

where G ¼ E=½2ð1 þ μÞ�.

well as the evolution of strength parameters et al.) at rðnÞ ¼

as

ðnÞ

ðnÞ

history-dependent material behavior (e.g., the motion and stress path, as ðnÞ rðjþ1Þ

ðnÞ r1

1Þ

displacement uRp after the n excavation stages can be obtained:

the total tangential stain increment Δεtðjþ1Þ and total elastic tangential stain increment Δεetðjþ1Þ during each excavation stage at the material

should

be tracked in the process of repetitive excavation disturbances, thus Δεtðjþ1Þ and Δεetðjþ1Þ can be then calculated by Eq. (23) and Eq. (17),

σ0

σ Rp

(42)

2G 0ðnÞ

From Eq. (8) and Eq. (42), the initial undisturbed location Rp

� σ0 ðnÞ ðnÞ R0ðnÞ 1þ ¼ RðnÞ p p þ uRp ¼ Rp

respectively. It is worth pointing out that for the nth excavation stage,

reads

σ Rp �

(43)

2G

The previous location of the material point at r1 ¼ RðnÞ p before the ðnÞ

Δεtðjþ1Þ and Δεrðjþ1Þ are determined directly by combining Eq. (31), Eq. ðnÞ

�η

(41)

ðnÞ c

ðnÞ

mined as follows

� ðiÞ ðnÞ Δηðjþ1Þ ¼ εtðjþ1Þ

ðnÞ

ωp

where ηðjþ1Þ is computed by Eq. (40) using spatial interpolation, and ηc

ðnÞ

n X

8 > < ωp

ωr

The total plastic internal variable increment at rðnÞ ¼ rðjþ1Þ is deter­

ηðnÞ ðjþ1Þ ¼

�

ω ηðnÞ ðjþ1Þ ¼

pðnÞ rðjþ1Þ

ðnÞ

ðnÞ

ðnÞ Rp0 .

nth excavation stage is denoted as It is apparent that Rp0 was located in the elastic zone in the previous excavation stage, so the

(23) and Eq. (17). Fig. 6 illustrates the method to inversion of history-dependent vari­ able based on the current material state, taking the estimation of ma­ terial point location in each excavation stage for example (red line in

ðnÞ

induced total displacement uRp0 of the material point at r1 ¼ RðnÞ p after the n-1 excavation stages can be computed by the Kirsch’s solution: ðnÞ

7

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

uRp0 ¼

σ0

σ Rp 2G

Rðn p

1Þ

!2 (44)

ðnÞ

Rp0

ðnÞ

Rp0

ðn 1Þ

where Rp is the plastic radius of rock mass undergoing the total of n-1 excavation disturbances. Therefore, taking account of Eq. (8) and Eq. (44), the initial undis­ 0ðnÞ

turbed location Rp calculated as

of the material point at r1 ¼ Rp ðnÞ

" ðnÞ

ðnÞ

R0ðnÞ ¼ Rp0 þ uRp0 ¼ Rp0 1 þ p

σ0

σ Rp

Rðn p

1Þ

ðnÞ

can be also

!2 # (45)

ðnÞ

2G

Rp0

According to Eq. (43) and Eq. (45), Rp0 is expressed with respect to ðnÞ

RðnÞ p

and (

1Þ Rðn p

as follows

ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi� ffiffiffiffiffiffiffiffiffiffiffiffiffiffi� ffiffiffiffi � � rffihffiffiffiffiffiffiffiffiffi�ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi�i 2 2 2ðσ0 σ Rp Þ σ 0 σ Rp σ 0 σ Rp ðnÞ ðn 1Þ RðnÞ 1 þ þ R 1 þ R p p p 2G 2G G

Fig. 7. Stress–strain curve in plastic softening zone and residual zone under different excavation stages.

)

ðnÞ

Rp0 ¼

2

Fig. 8. Calculation process of material behavior under repetitive excavation disturbances. 8

(46)

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

Table 1 Parameters for comparison example.17 Material parameter

σ ðnÞ rðkÞ ¼ σ 0

Value 1

Excavated tunnel radius in each stage a0 (m) ðnÞ

In-situ stress σ0 (MPa)

ðnÞ

Young’s modulus E (GPa) Poisson’s ratio μ cp (MPa) cr/cp ϕp (� ) ϕr/ϕp ψ p (� ) ψ r/ψ p

σ Rp

σðnÞ tðkÞ ¼ 2σ 0

25 0

Support pressure in each stage σa (MPa)

σ0

ðnÞ

uk ¼

1.25 0.3 0.8 0.4, 1 25 1 0, ϕp 1

� �2 � RðnÞ p ðnÞ

rk

(48)

σ ðnÞ rðkÞ

ð1 þ μÞ σ 0 E

σ Rp

��

�2 RðnÞ p ðnÞ rk ðnÞ rk

where E and μ denotes Young’s modulus and Poisson’s ratio, respectively. The initial material point location rk is calculated as follows ð0Þ

ð0Þ

ðnÞ

(49)

ðnÞ

rk ¼ rk þ uk

Since the initial location rk of material point is always fixed during ð0Þ

ðn 1Þ

If RðnÞ p is given, the previous location

ðnÞ Rp0

of the material point at

ðnÞ r1

repetitive excavations, the previous material point location rk ð0Þ rk

RðnÞ p can be determined by Eq. (46), thus leading to estimation of the material state (including location, stresses, strains, and strength pa­ rameters et al.) during each stage of excavation by adopting the spatial interpolation as shown in Fig. 6. Although Vrakas and Anagnostou17 derived a different expression of

ðnÞ

ðnÞ

of the plastic radius Rp after the nth excavation stage (n � 2). Only if ðnÞ

within the framework of Lagrangian formulation.

ðnÞ Rp

ðnÞ

guess value

by

RðnÞ p ðk ne

1Þ;

1 � k � ne þ 1

¼

ðnÞ Rp ini

ðnÞ

ðnÞ

with Rp

ðnÞ ini

> Rp

ð1Þ

and a tolerance Tol are ðn 1Þ

ðnÞ

material point at tunnel wall (i.e., rðnp þ1Þ ¼ aðnÞ ) before the nth exca­ ðnÞ

ðnÞ

Re

ðnÞ r1

~0 of the assumed in advance, and then the previous location rðnp þ1Þ ¼ a

Hence, the current location of material point at rk after the nth stage can be obtained as

ðnÞ

is given, can the material behavior in plastic zone be solved

sequentially from r1 ¼ Rp to rðnp þ1Þ ¼ aðnÞ (see Fig. 5). Thus, an initial

3.2.2. Derivation of rock behavior in elastic zone The elastic zone is composed of neþ1 equally spaced material points, which is similar to that in the first excavation stage (see Fig. 4(a)). The material point location in the outmost annulus in elastic zone is set as Re .

rk ¼ RðnÞ p þ

(50)

ðnÞ

rk

3.2.3. Numerical implementation The key point of numerical implementation lies in the determination

of the material coordinate, which violated the adopted Lagrangian approach. By contrast, we provide a rigorous and accurate theoretical derivation of

ðn 1Þ

Δuk ¼ rk

Rp0 , the calculation process was based on the spatial coordinate instead

ðnÞ Rp0

can be

interpolated based on and the locations of other existing material points after (n-1)th stage. Therefore, according to Eq. (10), the displacement increment of material point in elastic zone induced by the nth excavation stage equals

¼

~0 should vation stage can be solved iteratively via Eq. (31). In theory, a ðnÞ

ðnÞ a0 .

be equal to the designed tunnel radius If the convergence condition � ðnÞ ðnÞ � ðnÞ ~0 of �a a0 � � Tol is satisfied, the true value of Rp is achieved; otherwise, the initial guess value Rp

ðnÞ ini

(47)

should be adjusted and the above

calculation is repetitive. In order to accelerate convergence during the numerical implementation, the method of bisection is adopted, which

The stresses and total displacement of the material points are given

Rp

ðnÞ ini

is assigned to the mid-value of an given interval which varies with

~0 and a0 . The detail calculation the relative magnitudes between a procedure is illustrated in Fig. 8, while the input parameters including radius of research domain Re, peak and residual values of rock ðnÞ

Fig. 9. Variation of (a) plastic radius and (b) wall convergence during re-profiling works. 9

ðnÞ

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

Fig. 10. Errors of (a) plastic radius and (b) wall convergence induced by neglecting the elastic strain in plastic zone.

Fig. 11. Effects of the varying critical plastic internal variable on (a) plastic radius and (b) wall convergence.

Fig. 12. Effect of the varying critical plastic internal variable on the maximum volumetric strain in surrounding rock.

Fig. 13. The evolution of rock pressure acting on a rigid support to prevent deformation induced by the first excavation.

10

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

K. Guan et al.

Fig. 14. Under unsupported condition, prediction of rock pressure for different excavation stage, according to the tunnel convergence measured in strainsoftening rock mass during the current stage.

Fig. 16. Prediction of the required over-excavation size based on the measured initial excavation-induced wall convergence.

the increasing re-profiling stages, Rp tends to increase and then maintain constant, while the wall convergence, on the contrary, decrease gradu­ ally to zero, which indicates that the squeezing potential of tunnel is reduced due to the re-profiling works, and the desired profile can be

obtained (i.e., Δua ¼ 0 and aðnÞ ¼ a0 ) after several stages of the re­ petitive excavations. However, in the case of convergence rate exceeding 30% after the first excavation (see the red line and blue line in Fig. 9(b)), large deformation phenomenon may also take place despite of three re-excavation stages. Higher dilatancy of rock mass leads to a smaller plastic radius and a larger wall convergence after each excava­ tion stage. It appears evident that, the lower the residual strength of rock mass, the more are the re-profiling works to acquire the desired profile, reaching up to at least eight stages in extreme case of brittle-plastic model with cr/cp ¼ 0.4 (see the blue line in Fig. 9(b)). Comparative computations in Fig. 9 also shows that, in the whole process of repetitive excavation disturbances, the negligence of the elastic strains in all plastic zones overestimates the plastic radius Rp and underestimates the wall convergence Δua of tunnel for each excavation stage. As can be seen in Fig. 10, the induced error of Δua increases generally with the increasing excavation stages and Δua may be underestimated up to 60%. However, the overestimation of Rp is slight since the error of Rp is limited within 10%. It is worthwhile pointing out that, for the perfectly-plastic rock mass with high dilatancy, even if the discrepancy from the accurate solution can be ignored for the first excavation, the prediction of wall convergence Δua is underestimated significantly for the increasing stages of the subsequent re-profiling work (see the blue line in Fig. 10(b)). ðnÞ

Fig. 15. Evolution of over-excavation size along with the critical plastic the internal variable and the support pressure.

properties, total excavation stages Nc, excavated tunnel radius a0 and ðnÞ

support pressure as well as critical plastic internal variable ηc each stage, and tolerance Tol.

σ ðnÞ a

ðnÞ

of

ðnÞ

4. Rock responses to re-profiling works

4.2. Effect of the varying critical plastic internal variable

4.1. Errors induced by ignoring elastic strains in plastic zone of each excavation stage

During re-profiling works, the critical plastic internal variable ηc may be different for different excavation stage (see Fig. 7). The parameters in Table 1 are adopted except for cr ¼ 0.2 MPa, ϕr ¼ 20� , ψ p ¼ 5� and ψ r ¼ 0. The total stage of excavation is set to five. Fig. 11 compares the results of rock deformation under the identical

Vrakas and Anagnostou17 have analyzed the ground response curve by totally neglecting the elastic strains of plastic zone in all stages including the first excavation and the subsequent re-profiling works. The

critical plastic internal variable with ηc ¼ 0:5 and the varying one with ðnÞ

parameters are given in Table 1. Note that ne ¼ np ¼ 500, Re ¼ 5Rp and Tol ¼ 1e 5 are assumed to take into account the accuracy and conver­ gence speed of solutions in the following analysis, while the critical plastic internal variable is set to be zero and infinity to denote the brittle-plastic model and perfectly-plastic model adopted by Vrakas and Anagnostou,17 respectively. Fig. 9 plots the evolution of plastic radius and convergence rate of tunnel wall during repetitive excavations. It is observed that along with ð1Þ

ηð1:5Þ ¼ ½0:5;0:55;0:6;0:65;0:7�. It is apparent that, for the case of varying c ηðnÞ , all of the plastic radius Rp, the residual radius Rr and the wall c

convergence Δua are smaller than those predicted in the case of identical

ηðnÞ c , indicating that the tunnel maintenance and the re-profiling works

are quite effectively in this ground condition. Hence, it highlights the importance of investigating the post-peak behavior of rock mass not only in the first stage but also in the re-profiling stages. 11

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

Fig. 17. Monitoring of the rock behavior at four spatially-fixed locations: (a) location of monitoring point; (b) total plastic internal variable; (c) rock strength; (d) displacement increment; (e) total displacement.

Fig. 12 presents the variation of the maximum volumetric strain εv, in surrounding rock. It is shown that εv,max almost keep fixed regardless of the excavation stages when the critical plastic internal variable ηc remains constant. By contrast, for the case of the varying ηc, εv,max increases directly with the re-profiling stages, which manifests that the energy release induced by the stepwise repetitive excavation lead to the growing volume expansion of rock mass during re-profiling works.

commonly characterized by the resistance of a rigid support measure with sufficient stiffness to prevent the rock deformation to occur after tunnel excavation,9,24 which is acted at the future excavation boundary a0 (see Fig. 3).

max

Fig. 13 shows the developed rock pressure σs under different values ð1Þ

of support pressure σ a and the critical plastic internal variable ηc after the first excavation stage. As can be seen that, along with the increasing ð1Þ

ηc, σ ð1Þ s decreases initially until the residual (failure) zone of surrounding

rock is eliminated (i.e., Rr/a ¼ 1), and it then increases slightly to a

stable value. Moreover, the larger the acted support pressure σa after excavation, the higher load developing on the rigid support measure during the next re-profiling work. Vrakas and Anagnostou17 derived the relation between rock pressure

4.3. Rock pressure acting on a rigid support measure

ð1Þ

For the purpose of simplifying calculations, the identical critical plastic internal variable ηc of rock mass during repetitive excavation stages is assumed in the following analysis. The rock pressure is 12

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

Fig. 18. Tracking mechanical behaviors of three material points based on Lagrangian description: (a) location of material point; (b) total plastic internal variable; (c) rock strength; (d) displacement increment; (e) total displacement.

stain-softening rock mass after each excavation stage, based on the wall

ðnÞ σ ðnÞ s and tunnel convergence increment Δua for two borderline cases, i.

e., perfectly-plastic and brittle-plastic, which is expressed as

σ

ðnÞ s

¼

σ Dp mp

1

�� 1

ΔuðnÞ a a0

�1

mp

� 1

convergence increment Δua measured in the current stage. The support ðnÞ

pressures are all assumed to offer negligible resistance (i.e., σa ¼ 0). The required minimum bearing capacity of a rigid support can be esti­ mated by the chart. If the exerted temporary support (such as linings) ðnÞ

(51)

after each excavation stage is lower than σ s , the tunnel becomes deformed and the extensive re-profiling works may be required due to the insufficient bearing capacity of linings. For a given low tunnel ðnÞ

and

σ ðnÞ s ¼

σ Dr mr

1

�� 1

ΔuðnÞ a a0

�1

mr

� 1

(52)

convergence in stain-softening rock mass (A1B1 and A2B2), σ s is over­ estimated significantly by perfectly-plastic model but underestimated by brittle-plastic model. By contrast, when the tunnel convergence is rela­ tively high (B1C1 and B2C2), the rock behavior can be treated as brittleplastic, and the rock pressure can be calculated straightforwardly by Eq. ðnÞ

, respectively. It can be observed that, the curves plotted by Eqs.51 and 52 are fixed and independent on the excavation stages (see the dash lines in Fig. 14).

The solid lines in Fig. 14 allow to prediction of the rock pressure σ s in ðnÞ

(52). Compared with σs

ð1Þ

13

after the first excavation stage, the rock

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

4.5. The history-dependent material behaviors This part provides some useful insight into the history-dependent behavior and path of material points in surrounding rock in the pro­ cess of the repetitive excavation stage. The strength parameters are given as cp ¼ 4.5 MPa, cr ¼ 2 MPa, ϕp ¼ 25� , ϕr ¼ 8� , ψ p ¼ 5� , ψ r ¼ 0 and ηc ¼ 0.7, while the other parameters are listed in Table 1. Note that the designed excavated radius of each excavation stage is assumed to be identical (a0 ¼ a0 ¼ 1 m). Fig. 17 presents the evolution of total plastic internal variable, strength, displacement increment and total displacement at four spatially-fixed monitoring locations (i.e., P1 is located at the deformed tunnel wall a(1) after the first excavation, P2 is located at the desired profile a0, P3 is located at 2.5a0, and P4 is located at the elasto-plastic ðnÞ

interface Rp , as shown in Fig. 17(a)), based on the Eulerain descrip­ tion. As can be seen that, after the first excavation, the residual zone develops initially from P1 to P2 (Fig. 17(b)), and the total displacements of the four locations increases remarkably (Fig. 17(e)). The rock strength at P4 remains peak condition since no irreversible plastic deformation occurs, thus leading to the minimum (elastic) displacement. At the end of the next excavation stage, the rock mass at P1 is removed, while the location P4 becomes plastified and P3 reaches the residual state following P1 and P2, which indicates that the plastic zone is expanded after the re-profiling work. Accounting for the accumulated excavation disturbances to surrounding rock, the total plastic internal variables at all locations increase constantly with the excavation stage (Fig. 17(b)), which results in the progressive degradation of rock strength (Fig. 17(c)) and the increasing total displacement (Fig. 17(e)). However, their variation rates decrease along with the repetitive excavation stage since the energy release of rock mass induced by re-profiling works become less; hence, the displacement increment during each excavation is decreasing (Fig. 17(d)), which is contrary to the evolution of the total displacement. Moreover, in the case of heavily squeezing ground which requires repetitive repairs,2 the emphasis should be put on the effec­ tiveness of the first stage of re-profiling, which has a more prominent role on the profile maintenance and squeezing control than other further re-profiling stages. A significant improvement over the Eulerain description is that the history-dependent behavior of each material point can be accurately tracked when the material is deforming during the re-profiling works according to the Lagrangian method.38 The changing paths of several parameters for three material points are provided in Fig. 18. The initial undisturbed locations r(0) of C1, C2, and C3 are situated at a0, 2.5a0 and ð1Þ

Fig. 19. Stress path of a material point with r(0) ¼ 3a0 under different critical plastic internal variable during the repetitive excavation disturbance.

pressure σ s can be reduced effectively by re-profiling work (see the solid lines), thus leading to the mitigation of the loads acting on the rigid support and linings. ðnÞ

4.4. The amount of over-excavation to acquire the desired profile after the first excavation In order to acquire the desired profile and improve the repair effi­ ciency following the occurrence of large deformation involved in the first excavation, the tunnel can be over-excavated to a magnitude to accommodate the squeezing ground, which is consistent with the pas­ sive excavation approach.9 The amount of over-excavation Δa is given by ð2Þ

ð2Þ

að2Þ ¼ a0

Δa ¼ a0

(53)

ð1Þ

a0

where the deformed tunnel radius a(2) after over-excavation equals

justly the size of desired profile a0 (i.e., að2Þ ¼ a0 ). The evolution of the over-excavation size along with the critical plastic internal variable and support pressure is presented in Fig. 15. As can be seen that, in the presence of a new negligible support after the ð1Þ

over-excavation (i.e.,

σ ð2Þ a

ð1Þ

Rp , respectively (Fig. 18(a)). It is apparent that the material points are moving towards the tunnel wall due to the increasing total displace­ ments during the repetitive excavation stage (Fig. 18(e)), thus resulting in the material point C1 removed after the first excavation stage. Compared with C3, C2 enters into residual zone earlier (Fig. 18(b)). The evolution of material behavior based on Lagrangian description is similar to that in Eulerain description. More specifically, along with the increasing re-profiling stage, both of the total plastic internal variable and total displacement increase, but the values of strength and the displacement increment of material point decrease gradually. Fig. 19 plots the stress path for a material point with r(0) ¼ 3a0 under different critical plastic internal variable ηc during three stages of re­ petitive excavation. It can be observed that, the bigger the value of ηc, the shorter the stress path, and thus the closer to the peak envelope line will be the stress state. In particular, when ηc equals 0.1, the material point lying initially on the elastic path after the first excavation stage, is now located on a stress path in plastic zone after re-profiling works, with the values of stress component remaining high (see the green line). Furthermore, along with the increasing excavation stage, both of the tangential stress and radial stress decrease and the stress state of ma­ terial point is approaching to the residual envelope line, which manifests that the re-profiling works is also the process of rock de-stressing. ð1Þ

¼ 0), Δa tends to decreases inversely with the

value of ηc, while the state of the old temporary support σ a (e.g., totally destroyed or keep complete) has no effect on the over-excavation size. ð1Þ

When σa ¼ 1 MPa, the value of Δa is reduced remarkably, which in­ dicates that the construction procedure of over-excavation together with a new support measure can decrease effectively both of the excavated volume and economic cost involved in the re-profiling work. Fig. 16 plots the relation between the initial convergence of tunnel and the required size of over-excavation. It can be observed that the ð2Þ

higher the first excavation-induced convergence Δua , the larger volume of rock mass required to be over-excavated in the next stage. If a suffi­ cient support pressure is acted on tunnel after the over-excavation ð1Þ

(σ a ¼ 1 MPa), the over-excavation construction may become unnec­ essary due to the over-excavation size equals zero (i.e., Δa ¼ 0 and ð2Þ

a0 ​ ¼ a0 ), even for the cases of extremely high convergence rate after the first excavation reaching up to 60%. In this condition, the clearance profile can be achieved just by repetitive excavation of the deformed ð2Þ

ð1Þ

tunnel to the desired size a0 . ð1Þ

14

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

Fig. 20. Radial stress, plastic internal variable, rock strength of surrounding rock under different initial temporary support after repetitive excavation (Note: Dash line in this picture donates the desired profile.).

Fig. 20 gives cloud pictures of radial stress, plastic internal variable and rock strength after excavation under different initial support con­

be acquired only by a single stage of re-profiling when σa ¼ 0. When the first excavation-induced rock deformation is constrained by a high ð1Þ

dition σ a , with σ a ¼ 0 and σa ¼ 0, assuming that the re-profiling work after the first excavation stage is inevitable. If the initial tempo­ rary support is totally destroyed by high rock pressure or no support ð1Þ

ð2Þ

ð3Þ

support resistance with σa ¼ 4MPa, both of the stress and deformation are released considerably during the next excavation, which leads to the large contraction of tunnel and necessitates the further re-profiling works. In such case, at least two stages of re-profiling are required to obtain the desired profile. Therefore, in order to maintain tunnel effi­ ciently during the re-profiling works, the concentrated stress and squeezing potential of rock mass should be released in advance during the first excavation, which is beneficial to decrease the number of tunnel repair. ð1Þ

measure is provided (i.e., σa ¼ 0) after the first excavation, the reprofiling works can alleviate the stress concentration, and increase adversely the plastic deformation and strength degradation. However, the tunnel profile is enlarged with the further excavation stage, and radius of tunnel is quite approximated to the desired size a0 after the second excavation stage, which means that the clearance of tunnel can ð1Þ

15

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

controversial.17,24–26 This section investigates this problem by consid­

ering different post-peak behavior and pilot tunnel radius a0 . The ð1Þ

desired radius of main tunnel is assumed to be ð1Þ ð2Þ for the ratio a0 =a0

ð2Þ a0

¼ a0 ¼ 1m. Values

of 0.1, 0.3, 0.5, 0.7 and 0.9 are selected. To describe different post-peak behavior, the critical plastic internal variable ηc ranges from 1e-3 to 10 in this case. Fig. 21 presents the variation of the convergence ratio Δua0 =a0 at ð1Þ

the future excavation boundary tunnel. As can be seen,

ð1Þ ð2Þ Δua0 =a0

ð2Þ a0

ð2Þ

after the excavation of the pilot

decreases directly with the radius a0

ð1Þ

ð1Þ ð2Þ a0 =a0

of pilot tunnel, and its value as low as 2.5% if � 0:3. Therefore, it can be inferred that when the pilot tunnel with small radius is exca­ vated in practice, the stress release and deformation reduction of rock mass in the vicinity of the future main tunnel is not high. However, in the case of large size of pilot tunnel with a0 =a0 � 0:5, the value of ð1Þ

ð1Þ ð2Þ Δua0 =a0

ð2Þ

is more than 5%, which corresponds to the high critical stain

for very severe squeezing problem.3 In such condition, Δua0 increases ð1Þ

significantly along with the increasing value of

Fig. 21. Convergence ratio at the future enlarged boundary after excavation of the pilot tunnel, considering different values of critical plastic internal variable

ð1Þ a0 ,

even reaching up to

20% of the desired radius a0 of main tunnel. To investigate the effectiveness of advancing a pilot tunnel more clearly, comparative analysis of the construction methods without/with pilot tunnel is conducted. Fig. 22 provides the calculated convergence ratio of main tunnel after completing the entire excavation process under different construction method, post-peak behavior and pilot tunnel radius. It can be observed that in the moderate/minor squeezing problems with wall convergence ratio less than 5%,3 the effect of pilot tunnel on deformation reduction of main tunnel can be ignored regardless of the radius of pilot tunnel. However, for the severe squeezing ground with wall convergence ratio more than 5%, the effectiveness of the construction method depends closely on the size of ð2Þ

ηc and the radius að1Þ 0 of pilot tunnel.

pilot tunnel. More specifically, in the case of a0 =a0 � ð0:3e0:5Þ, the wall convergence increment by the construction method without pilot tunnel is almost consistent with the one with pilot tunnel In contrast, ð1Þ

ð2Þ

when a0 =a0 � ð0:3e0:5Þ, the squeezing potential and large deforma­ tion of the main tunnel are alleviated effectively by driving a pilot tunnel ð1Þ

ð2Þ

first, and the maximum reduction of Δua0 =a0 amounts to 17.56%– 54.85% as show in Table 2. Therefore, in the severe squeezing ground (i.e., convergence ratio is ð2Þ

ð2Þ

higher than 5%), if the pilot tunnel size is large (more than (0.3–0.5)a0 in the case of this study, and this range is closely related to the specific geological conditions), the construction method of advancing a pilot tunnel has a significant effect on reducing the squeezing potential and large deformation of the main tunnel. However, in the cases of moder­ ate/minor/no squeezing ground or small radius of tunnel pilot, preconstruction of a pilot tunnel does not seem to be manifested effec­ tively, which explained the conclusion by Vrakas and Anagnostou17 that the construction method is not economic with respect to the excavation ð2Þ

Fig. 22. Convergence ratio of main tunnel after completing the entire exca­ vation process under different construction method, post-peak behavior and

pilot tunnel radius (Note: in the case without pilot tunnel, Δua = a0 is calcu­ ð2Þ

lated as

ð2Þ

ð1Þ Δua =a0 ).

Table 2 In the severe squeezing problem, the maximum convergence reduction of main tunnel induced by driving a pilot tunnel. a0 =a0

0.1

0.3

0.5

0.7

0.9

Maximum convergence reduction

1.91%

6.55%

17.56%

34.50%

54.85%

ð1Þ

ð2Þ

of the future main tunnel with minor pilot tunnel size a0 =a0 ð1Þ

ð2Þ

¼ 1=3.

6. Conclusions This paper develops a rigorous theoretical model to characterize the repetitive tunnel re-profiling works in heavily squeezing ground, based on the finite strain elasto-plasticity and Lagrangian material description. The post-peak behavior of rock mass and elasto-plastic strains in all plastic zones during the entire excavation process are taken into ac­ count. In fact, the re-profiling works are essentially the process of re­ petitive disturbances to the plastified surrounding rock and accumulated degradation of material properties, thus necessitating the consideration of the influence of rock responses in previous excavation stages. The strain-softening behavior of each material point in plastic zone is treated as the brittle-plastic response of a small enough plastic annulus, and then

5. The effectiveness and applicability of an advancing pilot tunnel considering different post-peak behavior The primary role of pilot tunnel lies on ventilation, drainage and advanced geological exploration,39,40 while the construction method of advancing a pilot tunnel with respect to the subsequent enlargement to the full cross-section, is adopted in mining engineering and tunnelling. However, the effectiveness of this procedure in relation to the defor­ mation reduction of the future main tunnel still remain 16

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153

the incremental relations of some variables (including location, stresses, strains and strength parameters et al.) for the adjacent material points are obtained at the current excavation stage, which leads to the deter­ mination of both the current and previous material states, in combina­ tion with the inversion of material history-dependent path and property by spatial interpolation. The corresponding numerical procedure is proposed, which considers the whole unloading history in the process of repetitive excavation. Based on the numerical results, the following conclusions can be drawn:

is beneficial to reducing the squeezing potential and large deformation of the main tunnel. However, in the cases of mod­ erate/minor/no squeezing ground or small radius of tunnel pilot, pre-construction of a pilot tunnel does not seem to be manifested effectively. Declaration of competing interest None.

(1) The squeezing potential of tunnel is reduced due to the reprofiling works. Along with the increasing re-profiling stage, the plastic radius tends to increase and then maintain constant, while the wall convergence, on the contrary, decrease gradually to zero. For rock mass with higher dilatancy or lower residual strength, more re-profiling works are required to achieve the desired profile. (2) The negligence of elastic strains in all plastic zones overestimates slightly the plastic radius and underestimates the wall conver­ gence of tunnel after each excavation stage. For perfectly-plastic rock mass with high dilatancy, although the discrepancy from the accurate solution can be ignored for the first excavation, the underestimation of the wall convergence can become consider­ able for the subsequent stages of excavation. (3) It is necessary to investigate the post-peak behavior of rock mass not only in the first excavation stage but also in re-profiling stages, since the case with the increasing critical plastic internal variable ηc

ðnÞ

(4)

(5)

(6)

(7)

(8)

Acknowledgements This work is funded by the National Key Research and Development Program of China (Grant No. 2016YFC0801607), National Science Foundation of China (Grant Nos. 51525402, 51874069 and 51761135102), the Fundamental Research Funds for the Central Uni­ versities of China (Grant Nos. N170108028, N180101028, N180106003 and N180115009) and Postdoctoral Science Foundation of China (Grant No. 2018M641706). These supports are gratefully acknowledged. Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi. org/10.1016/j.ijrmms.2019.104153. References

predicts lower rock deformation than the one with

the identical value of thus leading to the tunnel repairs more effectively in such ground condition. The rock pressure developing on a rigid support installed after excavation of a squeezing tunnel can be decreased significantly after re-profiling work. For a low tunnel convergence in strainsoftening rock mass, the rock pressure is overestimated by perfectly-plastic model but underestimated by brittle-plastic model. However, when the tunnel convergence is relatively high, the rock pressure can be estimated directly by brittle-plastic model. The higher the first excavation-induced convergence, the larger additional over-excavation size is required to re-establish the desired profile in the next excavation stage, regardless of the state of the old temporary support (e.g., totally destroyed or keep complete). Furthermore, the construction procedure of overexcavation together with a slight new support can considerably reduce the excavated volume of rock mass involved in the reprofiling work. The history-dependent material behaviors are investigated based on Eulerain description and Lagrangian description. Along with the increasing re-profiling stage, both of the total plastic internal variable and total displacement increase, but the values of strength and the displacement increment of material point decrease gradually. In the case of heavily squeezing ground which requires repetitive repairs, the emphasis should be put on the effectiveness of the first stage of re-profiling. The bigger the value of critical plastic internal variable, the shorter the stress path, and thus the closer to the peak envelope line will be the stress state of material point. To stabilize tunnel efficiently during the re-profiling works, the concentrated stress and squeezing potential of rock mass should be released in advance during the first excavation, which is beneficial to decrease the number of tunnel repair. In severe squeezing ground (i.e., convergence ratio is higher than 5%), if the pilot tunnel size is large (more than 30%–50% of the main tunnel radius which depends on the specific geological conditions), the construction method of advancing a pilot tunnel

ηðnÞ c ,

€ Hasanpour R. Estimation of ground pressures on a shielded TBM in 1 Aydan O, tunneling through squeezing ground and its possibility of jamming. Bull Eng Geol Environ. 2019:1–15. 2 Yang S-Q, Chen M, Jing H-W, Chen K-F, Meng B. A case study on large deformation failure mechanism of deep soft rock roadway in Xin’An coal mine. China. Engineering Geology. 2017;217:89–101. 3 Hoek E. Big tunnels in bad rock. J Geotech Geoenviron Eng. 2001;127(9):726–740. 4 Hoek E, Guevara R. Overcoming squeezing in the Yacambú-Quibor tunnel. Venezuela. Rock Mechanics and Rock Engineering. 2009;42(2):389–418. 5 Gao F, Stead D, Kang H. Numerical simulation of squeezing failure in a coal mine roadway due to mining-induced stresses. Rock Mech Rock Eng. 2015;48(4): 1635–1645. 6 Bilgin N, Algan M. The performance of a TBM in a squeezing ground at Uluabat, Turkey. Tunnelling & Underground Space Technology Incorporating Trenchless Technology Research. 2012;32(32):58–65. 7 Barla G, Bonini M, Semeraro M. Analysis of the behaviour of a yield-control support system in squeezing rock. Tunnelling & Underground Space Technology Incorporating Trenchless Technology Research. 2011;26(1):146–154. 8 Cao C, Shi C, Lei M, Yang W, Liu J. Squeezing failure of tunnels: a case study. Tunn Undergr Space Technol. 2018;77:188–203. 9 Barla G. Tunnelling under squeezing rock conditions. Eurosummer-school in tunnel mechanics. Innsbruck. 2001:169–268. 10 Wang HN, Li Y, Ni Q, Utili S, Jiang MJ, Liu F. Analytical solutions for the construction of deeply buried circular tunnels with two liners in rheological rock. Rock Mech Rock Eng. 2013;46(6):1481–1498. 11 Wang HN, Utili S, Jiang MJ. An analytical approach for the sequential excavation of axisymmetric lined tunnels in viscoelastic rock. Int J Rock Mech Min Sci. 2014;68(6): 85–106. 12 Wang HN, Utili S, Jiang MJ, He P. Analytical solutions for tunnels of elliptical crosssection in rheological rock accounting for sequential excavation. Rock Mech Rock Eng. 2015;48(5):1997–2029. 13 Wang HN, Jiang MJ, Zhao T, Zeng GS. Viscoelastic solutions for stresses and displacements around non-circular tunnels sequentially excavated at great depths. Acta Geotechnica. 2018;(1):1–29. 14 Wang HN, Zeng GS, Utili S, Jiang MJ, Wu L. Analytical solutions of stresses and displacements for deeply buried twin tunnels in viscoelastic rock. Int J Rock Mech Min Sci. 2017;93:13–29. 15 Peng N, Zhang Y, Wang C. Study on mechanical properties of surrounding rock during tunnel construction with in-situ expansion step excavation method. DISASTER ADVANCES. 2010;3(4):194–199. 16 Belytschko T, Liu WK, Moran B, Elkhodary K. Nonlinear Finite Elements for Continua and Structures. John wiley & sons; 2013. 17 Vrakas A, Anagnostou G. Ground response to tunnel Re-profiling under heavily squeezing conditions. Rock Mech Rock Eng. 2016;49(7):2753–2762. 18 Cui L, Zheng J-j, Sheng Q, Pan Y. A simplified procedure for the interaction between fully-grouted bolts and rock mass for circular tunnels. Comput Geotech. 2019;106: 177–192. 19 Peng J, Cai M, Rong G, Yao MD, Jiang QH, Zhou CB. Determination of confinement and plastic strain dependent post-peak strength of intact rocks. Eng Geol. 2017:218.

17

K. Guan et al.

International Journal of Rock Mechanics and Mining Sciences 125 (2020) 104153 32 Schuerch R, Anagnostou G. The applicability of the ground response curve to tunnelling problems that violate rotational symmetry. Rock Mech Rock Eng. 2012;45 (1):1–10. 33 Alejano LR, Rodríguez-Dono A, Veiga M. Plastic radii and longitudinal deformation profiles of tunnels excavated in strain-softening rock masses. Tunnelling & Underground Space Technology Incorporating Trenchless Technology Research. 2012;30 (3):169–182. 34 Alejano LR, Alonso E, Rodríguez-Dono A, Fern� andez-Manín G. Application of the convergence-confinement method to tunnels in rock masses exhibiting Hoek–Brown strain-softening behaviour. Int J Rock Mech Min Sci. 2010;47(1):150–160. 35 Lee YK, Pietruszczak S. A new numerical procedure for elasto-plastic analysis of a circular opening excavated in a strain-softening rock mass. Tunnelling & Underground Space Technology Incorporating Trenchless Technology Research. 2008;23(5):588–599. 36 Alonso E, Alejano LR, Varas F, Fdez-Manin G, Carranza-Torres C. Ground response curves for rock masses exhibiting strain-softening behaviour. Int J Numer Anal Methods Geomech. 2003;27(13):1153–1185. 37 Vrakas A, Anagnostou G. A finite strain closed-form solution for the elastoplastic ground response curve in tunnelling. Int J Numer Anal Methods Geomech. 2014;38 (11):1131–1148. 38 Soga K, Alonso E, Yerro A, et al. Trends in large-deformation analysis of landslide mass movements with particular emphasis on the material point method. Geotechnique. 2017;66(3):248–273. 39 Steiner W. Tunnelling in squeezing rocks: case histories. Rock Mech Rock Eng. 1996; 29(4):211–246. 40 Zingg S, Anagnostou G. An investigation into efficient drainage layouts for the stabilization of tunnel faces in homogeneous ground. Tunn Undergr Space Technol. 2016;58:49–73. 41 Wei J, Zhu W, Guan K, Zhou J, Song J-J. An Acoustic Emission Data-Driven Model to Simulate Rock Failure Process. Rock Mechanics and Rock Engineering. 2019:1–17.

20 Walton G, Diederichs MS. A new model for the dilation of brittle rocks based on laboratory compression test data with separate treatment of dilatancy mobilization and decay. Geotech Geol Eng. 2015;33(3):661–679. 21 Schuerch R, Vrakas A, Anagnostou G. On manifestations of delayed failure and the effect of dilatancy in transient poro-elasto-plastic analyses of slopes and excavations. Geotechnique. 2017;67(11):939–952. 22 Cai W, Dou L, Ju Y, Cao W, Yuan S, Si G. A plastic strain-based damage model for heterogeneous coal using cohesion and dilation angle. Int J Rock Mech Min Sci. 2018; 110:151–160. 23 Guan K, Zhu W, Wei J, Liu X, Niu L, Wang X. A finite strain numerical procedure for a circular tunnel in strain-softening rock mass with large deformation. Int J Rock Mech Min Sci. 2018;112:266–280. 24 Kov� ari K. Tunnelling in squeezing rock. Tunnel. 1998;98(5):12–31. 25 Wu B, Gamage RP, Xi Z, Jeffrey RG, Mills K, Wang X. An analytical thermo-poroelasticity model for the mechanical responses of a wellbore and core during overcoring. Int J Rock Mech Min Sci. 2017;98:141–158. 26 Yan X, Ye YS. Research on test of reducing deformation of soft rock tunnel with high stress level by pilot heading. Appl Mech Mater. 2012;170–173:1565–1568. 27 Anagnostou G. Practical consequences of the time-dependency of ground behavior for tunneling. In: Rapid Excavation and Tunneling Conference: 2007 Proceedings. SME; 2007:255–265. 28 Vrakas A. Analysis of Ground Response and Ground-Support Interaction in Tunnelling Considering Large Deformations. ETH Zurich; 2016. 29 Vrakas A, Anagnostou G. Finite strain elastoplastic solutions for the undrained ground response curve in tunnelling. Int J Numer Anal Methods Geomech. 2015;39(7): 738–761. 30 Belytschko T, Liu WK, Moran B, Elkhodary KI. Nonlinear Finite Elements for Continua and Structures. 2014. 31 Vrakas A. Relationship between small and large strain solutions for general cavity expansion problems in elasto-plastic soils. Comput Geotech. 2016;76:147–153.

18