Dynamic soil-tunnel interaction in layered half-space for incident P- and SV-waves

The dynamic soil-tunnel interaction is studied by the model of a rigid tunnel embedded in layered half-space, which is simplified as a single soil layer on elastic bedrock to the excitation of P- and SV-waves. The indirect boundary element method is used, combined with the Green’s function of distributed loads acting on inclined lines. It is shown that the dynamic characteristics of soil-tunnel interaction in layered half-space are different much from that in homogeneous half-space, and that the mechanism of soil-tunnel interaction is also different much from that of soil-foundation-superstructure interaction. For oblique incidence, the tunnel response for in-plane incident SV-waves is completely different from that for incident SH-waves, while the tunnel response for vertically incident SV-wave is very similar to that of vertically incident SH-wave.


Introduction
Dynamic soil-structure interaction (SSI) is an interdisciplinary field involving in the knowledge of soil and structural dynamics, earthquake engineering, and geophysics. Most studies on this problem mainly focus on soil-superstructure interaction, using a model of a rigid foundation with or without a building on it. For example, the classic solutions of a semi-circular rigid foundation with a shear wall on it were obtained by a kind of analytical method (Luco 1969;Trifunac 1972). More recently, the influences of site dynamic characteristics on SSI were studied separately using the same model in elastic layered half-space (Liang et al. 2013a, b) by indirect boundary element method.
The scholars have already obtained the solutions of dynamic responses of underground tunnel by analytical methods (Lee and Trifunac 1979) or numerical methods De Barros and Luco 1994) for several decades. However, the interaction between soil and underground structure, although as an important part of soil-structure interaction, are rarely studied up to now, only Hatzigeorgiou and Beskos (2010) compared damage evolution between lined tunnel and soil cavity to study the interaction between an underground runnel and the surrounding soil. Parvanova et al. (2014) analyzed the surface displacement and stress distribution along tunnel circumference to study respectively the interaction between one tunnel or twin tunnels and local topography.
In a companion paper (Fu J, Liang J and Qin L (2015). Dynamic soil-tunnel interaction in layered half-space for incident plane SH waves. In review, cited as ' '(Fu et al. 2015)'' in the following for convenience), the influence of site dynamic characteristics on soil-tunnel interaction is already studied by the model of a rigid tunnel embedded in layered half-space to the excitation of SH-waves, and the main milestones and methods on dynamic responses of underground tunnel are reviewed in the introduction. In this paper, the problem is continuously discussed using the same model to the excitation of P-and SV-waves, by indirect boundary element method combined with Green's functions of distributed loads acting on inclined lines (Wolf 1985). The further study on soil-tunnel interaction to incident surface waves is our work in the future.

Model
In Fig. 1, an underground lined tunnel is completely rigid and infinitely long, with outer radius of a, the inner radius of b, the mass of M 0 , and mass density of q 0 , also its embedded depth from ground surface to the center is d. It is bonded tightly to the surrounding soil at interface C without slippage. The layered half-space is simplified to a single soil layer with thickness D over bedrock. Both the soil layer and bedrock are elastic, homogeneous, and isotropic medium. The material parameters of the bedrock are characterized by shear-wave velocity b R , mass density q R , Poisson ratio m R , and damping ratio n R , while the material parameters of the soil layer are characterized by shear-wave velocity b L , mass density q L , Poisson ratio m L , and damping ratio n L . The plane P-waves or SV-waves are incident from depth D 0 with horizontal incident angle h, circular frequency x, and unit amplitude.
In order to use indirect boundary element method (IBEM) to solve the problem, the layered half-space should be divided into sub-layers, with the boundary C of the tunnel divided into N elements of straight lines meanwhile, and it is better to make all elements the same length in order that IBEM can perform best. Also, as the Green's functions used in this paper are distributed line loads in horizontally layered half-space, the elements should be symmetrical about z-axis.

Impedance function
In order to apply IBEM, a set of fictitious horizontal loads q j e ixt and vertical loads r j e ixt (j = 1, 2, …, N) which compose the fictitious load vector P ¼ q 1 ; q 2 ; . . .; q N ; r 1 ; r 2 ; . . .; is imposed onto every element as Fig. 2, with time factor e ixt is omitted hereafter. The values of these loads are all unknowns which should be determined by boundary condition that the tunnel produces the rigid displacement to the excitation of these loads. For in-plane excitation, the rigid displacement of the tunnel is D = [D x , au, D z ] T , with D x and D z being the horizontal and vertical displacements, and u being the rotational angle about its center, respectively. So on boundary C, the no slippage assumption gives Symbol U(x, z) is a two-dimensional vector whose elements represent the horizontal and vertical displacement at the point (x, z), respectively.
Under the excitation of fictitious loads, the displacements at the point (x, z) belonging to lth element C l can also be represented by in which g h ðx; zÞ is a 2 9 2N matrix of displacement Green's functions with g lj h and g lj v are the horizontal and vertical displacements at the point (x, z) of lth element when a unit distributed load p j or r j is imposed onto jth element (Wolf 1985).
If it is assumed that the symbol K is a 2N 9 3 matrix, with its three columns being the values of fictitious loads when the tunnel moves with unit horizontal displacement D x , unit rotation arclength au, and unit vertical displacement D z , respectively. Introducing Eqs. (3) and (5) into (2) gives For each column of K and X(x, z), every point on boundary C determines a set of 2 9 2N equations like Eq. (6), and if N target points on boundary C are chosen (usually one target point from one element in order that IBEM can perform best), there comes a set of 2N 9 2N equations, from which this column of K can be solved.
Then, the traction at the point (x, z) of lth element is in which Tðx; zÞ is a two-dimensional vector whose elements represent horizontal and vertical tractions at the point (x, z), respectively, and g t is a 2 9 2N matrix of traction Green's functions with P lj and H lj are the horizontal and vertical traction at the point (x, z) when a unit load p j or r j is imposed onto jth element (Wolf 1985), and e xl and e zl are the unit normal vector in x-direction and z-direction of lth element. Finally, the force vector F = [F x , M/a, F z ] T , with F x , M, and F z being the total horizontal force, rotational moment, and vertical force imposed on the tunnel, is obtained by integral with respect to the tractions along C Introducing Eq. (7) into (9) gives the desired relationship between tunnel displacement and the force imposed on it So impedance function matrix of the tunnel is and its form is as follows in which K HH , K MM , K VV , K MH , and K HM being the horizontal, rotational, vertical, and two-coupling impedance functions, respectively, with K MH = K HM . Taking K HH for example, it is convenient to write the impedance function as

Tunnel response
Effective input motion D is the tunnel displacement under harmonic-wave excitation, and it can be decomposed into two parts (Luco and Wong 1987) in which D 1 corresponds to the tunnel displacement when its mass M 0 is not taken into account (Luco 1986) with U f ðx; zÞ and T f ðx; zÞ are two-dimensional vectors, corresponding to the displacements and tractions in two directions of free-field response, respectively. Symbol D 2 is the additional displacement caused by inertia force F = [F x , M/a, F z ] T , with F Tx , M T , and F Tz being the horizontal force, rotational moment, and vertical force caused by tunnel mass, respectively, and based on the Earthq Sci (2015) 28 (4):275-284  277 concept of impedance functions, the additional displacement is solved by Fig. 3 Spectrum of tunnel impedance functions (d/a = 2, q R = q L = 2000 kg/m 3 , m R = m L = 0.25, damping ratio n R = 0.05 and n L = 0.02 for layered half-space, and n R = n L = 0.05 for For a rigid body with I 0 being the rotational inertia with respect to tunnel center. Introducing Eqs. (16) and (17) into (14) gives the final solution of tunnel displacement As the amplitude of incident wave is assumed to be unit 1, the tunnel displacements D x , au, and D z are all dimensionless, in fact they represent the amplification factor of incident excitation.

Numerical results and analysis
3.1 Impedance function Figure 3 is the impedance function of tunnel in frequency domain of homogeneous half-space and layered half-space. The embedded depth of the tunnel is d/a = 2. The parameters of homogeneous half-space are q R = q L = 2000 kg/ m 3 , m R = m L = 0.25, and n R = n L = 0.05. While the parameters of layered half-space are q R = q L = 2000 kg/ m 3 , m R = m L = 0.25, n R = 0.05, n L = 0.02, with the shearwave velocity ratio of the soil layer to the bedrock (''bedrock stiffness'' for short) varying with four values b R /b L = 2, 3, 5, and ?, and the ratio of the soil-layer thickness to the tunnel radius (''soil-layer thickness'' for short in the following) varying with three values D/a = 4, 6, and 8, so the tunnel is completely embedded within the soil layer in order to be convenient to analyze the site effect on tunnel response.
The impedance functions of homogeneous half-space are vibrating functions and the impedance functions of layered half-space vibrate around that of homogeneous half-space, because the layered half-space involves the dynamic characteristics of site while homogeneous halfspace cannot reflect these characteristics. When the bedrock stiffness b R /b L increases, the influence of the site dynamic characteristics also increase in the way that the vibrating period of impedance functions keeps invariable, and the shape of the curves is such as that the curves are multiplied by a factor in y-axis. Figure 4a and b are the spectrum of tunnel horizontal displacement D x , rotational-arc au, and vertical displacement D z in homogeneous half-space and layered half-space with bedrock stiffness b R /b L = 2 for incident P-wave and SV-wave, respectively. The mass density of the tunnel is q 0 = 2500 kg/m 3 , and the dimensionless tunnel mass is M 0 /M s = 1/4 with M s being the mass of soil replaced by the tunnel. The incident P-and SV-wave comes from D 0 / a = 8 for all sites with four incident angle h = 5°, 30°, 60°, and 90°, and the parameters of the half-space are the same as that in Fig. 3. The tunnel displacement spectrum for incident SH-wave is also plotted in Fig. 4c (Fu et al. 2015) for comparison, there is only out-of-plane translational displacement D j j in this condition. For oblique incidence, the spectrum of in-plane displacements is more complicated than that of incident SHwave, especially for small incident angle (h = 5°and 30°). Also, the displacements in layered half-space increase with incident angle increasing for incident SH-waves, while this is not the fact for in-plane excitation. For vertical incidence (h = 90°), the symmetry of the tunnel gives that D P

Tunnel response
x ¼ au P ¼ D SV z ¼ 0, it is also noticed that the tunnel displacements in layered half-space are larger than those in homogeneous half-space.
The tunnel displacement spectrum of homogeneous half-space is much smoother than that of layered halfspace; also, there is an evident peak for layered half-space on the spectrum of translational displacements D x and D z for large incident angle (h = 60°and 90°), while the peak does not exist for homogeneous half-space. This is because the site dynamic characteristics introduce much influence on tunnel response of layered half-space, while the homogeneous half-space does not involve these characteristics. It is noticed that there is also a peak on the spectrum of D x for incident SV-wave with small angle (h = 5°and 30°), but there is no interest in this peak which does not reflect the site dynamic characteristics. Figure 5a and b are the spectra of tunnel translational displacement in homogeneous half-space and layered halfspace for vertically incident P-wave (D z ) and SV-wave (D x ), respectively. The dimensionless tunnel mass is M 0 / M s = 0, 1/4, 1/2, with other parameters the same as those in Figs. 3 and 4. For comparison, the tunnel displacement spectrum for vertical incident SH-wave is also plotted in Fig. 5c (Fu et al. 2015). It is noticed that for both incident P-wave and SV-wave, the peak value becomes larger with soil-layer thickness increasing, which is similar to incident b Fig. 5 Spectrum of tunnel displacement in homogeneous half-space and layered half-space for vertical incidence (h = 5°, d/a = 2, q R = q L = 2000 kg/m 3 , m R = m L = 0.25, q 0 = 2500 kg/m 3 , D 0 / a = 8, damping ratio n R = 0.05 and n L = 0.02 for layered halfspace, n R = n L = 0.05 for homogeneous half-space). a Tunnel displacement |D z | for vertically incident P-wave. b Tunnel displacement |D x | for vertically incident SV-wave. c Tunnel out-of-plane displacement for vertically incident SH-wave (Fu et al. 2015) SH-waves because the path the incident wave propagates and amplifies is longer in thicker soil layer (Fu et al. 2015). The tunnel mass have little influence on tunnel displacement spectrum for both incident P-wave and SV-wave as the condition of incident SH-wave because the tunnel mass itself is small. It can be concluded that the kinematic interaction also overwhelmingly dominates for in-plane excitation, and the inertia interaction can hardly have influence on soil-tunnel interaction.
For free-field ground motion to vertically incident Pwave, the frequencies for which interference produces maximum response of the soil layer are (''resonant frequencies'' for short) So D/a = 4 corresponds to x a a/b L = 0.68, 2.04, 3.40, …, D/a = 6 corresponds to x a a/b L = 0.45, 1.36, 2.27, …, D/a = 8 corresponds to x a a/b L = 0.34, 1.02, 1.70, … and so on. It is observed that the peak frequency of tunnel displacement evidently becomes lower with soil-layer thickness increasing, and it is lower than the first resonant frequency x a of free-field response for b R /b L = 2, while higher than x a for b R /b L = 3, 5 and ?. Nevertheless, the difference between the peak frequency of tunnel displacement and x a is not large. While for free-field ground motion to incident SV-wave, the resonant frequencies of the soil layer are So D/a = 4 corresponds to x b a/b L = 0.39, 1.18, 1.96, …, D/a = 6 corresponds to x b a/b L = 0.26, 0.79, 1.31, …, D/a = 8 corresponds to x b a/b L = 0.20, 0.59, 0.98, … and so on. The peak frequency of tunnel displacement also becomes lower with soil-layer thickness increasing as incident P-wave, but it is lower than the first resonant frequency x b of free-field response for all bedrock stiffness, and the difference between the peak frequency of tunnel response and x b is not large either. Moreover, it is noticed that the peak frequency of tunnel response for both incident P-wave and SV-wave becomes lower with bedrock stiffness decreasing (there exists abnormal case for b R /b L = ? to incident P-wave), but this phenomenon is not evident. While in the papers by Liang et al. (2013a, b) studying the soil-foundation-superstructure interaction, although the foundation is also assumed to be completely rigid, the difference between the peak frequency of foundation displacement spectrum and the resonant frequency of freefield response is much larger, and the variation of bedrock stiffness can have more evident influence on the peak frequency of foundation displacement. This is because the soil-tunnel interaction is dominated by kinematic interaction, which can be influenced only by site dynamic characteristics, so the dynamic characteristics of tunnel response are similar to site dynamic characteristics; while the dynamic characteristics of foundation response are also influenced strongly by superstructure dynamic characteristics, and the system mass is large with the inertia interaction also introducing much influence on foundation response, so the dynamic characteristics of foundation response are different much from site dynamic characteristics.
It is also noticed that although the spectra shapes of oblique incident SV-wave differ much from that of oblique incident SH-wave in Fig. 4, the two spectra of vertical incidence is very similar to each other, especially for peak value and peak frequency. This is because the dynamic characteristics of underground tunnel can be influenced only by the site characteristics. As the free-field response for vertically incident SV-wave is exactly identical to that of vertically incident SH-wave, the two spectra of tunnel displacement for vertical incidence are similar to each other; while as the free-field ground motion for obliquely incident SV-wave and SH-wave is essentially different, the two spectra of tunnel displacement for oblique incidence are also different much just as the spectrum of free-field response. Nevertheless, in Liang et al. (2013a, b), even for vertically incident SV-wave and SH-wave, the displacement spectrum of rigid foundation still holds much difference. This is because the foundation displacement spectrum can be influenced by both the site characteristics and the superstructure characteristics. As the dynamic characteristics of superstructure of in-plane direction are different from that of out-of-plane direction, the two spectra of foundation displacement holds little similarities although the site characteristics for vertically incident SVand SH-wave are identical. In conclusion, the mechanism of soil-tunnel interaction which is a rigid system, is totally different from that of soil-foundation-superstructure interaction which is a flexible system.

Analysis in time domain
The similarity of tunnel response to vertically incident SVwave and SH-wave can further be justified in time domain. Figure 6 is the tunnel response in time domain for b Fig. 6 Time history (left) and response spectrum (right) of tunnel acceleration for vertically incident El Centro wave with peak acceleration of 0.1 g in homogeneous half-space ( vertically incident El Centro wave with peak ground acceleration of 0.1 g as SV-wave and SH-wave (Fu et al. 2015). The left part of each sub-figure is the time history of tunnel acceleration with x-axis being the time history by interval 0.02 s and y-axis being the acceleration of 1 g; the right part is the response spectrum of tunnel acceleration with x-axis being the period and y-axis being the maximum acceleration of 1 g. In this section, the outer radius of the tunnel is a = 5 m, the inner radius is b = 4 m, the embedded depth is d = 10 m, and the mass density is q 0 = 2500 kg/m 3 . The vertical incident SV-wave or SHwave all come from depth D 0 = 40 m. For homogeneous half-space (a), the shear-wave velocity is 250 m/s, the mass density is 2000 kg/m 3 , and the damping ratio is 0.05. For layered half-space, the soil-layer thickness is D = 20 m (b), 30 m (c), and 40 m (d), which corresponds to D/a = 4, 6, and 8, respectively. The soil layer is of shear-wave velocity b L = 250 m/s, mass density q L = 2000 kg/m 3 , and damping ratio n L = 0.02; the bedrock is of b R = 500 m/s, q R = 2000 kg/m 3 , and n R = 0.05. It is observed that in time domain, the tunnel responses for vertically incident SV-wave and SH-wave are more similar than that in frequency domain-they are nearly identical.

Conclusions
The spectrum of tunnel impedance function of layered halfspace vibrates around that of homogeneous half-space; the mechanism of dynamic soil-tunnel interaction in layered half-space is different much from that in homogeneous half-space, and the former is larger than the latter. This is because the layered half-space involves the site dynamic characteristics while the homogeneous half-space cannot reflect these characteristics. The mechanism of dynamic soil-tunnel interaction is different much from that of dynamic soil-foundation-superstructure interaction, because the soil-tunnel interaction is dominated by kinematic interaction, so the dynamic characteristics of tunnel response are similar to site dynamic characteristics, while for soil-foundation-superstructure interaction, the foundation response can be influenced by both site dynamic characteristics, and superstructure dynamic characteristics which can be represented by inertia interaction, so the difference of the peak frequency of foundation response to the resonant frequency of the free-field response is much larger than that of tunnel response to the resonant frequency of the free-field response in soil-foundation-superstructure interaction.
For oblique incidence, the tunnel response for in-plane incident waves is completely different from that for incident SH-waves, especially for small incident angle; while the tunnel response for vertically incident SV-wave is very similar to that of vertically incident SH-wave, because the tunnel response is influenced strongly by the site dynamic characteristics which are identical for vertically incident SV-wave and SH-wave, while differ much for oblique incidence.