Analytical study on interaction between existing and new tunnels parallel excavated in semi-infinite viscoelastic ground

Computers and Geotechnics 120 (2020) 103385

Contents lists available at ScienceDirect

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

Research Paper

Analytical study on interaction between existing and new tunnels parallel excavated in semi-infinite viscoelastic ground

T

H.N. Wanga,b, , X. Gaoa, L. Wua, M.J. Jiangb,c ⁎

a

School of Aerospace Engineering and Applied Mechanics, Tongji University, Shanghai 200092, China State Key Laboratory of Disaster Reduction in Civil Engineering, Tongji University, Shanghai 200092, China c Department of Civil Engineering, Tianjin University, Tianjin 300072, China b

ARTICLE INFO

ABSTRACT

Keywords: Viscoelasticity Shallow tunnels Analytical solution Existing and new tunnels interaction

A new tunnel is usually excavated in close proximity to existing ones, which leads to significant negative influence on the existing tunnels. Moreover, the rheology of the ground probably induces quite different timedependent ground deformation if the new tunnel is excavated at different time. This study focuses on the interaction between the existing and the new tunnels which excavated in rheological ground at a shallow depth. Through the strict derivation, a new analytical solution is proposed for ground stress and displacement induced by the interaction of new and existing tunnels. The ground rheology, excavation delay of the new tunnel, tunnel size and various tunnel arrangement, are all taken into account. The complex variable theory combined with the extension of corresponding principle are employed in the derivation. By deriving the potentials in complex variable theory for the problems in all the excavation stages, the time-dependent stresses and displacements are finally addressed for the whole excavation process, where the ground is simulated by any linear viscoelastic models. To verify and validate the analytical solutions, the analytical solution is compared with numerical results under simplified and complex ground conditions, which shows good consistency except the solution for Case 1 (considering gravity gradient). A parametric study is finally preformed to find the influence of excavation time and location of the new tunnel, the tunnel spacing and relative size of the new tunnel, on stresses/displacements around tunnels and surface settlements. The results show that the excavation time of the new tunnel (t2) significantly influence the additional displacements around the existing tunnel which is a decrease exponential function of t2; when the distance from center to center is larger than 2.5R1 (2.0R1), the interaction between two tunnels can be neglected from perspective of displacement (stress).

1. Introduction New tunnel excavation in close proximity to the existing ones, as well as the sequential excavation of twin tunnels, are commonly encountered in underground constructions. The new tunnel excavation induces additional displacement and stress around existing tunnels, which probably results in instability of the existing tunnels [1–3]. Furthermore, the stress and displacement around the new tunnels are also quite different from those in single tunnel problem, due to the influence of the existing tunnel boundaries. Therefore, the interaction between existing and new tunnels should be exactly taken into account in assessment of tunnel stability. Rock or soil probably exhibits time-dependent behaviors (rheology), which induces gradual deformation over time even after the completion of the tunnel excavation [4,5]. Ground rheology induces quite different ⁎

ground responses when the new tunnel is excavated at different time. In urban areas, tunnels are generally constructed at a shallow depth, and greatly influenced due to the presence of ground surface. and This study focuses on the time-dependent ground responses induced by the interaction between existing and new tunnels at shallow depth. Numerical methods are usually employed in detailed design to predict the ground responses induced by a new tunnel excavation near the existing ones, with consideration on complex geological conditions and sequential excavation [6–10]. Numerical method provides very helpful results but requires longer running-times, especially when complete parametric analyses are performed. Empirical formulas with simple expressions are also widely used in engineering practice [11–13]. For example, superposition technique on the Gaussian equation suggested by Peck [11] is one of the most popular empirical formulas based upon field observations to predict the soft ground surface

Corresponding author at: School of Aerospace Engineering and Applied Mechanics, Tongji University, Shanghai 200092, China. E-mail address: [email protected] (H.N. Wang).

https://doi.org/10.1016/j.compgeo.2019.103385 Received 29 August 2019; Received in revised form 5 November 2019; Accepted 8 December 2019 0266-352X/ © 2019 Elsevier Ltd. All rights reserved.

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

movement induced by shallow twin tunneling [12,13]. However, the empirical equations fail to consider the time effect induced by the rheology properties of ground and the sequential excavation of two tunnels. Also the empirical formulas can only predict the surface settlement instead of the displacement and stress fields of the whole ground. In contrast, the simplified analytical models can be employed to efficiently obtain initial estimation of the design parameters in preliminary design. The analytical models provide the stress and displacement field of the whole ground through strict mathematical derivation, adopting all of the concerned parameters. They also provide an approach to gain deep insight into the mechanical mechanisms of engineering problems. Because of the difficulty in analytical derivation, the ground is assumed to be elastic in most of analytical research for two-tunnel (hole) problem in infinite or semi-infinite ground. Aiming at the stress around two holes in infinite elastic plane, the solutions were proposed by introducing the Airy stress functions in bipolar coordinates in many studies [14–18]. Based on the bipolar coordinates and Schwartz alternating method, the elastic analytical solutions were presented by Kooi and Verruijt [16] for displacement and stress caused by the interaction of twin tunnels in deep depths. Hoang and Abousleiman [17] presented the stress solution for an infinite plate containing two equal/unequal circular holes, subjected to general stresses at infinity and internal pressure along hole boundaries. For the holes in half-space, Spencer and Sinclair [19] employed a sequence of Airy stress function to derive an analytical solution for twin circular holes under gravitational load. These approaches are difficult in providing the displacement fields. The complex variable theory [20–25] is powerful in solving the elastic cases of non-circular openings or the cases with multiply connected regions. Based on complex variable theory, Manh et al. [20] presented the solutions for the stresses and displacements around two deep circular tunnels in an elastic ground. For non-circular holes, based on the Schwarz alternating method, the stress concentration factors were obtained by Ukadgaonker and Patil [22] for two unequal elliptical holes with different sizes, subjected to internal pressure and shear stresses along hole boundaries. Zhang et al. [23,24] provided the accurate stress solution for two/multiple elliptical holes in infinite region by using Schwarz alternating method and complex variable method. In recent analytical research [26], the stresses and displacements were provided around two closely located circular tunnels in deep ground, considering the rock viscoelasticity and the construction sequence. However, the solution was only applied to deep tunnel analysis. For twin tunnels in an elastic half space, complex variable method combined with Schwarz alternating method are usually employed: Fu et al. [25] presented the solutions for stress and displacement around shallow buried twin parallel tunnels, taking into account the prescribed uniform radial displacement along tunnel boundaries; Wang et al. [27] proposed the elastic solutions for shallow twin tunnels subjected to arbitrary distributed surcharge loads. In summary, the aforementioned analytical studies fail to consider the ground rheology and the sequential excavation of tunnels, or the tunnels were set in infinite viscoelastic ground. In this study, the analytical solution for ground displacement and stress after new tunnel excavation will be strictly derived, considering the ground viscoelasticity, the delay excavation of new tunnels and tunnels interaction. We expect that the solutions can help reveal the mechanical mechanism and the law of the interaction between the new and existing tunnels.

Fig. 1. Geometry of the problem, and reference coordinate systems.

medium/stiff clay, is in a state of small deformation. The ground rheology is accounted for by various linear viscoelastic models [28]. Assumed that the cross-section excavation is instantaneous. To tackle the problem as two dimensional plane-strain one, two advancement parameters 1 (t ) and 2 (t ) (0 < 1 ( 2 ) 1), are introduced to account for the progressive release of the stresses cause by the existing and second tunnel longitudinal advancement, respectively. i (t ) , i = 1 and 2, is a function of the variation of the radial displacement with the distance z between the tunnel face and the cross-section considered, that can be determined by in-situ measurements, or calculated numerically [29–31]. i (t ) equal to 1.0 in case of the distance z being sufficient large. According to the aforementioned assumptions, the equivalent planestrain problem in the plane of tunnel cross-section of the tunnels can be simplified, as illustrated in Fig. 1. Both the Cartesian coordinates (x1, y1) , (x2 , y2 ) and polar coordinates (r1, 1) , (r2, 2) are employed in the following derivation of analytical solutions. A tension-positive notation is used throughout this study. According to the start excavation time of two tunnels, the excavation process is divided into two stages in the derivation presented herein: in the first tunnel excavation stage spanning from t = t1 to t = t 2 , only the first tunnel is excavated, as shown in Fig. 2(a). In the second tunnel excavation stage spanning from t = t 2 onwards, both the two tunnels are excavated, as shown in Fig. 2(c). Note that t2 may be either much larger than or equal to t1 so that the sequential excavation process is accounted for including the excavation of both tunnels occurring at the same time. In the first tunnel excavation stage (Fig. 2a), x(2 1) and y(2 1) are the tractions acting on the future boundary of the second tunnel; whereas no tractions are exerted on the internal tunnel boundaries in the second tunnel excavation stage (Fig. 2c). Therefore, according to the superposition principle (viscoelasticity), the incremental mechanical responses after the second tunnel excavation can be addressed by the model in Fig. 2(b). 3. Derivation of analytical solutions 3.1. Formulation for the general viscoelastic problem The different stress-strain behavior of linear viscoelastic model can be schematized by a number of springs and dashpots connected either in series or in parallel. The methodology is provided in the literature [5] for solving a general viscoelastic problem involving time-dependent boundaries. In a sequential excavation, the state of stress is dependent of the excavation steps carried out over time, and the displacements in viscoelastic medium depend on the entire previous stress history. For a general problem, l loads are assumed to exert on the structure at different time before the generic time t , i.e., the k–th load (k = 1,2,…., l) is applied on the structure at time tbk and removed at tmk . According to the superposition principle for viscoelastic case along with complex variable theory [32], the total displacement at the generic time t can be calculated as follows:

2. Problem definition Prior to tunnel excavation, the semi-infinite ground without any tunnels is subjected to gravity induced initial stresses. Two parallel circular tunnels, i.e., the first (existing) tunnel and the second (new) one are then successively begin excavated at time t = t1 = 0 and t = t2 , respectively. In this study, the ground which composes of rock or 2

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 2. Schematic representation of the stress boundary conditions for calculation of (a) stresses and displacements before excavation of second tunnel; (b) stresses and displacements induced by excavation of second tunnel; (c) total stresses and displacements. The superposition relation is also shown.

u1vx (z , t ) + i u1vy (z , t )

uxv (z , t ) + iu yv (z , t ) =

1 2

l

Tk tbk

k=1

l k=1

Tk tbk

I (t

H (t

) ) z

(k ) 1

1 2

(z , ) d

(k )¯ 1 (z ,

z

)

= (k ) 1

+

I (t ) = L

1

1

0

(1)

v 1x (z , v 1y ( z ,

1 , sL [G (t )]

v xy (z ,

where

2

(k ) (z , t )

(k )

l

z¯

k=1

(z , t ) =

z¯

z

k=1

t ) = Im

) H (t

(k ) 1

0

2 (k ) (z , t ) z2

2 (k ) (z , t ) z2

(z , t )

+

+

(k ) (z , t )

z

t

[tbk , tmk ]

t

[tbk, tmk ]

(k ) (z , t )

z

.

1 z 2

)d

( ¯1)

(z )

z

¯ ( 1)

+

(z )

)d

(4)

p0x t) = p + Re 2 1 (t ) t) 0y

( 1) (z )

1 (t )

z

t ) = Im z 1 (t )

2 ( 1) (z ) z2

+

z

2 ( 1) (z ) z2

( 1) (z )

+

z

( 1) (z )

1 (t )

z

3.3. Solutions for the second tunnel excavation stage In order to obtain the time-dependent ground responses, the variation of boundary conditions during the whole excavation stage should be clarified. A model of two holes in half-plane is assumed with their boundaries loaded by various tractions during different excavation stages: (a) during first excavation stage, the ground is subjected to gravity, and the first tunnel boundary is stress-free. The future boundary of the second tunnel is subject to surface forces x(2 1) and (2 1) (see Fig. 2a) determined by the stress field in Eq. (5); (b) during y second tunnel excavation stage, the ground is subject to gravity and the two hole boundaries are both stress-free (see Fig. 2c). If considering the progressive stress release, this can be equivalent to exerting the trac(2 1) (2 1) tions, and on second tunnel boundary (see 2 (t ) x 2 (t ) y Fig. 2b) from t2 to generic time t on the bases of the boundary conditions in Fig. 2(a). According to Eq. (1), the additional displacements at time t (t [t2, ) ) occurred after the second tunnel excavation can be expressed as follows:

,

(3) and

) I (t

where p0x and p0x are horizontal and vertical initial stresses induced by gravity, respectively.

(k )

(k ) 1

(z , t ) t [tbk , tmk ] ; Re[·] and Im[·] denote the real and 0 t [tbk , tmk] imaginary component of a generic complex variable [·]. (z , t ) =

0

1(

(5)

(plane strain problem),

where L [f (t )]) with respect to variable s is defined in the Laplace transform of the function f (t ) , andL 1 [g (s )] indicates the inverse Laplace transformation of g(s); G(t) and K(t) represent the shear and bulk relaxation moduli of viscoelastic model, respectively. The expressions of H(t) and I(t) for generalized Kelvin and Poynting-Thomson viscoelastic models can be found in Table 3 [33]. The stresses are addressed by exploiting the principle of superposition as: l

t

v 1xy (z ,

3L [K (t )] + 7L [G (t )] 1 sL [G (t )] 3L [K (t )] + L [G (t )]

t) = Re t)

1(

(z )

According to Eq. (3) and considering the initial stress fields, the total stresses can be expressed as follows:

(2)

v x (z , v y (z ,

( 1)

t

¯(z, ) d

where Tk = min{tmk, t } . 1(k ) and 1(k ) are the potentials corresponding to the case subjected to k–th load in the loading period [tbk, tmk ]; z = x + iy 1 ; e (z¯, t ) are the conjugates of the complex function and i = e = e (z, t ) , and

H (t ) = L

1 2

u2vx (z , t ) + i u2vy (z, t )

3.2. Analytical solutions for the first tunnel excavation stage

v (1) v (2) v (2) v = u2(1) x (z , t ) + i u 2y (z , t ) + [u 2x (z , t ) + iu 2y (z , t )]

The additional stresses and displacements caused by excavation are crucial in engineering application. Prior to tunnel excavation, the gravity-induced initial stresses acted on the future boundaries of the tunnels, e.g., x(1 0) and y(1 0) on the boundary of the first tunnel. The additional ground responses induced by the first tunnel excavation can be determined by the model of the first tunnel subjected to tractions (1 0) (1 0) and along tunnel boundary. Based on the 1 (t ) x 1 (t ) y previous research [27], the corresponding elastic complex potentials of this problem in case of 1 = 1, ( 1) ( ( 1) ), can be addressed (detailed in Appendix A). According to Eq. (1), the additional displacements at generic time t (t [t1 = 0, t2]) of viscoelastic problem can be expressed as follows:

(6)

where the displacements are divided into two parts, i.e., u.(2) v , whose expressions are the follows: 1

v (1) v Part 1: u2(1) x + i u 2y = 2 [

1 t [ ( 2 0 1 t2 0

) H (t 1(

t ( ) I (t 0 1 t2 1 ( ) I (t2 0

u.(1) v

and

)d )d ]

( 1)

(z )

)d

) H (t2

)d ] z

( ¯1) (z )

z

+

¯ ( 1)

(z ) (7)

3

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al. v (2) v Part 2: u2(2) x (z , t ) + iu 2y (z , t ) = 1 2

t t2

1(

1 2

( 2)

) 2 ( ) H (t

(z )

t t2

)d

z

1(

) 2 ( ) I (t

( ¯2) (z )

z

+

)d

¯ ( 2)

(z ) (8)

and in part 2 denote the potentials for the twin tunnel model in Fig. 2(b) in cases of 1 = 2 = 1.0 . Determination of these potentials is detailed in Appendix B, where Schwartz alternating method is introduced to translate the two tunnel problem into a series of single tunnel problems. Noted that the part 2 displacements are induced by excavation of the second tunnel, and the part 1 displacements are the additional rheology displacements in this stage induced by the first tunnel excavation. Therefore, the displacements occurred after excavation of the second tunnel are significantly dependent of the whole sequential excavation process of the two tunnels. However for the elastic problem, the part 1 displacements in Eq. (7) are zero. The additional stresses occurred during this stage are: ( 2)

v 2x (z , v 2y (z ,

( 2)

t) = Re 2 1 (t ) 2 (t ) t) v 2xy (z ,

( 2) (z )

1 (t ) 2 ( t )

z

t ) = Im z 1 (t ) 2 (t )

2 ( 2) (z )

z2

+

z

2 ( 2) (z )

z2

1 (t ) 2 (t )

+

( 2) (z )

z

( 2) (z )

z

(9) Because the effects of longitudinal advancement are investigated only in case of single tunnels [30,31], the expressions for the advancement parameters, 1 and 2 , are not available in the literature. Therefore, in the following of the paper 1 = 2 = 1.0 is assumed. This means that we consider only cross-sections located at a distance from both tunnel faces such that three dimensional effects are not felt. 4. Verification and validation of analytical solutions Fig. 3. Geometry of the domain together with the boundary conditions in numerical simulation for verification of derivation of (a) Case 1 and (b) Case 2.

4.1. Verification of analytical solutions To verify the analytical solution, the derived solutions are compared with the results from numerical simulations carried out using the software ANSYS (Version 14.0, employing the “structure mechanics” module). To maintain the consistency with the analytical derivation, all the FEM analyses are performed under plane-stain conditions and the same assumptions in the analytical derivation. A semi-infinite medium with two circular holes are considered. In the simulation, the generalized Kelvin model with shear modulus G M = 20MPa (for Maxwell part), G K = 10MPa (for Kelvin part), and coefficient of viscosity K = 100MPa· day , is assumed for the ground, and the Poisson’s ratio v is 0.25. Meanwhile, the unit weight and lateral pressure coefficient K 0 are 18.5kN m3 and 0.27, respectively. The effect of longitudinal advancement is not considered, i.e., 1 = 1 = 1.0 . According to the conditions in practical engineering, the comparison for the following two cases is performed:

In the numerical simulation, the gravity induced the initial stresses are first generated in the half-plane without any holes. Afterwards, the elements inside the left hole are removed to simulate first tunnel excavation at t = 0, where the stiffness matrix of these elements was multiplied by a coefficient of 10 6 (to deactivate elements for the death capability in ANSYS). Consequently, the elements inside the right hole are removed after 30 days (t2 = 30 day) to simulate the excavation of second tunnel. The excavation induced incremental displacements and stresses can be obtained by subtracting the displacements and stresses that occurred before excavation from the total ones. The numerical model has been verified in our previous studies [26,35,36]. Figs. 5 and 6 show the comparison on additional displacements after the first tunnel excavation, as well as the total stresses versus time at points P1- P4 (point locations are shown in Fig. 4), predicted by the analytical solutions and FEM simulations. Noted that the results from analytical and numerical models exhibit a close agreement for both Cases 1 and 2. Thus, the analytical and numerical results are mutually verified.

Case 1:. the variation of vertical gravity loading across the height of excavation is considered in the initial stress state; Case 2:. the variation of vertical gravity loading across the height of excavation is neglected. The initial stresses along anticipated tunnel boundary are assumed approximately equal to the stresses at the center of the tunnel [34].

4.2. Validation of analytical solutions

The geometry and boundary conditions applied in FEM simulation for two cases are shown in Fig. 3: the vertical displacements along the bottom boundary and the horizontal displacements along left and right boundaries are restrained; the two tunnels have the same size and buried depth in verification of Case 1, whereas in Case 2 verification, the size and buried depth of two tunnels are quite different. Fig. 4 presents the mesh in the vicinity of the tunnels for two cases, where the mesh in the yellow area represents the excavated rock or soil.

In order to further validate the applicability of proposed analytical model, the comparisons on surface settlement and displacements along tunnel boundaries were performed between analytical solutions and results from finite differential method (FDM), where the more complex constitutive relations and geological condition of the ground, are adopted in the numerical model. Numerical analyses were performed by software FLAC3D for surface settlement in the Istanbul metro [37]. The numerical results were 4

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

balance out the resultant forces Y (see Eq. (A2) in Appendix A) induced by the variation of vertical gravity loading (Y is equal to the weight of excavated rock/soil, the direction is vertical upward), which makes the Case 2 model being more close to the actual situation than Case 1 model. Therefore, the analytical solutions of Case 2 will be employed in the following parametric investigation. By using software FLAC3D, the deformation along the first and the second tunnel boundary induced by the driving of new tunnel in the shallow depth was presented in the reference [8]. The geometry of the problem is shown in Fig. 9 and the parameters of mechanical properties are presented in Table 2. The comparisons are plotted in Fig. 10 on radial displacements induced by the second tunnel excavation between analytical solution of Case 2 and numerical results for various tunnel spacing. It is noted that the variation pattern of induced radial displacements for all the cases are very similar between analytical and numerical results, and the magnitude of analytical results is somewhat larger than numerical one due to the absent of lining in the analytical model. From the good agreement between analytical and numerical results for real engineering, we believe that the proposed analytical solution of Case 2 can provide reliable prediction of displacement field qualitatively. 5. Mechanical mechanism of time-dependency of displacement Different viscoelastic models can simulate different rheological characteristics of geo-material [38]. Generalized Kelvin and PoyntingThomson viscoelastic models have limited viscous deformation under constant stress, which are usually employed for simulating the rock with good mechanical properties or subjected to low stress. Different from the elastic problems where the additional displacements after the second tunnel excavation are completely comes from the excavation of the second tunnel, the additional displacements are both influenced by the excavation of the first and second tunnels when the rheology of the ground is considered (part 1 and part 2 displacements shown in Eqs. (7) and (8). According to the expressions of functions H(t) and I(t) provided in Table 3, the results of the four integrations about t in expressions of part 1 and part 2 displacements are presented in Table 4 for the two viscoelastic models ( 1 = 2 = 1.0 ). The functional forms of the integration expressions for the two models are similar. According to the expressions of J1 and J2, the influence of the first tunnel excavation (part 1) on the additional displacements exponential decreases against time, and is a constant C when t tends to infinity. The constant C is an exponential function of t2, i.e., exp( Ct2 ) , which shows the influence of excavation time of the second tunnel on additional displacement is a decrease function. The displacements induced by the second tunnel excavation (part 2) are the exponential increase functions of (t2-t). Noted the excavation time of the second tunnel has no influence on the final part 2 displacement.

Fig. 4. FEM mesh of the vicinity domain: (a) and (b) the mesh for Cases 1 and 2, respectively.

closely consistent with field measurements [37], which proved the applicability of the numerical model. Fig. 7 shows the geological geometrical conditions. The parallel twin tunnels with diameter of D = 6.16 m, both located at depth of 9.74 m. In the numerical simulation, the elasto-plastic analysis based on Mohr-Coulomb failure criterion for the ground was carried out with the shotcrete lining and soil nails modeled by different elements. The thickness and mechanical properties of soil layers are provided in Table 1 [37]. The equivalent material parameters employed in analytical solution are determined as follows:

Qe =

Q1 h1 + Q2 h2 + Q3 h3 + Q4 h4 h1 + h2 + h3 + h 4

(10)

where Qe represents the equivalent parameter (i.e. Young’s modulus, Poisson’s ratio and unit weight) in analytical solution; Qj and hj (j = 1, 2, 3, 4) represent the material parameters and thickness of soil layer j in Table 1, respectively. The effect of support, the construction process and ground plasticity are neglected in the analytical calculation. Based on these parameters, the final ground settlements obtained by the proposed analytical solutions (for Cases 1 and 2), as well as the corresponding FDM results, are plotted in Fig. 8 for various tunnel spacing. The variation patterns of surface settlement predicted by analytical model of Case 2 (without gravity variation) are closely consistent with the FDM results in tendency except for the larger values. This may be due to the simplification in the analytical model, e.g., the neglecting of lining support and ground plasticity which were considered in the numerical model. In contrast, the variation pattern and magnitude of settlements obtained by analytical model of Case 1 (with gravity variation) are both very different from FDM results. The reason probably is that: in the real engineering, the weight of liner and shotcrete-bolt support can partly

6. Parametric investigation For geo-materials with good mechanical properties and low stresses in application, limited viscosity is present thus the Kelvin, generalized Kelvin or Poynting-Thomson viscoelastic models are commonly employed. In this section, a parametric analysis is carried out considering the effects of the various parameters on the displacements and stresses around the tunnels. The ground rheology is described by the generalized Kelvin viscoelastic model with constitutive parameters GM GK = 2 . The displacement in this section is the incremental one occurred after the second tunnel excavation, and the stress is the total one. 1 = 2 = 1.0 in the following analysis (fast advancement). For the sake of generality, the displacements and stresses are normalized by ur and p0 , respectively, where p0 is the initial vertical stress at the depth of first tunnel center; ur is the radial displacement on a tunnel boundary for a circular tunnel excavated in an infinite elastic 5

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 5. Comparison on displacements between analytical and FEM results: (a) and (b) displacements for Case 1; (c) and (d) displacements for Case 2. P1 as shown in Fig. 5.

medium, subjected to hydrostatic far field stress p0 :

ur

p R1 = 0 2Gs

P4 are points

the first tunnel excavation, called non-simultaneous excavation; (c) t2 = 5.0TK , i.e., the second tunnel is excavated after a very long time of the first tunnel excavation, called stable state excavation. Meanwhile, the buried depths and radii of tunnels are H1 = H2 = 2R1 and R1 = R2 , respectively; the spacing between twin tunnels is c = R1. In the following analysis, the additional displacements occurred after excavation of the second tunnel, and the total stress of the ground, are provided by analytical model. The additional displacements at the given points versus time are plotted in Fig. 11 for stable state excavation case. The plots illustrate that the displacements increase with time and tend to be stable approximate at time t TK = 8 . Because the final stable displacements are the most significant in tunnel design, these displacements are adopted in the following parametric analysis. In Fig. 12 the final surface settlements occurred after the excavation of the second tunnel against x1 are presented for the aforementioned three cases. The figure shows that the settlements along the ground surface at the region over the first tunnel are quite different when the second tunnel excavated at different times, e.g., they obviously decrease

(14)

where R1 is radius of the first tunnel; Gs = GM GK (GM + GK ) is the permanent shear modulus. Other parameters indicated dimension are normalized by R1. The generic time t in viscoelastic cases is normalized by the retardation time TK (TK = K GK for the generalized Kelvin model). The sign convention adopted in this section is the same as that in Section 4. 6.1. Influence of the excavation time of the second tunnel In this sub-section, three cases with following excavation time of second tunnel are considered: (a) t2 = 0 , i.e., the second tunnel is excavated instantaneously after the first tunnel excavation (simultaneous excavation), thus only elastic deformation has been occurred before the second tunnel excavation; (b) t2 = TK , i.e., the second tunnel is excavated after a certain time of 6

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 6. Comparison on stresses between analytical and FEM results: (a) and (b) stresses for Case 1 and Case 2, respectively.

Fig. 7. Geological properties and geometry of twin tunnel problem in validation of analytical model: (a) geological profile of longitudinal cross section; (b) geometry of twin tunnels. Table 1 Thickness and properties of soil layers in Zone C of project area [37]. Layer

Layer Layer Layer Layer

Soil type

1 2 3 4

Artificial filling Sand Suleymaniye formation Trakya formation

Thickness (m)

2.95 10.6 11.32 Base

Unit weight (kN/m3)

18 17 18.9 25

Mohr-Coulomb shear strength parameters Cohesion (kPa)

Friction angle (deg)

1 1 20 80

10 25 14 25

7

Young’s modulus (MPa)

Poisson’s ratio

5 15 38 60

0.4 0.35 0.33 0.2

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 8. Comparison on surface settlement versus x1 between analytical solution for Case 1/Case 2, and FDM results [37]: (a) and (b) the settlements in case of tunnel spacing c = 5D and c = D , respectively.

that in Section 5, and the function form of displacement variation with time can be clearly found through the analytical solutions. The hoop stresses along first tunnel boundary versus 1 are plotted in Fig. 14 for twin tunnels and single tunnel problem. The interaction of the two tunnels causes the stresses around 1 = 0° in twin tunnels problem larger than that in single tunnel problem, and the maximum difference is up to 20% of the maximum stresses in single tunnel problem. Meanwhile, from Fig. 13(b) it is seen that the maximum difference of final displacements between twin tunnels and single tunnel problem is occurred around 2 = 135° and the difference is up to 22% of that in single tunnel problem.

Fig. 9. Geometry of the tunnels in FDM simulation [8].

6.2. Influence of location of the second tunnel

Table 2 Soil properties in the project area [8]. Unit weight (kN/m3)

Soil

18.5

Mohr-Coulomb shear strength parameters Cohesion (kPa)

Friction angle (deg)

0

37

Young’s modulus (MPa)

Poisson’s ratio

150

0.3

The new tunnel can be constructed at different position near the existing tunnel in practical engineering, which induces the quite different ground responses. According to the relative position of the two tunnels, the arrangement angle is defined in Fig. 1. In this sub-section, the first tunnel is fixed at a depth H1 = 2R1, and the following five cases with different location of the second tunnel are considered: (1) = 15° ; (2) = 0° ; (3) = 30° ; (4) and (5) = 60° = 90° . It is assumed that the size of the two tunnels are equal. The tunnel spacing is c = R1 and the second tunnel is excavated after a very long time of the first tunnel excavation. The final additional displacements along tunnel boundaries versus 1 or 2 are plotted in Fig. 15 for the five cases, where the variation and magnitude of displacements observed are both quite different, especially for those along the first tunnel boundary. Fig. 15(a) and (b) indicate that the displacements along the first tunnel boundary are a great fluctuation, and are notable larger in the cases with = 60° and 90° . The radial displacements show a peak value at 1 = 330° 300° 270°when the arrangement angle is 30° 60° 90° , respectively. It is emerged from Fig. 15(c) and (d) that the variations of displacements along the boundary of second tunnel are similar in the five cases, except for the magnitudes. The maximum radial displacement along second tunnel boundary occur at 2 = 90° in the case with = 15° , due to the closer distance from second tunnel center to ground surface. The maximum

with the delay time of second tunnel excavation, whereas the settlements are almost the same in the three cases at the region over the second tunnel. Fig. 13 presents the final additional displacements along tunnel boundaries versus 1 or 2 for various excavation times of the second tunnel. The differences of displacement among these cases are obvious along the first tunnel boundary, however little difference is observed of displacements along the second tunnel boundary. Therefore, it is concluded that the excavation time of the second tunnel only significantly influences the additional displacements around the first tunnel. The earlier the second tunnel excavated, the larger the displacements around the first tunnel are. The final displacements around the second tunnel are almost unaffected by the excavation time of the second tunnel. It can be noted that the above conclusions are consistent with

8

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 10. Comparison on deformations of tunnel boundaries between analytical and FDM results [8]: (a) and (b) radial displacements along the first tunnel boundary obtained by analytical and numerical results, respectively; (c) and (d) radial displacements along the second tunnel boundary obtained by analytical and numerical results, respectively.

9

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Table 3 Expressions of H(t) and I(t) for the two viscoelastic models (modified from [33]). H(t) Generalized Kelvin model

I(t)

1 GM

1

(t ) +

K

GK t K

exp

H (t ) +

(

exp Poynting-Thomson model

1 P

(

) exp 2

GP GP + GH

1 GP + GH

+

(

GP GH t (GP + GH) P

(t )

6 3K e + GM

)

(t ) +

3K e GK + GM (3K e + GK ) t K (3K e + GM )

6 3K e + GP + GH

+

1 GP + GH

GP (3K e + GH ) t P (3K e + GP + GH)

exp

2 6GM 2· K (3K e + GM )

+

1

) (t) + ( ) exp

P

6GP2 2· P (3K e + GP + GH) 2 GP GP GH t GP + GH (GP + GH) P

Table 4 Integration results for the expressions of displacements for Generalized Kelvin and Poynting-Thomson viscoelastic models. Displacement occurred after the second tunnel excavation = Part 1 + Part 2 displacements Part 1 displacements, induced by the first tunnel excavation ( ¯1) (z ) 1 1 v (1) v ¯ (z ) , with J (t ) = ( 1) u 2(1) (z ) J2 (t ) z + ( 1) 1 x + i u 2y = J1 (t ) 2

2

z

t 0

I (t

Part 2 displacements: induced by the second tunnel excavation ( ¯2) (z ) 1 1 ¯ (2) v v ( 2) u2(2) (z ) J4 (t ) z + ( 2) (z ) , with J3 (t ) = x (z , t ) + iu 2y (z , t ) = J3 (t ) 2

Integration J1 (t )

z

t t2

I (t

I (t2

) d ; J2 (t ) =

) d , J4 (t ) =

Generalized Kelvin model BK CK

J2 (t )

J2 (t )

1 GK

J3 (t )

J4 (t ) Generalized Kelvin: AK =

2

t2 0

)d

GK

exp

K

t

1 GK

1

exp

GK

1 GK

1

exp

GK

6 , 3K e + GM

Poynting-Thomson: AP1 =

exp( CK t )[1

BK =

1 , GP + GH

K

K

exp

t)

+

(t

t2)

2 6GM K (3K e + GM )

AP2 =

2,

H (t

)d

t2 0

H (t2

)d .

)d .

H (t

Poynting-Thomson model

exp( CK (t2

1

t t2

t 0

GK (t2 K

1 GM

+

CK =

6 , 3K e + GP + GH

+

t ))]

J2 (t ) 1 BP 1 CP

t)

BK [1 CK

exp( CK (t

1 BP 1 [1 CP 1 BP 1 [1 CP

t2))] + AK

1 GM

BP2 1 CP

exp( CP2 t )[1 exp( CP2 (t2

exp( CP1 t )[1 exp( CP1 (t2 exp( CP1 t )] + AP1 +

exp( CP1 (t

BP2 1 [1 CP

t ))]

t )] exp( CP2 (t

t2))] + AP2

t2))] + AP1

3K e GK + GM (3K e + GK ) ; K (3K e + GM )

BP1 =

1 P

(

GP GP + GH

Fig. 11. Normalized displacements versus time at points with horizontal and vertical displacements, respectively.

1

) ,B 2

2 P

= 0° 90° ,

6GP2

=

P (3K e + GP + GH)

2

2

, CP1 =

GP GH , (GP + GH) P

CP2 =

GP (3K e + GH ) P (3K e + GP + GH)

.

= 180° along tunnels boundaries and O1 for stable state excavation: (a) and (b) the

10

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 14. Normalized hoop stresses along the first tunnel boundary versus single tunnel and twin tunnels problem.

Fig. 12. Surface settlements versus x1 for the simultaneous, non-simultaneous and stable state excavations of the second tunnel.

1

for

tunnel boundary, especially when = 60° 90°. The maximum compressive stress occurs in the case with = 15° at the point with 1 = 0° , while the maximum tensile stress occurs in the case with = 60° at the point with 1 = 260° . Fig. 16(b) shows that the maximum stress concentration is greater along the second tunnel boundary than that along the first tunnel boundary. The maximum compressive stress along the second tunnel boundary occurs in the case with = 30° and 60° , while

circumferential displacement occur at 2 = 40° in the case with = 60°. Overall, the displacements along two tunnel boundaries are both small when the tunnels are horizontally arranged. Fig. 16 presents the hoop stresses versus 1 or 2 along tunnel boundaries for the five cases. It is seen from Fig. 16(a) that the variations of hoop stress in the five cases are much different along the first

Fig. 13. Normalized displacements along tunnels boundaries versus 1 or 2 for simultaneous, non-simultaneous and stable excavation of the second tunnel: (a) and (b) displacements along the first and second tunnels boundaries, respectively.

11

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 15. Normalized displacements versus 1or 2 for various locations of the second tunnel in stable state excavation: (a) and (b) the radial and circumferential displacements along the first tunnel boundary, respectively; (c) and (d) the radial and circumferential displacements along the second tunnel boundary, respectively.

the maximum tensile stress occurs in the case with = 15° . The figures also indicate that the horizontal arrangement of two tunnels leads to relative smaller hoop stresses. The ground surface settlements in the cases with different

arrangement angles are plotted in Fig. 17. Because of the shallower buried depth of the second tunnel with = 15° , the settlement in this case is obviously larger than those in other cases, while the settlement is much smaller when the second tunnel is located blow the first tunnel

12

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 16. Normalized hoop stresses in stable state excavation versus second tunnel boundaries, respectively.

1or

2

for various locations of the second tunnel: (a) and (b) the hoop stresses at the first and

over the tunnel center when the two tunnels are vertically arranged ( = 90° ), instead of a shape of “V” with flat bottom in the other cases (it shows as a flat bottom because of the small tunnel spacing). It is probably because the presence of first tunnel significantly reduce the vertical ground settlement induced by the second tunnel excavation. 6.3. Influence of tunnel spacing and the relative size of the second tunnel In this sub-section, the buried depth of tunnels is assumedH1 = H2 = 3R1, and the stable state excavation of second tunnel is adopted. Fig. 18 presents the final displacements at the given points versus the tunnel spacing. It is noted that the absolute value of displacements all decrease gradually with the increase of tunnel spacing, and a similar decreasing rate can be observed. The horizontal displacements at points with 1 = 90° and 2 = 90° tend to zero when c R1 larger than 2.5, which indicates that the interaction between tunnels can be neglected. Fig. 19 presents the hoop stresses at the given points versus the tunnel spacing. Compared with the variation of displacement in Fig. 18, the stresses also decrease with the increase of tunnel spacing, however, the decreasing rate is much larger than that of displacement. The stresses are almost unchanged with spacing when the spacing is larger than 2.0R1. In order to investigate the influence of relative size of second tunnel, the displacements versus the relative radius of second tunnel are plotted in Fig. 20. The displacements nonlinearly increase with the increase of radius of second tunnel. The greatest increasing rate of horizontal displacement appears at first tunnel boundary with 1 = 0° , while that of vertical displacement appears at second tunnel boundary with 2 = 90°.

Fig. 17. Final surface settlement in stable excavation for various locations of the second tunnel.

( = 90° ). It can be also observed that the position corresponding to maximum settlement gradually moves to left with the arrangement angle change from 15° to 90° . It is interesting to see that the variation of settlement against x1 is taken on a shape of “W” with a slight heave

13

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 18. Normalized final displacements versus tunnel spacing at points with horizontal and vertical displacements, respectively.

1

= 0° 90° ,

2

= 180° 90° along tunnels boundaries, and points O1 and O2: (a) and (b) the

unchanged. The increase of hoop stresses at point with most obvious.

2

= 90° is the

7. Conclusions In this study, the analytical approaches were provided for obtaining the time-dependent stress and displacement induced by excavation of new tunnels in close proximity to the existing ones in semi-infinite viscoelastic ground, by complex variable theory and the methodology for solving a generalized viscoelastic problem. Two types of initial stress field due to the gravity were introduced. A close agreement was shown between the analytical solutions and numerical results under the same simplification. The applicability of analytical solution in real engineering was validated through the comparison between analytical and FDM results under complicated geometry and ground conditions. Finally, the parametric investigation was performed to investigate the influence of some core parameters. From the investigation, the following conclusions are drawn:

Fig. 19. Normalized hoop stresses versus tunnel spacing at points with ° ° ° ° 1 = 0 90 and 2 = 180 90 along tunnels boundaries.

(a) The variation patterns of surface settlement predicted by analytical model of Case 2 (without gravity variation across the height of excavation) are close consistent with the FDM results for real engineering. However, the variation pattern and magnitude of settlements obtained by analytical model of Case 1 (with considering the gravity gradient) are both very different from FDM results. (b) The additional displacements around the first tunnel, which is an exponential function of t2 (exp( Ct2 )), is significantly influenced by the excavation time of the second tunnel. The final displacements around the second tunnel are almost unaffected by the excavation time of the second tunnel. (c) The displacements and hoop stresses along two tunnel boundaries are both relatively small when the tunnels are horizontally arranged. When the tunnels are vertically arranged, the presence of the first tunnel significantly reduce the vertical ground settlement

It is also emerged from Fig. 20(b) that the vertical displacements at ° 1 = 90 and O1 are almost unchanged with the increase of the second tunnel radius when R2 R1 is larger than 1.5. The transition from negative to positive value for horizontal displacements at point with ° 2 = 180 is observed. The horizontal displacements along first tunnel boundary and vertical displacements along the second tunnel boundary both increase obviously with the increase of the second tunnel radius. Fig. 21 presents the hoop stress at the given points versus the relative second tunnel radius. It is noted that the absolute value of hoop stresses increase with the increase of second tunnel radius at the second tunnel boundary with 2 = 180° 90° , whereas it decreases at the first tunnel boundary ( 1 = 90° ). The stresses at point with 1 = 0° are almost

14

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

Fig. 20. Normalized displacements versus the relative radius of second tunnel, at points with (a) and (b) the horizontal and vertical displacements, respectively.

1

= 0° 90° ,

2

= 180° 90° along tunnels boundaries, and points O1 and O2 :

The analytical solution proposed in this study can be used as an alternative approach for the preliminary design of future shallow tunnels. However, the solution is under the simplifying assumption of small deformation and unsupported condition, due to the difficulty in mathematical and mechanical modelling if complex conditions are considered. In the future study, the excavation of multi tunnels and the effect of support deserve attraction. CRediT authorship contribution statement H.N. Wang: Conceptualization, Methodology. X. Gao: Methodology, Writing - review & editing. L. Wu: Software, Validation, Writing - original draft. M.J. Jiang: Supervision, Project administration. Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 21. Normalized hoop stresses versus the relative radius of the second tunnel, at pints with 1 = 0° 90° and 2 = 180° 90° along tunnels boundaries.

Acknowledgements

induced by the second tunnel excavation. (d) When the distance from center to center is larger than 2.5R1 (2.0R1), the interaction between two tunnels can be neglected from perspective of displacement (stress).

This study was supported by the National Natural Science Foundation of China (Grant Nos. 11572228, 11872281, 51639008); the State Key Laboratory of Disaster Reduction in Civil Engineering (SLDRCE19-A-06). These supports were greatly appreciated.

Appendix A: Determination of elastic potentials for the first tunnel excavation stage [27,39] In the first tunnel excavation stage, only one circular tunnel (the first tunnel) contains in the ground. The region of a half-plane containing a circular hole can be mapped conformably onto an annulus in the ξ-plane, with the outer and inner radii of 1 and , respectively [34]. The mapping equation can be expressed as follows:

15

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

z=

( )=

where

=

ih h R0

1 1+

1

2 2

1+ 1

(A1)

R02 h2

1

; and h is the depth of the tunnel center,

=

+ i = ei . The following presents the outline of the derivation of

(1 0) (1 0) potentials for single tunnel problem, subject to tractions and along tunnel boundary. x y The two potentials in complex variable theory for generalized shallow tunnel problem can be expressed as follows [39]: ( 1)

(z ) =

( 1)

(z ) =

X + iY [ 2 (1 + )

ln(z

z¯0 ) + ln(z

z 0 )] +

k

ak

+

k=0

X iY [ln(z 2 (1 + )

z¯0 ) + ln(z

z 0 )] +

k

bk k=1

ck

k

+

k=0

dk

k

(A2)

k=1

1

2

where (X , Y ) represents the resultant vector of the tractions acting on the tunnel boundary; z 0 = ih 1 + 2 ; = 3 4 for plane strain problem with being the Poisson's ratio; the coefficients ak , bk , ck , and dk are determined according to the stress boundary conditions. According to the real engineering practice, two cases with different initial stress states are considered: Case 1: the variation of vertical gravity loading across the height of excavation is considered in the initial stress state; Case 2: the variation of vertical gravity loading across the height of excavation is neglected. The initial stresses along anticipated tunnel boundary are assumed approximately equal to the stresses at the center of the tunnel [34]. In Eq. (A2), X = 0 and Y = R2 in Case 1, whereas X = Y = 0 in Case 2. The boundary conditions, in terms of the potentials for the ground surface (zero stress), are as follows: ( 1)

( , t) +

( ) ¯( )

( 1)¯ (

, t)

+

( 1)¯(

, t)

=0

(A3)

= ei

The stress condition on the tunnel boundary is as follows:

(1

ei )

( 1)

where C3 is a constant, dy

( ) ¯( )

( , t) +

dx

b x

and

( 1)¯ (

, t)

+

( 1)¯(

, t)

= (1 = e

b y

ei ) i

(

b x

+i

b y )ds

+ C3 ;

(A4)

i

are projections along the x and y of boundary stresses, respectively.

b x

=

K 0 y cos ,

b y

= y sin

in Case 1, and

= K 0 h ds , yb = h ds in Case 2, where K 0 is the lateral pressure coefficient and is the angle between the outward normal direction of the tunnel boundary and the x-axis. According to the boundary conditions, the linear system of the equations with regard to the coefficients in the potentials can be finally determined [39]. b x

Appendix B: Determination of elastic potentials for twin tunnel model [27] Fig. A1 presents the alternating iterative processes for solving the problem of semi-infinite ground containing the two tunnels loaded by tractions (2 1) and along the second tunnel boundary (Fig. A1a) ( 1 = 2 = 1.0 ). Analysis for single-tunnel problems subject to redundant surface y tractions is repeated in the alternating iterations to eliminate the non-zero tractions on the tunnel boundaries. (2 1) (2 1) In the first iteration, as illustrated in Fig. A1(b), only the second tunnel contained in the half plane is subject to tractions and x y along the second tunnel boundary. By the derivation in Appendix A, the potentials 11 ( ) and 11 ( ) can be addressed. Subsequently, the non-zero redundant tractions S 1 1 (Sx1 1 and Sy1 1 are the horizontal and vertical components), induced along the first tunnel boundary will be determined by (2 1) x

Fig. A1. Schwartz alternating iterative steps for solving the problem of semi-infinite ground containing the two tunnels loaded by tractions the second tunnel boundary ( 1 = 2 = 1.0 ). 16

(2 1) x

and

(2 1) y

along

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al.

the following equation. zb1 zb01

f1 (z b1) = i

[Sx1

1

(z b1) + iSy1

1

(z b1)] ds =

1 1 (z b1)

+

1 (z b1)

¯

1 (z b1)

d 11¯ ( ) d

+ = z b1

1¯ 1 (z b1)

(A5)

where and z b1 represent the certain and generic points on boundary of the first tunnel in mapping planes; 1 (z ) is the mapping function of the problem only the first tunnel containing in the ground. In order to eliminate the redundant tractions to satisfy the stress boundary condition in Fig. A1(a), the ground only containing the first tunnel is then investigated (see Fig. A1c), where the reverse redundant tractions S 1 1 (z b1) are exerted on tunnel boundary. Therefore, the boundary condition yields as follows: 0 z b1

1 2 (z b1)

1 (z b1)

+

¯ 1 (z b1)

d 21¯ ( ) d

+

1¯ 2 (z b1)

=

f1 (z b1)

(A6)

= zb1

In order to easily obtain the analytical solutions through Cauchy’s Integral Formulas in the next iteration, the functions f1 (z b1) is expanded as Fourier series. The potentials 21 ( j ) and 21 ( j ) are addressed via the procedure in Appendix A. Through the follow-up repeating iterations, the redundant tractions along tunnel boundaries will approach zero with the increase of the iterations. Assuming that ( 2) and ( 2) are the potentials of the twin-tunnel problem in Fig. 2(b), they are finally derived by superposing the potentials in all iteration steps, as the following: ( 2)

( 2)

n

( )= ( )=

k=1 n k=1

[

k 1 (

)+

k 2 (

)],

[

k 1 (

)+

k 2 (

)].

where n is the final iteration number; th iteration.

(A7) k j (

) and

k j (

) denote the potentials in the case of half-plane medium only containing the j-th tunnel in the k-

[17] Hoang SK, Abousleiman YN. Extended Green's solution for the stresses in an infinite plate with two equal or unequal circular holes. J Appl Mech 2008;75(3). [18] Radi E. Path-independent integrals around two circular holes in an infinite plate under biaxial loading conditions. Int J Eng Sci 2011;49(9):893–914. [19] Spencer SE, Sinclair GB. Stress analysis of an elastic half-space perforated by twin circular holes and under gravitational loading. Int J Mech Sci 1982;24(1):27–35. [20] Manh HT, Sulem J, Subrin D. Interaction of circular tunnels in anisotropic elastic ground. Géotechnique 2015;2(4):287–95. [21] Guo ZH, Liu XR, Zhu ZY. Elastic solution for a deep twin tunnel’s stress based on complex variable theory and the superpostion principle. J Eng Res 2017;5(2):68–86. [22] Ukadgaonker VG, Patil DB. Stress analysis of a plate containing two elliptical holes subjected to uniform pressures and tangential stresses on hole boundaries. J Eng Indust 1993;115(1):93–101. [23] Zhang LQ, Yue ZQ, Lee CF, Tham LG, Yang ZF. Stress solution of multiple elliptic hole problem in plane elasticity. J Eng Mech 2003;129(12):1394–407. [24] Zhang LQ, Lu AZ, Yue ZQ, Yang ZF. An efficient and accurate iterative stress solution for an infinite elastic plate around two elliptic holes. subjected to uniform loads on the hole boundaries and at infinity. Eur J Mech - A/Solids 2009;28(1):189–93. [25] Fu J, Yang JS, Yan L, Abbas SM. An analytical solution for deforming twin-parallel tunnels in an elastic half plane. Int J Numer Anal Meth Geomech 2014;39(5):524–38. [26] Wang HN, Zeng GS, Utili S, Jiang MJ, Wu L. Analtycial soltuions of stresses and displacements for deeply buried twin tunnels in viscoelastic rock. Int J Rock Mech Min Sci 2017;93:13–29. [27] Wang HN, Wu L, Jiang MJ, Song F. Analytical stress and displacement due to twin tunneling in an elastic semi-infinite ground subjected to surcharge loads. Int J Numer Anal Meth Geomech 2018;42(6):809–28. [28] Jaeger JC, Cook NGW, Zimmerman RW. Fundamentals of rock mechanics. 4th ed. Oxford: Blackwell Publishling; 2007. [29] Unlu T, Gercek H. Effect of Poisson's ratio on the normalized radial displacements occurring around the face of a circular tunnel. Tunn Undergr Space Technol 2003;18(5):547–53. [30] Vlachopoulos N, Diederichs MS. Improved longitudinal displacement profiles for convergence confinement analysis of deep tunnels. Rock Mech Rock Eng 2009;42:131–46. [31] Paraskevopoulou C, Diederichs MS. Analysis of time-dependent deformation in tunnels using the Convergnece-Confinement Method[J]. Tunn Undergr Space Technol 2017;71:62–80. [32] Muskhelishvili NI. Some basic problems of the mathematical theory of elasticity, 2nd ed. Groningen: Noordhoff International Publishing; 1977.

References [1] Chang CT, Sun CW, Duann SW, Hwang RN. Response of a Taipei Rapid Transit System (TRTS) tunnel to adjacent excavation. Tunnel Undergr Space Technol Incorporat Trenchless Technol Res 2001;16(3):151–8. [2] Chakeri H, Ozcelik Y, Unver B. Investigation of ground surface settlement in twin tunnels driven with EPBM in urban area. Arab J Geosci 2015;8(9):7655–66. [3] Liang RZ, Xia TD, Hong Y, Yu F. Effects of above-crossing tunnelling on the existing shield tunnels. Tunn Undergr Space Technol 2016;58:159–76. [4] Sulem J, Panet M, Guenot A. Closure analysis in deep tunnels. Int J Rock Mech Min Sci Geomech Abstracts 1987;24(3):145–54. [5] 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. [6] Chehade FH, Shahrour I. Numerical analysis of the interaction between twin-tunnels: influence of the relative position and construction procedure. Tunn Undergr Space Technol 2008;23(2):210–4. [7] Liang QG, Li J, Li DW, Qu EF. Effect of blast-induced vibration from new railway tunnel on existing adjacent railway tunnel in Xinjiang, China. Rock Mech Rock Eng 2013;46(1):19–39. [8] Do NA, Dias D, Oreste P. 3D numerical investigation on the interaction between mechanized twin tunnels in soft ground. Environ Earth Sci 2015;73(5):2101–13. [9] Das R, Singh PK, Kainthola A, Panthee S, Singh TN. Numerical analysis of surface subsidence in asymmetric parallel highway tunnels. J Rock Mech Geotech Eng 2017;9:170–9. [10] Dang VK, Dias D, Do NA, Vo TH. Impact of blasting at tunnel face on an existing adjacent tunnel. Int J Geomate 2018;15(47):22–31. [11] Peck RB. Deep excavations and tunnelling in soft ground. In: Proc 7th int conf on SMFE; 1969. [12] Mair RJ, Taylor RN, Burland JB. Prediction of ground movements and assessment of risk of building damage due to bored tunneling. Geotech Aspects Undergr Constr Soft Ground; 1996. [13] Suwansawat S, Einstein HH. Describing settlement troughs over twin tunnels using a superposition technique. J Geotech Geoenviron Eng 2007;133(4):445–68. [14] Ling CB. On the stresses in a plate containing two circular holes. J Appl Phys 1948;19(1):77–82. [15] Zimmerman RW. Second-Order approximation for the compression of an elastic plate containing a pair of circular holes. ZAMM-J Appl Math Mech: Zeitschrift für Angewandte Mathematik und Mechanik 1988;68(11):575–7. [16] Kooi CB, Verruijt A. Interaction of circular holes in an infinite elastic medium. Tunn Undergr Space Technol 2001;16(1):59–62.

17

Computers and Geotechnics 120 (2020) 103385

H.N. Wang, et al. [33] Wang HN, Chen XP, Jiang MJ, Song F, Wu L. The analytical predictions on displacement and stress around shallow tunnels subjected to surcharge loadings. Tunn Undergr Space Technol 2018;71:403–27. [34] Verruijt A. Deformations of an elastic half plane with a circular cavity. Int J Solids Struct 1998;35(21):2795–804. [35] Song F, Wang HN, Jiang MJ. Analytical solutions for lined circular tunnels in viscoelastic rock considering various interface conditions. Appl Math Model 2018;55:109–30.

[36] Wang HN, Zeng GS, Jiang MJ. Analytical stress and displacement around non-circular tunnels in semi-infinite ground. Appl Math Model 2018;63:303–28. [37] Hasanpour R, Chakeri H, Ozcelik Y, Denek H. Evaluation of surface settlements in the Istanbul metro in terms of analytical, numerical and direct measurements. Bull Eng Geol Environ 2012;71(3):499–510. [38] Cristescu ND. Rock rheology. Dordrecht: Kluwer Academic Publishers; 1989. [39] Lu AZ, Zeng XT, Xu Z. Solution for a circular cavity in an elastic half plane under gravity and arbitrary lateral stress. Int J Rock Mech Min Sci 2016;89:34–42.

18