Elastic-slip interface effect on dynamic stress around twin tunnels in soil medium subjected to blast waves

Elastic-slip interface effect on dynamic stress around twin tunnels in soil medium subjected to blast waves

Computers and Geotechnics xxx (xxxx) xxxx Contents lists available at ScienceDirect Computers and Geotechnics journal homepage: www.elsevier.com/loc...

1MB Sizes 0 Downloads 25 Views

Computers and Geotechnics xxx (xxxx) xxxx

Contents lists available at ScienceDirect

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

Research Paper

Elastic-slip interface effect on dynamic stress around twin tunnels in soil medium subjected to blast waves ⁎

Xue-Qian Fanga,b,c, , Teng-Fei Zhanga, Bai-Lin Lid, Rui-Jia Yuane a

Department of Engineering Mechanics, Shijiazhuang Tiedao University, Shijiazhuang 050043, China Hebei Key Laboratory of Intelligent Materials and Structures Mechanics, Shijiazhuang 050043, China c State Key Laboratory of Mechanical Behavior and System Safety of Traffic Engineering Structures, Shijiazhuang 050043, China d School of Civil Engineering, Shijiazhuang Tiedao University, Shijiazhuang 050043, China e College of Aerospace Engineering, Nanjing University of Aeronautics and Astronautics, Nanjing, Jiangsu 21106, China b

A R T I C LE I N FO

A B S T R A C T

Keywords: Twin tunnels Elastic-slip interface Blast waves Dynamic response

Blasting excavation of deep-buried structures is a complex dynamic process which can cause blast loading and insitu stress redistribution. An elastic-slip interface model is introduced to analyze the dynamic response of twin tunnels subjected to blast waves. The blast waves is expressed by Fourier-Bessel expansion method of wave function, and the wave field is superimposed by Graf coordinate transformation method. To analyze the interface effect, the elastic and slip coefficients of interface are introduced to describe the interface properties. Through numerical examples, it is found that the dynamic stress at the upper points of tunnels is larger than that at the sides of tunnels. However, the interface effect on the dynamic stress at the sides of tunnels is larger than that at the upper points, especially a small distance between twin tunnels is selected. Due to the interface effect, not only compressive but also tensile deformations can be seen, and the large oscillation of dynamic stress at the sides of tunnels can be observed. Excellent agreement with existing numerical results validates this dynamic model.

1. Introduction

stresses and liquefaction of the soil, and it is found that the gas in the tunnel makes the soil stiffness stronger [4]. The dynamic response of a long cavity subjected to internal transient loads was studied by Wang et al., and the surrounding soil was assumed to be a nearly-saturated porous medium. The exact solutions for three types of transient loads were obtained in the Laplace transform domain [5]. Gao et al. [6] presented the three-dimensional dynamic response of a cylindrical lined tunnel in full-space saturated soil due to an internal blast loading. The dynamic deformation of the tunnel lining caused by a blast-induced pressure pulse was investigated by Osinov et al., and it is concluded that the explosion in a tunnel can produce not only compressive but also tensile deformations around the tunnel [7]. Li and Li derived the theoretical formulas for the dynamic response of underground circular tunnel induced by blasting loads, and a two-dimensional numerical model was established by the discrete element code [8]. The three-dimensional dynamic response of a tunnel embedded in soil medium subjected to blasting waves was studied analytically. By simplification, the solutions for the displacement, stress and pore water pressure were derived in time domain by means of the Fourier and Laplace transforms [9]. Li et al. established a theoretical formulation to investigate the

Deep-buried structures are excavated by drilling and blasting techniques which have intensive adaptability to many kinds of geological conditions. The rock mass is excited by blasting waves, and the in-situ stress exerted on the surrounding medium is suddenly released. The released stress in the rock mass during the blasting excavation will propagate in the form of stress waves. The stress waves can cause severe damage to the adjacent structures such as the existing tunnels and pipelines, and ultimately worsen the overall serving behavior of the underground structures [1–3]. The tunnels in the rock mass can result in more complicated propagation of blast waves. The blast-induced stress waves can change the directions of the maximum and minimum principal stresses around the tunnels. The blast-created tensile failure will bring great damage to the adjacent engineering structures. To reduce the blast-induced damage or failure of adjacent tunnels, lots of theoretical and numerical methods dealing with the excavation-induced damage have been proposed. A theoretical modelling of waves in saturated granular soil embedded with a cylindrical tunnel was constructed to analyze the effective



Corresponding author. E-mail address: [email protected] (X.-Q. Fang).

https://doi.org/10.1016/j.compgeo.2019.103301 Received 6 July 2019; Received in revised form 15 September 2019; Accepted 13 October 2019 0266-352X/ © 2019 Elsevier Ltd. All rights reserved.

Please cite this article as: Xue-Qian Fang, et al., Computers and Geotechnics, https://doi.org/10.1016/j.compgeo.2019.103301

Computers and Geotechnics xxx (xxxx) xxxx

X.-Q. Fang, et al.

discussed [25]. By using the spring elastic model, the anti-plane dynamic response of a non-circular tunnel embedded in the anisotropic medium was investigated, and the imperfect interface was simulated by the elastic spring [26]. However, the imperfect interface around twin lining tunnels subjected to blast waves has not addressed. Actually, in excavating deep-buried structures, the imperfect interface effect on the dynamic stress and rock damage play an important role in predicting the dynamic response of underground tunnels. The purpose of this paper is to investigate the dynamic response around twin lining tunnels with consideration of imperfect interface effect. The transient response of blasting loadings is expressed by using Fourier transformation. To facilitate the theoretical calculations, the wave function expansion method and Graf coordinate transform method are combined to obtain the analytical solutions. The dynamic stress under different interface coefficients is discussed in detail. The interface effect with different distances between twin tunnels is also examined.

dynamic stress concentration factor around a circular tunnel subjected to blasting stress, and a two-dimensional numerical model was also presented [10]. The experimental works including laboratory test and field measurements have been attracting lots of attention. By simulating the field tests and modifying the empirical detonation pressure equations, the effect of tunnel depth and explosive amount on the lining stresses was analyzed [11]. Yang et al. performed the physical tests on the dynamic response of tunnels, and the effect of long-term train loads simulated by an electromagnetic shaker was analyzed [12]. Kristoffersena et al. used the digital image correlation to measure deformations in submerged floating tunnels in concrete with circular and rectangular cross-sections subjected to blast loading [13]. Ma et al. applied a qualitative method to evaluate tunnel ground vibrations induced by fault-slip in underground mine, and the peak ground velocity and the peak ground acceleration were analyzed [14]. Some numerical methods were also introduced to investigate the dynamic behavior of underground structure subjected to blast waves. Buonsanti and Leonardi presented the finite element analysis on the thermal stress distribution of an underground structure subjected to a fire generation and a subsequent explosion [15]. Li et al. studied the influence of two empty holes on the cracks propagating behavior in rock mass by using blasting experiments, and the numerical models were established by using AUTODYN code [16]. Most recently, discrete element method was used to investigate the fracture extension in jointed rock masses subjected to explosion loading [17]. The excavated section of the tunnel expresses significant effect on the structural integrity of medium and leads to more complicated transmissions and reflections of blast waves. To exercise control over the blast vibrations, it is requisite to make accurate predictions of PPV (peak particle velocity) locations. The predictions of PPV resulting from tunnel blasts have been attracting more and more attention in past years. Singh conducted the field tests at seven coal mines in India to investigate the influence of heavy blasting on the adjacent underground coal mine workings, and the PPV threshold for the safety of underground structures was obtained. It is found that the explosion-induced vibrations can result in severe damage to the roof [18]. Xia et al. investigated the effect of underground water tunnel on the ground vibration resulting from a series of blasts, and an empirical formula was proposed to calculate the PPV [19]. Salum and Murthy proposed a methodology to characterize the rock mass quality subjected to drilling and blasting, and various theoretical approaches for controlling the PPV for the various types of rock mass were given [20]. Lak et al. derived a two dimensional elasto-dynamic Green's function in terms of displacement, the PPV related to the displacement Green's function was analyzed and discussed [21]. Although some progress has been made in predicting the dynamic response of tunnels subjected to blast waves, there are still some issues for increasing the strength of underground tunnels which attract little attention. The interface properties play an important role in predicting the dynamic response of underground tunnels. Most recently, some interface models were proposed, and the interface effect on the dynamic strength of tunnels was analyzed [22–26]. Yi et al. investigated the dynamic stress concentration around a circular lined tunnel with an imperfect interface subjected to cylindrical P-waves, and the imperfect interface was modeled with a spring model [22]. Fang et al. proposed an elastic-slip interface model to study the dynamic stress distribution around the circular tunnel subjected to SH waves, and the displacement contours and dynamic stresses with different interface coefficients were analyzed [23]. The dynamic response of a circular tunnel with an imperfect interface subjected to plane SV-waves was investigated, and the linear spring model was presented to depict the imperfect interface between the rock mass and the lining [24]. Huang et al. employed a viscous-slip interface model to simulate the imperfect contact between the tunnel lining and the surrounding rock mass, and the dynamic stress concentration around twin shallowly buried lining tunnels was

2. Problem formulation Twin lining tunnels embedded in the rock mass are considered, as depicted in Fig. 1. The inner and outer radii of tunnels are denoted by a1 and a2 . The distance between the two centers of tunnels is d. It is assumed that the tunnels are so long. The blast waves with wave frequency of ω propagate in the rock mass. The rock mass and tunnel lining are supposed to be elastic and isotropic. According to the characters of tunnels and waves, the plane strain problem is employed. Around the twin lining tunnels, the imperfect interface is assumed, as shown in Fig. 1, and the interface properties are optimized to increase the strength of tunnels. During the blast process, the pressure variation with the time includes three steps: the increase of loadings, the initial expansion of gas and the blast of gas. The pressure resulting from blast impinges the tunnel in very short time, and the rock mass near the blast point is pulverized. The blast gas expands further and impacts the surrounding rock mass. The crack in the impacted rock mass makes the pressure decrease. When the pressure decreases to the original stress, the stress begins to relax. For simplification, the blast loading is substituted by a triangular load with one peak. The blast process includes the linear loading and unloading. To show how the formulations and the boundary conditions are used, as shown in Fig. 2. According to Fig. 2, the following formulations and the boundary conditions are used to solve this problem. The blast wave is a kind of compressional wave. From the beginning of reaching the boundary of tunnel, the propagating time of incident waves is normalized as [10]

τ=

cp t a2

(1)

where cp is the wave speed of compressional wave, t is the time and τ is the normalized time.

Fig. 1. Twin tunnels with elastic-slip interface subjected to blast loading. 2

Computers and Geotechnics xxx (xxxx) xxxx

X.-Q. Fang, et al.

g (x i , t ) =

∫0 ∫0

=

2 πtr

∫0

∞ R (ω) sin ω (t − τ ) 1 2 dτ dω + tr π 0 ω R (ω) sin ω (t − τ ) dω ω

tr





∫t

t

r

2 −1 dτ ts − tr π

2 R (ω)[cos ω (t − tr ) − cos ωt ] dω − ω2 π (ts − tr ) R (ω)[1 − cos ω (t − tr )] dω ω2

∫0





(9)

When t ⩾ ts , Fig. 2. Flowchart explaining the algorithm.

g (x i , t ) =

To connect the steady response with transient response, the Fourier transformation technique is introduced. For the arbitrary aperiodic vibration of f (t ) , the transient response of system can be given by the following equation [10]

1 2π

g (x i , t ) =

∫0



∫−∞ χ (xi , ω) F (ω) e−iωtdω

(2)

where χ (x i , ω) is the admittance of system, F (ω) is the Fourier transformation of aperiodic vibration f (t ) . To obtain the transient response, the Fourier transformation of f (t ) is used. The Fourier transformation of Heaviside unit step function is given as

i , Imζ > 0, 2π ζ

F (ζ ) =

where Re (ζ ) = αa2 and α = ω cp is the wave number of compressional wave. The admittance of system is simplified as

2 π

∫0



R (ω) cos ωtdω =

∫0



I (ω) sin ωtdω.

2 π

∫0

t



∫0



R (ω) cos ωτdτ =

2 π

∫0



R (ω) sin ωt dω. ω

∞ R (ω) sin ω (t − τ ) 1 2 dτ dω tr π 0 ω ∞ R (ω)(1 − cos ωt ) 2 dω = πtr 0 ω2 t

∫0



R (ω)[cos ω (t − ts ) − cos ω (t − tr )] dω. ω2

(10)





εn1 in1 Jn1 (αr1) cos(n1 θ1 ) e−iωt

(11)

n1= 0

where φ0 is the vibration amplitude and Jn ( ∗ ) is the Bessel function of the first kind. In the coordinate system o2 x2 y2 , the incident wave can be expressed as [27] ∞

φi2 = φ0

(7)



εn2 in2 Jn2 (αr2) cos(n2 θ2 ) e−iωt

(12)

n2 = 0

{

1 n=0 . Superscripts 1 and 2 denote the coordinate where εn = 2 n⩾1 systems of o1 x1 y1 and o2 x2 y2 . When the incident wave impacts the tunnel lining, the scattered P and SV waves come into being. The scattered P and SV waves can be expressed as





2 π (ts − tr )

φi1 = φ0

(6)

When 0 ⩽ t < tr ,

∫0



In the following, the wave function expansion method is introduced to express the wave field in the tunnel and soil medium. By using wave function expanded method, the incident wave can be expanded as [27]

where tr and ts are the standardized rise time and total time of blasting loading, respectively.

g (x i , t ) =

2 −1 dτ ts − tr π

(5)

For the transient response under blast loadings, the input function can be expressed as

t<0 ⎧ 0, t ⎪ ⩽ , 0 t < tr ⎪ tr f (t ) = ts − t . ⎨ , tr ⩽ t < ts ⎪ ts − tr ⎪ 0, t ⩾ ts ⎩

2 πtr



ts

R (ω)[cos ω (t − tr ) − cos ωt ] dω ω2

=

∫0

r

3.1. Wave fields in the soil medium

By using the response of Heaviside step function, Eq. (5) can be written as

gh (x i , t ) =



∫t

Blast waves is a kind of compressional wave. The transient loading and unloading can be expressed by superimposing the harmonic wave of all frequencies. So, the dynamic response under the compressional wave should be given firstly.

(4)

2 π



3. The wave fields under compressional waves

The transient response can be expressed by using sine and cosine functions. When the input signal is impulse function, the impulse response can be expressed as [10]

gδ (x i , t ) =

∞ R (ω) sin ω (t − τ ) 1 2 dτ dω + tr π 0 ω R (ω) sin ω (t − τ ) dω ω

tr

Then, the transient response of tunnels under blast loadings can be obtained by using Eqs. (8)–(10). After obtaining the numerical result of R (ω) under all wave numbers, the transient response is substituted by trapezoid function.

(3)

χ (x i , ω) = R (ω) + iI (ω).

∫0



φs1 =

(8)



As1, n1 Hn(1) (ksα r1) cos(n1 θ1) 1

(13)

n1= 0 ∞

When tr ⩽ t < ts ,

ψs1 =



Bs1, n1 Hn(1) (ksβ r1) sin(n1 θ1) 1

(14)

n1= 0

Hn(1) 1

( ∗ ) denotes where As1, n1 and Bs1, n1are the expanded coefficients, the n1 order Hankel function of the first kind, ksα and ksβ are the wave 3

Computers and Geotechnics xxx (xxxx) xxxx

X.-Q. Fang, et al.

numbers of P and SV waves. Similarly, the scattered P and SV waves can be written as ∞

φs2 =

As2, n2 Hn(1) 2 (ksα r2 ) cos(n2 θ2 )

∑ n2 = 0

(15)



ψs2 =



Bs2, n2 Hn(1) 2 (ksβ r2 ) sin(n2 θ2 )

n2 = 0

(16)

where As2, n2 and Bs2, n2 are the expanded coefficients. Then, the potential functions of waves in the soil medium are written as

φs = φi + φs1 + φs2

(17)

ψs = ψs1 + ψs2 .

(18)

Fig. 3. Comparison with existing numerical results.

4. Boundary conditions and solving the expanded coefficients For perfect contact surface, the displacements and stresses are continuous at the boundary. However, the ideal boundary conditions may result in large stresses around the tunnel. To optimize the interface properties, an elastic-slip interface model is proposed. In this model, the elastic and slip characters which are very common in practical engineering, are considered.

3.2. Wave fields in the lining Around the lining, there exist two interfaces. At the outer interface, the ingoing P and SV waves come into being. At inner interface, the outgoing P and SV waves around the lining can be seen. The ingoing and outgoing P and SV waves are expanded, in the coordinate system of o1 x1 y1, as

1 1 σrr1 ,(s) = σrr1 ,(l) , σrθ ,(s ) = σrθ,(l) ,

ur1,(s) = ur1,(l) +



φl11 =



A11, n1 Jn1 (klα r1) cos(n1 θ1)

n1= 0

=

ur2,(s) = ur2,(l) +



B11, n1 Jn1 (klβ r1) sin(n1 θ1)

n1= 0

(20)

σrr1 ,(l)



φl12 =



A12, n1 Hn(1) (klα r1) cos(n1 θ1) 1

n1= 0



B12, n1 Hn(1) (klβ r1) sin(n1 θ1). 1

n1= 0

(22)





A21, n2 Jn2 (klα r2) cos(n2 θ2)

n2 = 0

(23)



B21, n2 Jn2 (klβ r2) sin(n2 θ2)

n2 = 0

(24)



A22, n2 Hn(1) (klα r2) cos(n2 θ2) 2

n2 = 0

(25)





B22, n2 Hn(1) (klβ r2) sin(n2 θ2) 2

n2 = 0

(32)

kr

(33)

, uθ2,(s) = uθ2,(l) +

σrθ2 ,(s) kθ

+ Dσrθ2 ,(s) , (r2 = a2),

(34)

= 0, (r1 = a1)

(35)

S 2(3) 4 × 2, n2 × (26)

(27)

ψl1 = ψl11 + ψl12 .

∑ n1= 0

where klα and klβ are the wave numbers of P and SV waves in the lining. Then, the total potential function in the lining structure can be written, in the coordinate system of o1 x1 y1, as

φl1 = φl11 + φl12

(36)



∑ n2 = 0

+(3) ⎡T 1α As2, n2 ⎤ ⎥ ⎢ −(3) T 1 Bs2, n2 ⎦ ⎣ β n



ψl22 =

1 + Dσrθ ,(s ) , (r1 = a2 ),

1 ⎡ A11, n1 ⎤ + L1(3) ⎡ A12, n1 ⎤ = −S1(1) ⎡ φ0 εn1 i ⎤ + L1(1) 4 × 2, n1 × 4 × 2, n1 × 4 × 2, n1 × ⎢ B11, n1 ⎥ ⎢ B12, n1 ⎥ 0 ⎦ ⎣ ⎣ ⎦ ⎣ ⎦ (37)



φl22 =

= 0,

1 σrθ ,(l)

σrr2 ,(s)

⎡ As1, n1 ⎤ + S1(3) S1(3) 4 × 2, n1 × 4 × 2, n1 × ⎢ Bs1, n1 ⎥ ⎣ ⎦



ψl21 =



where ur ,(s) and uθ,(s) are the displacements in the soil medium, σrr ,(s) and σrθ,(s) are the displacements in the soil medium. σrr ,(l) and σrθ,(l) denote the stresses in the tunnel lining. Note that superscripts 1 and 2 represent the tunnels 1 and 2. kr and kθ are the elastic coefficients of interface in the radial and circumferential directions. D is the circumferential slip coefficient of interface. Substituting Eqs. (11)–(29) into Eqs. (31)–(36), the following system of linear algebraic equations can be obtained

In the coordinate system of o2 x2 y2 , the ingoing and outgoing P and SV waves are expressed as

φl21 =

, uθ1,(s) = uθ1,(l) +

σrr2 ,(l) = 0, σrθ2 ,(l) = 0, (r2 = a1), (21)



ψl12 =

kr

(31) 1 σrθ ,(s )

σrr2 ,(s) = σrr2 ,(l) , σrθ2 ,(s) = σrθ2 ,(l) , (19)



ψl11

σrr1 ,(s)

+(3) ⎡T 2α As1, n1 ⎤ (3) ⎡ As2, n2 ⎤ ⎢ −(3) ⎥ + S 24 × 2, n2 × ⎢ B ⎥ T 2 Bs1, n1 ⎣ s2, n2 ⎦ ⎣ β ⎦

2 ⎡ A21, n2 ⎤ + L1(3) ⎡ A22, n2 ⎤ = −S 2(1) ⎡ φ0 εn2 i ⎤ + L1(1) 4 × 2, n2 × 4 × 2, n2 × 4 × 2, n2 × ⎢ B21, n2 ⎥ ⎢ B22, n2 ⎥ 0 ⎦ ⎣ ⎣ ⎦ ⎣ ⎦ (38)

n

(28)

A11, n1 ⎤ A12, n1 ⎤ L22(1)× 2, n1 × ⎡ + L22(3)× 2, n1 × ⎡ =0 ⎢ B11, n1 ⎥ ⎢ B12, n1 ⎥ ⎣ ⎦ ⎣ ⎦

(39)

In the coordinate system of o2 x2 y2 , the total potential function in the lining structure can be expressed as

A21, n2 ⎤ A22, n2 ⎤ L22(1)× 2, n2 × ⎡ + L22(3)× 2, n2 × ⎡ =0 ⎢ B21, n2 ⎥ ⎢ B22, n2 ⎥ ⎣ ⎦ ⎣ ⎦

(40)

φl2

(29)

where

(30)

T 1±q (3) =

=

φl21

+

φl22

ψl2 = ψl21 + ψl22 .

4

εn1 [± (−1)n2 Jn(3) (ksq d ) + Jn(3) (ksq d )] 1+ n2 1− n2 2

Computers and Geotechnics xxx (xxxx) xxxx

X.-Q. Fang, et al.

Fig. 4. DSCF variation at the positions of θ = 0 and θ = π 2 with different elastic coefficients of interface (kθ = 1011 N m3 and D = 10−11m3 N ).

T 2±q (3) =

εn2 [± (−1)n2 Jn(3) (ksq d ) + Jn(3) (ksq d )] 1+ n2 1− n2 2

S1±4 ×(j2,) n

(j ) (j ) ⎡ S11σrr (n, a2) ± S12σrr (n, a2) ⎤ ⎥ ⎢ (j ) (j ) ⎢ ∓ S11σrθ (n, a2 ) S12σrθ (n, a2 ) ⎥ =⎢ ⎥ (j ) (j ) ⎢ S11ur (n, a2) ± S12ur (n, a2) ⎥ ⎢ ∓ S11(j) (n, a ) S12(j) (n, a ) ⎥ 2 2 uθ uθ ⎦ ⎣

L1±4 ×(j2,) n

− L11(σjrr) (n, a2) ⎡ ⎢ ± L11(σjrθ) (n, a2) ⎢ ⎢ 1 − L11u(jr) (n, a2) + k ⎢ r =⎢ (j ) ( 11 ( , )) L n a − ⎢ 2 σrr ⎢ 1 ⎢ ± L11u(jθ) (n, a2) + k + D θ ⎢ ⎢ (± L11(j) (n, a2)) σrθ ⎣

(

DSCF is the ratio of the circumferential stress around the tunnel to the maximum dynamic stress. Thus, the DSCF around the twin tunnels is defined as ∗ DSCF = τθz = τθz τ0

where τ0 is the maximum magnitude of the stress in the incident di2 rection, and τ0 = μksβ φ0 . The following geometrical and material parameters used in the numerical examples are defined. The inner radius of tunnels a1 = 5m , the outer radius of tunnels is a2 = 1.1a1. The elastic modulus of soil medium is Es = 18.73Gpa and the Poisson ratio of soil medium is υs = 0.206. The elastic modulus of tunnel lining is El = 9.3Gpa and the Poisson ratio of tunnel lining is υl = 0.18 [10]. The standardized rise time is tr = 5 and the total time of blasting loading is ts tr = 5. To validate this present model, comparison with existing results is given in Fig. 3. In Ref. [10], only one tunnel subjected to blast waves is considered. If the distance between twin tunnels is defined as d = 50mm , the interacting effect between the tunnels can be ignored. If kr , kr and D are defined as kr = kr = 1015N mm3 and D = 0, the interface effect vanishes. Excellent agreement with the result in Ref. [10] can be seen in Fig. 2. Under the explosion loading, the dynamic stresses at the positions of θ = 0 (sides of tunnels) and θ = π 2 (upper points of tunnels) are both compress stress concentration. The dynamic stress at the position of θ = π 2 is larger than that at θ = 0 .

∓ L12(σjrr) (n, a2)

⎤ ⎥ ⎥ ⎥ 1 ∓ L12u(jr) (n, a2) + k ⎥ r ⎥ (∓ L12(σjrr) (n, a2)) ⎥ ⎥ 1 − L12u(jθ) (n, a2) + k + D ⎥ θ ⎥ ⎥ (−L12(σjrθ) (n, a2)) ⎦ − L12(σjrθ) (n, a2)

)

(

(41)

)

(j ) (j ) ⎡ L11σrr (n, a1) ± L12σrr (n, a1) ⎤ L2±2 ×(j2,) n = ⎢ ⎥. ∓ L11(σjrθ) (n, a1 ) L12(σjrθ) (n, a1) ⎦ ⎣

The elements of matrix S1, L1 and L2 are shown in Appendix. 5. Numerical examples and discussion Dynamics stress is an important character in predicting the dynamic response of underground tunnels. To analyze the imperfect interface effect on the dynamic stress of twin tunnels under blast loadings, the dimensionless dynamic stress is introduced. Following the work in Ref. [27], the dynamic stress concentration factor (DSCF) is defined. The

5.1. Effect of elastic coefficient of interface in the radial direction Fig. 4 illustrates the DSCF at the positions of θ = 0 and θ = π 2 with different elastic coefficients of interface (kθ = 1011 N m3 and D = 10−11m3 N ). The distance between twin tunnels is defined as 5

Computers and Geotechnics xxx (xxxx) xxxx

X.-Q. Fang, et al.

Fig. 5. DSCF variation at the positions of θ = 0 and θ = π 2 with different elastic coefficients of interface (kr = 1011 N m3 and D = 10−11m3 N ).

result in larger negative stress.

d = 15m , d = 20m and d = 25m . It can be seen that the effect of elastic coefficient of interface in the radial direction on the dynamic stress at the position of θ = 0 is greater than that at the position of θ = π 2 . This phenomenon can be interpreted as follows. The elastic property of interface buffers the impact of blast loadings, and the strength of blasting loadings becomes weaker. In the unloading process, the interface effect at the position of θ = π 2 is very small, so this process is stable. However, the large variation of stress at the position of θ = 0 can be observed. This phenomenon results from the repeated energy charge and release at the interface. Therefore, a large elastic coefficient of interface is preferable in the earthquake region. If the distance between twin tunnels is large, the interface effect becomes small. When the interacting effect between twin tunnels vanishes, only the dynamic stress near the position of θ = 0 shows variation with the elastic coefficient of interface, which results from the weak multiple scattering between twin tunnels.

5.3. Effect of slip coefficient of interface Fig. 6 illustrates the DSCF at the positions of θ = 0 andθ = π 2 with different slip coefficients of interface in the circumferential direction (kr = 1011 N m3 and kθ = 1011 N m3 ). The distance between twin tunnels is defined as d = 15m , d = 20m and d = 25m . Due to the slip of interface, the large oscillation of dynamic stress at the position of θ = 0 can be observed, especially a small distance between twin tunnels is selected. However, the effect of slip interface on the dynamic stress at the position of θ = π 2 is nearly unrelated with the distance between twin tunnels. When a large slip coefficient of interface is selected, the absorptive capacity of tunnels under blast waves decreases, and makes the waves diffuse around the tunnels. 6. Conclusion

5.2. Effect of elastic coefficient of interface in the circumferential direction

An elastic-slip interface is proposed to investigate the dynamic stress around twin tunnels under blast loading. The transient response is described by introducing the arbitrary aperiodic vibration. The wave fields in the soil medium and tunnel lining are expressed by using wave function expansion method. In numerical examples, some important findings are found.

Fig. 5 illustrates the DSCF at the positions of θ = 0 and θ = π 2 with different elastic coefficients of interface in the circumferential direction (kr = 1011 N m3 and D = 10−11m3 N ). The distance between twin tunnels is defined as d = 15m , d = 20m and d = 25m . Due to the interface effect, the large oscillation of dynamic stress at the position of θ = 0 can be observed, and the change of stress direction is also seen. If the distance between twin tunnels is smaller, the oscillation time of dynamic stress is longer. However, there is no oscillation of dynamic stress at the position of θ = π 2 . The interface effect at the position of θ = π 2 shows little variation with the distance between twin tunnels. The phenomenon results from the elastic property of interface. When the interacting effect between twin tunnels vanishes, the interface can

a. The effect of elastic coefficient of interface on the dynamic stress at the sides is much greater than that at the upper positions. b. Due to the effect of elastic coefficient of interface in the circumferential direction, not only compressive but also tensile deformations at the sides of tunnels can be seen. c. The dynamic stress decreases with decreasing slip coefficients of 6

Computers and Geotechnics xxx (xxxx) xxxx

X.-Q. Fang, et al.

Fig. 6. DSCF variation at the positions of θ = 0 and θ = π 2 with different slip coefficients of interface (kr = 1011 N m3 and kθ = 1011 N m3 ).

Acknowledgement

interface. Due to the slip of interface, the large oscillation of dynamic stress at the sides of tunnels can be observed, especially a small distance between twin tunnels is selected.

This work is supported by the project of Natural Science Foundations of Hebei Province (No. A2019210037), and the funded project for Innovative Graduates in Hebei, China (No. CXZZBS2018149).

The dynamic interface model and the developed formulations can provide an important solving method for underground tunnel design in earthquake region and significant reference for the interface design of tunnels in engineering construction. Appendix A

S11(σjrr) (n, r ) = 2μs

ksα (j) n2 − n 2 ⎤ (j ) Jn + 1 (ksα r ) + ⎡2μs − (λs + 2μs ) ksα Jn (ksα r ) ⎢ ⎥ r r2 ⎣ ⎦

(A1)

S12σ(jrr) (n, r ) = 2μs

n n − 1 (j) ⎡−ksβ Jn(j+) 1 (ksβ r ) + Jn (ksβ r ) ⎤ r⎣ r ⎦

(A2)

S11σ(jrθ) (n, r ) = 2μs

n n − 1 (j) ⎡−ksα Jn(j+) 1 (ksα r ) + Jn (ksα r ) ⎤ r⎣ r ⎦

(A3)

2 n2 − n ⎞ (j) ⎡ ksβ ⎤ ⎛ ksβ S12(σjrθ) (n, r ) = 2μs ⎢− Jn(j+) 1 (ksβ r ) + ⎜ − ⎟ Jn (ksβ r ) ⎥ 2 r 2 r ⎝ ⎠ ⎣ ⎦

S11u(jr) (n, r ) = −ksα Jn(j+) 1 (ksα r ) +

(A4)

n (j) Jn (ksα r ) r

(A5)

S12u(jr) (n, r ) =

n (j) Jn (ksβ r ) r

(A6)

S11u(jθ) (n, r ) =

n (j) Jn (ksα r ) r

(A7)

S12u(jθ) (n, r ) = ksβ Jn(j+) 1 (ksβ r ) −

n (j) Jn (ksβ r ) r

(A8) 7

Computers and Geotechnics xxx (xxxx) xxxx

X.-Q. Fang, et al.

L11σ(jrr) (n, r ) = 2μl

klα (j) n2 − n Jn + 1 (klα r ) + ⎡2μl − (λl + 2μl ) klα2 ⎤ Jn(j) (klα r ) ⎢ ⎥ r r2 ⎣ ⎦ (A9)

L12(σjrr) (n, r ) = 2μl

n n − 1 (j) ⎡−klβ Jn(j+) 1 (klβ r ) + Jn (klβ r ) ⎤ r⎣ r ⎦

(A10)

L11σ(jrθ) (n, r ) = 2μl

n n − 1 (j) ⎡−klα Jn(j+) 1 (klα r ) + Jn (klα r ) ⎤ r⎣ r ⎦

(A11)

2 n2 − n ⎞ (j) ⎡ klβ ⎤ ⎛ klβ L12(σjrθ) (n, r ) = 2μl ⎢− Jn(j+) 1 (klβ r ) + ⎜ − ⎟ Jn (klβ r ) ⎥ 2 r 2 r ⎝ ⎠ ⎣ ⎦

(A12)

n L11u(jr) (n, r ) = −klα Jn(j+) 1 (klα r ) + Jn(j) (klα r ) r

(A13)

L12u(jr) (n, r ) =

n (j) Jn (klβ r ) r

(A14)

L11u(jθ) (n, r ) =

n (j) Jn (klα r ) r

(A15)

L12u(jθ) (n, r ) = klβ Jn(j+) 1 (klβ r ) −

n (j) Jn (klβ r ) r

(A16)

[14] Ma J, Dong L, Zhao G, Li X. Qualitative method and case study for ground vibration of tunnels induced by fault-slip in underground mine. Rock Mech Rock Eng 2019;52(6):1887–901. [15] Buonsanti M, Leonardi G. 3-D simulation of tunnel structures under blast loading. Arch Civ Mech Eng 2013;13(1):128–34. [16] Li M, Zhu ZM, Liu RF, Liu B, Zhou L, Dong YQ. Study of the effect of empty holes on propagating cracks under blasting loads. Int J Rock Mech Min 2018;103:186–94. [17] Lak M, Marji MF, Bafghi AY. Discrete element modeling of explosion-induced fracture extension in jointed rock masses. J Min Environ 2019;10(1):125138. [18] Singh PK. Blast vibration damage to underground coal mines from adjacent open-pit blasting. Int J Rock Mech Min 2002;39(8):959–73. [19] Xia X, Li HB, Liu YQ, Yu C. A case study on the cavity effect of a water tunnel on the ground vibrations induced by excavating blasts. Tunn Undergr Sp Tech 2018;71:292–7. [20] Salum AH, Murthy V. Optimising blast pulls and controlling blast-induced excavation damage zone in tunnelling through varied rock classes. Tunn Undergr Sp Tech 2019;85:307–18. [21] Lak M, Marji MF, Bafghi AY, Abdollahipour A. Analytical and numerical modeling of rock blasting operations using a two-dimensional elasto-dynamic Green's function. Int J Rock Mech Min 2019;114:208–17. [22] Yi C, Zhang P, Johansson D, Nyberg U. Dynamic response of a circular lined tunnel with an imperfect interface subjected to cylindrical P-waves. Comput Geotech 2014;55:165–71. [23] Fang XQ, Zhang TF, Li HY. Elastic-slip interface effect on dynamic response of lined tunnel in a semi-infinite alluvial valley under SH waves. Tunn Undergr Sp Tech 2018;74:96–106. [24] Fan ZF, Zhang JC, Xu H. Theoretical study of the dynamic response of a circular lined tunnel with an imperfect interface subjected to incident SV-waves. Comput Geotech 2019;110:308–18. [25] Huang L, Liu ZX, Wu CQ, Liang JW. The scattering of plane P, SV waves by twin lining tunnels with imperfect interfaces embedded in an elastic half-space. Tunn Undergr Sp Tech 2019;85:319–30. [26] Zhang XP, Jiang YJ, Sugimoto S. Anti-plane dynamic response of a non-circular tunnel with imperfect interface in anisotropic rock mass. Tunn Undergr Sp Tech 2019;87:134–44. [27] Pao YH, Mow CC. Diffraction of elastic waves and dynamic stress concentrations. Russak (New York): Crane; 1973.

References [1] Dong LJ, Sun DY, Li XB, Ma J, Zhang LY, Tong XJ. Interval non-probabilistic reliability of surrounding jointed rock mass considering microseismic loads in mining tunnels. Tunn Undergr Sp Tech 2018;81:326–35. [2] Vosoughifar H, Madadi F, Rabiefar A. Modified dynamic stress concentration factor for twin tunnels using a novel approach of FEM-scattering. Tunn Undergr Sp Tech 2017;70:30–41. [3] Ma J, Dong L, Zhao G, Li X. Ground motions induced by mining seismic events with different focal mechanisms. Int J Rock Mech Min 2019;116:99–110. [4] Osinov VA. Blast-induced waves in soil around a tunnel. Arch Appl Mech 2011;81(5):543–59. [5] Wang Y, Gao GY, Yang J, Song J. The influence of the degree of saturation on dynamic response of a cylindrical lined cavity in a nearly saturated medium. Soil Dyn Earthq Eng 2015;71(8):27–30. [6] Gao M, Zhang JY, Chen QS, Gao GY, Yang J, Li DY. An exact solution for threedimensional (3D) dynamic response of a cylindrical lined tunnel in saturated soil to an internal blast load. Soil Dyn Earthq Eng 2016;90:32–7. [7] Osinov VA, Chrisopoulos S, Triantafyllidis T. Numerical analysis of the tunnel-soil interaction caused by an explosion in the tunnel. Soil Dyn Earthq Eng 2019;122:318–26. [8] Li CJ, Li XB. Influence of wavelength-to-tunnel-diameter ratio on dynamic response of underground tunnels subjected to blasting loads. Int J Rock Mech Min 2018;112:323–38. [9] Wang Y, Gao GY, Yang J. Three-dimensional dynamic response of a lined tunnel in a half-space of saturated soil under internal explosive loading. Soil Dyn Earthq Eng 2017;101:157–61. [10] Li XB, Li CJ, Cao WZ, Tao M. Dynamic stress concentration and energy evolution of deep-buried tunnels under blasting loads. Int J Rock Mech Min 2018;104:131–46. [11] Ramulu M, Chakraborty AK, Sitharam TG. Damage assessment of basaltic rock mass due to repeated blasting in a railway tunnelling project-a case study. Tunn Undergr Sp Tech 2009;24:208–21. [12] Yang W, Li L, Shang Y, Yan Q, Fang Y, He C, et al. An experimental study of the dynamic response of shield tunnels under long-term train loads. Tunn Undergr Sp Tech 2018;79:67–75. [13] Kristoffersena M, Minorettib A, Børvik T. On the internal blast loading of submerged floating tunnels in concrete with circular and rectangular cross-sections. Eng Fail Anal 2019;103:462–80.

8