Influence of the Seismic Wave Velocity of the Damage Zone on Near-Field Ground Motions

In burst-prone deep underground engineering, seismic waves generated from a near-field ground motion event may play a critical role in causing localized rockburst damage. Accurate estimation of near-field ground motions around excavations is important for seismic hazard risk assessment and dynamic rock support design in underground engineering. During the excavation of an underground cavern, stress redistribution in the surrounding rock leads to the formation of damage zones, including the excavation damage zone (EDZ) and excavation fracture zone (EFZ). The poor properties of the rock in the damage zones cause the wave velocities of the rock mass to decrease and the dynamic wave interaction to change, thereby affecting the ground motions around the excavation. This paper studies the near-field ground motion behavior and reveals the control effect of the seismic wave velocity in the damage zones on the near-field ground motions by the aid of the finite fracturing source model (FFSSM). The research results provide a new knowledge of the influence of excavation disturbance on the ground motion distribution around the excavation, and provide new ideas for the seismic hazard risk assessment and prevention in underground engineering.


INTRODUCTION
With the flourishing development of China's economy, infrastructure construction under complex geological conditions has gradually increased, especially deep hydropower engineering and energy mining projects in the western region. The outstanding features of these engineering construction projects are the large burial depth, complex geological conditions, and high in situ stresses. As a result, rock mass seismic hazards induced by highly stressed environments, such as rockbursts, pressure bumping, and mine earthquakes, have become increasingly severe and greatly challenged project construction. Therefore, the accurate estimation of ground motions around excavations is important for seismic hazard risk assessment and dynamic rock support design in underground engineering.
During the excavation of deep underground engineering, due to the newly created stress state exceeding the rock mass strength, a low-damage excavation damage zone (EDZ) and a high-damage excavation fracture zone (EFZ) form in the surrounding rock (Li et al., 2016;Xue et al., 2017). The physical, mechanical and hydraulic properties of the rock mass in the EDZ and EFZ deteriorate, and many ruptures and microcracks form in the rock mass, which cause the wave velocities to decrease and the dynamic wave interactions to change, thereby affecting the ground motions around excavations (Backblom and Martin, 1999;Tsang et al., 2005;Ge, 2008;Chen et al., 2018;Jia et al., 2021). When studying ground motion behavior around excavation under dynamic loading conditions, the damage zone must be taken into consideration.
According to a series of field tests, in-situ tests, laboratory tests, analysis and research carried out in the Underground Research Laboratory in Canada, after the excavation of deep tunnels, when the stress redistribution meets certain conditions, damage zones form in the surrounding rock within a certain range of the tunnel (Read, 1994(Read, , 2004Martin et al., 1997;Daraei and Zare 2018;Wang and Zhao, 2021). Siren et al. (2015) divided the damage zones in more detail according to the excavation method used, including the highly damage zone (HDZ), construction-induced excavation damage zone (EDZCI), stress-induced excavation damage zone (EDZSI) and excavation disturbed zone (EdZ). Li et al. (2016) divided the surrounding rock of an excavated cavern into the EFZ, construction damage zone (CDZ), stress-induced damage zone (SIDZ) and excavation influence zone (EIZ) on the basis of the work of Harrison and Hudson (2000), Eberhardt and Diederichs (2012), and Siren et al. (2015). Among them, the EFZ corresponds to the HDZ and the EIZ corresponds to the EdZ defined by Siren et al. (2015). According to the research results of the above-mentioned scholars, the author divides the damage zones into the EFZ, EDZ and EIZ, as shown in Figure 1. Figure 1 shows the surrounding rock around an excavated cavern can be divided into an EFZ, EDZ and EIZ from the surface of the tunnel to the surrounding rock. As the distance from the excavation boundaries increases, the influence of the surrounding rock excavation disturbance and stress redistribution gradually decreases. The EFZ rock mass is characterized by macrofractures caused by excavation or spalling due to high in situ stress, which is only distributed within a small range from the surface of the tunnel, and the physical, mechanical and hydraulic properties of the rock mass decrease significantly. The EDZ rock not only contains some macrofractures but also includes mesoscale damage not visible to the naked eye, such as microcracks. In addition, the EDZ is distributed more widely around the tunnel, and the properties of the rock mass are slightly reduced. The excavation disturbance is small in the EIZ, so no new cracks are generated in the corresponding rock mass; it is basically in the initial stress state, and the rock mass properties are almost unchanged (Li et al., 2016).
In rock, the wave velocity decreases with the development of cracks and increases with the increase in rock integrity and stress. The surrounding rock of a cave wall is affected by excavation unloading, causing the wave velocity to decrease to different degrees. In general, as the depth from the cave wall increases, the influence of external disturbance on the surrounding rock gradually weakens until it disappears (there is a certain degree of local jump). The quality of the EFZ rock mass is the worst, and the wave velocity of this rock mass has dropped significantly, which is the lowest among the three zones. The rock mass in the EIZ is relatively intact, so the decrease in mechanical properties is small, and the wave velocity in the EIZ rock mass is not much different from that in the original rock mass. Therefore, the highest wave velocity among the three zones is observed in the EIZ. The EDZ is a transition zone between the EFZ and EIZ, and the wave velocity in the EDZ is between those of the EFZ and EIZ. The existence of the damage zones changes not only the physical, mechanical and hydraulic properties of the rock but also the dynamic wave interaction between the excavation-damage structure (including the excavation and the damage zones) and the seismic waves, such as the reflection, refraction and diffraction of the seismic waves. These two factors have led to changes in site effects and affected ground motions around excavations.
When the damage degree of underground excavation is different, it has different characteristics of the damage zones. The properties of the rock mass in the damage zones gradually deteriorate with the increasing damage degree, and the most obvious change is the decrease in the wave velocity of the rock mass mentioned above. Factors affecting ground motion around underground excavation mainly include excavation shape, excavation span, seismic magnitude, wavelength, rock mass characteristics, and supporting conditions (Cai et al., 2012;Zhu et al., 2012;Cai, 2013;Deng et al., 2014;Lei and Omer, 2015;Wang and Cai, 2015). Among them, the wavelength and span have the most significant impact on the ground motion amplification effect (Wang and Cai, 2015). The wave velocity in the damage zones will affect the incident wave wavelength ( λ v/f, where λ, v and f represent the wavelength, wave velocity and frequency, respectively). Therefore, study of the influence of seismic wave velocity in the damage zones on ground motion is of great significance to the seismic risk assessment around the tunnel. Ma et al. (2015) found that a seismic event concentration area is usually no more than three times or even the same as the size of Frontiers in Earth Science | www.frontiersin.org December 2021 | Volume 9 | Article 803522 the seismic source from the excavation. Research on the distribution of mine earthquakes in the Upper-Silesian coal mine in Poland found that most of the mine earthquake events that induced rockburst damage were within 100 m of the location of the rockburst damage. In view of the above research results, this paper focuses on near-field ground motion events. With the help of the finite fracturing source model (FFSSM) proposed by the authors and the Specfem2D based on spectral element method (SEM), the influence of damage zones on the dynamic wave interaction and the impact of seismic wave velocity in the damage zones on the near-field ground motion is studied.

FFSSM Description
In the simulation of far-field ground motion events, the size of the seismic source is much smaller than the distance between the seismic source and the excavation, so the size of the seismic source can be ignored. The size and the dynamic fracturing process of the seismic source have very little influence on the distribution of ground motion around the tunnel. Therefore, the source is generally expressed by a point seismic source model (PSSM, as shown in Figure 2A). However, the near-field seismic source is relatively close to the excavation, and the seismic source size and dynamic fracturing process have an important control effect on the dynamic wave interaction and characteristics of ground motion. Therefore, the near-field seismic source cannot be simplified as a point source model. The authors propose the FFSSM ( Figure 2B) under the premise of considering the size of the near-field seismic source and the dynamic fracturing process . The basic idea of the FFSSM is to divide the main seismic source into a finite number of sub-sources with the same size (L e , which represents the size of the seismic sub-sources). The number of seismic sub-sources (n) is determined based on the scale of the near-field ground motion event. The sub-sources have two characteristics: delayed fracture and bilateral fracture. The delayed fracture mode means that the sub-sources fracture in sequence. The first one that fractures (t 1 0 s) is the initially activated sub-source (SS#1), and then the follow-up seismic subsources (SS#2 − SS#5) are activated in turn ( SS#2 − #5 are activated at t 2 − t 5 , and t 1 ＜t 2 t 3 ＜t 4 t 5 ) when the fracture spreads to follow-up seismic sub-sources. When all the seismic sub-sources are excited and fractured, the dynamic fracturing process of the near-field seismic event is completed. Bilateral fracturing refers to the mode that fracture propagates from the starting point of the fracture proceeds to both sides simultaneously, i.e., the sequence of fracture is that SS#1 fractures first, then SS#2 and SS#3 are excited at t 2 t 3 , and finally fracture propagate to SS#4 and SS#5 at t 4 t 5 . Moreover, the seismic subsource adopts the dynamic corner frequency, and the seismic subsource frequency decreases with the progress of the seismic subsource dynamic fracturing process, that is, the frequency of the seismic sub-source that fractures first is higher, and the frequency of the seismic sub-source that ruptures later is higher.  Figure 2B show the propagation direction of the seismic sub-source fracture.
Frontiers in Earth Science | www.frontiersin.org December 2021 | Volume 9 | Article 803522 The FFSSM includes two sets of calculation processes, forward simulation and back analysis, which have been implanted in the SEM software. This paper focuses on parameter analysis, that is, adjusting parameters to simulate different cases, so as to study the distribution of ground motion around excavation under different conditions and reveal the control effect of seismic wave velocity in the damage zones on the near-field motion. Therefore, this paper follows the forward simulation calculation process of FFSSM ( Figure 3). There are five parameters involved in forward simulation calculation process, including the size of the seismic sub-source (L e ), the time shift (t i ), the number of seismic sub-sources (n), the sub-seismic moment (M e ), and the dynamic corner frequency (f 0 ). For more information about FFSSM, please refer to the author's FFSSM-related paper .

The Validity of the FFSSM
To validate the reliability of the FFSSM, taking the Kiirunavaara mine as an example, the simulation results of FFSMM were compared with the results simulated by Chen et al. (2017) with the three dimension distinct element code (3DEC) and the results obtained through the empirical scaling law proposed by Potvin and Wesseloo (2013). Potvin and Wesseloo (2013) combined the original far-field relationship and near-field saturation, and proposed a scaling law that is applicable to both the far-field and the nearfield without considering site effects, as shown in Eq. (1).
where, M L is the local magnitude of a seismic event; R is the hypocentral distance in meters; R 0 is the source radius estimated as R 0 α · 10 (1/3)(ML+1.5) ; α and C are empirical constants that need to be determined according to site situation, and α 0.53 and C 0.2 ∼ 0.3 are recommended for prefeasibility type studies. In this paper, C takes 0.2 and 0.3, respectively. Frontiers in Earth Science | www.frontiersin.org December 2021 | Volume 9 | Article 803522 To eliminate the influence of the site effect, the numerical simulation did not carry out the excavation ( Figure 4A), and the comparison of the peak particle velocity (PPV) curves obtained by different methods is shown in Figure 4B. The simulation result curve of the FFSSM (curve c2) is roughly consistent with the result curve simulated by Chen et al. (2017) with 3DEC (curve c1), and it is between the scaling law curves where C takes 0.2 and 0.3 (curves c3 and c4), indicating that FFSSM has high reliability and applicability, and can be used to simulate seismic wave propagation and ground motion distribution.
In addition, Figure 4B also shows the limitations of Eq. 1. This empirical scaling law does not consider the dynamic fracturing process and heterogeneous radiation of near-field ground motion events, which leads to deviations in the evaluation of near-field ground motions. Figure 4B shows that the curves obtained by FFSSM (curve c2) and 3DEC (curve c1) are in the high coincidence degree, while the other two curves (curves c3 and c4) obtained by Eq. 1 are significantly different from the curves c1 and c2 in terms of changing trend and gradient. As the distance between the seismic source and the receiver decreases, curves c1 and c2 gradually approach curve c3. It can be predicted that when the distance is smaller (40 m or less), curves c1 and c2 will exceed the range defined by curves c3 and c4, resulting in higher values obtained by the PPV empirical scaling law. This is because the smaller the distance between the tunnel and the seismic source, the more obvious the size effect of the seismic source on the tunnel, and the more it cannot be ignored. Moreover, because the vibration risk assessment is often carried out before tunnel construction, the insufficiency of engineering data and geological information makes the estimation of empirical constants α and C often empirical and subjective, and errors are prone to occur. Therefore, compared to Eq. 1, the FFSSM is more suitable for near-field ground motion simulation.

Numerical Model Description
In the SPECFEM2D models, heterogeneous media with different seismic wave velocities in the EDZ and EFZ are considered to examine the influence of the EDZ and EFZ on near-field ground motions. The models are shown in Figure 5.
The calculation model is based on the Jinping II Hydropower Station ( Figure 5A). The rock mass in the numerical model of the P-SV system is set as Jinping marble with a density of 2,780 kg/ m 3 , a Young's modulus of 56 GPa, a Poisson's ratio of 0.27 (Fan et al., 2021; and a size of 104 × 104 m 2 . The four sides of the numerical model (with a thickness of 2 m) are set as perfect matching layers (PML) to eliminate unrealistic seismic wave reflections, and excavate the center (X 50, Z -50) of the numerical model.
Three different sets (Set A, Set A and Set C) of receivers are arranged in the numerical model to monitor the ground motion around the excavation, as shown in Figure 5B. Set A are placed in a range of 25 × 25 m 2 around the tunnel with an interval of 0.5 m. Set B are placed on the excavation surface, considering the strong ground motion on the excavation surface, the interval is set to 0.25 m. Set C are placed in a range of 100 × 100 m 2 with an interval of 1m, and are arranged around Set A.

Simulation Scheme
The excavation section of the tunnel adopts the horseshoe shape of diversion tunnel #2, as shown in Figure 6A. The section   The influence of the damage zones (including the EDZ and EFZ) on the amplification effect is controlled by the damage degree of the damage zones, and there are two main controlling factors: the range of damage zones and the wave velocity in the damage zones. When the degree of damage increases, the depth of the EDZ and EFZ gradually increases, and the wave velocity in the EDZ and EFZ also gradually decreases. This paper focuses on studying the influence of the seismic wave velocity in the damage zones on the amplification effect. On this basis, the following scheme is set up, as shown in Figure 6C below.
The seismic source is represented as the FFSSM and is located in the lower left corner of the numerical model. There are two reasons for setting the seismic source in the lower left corner: 1) The diagonal distance of the rectangular model is greater than the side length, and the seismic source is set in the lower left corner to maximize the use of the numerical model, in other words, reduce the computational domain range and reduce the computing time; 2) The dip of the seismic source is 45°. When it is located at the lower left corner, the largest wavefront will cover the tunnel area, and the vibration effect is the most obvious. According to previous research results , when the main seismic source moment magnitude M w is less than 2, n can be set to three. Therefore, the FFSSM in the numerical simulation consists of three seismic sub-sources, namely seismic subsources 1, 2 and 3, and seismic sub-sources 2 and 3 are located on both sides of seismic sub-source 1 ( Figure 6C). As mentioned in FFSSM Description, the FFSSM has the bidirectional fracture mode and the delayed fracture mode, i.e., seismic sub-source 1 is activated first and begins to fracture, and after the delay time (t 2 t 3 ), seismic sub-source 2 and seismic sub-source 3 were excited simultaneously. The seismic source is represented as the FFSSM with a strike (ϕ s ) of −90°, a rake (c r ) of −90°and a dip (δ d ) of 45°. The widely accepted Ricker time wavelet was adopted as a time function, and 5 MPa is set as stress drop based on comprehensive consideration of engineering background and related research (Hanks, 1976;Abercrombie, 1995;Ide and Beroza, 2001;Atkinson and Pierre, 2004;Allmann and Shearer, 2009;Hardebeck and Aron, 2009;Wang et al., 2013;Deng et al., 2014).
To study the influence of the seismic wave velocity in the EDZ and EFZ on near-field ground motions around excavation due to a fault slip-induced seismic event, three sets (Set 1 to Set 3) of cases with different seismic sub-source length (L e is 2, 4 and 6 m, respectively) are considered, and each set has four cases (Case 1 to Case 4), as shown in Figure 6C. In the same set, the seismic source parameters of the four simulation cases remain unchanged ( Table 1) but the seismic wave velocity in the EDZ and EFZ changes ( Table 2). In addition to the seismic wave velocity in the EDZ and EFZ, other material properties of the EIZ, EDZ, EFZ and tunnel (such as the density, Young's modulus, and Poisson's ratio of the rock mass) in all cases remain unchanged. Moreover, all the 12 cases have the same EDZ and EFZ range as shown in Figure 6C, which is based on the actual range of the EDZ and EFZ at 11 + 205 in diversion tunnel #2 of the Jinping II Hydropower Station (extremely strong rock burst event, with extremely high damage degree to the surrounding rock around the tunnel).
A non-uniform velocity distribution will lead to high modeling difficulty and complex parameter input, which is difficult to achieve, so we simplified the numerical model. Different seismic wave velocities are assigned to EDZ, EFZ and EIZ, while the seismic wave velocities in the same zone remain unchanged. The rock mass wave velocities in the damage zones are selected based on the long-term monitoring data of Jinping II Hydropower Station. As the damage degree increases, the wave velocity in the surrounding rock decreases. The wave velocities in the EDZ of Case 1 to Case 4 are 0.95, 0.9, 0.8, and 0.7 times of the wave velocity in the EIZ, and the wave velocity in the EFZ is 0.9, 0.8, 0.6, and 0.4 times the wave velocity in the EIZ, respectively.

RESULTS OF THE NUMERICAL SIMULATIONS
In the numerical modeling, propagation of seismic waves excited by the FFSSMs in models with different seismic wave velocity in the EDZ and EFZ, including background models (without an excavation) and excavation models, is conducted first. Then, the PPV distribution around the excavation is obtained based on the monitoring data recorded by the receivers. The receiver can set different seismotypes such as velocity, displacement, acceleration, pressure and fluid potential, and the PPV distribution can be obtained by processing the velocity-time curve of each receiver. Finally, an amplification factor ( α, obtained by the ratio between the PPVs monitored by the receivers in the excavation model and the PPVs in background model) at each receiver is calculated to explore the influence of the seismic wave velocity in the EDZ and EFZ on the near-field ground motions. Figure 7 presents snapshots of the vertical velocity component at 0.0096, 0.016, and 0.0224 s for the four cases with different seismic wave velocity in the EDZ and EFZ (Set 2 of simulation scheme). In the simulation of this paper, because the selected seismic sub-source represents a normal fault with a dip of 45°, the strong ground motion of the S-waves can be observed at 45°/135°t o the X axis in the background model (to facilitate the analysis, the results are drawn in the excavation model, as the white dashed line shows in Figure 7A), and the ground motion at other locations is relatively weak (Stein and Wysession, 2009). The following conclusions can be drawn from Figure 7.
When S-waves propagate to the excavation structure, reflected S-waves will be generated at the upper-left side and the lowerright side of the excavation. Reflected S-waves will form a seismic wave superposition area with both S-waves and diffracted S-waves, where the ground motion is enhanced, and the wave superimposed area will gradually move farther away over time, as shown by the red dotted box in Figures 7B,C. Diffracted S-waves will also be generated when the S-wave propagates to the excavation structure. Because of the change in the propagation direction and the generation of diffraction loss, the ground motion at the upper-right side of the excavation, which is the largest wavefront coverage area in the background model, is weakened, forming a wave weakened area that gradually moved farther away, as shown by the green dotted box in Figures 7B,C.
As the seismic wave velocity in the damage zones decreases, the properties of surrounding rock in different zones differ more obviously. The multiple refraction and reflection effects of the EDZ, EFZ and tunnel cause the seismic wave components to be diversified and complicated (such as the incident waves, reflected waves and refracted waves caused by different interfaces), so the degree of wave field disturbance around the tunnel gradually increases, and the site effect of the excavation-damage structure on near-field ground motion gradually increases. Moreover, changes in the seismic wave velocity in the damage zones also cause changes in the shielding effect of the excavation-damage structure. As the seismic wave velocity in the damage zones decreases, the ground motion on the upper-right side of the excavation weakens, which means that the shielding effect will gradually increase. Figure 8 presents the PPV and amplification factor distributions of four cases with different rock mass wave velocities in the EDZ and EFZ, corresponding to Set 3. As shown in Figure 8A, the reflected and the diffracted S-waves are generated and started to propagate under the dynamic wave interaction, resulting in the appearance of the Zones A and B with high-PPV. The seismic sub-source 1 contacts with the excavation-damage structure first and generates reflected S-waves, and then interacted with the S-waves of the seismic sub-sources 2 and 3 arrived later in the left floor and left sidewall areas of the excavation to form a high-PPV area Zone A. The S-waves superimpose with reflected and diffracted S-waves in the lower-right and upper-left sides of the excavation to form a high-PPV area Zone B. In addition, there is a low-PPV area Zone C ( Figure 8A), which is formed by the shielding effect of the excavation-damage structure and the diffraction loss, and is located in the upper-right side area of the tunnel.
Table 2 and Figure 8A show that with the decrease in the seismic wave velocity in the EDZ and EFZ, the overall trend of the ground motion around excavation gradually increases. From Case 1 to Case 4, the seismic wave velocity gradually decreases, the values of PPV m are 7, 6.4, 7.4, and 10 cm/s, and the PPV m area is initially located in the left floor area, and gradually moved to the left side wall area. In addition, with the decrease in seismic wave velocity in the EDZ and EFZ, several small areas of high-PPV gradually form in the right sidewall and right haunch areas of the excavation, and the value gradually increases, reaching 6.4 cm/s in Case 4.
The effect of the seismic wave velocity on the PPV contours is shown in Figures 8B,C. Due to the site effect of the excavation-damage structure, two types of representative areas are formed: one is the PPV amplification area (A1 and A2) with α ＞ 1, and the other is the PPV weakening area with α ＜ 1, which shows that the site effect of the excavation-damage structure on the ground motion includes both amplification and shielding. With the decrease in the seismic wave velocity in the EDZ and EFZ, the amplification effect of the tunnel excavation-damage structure increases, and the ground motion around the tunnel gradually increases. From Case 1 to Case 4, the seismic wave velocity gradually decreases, the maximum amplification factor α m in PPV amplification area A1 gradually increases from 2.1 to 3.1, and the amplification factor in the right haunch area and right sidewall area of the tunnel gradually increases to 2.3. In contrast, PPV amplification area A2 gradually decreases in both range and intensity, and the maximum amplification factor α m in A2 gradually decreases from 4.3 to 2.4.
Combining Figure 7 and Figure 8, the following conclusions can be drawn: the seismic wave velocity in the EDZ and EFZ greatly influences the site effect of the tunnel. As the seismic wave velocity decreases, on the one hand, the amplification effect in the near-excavation area increases, which is presented as an increase in PPV m and α m ; on the other hand, the amplification effect of the far excavation area decreases, and the range shrinks. Figure 9 presents the maximum PPV in the excavation models and the maximum amplification factor in the nearexcavation area (A1) with different seismic wave velocity in the EDZ and EFZ, corresponding Table 2. There are two reasons for selecting the near-excavation area A1 for analysis: on the one hand, the near-excavation area has a high PPV value, which can be seen from Figure 8, and on the other hand, the rockburst mainly occurs in the shallow rock mass of the tunnel. The following conclusions can be drawn from Figure 9. When the quality of the surrounding rock decreases, as the seismic wave velocity decreases, both PPV m and α m show a gradual increasing trend overall, which shows that the changes in the quality of the surrounding rock masses in different damage zones have an amplification effect on the ground motion in the near-excavation area, and the greater a change in the surrounding rock quality is, the more significant the corresponding amplification effect. Of course, the above conclusion is an overall trend based on comprehensive consideration of all data, with some exceptions, e.g., the PPV m and α m of Case 2 in Set 2 and Set 3 are slightly lower than those of Case 1. There are two reasons for these exceptions: 1) The wave velocity in Case 2 is very close to that in Case 1. The wave velocities in the EDZ and EFZ of Case 2 are 0.95 and 0.9 times those of Case 1, respectively, which results in a small impact on ground motion. 2) The FFSSM is used in the numerical simulation. The seismic source is composed of three sub-sources at different locations, so the seismic response around the tunnel is determined by the composite wave of the three sub-sources at this location and the excavation damage structure. When the seismic wave velocity changes, the wavelength and frequency of the composite waves at different locations will also change, and this change is different from the change of the a single seismic sub-source and more complicated. It can be found from Figure 9 that there are no exceptions to the trend line of Set 1 ( L e 2 m), and the trend lines of Set 2 ( L e 4 m) and Set 3 ( L e 6 m) have exceptions, which shows that as the size of the sub-sources increases, its influence on the ground motion also increases. When the size of the sub-sources is large and the difference between the wave velocity in the damage zones and the EIZ is small, e.g., Case 2 in Set 2 and Set 3, the influence of the sub-source size exceeds the influence of the wave speed change, and some exceptions appear. These special exceptions do not affect the overall trend. When the difference between the wave velocity in the damage zones and the wave velocity in the EIZ is large, the influence of the sub-source size does not play a leading role, as shown in Figure 9.

SUMMARY AND CONCLUSIONS
The excavation unloading of the tunnel forms a high-stress zone, leading to the formation of the damage zones. The quality of the rock mass and the seismic wave velocity in these zones decrease, so the ground motion behavior also changes. This paper studies the near-field ground motion behavior in the damage zones and reveales the influence of the seismic wave velocity in the damage zones on the near-field ground motion by the aid of the finite fracturing source model (FFSSM). The ground motion behavior around the excavation is affected by the seismic wave velocity in the EDZ and EFZ. When the seismic wave velocity in the damage zones changes, the intensity of ground motions around the tunnel also changes. As the seismic wave velocity in the damage zones decreases, PPV m and α m gradually increase.
The numerical results provide additional insights into the ground motion behavior of underground excavations under dynamic load conditions, such as hydropower station construction (powerhouses and tunnels), mining (shafts, roadways, and stopes), and enrich the theoretical research of seismic hazard risk assessment and dynamic rock support design. In addition, the next step of this research is to extend the numerical model to three dimensions, couple software to take the in-situ stress and disturbance stress into consideration, and further analyze the impact of underground excavation damage zones on the risk of rock mass seismic hazards around deep tunnels. Frontiers in Earth Science | www.frontiersin.org December 2021 | Volume 9 | Article 803522