CN105022031B - A kind of layered velocity localization method of region rock mass microseism focus - Google Patents

Abstract

The invention belongs to the field of geotechnical engineering, and provides a layered velocity positioning method for regional rock mass microseismic sources. Obtain the analytical equation of the interface of each wave velocity layer; ② Install sensors in the rock mass area and measure the spatial coordinates of each sensor; ③ Carry out blasting experiments in the rock mass area, and calculate the rock mass wave velocity of each wave velocity layer; ④ Calculate the microseismic source The space coordinates (x ₀ ,y ₀ ,z ₀ ). The method of the invention has high positioning accuracy for the microseismic source, and can promote the microseismic monitoring technology to better play the role of prediction and early warning in engineering practice.

Description

A Layered Velocity Location Method for Microseismic Sources in Regional Rock Mass

technical field

The invention belongs to the field of geotechnical engineering, in particular to a layered velocity positioning method of regional rock microseismic sources.

Background technique

Under the action of external stress, local elastic-plastic energy concentration will occur inside the rock mass. When the energy accumulates to a certain critical value, it will cause the generation and expansion of micro-cracks in the rock mass. The generation and expansion of micro-cracks are accompanied by elastic wave Or the release of stress waves and rapid propagation in the surrounding rock mass, this elastic wave is called microseismic in geology. Microseismic monitoring technology is to collect microseismic signals by arranging sensors in the rock mass of the monitoring area, and analyze the collected signals to locate the microseismic source and delineate the damaged area of the rock mass, providing an effective basis for the stability evaluation of the regional rock mass. Accurate positioning of the microseismic source is the basis for explaining the focal mechanism, calculating the source rupture scale, and delineating the rock mass damage area. The low accuracy of the source positioning will lead to inaccurate calculation of the source rupture scale and incorrect delineation of the rock mass damage area. This leads to misjudgment of the instability of the rock mass in the construction site area, which cannot provide effective guidance for the construction organization in engineering practice. Therefore, in order to make microseismic monitoring technology play a better role in prediction and early warning in engineering practice, it is very important to improve the accuracy of microseismic source location.

At present, the microseismic source location methods include geometric mapping method, relative location method, spatial domain location method, linear location method and nonlinear location method, among which the methods that receive more attention are linear location method and nonlinear location method. During the positioning process of the microseismic source, the waveform signal received by the sensor contains background noise such as current interference and construction machinery impact, so it is not easy to accurately pick up the microseismic signal waveform at the take-off time, and the pick-up accuracy at the take-off time is not high. , the positioning effect of the linear positioning method is poor. Because the nonlinear positioning method can obtain the best approximation solution under the condition of low picking accuracy at the take-off time, its stability and reliability are higher, and it has become the focus of current researchers.

The existing nonlinear localization method is to locate the microseismic source on the basis of simplifying the rock mass in the microseismic monitoring area to a single and uniform rock mass wave velocity. For this simplified single velocity model, the i-th sensor receives the P wave The relationship between the first arrival time t _i and source parameters (x ₀ , y ₀ , z ₀ , t ₀ ) is shown in formula (11):

In formula (11), l _i is the spatial distance between the microseismic source and the i-th sensor, ( _xi , y _i , z _i ) is the spatial coordinates of the i-th sensor, and m is the sensor that receives the P wave in the monitoring area Quantity, V is the P-wave velocity in the monitoring area, i=1,2,...,m.

Theoretically, the sum of the microseismic event occurrence time and travel time t ₀ +Δt _i should be equal to t _i , as shown in formula (12):

ξ _i =t _i -(Δt _i +t ₀ )=0 (12)

However, in actual engineering, due to the influence of microseismic monitoring instruments and human factors, there is an error in picking up the microseismic signal when it jumps, that is, the value of _ξi is not zero. In this case, in order to realize microseismic source location, it is usually necessary to construct an objective function whose dependent variable is the arrival time or the dependent variable is the arrival time difference or the dependent variable is the arrival time difference quotient, and use nonlinear optimization tools to solve the above objective function, with the The source parameter approximation solution whose objective function is as small as possible is taken as the final calculation result of the microseismic source parameters.

In fact, many rock masses in nature are formed by geological structures in different geological periods. 75% of the earth's land area is composed of sedimentary rocks, which are characterized by layered distribution of rock mass media, such as layered coal seams, layered Slopes, etc., the rock mass wave velocity of each layer of rock mass medium is not the same. In addition, due to the influence of on-site construction, cracks or fissures are produced in the rock mass, which causes the wave velocity of the rock mass to change with the advancement of on-site construction, that is, the wave velocity of the rock mass in the microseismic monitoring area is not always a constant value. Therefore, the existing non-linear positioning method locates the microseismic source on the basis of uniform and single rock mass wave velocity, which will inevitably be affected by the layering of the rock mass medium and on-site construction, and the accuracy of its source positioning will inevitably be limited.

Contents of the invention

The purpose of the present invention is to overcome the deficiencies of the prior art and provide a layered velocity positioning method for regional rock mass microseismic sources to improve the positioning accuracy of microseismic sources and promote microseismic monitoring technology to better play the role of prediction and early warning in engineering practice .

The layered velocity positioning method of the regional rock mass microseismic source of the present invention, the steps are as follows:

① According to the rock mass wave velocity of the rock mass area of the seismic source to be measured, the rock mass area is divided into different wave velocity layers. The wave velocity layer is recorded as the first wave velocity layer, the second wave velocity layer, ..., the kth wave velocity layer from one side to the other in the direction perpendicular to the wave velocity layer, and the rock mass wave velocity of each wave velocity layer is recorded as V ₁ , V ₂ ,...,V _k , establishing a three-dimensional Cartesian coordinate system with any point in the rock mass region as the coordinate origin, and calculating the analytical equation of the interface of each wave velocity layer according to the thickness and occurrence of each wave velocity layer;

② Install m sensors in the rock mass area, m≥5, each sensor forms a spatial network structure distribution, and sensors must be installed in the first wave velocity layer and the kth wave velocity layer, measure the spatial coordinates of each sensor, and place the first The spatial coordinates of the _i sensors are denoted as ( _xi , y, _zi );

③ Carry out blasting experiments in the rock mass area, so that microseisms with known source locations occur in the rock mass area. The spatial coordinates of the source are (x _b , y _b , z _b ), and the i-th sensor collects The take-off time of the P wave generated by the microseism is recorded as t _i , and the wave velocity layers that the P wave generated by the microseism propagates from the source to the i-th sensor pass through are numbered 1, 2,..., n. The wave velocity of the rock mass is recorded as V _i1 , V _i2 ,...,V _in , and the constant items of the analytical equations of the wave velocity layer interface passing through in turn are recorded as D _i(1,2) , D _i(2,3) ,..., D _i(n-1,n) , D _i(1,2) represents the constant term of the analytical equation of the P-wave generated by the microseism that propagates from the source to the interface between the No. 1 wave velocity layer and the No. 2 wave velocity layer that the i-th sensor passes through , D _i(2,3) ,..., the meanings represented by D _i(n-1,n) and so on;

If the above-mentioned seismic source and the i-th sensor are located in different wave velocity layers, the theoretical value Δt _i of the elapsed time for the P-wave generated by the microseismic source at the known source position to propagate from the seismic source to the i-th sensor is shown in formula (1) :

In formula (1), α _i and β _i are spatial distance parameters, α _i = |Axi +By _i +Cz _{i -Ax b} -By _b -Cz _b |, A, B, and C are the coefficients of x, y, and z in the analytical equation of the interface of each wave velocity layer;

If the above-mentioned seismic source and the i-th sensor are located in the same wave velocity layer, the theoretical value Δt _i of the elapsed time for the P-wave generated by the microseismic source at the above-mentioned known source position to propagate from the seismic source to the i-th sensor is shown in formula (2) :

Establish the objective function f for solving V ₁ , V ₂ ,..., V _k described in formula (3),

In formula (3), t _i and t _j are the take-off moments when the i-th and j-th sensors receive the P waves generated by microseisms respectively; in the process of solving V ₁ , V ₂ ,...,V _k , it is necessary The corresponding relationship between V _i1 , V _i2 ,…,V _in and the wave velocity of rock mass in each wave velocity layer, replace V _i1 ,V _i2 ,…,V _{in in} formula (1) and formula (2) with V ₁ , V ₂ ,...,V _k , use nonlinear optimization tools to solve the minimum value of the objective function f described in formula (3), and use the rock mass wave velocities V ₁ , V ₂ ,...,V _k corresponding to the minimum value of f as step ④ to solve the unknown seismic source wave velocity of rock mass;

④ When a microseism with an unknown source location occurs in the rock mass area, assume that the spatial coordinates of the source are (x ₀ , y ₀ , z ₀ ),

If the seismic source and the i-th sensor are located in different wave velocity layers, the theoretical value Δt _i of the elapsed time for the P-wave generated by the microseismic source at the unknown source position to propagate from the seismic source to the i-th sensor is shown in formula (4):

In formula (4), α _i and β _i are spatial distance parameters, α _i =|Axi +By _i +Cz _{i -Ax 0} -By ₀ -Cz ₀ |, A, B, and C are the coefficients of x, y, and z in the analytical equation of the interface of each wave velocity layer;

If the above-mentioned source and the i-th sensor are located in the same wave velocity layer, the theoretical value Δt _i of the elapsed time for the P-wave generated by the microseismic source at the above-mentioned unknown source position to propagate from the source to the i-th sensor is shown in formula (5):

Establish the sub-velocity layer described in formula (6) to solve the objective function f _η of the unknown source space coordinates (x ₀ , y ₀ , z ₀ ),

In formula (6), t _i , t _j are the take-off moments when the i-th and j-th sensors receive the P waves generated by microseisms, η=1,2,...,k, and k is the The total number of wave velocity layers. _In the process of solving the unknown source space coordinates (x ₀ , y ₀ , _{z 0} ), it is necessary to _formulate (4) and V _i1 , V _i2 ,...,V _{in in} formula (5) are replaced by V ₁ , V ₂ ,...,V _k , and the objective function f _η described in formula (6) is solved by nonlinear optimization tools The minimum value of η takes 1, 2, ..., k, and the minimum value of f _η obtained from the solution is respectively recorded as f _1min , f _2min , ..., f _kmin ; then solve f _1min , f _2min , …, the minimum value f _θ in f _kmin ,

f _θ ＝min{f _1min ,f _2min ,...,f _kmin } (7)

The (x ₀ , y ₀ , z ₀ ) corresponding to f _θ is the spatial coordinate of the unknown source;

When using nonlinear optimization tools to solve the minimum value of f _η , the search for the feasible region Ω _η should satisfy the linear programming inequality described in formula (8):

R _η ·[x ₀ ,y ₀ ,z ₀ ] ^T ≤b _η (8)

In formula (8),

A _η1 ,…,A _ηy ,B _η1 ,…,B _ηy ,C _η1 ,…,C _ηy ,D _η1 ,…,D _ηy are the parameters of the analytical equations of the boundary surfaces of the nth wave velocity layer respectively, and y is the The number of each boundary surface of the η wave velocity layer.

Because in engineering practice, the interface of each wave velocity layer of the regional rock mass is not absolutely parallel, therefore, the wave velocity of the rock mass in the same wave velocity layer described in step ① of the above method is equal and the interface of each wave velocity layer is regarded as parallel to each other It is an approximate treatment of each wave velocity layer of the regional rock mass in engineering practice, in order to simplify the operation of solving the unknown seismic source space coordinates (x ₀ , y ₀ , z ₀ ) and reduce the calculation amount when solving the objective function.

In step ② of the above method, in order to achieve a good microseismic monitoring effect, it is best to arrange the sensors in different wave velocity layers.

In the above method, the nonlinear optimization tool is simplex method, particle swarm optimization algorithm or genetic algorithm.

In the above method, one burst test is carried out in step ③, or at least two burst tests are carried out intermittently. Considering that the on-site construction will cause cracks or cracks in the rock mass, which will cause the wave velocity of the rock mass to change with the advancement of the on-site construction, therefore, as the on-site construction advances, it is preferable to conduct at least two blasting tests intermittently in step ③ to ensure Updating the wave velocity values of rock mass in each wave velocity layer in real time, and using it as a basis to locate the microseismic source will help improve the positioning accuracy of the source.

Taking the source and sensor locations shown in Figure 2 as an example, the derivation process of formula (4) is explained below:

The rock mass area in Figure 2 is divided into four wave velocity layers, which are recorded as the first, second, third, and fourth wave velocity layers from bottom to top, and the rock mass wave velocities of each wave velocity layer are recorded as V ₁ , V ₂ , V ₃ , V ₄ , the top of the fourth wave velocity layer is the ground line, the microseismic source is located in the first wave velocity layer, the position of the microseismic source is recorded as point P, and its spatial coordinates are (x ₀ , y ₀ , z ₀ ), The sensor numbered 2-1 is located in the fourth wave velocity layer, and the position of the sensor is recorded as point S, and its space coordinates are (x ₁ , y ₁ , z ₁ ).

When the actual seismic wave propagates between different media, it follows Feynman's principle and propagates along a broken line (the path shown by the solid arrow in Figure 2), that is, it propagates along the path of minimum travel time. But because the propagation speed of P wave is fast, and the microseismic monitoring area size is only a few hundred meters in a circle usually, the thickness of each wave velocity layer is not big, so the present invention approximates the propagation path of the P wave that microseismic produces as straight line (in Fig. 2 The path shown by the dotted arrow), based on this approximate assumption, the P wave generated by the microseismic will pass through different wave velocity layers during the process of propagating from the microseismic source to the sensor, and the corresponding travel time in each wave velocity layer can be determined by the corresponding spatial distance and wave velocity value show.

Taking the P wave generated by the microseismic source received by the sensor in Figure 2 as an example, the derivation (4) is as follows:

Since the interfaces of the wave velocity layers are parallel to each other, the interface between the first wave velocity layer and the second wave velocity layer, the interface between the second wave velocity layer and the third wave velocity layer, and the interface between the third wave velocity layer and the fourth wave velocity layer The analytical equation is shown in formula (13):

The intersection point of the straight line PS with the interface between the first wave velocity layer and the second wave velocity layer is L, the intersection point with the interface between the second wave velocity layer and the third wave velocity layer is M, and the intersection point with the interface between the third wave velocity layer and the fourth wave velocity layer is N, the lengths of the line segments PL, LM, MN and NS are shown in formula (14):

In formula (14), α ₁ =|Ax ₁ +By ₁ +Cz ₁ -Ax ₀ -By ₀ -Cz ₀ |,

According to the length of the line segments PL, LM, MN and NS, combined with the rock mass wave velocity of each wave velocity layer, the time Δt ₁ for the P wave generated by the microseismic source to propagate from the source to the sensor can be obtained,

According to the analogy of the derivation process of Δt ₁ , it can be known that when microseisms occur in the rock mass area, if the seismic source and the i-th sensor are located in different wave velocity layers, if the spatial coordinates of the i-th sensor are denoted as (x _i , y _i , z _i ), then the theoretical value Δt _i of the elapsed time for the P wave generated by the microseism to propagate from the source to the i-th sensor can be expressed as:

In formula (4), α _i and β _i are spatial distance parameters, α _i =|Axi +By _i +Cz _{i -Ax 0} -By ₀ -Cz ₀ |,

The derivation process of formula (1) is similar to formula (4).

Compared with the prior art, the present invention has the following beneficial effects:

1. The present invention provides a new method for regional rock mass microseismic source location. Because the method divides the rock mass region into different wave velocity layers according to the actual geological structure of the microseismic source rock mass area to be measured, the non-linear The optimization tool divides the wave velocity layer to solve the spatial coordinates of the seismic source, and then takes the spatial coordinates corresponding to the minimum value of the objective function f _η as the final spatial coordinates of the seismic source. Compared with the existing nonlinear positioning method based on a single velocity model to determine the source position, this The method described in the invention is more in line with the actual geological structure of the regional rock mass, and is especially suitable for the location of the microseismic source of the regional rock mass in layered distribution. Experiments show that the method of the present invention has effectively improved the accuracy of source location (see implementation example and comparative example), which is conducive to microseismic monitoring technology to play a better role in prediction and early warning in engineering practice.

2. Since the method of the present invention has fully considered the impact of on-site construction on the rock mass wave velocity, the blasting test of the known microseismic source position is combined with the use of nonlinear optimization tools to invert the rock mass wave velocity of each wave velocity layer. The wave velocity of the rock mass in the wave velocity layer can more truly and accurately reflect the real-time rock mass wave velocity, which is conducive to the further improvement of the accuracy of the microseismic source.

Description of drawings

Fig. 1 is a schematic diagram of the wave velocity layer of the regional rock mass, the installation position of the sensor and the position of the preset blast hole in the embodiment;

Fig. 2 is a schematic diagram of deriving Δt _i . In Fig. 2, the solid line arrow indicates the minimum travel time path of the P wave generated by the microseism, and the dotted line arrow represents the shortest distance path of the P wave generated by the microseism;

In the figure, 1—rock mass, 1-1—first wave velocity layer, 1-2—second wave velocity layer, 1-3—third wave velocity layer, 1-4—fourth wave velocity layer, 2-1, 2- 2. 2-3, 2-4, 2-5, 2-6, 2-7, 2-8 are sensor numbers, 3-1—the first blast hole, 3-2—the second blast hole, 4—microseismic source.

detailed description

The layered velocity positioning method of the regional rock microseismic source of the present invention will be further described below in conjunction with the accompanying drawings. It must be pointed out that the following examples are only used to further illustrate the present invention, and should not be interpreted as limiting the protection scope of the present invention. Those skilled in the art can make some non-essential improvements and adjustments to the present invention according to the above-mentioned content of the invention. Implementation still belongs to the scope of protection of the invention.

Example

In this embodiment, the size of the rock mass area of the microseismic source to be tested is 200m × 200m × 200m, and the positioning steps of the microseismic source are as follows:

①According to the rock mass wave velocity in the rock mass area of the microseismic source to be tested recorded in the geological exploration data, the above rock mass area is divided into three wave velocity layers, the rock mass wave velocity in the same wave velocity layer is equal and the interface of each wave velocity layer In order to be parallel to each other, each wave velocity layer is recorded as the first wave velocity layer 1-1, the second wave velocity layer 1-2, and the third wave velocity layer 1-3 from the upper side to the lower side in the direction perpendicular to the wave velocity layer. The wave velocity of the rock mass in the layer is recorded as V ₁ , V ₂ , V ₃ in turn; the three-dimensional Cartesian coordinate system is established with the O point in the rock mass area as the coordinate origin, and the X, Y, and Z directions in the three-dimensional Cartesian coordinate system correspond to Geology E (east), N (north) and U (elevation), O point is 60m away from the intersection of the interface between the second wave velocity layer and the third wave velocity layer and the Z axis, and O point is 60m away from the first wave velocity layer and the third wave velocity layer. The distance between the interface of the two wave velocity layers and the intersection point of the Z axis is 150m. The occurrences of the above three wave velocity layers are all N30°E, SE∠15°. According to the thickness and occurrence of each wave velocity layer, the analytical equation Ax+By+Cz+D=0 of the interface of each wave velocity layer is calculated , respectively as follows:

The analytical equation of the interface between the first wave velocity layer and the second wave velocity layer is: 2241x-1294y+9659z-1448850=0

The analytical equation of the interface between the second wave velocity layer and the third wave velocity layer is: 2241x-1294y+9659z-579540=0

It can be seen from the analytical equations of the interfaces of the above-mentioned wave velocity layers that A=2241, B=1294, and C=9659.

②As shown in Figure 1, install eight acceleration sensors in the above rock mass area, and number each sensor as 2-1, 2-2, 2-3, 2-4, 2-5, 2-6, 2- 7 and 2-8, these eight sensors form a spatial network structure distribution, measure the spatial coordinates of each sensor, and denote the spatial coordinates of the sensor numbered 2-i as ( _xi , y _i , z _i ), each The spatial coordinates of the sensors and the wave velocity layer where each sensor is located are shown in Table 1.

Table 1 The spatial coordinates and wave velocity layers of each sensor

Drill two blast holes in the above-mentioned rock mass area, respectively denoted as the first blast hole 3-1 and the second blast hole 3-2, respectively install 200g of emulsion explosive at the bottom of each blast hole, and connect the detonating wire And the high-voltage electrostatic detonator, the orifice of each blast hole is sealed with loose soil particles on site to reduce the energy loss during blasting. The first blast hole 3-1 is located in the second wave velocity layer 1-2, the coordinates of the center of the hole bottom are (94,98,133), and the coordinates of the center of the second blast hole 3-2 are (61,73,101) .

In the subsequent steps, a blast test is carried out in the first blast hole to cause microseisms with a known source location in the rock mass region; a blast test is carried out in the second blast hole for simulating unknown Microseisms at the source location.

③Detonate the emulsion explosive in the first blast hole 2-1 to carry out the blasting experiment to cause a microseismic earthquake with a known source location in the rock mass area, and the spatial coordinates (x _b , y _b , z _b ) of the source are the first The coordinates at the center of the hole bottom of blast hole 3-1, namely (x _b , y _b , z _b )=(94,98,133), start from the detonation moment of the emulsion explosive in the first blast hole (that is, the moment when microseisms occur) For timing, the take-off time of the P wave generated by the microseism collected by the sensor numbered 2-i is recorded as t _i , as shown in Table 2.

Table 3 The take-off time t _i of the P wave received by each sensor

Sensor number 2-1 2-2 2-3 2-4 2-5 2-6 2-7 2-8 25.5 27.2 27.4 25.8 16.6 18.9 19.9 17.2

The P waves generated by the microseisms propagate from the seismic source to the i-th sensor and pass through the wave velocity layers in turn. They are numbered as 1, 2,..., n wave velocity layers, and the rock mass wave velocities of the wave velocity layers passing through in turn are recorded as V _i1 , V _i2 ,... ,V _in , the constant items of the analytical equations of the wave velocity layer interface passing through in turn are recorded as D _i(1,2) , D _i(2,3) ,…, D _i(n-1,n) , D _{i( 1,2)} represents the constant term of the analytical equation of the P-wave generated by the microseism that propagates from the seismic source to the interface between the No. 1 wave velocity layer and the No. 2 wave velocity layer that the i-th sensor passes through, D _i(2,3) ,...,D _{i (n-1,n)} represents the meaning and so on;

If the above-mentioned seismic source and the sensor numbered 2-i are located in different wave velocity layers, the theoretical value Δt _i of the elapsed time for the P wave generated by the microseismic source at the known source position to propagate from the source to the sensor numbered 2-i is as follows Formula (1) shows:

In formula (1), α _i and β _i are spatial distance parameters, α _i = |Axi +By _i +Cz _{i -Ax b} -By _b -Cz _b |, A, B, and C are the coefficients of x, y, and z in the analytical equation of the interface of each wave velocity layer;

If the above-mentioned seismic source and the sensor numbered 2-i are located in the same wave velocity layer, then the theoretical value Δt _i of the elapsed time for the P wave generated by the microseismic source at the known source position to propagate from the source to the sensor numbered 2-i is as follows Formula (2) shows:

Specifically, Δt ₁ to Δt ₈ are as follows:

Taking the sensor numbered 2-1 as an example, explain the corresponding relationship between V ₁₁ , V ₁₂ and V ₁ , V ₂ , V ₃ in Δt ₁ and the method for determining the value of D _1(1,2) . Since the sensor numbered 2-1 is located in the third wave velocity layer 1-3, and the first blasthole 3-1 is located in the second wave velocity layer 1-2, the P The propagation of the wave into the seismic source to the sensor numbered 2-1 needs to pass through the second wave velocity layer, the interface between the second wave velocity layer and the third wave velocity layer, and the third wave velocity layer, that is, V ₁₁ =V ₂ , V ₁₂ =V ₃ , The value of D _1(1,2) is a constant term of the analytical equation of the interface between the second wave velocity layer and the third wave velocity layer, D _1(1,2) =-579540. Taking the sensor numbered 2-2 as an example, explain the corresponding relationship between V ₂₁ in Δt ₂ and V ₁ , V ₂ , V ₃ , since the sensor numbered 2-2 and the first blast hole 3-1 are located at the first In the second wave velocity layer 1-2, the P wave generated by the microseism caused by the blasting test of the first blast hole only needs to pass through the second wave velocity layer, ie, V ₂₁ =V ₂ , to propagate to the sensor numbered 2-2.

Refer to the corresponding relationship between V _i1 , V _i2 ,...,V _in and V ₁ , V ₂ , V ₃ in Δt ₁ and Δt ₂ and D _i(1,2) , D _i(2,3) ,...,D The method of determining the value of _i(n-1,n) is to obtain the corresponding relationship between V _i1 , V _i2 ,...,V _in and V ₁ , V ₂ , V ₃ in the above-mentioned Δt ₁ to Δt ₈ and D _{i( 1,2)} , D _i(2,3) ,..., the values of D _i(n-1,n) are shown in Table 4 and Table 5 respectively.

Table 4 Correspondence between V _i1 , V _i2 ,...,V _in and V ₁ , V ₂ , V ₃ in Δt ₁ ～Δt ₈

Table 5 Δt ₁ ~ Δt ₈ , D _i(1,2) , D _i(2,3) , ..., D _i(n-1,n) values

corresponding value -579540 -579540 -1448850 -1448850 -1448850

Set up the described solution V ₁ of formula (3), V ₂ , the objective function f of V ₃ ,

In formula (3), t _i and t _j are the take-off moments when the sensors numbered 2-i and 2-j respectively receive the P waves generated by microseisms; in the process of solving V ₁ , V ₂ , V ₃ , replace V _i1 , V _i2 ,...,V _{in in} formula (1) and formula (2) with V ₁ , V ₂ , V ₃ correspondingly according to Table 4, and assign V in formula (1) according to Table 5 D _i(1,2) , D _i(2,3) ,..., the corresponding value of D _i(n-1,n) , will (x _b ,y _b ,z _b ) and ( _xi ,y _i , z _i ), and the values of A, B, and C are brought into formula (1), and then this embodiment adopts the genetic algorithm to solve the minimum value of the objective function f described in formula (3), and the rock mass corresponding to the minimum value of f The wave velocities V ₁ , V ₂ , and V ₃ are used as step ④ to solve the rock mass wave velocities of the unknown source, and finally get V ₁ =3891m/s, V ₂ =4084m/s, and V ₃ =4195m/s.

④ Detonate the emulsion explosive in the second blast hole 2-2 to conduct a blasting experiment to simulate a microseism in the rock mass area with an unknown source location, assuming that the spatial coordinates of the source are (x ₀ , y ₀ , z ₀ ), from The detonation moment of the emulsion explosive in the second blasting hole (that is, the moment when the microseism occurs) starts counting, and the take-off time when the sensor numbered 2-i collects the P wave generated by the microseism is recorded as t _i , as shown in Table 6 .

Table 6 The take-off time t _i of the P wave received by each sensor

Sensor number 2-1 2-2 2-3 2-4 2-5 2-6 2-7 2-8 13.7 25.9 31.1 22.5 13.3 25.7 31.2 22.3

The P wave generated by the microseism is propagated from the source to the sensor numbered 2- _i . The wave velocity layers that pass through in turn are numbered 1, 2,..., n. _i2 ,…,V _in , the constant items of the analytical equations passing through the wave velocity layer interface in turn are denoted as D _i(1,2) , D _i(2,3) ,…,D _i(n-1,n) , D _i(1,2) represents the constant term of the analytical equation of the interface between the No. 1 wave velocity layer and the No. 2 wave velocity layer where the P wave generated by the microseism propagates from the source to the sensor numbered 2-i, and D _{i(2,3 )} ,..., the meaning represented by D _i(n-1,n) and so on;

If the above-mentioned seismic source and the sensor numbered 2-i are located in different wave velocity layers, then the theoretical value Δt _i of the elapsed time for the P wave generated by the microseismic source at the unknown source position to propagate from the source to the sensor numbered 2-i is as follows: (4) as shown:

In formula (4), α _i and β _i are spatial distance parameters, α _i =|Axi +By _i +Cz _{i -Ax 0} -By ₀ -Cz ₀ |, A, B, and C are the coefficients of x, y, and z in the analytical equation of the interface of each wave velocity layer;

If the above-mentioned seismic source and the sensor numbered 2-i are located in the same wave velocity layer, then the theoretical value Δt _i of the elapsed time for the P wave generated by the microseismic source of the above-mentioned unknown source to propagate from the source to the sensor numbered 2-i is as follows: 5) As shown:

Establish the sub-velocity layer described in formula (6) to solve the objective function f _η of the unknown source space coordinates (x ₀ , y ₀ , z ₀ ),

In formula (6), t _i and t _j are respectively the take-off time when the sensors numbered 2-i and 2-j receive the P wave generated by the microseism, η=1, 2, 3, and the genetic algorithm is used to solve the formula The minimum value of the objective function f _η shown in (6), η takes times 1, 2, 3, and the minimum value of f _η obtained from the solution is recorded as f _1min , f _2min , f _3min respectively; In the process of (x ₀ , y ₀ , z ₀ ), replace V _i1 , V _i2 ,…,V _{in in} formula (4) and formula (5) with V ₁ , V ₂ , V ₃ , respectively Assign corresponding values to D _i(1,2) , D _i(2,3) ,..., D _i(n-1,n) in formula (4) (V _i1 , V _i2 ,...,V _in and V ₁ , V ₂ , V ₃ and the method for determining the value of D _i(1,2) , D _i(2,3) ,..., D _i(n-1,n) refer to step ③), and the ( The values of x _i , y _i , z _i ) and the values of A, B, and C are brought into formula (4), and then the minimum value of the objective function f _η described in formula (6) is solved by genetic algorithm.

When the present embodiment utilizes the genetic algorithm to _solve the minimum value of f, the search feasible region _Ω should satisfy the linear programming inequality described in formula (8):

R _η ·[x ₀ ,y ₀ ,z ₀ ] ^T ≤b _η (8)

In formula (8),

A _η1 ,…,A _ηy ,B _η1 ,…,B _ηy ,C _η1 ,…,C _ηy ,D _η1 ,…,D _ηy are the parameters of the analytical equations of the boundary surfaces of the nth wave velocity layer respectively, and y is the The number of each boundary surface of the η wave velocity layer.

Specifically, since the first wave velocity layer, the second wave velocity layer, and the third wave velocity layer are all hexahedral, that is, each wave velocity layer has six boundary surfaces, therefore, in formula (8), y=6. In the process of solving f _1min , f _2min , and f _3min , searching for the feasible regions Ω ₁ , Ω ₂ , and Ω ₃ should satisfy the following linear programming inequality:

η=1, which means that the hypocenter space coordinates (x ₀ , y ₀ , z ₀ ) are assumed to be located in the first wave velocity layer,

The search feasible region Ω ₁ is:

η=2, which means that it is assumed that the source space coordinates (x ₀ , y ₀ , z ₀ ) are located in the second wave velocity layer,

The search feasible region Ω ₂ is:

η=3, which means that the hypothetical source space coordinates (x ₀ , y ₀ , z ₀ ) are located in the third wave velocity layer,

The search feasible region Ω ₃ is:

According to the above-mentioned searching feasible domain Ω ₁ , Ω ₂ , Ω ₃ , the minimum values of f _η obtained by using the genetic algorithm are respectively:

f _1min = 0.016s, the corresponding spatial coordinates of the seismic source (x ₀ , y ₀ , z ₀ ) = (52.6, 84.1, 113.8),

f _2min = 0.004s, the corresponding spatial coordinates of the seismic source (x ₀ , y ₀ , z ₀ ) = (60.4, 76.1, 99.5),

f _3min =0.027s, corresponding to the spatial coordinates of the seismic source (x ₀ , y ₀ , z ₀ )=(64.5, 70.2, 95.7).

Solve the minimum value f _θ among f _1min , f _2min , f _3min according to formula (7),

f _θ ＝min{f _1min ,f _2min ,f _3min }＝f _2min ＝0.004s (7)

(x ₀ , y ₀ , z ₀ )=(60.4, 76.1, 99.5) corresponding to f _θ is the spatial coordinate of the unknown source.

And the real space coordinates of this seismic source are the coordinates (61,73,101) at the center of the hole bottom of the second blast hole 3-2, and the distance between this real seismic source and the seismic source located in this embodiment is:

comparative example

In this comparison, the existing nonlinear positioning method is used to locate the microseismic source, that is, the rock mass area with a size of 200m×200m×200m as shown in Figure 1 is regarded as a medium with uniform wave velocity. According to the exploration data, it is known that the rock mass area The average rock mass wave velocity V in the body area is 4000m/s.

Since this comparative example does not need to divide the rock mass region into different wave velocity layers, after detonating the emulsion explosive in the second blast hole 2-2 in step ④ of the embodiment in this comparative example, the number is 2-i The spatial coordinates (x ₀ , y ₀ , z ₀ ) of the location seismic source are determined by the take-off time t _i (see Table 6) of the P wave collected by the sensors and the average wave velocity V=4000m/s in the above rock mass region, the steps are as follows :

The theoretical value Δt _i of the elapsed time for the P wave generated by the microseismic source with unknown source position to propagate from the source to the sensor numbered 2-i is shown in formula (8):

Establish the objective function F for solving the unknown seismic source space coordinates (x ₀ , y ₀ , z ₀ ) described in formula (9),

In formula (10), t _i , t _j are respectively the take-off time when the sensors numbered 2-i and 2-j receive the P wave generated by the microseism;

Put the average rock mass wave velocity value V=4000m/s and the value of ( _xi , y _i , _zi ) in Table 6 into formula (10), and then use the nonlinear optimization tool genetic algorithm to solve the formula (10) The minimum value F _min of the above-mentioned objective function F is F _min =0.031s.

(x ₀ , y ₀ , z ₀ )=(57.3, 65.1, 112.5) corresponding to F _min is the spatial coordinate of the unknown source, and the real spatial coordinate of the source is the center of the bottom of the second blast hole 3-2 Coordinates (61,73,101), the distance between the real seismic source and the seismic source located in this comparative example is:

The distance between the microseismic source and the real seismic source located by the embodiment and the comparative example can be seen that the position of the seismic source positioned by the method of the present invention is closer to the real position of the seismic source, illustrating the seismic source positioning accuracy of the method of the present invention higher.

Claims (3)

1. A layered velocity location method of regional rock mass microseismic source, characterized in that the steps are as follows: ① According to the rock mass wave velocity of the rock mass area of the seismic source to be measured, the rock mass area is divided into different wave velocity layers. The wave velocity layer is recorded as the first wave velocity layer, the second wave velocity layer, ..., the kth wave velocity layer from one side to the other in the direction perpendicular to the wave velocity layer, and the rock mass wave velocity of each wave velocity layer is recorded as V ₁ , V ₂ ,...,V _k , establishing a three-dimensional Cartesian coordinate system with any point in the rock mass region as the coordinate origin, and calculating the analytical equation of the interface of each wave velocity layer according to the thickness and occurrence of each wave velocity layer; ② Install m sensors in the rock mass area, m≥5, each sensor forms a spatial network structure distribution, and sensors must be installed in the first wave velocity layer and the kth wave velocity layer, measure the spatial coordinates of each sensor, and place the first The spatial coordinates of the _i sensors are denoted as ( _xi , y, _zi ); ③ Carry out blasting experiments in the rock mass area, so that microseisms with known source locations occur in the rock mass area. The spatial coordinates of the known source are (x _b , y _b , z _b ), and the i-th sensor The take-off time of the P wave generated by the microseism is recorded as t _i , and the wave velocity layers through which the P wave generated by the microseism propagates from the source to the i-th sensor are numbered 1, 2,..., n. The wave velocity of the rock mass in the wave velocity layer is recorded as V _i1 , V _i2 ,...,V _in , and the constant items of the analytical equations of the interface of the wave velocity layers passing through in turn are recorded as D _i(1,2) , D _i(2,3) , ..., D _i(n-1,n) , D _i(1,2) represent the analytical equation of the P-wave generated by the microseism that propagates from the source to the interface between the No. 1 wave velocity layer and the No. 2 wave velocity layer that the i-th sensor passes through Constant term, D _i(2,3) ,..., the meaning represented by D _i(n-1,n) and so on; If the above-mentioned known seismic source and the i-th sensor are located in different wave velocity layers, then the theoretical value Δt _i of the elapsed time for the P-wave generated by the microseismic source at the above-mentioned known seismic source to propagate from the known seismic source to the i-th sensor is as follows: 1) as shown: ${\mathrm{<span\; class="notranslate">\; \&Delta;t</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&beta;</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}}{{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}}<span\; class="notranslate">\; (</span>\frac{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; Ax</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; By</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; Cz</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; |</span>}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; 1</span>}}}<span\; class="notranslate">\; +</span>\frac{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; 3</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; |</span>}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; 2</span>}}}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; ...</span>}<span\; class="notranslate">\; +</span>\frac{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; |</span>}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>}}<span\; class="notranslate">\; +</span>\frac{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; Ax</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; By</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; Cz</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; |</span>}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; no</span>}}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>$ In formula (1), α _i and β _i are spatial distance parameters, α _i = |Axi +By _i +Cz _{i -Ax b} -By _b -Cz _b |, A, B, and C are the coefficients of x, y, and z in the analytical equation of the interface of each wave velocity layer; If the above-mentioned known source and the i-th sensor are located in the same wave velocity layer, then the theoretical value Δt _i of the elapsed time for the P-wave generated by the microseismic source at the above-mentioned known source to propagate from the known source to the i-th sensor is as follows: 2) as shown: ${\mathrm{<span\; class="notranslate">\; \&Delta;t</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; =</span>\frac{\sqrt{{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; the\; y</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; the\; y</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; z</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; z</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; 1</span>}}}<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; )</span>$ Establish the objective function f for solving V ₁ , V ₂ ,..., V _k described in formula (3), where i≠j (3) In formula (3), t _i , t _j are the take-off moments of the P waves generated by the i-th and j-th sensors receiving microseisms with known source locations respectively; when solving V ₁ , V ₂ ,...,V _k In the process, it is necessary to _{replace V i1 ,} V _i2 ,...,V _{in in} formula (1) and formula ₍ 2) with V ₁ , V ₂ ,…,V _k , use nonlinear optimization tools to solve the minimum value of the objective function f described in formula (3), and take the rock mass wave velocity V ₁ , V ₂ ,…,V _k corresponding to the minimum value of f as Step 4. Solve the rock mass wave velocity of unknown source; ④ When a microseism with an unknown source location occurs in the rock mass area, it is assumed that the spatial coordinates of the unknown source are (x ₀ , y ₀ , z ₀ ), If the unknown source and the i-th sensor are located in different wave velocity layers, the theoretical value Δt _i of the elapsed time for the P-wave generated by the microseismic source at the unknown source to propagate from the unknown source to the i-th sensor is given by Eq. (4) Show: ${\mathrm{<span\; class="notranslate">\; \&Delta;t</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&beta;</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}}{{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}}<span\; class="notranslate">\; (</span>\frac{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; Ax</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; By</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; Cz</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; |</span>}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; 1</span>}}}<span\; class="notranslate">\; +</span>\frac{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; 3</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; |</span>}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; 2</span>}}}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; ...</span>}<span\; class="notranslate">\; +</span>\frac{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; |</span>}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>}}<span\; class="notranslate">\; +</span>\frac{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; Ax</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; By</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; Cz</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; |</span>}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; no</span>}}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 4</span>}<span\; class="notranslate">\; )</span>$ In formula (4), α _i and β _i are spatial distance parameters, α _i =|Axi +By _i +Cz _{i -Ax 0} -By ₀ -Cz ₀ |, A, B, and C are the coefficients of x, y, and z in the analytical equation of the interface of each wave velocity layer; If the above-mentioned unknown source and the i-th sensor are located in the same wave velocity layer, then the theoretical value Δt _i of the elapsed time for the P-wave generated by the microseismic source at the above-mentioned unknown source to propagate from the unknown source to the i-th sensor is given by Equation (5) Show: ${\mathrm{<span\; class="notranslate">\; \&Delta;t</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; =</span>\frac{\sqrt{{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; the\; y</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; the\; y</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; z</span>}}_{\mathrm{<span\; class="notranslate">\; b</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; z</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}}{{\mathrm{<span\; class="notranslate">\; V</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; 1</span>}}}<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 5</span>}<span\; class="notranslate">\; )</span>$ Establish the sub-velocity layer described in formula (6) to solve the objective function f _η of the unknown source space coordinates (x ₀ , y ₀ , z ₀ ), where i≠j (6) In formula (6), t _i , t _j are respectively the take-off moments of the P waves generated by the i-th and j-th sensors receiving the microseisms with unknown source locations, and η=1, 2,..., k, where k is the The total number of wave velocity layers in the rock mass area, in the process of solving the unknown source space coordinates (x ₀ , y ₀ , z ₀ ), it is necessary to base the correspondence between V _i1 , V _i2 ,…,V _in and the rock mass wave velocities relationship, replace V _i1 , V _i2 ,...,V _{in in} formula (4) and formula (5) with V ₁ , V ₂ ,...,V _k , and use nonlinear optimization tools to solve the objective described in formula (6) The minimum value of the function f _η , η takes 1, 2, ..., k, and the minimum value of the obtained f _η is recorded as f _1min , f _2min , ..., f _kmin respectively; then solve f _1min according to formula (7) ,f _2min ,…,f _kmin the minimum value f _θ , f _θ ＝min{f _1min ,f _2min ,...,f _kmin } (7) The (x ₀ , y ₀ , z ₀ ) corresponding to f _θ is the spatial coordinate of the unknown source; When using nonlinear optimization tools to solve the minimum value of f _η , the search for the feasible region Ω _η should satisfy the linear programming inequality described in formula (8): R _η ·[x ₀ ,y ₀ ,z ₀ ] ^T ≤b _η (8) In formula (8), A _η1 ,…,A _ηy ,B _η1 ,…,B _ηy ,C _η1 ,…,C _ηy ,D _η1 ,…,D _ηy are the parameters of the analytical equations of the boundary surfaces of the nth wave velocity layer respectively, and y is the The number of each boundary surface of the η wave velocity layer. 2. according to the layered velocity localization method of regional rock mass microseismic source of claim 1, it is characterized in that described nonlinear optimization tool is simplex method, particle swarm algorithm or genetic algorithm. 3. according to the layered velocity location method of regional rock mass microseismic source according to claim 1 and 2, it is characterized in that in step 3. a blasting test is carried out, or at least two blasting tests are carried out intermittently.