Tunnel earthquake three-dimensional full-waveform inversion method

Abstract

The invention relates to the technical field of tunnel engineering, and discloses a tunnel earthquake three-dimensional full-waveform inversion method, which comprises the following steps: S1, establishing a three-dimensional initial velocity model: establishing the three-dimensional initial velocity model based on a common reflection point migration method of ellipsoid expansion; s2, forward modeling three-dimensional full-space wave field characteristic analysis: performing forward modeling in the established typical unfavorable geologic body model, setting a seismic source point on a tunnel wall or a tunnel face, setting a detection point on the tunnel wall, and analyzing three-dimensional full-space wave field signal characteristics received at the positions of the tunnel wall and the tunnel face; s3, full-waveform inversion numerical simulation: realizing a full-waveform inversion numerical model experiment on the basis of the established initial velocity model to obtain simulation data; s5, performing three-dimensional full-waveform inversion imaging; and S6, result verification. The tunnel earthquake three-dimensional full-waveform inversion precision can be improved, and the tunnel earthquake three-dimensional advanced prediction method which is low in cost, safe, efficient and high in precision is achieved.

Description

Technical Field

[0001] The invention relates to the technical field of tunnel engineering, and in particular to a tunnel seismic three-dimensional full waveform inversion method. Background Art

[0002] The construction of tunnels has shown significant advantages in solving terrain restrictions, improving road alignment and shortening travel distances, and has become a key measure to improve the quality of road transportation and service quality. Tunnel projects are buried deep underground and face complex and changeable engineering geological and hydrogeological conditions. During the construction process, unfavorable geological bodies such as fault fracture zones, karst development zones, and water-rich zones are often encountered. Due to the unclear geological conditions in front of the face, if effective measures are not taken in time, it is very easy to cause large deformations, landslides, water and mud gushing, and other disasters, which may even cause major casualties and property losses in severe cases. This situation is particularly prominent in areas with complex geological conditions such as high mountain canyons with intense tectonic movements and developed surface water systems.

[0003] In order to ensure the safety of the project and personnel, ensure the smooth construction of the tunnel, and reduce economic losses, advanced geological forecasting is particularly important. Through advanced geological forecasting, we can accurately understand the changes in the surrounding rock structure and geological disasters in front of the face, so as to timely and reasonably modify the construction plan, arrange protective measures, and avoid or reduce the harm of dangerous situations. Therefore, how to accurately and timely predict the adverse geological factors in front of the face has become the key to determining whether the tunnel construction is smooth or not.

[0004] Elastic wave detection is one of the important means of tunnel geological advance detection, which can accurately predict the spatial distribution characteristics of unfavorable geological bodies within a certain range in front of the face. Although predecessors have made breakthrough progress in tunnel accelerometers, tunnel seismic sources and data processing, the accuracy of geological prediction still needs to be improved due to the low velocity accuracy of geological bodies in the advance direction. As the highest accuracy velocity analysis method at present, three-dimensional elastic wave full waveform inversion faces the following difficulties: the existing tunnel earthquake prediction methods lack data on the positions of the face and both sides of the tunnel wall, and the seismic sources mostly use explosives and mechanical impact seismic sources. The acquisition cost of explosive seismic sources is high, but the data acquisition efficiency is low and the data volume is small. The prediction results are mainly two-dimensional or pseudo-three-dimensional, and true three-dimensional detection imaging cannot be achieved.

[0005] The explosive source is excited in the borehole to avoid the influence of the excavation damage zone, and the data signal-to-noise ratio is high. The mechanical impact source is excited on the surface of the tunnel wall, and the data quality is low due to the influence of the excavation damage zone. A small amount of high signal-to-noise ratio data reduces the accuracy of full waveform inversion due to the small amount of data, and a large amount of low signal-to-noise ratio data also reduces the accuracy of full waveform inversion due to the low signal-to-noise ratio.

[0006] The accuracy of the initial model is a key parameter of full waveform inversion, which directly affects the calculation stability, number of iterations and calculation accuracy. A high-precision initial model can make the inversion calculation process converge quickly and save calculation time. The existing initial velocity model established by direct wave fitting and tomographic inversion method has low accuracy and is quite different from the actual surrounding rock velocity, making the full waveform inversion result unstable.

[0007] The three problems of small data volume, low data signal-to-noise ratio and low initial velocity model accuracy seriously restrict the accuracy of tunnel seismic full waveform inversion, resulting in limited prediction accuracy of poor geological bodies. Summary of the invention

[0008] In view of the above-mentioned deficiencies in the prior art, a three-dimensional full waveform inversion method for tunnel seismic is provided by the present invention, which solves the problems of high cost of existing tunnel seismic detection, low efficiency of explosive source excitation, low signal-to-noise ratio of mechanical impact source data and small amount of data.

[0009] In order to achieve the above-mentioned purpose, the technical solution adopted by the present invention is: a tunnel seismic three-dimensional full waveform inversion method, comprising the following steps:

[0010] S1, establishing a three-dimensional initial velocity model: based on the common reflection point migration method of ellipsoid expansion, three-dimensional migration imaging is performed to obtain a velocity spectrum, the velocity corresponding to the maximum amplitude is selected as the optimal velocity, and a three-dimensional initial velocity model is established;

[0011] S2, forward simulation wave field characteristics analysis: establish a typical bad geological body model, the typical bad geological body includes three-dimensional karst caves and fault fracture zones; conduct forward simulation in the established typical bad geological body model, set source points on the tunnel wall or the tunnel face, set detection points on the tunnel wall, and analyze the characteristics of seismic wave signals received at the tunnel wall and the tunnel face;

[0012] S3, full waveform inversion numerical simulation: using the seismic wave signal data received in S2 and based on the three-dimensional initial velocity model in S1, a full waveform inversion numerical model experiment is performed to obtain simulation data;

[0013] S4, deployment of observation system and data acquisition and processing: a symmetrical joint observation system is deployed in the actual tunnel, and the data of the explosive source and the mechanical impact source are jointly processed to obtain the actual data;

[0014] S5, 3D full waveform inversion imaging: by comparing the simulation data obtained in S3 with the actual data obtained in S4, the residual is calculated, the gradient is calculated using the residual, and the gradient is used to update the parameters of the 3D initial velocity model to obtain a reliable 3D velocity model, realize the 3D full waveform inversion of the tunnel, and accurately calculate the physical quantities of the surrounding rock, including the longitudinal wave velocity, the shear wave velocity, the Poisson's ratio, the Young's modulus and the volume modulus;

[0015] S6, result verification: reversely verify the accuracy of the calculation results, obtain the actual correction coefficients, and finally determine a reliable three-dimensional initial velocity model.

[0016] Furthermore, in the above tunnel seismic 3D full waveform inversion method, the forward simulation in S2 is the process of solving the wave equation. For 3D multi-component, the wave equation is expressed as formula (1):

[0017]

[0018] Among them, x, y, z represent the coordinates of underground space, m is the longitudinal wave velocity of the formation, It represents the gradient of the objective function with respect to the parameter m during the kth iteration, α k is the update step size.

[0019] Furthermore, in the above tunnel seismic 3D full waveform inversion method, the cross-correlation calculation formula of the forward extension wave field from the shot point and the reverse time-delay extension wave field from the detection point is as follows:

[0020]

[0021] Among them, S(x,y,z,t) and R(x,y,z,t) are the forward extended wave field and the reverse time-delayed wave field respectively.

[0022] Furthermore, in the above-mentioned tunnel seismic three-dimensional full waveform inversion method, the deployment of the observation system and data acquisition processing in S4 are specifically as follows: detectors and source points are deployed on both sides of the tunnel wall, and the source points include explosive sources and mechanical impact sources; the two sources are stimulated, and elastic wave data at different positions in the tunnel are obtained through the detectors; the longitudinal and transverse wave signals in the elastic wave data are extracted, and the surrounding rock velocity is calculated using the direct wave velocity fitting.

[0023] Furthermore, in the above-mentioned tunnel earthquake 3D full waveform inversion method, the combined processing of explosive source and mechanical impact source adopts the empirical Green function calculation method, which is expressed as formula (3):

[0024]

[0025] Among them, x belongs to the variable of spatial range A, x∈A, r0 is the distance in space A, u(x,t) is the observation record at point x, and v(x+r0,t) is the observation record at point x+r0.

[0026] Furthermore, in the above-mentioned tunnel seismic three-dimensional full waveform inversion method, the simulation data in S3 includes rock physical quantities, and the physical quantities include longitudinal wave velocity, shear wave velocity, Poisson's ratio, Young's modulus and volume modulus.

[0027] The beneficial effects of the present invention are as follows: the dual-source joint symmetrical observation system of the present invention takes into account the influence of the tunnel excavation damage zone on the data signal-to-noise ratio and resolution, selects the optimal excitation point and detection point, effectively increases the amount of data on both sides of the tunnel wall, and provides more abundant basic information for subsequent data processing and interpretation. The dual-source data joint processing method improves the signal-to-noise ratio of a large amount of mechanical impact source data with low signal-to-noise ratio by constraining a small amount of source data with high signal-to-noise ratio, improves the overall quality of the data, and provides more accurate data support for subsequent three-dimensional full waveform inversion. The three-dimensional initial velocity modeling method based on the common reflection point offset of the ellipsoid expansion sets the velocity scanning range, and selects the velocity corresponding to the maximum amplitude value in the three-dimensional offset scanning diagram as the optimal velocity of the spatial point. Combined with the calculation of massive data, the relationship between the velocity corresponding to the maximum amplitude energy in the three-dimensional offset data and the true velocity of the surrounding rock is established, and the empirical correction coefficients of different lithologies are summarized. Based on these correction coefficients, a reliable three-dimensional initial velocity model is established, thereby improving the accuracy of the three-dimensional full waveform inversion, and then improving the prediction accuracy. BRIEF DESCRIPTION OF THE DRAWINGS

[0028] Figure 1 A flowchart of the method is shown in FIG.

[0029] Figure 2 This is the principle diagram based on the common reflection point migration of the ellipsoid;

[0030] Figure 3 This is a diagram of the dual-source symmetrical joint observation system. DETAILED DESCRIPTION

[0031] The specific implementation modes of the present invention are described below so that those skilled in the art can understand the present invention. However, it should be clear that the present invention is not limited to the scope of the specific implementation modes. For those of ordinary skill in the art, as long as various changes are within the spirit and scope of the present invention as defined and determined by the attached claims, these changes are obvious, and all inventions and creations utilizing the concept of the present invention are protected.

[0032] like Figure 1 As shown, this embodiment provides a three-dimensional full waveform inversion method for tunnel earthquakes. This embodiment studies a dual-source symmetrical joint observation system, a dual-source seismic data joint processing, and an initial velocity modeling method of common reflection point offset based on ellipsoidal expansion to improve the accuracy of three-dimensional full waveform inversion of tunnel earthquakes and realize a low-cost, safe, efficient and high-precision three-dimensional advance prediction method for tunnel earthquakes.

[0033] The common reflection point migration method based on ellipsoid expansion can be used to realize the migration imaging of a certain point in space. It is believed that the optimal amplitude value of a certain point in space can be obtained when the correct velocity is given. Conversely, the velocity corresponding to the strongest amplitude value of a certain point in space is the optimal velocity of the space point. Based on this basic theory, a certain velocity interval is set for all points in space, and three-dimensional migration imaging is performed with a fixed velocity step to obtain the velocity spectrum. The velocity corresponding to the maximum amplitude value at the space point is selected as the optimal velocity, and the three-dimensional initial velocity model is obtained. On this basis, three-dimensional full waveform inversion imaging can be implemented.

[0034] like Figure 2 As shown, P0 to P0 are the common reflection points of curves C1 and C2 on the ellipsoid. G1 and G2 are the coordinates of the detection points, and S is the coordinates of the shot point. a and b are the major and minor axes of the ellipsoid curve C2. Given the speed v, the amplitude value of the point P0 (x, y, z) in the tunnel advance direction is obtained. Conversely, the speed corresponding to the maximum amplitude value at point P0 is the optimal speed.

[0035] On the basis of the above, a three-dimensional elastic wave numerical simulation is carried out to establish a three-dimensional velocity model of poor geological bodies such as caves and weak interlayers. The numerical simulation obtains three-dimensional seismic records. The model data is used to implement ellipsoidal common reflection point migration imaging to obtain three-dimensional migration imaging, which is compared with the three-dimensional velocity model to verify the accuracy of the migration method. Given the velocity range and velocity step, a constant velocity three-dimensional scan is performed to obtain the three-dimensional velocity initial model, which is compared with the correct three-dimensional velocity model, and the empirical correction coefficient is statistically calculated.

[0036] On the basis of the above, multiple tunnel experiments were carried out, and the measured seismic data of the tunnel were processed by using methods such as automatic gain, bandpass filtering, longitudinal and transverse wave separation, FK wave field separation, and velocity spectrum scanning. According to the simulation data analysis process, the empirical correction coefficients of the measured data were obtained, and finally a reliable three-dimensional initial velocity model was established.

[0037] Forward simulation is the basis and premise of full waveform inversion. It derives the three-dimensional multi-wave and multi-component forward simulation equations of linear waves, the boundary conditions of the perfectly matched layer, and the stress equations of source loading. The partial derivative equations are transformed into second-order time and twelfth-order space differential equations, and the forward simulation is implemented by programming. The adjoint equations, target functionals, stability factors, etc. of full waveform inversion are derived and implemented by Fortran programming. The calculation formula for model update during inversion calculation is: Where x, y, z represent the coordinates of underground space, m is the longitudinal wave velocity of the formation, It represents the gradient of the objective function with respect to the parameter m during the kth iteration. kis the update step size, which is usually a positive constant. Full waveform inversion uses the adjoint state method to calculate the gradient, avoiding the huge amount of calculation required to calculate the Frechet derivative. The model gradient can be obtained by cross-correlating the forward extended wave field from the shot point and the reverse time-delayed wave field from the detector point. The calculation formula is: S(x,y,z,t) and R(x,y,z,t) are the forward extension wave field and the reverse extension wave field respectively. In the forward extension, the excitation function is the wavelet estimated by the inverse bispectrum method. In the reverse extension, the excitation function becomes the residual of the simulated data and the observed data.

[0038] During implementation, a symmetrical joint observation system is deployed in the actual tunnel and dual-source data is processed jointly. Figure 3 As shown in the figure, the arrangement of observation system and data acquisition and processing are as follows: 1 detector and 24 source points are installed on each side of the tunnel wall, 2 detectors and 48 source points are used in common, and a total of 2 detector holes and 2 blast holes need to be drilled. The depth of the 2 blast holes is 1.5m, and the depth of the 2 detector holes is 2m. During data acquisition, 2 shots are fired in sequence at the blast hole positions using explosive sources, and 48 shots are fired at other source points using mechanical impact sources. The elastic wave data at different positions of the tunnel are obtained through the detectors, and the longitudinal and transverse wave signals in the elastic wave data are extracted. The surrounding rock velocity is calculated using the direct wave velocity fitting to obtain seismic data with high signal-to-noise ratio on both sides of the tunnel wall; the joint processing of the explosive source and the mechanical impact source adopts the empirical Green's function calculation method, and its expression is: Among them, x belongs to the variable of spatial range A, x∈A, r0 is the distance in space A, u(x,t) is the observation record at point x, and v(x+r0,t) is the observation record at point x+r0. Combined with the symmetrical observation system, it can be considered that u(x,t) is the mechanical-driven impact source record of the 1st to 24th shots, and v(x+r0,t) is the explosive source record of the 24th shot. If the explosive source record is missing, v(x+r0,t) is the mechanical-driven impact source record of the 24th shot. On this basis, the common reflection point migration method based on ellipsoid expansion is used to establish a reliable three-dimensional initial velocity model; finally, the tunnel seismic true three-dimensional full waveform inversion imaging is implemented to accurately calculate the longitudinal and transverse wave velocities and Poisson's ratio parameters of the surrounding rock. The later excavation results are used to reversely verify the accuracy of the calculation results, obtain the actual correction coefficient, and finally determine a reliable three-dimensional initial velocity model.

Claims

1. A tunnel seismic three-dimensional full waveform inversion method, characterized in that: The following steps are involved: S1, establishing a three-dimensional initial velocity model: based on the common reflection point migration method of ellipsoid expansion, three-dimensional migration imaging is performed to obtain a velocity spectrum, the velocity corresponding to the maximum amplitude is selected as the optimal velocity, and a three-dimensional initial velocity model is established; S2, forward simulation wave field characteristics analysis: establish a typical bad geological body model, the typical bad geological body includes three-dimensional karst caves and fault fracture zones; conduct forward simulation in the established typical bad geological body model, set source points on the tunnel wall or the tunnel face, set detection points on the tunnel wall, and analyze the characteristics of seismic wave signals received at the tunnel wall and the tunnel face; S3, full waveform inversion numerical simulation: using the seismic wave signal data received in S2 and based on the three-dimensional initial velocity model in S1, a full waveform inversion numerical model experiment is performed to obtain simulation data; S4, deployment of observation system and data acquisition and processing: a symmetrical joint observation system is deployed in the actual tunnel, and the data of the explosive source and the mechanical impact source are jointly processed to obtain the actual data; S5, 3D full waveform inversion imaging: by comparing the simulation data obtained in S3 with the actual data obtained in S4, the residual is calculated, the gradient is calculated using the residual, and the gradient is used to update the parameters of the 3D initial velocity model to obtain a reliable 3D velocity model, realize the 3D full waveform inversion of the tunnel, and accurately calculate the physical quantities of the surrounding rock, including the longitudinal wave velocity, the shear wave velocity, the Poisson's ratio, the Young's modulus and the volume modulus; S6, result verification: reversely verify the accuracy of the calculation results, obtain the actual correction coefficients, and finally determine a reliable three-dimensional initial velocity model.

2. The tunnel seismic 3D full waveform inversion method according to claim 1 is characterized in that: The forward modeling in S2 is the process of solving the wave equation. For three-dimensional multi-component, the wave equation is expressed as formula (1): Among them, x, y, z represent the coordinates of underground space, m is the longitudinal wave velocity of the formation, It represents the gradient of the objective function with respect to the parameter m during the kth iteration, α k is the update step size.

3. The tunnel seismic 3D full waveform inversion method according to claim 2 is characterized in that: The cross-correlation calculation formula of the forward extended wave field from the shot point and the reverse time-delayed wave field from the detection point is as follows: Among them, S(x,y,z,t) and R(x,y,z,t) are the forward extended wave field and the reverse time-delayed wave field respectively.

4. The tunnel seismic 3D full waveform inversion method according to claim 1, characterized in that: The deployment of the observation system and data acquisition processing in S4 is specifically as follows: detectors and source points are deployed on both sides of the tunnel wall, and the source points include explosive sources and mechanical impact sources; the two sources are stimulated, and elastic wave data at different positions of the tunnel are obtained through the detectors; the longitudinal and transverse wave signals in the elastic wave data are extracted, and the surrounding rock velocity is calculated using the direct wave velocity fitting.

5. The tunnel seismic 3D full waveform inversion method according to claim 4 is characterized in that: The combined processing of the explosive source and the mechanical impact source adopts the empirical Green's function calculation method, which is expressed as formula (3): Among them, x belongs to the variable of spatial range A, x∈A, r0 is the distance in space A, u(x,t) is the observation record at point x, and v(x+r0,t) is the observation record at point x+r0.

6. The tunnel seismic 3D full waveform inversion method according to claim 1, characterized in that: The simulation data in S3 include rock physical quantities, including longitudinal wave velocity, transverse wave velocity, Poisson's ratio, Young's modulus and volume modulus.

Images