Coal mining whole process advanced detection method and system

Abstract

The invention provides a coal mining whole process advanced detection method and system, and belongs to the technical field of coal mining advanced detection. The method comprises the steps that a seismic source used for generating seismic signals is arranged in a coal mine area, and a detector is arranged in a coal seam; obtaining a seismic signal collected by the detector, and performing denoising processing on the seismic signal by using a wavelet threshold denoising method; constructing a forward modeling model by using a viscoelasticity forward modeling method based on the seismic signals subjected to denoising processing; performing inversion calculation on the forward model to generate an inversion speed model, and performing reverse time migration imaging on the inversion speed model to generate an imaging result; and analyzing and judging geological conditions of the coal mining area according to an imaging result, and carrying out early warning according to an analysis and judgment result. According to the invention, a high-precision imaging result of mining area geology is constructed by using a reverse time migration imaging technology, so that the accuracy and effectiveness of geological detection of coal mine engineering are comprehensively improved.

Description

Technical Field

[0001] The present invention relates to the technical field of advance detection in coal mining, and more particularly to a method and system for advance detection in the entire coal mining process. Background Art

[0002] Advance detection of the entire coal mining process is the core of coal mining, and its purpose is to ensure the safety and efficiency of coal mining. Traditional coal mining detection methods include slot wave seismic method, transient electromagnetic method, borehole inclinometer technology and other methods to ensure smooth construction. In terms of slot wave seismic exploration, it is necessary to further clarify the three-dimensional full wave field characteristics of advance detection in coal mines, carry out numerical simulation and physical simulation of viscoelastic wave fields based on half-space-full-space propagation, and establish a reliable basis for the inversion and interpretation of slot wave exploration data in coal mines and form specifications.

[0003] There are many uncertainties in the coal mining process, such as the hydrogeological conditions and coal seam complexity of the mine mining site, which poses many risks to mining. The full-process advance detection can reflect the working conditions with the help of detection parameters so that corresponding adjustments can be made. The main purpose is to detect coal seams and their structural characteristics, gas safety geological conditions, and hydrogeological and engineering geology. Slot wave seismic exploration is based on the classical seismic wave equation, but a lot of actual data evidence shows that the propagation of waves in actual media will show waveform attenuation, and the classical non-attenuated acoustic wave equation cannot explain this phenomenon. Summary of the invention

[0004] In view of the above problems, the purpose of the present invention is to provide a method and system for advanced detection of the entire coal mining process. Based on the viscoelastic acoustic wave equation, it can better fit the real physical process of wave propagation, thereby improving the inversion imaging accuracy in geological detection. At the same time, an early warning system is established to timely report the geological conditions of coal seams and reduce the impact of adverse construction factors.

[0005] In order to achieve the above object, the present invention is implemented through the following technical solutions: In a first aspect, the present invention discloses a method for advanced detection of the entire coal mining process, comprising: Arrange a seismic source for generating seismic signals in the coal mine area and arrange geophones in the coal seam; Acquire the seismic signal collected by the geophone, and use the wavelet threshold denoising method to denoise the seismic signal; Based on the de-noised seismic signal, a forward model is constructed using the viscoelastic forward modeling method. Perform inversion calculation on the forward model to generate an inversion velocity model, perform reverse time migration imaging on the inversion velocity model to generate imaging results; The geological conditions of the coal mining area are analyzed and judged based on the imaging results, and early warnings are issued based on the analysis and judgment results.

[0006] Further, the method of arranging a seismic source for generating seismic signals in the coal mine area and arranging a detector in the coal seam includes: A seismic source is arranged every 1.5m on the side walls of the coal mine tunnel, and three seismic sources are arranged on each side wall; The detectors are arranged at 3m away from the working face of the coal seam, and 6 detectors are arranged on each side wall; the detectors are installed in the boreholes, the borehole diameter is 50mm, the borehole depth is 1.2m, and the borehole interval is 1.5m.

[0007] Further, the denoising process of the seismic signal using the wavelet threshold denoising method includes: The wavelet decomposition formula is used to perform wavelet decomposition, and the seismic signal is divided into a low-frequency part and a high-frequency part, wherein the high-frequency part contains noise information; Calculate the soft threshold according to the number of wavelet decomposition layers; Use the soft threshold denoising formula to compress the high-frequency part of the seismic signal; The denoised seismic signal is generated using the wavelet reconstruction formula.

[0008] Furthermore, the forward model is constructed based on the de-noised seismic signal using a viscoelastic forward modeling method, including: Construct the first-order stress-velocity equations for two-dimensional fractional viscoelastic acoustic waves; The short memory algorithm is used to transform the first-order stress-velocity equations of two-dimensional fractional viscoelastic acoustic waves into an equivalent problem that is non-local in time and local in space. The transformed equivalent problem is numerically solved using a finite difference method to generate a forward model; The boundary of the calculation region is processed using the auxiliary differential equation perfectly matched layer method to obtain the absorbing boundary region equation group.

[0009] Furthermore, the first-order stress-velocity equations of the two-dimensional fractional viscoelastic acoustic wave are specifically as follows:

[0010]

[0011]

[0012] The stress-strain relationship can be expressed using fractional order operators as:

[0013] in, , , is the fractional operator form of the elastic modulus, For frequency , , and is the strain component, , and is the stress component, is the bulk elastic modulus, The order is The time fractional operator of is the fractional order, is the quality factor, is the Euler gamma function, is the density, , is the Lamé constant; The method uses a short memory algorithm to transform the two-dimensional fractional viscoelastic acoustic wave first-order stress-velocity equations into an equivalent problem that is non-local in time and local in space, including: Will Flat exotic core Transformed into relatively localized exotic nuclei , that is,

[0014] The first-order stress-velocity equations of two-dimensional fractional viscoelastic acoustic waves are transformed into an equivalent problem that is non-local in time and local in space. Furthermore, the finite difference method is used to numerically solve the transformed equivalent problem to generate a forward model, including: The spatial derivative discretization is based on the second-order central finite difference method, including: Assume that at the grid point The time steps are , , ;

[0015]

[0016]

[0017] The time fractional derivative is discretized based on the generalized Gauss-Laguerre quadrature formula, and the time step iteration is performed, including: Define the helper function:

[0018] And satisfy:

[0019] Further we get:

[0020] Based on the generalized Gauss-Laguerre numerical quadrature formula, the discrete format of the fractional-order operator is obtained: So for the original integral we have:

[0021] in, , , is the weight, is the attenuation parameter; This results in the time iteration format of the viscoelastic acoustic wave equation:

[0022] Based on the time iteration format, a forward model is generated to obtain wavelength data.

[0023] Furthermore, the auxiliary differential equation perfectly matched layer method is used to process the boundary of the calculation area to obtain the absorbing boundary area equation group, including: transforming the time iteration format of the viscoelastic acoustic wave equation into the frequency domain; Using the calculation formula of the stretched coordinate system, the wave propagates outside the calculation area and automatically decays, and then transforms back to the time domain; The exponential integrator is used to solve the stiff auxiliary differential equation, and finally the absorbing boundary region equations are obtained. Further, the inversion calculation is performed on the forward model to generate an inversion velocity model, and reverse time migration imaging is performed on the inversion velocity model to generate an imaging result, including: The forward model is inverted using full waveform inversion based on quasi-Newton method optimization iteration to generate an inversion velocity model. The full waveform inversion using the optimization iteration based on the quasi-Newton method includes: The objective function is established based on the data obtained from simulation and observation: Construct the relationship between the objective function and the model parameters: The gradient of the objective function is constructed by convolution of the forward and reverse wave fields; The velocity model is updated through the full optimization iteration based on the local optimal algorithm of the derivative; After the iterative update of the velocity model is completed, the inversion velocity model is generated.

[0024] Further, performing reverse time migration imaging on the inversion velocity model includes: Prestack wavefield separation and reverse time migration are used to optimize imaging, specifically: By calculating the zero-delay cross-correlation between the source wavefield and the receiving wavefield, the following imaging formula is obtained:

[0025]

[0026]

[0027] in, is the cross-correlation value between the downgoing source wavefield and the upgoing receiving wavefield; is the cross-correlation value between the upgoing source wavefield and the downgoing receiving wavefield; , are the cross-correlation values ​​of the source wavefield and the received wavefield along the same downlink and uplink propagation directions, respectively; , They are the downgoing and upgoing components of the source wavefield, respectively; , They are the downlink and uplink components of the receiving wave field respectively; , , denote the space and time coordinates respectively, It is the maximum time of TRON extrapolation.

[0028] In a second aspect, the present invention also discloses a system for advanced detection of the entire coal mining process, comprising: A seismic source arrangement module configured to arrange seismic sources for generating seismic signals in a coal mine area and to arrange geophones in a coal seam; A signal acquisition and preprocessing module is configured to acquire seismic signals acquired by the geophone and denoise the seismic signals using a wavelet threshold denoising method; A forward modeling module configured to construct a forward model using a viscoelastic forward modeling method based on the de-noised seismic signal; An inversion imaging module is configured to perform inversion calculation on the forward model, generate an inversion velocity model, perform reverse time migration imaging on the inversion velocity model, and generate imaging results; The early warning module is configured to analyze and identify the geological conditions of the coal mining area according to the imaging results, and to issue an early warning according to the analysis and identification results.

[0029] Compared with the prior art, the beneficial effects of the present invention are as follows: the present invention discloses a method and system for advanced detection of the entire coal mining process. First, the layout of the seismic source and sensors is determined according to the construction environment of the mining area; the data collected by the sensor is transmitted to the host and denoised. Subsequently, a forward simulation scheme for viscoelastic waves is designed, and full waveform inversion and reverse time migration imaging techniques are used emphatically. Through a full waveform inversion method, the data is deeply analyzed, the propagation process of the wave field in the medium is accurately restored, and static geological conditions such as partial faults, collapse columns, and coal seam thinning areas are finely detected, and the roof fracture zone, stress concentration area, and protruding danger area are accurately identified. A reverse time migration imaging technology is used to construct high-precision imaging results of the mining area geology, thereby comprehensively improving the accuracy and effectiveness of coal mine engineering geological detection.

[0030] It can be seen that compared with the prior art, the present invention has outstanding substantive features and significant progress, and the beneficial effects of its implementation are also obvious. BRIEF DESCRIPTION OF THE DRAWINGS

[0031] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on the provided drawings without paying creative work.

[0032] Figure 1 It is a method flow chart of a specific implementation mode of the present invention.

[0033] Figure 2 It is a schematic diagram of the arrangement of the seismic source and the detector according to a specific embodiment of the present invention.

[0034] Figure 3 It is a system structure diagram of a specific implementation mode of the present invention. DETAILED DESCRIPTION

[0035] In order to enable those skilled in the art to better understand the scheme of the present invention, the present invention is further described in detail below in conjunction with the accompanying drawings and specific implementation methods. Obviously, the described embodiments are only part of the embodiments of the present invention, rather than all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention.

[0036] See also Figure 1 As shown, this embodiment provides a method for advanced detection of the entire coal mining process, comprising the following steps: S1: Arrange a seismic source for generating seismic signals in the coal mine area and arrange detectors in the coal seam.

[0037] In a specific embodiment, Figure 2 As shown, seismic sources 2 are arranged on the side walls of the coal mine tunnel 1, one seismic source 2 is arranged every 1.5m, and three seismic sources 2 are arranged on each side wall to ensure signal coverage and detection depth. Geophones 3 are arranged about 3m away from the coal seam working face, the borehole interval is 1.5m, and six geophones 3 are arranged on each side wall to increase the density of signal reception.

[0038] The borehole diameter was kept at 50 mm and the depth was adjusted to 1.2 m to suit the coal seam thickness and environment. The minimum distance between the seismic source and the adjacent detector was 10 m and the width of the coal mine working face was 12 m to cover the entire detection area.

[0039] S2: Obtain the seismic signal collected by the detector and denoise the seismic signal using the wavelet threshold denoising method.

[0040] In a specific implementation, when acquiring the seismic signal collected by the detector, the seismic signal reception time is set to 0.15s to capture more echo information, the sampling frequency is maintained at 8000Hz to improve the resolution, and each collected data is automatically stored until all seismic sources are excited.

[0041] After the data is stored, the wavelet threshold denoising method is used for denoising, which includes: First, the wavelet decomposition formula is used to perform wavelet decomposition, dividing the seismic signal into a low-frequency part and a high-frequency part, where the high-frequency part contains noise information. Then, the soft threshold is calculated according to the number of wavelet decomposition layers, and the soft threshold denoising formula is used to compress the seismic signal in the high-frequency part. Finally, the wavelet reconstruction formula is used to generate the denoised seismic signal.

[0042] S3: Based on the denoised seismic signal, a forward model is constructed using the viscoelastic forward modeling method.

[0043] In a specific implementation, first, a two-dimensional fractional viscoelastic acoustic wave first-order stress velocity equation group is constructed. Then, the two-dimensional fractional viscoelastic acoustic wave first-order stress velocity equation group is transformed into a time non-local and space local equivalent problem using a short memory algorithm. At this time, the transformed equivalent problem is numerically solved using a finite difference method to generate a forward model. Finally, the boundary of the calculation area is processed using an auxiliary differential equation perfectly matched layer method to obtain an absorbing boundary area equation group.

[0044] Among them, the auxiliary differential equation perfectly matched layer method is used to process the boundary of the calculation area to obtain the absorbing boundary area equation group, which specifically includes: transforming the time iteration format of the viscoelastic acoustic wave equation into the frequency domain; Using the calculation formula of the stretched coordinate system, the wave propagates outside the calculation area and automatically decays, and then transforms back to the time domain; The exponential integrator is used to solve the stiff auxiliary differential equation, and finally the absorbing boundary region equations are obtained.

[0045] S4: Perform inversion calculation on the forward model to generate an inversion velocity model, perform reverse time migration imaging on the inversion velocity model, and generate imaging results.

[0046] In a specific implementation, the forward model is inverted using full waveform inversion based on quasi-Newton optimization iteration to generate an inversion velocity model. First, an objective function is established based on the data obtained by simulation and observation. Then, the relationship between the objective function and the model parameters is constructed, and the gradient of the objective function is constructed by the convolution of the forward and reverse wave fields. At this time, the velocity model is updated based on the local optimal algorithm of the derivative through full optimization iteration; after the iterative update of the velocity model is completed, the inversion velocity model is generated.

[0047] S5: Analyze and identify the geological conditions of the coal mining area based on the imaging results, and issue early warnings based on the analysis and identification results.

[0048] In a specific implementation method, the geological conditions of the coal mining area are systematically and meticulously analyzed and identified based on the imaging results, and various key geological features and unfavorable factors in the strata are accurately displayed. In particular, potential geological structures such as hidden faults, drastic changes in coal seam thickness, gas-rich zones, karst collapse zones, and weak surrounding rocks can be warned in a timely manner. Based on the results of model identification, coal mine management departments and engineering personnel can take preventive measures in a timely manner, including laying out support and reinforcement facilities in dangerous areas, taking gas extraction or ventilation measures, and strengthening drainage systems, so as to effectively avoid and control the risks of geological disasters.

[0049] In an embodiment of the present invention, based on step S2, a possible embodiment is given below to illustrate its specific implementation scheme in a non-limiting manner.

[0050] S201: Use the wavelet decomposition formula above to perform wavelet decomposition. Select four wavelet functions, Haar, Daubechies, Coiflets, and Symlets, to perform 32-layer wavelet decomposition on the seismic data. First, use the step signal to approximate the original signal:

[0051]

[0052]

[0053] In the formula, is the approximate coefficient, is a scaling function that roughly represents the original signal, is the set of scaling functions The space created; Will This breaks down into:

[0054] in:

[0055]

[0056]

[0057]

[0058]

[0059]

[0060] Finally, we can get:

[0061] In the formula, is the fine coefficient, To correct the difference between the scale function representation and the original signal, Represents a set of wavelet functions The space created.

[0062] S202: Calculate the threshold value, the calculation formula is:

[0063] Where N is the number of decomposition levels.

[0064] Calculated .

[0065] S203: Perform denoising using a soft threshold denoising formula:

[0066] in is the wavelet coefficient.

[0067] S204: Perform wavelet reconstruction:

[0068]

[0069] Finally, the denoised signal is obtained using the wavelet reconstruction formula mentioned above.

[0070] In an embodiment of the present invention, based on step S3, a possible embodiment is given below to illustrate its specific implementation scheme in a non-limiting manner.

[0071] In step S3, the detection area and the locations of the source and the detector are constructed in the numerical model, and the source is set as the Ricker wavelet, and the formula is:

[0072] in, is the amplitude of the Ricker wavelet signal, The main frequency of the earthquake source is According to the sampling frequency, the time interval dt is:

[0073] The constructed viscoelastic acoustic wave equations are numerically solved using the short-memory operator format and the finite difference method, and the auxiliary differential equation perfectly matched layer (ADE-Pml) is used to solve the absorbing boundary.

[0074] The specific process of this embodiment includes: S301: Two-dimensional fractional viscoelastic acoustic wave first-order stress-velocity equations:

[0075]

[0076]

[0077] The stress-strain relationship can be expressed using fractional-order operators as follows:

[0078]

[0079]

[0080]

[0081]

[0082] Where C is the fractional operator form of the elastic modulus, is the frequency, , and is the strain component, , and is the stress component, is the bulk elastic modulus, The order is The time fractional operator of is the fractional order, is the quality factor, is the Euler gamma function, is the density, , is the Lamé constant.

[0083] S302: Flat exotic core Transformed into relatively localized exotic nuclei , that is,

[0084] Convert the equation into an equivalent problem that is non-local in time and local in space.

[0085] S303: Solve the above problem by using the finite difference method: The spatial derivatives are discretized based on the second-order central finite difference method. Approximate the spatial derivative using the central difference method in space: Assume that at the grid point The time steps are , , .

[0086]

[0087]

[0088]

[0089] In terms of time, the time fractional derivative is discretized based on the generalized Gauss-Laguerre quadrature formula, and the time step iteration is performed: Define the helper function:

[0090] And satisfy:

[0091] Further we get:

[0092] Based on the generalized Gauss-Laguerre numerical quadrature formula, the discrete format of the fractional-order operator is obtained: So for the original integral we have:

[0093] in , , is the weight, is the attenuation parameter.

[0094] This results in the time iteration format of the viscoelastic acoustic wave equation:

[0095] in, is the second-order format Strong split differential operator.

[0096] Based on the time iteration format, a forward model is generated to obtain wavelength data.

[0097] S304: Use the Perfect matched layer with auxiliary differential equation (ADE-Pml) method to process the boundary and obtain the absorbing boundary region equations. The steps are as follows: The viscoelastic acoustic wave equation is transformed into the frequency domain, and then the stretching coordinates are given so that the wave propagates outside the calculation area and automatically attenuates, and then transformed back to the time domain, thereby avoiding the contamination of the wave field data in the calculation area by pseudo wave reflection. The relevant calculation formula for the given stretching coordinate system is:

[0098] In the above formula, is the stretching coordinate conversion factor is the frequency shift term of the coordinate stretching factor is the proportional term of the coordinate stretching factor, is the angular frequency of the sound wave signal, is the loss term of the coordinate stretching factor; in, represents the stretched coordinate system, and in addition:

[0099]

[0100]

[0101] Where L is the thickness of the PML layer, , i is the distance from the current calculation point to the boundary, is the maximum loss coefficient, is the velocity of the pressure wave, R is the target theoretical reflection coefficient, which is selected as 0.05% here. Take p=1, is the main frequency of the earthquake source. Usually between 1 and 20, we take m=3.

[0102] In an embodiment of the present invention, based on step S4, a possible embodiment is given below to illustrate its specific implementation scheme in a non-limiting manner.

[0103] In this embodiment, the initial velocity model matrix is ​​first imported, and a forward simulation is performed in this model. The total duration of the forward simulation is 0.2 seconds, and the time nodes are from 0 seconds to 0.2 seconds. Wave field data and seismic records are recorded during the simulation. The generated simulated seismic records are interpolated with the actual geophone data stored in advance to obtain residual data.

[0104] Next, the residual data is used to perform the reverse wavefield simulation in the wavefield. The reverse simulation duration is also 0.2 seconds, and the time nodes are 0.2 seconds to 0 seconds. The wavefield data of the reverse wavefield is recorded for subsequent processing. Finally, the forward wavefield data and the reverse wavefield data are correspondingly multiplied in two-dimensional space at each time point, and the multiplication results of all time nodes are gradually accumulated to complete the cross-correlation operation.

[0105] The specific process of this embodiment is as follows: S401: Establish an objective function based on the data obtained from simulation and observation:

[0106] In the formula, the objective function is a scalar with units of , and denote the simulated and observed Seismogram on the component, and represent the source and the detector respectively, is the model parameter vector, such as velocity, density, anisotropy, absorption attenuation parameters, etc. , , denote the number of earthquakes, geophones and time windows respectively, S402: Construct the relationship between the objective function and the model parameters as follows:

[0107] In the formula, is the change of the objective function, , and represent the speed and density of P-wave and S-wave respectively, , and are the gradients of the objective function relative to the above three model parameters.

[0108] S403: Construct the gradient of the objective function through the convolution of the forward and reverse wave fields, and select the density , bulk modulus and shear modulus As model parameters, their gradient formulas are as follows:

[0109]

[0110]

[0111] In the formula, and They are the main and accompanying wave fields. and is the normal and adjoint strain deviator, It’s record time.

[0112] S404: The velocity model is updated through full optimization iteration based on the local optimal algorithm of the derivative, so that the error calculation between the continuously updated velocity model and the actual model is less than the allowable error. , obtain the local optimal solution near the initial model, the calculation formula is:

[0113] in, is the initial model, For the actual model.

[0114] S405: The full optimization iteration uses the quasi-Newton method. The error functional satisfies the quadratic form of the model change near the initial model and is expanded using the Taylor formula:

[0115] In the formula, , is the error functional derived from the model in the initial model The value at Let the error functional be The second-order derivative of the model, namely the Hessian operator

[0116] In the formula, is its nonlinear term.

[0117] Create an iterative formula:

[0118] Since the Hessian operator requires a lot of computation, the quasi-Newton method uses an approximate Hessian matrix inverse method to reduce the computation.

[0119]

[0120] In the formula, is the iteration step length, is the inverse of the approximate Hessian matrix.

[0121] Given allowable error =0.01, if Then continue to calculate the gradient J and update the model parameters until , that is, the inversion velocity model is obtained .

[0122] S406: Using the above obtained inversion velocity model For reverse time migration imaging, the present invention preferably uses prestack wavefield separation reverse time migration to suppress noise in low-frequency imaging and make imaging more accurate. The separation reverse time migration imaging condition is a cross-correlation imaging condition, and the source wavefield and the receiving wavefield include:

[0123]

[0124] In the formula, , They are the downgoing and upgoing components of the source wavefield, respectively; , They are respectively the downlink and uplink components of the receiving wave field.

[0125] By calculating the zero-delay cross-correlation between the source wavefield and the receiving wavefield, we can obtain:

[0126]

[0127]

[0128] In the formula, It is the cross-correlation value between the downgoing source wavefield and the upgoing receiving wavefield, i.e., the cross-correlation imaging of the one-way wave equation migration; is the cross-correlation value between the upgoing source wavefield and the downgoing received wavefield, i.e., the prismatic wave ray path calculated by ray tracing; , are the cross-correlation values ​​of the source wavefield and the received wavefield along the same downlink and uplink propagation directions, respectively; , , denote the space and time coordinates respectively, It is the maximum time of TRON extrapolation.

[0129] Since the latter two cross-correlations are performed on the entire propagation path, they are low-frequency imaging noise in the profile, that is, the first two items can be imaged correctly:

[0130] The present invention provides a method for advanced detection of the entire coal mining process, which utilizes a reverse time migration imaging technology to construct high-precision imaging results of the mining area geology, thereby comprehensively improving the accuracy and effectiveness of coal mine engineering geological detection.

[0131] See also Figure 3 As shown, the present invention also discloses a whole-process advance detection system for coal mining, including: a source arrangement module 101, a signal acquisition and preprocessing module 102, a forward modeling module 103, an inversion imaging module 104 and an early warning module 105. The seismic source arrangement module 101 is configured to arrange seismic sources for generating seismic signals in the coal mine area and to arrange detectors in the coal seam.

[0132] The signal acquisition and preprocessing module 102 is configured to acquire the seismic signal acquired by the detector and perform denoising on the seismic signal using a wavelet threshold denoising method.

[0133] The forward modeling module 103 is configured to construct a forward model using a viscoelastic forward modeling method based on the de-noised seismic signal.

[0134] The inversion imaging module 104 is configured to perform inversion calculation on the forward model to generate an inversion velocity model, perform reverse time migration imaging on the inversion velocity model, and generate an imaging result.

[0135] The early warning module 105 is configured to analyze and identify the geological conditions of the coal mining area according to the imaging results, and issue an early warning according to the analysis and identification results.

[0136] The specific implementation of the system for advanced detection of the entire coal mining process of this embodiment is basically consistent with the specific implementation of the method for advanced detection of the entire coal mining process described above, and will not be repeated here.

[0137] In this specification, each embodiment is described in a progressive manner, and each embodiment focuses on the differences from other embodiments. The same or similar parts between the embodiments can be referred to each other. As for the method disclosed in the embodiment, since it corresponds to the system disclosed in the embodiment, the description is relatively simple, and the relevant parts can be referred to the method part.

[0138] Professionals may further appreciate that the units and algorithm steps of each example described in conjunction with the embodiments disclosed herein can be implemented in electronic hardware, computer software, or a combination of the two. In order to clearly illustrate the interchangeability of hardware and software, the composition and steps of each example have been generally described in the above description according to function. Whether these functions are performed in hardware or software depends on the specific application and design constraints of the technical solution. Professionals and technicians may use different methods to implement the described functions for each specific application, but such implementation should not be considered to be beyond the scope of the present invention.

[0139] In the several embodiments provided by the present invention, it should be understood that the disclosed systems, systems and methods can be implemented in other ways. For example, the system embodiments described above are only schematic. For example, the division of the units is only a logical function division. There may be other division methods in actual implementation, such as multiple units or components can be combined or integrated into another system, or some features can be ignored or not executed. Another point is that the mutual coupling or direct coupling or communication connection shown or discussed can be through some interfaces, indirect coupling or communication connection of systems or units, which can be electrical, mechanical or other forms.

[0140] The units described as separate components may or may not be physically separated, and the components shown as units may or may not be physical units, that is, they may be located in one place or distributed on multiple network units. Some or all of the units may be selected according to actual needs to achieve the purpose of the solution of this embodiment.

[0141] In addition, each functional module in each embodiment of the present invention may be integrated into one processing unit, or each module may exist physically separately, or two or more modules may be integrated into one unit.

[0142] Similarly, each processing unit in each embodiment of the present invention may be integrated into one functional module, or each processing unit may exist physically, or two or more processing units may be integrated into one functional module.

[0143] The steps of the method or algorithm described in conjunction with the embodiments disclosed herein may be implemented directly using hardware, a software module executed by a processor, or a combination of the two. The software module may be placed in a random access memory (RAM), a memory, a read-only memory (ROM), an electrically programmable ROM, an electrically erasable programmable ROM, a register, a hard disk, a removable disk, a CD-ROM, or any other form of storage medium known in the art.

[0144] Finally, it should be noted that, in this article, relational terms such as first and second, etc. are only used to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply any such actual relationship or order between these entities or operations. Moreover, the terms "include", "comprise" or any other variants thereof are intended to cover non-exclusive inclusion, so that a process, method, article or device including a series of elements includes not only those elements, but also other elements not explicitly listed, or also includes elements inherent to such process, method, article or device. In the absence of further restrictions, the elements defined by the sentence "comprise a ..." do not exclude the presence of other identical elements in the process, method, article or device including the elements.

[0145] The above is a detailed introduction to the method for advanced detection of the entire coal mining process provided by the present invention. Specific examples are used herein to illustrate the principles and implementation methods of the present invention. The description of the above embodiments is only used to help understand the method of the present invention and its core idea. It should be pointed out that for ordinary technicians in this technical field, without departing from the principles of the present invention, the present invention can also be improved and modified, and these improvements and modifications also fall within the scope of protection of the present invention.

Claims

1. A method for advanced detection of the entire coal mining process, characterized in that: include: Arrange a seismic source for generating seismic signals in the coal mine area and arrange geophones in the coal seam; Acquire the seismic signal collected by the geophone, and use the wavelet threshold denoising method to denoise the seismic signal; Based on the de-noised seismic signal, a forward model is constructed using the viscoelastic forward modeling method. Perform inversion calculation on the forward model to generate an inversion velocity model, perform reverse time migration imaging on the inversion velocity model to generate imaging results; The geological conditions of the coal mining area are analyzed and judged based on the imaging results, and early warnings are issued based on the analysis and judgment results.

2. The method for advanced detection of the entire coal mining process according to claim 1 is characterized in that: The method of arranging a seismic source for generating seismic signals in a coal mine area and arranging a geophone in a coal seam includes: A seismic source is arranged every 1.5m on the side walls of the coal mine tunnel, and three seismic sources are arranged on each side wall; The detectors are arranged at 3m away from the working face of the coal seam, and 6 detectors are arranged on each side wall; the detectors are installed in the boreholes, the borehole diameter is 50mm, the borehole depth is 1.2m, and the borehole interval is 1.5m.

3. The method for advanced detection of the entire coal mining process according to claim 1 is characterized in that: The method of using the wavelet threshold denoising method to denoise the seismic signal comprises: The wavelet decomposition formula is used to perform wavelet decomposition, and the seismic signal is divided into a low-frequency part and a high-frequency part, wherein the high-frequency part contains noise information; Calculate the soft threshold according to the number of wavelet decomposition layers; Use the soft threshold denoising formula to compress the high-frequency part of the seismic signal; The denoised seismic signal is generated using the wavelet reconstruction formula.

4. The method for advanced detection of the entire coal mining process according to claim 3 is characterized in that: The method of constructing a forward model based on the de-noised seismic signal using a viscoelastic forward modeling method comprises: Construct the first-order stress-velocity equations for two-dimensional fractional viscoelastic acoustic waves; The short memory algorithm is used to transform the first-order stress-velocity equations of two-dimensional fractional viscoelastic acoustic waves into an equivalent problem that is non-local in time and local in space. The transformed equivalent problem is numerically solved using a finite difference method to generate a forward model; The boundary of the calculation region is processed using the auxiliary differential equation perfectly matched layer method to obtain the absorbing boundary region equation group.

5. The method for advanced detection of the entire coal mining process according to claim 4 is characterized in that: Said The first-order stress-velocity equations of two-dimensional fractional viscoelastic acoustic waves are as follows: Among them, The stress-strain relationship can be expressed using fractional order operators as: in, , , is the fractional operator form of the elastic modulus, For frequency , , and is the strain component, , and is the stress component, is the bulk elastic modulus, The order is The time fractional operator of is the fractional order, is the quality factor, is the Euler gamma function, is the density, , is the Lamé constant; The method uses a short memory algorithm to transform the two-dimensional fractional viscoelastic acoustic wave first-order stress-velocity equations into an equivalent problem that is non-local in time and local in space, including: Will Flat exotic core Transformed into relatively localized exotic nuclei , that is, The first-order stress-velocity equations of two-dimensional fractional viscoelastic acoustic waves are transformed into an equivalent problem that is non-local in time and local in space.

6. The method for advanced detection of the entire coal mining process according to claim 5 is characterized in that: The method of using the finite difference method to numerically solve the converted equivalent problem to generate a forward model includes: The spatial derivative discretization is based on the second-order central finite difference method, including: Assume that at the grid point The time steps are , , ; The time fractional derivative is discretized based on the generalized Gauss-Laguerre quadrature formula, and the time step iteration is performed, including: Define the helper function: And satisfy: Further we get: Based on the generalized Gauss-Laguerre numerical quadrature formula, the discrete format of the fractional-order operator is obtained: So for the original integral we have: in, , , is the weight, is the attenuation parameter; This results in the time iteration format of the viscoelastic acoustic wave equation: in, is the second-order format Strong split differential operator; Based on the time iteration format, a forward model is generated to obtain wavelength data.

7. The method for advanced detection of the entire coal mining process according to claim 6, characterized in that: The auxiliary differential equation perfectly matched layer method is used to process the boundary of the calculation area to obtain the absorbing boundary area equation group, including: transforming the time iteration format of the viscoelastic acoustic wave equation into the frequency domain; Using the calculation formula of the stretched coordinate system, the wave propagates outside the calculation area and automatically decays, and then transforms back to the time domain; The exponential integrator is used to solve the stiff auxiliary differential equation, and finally the absorbing boundary region equations are obtained.

8. The method for advanced detection of the entire coal mining process according to claim 7, characterized in that: The method of performing inversion calculation on the forward model to generate an inversion velocity model, performing reverse time migration imaging on the inversion velocity model to generate an imaging result includes: The forward model is inverted using full waveform inversion based on quasi-Newton method optimization iteration to generate an inversion velocity model. The full waveform inversion using the optimization iteration based on the quasi-Newton method includes: The objective function is established based on the data obtained from simulation and observation: Construct the relationship between the objective function and the model parameters: The gradient of the objective function is constructed by convolution of the forward and reverse wave fields; The velocity model is updated through the full optimization iteration based on the local optimal algorithm of the derivative; After the iterative update of the velocity model is completed, the inversion velocity model is generated.

9. The method for advanced detection of the entire coal mining process according to claim 8, characterized in that: The performing reverse time migration imaging on the inversion velocity model comprises: Prestack wavefield separation and reverse time migration are used to optimize imaging, specifically: By calculating the zero-delay cross-correlation between the source wavefield and the receiving wavefield, the following imaging formula is obtained: in, is the cross-correlation value between the downgoing source wavefield and the upgoing receiving wavefield; is the cross-correlation value between the upgoing source wavefield and the downgoing receiving wavefield; , are the cross-correlation values ​​of the source wavefield and the received wavefield along the same downlink and uplink propagation directions, respectively; , They are the downgoing and upgoing components of the source wavefield, respectively; , They are the downlink and uplink components of the receiving wave field respectively; , , denote the space and time coordinates respectively, It is the maximum time of TRON extrapolation.

10. A coal mining whole process advance detection system, characterized in that: The system adopts the method for advanced detection of the entire coal mining process as described in any one of claims 1 to 9; The system comprises: A seismic source arrangement module configured to arrange seismic sources for generating seismic signals in a coal mine area and to arrange geophones in a coal seam; A signal acquisition and preprocessing module is configured to acquire seismic signals acquired by the geophone and denoise the seismic signals using a wavelet threshold denoising method; A forward modeling module configured to construct a forward model using a viscoelastic forward modeling method based on the de-noised seismic signal; An inversion imaging module is configured to perform inversion calculation on the forward model, generate an inversion velocity model, perform reverse time migration imaging on the inversion velocity model, and generate imaging results; The early warning module is configured to analyze and identify the geological conditions of the coal mining area according to the imaging results, and to issue an early warning according to the analysis and identification results.

Images