Prediction of Bench Blasting Vibration on Slope and Safety Threshold of Blasting Vibration Velocity to Undercrossing Tunnel

,


Introduction
Strategic petroleum reserves are the most important links in a national energy security system. As of 2015, China has built eight national petroleum reserve bases with a total reserve capacity of 28.6 million cubic meters, including seven surface and one underground oil storage depots. With the undertaking of large-scale construction projects related to these petroleum reserves, large numbers of adjacent petroleum pipeline tunnels passing through open high-steep slopes continue to emerge. During the excavation of the high-steep slope and petroleum pipeline tunnel, blasting, as a fast and efficient excavation method for hard rock masses, has been widely used. Reasonable and proper evaluation of the impact of the blasting vibration and control of the detrimental effects of the vibration are the key technical issues to ensure the safety and stability of the slope and tunnel during blasting excavation, thus realizing the safe and efficient construction of national strategic petroleum reserve projects.
ere are many extensive studies on the impact of rock blasting excavations on open-pit slopes and tunnels [1][2][3]. For example, Jiang et al. established a mathematical model to describe the attenuation of the PPV on open-pit slopes subjected to underground mining blasting [4]. Li et al. proposed a method to predict the time history of blast vibration on high slopes, from which a blast vibration spectral control scheme was presented [5]. Mohammadi Azizabadi et al. modeled the effect of blast vibration on slope stability in jointed rock masses by coupling waveform superposition and numerical methods [6]. Ma et al. proposed an integrated method of microseism energy density and the magnitude-frequency relation based on microseismic monitoring to analyse the slope stability [7]. Huang et al. discussed the time-frequency characteristics and multifrequency band energy distribution characteristics of tunnel blasting vibration signals by using Fourier transform and wavelet packet transform [8]. Lu et al. derived analytical solutions for the particle velocity response of the surrounding rock of a circular tunnel subjected to cylindrical P-waves [9]. Li and Li presented a mathematical method to calculate the relative velocity around a circular tunnel induced by blasting loads and theoretically analyzed the influence of the wavelength-totunnel-diameter ratio on the dynamic response of underground tunnels [10]. Jiang and Zhou proposed an approach to mathematically model the influences of blasting vibration on the tunnel structure and finally provided blasting vibration safety criteria [11]. In recent years, with the development of computer technology, numerical simulation methods are frequently adopted to study the influence of blasting on slopes and tunnels [12][13][14][15][16]. However, construction projects related to the national petroleum reserve bases usually involve cooperative blasting excavation of high-steep slopes and undercrossing tunnels, which is the key constraint on the construction period and engineering safety of the project [17]. e unreasonable schedule of the slope and tunnel construction directly affects the progress of the overall construction project. erefore, it is necessary to conduct an in-depth study on the schedule design of the construction of the upper slope and undercrossing tunnel to guide the overall project construction.
In this paper, we focus on the expansion project of the national petroleum reserve base of Aoshan Island in Zhoushan City (Section 2). First, based on the field monitoring and analysis of the blasting vibration, a mathematical model is established to describe the attenuation law of PPV on the slope surface affected by bench blasting (Section 3 and Section 4). Second, a 3D numerical model is established to analyze the vibration response of the slope subjected to bench blasting by using the dynamic finite element software ANSYS/LS-DYNA. In addition, the reliability of the numerical simulation results is verified using the field monitoring data (Section 5 and Section 6). en, the dynamic response of the undercrossing tunnel subjected to bench blasting was analyzed by numerical simulations, and a PPV safety criterion of the tunnel surrounding rock is proposed (Section 7 and Section 8). Finally, based on these primary results, we calculated key nodes for the cooperative blasting construction of the high-steep slope and undercrossing tunnel and proposed an excavation timing arrangement of the upper slope and undercrossing tunnel (Section 9).

Engineering Background
e Zhoushan National Petroleum Reserve Base Expansion Project is located in the west of Aoshan Island, Zhoushan City, Zhejiang Province. e project is located at 29°56′42″-29°59′00″N in latitude and 122°8′12″-122°9′20″E in longitude. e estimated rock excavation volume of the project is 3.05 million cubic meters, and the designed reserve of petroleum is 2.4 million cubic meters. e project construction involves cooperative blasting excavation of a highsteep rock slope and undercrossing oil pipeline tunnel, as shown in Figure 1. e slope to be excavated is steep in height with a slope angle of 35−45 degrees and gentle at the lower part with a slope angle of 20−30 degrees. e highest excavation elevation is +128 m, and the lowest excavation elevation is +8.0 m. In order to ensure the stability of the undercrossing tunnel surrounding rock, tunneling is stopped when the tunnel working face is 30−50 m away from slope surface. At this time, the blasting excavation elevation of the upper high-steep slope is approximately 40 m. Moreover, there are enlarged tunnel sections with a length of 14 m, at 94 m and 294 m from the designed tunnel entrance, as shown in Figure 2. e surface of the slope body is mostly silty clay, and downwards are strong weathering, moderate weathering, and slightly weathering crystal debris vitreous tuff, with medium hardness in most rock masses. e width of each stage platform is 3 m (Figure 2). e slope is excavated by the bench presplitting blasting method. Table 1 lists the parameters of bench blasting. Figure 3 shows the arrangement of the blast holes and the design of delay detonators. e tunnel is excavated by the drilling-and-blasting method. Figure 4 shows the arrangement of the blast holes and the design of delay detonators, and Table 2 lists the parameters of the blast holes.

Layout of Blasting Vibration Monitoring Points.
To accurately evaluate the dynamic response of the upper slope and undercrossing tunnel, 8 locations are designated as the blasting vibration monitoring points according to the blasting characteristics of the slope excavation, to monitor and analyze the vibration of the blasting construction, as shown in Figure 2. In order to ensure the authenticity of the vibration date, considering the influence of the loosening of the tunnel surrounding rock on the particle vibration, a 2-2.5 m deep hole is drilled at the measuring point and a threaded steel bar with a diameter of 18 mm is inserted and densified by grouting. A steel plate of 10 cm × 10 cm × 1.5 cm in size is welded horizontally on the exposed end of the steel bar, as shown in Figure 5(a). For the measuring point on the slope surface, the ground is smoothed with cement mortar after the surface slag and gravel were removed, and the steel plate with prefabricated threaded holes is wedged horizontally. e vibration velocity transducer is fixed on the steel plate with screws tightened to establish a rigid      Shock and Vibration connection between the transducer and the slope surface, as shown in Figure 5(b).

Vibration Monitoring Results and Analysis.
During the slope excavation, 10 sets of effective field blasting vibration tests on the slope surface and 10 sets of effective field blasting vibration tests on the tunnel surrounding rock were carried out. In order to better study the blasting vibration response of the upper slope and tunnel surrounding rock, the formula for the peak velocity vector sum (PVS) is defined as follows [18]: where V x , V y , and V z are time-history functions of vibration velocity in x, y, and z directions, respectively. PVS comprehensively covers the time history information on the particle vibration velocity in three directions, although the occurrence moment of PVS is not strictly consistent with the moment of PPV in all directions (x, y, and z).
PPV and PVS data are listed in Table 3 and 4. It can be seen from Table 3 and 4 that PPV is the largest in the z direction at most measuring points on the slope surface. However, there is no obvious rule for the maximum PPV direction of the measured points in the tunnel. erefore, when evaluating the influence of the blasting vibration velocity on the tunnel surrounding rock, the maximum PPV should be selected.
Frequency analyses show that approximately 94% of the monitoring data contain vibration frequencies of 25−100 Hz for the upper slope, with main vibration frequencies concentrated in the range of 25-60 Hz and 50−150 Hz for the three directions and in the range of 50-100 Hz for the undercrossing tunnel. e blasting vibration frequency of the slope excavation is relatively high compared with the natural frequency of the upper slope and undercrossing tunnel. erefore, resonance of the blasting vibration with the upper slope and undercrossing tunnel is difficult to achieve.

Attenuation Rule and Prediction Model of PPV
Previous studies have shown that, during the propagation of the blasting seismic waves along the slope surface, the attenuation of seismic waves is affected by many factors such as explosion source, propagation medium, rock properties, and distance from the explosion source [12,19,20]. Multiple PPV prediction models are proposed as follows: Note. e cross-sectional area is 29 m 2 , the total number of blast holes is 86, the total charge is 145 kg, the maximum charge per delay is 18.9 kg, and the specific consumption of the explosives is 1.6 kg/m 3 .

Shock and Vibration
where K is the field influence coefficient, Q is the maximum charge per delay, D is the horizontal distance, H is the vertical distance, α is the PPV attenuation coefficient, β is the influence coefficient of elevation difference, and R is the distance from the explosion source. e relative slope gradient (H/D) affects the propagation of the blasting seismic waves in the rock and soil because of the propagation path of the blasting seismic waves along the slope surface. Table 5 summarizes the main variables involved in the propagation attenuation of the blasting seismic waves on the slope surface. e functional relationship between the variables can be expressed as where variables D, Q, and C P are independent and satisfy the following formulas according to Buckingham's Pi theorem in dimensional analysis [21]: where π n denotes the dimensionless form of the dependent variables, n � 1, 2, . . ., 6, and the exponents α n , β n , and c n are dimensional exponents. π n can be calculated as follows: Substituting (5) into (3) gives For the excavation blasting operation of the same site, E, ρ, μ, and C P of the propagation medium can be regarded as constants. erefore, (6) can be simplified to erefore, the similarity criterion equation of PPV can be written as follows: where α 7 , β 7 , and c 7 are the coefficient related to engineering geological conditions. e frequency of PPV caused by the blasting seismic wave is influenced by the properties of the propagation medium, propagation distance, and mass of the explosives. Assuming that the rock property and mass of the explosives are constant, the frequency can be calculated as follows [22]: where k f is the frequency coefficient, k f � 0.01-0.03, and C s is the shear wave velocity. Substituting (9) into (8) gives where α 8 and β 8 are the PPV attenuation coefficient and c 8 is the influence coefficient of the relative slope gradient.
Regression analyses of the test results listed in Table 3 using (10) lead to the prediction model of PPV and PVS on the slope surface of bench blasting. In addition, to validate the established prediction model, its prediction accuracy is compared with that obtained using (2), based on the fitting coefficients of the fitting curves, as shown in Table 6.
(a) e fitting correlation coefficients of PPV and PVS at the monitoring points on the slope surface obtained using the established mathematical model, (10), are all larger than those obtained using previous classical formulas, (2), indicating that the propagation attenuation of PPV and PVS on the slope surface is more complicated and is affected by the relative slope gradient. e mathematical prediction model established by considering the influence of the relative slope gradient can better describe the propagation attenuation of the bench blasting vibration on the slope surface. (b) e regression analysis results show that the prediction formula of PVS obtained by using the established mathematical prediction model has a high prediction accuracy with a correlation coefficient of 0.953. Additionally, the PVS comprehensively considers the vibration effect in all directions. us, it can better reflect the vibration response characteristics of the PPV of the slope surface. erefore, the following prediction formula can be used to predict the PVS of the slope surface particles when the bench blasting vibration propagates to the slope surface:

Numerical Modeling and Parameter Selection
e dynamic finite element software ANSYS/LS-DYNA is used to establish the numerical model of the cooperative blasting construction of the upper slope and undercrossing tunnel according to the actual engineering of the blasting excavation, as shown in Figure 2. e coordinates are defined as follows: the radial direction of the undercrossing Note. M is the mass, L is the length, and T is the time.   Table 6.

Shock and Vibration
tunnel is the x-axis, the axial direction of the undercrossing tunnel is the y-axis, and the vertical direction is the z-axis. In order to avoid the boundary effect of the model, the lengths along the tunnel axis direction (y direction) and horizontal radial direction (x direction) are set to 250 m and 200 m, respectively, and the vertical direction (z direction) is set to 157 m. e 8-node SOLID164 solid element is used to establish the model, including 307708 units and 323265 nodes. According to the characteristics of the project site, the top surface of the numerical model is considered as a free constrained boundary and the other surfaces are nonreflecting boundaries. Moreover, displacement constraints are applied to the bottom of the model. e numerical model is shown in Figure 6. e diameter of the blasting holes used in the blasting excavation of the upper open-pit slope is 70−115 mm, which is very small compared with the size of the numerical calculation model. erefore, if the numerical calculation model of the blasting hole is established according to the actual blasting parameters and the explosive material model or the equivalent load is applied on the wall of the blasting hole, the number of numerical calculation model elements will be large, which may prevent the calculation process from completing. In this paper, the equivalent approximate blasting load is applied to the slope surface above the undercrossing tunnel with an elevation of +40 m. e numerical calculation is based on the following basic assumptions [23]: (a) e blasting impact load is equivalent to the triangular load (b) e boosting time of the triangular load is 100 μs, and the positive pressure acting time of the triangular load is 600 μs (c) e blasting impact load acts on the slope in the form of uniform vertical pressure e peak value of the equivalent blasting load is calculated through the following formula [24]: where ρ e is the charge density, D is the explosive detonation velocity, k d is the decoupling charge coefficient, and d c are the blast hole diameter and equivalent charge diameter, respectively, and ƞ is the increasing multiples of detonation pressure, ƞ � 8∼11. e peak value of the blasting load is calculated by using equation (12), and the bench blasting parameters is 730 MPa. Figure 7 shows the time-history curve of the blasting load on the slope surface [25,26]. e numerical calculation parameters are selected based on the results of indoor mechanical tests. e surrounding rock in the research area is simplified to be homogenous without considering the influence of cracks and weak planes.
e * MAT_PLASTIC_KINEMATIC material model is used for the surrounding rock [27]. e physical and mechanical parameters of the surrounding rock are listed in Table 7.

Numerical Calculations and Verification
e Rayleigh damping parameters in the numerical calculation are determined by repeated trial calculation and comparison. In this paper, the numerical calculation results of PPV on the slope surface caused by bench blasting are compared with the field measured data, and then, the Rayleigh damping parameters are adjusted manually step by step according to the comparative analysis results until the relative error between the numerical calculation results and the field measured results is less than 5%. Finally, the reasonable Rayleigh damping parameters are determined to be α 0 � 0.4 and β 0 � 0.0003.
To validate the numerical calculation results, 5 field monitoring points are set up on the slope surface, as shown in Figure 2, and the monitoring points in the numerical model are set up at the same locations. e time-history curves of monitoring points in the numerical model are obtained; one typical curve is shown in Figure 8. e numerically simulated and measured PPVs at each monitoring point are listed in Table 8.
From Table 8, it is evident that the simulated PPVs are slightly higher than but similar to those obtained by field monitoring, and the maximum relative error of PPV is Shock and Vibration 9 10.16%. e particle vibration frequencies of 45.25-120.5 Hz at each monitoring point obtained by numerical simulations are slightly higher than those obtained by in situ monitoring points. e main reason for the above phenomenon is that the particle vibration frequencies of the numerical simulations are obtained without considering the attenuation and dissipation of the blasting seismic waves caused by joints and weakening faces in the rock mass [28,29]. In summary, the comparison between numerical simulation and field monitoring data shows that the three-dimensional model and model parameters adequately describe the field data. erefore, a study of the dynamic response and safety effects of the undercrossing tunnel subjected to upper slope excavation is feasible by numerical simulations.

Tunnel Dynamic Response Caused by
Bench Blasting e reliability of the numerical model was verified using the field monitoring data. On this basis, the dynamic response characteristics of the undercrossing tunnel subjected to bench blasting can be analyzed based on numerical calculations. e rock surrounding the undercrossing tunnel in the expansion project is of high grade and good stability. Initial shotcrete and anchor support were carried out for the undercrossing tunnel during blasting construction on the slope. In order to avoid increases in the numerical calculation time due to excessive numbers of numerical model elements, the initial shotcrete and anchor support of the undercrossing tunnel is not modeled separately in the modeling process. erefore, considering this actual situation, the safety control standards of the surrounding rock are used to analyze the influence of blasting vibration on the tunnel.
According to the theory of seismic wave propagation, the surrounding rock and its enlarged section are most affected when the blasting area is closest to the tunnel. erefore, the specific calculation condition in this section is the blasting area which is located on the bench directly above the tunnel, and the distance from the explosive source is 25.8 m.  x direction y direction z direction To evaluate the influence of bench blasting on the surrounding rock of the tunnel, the safety measure stipulated in the initial scheme of the project blasting design is used, i.e., tunnel blasting is suspended when the tunnel working face is 21 m from the slope surface. erefore, four featured points A, B, C, and D on the tunnel section 21 m from the slope surface are selected, respectively corresponding to the vault crown, arch lumbar, midpoint of the side wall, and bottom of the tunnel wall, as shown in Figure 4. e distribution of PPVs and the effective tensile stress (ETS) at each inspection point of the tunnel section are listed in Table 9.
e data in Table 9 show that PPVz and PPVy are the largest at the vault crown and gradually decrease moving downward along the tunnel section. PPVx is small at the vault crown because of the radial symmetry constraint, followed by an increase, and then decreases moving downward along the tunnel section. e attenuation speed of PPVz is evidently faster than those in the other two horizontal directions. It can be seen that the stress wave generated by the blasting load initially propagates downward with a high-strength longitudinal wave and gradually attenuates with an increase in distance. e surface wave formed at the free surface has a slower attenuation speed than that of the longitudinal wave. e maximum PPV of the tunnel section appears at the vault crown, which is 13.92 cm/s. e distribution of ETS at each inspection point shows that ETS gradually decreases along the tunnel section, and the attenuation speed gradually decreases. PPV and ETS both produce maximum values at the vault crown of the tunnel, which are 13.92 cm/s and 4.78 MPa, respectively. erefore, the safety and stability of the tunnel are evaluated according to the dynamic response of the surrounding rock at the vault crown. e maximum ETS at the vault crown reaches 4.78 MPa, exceeding the tensile strength standard, which is the dynamic tensile strength (3.7 MPa) shown in Table 7. erefore, the surrounding area of the tunnel is not safe according to the maximum tensile strength theory.

PPV Safety Criterion of Tunnel Surrounding Rock
To further study the dynamic response of the tunnel, ETS and PPVs at different measuring points are calculated by numerical simulation and listed in Table 10. e numerical simulation results indicate that ETS are different from PPVs. A statistical relationship model between PPV and ETS is established and shown in Figure 9. e linear statistical relationship between PPV and ETS shown in Figure 9 is established as follows: Equation (13) indicates that a linear relationship exists between PPV and ETS [23,30]. Based on the first strength theory, ETS exceeds the dynamic tensile strength (3.7 MPa) of the tunnel surrounding rock when PPV reaches 10.38 cm/s. To be more practical, the PPV safety criterion of the tunnel surrounding rock is determined to be 10 cm/s.

Excavation Timing Arrangement of the Upper Slope and Undercrossing Tunnel
Based on the simulation calculations, the relationship between PPV at the vault crown and the length of retained rock mass during bench blasting directly above the tunnel entrance is obtained, as shown in Figure 10. It can be seen from Figure 10 that PPVz decrease with an increase in the axial distance, and their attenuation speed is faster than that of PPVx and PPVy. e attenuation speed in the excavated direction of the tunnel is less than that in the direction of no excavation. e horizontal vibration velocity (PPVx and PPVy) increases first, then decreases in the excavated direction, and decreases rapidly in the direction of no excavation. PPVs in the three directions increased slightly adjacent to the slope surface and then tended to become stable. e reason for that is mainly due to the coupling effect of the increase in the propagation distance and the decrease in the original rock stress constraint.
To ensure construction progress and the stability of the tunnel surrounding rock, PPV of the tunnel vault crown shall not exceed 10 cm/s and the length of the retained rock mass at the tunnel entrance shall not be less than 36 m; these values can be obtained from Figure 10. erefore, tunnel blasting is suspended when the distance from the tunnel working face to the slope surface is less than 36 m, and all the   working faces are transferred to slope excavation at this time. Subsequently, the tunnel is excavated from the slope surface to the slope body until tunnel hole-through is reached, after the excavation and reinforcement of the slope at the tunnel entrance have been completed. e excavation timing arrangement of the upper slope and undercrossing tunnel is shown in Figure 2. e research results have been successfully applied to the blasting project of the expansion project of the Zhoushan national petroleum reserve base. e cooperative blasting construction of the open-pit slope and undercrossing tunnel has been effectively guided with good results, and the construction period is shortened by 45 days while ensuring construction safety. e project site is shown in Figure 11.

Conclusions
(b) A 3D numerical calculation model is established to analyze the dynamic response of the tunnel surrounding rock subjected to bench blasting. A comparison of the numerical simulation data and the field measurement data shows that the numerical model and the selected parameters adequately describe the field data, and therefore, the numerical calculation model is feasible for studies on the dynamic response of the tunnel surrounding rock to bench blasting vibration. (c) A linear statistical relationship between PPV and ETS for the tunnel surrounding rock is established. Based on the first strength theory, the PPV safety criterion of the tunnel surrounding rock is determined to be 10 cm/s. e minimum safe distance from the tunnel working face to the slope surface, obtained by simulation calculations, is 36 m. (d) e excavation timing arrangement of the upper slope and undercrossing tunnel is proposed, which has been successfully applied to the expansion project of the Zhoushan national petroleum reserve base. e engineering application results show that the construction period is shortened by 45 days while ensuring construction safety. e research results have great guiding significance to similar cooperative blasting excavation engineering for high slope and adjacent tunnel with safety and efficiency.

Data Availability
No data were used to support the findings this study.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.