Abstract
In the tunnel and underground space engineering, the blasting wave will attenuate from shock wave to stress wave to elastic seismic wave in the host rock. Also, the host rock will form crushed zone, fractured zone, and elastic seismic zone under the blasting loading and waves. In this paper, an accurate mathematical dynamic loading model was built. And the crushed zone as well as fractured zone was considered as the blasting vibration source thus deducting the partial energy for cutting host rock. So this complicated dynamic problem of segmented differential blasting was regarded as an equivalent elastic boundary problem by taking advantage of Saint-Venant’s Theorem. At last, a 3D model in finite element software FLAC3D accepted the constitutive parameters, uniformly distributed mutative loading, and the cylindrical attenuation law to predict the velocity curves and effective tensile curves for calculating safety criterion formulas of surrounding rock and tunnel liner after verifying well with the in situ monitoring data.
1. Introduction
With the rapid infrastructure development in China, the improvement of control level on the blasting technique, and the well adaptability in various engineering, the blasting construction are used more and more widely; also, the drill and blast method has been becoming the most primary method in tunneling and underground space construction; in China, there are more than 95% mountain tunnel constructions which adopt the drill and blast method [1]. Also, the blasting is producing huge help on the engineering yet generating many adverse effects on existing structures. Also, many scholars have done a lot of researches on adverse effects [2–4]. Hence, predicting the vibration safety criterion has become an important and effective method for analyzing the adverse effects. Jiang and Zhou used the safety criterion method to analyze tunnel liner structure [5]. Karadogan et al. gave a damage criteria norm for blast-induced ground vibrations in Turkey [6]. Yang et al. studied the safety distance for secondary shotcrete in Jinping-II deep-buried tunnels [7]. Though many works have been done by all the world scholars, the most attempts were based on the monitoring data to summarize the empirical equation [8–18], which is restricted for the real engineering blasting vibration. At the same time, some new methods are used to analyze the blasting problem, such as the soft computing method [19], artificial neural networks method [20], MEMS-based commutation module, [21] and RES-based model [22]. But it should be based on large number of in situ monitoring data, not considering the effect of vibration frequency and duration.
The numerical simulation has become an important prediction method in the engineering, due to its high accuracy and less cost, as well as the fast development computer technique, to make many scholar use such ads ANSYS/LS-DYNA finite element software to solve the blasting vibration problems [23–29]. However, most models accepted the simplified triangle dynamic load, trapezoid dynamic load, or simplified data table which are not accurate enough for the blasting loading.
In blasting, the gas pressure in hole reached 50 to 100 thousand times atmosphere. The high gas pressure could crush the host rock around blast hole to form the crushed zone. And because of higher dynamic compressive strength, the crushed zone will consume large partial energy to make the shock wave attenuate very fast. After propagating a distance, the pressure cannot crush the host rock any more, to make the shock wave attenuate as stress wave. Then the surrounding rock would suffer radial compression to generate a strong radial movement to lead the rock shell expansion to cause the circular tensile stress. Also, the tensile stress tends to be greater than critical dynamic tensile strength which is usually less than 1/10 compression strength of rock masses, to lead to the radial fractures formed. After those radial fractures connected with crushed zone, the blasting gas which looked like a wedge entered the radial fractures to make them enlarge. Hence, the fractured zone would be formed which also could consume some energy to make the vibration wave attenuate as elastic seismic wave whose impact region was the elastic seismic zone.
In this paper, the accurate mathematical model of blasting loading was established on the analysis of the blasting pressure change, the blast hole volume expansion, the fracture development, and the blasting gas motion. All the blasting energy transformed to the kinetic energy of the host rock without considering energy loss, this energy was used to form the crushed zone and fractured zone, and the residual energy was used to produce the vibration in elastic seismic zone to propagate out. So the crushed zone and fractured zone were also regarded as the blasting vibration source by taking advantage of the unified mechanical and continuum damage mechanics thus deducting this partial energy for cutting host rock, and regarding the vibration wave on the elastic seismic zone as the cylindrical elastic seismic wave to attenuate along with radius direction [30].
For the far field vibration analysis, the separate blasting loading of each hole could meet Saint-Venant’s theorem to be regarded as the uniformly distributed loading on its own segmentation cutting zone. Hence the complicated dynamic problem in segmented differential blasting turned into an initial equivalent elastic boundary problem alongside the uniformly distributed and the attenuation law. At last, a 3D model in finite element software FLAC3D accepted the parameters to predict the velocity curve and effective tensile curve for confirming safety criterion after verifying well with the in situ monitoring data.
2. Equivalent Boundary of Elastic Seismic Loading
Based on above analysis, for the far field vibration, after a distance propagating, the blasting wave would attenuate as the elastic seismic wave, and this elastic seismic wave followed a circular attenuation law along the radius direction. Also, the crushed zone and fractured zone were considered as blasting vibration source boundary by taking advantage of continuum damage mechanics theory and unified constitutive theory [30].
2.1. Seismic Elastic Velocity
During the blasting, assuming the host rock as non-compressible linear media, and all the blasting energy would transform to kinetic energy of host rock without any energy loss. So, according to dynamic gradient theory in the semi-infinite medium, the formula of the velocity can be shown as [31]where where is the vibration rounded velocity, is the explosive density, is the blasting energy, is the charging length, is radius of blast hole, and is the propagated distance.
For a single blast hole in semi-infinite medium, the mainly destructional form of the crushed zone with the radius and for the fractured zone with radius , so the formula can be expressed as [32] where is Poisson’s ratio, is the P-wave of rock mass, and are, respectively, the dynamic compressive strength and dynamic tensile strength of rock, is the static compressive strength, is the blasting loading on the blast hole wall, and is the attenuation coefficient of blasting vibration wave.
In above,where is the shear modulus and is the density of rock mass.
Experimental study showed that the radius of crushed zone is 3 to 5 times blast hole radius and 10 to 15 times for fractured zone [33].
2.2. Seismic Elastic Loading on Equivalent Boundary
For the elastic seismic loading, the attenuation functions of shock wave and stress wave with the time and displacement, to make the blasting loading worked on the blast hole wall was shown as [34]where is the blasting loading after attenuating along with distance in radius and is the attenuation coefficient of vibration wave; in above formula, for the shock wave For the underground space and tunneling blasting construction, there are, respectively, the initial cut hole, destruction rock hole, peripheral hole, and bottom hole arranged on the cutting surface, among, the initial cut holes were the inclined hole located in the center of cut surface, and the other holes were vertical hole surrounded the initial cut holes. Hence the blasting vibration must be a multiple and superimposed waves on cutting surface.
(1) The Equivalent Blasting Loading of Initial Cut Zone. According to the above formula, for elastic dynamic loading of single blast hole, for the initial cut hole with ignoring the overlapping influence and energy loss, the equivalent elastic loading can be shown aswhere is dynamic load of single initial cut hole.
Assuming the angle between initial cut hole and cut surface is , so the projection of blast hole on the cut surface is ellipse with long radius and the short radius (the radius of blast hole), as shown in Figure 1.
From the above relationship, we can know . Usually, the initial cut hole was the first blasting segment, so the equivalent elastic dynamic load of initial cut section with considering the crushed zone and fractured zone as dynamic loading boundary was shown in Figure 2.
The equivalent elastic loading of all the initial cut holes by taking advantage of Saint-Venant’s principle was shown as where is equivalent elastic loading of all the initial cut holes, is the number of holes, is the area of a single elliptic fractured zone with taking the 10 to 15 times blast hole’s projection radius, and is the area of whole initial cut section.
(2) The Equivalent Blasting Loading of Other Segments. Except the initial cut holes, there was a free face due to last segmented blasting for other segments. And the free face make this segmented blasting could cut the whole rock mass between this segment and last segment, as shown in Figure 3. So taking , the formula is carried out as below:where is the area of a segmented blasting, is the area of whole cut surface, and is width of th segmented blasting.
Based on the time interval of two segments, the total equivalent elastic loading of other segments was shown aswhere is the number of segments and is the time interval between two segments.
3. Continuum Damage Mechanics
For the host rock of tunnel, there was a damage zone due to internal stress redistribution and fracture extension. This damage made constitutive parameter become lower, which also generated impact on vibration wave propagation. The Canadian scholars measured that the main fracture extension depth is 0.5 m, and the maximum depth could reach 1.0 m with using in situ ultrasonic wave velocity testing technology [35]. So according to continuum damage mechanics theory, the effective stress under damage condition can be expressed as [36]where is the effective stress under damage condition, is effective stress of protolith, is the scalar of damage condition, which can be carried out by Mazar damage model with two scalars being, respectively, tensile damage and compression damage , shown as [37]where and are, respectively, compression and tensile coefficients, which depended on the material properties, for the rock, their value is usually 0.5, and and are, respectively, uniaxial compression and uniaxial tensile critical strain value. As shown below:, are, respectively, the negative and positive strains.
So, the damage scalar is the superposition of compression and tensile:where is the effective strain, and when , , and when , .
Based on above analysis, in the damage host rock, the formula of Young modulus between damage rock and original rock iswhere is the Young modulus under damage rock and is Young modulus under original rock.
4. Equivalent Process of Blasting Loading
After blasting, the dynamic loading on the blast hole wall can make the volume of blast hole enlarged and the fracture expanded; also, the gas pressure and dynamic loading must be reduced with volume enlargement. At last, the explosion gas rapidly overflowed and the applied force decayed to zero with fracture development to connect each other.
Based on blasting mechanism, the process of blasting loading can be divided into four stages.(1)The first stage: the dynamic loading will increase with time till reaching the peak intensity of blasting.(2)The second stage: the blasting pressure will be reduced by the fracture expending, the fillings moving, and the volume increasing before the filling was ejected from the hole.(3)The third stage: the explosive gas erupts quickly from the blast hole to lead to lower pressure after the filling was ejected.(4)The fourth stage: the explosive gas rapidly overflowed and the applied force decayed to zero when fractures develop to connect together.
4.1. The First Stage of Blasting
The dynamic loading will increase with time till reaching the peak intensity of blasting when the detonation gas wave propagated to the bottom section of blast hole after exploding. Many research works show that the initial peak blasting loading must have relationship to detonation wave pressure, and according to the Chapman-Jouguet model, the detonation wave pressure in an explosion can be guided by the widely known equation [34]where is the largest detonation wave pressure, is the density of explosive, is detonation velocity, and is the ratio of the specific heats for the detonation gases; in this formula .
The initial explosion pressure which was the explosion gas acted on the blast hole wall just after detonation is approximately the half of the detonation pressure for the coupled charge: For decoupled charges, the initial explosion pressure also had the relationship to the proportion between the blast hole diameter and the charge diameter, so the formula iswhere is the charge diameter and is the blast hole diameter.
The rising time of loading was shown aswhere is the length of blast hole.
4.2. The Second Stage of Blasting
The blasting pressure will reduce with the fracture expending, the filling moving, and the volume increasing of blast hole before the filling was ejected from the hole. So the volume of detonation gas with time can be expressed aswhere is the blast hole radius with time, is the expansion velocity of blast hole wall, is the width of fracture, and is displacement of filling with time.
According to gas law, gas pressure with volume changed in detonation cavity can be shown as [5] Put the formula (21) into formula (22):where is the initial volume of blast hole, , , , , and are all the explosive material parameter, and is initial energy of explosive.
4.3. The Third Stage of Blasting
The explosive gas erupts quickly out from the blast hole to lead to pressure getting lower after the filling was ejected or no filling blasting; according to gas dynamics theory, the blast hole in this period can be looked at as a bottle structure which follows the next assumptions, (a) because of the volume of blast hole was enlarged in last stage, which make the section of exit is smaller than inside, the shape of blast hole looks like a bottle; (b) because of very short time of reaction and better heat-insulating property of rock, the whole process are assumed as adiabatic process; (c) because of the very slow speed of volume expansion during last stage, the initial velocity of gas is zero at the beginning of this stage, so , as shown in Figure 4.
In Figure 4, , , , , and are, respectively, the initial volume, pressure, density, temperature, and velocity. Also, , , , and are, respectively, the pressure, density, temperature, and velocity at the section of exit.
For this bottle structure, the inner pressure change of blast hole had the relationship with gas erupted statement; when , the exit statement was critical condition, and the velocity of gas erupted was sound velocity; when , the exit statement was subcritical condition, and the pressure at the exit was equal to air pressure ; when , the exit statement was supercritical statement, and there was congested phenomenon at exit to make the velocity of gas also the sound velocity, , where , are, respectively, critical statement of pressure and velocity.
According to isentropic gas formula, to determine critical pressure value , where is the critical temperature, is the specific heat ratio of the gas, and J/(kg·K) is the gas constant, so Usually, the gas pressure in blast hole can reach , so the critical pressure , which was supercritical statement to make the velocity of erupted gas be sound velocity.
Based on above analysis, the gas flow formula at the supercritical statement was From time to , based on the first law of thermodynamics with taking advantage of adiabatic process of gas in blast hole,So, Put the formula (25) into formula (27):So, So, the formula of gas pressure change was
4.4. The Fourth Stage of Blasting
The explosive gas rapidly overflows and the applied force decays to zero fast with fractures development to connect.
4.5. The Discipline of Loading
According the above separated analysis of blasting process, the blasting loading curve was shown as in Figure 5.
In Figure 5, is the beginning of blasting time, so ; the initial dynamic loading ; is the maximum blasting loading time, is the maximum blasting loading, is the time of fillings ejected, is the loading on the blast hole wall at that time, is the time just before fractures connect each other, is the loading on the blast hole wall at that time, is the end time of blasting, and is the gas pressure.
5. Project of Case Study
5.1. Engineering Introduction
The Chen-Chi West Connection Line Highway Tunnel, which was designed as a separate double tunnel each containing two lanes, was chosen as the case to be studied. The lengths of two tunnels are 1075 m and 1185 m, and the distance between them is between 18 m and 20 m. The cavern on the foothill region is 14 m wide and 12 m high. In this tunnel blasting, No. 2 rock emulsion explosive with diameter is 32 mm, length is 200 mm, and weight is 150 g was used. And positive multistage differential blasting method was taken by initial cut holes, destruction rock holes, peripheral holes, and bottom holes on the cutting surface. In working surface the distance between the holes was from 0.5 m to 0.7 m, the depth of the blast hole was 3.0 m, the charge length was 2.5 m without fillings, and the detonating cord across through the whole hole from opening to bottom. The blasting parameters were shown in Table 1, and the form of segments was shown in Figure 6 and interval time was shown in Table 2.
5.2. Constitutive Parameter
The formation lithology rock in the tunnel area is conglomerate, which belongs to the late Jurassic conglomerate, color is purple and structure is gravel with the diameter from 2 mm to 60 mm, and the gravels are cemented together with fillings to form the thick block layer. The compositions of the gravels are andesite, welded tuff, rhyolite, and siliceous rocks, and for the fillings they are debris of rock, fine sand, quartz, feldspar, and some other stable mineral sand. The constitutive parameters are shown in Table 3.
6. Numerical Simulation Analysis
6.1. Simulation Model
Based on the actual engineering geology and structures, a 3D model was built in the finite element software of FLAC3D to accept the constitutive parameters, blasting load, and wave propagation parameters, as shown in Figure 7. The geometrical dimensions are length 228 m ( direction), width 203 m ( direction), and height 159 m ( direction). The finite element mesh consists of 386527 nodes and 372200 elements. And for the blasting zone, the group was between 100 m and 103 m ( direction) after the cutting length 100 m, the initial cut hole zone was meshed as 192 elements (brownish zone) and other worked zone was meshed as 606 elements (orange zone), also, as shown in Figure 7.
6.2. Dynamic Loading in Simulation
For the decoupled charge blasting, every stage pressure and arrival time were calculated by using above formulas. The detail values were shown in Table 4.
In the numerical simulation, it is difficult to determine the dumping of the rock mass directly, so, the field monitoring and calculation is necessary; in this paper, the dumping which was input into model was 0.015 based on the field monitoring and calculation repeatedly.
6.3. Verification of Simulation
In this simulation, the solve time was set to 2 s, and after 1284933 calculation steps, the curves of velocity were predicted from the model compared with 9 times field monitoring data at the same location. The simulation data and monitoring data were shown in Table 5, and the typical velocity curves were shown in Figures 8 and 9.
As shown in Figure 8, the predicted peak velocity was 2.505 cm/s compared with monitored peak velocity 2.51 cm/s which was obtained from Figure 9, so, the error rate is 2%; also, for those two curves, they had similar arrival time from 0.3 s to 0.7 s. From the error rate and similar arrival time, the simulation result in the surrounding rock was accurate enough for analysis. Also, for the other group data in Table 4, the error rates are from 0.5% to 23.1%, which all show that the simulation had the higher accurate results.
6.4. Predicted Results Analysis
After the simulation, the velocity curves and effective tensile stress curves with the distances 10 m, 20 m, 30 m, 50 m, and 80 m from the blasting point in surrounding rock and tunnel liner were obtained. The typical curves were shown in Figures 10–13.
For the surrounding rock, as shown in Figure 10, the peak velocity is 7.099. And the maximum effective tensile stress is 0.4502, which was shown in Figure 11. Also, for the tunnel liner, it is easy to learn that the peak velocity is 2.719 and maximum effective tensile stress is 0.7434 from Figures 12 and 13.
6.5. Conforming with Safety Criterion
The peak velocity and maximum effective stress can be obtained from the simulated model, which were shown in Table 6.
Summarizing data in Table 6, the relationships between the effective tensile stress and peak velocity curves are, respectively, the safety criterion for the surrounding rock and tunnel liner, which were shown in Figures 14 and 15.
As shown in Figure 14, the safety criterion formula of the surrounding rock was In this case study, the rock dynamic tensile strength is 4 MPa, so after calculating, when the PPV in the surrounding rock at the junction of tunnel arch and wall reaches 13.865 cm/s, the effective tensile stress will approach the maximum tensile strength.
As shown in Figure 15, the safety criterion formula was Also, for the tunnel liner, because there was a damage zone caused by excavated stress resilience and before blasts, the damage scale was which was calculated by Mazar damage model based on the field fracture surveyed, so the safety criterion formula of the damage tunnel liner is The damage dynamic tensile strength was 3.2 MPa with the damage scale of , so, when the PPV in the tunnel linear at the junction of tunnel arch and wall reaches 6.213 cm/s, the effective tensile stress will approach the maximum tensile strength.
7. Conclusion
For the far field vibration analysis, based on the unified theory, the crushed zone and fractured zone can be regarded as blasting vibration source. Also, the wave propagation was essentially elastic. It was proved to be very accurate after comparing between the simulated velocity curves and the field monitoring data.
The accurate mathematic loading model was carried out by the analysis of the blasting pressure change, blast hole volume expansion, the fracture development, and the blasting gas motion; after well verification, this dynamic loading curve was more realistic and more accurate to the real loading, than other simplified loadings.
For the far field vibration analysis, the separated blasting loadings of each hole were regarded as uniformly distributed loading on its own segmented zone by taking advantage of Saint-Venant theorem; it also met the requirement of analysis by comparing the monitoring data.
Based on the predicted velocity curves and effective tensile stress curves in surrounding rock and tunnel linear, the safety criterion field can be carried out to analyze the safety and stability of every point in field. And through adopting the continuum damage mechanics theory, the safety criterion of tunnel linear was more corresponding to the real situation.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This work has been supported by the research on the characteristics of rock high strain rate under dynamic load, NSFC (Natural Science Foundation of China, no. 51474016). Also, thanks are to Dr. Ming Cui and Dr. Miao Li for help in gathering the monitoring data.