Dynamic Response of Underground Circular Lining Tunnels Subjected to Incident P Waves

1 Key Laboratory of Transportation Tunnel Engineering, Ministry of Education, School of Civil Engineering, Southwest Jiaotong University, Chengdu, Sichuan 610031, China 2 State Key Laboratory of Geohazard Prevention and Geoenvironment Protection, Chengdu University of Technology, Chengdu, Sichuan 610059, China 3Department of Structural Engineering, University of California-San Diego, La Jolla, CA 92093, USA


Introduction
Large earthquakes (i.e., the Kobe earthquake in Japan, 1995, the Chi-Chi earthquake in Taiwan, 1999, the Wenchuan earthquake of magnitude Ms = 8.0 in China, 2008) had caused a large number of seismic damages of tunnels and underground structures [1].This is inconsistent with the traditional concept [2] that underground structures were less susceptible to seismic damage during earthquakes because the seismic acceleration of underground structures was less than that of structures on the ground surface.The dynamic response of tunnels and underground structures during earthquakes is a process of stress wave propagation, reflection, and interaction.For rock lining tunnels, due to the similar stiffness between the tunnel lining and the surrounding rock medium, it is generally thought that the motion of tunnel lining will follow that of the surrounding rock, which means that the interaction between the surrounding rock and tunnel lining can be ignored.Built on this assumption, the fluctuation field and stress field of the surrounding rock and tunnel lining are solved by the wave equation of seismic wave propagation in an elastic half space rock medium [3].
Theoretical studies on dynamic response and aseismic design of tunnels and underground structures may be categorized into numerical methods and analytical methods [4].The numerical methods include finite difference, finite element, and boundary element methods and the analytical methods are mainly referred to as wave function expansion method, which are based on the theoretical analysis of elastic wave scattering and dynamic stress concentration.Although numerical methods can be used for arbitrary-shaped tunnels, analytical solutions are still valuable for providing insight into the formation mechanism of the problems and for checking accuracies of numerical methods.

Mathematical Problems in Engineering
In the 1970s, using the wave function expansion method, Pao and Mow [5] initiatively studied dynamic stress concentration of single cavity in the whole space for incident elastic plane P waves.Later, Lee et al. [6][7][8][9] extended the solution to half space based on the assumption of large circular arc and provided the analytical solutions considering single cavity in half space for incident plane P and SV waves, which is much more complicated due to the wave mode conversion.Manoogian and Lee [10] used the weighted residual method to study the problem of scattering and diffraction of plane SH-waves by including a different elastic medium of arbitrary shape in a half plane.Davis et al. [11] studied the transversal response of underground cylindrical cavities to incident SV waves and derived analytical solutions to evaluate the dynamic response of a flexible buried pipe during the Northridge earthquake.Just recently, Liang et al. [12][13][14], Ji et al. [15], and You and Liang [16] investigated the dynamic stress concentration of a cylindrical lined cavity in an elastic half space for incident plane P and SV waves and derived the series solution to study the amplification of ground surface motion due to underground group cavities for incident plane P waves.Kouretzis et al. [17] employed the 3-D shell theory in order to derive analytical expressions for the distribution along the cross section of axial, hoop, and shear strains for long cylindrical underground structures (buried pipelines and tunnels) subjected to seismic shear wave excitation.Esmaeili et al. [18] used hybrid boundary and finite element method (FEM) to study the dynamic response of lined circular tunnel subjected to plane P and SV harmonic seismic waves.Xu et al. [19] used Fourier-Bessel series expansion method to deduce a series solution of wave equation for dynamic response of underground circular lining tunnels and approximately studied dynamic stress concentration of lining with some influence factors.
Based on the assumption of a large circular arc, a series solution of dynamic stress is deduced by using a wave function expansion approach for a circular lining tunnel in an elastic half space rock medium subjected to incident plane P waves.Compared to previous studies, the boundary conditions of displacement and stress of circular lining tunnel are updated to adapt different surrounding rock conditions, computational parameters are reconfigured, and some cases are analyzed again by using MATLAB program to obtain more precise results.This paper focuses on the effect of the dynamic stress response on circular lining rock tunnels with different factors, including incident frequency, incident angle, buried depth, rock conditions, and lining stiffness.The aim and novelty in this paper mainly include that (1) the effect of low-frequency and high-frequency contents of seismic waves on the lining is estimated; (2) the most unfavorable seismic wave incident angle to tunnels is obtained by comparing the dynamic stress response of the lining subjected to seismic waves with different incident angle; (3) the effect of buried depth on the dynamic stress response of the lining in different surrounding rock conditions is investigated and the depth boundary where the dynamic stress response of the lining is significantly reduced is put forward; (4) the effect of the elastic modulus ratio between the lining and the surrounding rock is estimated to obtain the optimal elastic modulus ratio.These results may provide a good reference for studies on the dynamic response and aseismic design of tunnels and underground structures.

Model for Series Solution of Wave Equation
The schematic model is shown in Figure 1, consisting of half space rock medium and a tunnel lining, with buried depth, internal diameter, and external diameter of the tunnel lining being denoted by ,  1 , and  2 , respectively.The half plane rock medium and tunnel lining are assumed to be elastic, isotropic, and homogeneous.The material properties are characterized by the Lame constants   ,   ,   , and   , and the longitudinal and transverse wave velocities for elastic half space medium and tunnel lining are denoted by   ,   ,   , and   , where subscript  indicates half plane medium,  indicates tunnel lining,  indicates longitudinal wave (i.e., compressive wave), and  indicates transverse wave (i.e., shear wave).The large circular model is used here to simulate the surface of the half space (Figure 1).It was shown that when the circular becomes larger and larger, the solution converges to the exact solution [8,9].

Analyzing Wave Function of Free-Field in Half Space.
With the circular frequency of  and the incident angle of   , the wave potential functions of the incident plane P waves and reflected plane P and SV waves can be expressed in Fourier-Bessel series in the coordinates of  1 as [12] where   = /  is the longitudinal wave number,   = /  is the transverse wave number, and   () is the first type of Bessel function where the time factor exp(−) has been omitted.The coefficients in formula (1) are defined as follows: where  is the imaginary unit;   = 1 for  = 0 and   = 2 for  ≥ 1. Symbol   is the reflection angle of reflected SV waves, which satisfies the equation sin   /  = sin   /  .The reflection coefficients  1 and  2 are given by In O 2 coordinate system, the wave potential functions of the incident plane P waves and reflected plane P and SV waves are similar to O 1 coordinate system but are not presented here.

Analyzing Scattering Field at the Interface between Tunnel
Lining and the Half Space.The presence of the interface between the circular tunnel lining and the half space leads to scattering P waves  11 ( 1 ,  1 ) and SV waves  11 ( 1 ,  1 ), which can be expressed as series forms as below [15]: (1)  11, cos  1 +  (1)  11, sin  1 ) , 12, cos  1 +  (1)  12, sin  1 ) , where   () and  (1)   () are the first class Bessel function and the Hankel function, respectively.The wave numbers of P waves and SV waves of tunnel lining are  1 = / 1 and  1 = / 1 , respectively.In the lining medium, the expressions of P waves and SV waves are complete.Through the coordinate transformation,  2 ( 2 ,  2 ) and  2 ( 2 ,  2 ) can be converted into the Bessel series form in the coordinate system of  1 −  1 and   (   1 ) and  (1)   (   1 ) in the coordinate system  1 −  1 are a set of basis functions of the wave equation.
Similarly, for the presence of circular tunnel lining, there are scattering P waves  1 ( 1 ,  1 ) and SV waves  1 ( 1 ,  1 ), which are caused at the interface between the tunnel lining and the half space.Built on the assumption of large circular arc, the scattering P waves and SV waves are denoted as  12 ( 1 ,  1 ) and  2 ( 2 ,  2 ).They can be expressed as series forms: (1)   (   1 ) ( (1)  1, cos  1 +  (1)  1, sin  1 ) , where   = /  and   = /  are the wave numbers of P waves and SV waves in the half space.From a physical point of view,  1 ( 1 ,  1 ) and  1 ( 1 ,  1 ) are represented by the outward propagation waves from  1 point, while  2 ( 2 ,  2 ) and  2 ( 2 ,  2 ) are generated by the surface of the half space in the large circular arc.Therefore, the expressions of P waves and SV waves are complete in the half space medium.
The wave potential function in the lining medium can be expressed as The wave potential function in the half space medium can be written as

Solving Series Solution after Imposing Boundary Conditions.
Under a plane strain assumption, there are two families of boundary condition: (1) the stress at the surface of the half space and the tunnel lining inner surface is zero; (2) the stress and displacement are continuous at the interface of the tunnel lining and the half space.So the boundary conditions of the stress and displacement of the tunnel lining are summarized as Based on the boundary conditions in ( 8), the coefficients in ( 4) and ( 5) can be obtained.The wave potential functions  1 and  1 , in the tunnel lining, and   and   , in the half space, can be solved by ( 6) and (7).The stress and displacement expressions under a plane strain assumption for incident plane P waves can be described by After imposing the boundary conditions of    =    = 0 and  1  =  1  = 0, (10) and ( 11) can be obtained by combining (1), (2), and (9) [5]: (1)  11 (,  1 )  1(1) 12 (,  1 ) (1)   11, (1)   11, (1)   12,  (1)   12, (1)   11, (1)   11, (1)   12,  (1)   12, 1, where When  = 1,   () is a function of   (); when  = 2,   () is the function of  (1)   ().For the half plane, referred to as ", "   and   can be replaced by   and   , respectively.For the tunnel lining, referred to as "1, "   and   can be replaced by  1 and  1 , respectively.Potential functions  (+) ,  () ,  () , and  (+) satisfy the stress boundary in the half space.

11,𝑛
(1)   11, (1)   12,  lining stiffness  is represented by the shear wave velocity ratio of the lining and the rock half space; (5) the Poisson's ratio of rock medium and lining is assumed to be 0.25; (6) the longitudinal wave velocity of these two medium is 1.732 times of the shear wave velocity; (7) the radius of the large circular arc is  =  12 + , and the distance from the centers O 2 to O 1 is taken as  12 = 10 4  1 in the example, which is large enough to simulate the half space to ensure that the solution converges to the exact solution.

Effect of Incident Frequency on the Dynamic Stress Concentration
Coefficient for the Tunnel Lining.The dimensionless frequency  is defined as the ratio of inner diameter of the tunnel and the incident wave length [15]: where   is the shear wave length of the half space rock.
Considering the effect of different incident frequency contents of P waves on dynamic response of tunnel lining, three typical cases ( = 0.5, 1.0, and 2.0) are carried out for incidence P waves from vertical upward direction.
Figure 2 shows how the toroidal dynamic stress concentration coefficient for the tunnel lining ( *  ) changes with different incident frequency of P waves for two buried depths ℎ 1 = 20 m and ℎ 1 = 100 m.It is shown that, with increasing incident frequency,  *  decreases gradually and  *  at lowfrequency P waves is about 10 times larger than that at high-frequency.Therefore, low-frequency seismic waves are more detrimental to the tunnel lining than high-frequency.Also the increase of incident frequency leads to increasingly complex distribution of the toroidal dynamic stress in the lining, and the less principal stress concentration direction at low-frequency conditions is gradually transformed into more stress concentration directions at high-frequency conditions.

Effect of Incident Angle on the Dynamic Stress Concentration Coefficient of Lining.
The stress concentration coefficients of tunnels and underground structures are very sensitive to incident angle of seismic waves.Small changes of the incident angle may lead to a significant change of the stress and deformation at various locations of the tunnel structures.Herein, studies of three representative cases (i.e., vertical incidence with incident angle  = 0 ∘ , a small incident angle  = 30 ∘ , and large incident angle  = 60 ∘ ) were carried out for frequency coefficients of  = 1 and  = 2. Figure 3 displays the distribution of toroidal dynamic stress of in the lining for three different incident angles.
It is observed that the coefficient of toroidal dynamic stress concentration decreases gradually with increasing incident angle, but the change in the amplitude is small in both dimensionless frequency values,  = 1 and  = 2.It is worth noting that, for vertically incident P waves or a small incident angle ( = 0 ∘ ∼30 ∘ ), the coefficients of dynamic stress concentration of lining are very complex and uneven.The conclusion here agrees well with the results of shaking table tests [20].

Effect of Buried Depth on the Dynamic Stress Concentration Coefficient of Lining.
To demonstrate the effect of the buried depth and the material properties, two different types of rock mass (soft and hard rock mass) are chosen to study the seismic response of the tunnel lining within a buried depth range from 10 m to 150 m.
Case 1 (soft rock as surrounding rock medium of tunnel with elastic modulus   = 2 GPa). Figure 4 shows the distribution  of toroidal dynamic stress in the lining with different buried depth in a range of 10 m-150 m.With the increase of the buried depths, the dynamic stress concentration coefficient of the lining does not converge to the situation of single round hole in whole space.The distributions of lining toroidal dynamic stress with different values of buried depth are very complex and uneven and no directionality can be obviously observed.
As Figure 5 shows, when the buried depth in the soft surrounding rock is within 100 m, the maximum dynamic stress response of the lining varies with buried depth, and two peak values can be observed at the buried depth of ℎ 1 = 30m and 100 m.When the buried depth is greater than 100 m (ℎ 1 > 100 m), the coefficient of dynamic stress concentration of lining does not vary significantly with a high level of about 6.5.It appears that, with an increasing buried depth, the decreasing trend of dynamic stress response is not pronounced in the soft surrounding rock and the dynamic stress of tunnel lining may still be high in the case of a large buried depth.There are three main reasons: (1) the difference in stiffness is large between tunnel lining and soft surrounding rock; (2) the degree of stress concentration of tunnel lining is high and dynamic stress of the lining does not decrease obviously with increase of the buried depth; and (3) the pressure of surrounding rock increases gradually with increase of the buried depth.In this case, the soft surrounding rock has a low strength, large porosity, many joints and fractures, and low capability of self-arching and these rock properties lead to large tunnel dynamic response.
Case 2 (hard surrounding rock medium with elastic modulus   = 20GPa).Regarding the case with hard surrounding rock medium, no obvious directionality can be observed in terms of the distributions of toroidal dynamic stress caused by P waves, but the dynamic stress response of lining is more pronounced than the case of soft rock as illustrated in Figure 6.
Figure 7 shows that when the tunnel buried depth is less than 100 m, the coefficient of maximum dynamic stress concentration in the tunnel lining is generally high, ranging from 6.0 to 8.5 which indicates that the response of dynamic stress of the lining is large at the portal segment of hard rock tunnel and the tunnel lining may be damaged easily during a large earthquake.When the tunnel buried depth is greater than 75 m, the coefficient of dynamic stress concentration is reduced significantly.The coefficient of dynamic stress concentration starts to reach a constant level after the buried depth is greater than 100 m.It appears that the lining dynamic stress is small and tends to stabilize at the large buried depth segment of hard rock tunnel and the tunnel lining may be less susceptible to damage during a large earthquake.Therefore, it is concluded that, under the hard rock conditions, the tunnel buried depth (i.e., thickness of the overlying rock) is not a major factor that affects the dynamic stress and seismic damages of tunnel lining when the tunnel buried depth is greater than 100 m.This conclusion is consistent with the findings of Gao et al. [21], which states that the dynamic stress concentration coefficient tends to be constant when the buried depth is greater than 25 1 and the buried depth of 25 1 approaches 100 m.

Effect of Elastic Modulus on the Dynamic Stress Concentration Coefficient of Lining.
The elastic modulus ( 1 ) of the lining can be considered an index for the elastic deformation of tunnel lining, and the stresses caused by elastic deformation increase with increasing elastic modulus.This also means that, at a given level of stresses, the lining deformation decreases with increasing lining stiffness.We focus on the effect of lining elastic modulus on the dynamic stress concentration coefficient of the lining in a soft surrounding rock.Figure 8 shows the coefficients of dynamic stress concentration ( *  ) with different lining elastic modulus in a soft surrounding rock (  = 2 GPa).
As shown in Figure 8, the coefficient of dynamic stress concentration of the lining increases gradually with increasing lining elastic modulus, and the greater the elastic modulus is, the more uneven the distribution of the dynamic stress coefficient is.The different adaptabilities of lining with different stiffness to dynamic stresses may lead to the variation of the dynamic stress of lining.Herein the following example is adopted to illustrate this aspect.When the elastic modulus of the lining is  1 = 10 Gpa and soft surrounding rock is   = 2 Gpa (i.e.,  1 = 5  ), the maximum value of dynamic stress concentration coefficient of the lining is more than 100, and the variation of dynamic stress of lining across all directions is very complex.Therefore, under the premise of meeting the bearing capacity and deformation of the lining, appropriately soft lining is recommended to be adopted in the tunnel structure to mitigate the seismic damage of tunnel lining and the ratio between the elastic modulus of the lining and soft surrounding rock should be less than 5.

Conclusions
A series solution of wave equation for dynamic response of underground circular lining tunnels subjected to incident plane P waves is presented by Fourier-Bessel series expansion method in this paper.The effects on circular lining rock tunnels of different factors, including incident frequency, incident angle, buried depth, rock conditions, and lining stiffness, are studied through the proposed series solution.In this study, the following conclusions are drawn.
(1) In circular lining rock tunnels, the coefficient of toroidal dynamic stress concentration and dynamic stress decrease gradually with the increase of incident frequency of P waves. *  for low-frequency P waves is about 10 times larger than that for high-frequency, indicating that low-frequency contents of seismic waves are more detrimental to the tunnel linings than high-frequency contents.
(2) The incident angle of P waves has a significant influence on dynamic response of the tunnel lining.When seismic waves propagate vertically or with a small incident angle ( = 0 ∘ ∼30 ∘ ), the distribution of the coefficients of dynamic stress concentration in the tunnel lining is very complex and uneven, which is detrimental to the safety of the tunnel lining.
(3) When the buried depth of a tunnel in the soft surrounding rock is in a range of 0∼100 m, the decreasing tendency of dynamic stress response of the lining with increasing tunnel buried depth is not obvious; even if the buried depth is greater than 100 m, the coefficient of dynamic stress concentration of the lining still remains at a high level of about 6.5.In this situation, tunnels crossing the fault zone or the high stress segment may get damaged seriously during large earthquake.(4) When the tunnel buried depth in hard surrounding rock is less than 100 m, the coefficient of dynamic stress concentration in the tunnel lining is high, ranging from 6.0 to 8.5.While the buried depth is greater than 100 m, the dynamic stress of the lining is small and tends to stabilize.When the tunnel buried depth is more than 100 m, the tunnel buried depth is not a major factor that affects the dynamic stress.(5) The coefficient of dynamic stress concentration in the lining increases gradually with increasing elastic modulus of the lining, and the greater the elastic modulus is, the more uneven the dynamic stress coefficient is.Therefore, under the premise of meeting the bearing capacity and deformation of lining, appropriate soft lining is recommended to be adopted in the tunnel structures to mitigate seismic damage of tunnel lining, and the ratio between the elastic modulus of the lining and soft surrounding rock should be less than 5.

Figure 1 :
Figure 1: The model of half space rock medium and tunnel lining.

Figure 2 :
Figure 2: Dynamic stress concentration coefficient of tunnel lining for different incident frequency and buried depth.

Figure 3 :
Figure 3: Dynamic stress concentration coefficient of tunnel lining with different incident angles.

Figure 4 :Figure 5 :
Figure 4: Dynamic stress concentration coefficient of tunnel lining with different values of buried depth in the soft surrounding rock.

Figure 6 :Figure 7 :
Figure 6: Dynamic stress concentration coefficient of tunnel lining with different buried depths in hard surrounding rock.

Figure 8 :
Figure 8: Dynamic stress concentration coefficient of tunnel lining with different lining elastic modulus in the soft surrounding rock.