A multi-method collaborative detection method and system for karst collapse hidden danger area

Abstract

The application discloses a karst collapse hidden danger area multi-method cooperative detection method and system, relates to the technical field of geographic detection, and comprises the following steps: arranging an observation network in a target hidden danger area, respectively implementing an electric excitation sequence and a mechanical vibration excitation sequence, collecting potential difference data and vibration waveform data, inverting the potential difference data and the vibration waveform data at each collection time, analyzing the phase relationship and the change form of the resistivity sequence and the elastic wave velocity sequence with respect to the corresponding excitation time sequence, and dividing the area into regions with different response characteristics according to the characteristic differences of the multi-physical field dynamic response modes in space. The application can extract the multi-physical field dynamic response mode representing the dynamics process such as groundwater seepage and particle migration, realize direct description of the dynamics process activity and evolution time sequence in the hidden danger area, and improve the ability of identifying the hidden danger development stage from the physical response.

Description

Technical Field

[0001] This invention relates to the field of geographic exploration technology, and in particular to a multi-method collaborative detection method and system for karst collapse hazard areas. Background Technology

[0002] In the field of geological hazard investigation in karst areas, detection operations mainly rely on single physical field methods. However, the occurrence and development of karst collapse is a complex dynamic problem involving the coupling of multiple processes such as groundwater seepage, chemical dissolution, particle migration, and mechanical instability. Existing detection technologies based on single or static multi-method fusion lack dynamic correlation between the acquired data under a unified time reference, resulting in the inability to effectively capture and distinguish the co-evolutionary characteristics of geophysical fields stimulated by different dynamic processes.

[0003] Therefore, existing technologies are unable to effectively identify dynamic information directly related to the activity level and dominant causal mechanism of hidden dangers from the detection response, resulting in the assessment of collapse risk mostly remaining at the level of structural identification, while the ability to judge the development stage and evolution trend of hidden dangers is insufficient. Summary of the Invention

[0004] This invention provides a multi-method collaborative detection method and system for karst collapse hazard areas to solve the problems mentioned in the background art.

[0005] To achieve the above objectives, this invention provides a multi-method collaborative detection method for karst collapse hazard areas, comprising:

[0006] S1. Set up an observation network in the target hazard area and implement electrical excitation sequence and mechanical vibration excitation sequence respectively, wherein each excitation sequence contains at least two excitations with different parameter settings;

[0007] S2. For each excitation in the electrical excitation sequence, collect potential difference data at a series of moments before, during, and after its application; for each excitation in the mechanical vibration excitation sequence, collect vibration waveform data at a series of moments before, during, and after its application.

[0008] S3. Invert the potential difference data at each acquisition time to obtain the resistivity change sequence over time; invert the vibration waveform data at each acquisition time to obtain the elastic wave velocity change sequence over time.

[0009] S4. Analyze the phase relationship and change pattern of the resistivity sequence and elastic wave velocity sequence relative to their corresponding excitation time sequence, and extract the multiphysics dynamic response mode;

[0010] S5. Based on the spatial differences in the dynamic response modes of multiphysics fields, regions with different response characteristics are divided; and based on the characteristic combinations of response modes in each region, geological genetic correlation and hazard type determination are carried out to generate a hazard zoning map.

[0011] Preferably, the implementation of the electrical excitation sequence includes:

[0012] Based on the geological information of the target hazard area, two current intensity parameters were set for the power supply electrodes, one for shallow large pore response and the other for deep small fracture response.

[0013] Based on geological information, two pulse duty cycle parameters were set for the power supply electrodes to obtain steady-state background field and rapid dynamic transition information, respectively.

[0014] According to a preset timing sequence, at least two different excitations consisting of two current intensity parameters and two pulse duty cycle parameters are alternately applied to the power supply electrodes to obtain an electrical excitation sequence.

[0015] Preferably, the implementation of the mechanical vibration excitation sequence includes:

[0016] Based on the geological information of the target hazard area, two frequency parameters were set for the excitation source to enhance the sensitivity to low-density filling media and high-density intact bedrock, respectively.

[0017] Based on geological information, two amplitude parameters were set for the excitation source to detect the structures in the near-field strong attenuation region and the far-field weak attenuation region, respectively.

[0018] According to a preset timing sequence, at least two different excitations consisting of two frequency parameters and two amplitude parameters are alternately applied to the excitation source to obtain a mechanical vibration excitation sequence.

[0019] Preferably, step S1 further includes:

[0020] The excitations in the electrical excitation sequence and the excitations in the mechanical vibration excitation sequence are arranged in time according to a preset alternation order;

[0021] This ensures that both electrical and mechanical vibration excitations act on the same spatial region covered by the observation network.

[0022] Through arrangement and interaction, a coordinated excitation sequence that is temporally interleaved and spatially unified is obtained.

[0023] Preferably, the step of inverting the potential difference data at each acquisition time to obtain a resistivity-time sequence includes:

[0024] Differential processing is performed on the potential difference data before, during and after the same excitation in the electrical excitation sequence to obtain the potential difference change data caused by the excitation.

[0025] Based on potential difference change data, power supply current parameters and electrode position geometric relationship, and according to the principle of electric field superposition, a set of linear constraint relationships on the conductivity distribution of underground medium is constructed.

[0026] The linear constraint relationship is solved iteratively, and each iteration adjusts the equivalent characterization of the dielectric polarization effect based on the current solution until the solution converges, thereby outputting a sequence of resistivity changes over time.

[0027] Preferably, the step of inverting the vibration waveform data at each acquisition time to obtain the sequence of elastic wave velocity changing with time includes:

[0028] From the vibration waveform data corresponding to the same excitation in the mechanical vibration excitation sequence, extract the arrival time and waveform amplitude attenuation information of the direct wave and specific reflected wave;

[0029] Based on the first arrival time and the spatial location of the observation point, and according to the integral relationship between travel time and slowness along the ray path, a set of travel time constraints on the slowness distribution of the medium is constructed.

[0030] Based on waveform amplitude attenuation information and the exponential decay relationship between wave energy and quality factor, a set of attenuation constraints on the distribution of medium quality factor is constructed.

[0031] The travel time constraint and the decay constraint are solved jointly. The distribution of slowness and quality factor is updated alternately and iteratively until synchronous convergence, thereby outputting the sequence of elastic wave velocity changing with time.

[0032] Preferably, the analysis of the phase relationship and change pattern of the resistivity sequence and elastic wave velocity sequence relative to their corresponding excitation time sequence, and the extraction of multi-physics dynamic response modes, includes:

[0033] By correlating and comparing the resistivity sequence and the elastic wave velocity sequence with their corresponding excitation timing, the characteristic moments when the trend of resistivity sequence change turns, and the characteristic moments when the trend of elastic wave velocity sequence change turns, are identified.

[0034] By comparing the chronological relationship between the characteristic moments of the resistivity sequence and the characteristic moments of the elastic wave velocity sequence, the temporal order relationship of the electro-elastic response is obtained.

[0035] By comparing the relative proportions of the resistivity sequence changes before and after the characteristic time with the elastic wave velocity sequence changes before and after the characteristic time, the intensity coupling relationship of the electro-elastic response is obtained.

[0036] Based on the combination of time sequence relationship and intensity coupling relationship, a multiphysics dynamic response mode is defined.

[0037] Preferably, the step of dividing regions with different response characteristics based on the spatial differences in the dynamic response modes of multiphysics includes:

[0038] The synchronicity of the electrical response characteristics and elastic response characteristics of the multiphysics dynamic response mode at various locations in space is compared to obtain the inter-field synchronicity relationship at each location.

[0039] The statistical distribution characteristics of inter-field synchronicity relationships are analyzed throughout the entire space. Based on these statistical distribution characteristics, the main distribution interval of inter-field synchronicity values ​​is defined, and positions where the synchronicity values ​​are outside the main distribution interval and are spatially continuous are selected to form preliminary anomaly zones.

[0040] The spatial gradient changes of the field synchronicity relationship within the initial anomaly zone and between it and the surrounding area are analyzed. The boundary lines of the region are determined based on the extreme values ​​of the gradient changes, thus completing the division of the region.

[0041] Preferably, the step of performing geological genetic correlation and hazard type determination based on the feature combination of response patterns in each region to generate a hazard zoning map includes:

[0042] For each defined region, the dominant multiphysics dynamic response mode within it is extracted, and the specific temporal sequence relationship and intensity coupling relationship contained in the mode are analyzed.

[0043] The specific time sequence relationship and intensity coupling relationship obtained by analysis are compared with the dynamic correlation between resistivity and elastic wave velocity expected by different dominant karst geological processes.

[0044] Based on the consistency comparison results, the dominant geological process most consistent with the response pattern of the area is determined and associated with the corresponding karst collapse hazard type;

[0045] The spatial extent of the region, its associated dominant geological processes, and the identified types of hazards are comprehensively characterized to generate a hazard zoning map.

[0046] To address the aforementioned problems, this invention also provides a multi-method collaborative detection system for karst collapse hazard areas, the system comprising:

[0047] The network deployment and excitation setting module is used to deploy an observation network in the target hazard area and implement electrical excitation sequences and mechanical vibration excitation sequences respectively, wherein each excitation sequence contains at least two excitations with different parameter settings;

[0048] The data acquisition module is used to collect potential difference data at a series of moments before, during, and after each excitation in the electrical excitation sequence; and to collect vibration waveform data at a series of moments before, during, and after each excitation in the mechanical vibration excitation sequence.

[0049] The time-series parameter inversion module is used to invert the potential difference data at each acquisition time to obtain the resistivity change sequence over time; and to invert the vibration waveform data at each acquisition time to obtain the elastic wave velocity change sequence over time.

[0050] The sequence analysis module is used to analyze the phase relationship and change pattern of resistivity sequence and elastic wave velocity sequence relative to their corresponding excitation time sequence, and extract the multi-physics dynamic response mode.

[0051] The hazard zoning module is used to divide regions with different response characteristics based on the spatial differences in the dynamic response modes of multiphysics fields; and to perform geological genetic correlation and hazard type determination based on the characteristic combination of response modes in each region, generating a hazard zoning map.

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

[0053] 1. By simultaneously implementing electrical excitation sequences and mechanical vibration excitation sequences, and collecting dynamic response data throughout the entire process before, during, and after excitation application, the co-acquisition of two types of physical parameters, resistivity and elastic wave velocity, under a unified time reference was achieved. By analyzing the phase relationship and change pattern of the two parameters relative to the excitation sequence, multi-physics field dynamic response modes characterizing dynamic processes such as groundwater seepage and particle migration can be extracted. This enables direct characterization of the activity and evolution sequence of dynamic processes within the hazard area, enhancing the ability to identify the development stage of the hazard from the physical response.

[0054] 2. By constructing constraint relationships, the response channels of different physical mechanisms are isolated in principle, avoiding mechanism aliasing. In the iterative solution, the equivalent representation of secondary effects such as polarization and attenuation is dynamically adjusted according to the current solution, thereby realizing adaptive learning and accurate separation of mixed responses. This significantly improves the signal-to-noise ratio and fidelity of extracting weak real dynamic signals from strong background noise, and ensures the physical consistency and reliability of the basic data on which subsequent pattern extraction depends.

[0055] 3. By establishing inter-field synchronicity relationships to divide regions with different response characteristics, and comparing the dynamic patterns extracted within the regions with the expected response relationships of specific karst geological processes, it is possible to achieve direct correlation and judgment from the dynamic characteristics of physical fields to specific geological genesis mechanisms and hazard types. Attached Figure Description

[0056] Figure 1 A flowchart of a multi-method collaborative detection method for karst collapse hazard areas provided by the present invention;

[0057] Figure 2 The present invention provides a modular structure diagram of a multi-method collaborative detection system for karst collapse hazard areas. Detailed Implementation

[0058] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative effort are within the scope of protection of the present invention.

[0059] Example 1, referring to Figure 1 The diagram shown is a flowchart illustrating a multi-method collaborative detection method for karst collapse hazard areas according to an embodiment of the present invention. In this embodiment, the multi-method collaborative detection method for karst collapse hazard areas includes:

[0060] S1. Set up an observation network in the target hazard area and implement electrical excitation sequence and mechanical vibration excitation sequence respectively, wherein each excitation sequence contains at least two excitations with different parameter settings;

[0061] S2. For each excitation in the electrical excitation sequence, collect potential difference data at a series of moments before, during, and after its application; for each excitation in the mechanical vibration excitation sequence, collect vibration waveform data at a series of moments before, during, and after its application.

[0062] S3. Invert the potential difference data at each acquisition time to obtain the resistivity change sequence over time; invert the vibration waveform data at each acquisition time to obtain the elastic wave velocity change sequence over time.

[0063] S4. Analyze the phase relationship and change pattern of the resistivity sequence and elastic wave velocity sequence relative to their corresponding excitation time sequence, and extract the multiphysics dynamic response mode;

[0064] S5. Based on the spatial differences in the dynamic response modes of multiphysics fields, regions with different response characteristics are divided; and based on the characteristic combinations of response modes in each region, geological genetic correlation and hazard type determination are carried out to generate a hazard zoning map.

[0065] In this embodiment of the invention, implementing the electrical excitation sequence includes:

[0066] Based on the geological information of the target hazard area, two current intensity parameters were set for the power supply electrodes, one for shallow large pore response and the other for deep small fracture response.

[0067] Based on geological information, two pulse duty cycle parameters were set for the power supply electrodes to obtain steady-state background field and rapid dynamic transition information, respectively.

[0068] According to a preset timing sequence, at least two different excitations consisting of two current intensity parameters and two pulse duty cycle parameters are alternately applied to the power supply electrodes to obtain an electrical excitation sequence.

[0069] It should be noted that the selection of the above current intensity parameters is based on the difference in depth sensitivity of pore / fracture structures of different scales to the current field response:

[0070] Lower current intensities, such as 2A, generate a current field that is mainly concentrated in the shallow layer and is more sensitive to electrical changes in the shallow macropore development zone.

[0071] Higher current intensity, such as 5A, can penetrate to deeper layers, enhancing the ability to detect areas with deep micro-fractures.

[0072] In this embodiment of the invention, implementing a mechanical vibration excitation sequence includes:

[0073] Based on the geological information of the target hazard area, two frequency parameters were set for the excitation source to enhance the sensitivity to low-density filling media and high-density intact bedrock, respectively.

[0074] Based on geological information, two amplitude parameters were set for the excitation source to detect the structures in the near-field strong attenuation region and the far-field weak attenuation region, respectively.

[0075] According to a preset timing sequence, at least two different excitations consisting of two frequency parameters and two amplitude parameters are alternately applied to the excitation source to obtain a mechanical vibration excitation sequence.

[0076] Specifically, the implementation scenario is a karst development area. The shallow layer of this area is a Quaternary alluvial silty clay interbedded with sand, containing a large number of large pores, with a burial depth of 0 to 15 meters. The deep layer is a carbonate rock stratum, containing fine fissures, with a burial depth of less than 15 meters.

[0077] First, the typical karst development area of ​​the target area is macroscopically delineated using remote sensing image recognition technology. Then, preprocessing operations such as bold radiometric correction, geometric correction, and cloud and shadow removal are performed on the acquired remote sensing images of the target area.

[0078] Then, using object-oriented remote sensing interpretation methods, information on karst landform features such as karst depressions, dissolution fissures, and sinkholes is extracted. Combined with topographic elevation data, macroscopic anomaly zones of potential karst collapse hazards are delineated. Based on this, the layout range and density of the observation network are optimized.

[0079] Specifically, a high-density electrical resistivity tomography (ERT) network was set up, with 50 copper-clad steel power supply electrodes arranged in a cross pattern at 5-meter intervals, vibration sensors arranged at 10-meter intervals, and hydraulic exciters placed at the center of the observation network and at the surrounding karst fissure development points identified by remote sensing interpretation. Then, electrical excitation sequences and mechanical vibration excitation sequences with two different parameter settings were implemented.

[0080] Then, the power supply signal during the application of the electrical excitation sequence is subjected to constant current source amplitude stabilization processing and 50Hz power frequency notch filtering to avoid potential difference data distortion caused by power supply fluctuations and power frequency interference; the excitation signal during the application of the mechanical vibration excitation sequence is subjected to time-frequency modulation and wave field shaping processing to improve the directionality and penetration of elastic waves in the underground medium.

[0081] It should be noted that the core physical mechanism of the electrical excitation sequence implemented in this scheme is the conduction and polarization effect of the current field in the underground medium. The resulting Joule heating and electrostriction effects are extremely weak and insufficient to excite elastic wave signals that can be effectively detected by vibration sensors. Therefore, the implementation of electrical excitation will not cause perceptible mechanical interference to the subsequent elastic wave detection of mechanical vibration excitation.

[0082] Furthermore, the power supply electrode parameters are set using the electrode parameter adjustment panel of the high-density electrical resistivity meter. Based on the geological information that the shallow silty clay interbedded with sand layers in the karst development area have large pores and strong conductivity, while the deep carbonate rock strata have small fractures and weak conductivity, the power supply electrode current intensity can be adjusted to 2A for the response of shallow large pores and to 5A for the response of deep small fractures.

[0083] Meanwhile, to acquire the steady-state background field, the pulse duty cycle is adjusted to 80%, and to acquire the rapid dynamic transition information, it is adjusted to 30%. Then, according to the timing sequence of 2A current plus 80% duty cycle to 5A current plus 30% duty cycle, the power supply electrodes are controlled by the electrical resistivity meter controller to alternately apply the two excitations, thus obtaining an electrical excitation sequence suitable for the investigation of karst development areas.

[0084] Furthermore, the frequency and amplitude adjustment panel of the hydraulic vibrator is used to set parameters. Based on the geological information that the shallow silty clay interbedded with sand in the karst development area is a low-density filling medium that is easy to vibrate, and the deep carbonate rock strata are high-density intact bedrock with fast vibration attenuation, the excitation source frequency is adjusted to 10Hz to enhance the sensitivity to the low-density filling medium and to 50Hz to enhance the sensitivity to the high-density intact bedrock.

[0085] Meanwhile, the amplitude of the excitation source is adjusted to 80% for the near-field strong attenuation zone and to 40% for the far-field weak attenuation zone. Then, according to the timing sequence of 10Hz frequency plus 80% amplitude to 50Hz frequency plus 40% amplitude, the excitation source is controlled by the exciter controller to alternately apply the two excitations to obtain a mechanical vibration excitation sequence suitable for the investigation of karst development areas.

[0086] In this embodiment of the invention, step S1 further includes:

[0087] The excitations in the electrical excitation sequence and the excitations in the mechanical vibration excitation sequence are arranged in time according to a preset alternation order;

[0088] This ensures that both electrical and mechanical vibration excitations act on the same spatial region covered by the observation network.

[0089] Through arrangement and interaction, a coordinated excitation sequence that is temporally interleaved and spatially unified is obtained.

[0090] Specifically, relying on the synchronous control terminal of the high-density electrical resistivity meter and the hydraulic vibrator, the execution duration of each excitation sequence of the electrical excitation sequence and the mechanical vibration excitation sequence adapted to the karst development area survey is extracted. The two types of excitations are arranged sequentially on the time axis according to the order of the single execution of electrical excitation followed by the single execution of mechanical vibration excitation, so that the execution cycle of the next type of excitation starts immediately after the execution cycle of the previous type of excitation ends, thus completing the time arrangement of all excitations without time overlap.

[0091] Specifically, a high-precision clock synchronization module can be introduced here, with a timing accuracy of ≤1ms, to unify the time base of electrical excitation and mechanical vibration excitation, and avoid the deviation of multi-physics data correlation caused by timing misalignment.

[0092] Furthermore, the electrical excitation applied by the copper-clad steel power supply electrode and the mechanical vibration excitation applied by the hydraulic vibrator are both applied to the distribution area of ​​silty clay interbedded sand layer and carbonate rock strata covered by the cross-shaped observation network in the karst development area survey. This ensures that the electrical effect range of the electrical excitation and the elastic wave propagation range of the mechanical vibration excitation completely cover the physical boundary of the observation network, ensuring that the action space of the two types of excitation is not offset and completely overlaps.

[0093] Finally, through the above-mentioned time arrangement and spatial action setting, the electrical excitation sequence and the mechanical vibration excitation sequence are made to form an excitation form that is executed alternately in time and completely overlapped in space within the observation network coverage area of ​​the karst development area. In the end, a time-domain staggered and spatially unified synergistic excitation sequence suitable for the detection of karst collapse hazards in karst development areas is obtained.

[0094] In summary, this scheme deploys an observation network in the target hazard area and simultaneously implements electrical excitation sequences and mechanical vibration excitation sequences with at least two different parameter settings. These sequences are then time-interleaved and spatially unified to form a coordinated excitation. The different current intensities and pulse duty cycles of the electrical excitation can be adapted to the detection of shallow large pores and deep small cracks, respectively, taking into account both steady-state background field and dynamic transition information. The different frequencies and amplitude parameters of the mechanical vibration excitation can enhance the sensitivity to media of different densities, covering the structural detection needs of near-field and far-field.

[0095] Furthermore, this solution integrates remote sensing image recognition technology to delineate macroscopic hidden danger areas, solving the problem of blind deployment of traditional in-situ detection and observation networks. At the same time, it improves the quality of excitation signals through multi-dimensional signal processing, effectively reducing environmental and equipment interference.

[0096] Overall, the two work together in the same spatial region to achieve multi-parameter and all-round detection coverage. This not only ensures the breadth and depth of detection, but also improves the targeting of different geological features through parameter differentiation design. This lays a solid foundation for subsequent acquisition of dynamic sequences of resistivity and elastic wave velocity and extraction of multi-physics response modes, effectively solving the problem of the limitations of detection by a single excitation method.

[0097] Specifically, when collecting potential difference data and vibration waveform data, a high-density electrical resistivity meter can be connected to the copper-clad steel measuring electrode of the observation network, and a vibration waveform analyzer can be connected to the vibration sensor. The two devices can be synchronously calibrated with the execution timing of the electrical excitation and mechanical vibration excitation to ensure that the acquisition time is accurately matched with the excitation action.

[0098] Specifically, the acquired raw potential difference data can be preprocessed with baseline correction, bad pixel removal and interpolation completion, and the raw vibration waveform data can be preprocessed with detrending and 5-100Hz bandpass filtering to remove environmental interference and equipment errors during the data acquisition process, ensuring the integrity and validity of the raw data.

[0099] It should be noted that the vibration waveform data here refers to the original time-series record reflecting the propagation characteristics of underground elastic waves.

[0100] Furthermore, for each excitation of the electrical excitation sequence, the potential difference data between electrodes in the silty clay interbedded with sand layer and the carbonate rock strata distribution area in the karst development area are collected moment by moment using a high-density electrical resistivity tomography instrument at a series of moments before, during and after the excitation.

[0101] For each excitation in the mechanical vibration excitation sequence, the vibration waveform data of the elastic wave generated by the excitation propagating in the underground medium of the karst development area is collected by a vibration waveform analyzer at the same time. The excitation batch and time node corresponding to each data collection are recorded synchronously throughout the process.

[0102] In this embodiment of the invention, the potential difference data at each acquisition time are inverted to obtain a sequence of resistivity changes over time, including:

[0103] Differential processing is performed on the potential difference data before, during and after the same excitation in the electrical excitation sequence to obtain the potential difference change data caused by the excitation.

[0104] Based on potential difference change data, power supply current parameters and electrode position geometric relationship, and according to the principle of electric field superposition, a set of linear constraint relationships on the conductivity distribution of underground medium is constructed.

[0105] The linear constraint relationship is solved iteratively, and each iteration adjusts the equivalent characterization of the dielectric polarization effect based on the current solution until the solution converges, thereby outputting a sequence of resistivity changes over time.

[0106] The formula for the linear constraint relationship is as follows:

[0107]

[0108] In the formula, Indicates the source of motivation At the observation point The potential difference change caused by the location Represents the underground three-dimensional spatial domain. This represents the source current density distribution injected by the power supply electrode. Indicates spatial location At the excitation frequency The equivalent complex conductivity is below. Indicates the first The polarization parameter vector obtained after the iteration. Represents the Green's function, describing the point The unit current source at the observation point The generated potential, This represents the gradient operator.

[0109] Specifically, taking karst development areas as the implementation scenario, a high-density electrical resistivity meter and a vibration waveform analyzer can be used to retrieve potential difference data and vibration waveform data at various collection times in the karst collapse hazard area of ​​the karst development area. The resistivity value of the underground medium is calculated by inverting the potential difference data moment by moment according to the electrical excitation time sequence, and the elastic wave velocity value of the underground medium is calculated by inverting the vibration waveform data moment by moment according to the mechanical vibration excitation time sequence.

[0110] It should be noted that before the inversion calculation, the preprocessed potential difference data is regularized to reduce the ambiguity in the inversion process and improve the inversion accuracy of the resistivity sequence.

[0111] Then, all resistivity calculation results and elastic wave velocity calculation results are arranged in chronological order according to the acquisition time, respectively, to obtain the resistivity variation over time sequence and the elastic wave velocity variation over time sequence for detection in karst development areas.

[0112] Furthermore, using the data analysis terminal of a high-density electrical resistivity meter, potential difference data for the same excitation sequence in the karst development area were retrieved before, during, and after the excitation. The baseline potential difference data before the excitation was subtracted from the potential difference data at each time point during and after the excitation, and the data difference was calculated moment by moment to obtain the potential difference change data caused by the excitation in the distribution area of ​​silty clay interbedded with sand and carbonate rock strata in the karst development area.

[0113] Furthermore, based on the potential difference variation data of the karst development area, combined with the power supply current parameters set by the high-density electrical resistivity meter and the actual spatial geometric relationship between the copper-clad steel power supply electrode and the measuring electrode, the underground three-dimensional spatial domain is divided into calculation units that match the electrode arrangement according to the principle of electric field superposition.

[0114] For example, the midpoint of the line connecting two adjacent electrodes is used as the boundary to divide several rectangular calculation units. Each unit corresponds to a silty clay sand layer or carbonate rock stratum at a specific depth and range in the karst development area. An initial conductivity value is assigned to each calculation unit. The relationship between the potential difference change data and the current density distribution is that the potential difference change increases linearly with the increase of current density and decreases with the increase of conductivity of the calculation unit.

[0115] When establishing the correspondence, the actual current density of each calculation unit can be substituted into the relationship. Combined with the measured potential difference change data of the karst development area, the numerical ratio of the unit conductivity to the observed potential difference can be determined. For example, if the numerical ratio is 1:1.2, the numerical correspondence between the conductivity of each calculation unit and the observed potential difference can be established, thus constructing a set of linear constraint relationships for the conductivity distribution of the underground medium in the karst development area.

[0116] Finally, using the data analysis end of a high-density electrical resistivity tomography (EDT) instrument, the established linear constraint relationship was iteratively solved. First, based on the dielectric characteristics of the silty clay interbedded with sand and carbonate rock strata in the karst development area, and combined with the historical iterative calibration data of karst exploration in the same area, the convergence threshold of the conductivity solution was set to a change of less than 0.01 S / m. This value is the practical convergence standard for the inversion of the electrical properties of the underground medium in the karst development area, indicating that the change in the conductivity solution is within a negligible error range. At this point, the solution can be determined to be completely converged. Each iteration adjusts the equivalent characterization parameters of the polarization effect of the underground medium in the karst development area based on the calculation results of the current conductivity solution.

[0117] For example, if the current solution shows that the conductivity of a certain calculation unit is low, corresponding to the well-developed pores in the silty clay interbedded with sand, the polarization equivalent parameter of that unit is increased to match the equivalent characterization parameter with the current electrical response state of the medium. The calculation is continued iteratively until the change amplitude of the conductivity solution reaches 0.01, and the solution is completely converged. Then, the resistivity value is calculated time by time through the conversion relationship between conductivity and resistivity, and the resistivity is output in time sequence after being arranged in time order to adapt to the resistivity change over time in the karst development area.

[0118] In summary, this inversion process uses the finite element method for numerical simulation, which more accurately fits the electrical distribution of complex underground geological structures. Compared with the traditional finite difference method, the spatial resolution of the inversion results is significantly improved.

[0119] In addition, to further ensure data reliability, interference can be checked on the collected vibration waveform data before inversion: by using time-frequency analysis technology, specific frequency band signals that are strongly correlated with the electrical excitation timing can be identified and filtered out; if potential electromagnetic-mechanical coupling noise is detected, noise term constraints can be introduced in the joint inversion, and possible residual interference can be used as one of the parameters to be solved for iterative correction, thereby eliminating the influence of any potential interference on the inversion results.

[0120] In this embodiment of the invention, the vibration waveform data at each acquisition time are inverted to obtain a sequence of elastic wave velocity changing with time, including:

[0121] From the vibration waveform data corresponding to the same excitation in the mechanical vibration excitation sequence, extract the arrival time and waveform amplitude attenuation information of the direct wave and specific reflected wave;

[0122] Based on the first arrival time and the spatial location of the observation point, and according to the integral relationship between travel time and slowness along the ray path, a set of travel time constraints on the slowness distribution of the medium is constructed.

[0123] Based on waveform amplitude attenuation information and the exponential decay relationship between wave energy and quality factor, a set of attenuation constraints on the distribution of medium quality factor is constructed.

[0124] The travel time constraint and the decay constraint are solved jointly. The distribution of slowness and quality factor is updated alternately and iteratively until synchronous convergence, thereby outputting the sequence of elastic wave velocity changing with time.

[0125] The time-travel constraints are expressed by the following formula:

[0126]

[0127] In the formula, Indicates from the source of motivation To the observation point When the elastic wave first arrives, Indicates that the elastic wave originates from point Time The propagation path of the rays, Indicates spatial location The slowness of the medium at that point is the elastic wave velocity. The reciprocal, Indicates along the ray path Line integral infinitesimal element.

[0128] The attenuation constraint relationship formula is as follows:

[0129]

[0130] In the formula, and They represent frequencies respectively. The wave at the observation point The actual amplitude at the excitation source point and the amplitude at the excitation source point The initial amplitude at that point, Indicates the frequency of the elastic wave. Indicates spatial location Elastic wave velocity at that location, Indicates spatial location Medium quality factor at the point, integral term Describes the wave energy along the path The cumulative decay effect.

[0131] Specifically, taking karst development areas as the implementation scenario, a vibration waveform analyzer can be used to identify the specific reflected waves generated by the direct waves reaching the carbonate rock strata interface from the vibration waveform data corresponding to the same excitation in the mechanical vibration excitation sequence of the karst development area survey. The arrival time of each waveform at the receiving end of the vibration sensor is extracted, and the amplitude attenuation value of the waveform from the excitation source to the observation point is calculated to obtain complete arrival time and waveform amplitude attenuation information.

[0132] Specifically, short-time averaging or long-time averaging algorithms can be used to extract the first arrival wave. Compared with manual picking, this can improve the accuracy of first arrival time identification, thereby effectively identifying the weak signal first arrival wave in the far-field weak attenuation region.

[0133] Furthermore, based on the extracted arrival time and the actual spatial position of the vibration sensor, and according to the integral relationship between travel time and slowness along the ray path, the elastic wave propagation ray path from the hydraulic exciter to each vibration sensor is determined.

[0134] For example, starting from the hydraulic vibrator, the ray path extends to each vibration sensor along a direction perpendicular to the carbonate rock stratum interface, dividing the ray path into continuous differential segments. Each differential segment corresponds to a segment of underground medium. A numerical correspondence between the slowness of each differential segment and the total travel time is established, that is, the total travel time is equal to the sum of the products of the slowness of each differential segment and the segment length, thus constructing a set of travel time constraints on the distribution of slowness of underground medium in karst development areas.

[0135] Furthermore, based on the extracted waveform amplitude attenuation information, and according to the exponential decay relationship between wave energy and quality factor, combined with the actual propagation frequency of the elastic wave, a numerical correspondence between the quality factor and amplitude attenuation of each propagation path differential segment is established.

[0136] For example, the greater the amplitude attenuation, the smaller the corresponding differential segment medium quality factor, which indicates that the medium is more loose. Combining the geological characteristics of loose silty clay interbedded with sand layers and intact carbonate rock strata in karst development areas, specific numerical ratios are determined to construct a set of attenuation constraint relationships for the distribution of underground medium quality factors in karst development areas.

[0137] Furthermore, a vibration waveform analyzer can be used to incorporate the established travel time constraint relationship and attenuation constraint relationship into the same calculation system for joint solution. First, the slowness distribution of the medium is updated according to the travel time constraint relationship, and then the quality factor distribution is updated according to the updated slowness distribution and attenuation constraint relationship. The distribution states of the two types of parameters are updated alternately and iteratively, and the calculation is continued until the solutions of slowness and quality factor reach a stable state and converge synchronously.

[0138] It should be noted that this joint solution process introduces a multi-objective optimization method to achieve simultaneous optimal inversion of slowness and quality factor, avoiding parameter distortion caused by single-constraint inversion and reducing the inversion error of elastic wave velocity sequence.

[0139] Finally, the elastic wave velocity values ​​are calculated hourly using the reciprocal relationship between slowness and elastic wave velocity, and then output as a sequence of elastic wave velocity changes over time that is suitable for detection in karst development areas.

[0140] In summary, this scheme isolates response channels of different physical mechanisms in principle, avoids mechanism aliasing, and dynamically adjusts secondary effects such as polarization and attenuation during iteration. It achieves adaptive learning and precise separation of mixed responses, significantly improves the signal-to-noise ratio and fidelity of extracting weak real dynamic signals from strong background noise, and ensures the physical consistency and reliability of the basic data.

[0141] Furthermore, this solution addresses the technical challenges of traditional detection methods, such as signal interference, multiple solutions in inversion, and parameter distortion, through a full-process approach including excitation signal processing, raw data preprocessing, and inversion algorithm optimization, thereby significantly improving data quality and inversion accuracy.

[0142] The simultaneously acquired sequences of resistivity and elastic wave velocity over time provide accurate and comprehensive data support for subsequent analysis of their correlation with excitation timing and extraction of multi-physics dynamic response modes, helping to more accurately delineate response characteristic regions and determine the types of potential hazards.

[0143] In this embodiment of the invention, the phase relationship and change pattern of the resistivity sequence and the elastic wave velocity sequence relative to their corresponding excitation time sequence are analyzed to extract the multiphysics dynamic response mode, including:

[0144] By correlating and comparing the resistivity sequence and the elastic wave velocity sequence with their corresponding excitation timing, the characteristic moments when the trend of resistivity sequence change turns, and the characteristic moments when the trend of elastic wave velocity sequence change turns, are identified.

[0145] By comparing the chronological relationship between the characteristic moments of the resistivity sequence and the characteristic moments of the elastic wave velocity sequence, the temporal order relationship of the electro-elastic response is obtained.

[0146] By comparing the relative proportions of the resistivity sequence changes before and after the characteristic time with the elastic wave velocity sequence changes before and after the characteristic time, the intensity coupling relationship of the electro-elastic response is obtained.

[0147] Based on the combination of time sequence relationship and intensity coupling relationship, a multiphysics dynamic response mode is defined.

[0148] It should be noted that the characteristic moment in this scheme refers to the time node when the trend of the resistivity sequence or elastic wave velocity sequence changes significantly on the time axis. This node corresponds to a key event in the karst dynamics process, such as the arrival of the seepage front or the occurrence of a sudden change in pore pressure.

[0149] Specifically, the following method can be used to identify the characteristic moments: perform first-order difference processing on the sequence to obtain the instantaneous rate of change sequence, perform sliding window fitting on the difference sequence, identify the extreme points where the slope sign changes within the window, and take the time corresponding to the identified extreme points as the characteristic moments of the sequence.

[0150] Compared to simple threshold comparison, this method can more effectively eliminate random noise interference and accurately capture the true dynamic response inflection points.

[0151] Specifically, the resistivity-time sequence, elastic wave velocity-time sequence, and corresponding electrical excitation timing sequence and mechanical vibration excitation timing sequence of the resistivity-time sequence and elastic wave velocity sequence of the resistivity-time sequence and elastic wave velocity sequence of the resistivity-time sequence are retrieved to be suitable for the resistivity-time sequence of the resistivity-time sequence and the elastic wave velocity sequence of the elastic wave velocity sequence. Based on the measured data of the medium response of the silty clay interbedded sand layer and carbonate rock strata in the resistivity-time sequence, the reasonable range of the electrical-elastic response time difference is set to 0.5s-2s and the reference ratio of the electrical-elastic response intensity coupling is 3:1. This range and ratio are the practical judgment criteria for the medium response of the resistivity-time sequence in the resistivity-time sequence.

[0152] Furthermore, the resistivity sequence is correlated with the electrical excitation time sequence moment by moment, and the elastic wave velocity sequence is correlated with the mechanical vibration excitation time sequence moment by moment. The correlation between the sequence values ​​and the excitation application state is compared point by point to identify the characteristic moment when the resistivity sequence changes from stable to rising, and the characteristic moment when the elastic wave velocity sequence changes from stable to falling, thus clarifying the specific time nodes when the change trends of the two types of sequences change.

[0153] It should be noted that the first-order difference method combined with sliding window fitting is used here to identify the turning point of the sequence change trend, effectively eliminating small fluctuations in the sequence caused by random noise and accurately capturing the real dynamic response turning point.

[0154] Furthermore, the identified resistivity sequence characteristic moments and elastic wave velocity sequence characteristic moments are compared with the alternating timing of electrical excitation and mechanical vibration excitation in the karst development area, and the time difference between the two types of characteristic moments is recorded.

[0155] If the characteristic moment of the resistivity sequence is earlier than the characteristic moment of the elastic wave velocity sequence, and the time difference is within a reasonable range of the medium response in the karst development area, then the temporal order of the electro-elastic response is determined to be that the resistivity response precedes the elastic wave velocity response.

[0156] Furthermore, the numerical difference of the resistivity sequence before and after its characteristic moment is calculated to obtain the resistivity change amplitude. At the same time, the numerical difference of the elastic wave velocity sequence before and after its characteristic moment is calculated to obtain the elastic wave velocity change amplitude. The magnitudes of the two change amplitudes are then compared.

[0157] For example, if the resistivity change is 0.3 Ω·m and the elastic wave velocity change is 0.1 km / s, and the resistivity change is greater than the elastic wave velocity change, the relative ratio between the two is determined to be 3:1, thus obtaining the intensity coupling relationship of the electro-elastic response of the karst medium in the karst development area.

[0158] Finally, based on the obtained time sequence relationship and intensity coupling relationship of the electro-elastic response, and combined with the medium characteristics of silty clay interbedded with sand and carbonate rock strata in karst development areas, the combination of resistivity response preceding elastic wave velocity response and resistivity change amplitude being 3:1 is defined as a multi-physics dynamic response mode adapted to the detection of karst development areas. This mode clearly includes the electro-elastic response characteristics of the underground medium in karst development areas.

[0159] In summary, this scheme first correlates and compares the resistivity sequence and elastic wave velocity sequence with the corresponding excitation timing sequence to identify the characteristic moments when the changing trends of the two types of sequences change. Then, by comparing the order of the characteristic moments, the electro-elastic response timing sequence relationship is obtained. By comparing the relative proportions of the change amplitudes before and after the characteristic moments, the intensity coupling relationship is obtained. Finally, based on the combination of these two relationships, a multiphysics dynamic response mode is defined.

[0160] Overall, this method breaks through the limitations of traditional single-physical-field detection, enabling in-depth mining of the dynamic correlation characteristics of two core physical parameters. It accurately captures key information characterizing dynamic processes related to karst collapse, such as groundwater seepage and particle migration, and directly depicts the activity and evolution sequence of dynamic processes within the hazard area. This effectively compensates for the shortcomings of existing technologies in identifying the development stage and evolution trend of hazards, providing a core basis for subsequent division of response areas based on spatial feature differences, correlation of geological causes, and determination of hazard types. It significantly improves the pertinence and effectiveness of karst collapse hazard detection.

[0161] In this embodiment of the invention, regions with different response characteristics are divided based on the spatial differences in the dynamic response modes of multiphysics, including:

[0162] The synchronicity of the electrical response characteristics and elastic response characteristics of the multiphysics dynamic response mode at various locations in space is compared to obtain the inter-field synchronicity relationship at each location.

[0163] The statistical distribution characteristics of inter-field synchronicity relationships are analyzed throughout the entire space. Based on these statistical distribution characteristics, the main distribution interval of inter-field synchronicity values ​​is defined, and positions where the synchronicity values ​​are outside the main distribution interval and are spatially continuous are selected to form preliminary anomaly zones.

[0164] The spatial gradient changes of the field synchronicity relationship within the initial anomaly zone and between it and the surrounding area are analyzed. The boundary lines of the region are determined based on the extreme values ​​of the gradient changes, thus completing the division of the region.

[0165] Specifically, data visualization and analysis equipment can be used to retrieve the multi-physics dynamic response mode adapted to the detection of karst development areas and the response data of each spatial location of the observation network, and compare the changes in the trigger time and amplitude of electrical response characteristics and elastic response characteristics at each location.

[0166] For example, by comparing the response start time of resistivity and elastic wave velocity at a certain location in the distribution area of ​​silty clay interbedded with sand in karst development area, and then comparing the amplitude change trend of the two at a certain location in the distribution area of ​​carbonate rock strata, it is possible to directly determine whether the electrical properties and elastic response at each location change synchronously, and obtain the inter-field synchronization relationship corresponding to each spatial location in karst development area.

[0167] It should be noted that spatial interpolation and rasterization are performed on the inter-field synchronicity relationship here, which can generate a spatial distribution map of the inter-field synchronicity value with a resolution of 1m×1m, intuitively presenting the spatial differences in response characteristics.

[0168] Furthermore, using statistical analysis equipment, the inter-field synchronicity relationships within the entire observation network of the karst development area are transformed into quantifiable synchronicity values. The frequency and distribution range of all synchronicity values ​​are statistically analyzed, and the numerical range with the largest proportion is determined according to the frequency from high to low. For example, the largest numerical range is 0.6-0.8. This range is defined as the main distribution range of inter-field synchronicity values. Then, synchronicity values ​​are screened location by location, and locations where the values ​​are outside the main distribution range and are distributed in a continuous sheet in the spatial distribution of the karst development area are selected. These continuous locations are integrated and delineated as the preliminary anomaly zone of karst collapse hazard in the karst development area.

[0169] It should be noted that this screening process combines the macroscopic anomaly zone boundary of remote sensing image recognition with constraints to avoid misjudgment of anomaly zones caused by local anomalies in in-situ detection data, thereby improving the accuracy of anomaly zone delineation.

[0170] Furthermore, using spatial gradient analysis equipment, the difference in field synchronicity values ​​between adjacent locations within the initial anomaly zone of the karst development area is calculated, and then the difference in synchronicity values ​​at the boundary between the initial anomaly zone and the surrounding normal zone is calculated to identify the locations where abrupt changes in the difference occur.

[0171] For example, at the boundary between silty clay and sand layers and carbonate rock strata in karst development areas, the synchronicity value at this location suddenly drops from 0.7 in the normal range to 0.2, with a difference of 0.5. This value forms a significant extreme value compared to the average difference of 0.1 between adjacent locations. Connecting these boundary locations where the difference reaches an extreme value sequentially forms a closed boundary line. This boundary line serves as the dividing line between different response characteristic areas in the karst development area, thus completing the division of the karst development area into regions with different response characteristics.

[0172] It should be noted that the variational method is used here to find the extreme locations of gradient changes. Compared with the traditional threshold method, it can more accurately determine the boundary of the region and fit the actual distribution of the underground geological structure.

[0173] In this embodiment of the invention, based on the characteristic combinations of response patterns in each region, geological genetic correlation and hazard type determination are performed to generate a hazard zoning map, including:

[0174] For each defined region, the dominant multiphysics dynamic response mode within it is extracted, and the specific temporal sequence relationship and intensity coupling relationship contained in the mode are analyzed.

[0175] The specific time sequence relationship and intensity coupling relationship obtained by analysis are compared with the dynamic correlation between resistivity and elastic wave velocity expected by different dominant karst geological processes.

[0176] Based on the consistency comparison results, the dominant geological process most consistent with the response pattern of the area is determined and associated with the corresponding karst collapse hazard type;

[0177] The spatial extent of the region, its associated dominant geological processes, and the identified types of hazards are comprehensively characterized to generate a hazard zoning map.

[0178] Specifically, for each subdivided area of ​​the karst development region, a geological data processing workstation can be used to sort out the response data region by region, determine the multi-physics dynamic response mode with the highest proportion in each region as the dominant mode, decompose the corresponding electro-elastic response time sequence relationship and intensity coupling relationship from the dominant mode, and clarify the response time sequence and amplitude ratio characteristics of each region.

[0179] Furthermore, data on various dominant geological processes from karst geological surveys in karst development areas can be retrieved, such as groundwater seepage processes, soil and rock loosening processes, and fissure development and propagation processes. The dynamic correlation between resistivity and elastic wave velocity corresponding to each process can be analyzed, and the specific time sequence relationship and intensity coupling relationship obtained from the analysis of each region can be compared with the expected correlation of each dominant geological process to determine whether the response time characteristics and amplitude variation characteristics of the two are consistent.

[0180] Furthermore, based on the results of consistency comparison, if the response pattern of a certain area completely matches the expected correlation of the groundwater seepage process, then the dominant geological process of that area is determined to be groundwater seepage, and it is associated with the soil seepage type of karst collapse hazard. If the response pattern of a certain area completely matches the expected correlation of the soil and rock loosening process, then the dominant geological process of that area is determined to be soil and rock loosening, and it is associated with the bedrock collapse type of karst collapse hazard. This completes the determination of the geological genesis correlation and hazard type of each area in the karst development region.

[0181] It should be noted that the results of the hazard type determination are jointly verified with the karst landform features identified by remote sensing images. If the hazard area determined by in-situ detection coincides with the solution fissure area and sinkhole area interpreted by remote sensing, the accuracy of the hazard type is further confirmed, achieving dual determination by in-situ detection and remote sensing verification.

[0182] Finally, a geographic information plotter can be used to label the planar spatial range, boundary lines, corresponding dominant geological processes, and types of karst collapse hazards of each area in the karst development region in layers.

[0183] For example, in the preliminary anomaly areas where silty clay interbedded with sand layers are distributed in karst development areas, the groundwater seepage process and soil permeability-type hazards are marked. In the anomaly areas where carbonate rock strata are distributed, the fracture development and expansion process and bedrock collapse-type hazards are marked. Then, the location information of copper-clad steel power supply electrodes and vibration sensors of the observation network are superimposed to form a hazard zoning map that contains the unique geological features and hazard information of karst development areas.

[0184] It should be noted that this method integrates and overlays the hazard zoning map with the remote sensing image base map to achieve a combined macro and micro display of hazard information. It also supports layer editing, attribute querying, and spatial analysis of the map, providing visualized and analyzable technical support for the subsequent treatment and risk prevention of karst collapse hazards.

[0185] In summary, this scheme first compares the synchronicity of electrical and elastic response characteristics of the multiphysics dynamic response modes at different locations in space to obtain the synchronicity relationship between fields. Then, it screens out preliminary anomaly areas by analyzing their statistical distribution characteristics and determines the regional boundaries based on the extreme values ​​of spatial gradient changes, thus completing the division of regions with different response characteristics. Subsequently, it extracts the dominant response modes of each region, analyzes the coupling relationship between time sequence and intensity, compares them with the expected correlation with different dominant karst geological processes, determines the geological causes and hazard types, and finally generates a hazard zoning map.

[0186] By integrating remote sensing image recognition technology, the entire detection process was optimized and verified. This supplemented multi-dimensional signal processing techniques and clarified their technical effects. It achieved an organic combination of traditional in-situ detection and cutting-edge remote sensing technology, thereby enhancing the technology's advancement and practicality.

[0187] Overall, it has achieved a precise correlation between the dynamic characteristics of the physical field and the geological formation mechanism and the type of hidden danger, breaking the limitation of traditional detection that can only identify structures. It accurately locks the spatial range and core causes of the hidden danger area, and the generated zoning map intuitively presents the distribution and type information of the hidden danger. It not only improves the pertinence and reliability of karst collapse hidden danger detection, but also provides a clear and scientific basis for subsequent hidden danger management and risk prevention and control. It effectively makes up for the shortcomings of existing technologies in judging the development stage and evolution trend of hidden dangers.

[0188] Example 2, as Figure 2 The diagram shown is a modular structure diagram of a multi-method collaborative detection system for karst collapse hazard areas provided by the present invention, which includes:

[0189] The network deployment and excitation setting module 101 is used to deploy an observation network in the target hidden danger area and implement electrical excitation sequence and mechanical vibration excitation sequence respectively, wherein each excitation sequence includes at least two excitations with different parameter settings.

[0190] The data acquisition module 102 is used to acquire potential difference data at a series of moments before, during and after each excitation in the electrical excitation sequence; and to acquire vibration waveform data at a series of moments before, during and after each excitation in the mechanical vibration excitation sequence.

[0191] The time-series parameter inversion module 103 is used to invert the potential difference data at each acquisition time to obtain the resistivity change sequence over time; and to invert the vibration waveform data at each acquisition time to obtain the elastic wave velocity change sequence over time.

[0192] The sequence analysis module 104 is used to analyze the phase relationship and change pattern of the resistivity sequence and the elastic wave velocity sequence relative to their corresponding excitation time sequence, and extract the multi-physics dynamic response mode.

[0193] The hazard zoning module 105 is used to divide regions with different response characteristics based on the spatial differences in the dynamic response modes of multiphysics fields; and to perform geological genetic correlation and hazard type determination based on the characteristic combination of response modes in each region, thereby generating a hazard zoning map.

[0194] In the several embodiments provided by this invention, it should be understood that the disclosed methods and systems can be implemented in other ways. For example, the system embodiments described above are merely illustrative; for instance, the division of modules is only a logical functional division, and other division methods may be used in actual implementation.

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

[0196] Furthermore, the functional modules in the various embodiments of the present invention can be integrated into one processing unit, or each unit can exist physically separately, or two or more units can be integrated into one unit. The integrated unit can be implemented in hardware or in the form of hardware plus software functional modules.

[0197] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the present invention can be implemented in other specific forms without departing from the spirit or essential characteristics of the present invention.

[0198] The embodiments of this application can acquire and process relevant data based on artificial intelligence technology. Artificial intelligence is the theory, method, technology, and application system that uses digital computers or machines controlled by digital computers to simulate, extend, and expand human intelligence, perceive the environment, acquire knowledge, and use that knowledge to obtain optimal results.

[0199] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit it. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical solutions of the present invention without departing from the spirit and scope of the technical solutions of the present invention.

Claims

1. A multi-method collaborative detection method for karst collapse hazard areas, characterized in that, The method includes: S1. Set up an observation network in the target hazard area and implement electrical excitation sequence and mechanical vibration excitation sequence respectively, wherein each excitation sequence contains at least two excitations with different parameter settings; S2. For each excitation in the electrical excitation sequence, collect potential difference data at a series of moments before, during, and after its application; for each excitation in the mechanical vibration excitation sequence, collect vibration waveform data at a series of moments before, during, and after its application. S3. Invert the potential difference data at each acquisition time to obtain the resistivity change sequence over time; invert the vibration waveform data at each acquisition time to obtain the elastic wave velocity change sequence over time. S4. Analyze the phase relationship and change pattern of the resistivity sequence and elastic wave velocity sequence relative to their corresponding excitation time sequence, and extract the multiphysics dynamic response mode; S5. Based on the spatial differences in the dynamic response modes of multiphysics fields, regions with different response characteristics are divided; and based on the characteristic combinations of response modes in each region, geological genesis correlation and hazard type determination are carried out to generate a hazard zoning map. The analysis of the phase relationship and variation of the resistivity sequence and elastic wave velocity sequence relative to their corresponding excitation time sequence extracts multi-physics dynamic response modes, including: By correlating and comparing the resistivity sequence and the elastic wave velocity sequence with their corresponding excitation timing, the characteristic moments when the trend of resistivity sequence change turns, and the characteristic moments when the trend of elastic wave velocity sequence change turns, are identified. By comparing the chronological relationship between the characteristic moments of the resistivity sequence and the characteristic moments of the elastic wave velocity sequence, the temporal order relationship of the electro-elastic response is obtained. By comparing the relative proportions of the resistivity sequence changes before and after the characteristic time with the elastic wave velocity sequence changes before and after the characteristic time, the intensity coupling relationship of the electro-elastic response is obtained. Based on the combination of time sequence relationship and intensity coupling relationship, a multiphysics dynamic response mode is defined; The process of generating a hazard zoning map by performing geological genetic correlation and hazard type determination based on the feature combination of response patterns in each region includes: For each defined region, the dominant multiphysics dynamic response mode within it is extracted, and the specific temporal sequence relationship and intensity coupling relationship contained in the mode are analyzed. The specific time sequence relationship and intensity coupling relationship obtained by analysis are compared with the dynamic correlation between resistivity and elastic wave velocity expected by different dominant karst geological processes. Based on the consistency comparison results, the dominant geological process most consistent with the response pattern of the area is determined and associated with the corresponding karst collapse hazard type; The spatial extent of the region, its associated dominant geological processes, and the identified types of hazards are comprehensively characterized to generate a hazard zoning map.

2. The multi-method collaborative detection method for karst collapse hazard areas as described in claim 1, characterized in that, The implementation of the electrical excitation sequence includes: Based on the geological information of the target hazard area, two current intensity parameters were set for the power supply electrodes, one for shallow large pore response and the other for deep small fracture response. Based on geological information, two pulse duty cycle parameters were set for the power supply electrodes to obtain steady-state background field and rapid dynamic transition information, respectively. According to a preset timing sequence, at least two different excitations consisting of two current intensity parameters and two pulse duty cycle parameters are alternately applied to the power supply electrodes to obtain an electrical excitation sequence.

3. The multi-method collaborative detection method for karst collapse hazard areas as described in claim 1, characterized in that, The implementation of the mechanical vibration excitation sequence includes: Based on the geological information of the target hazard area, two frequency parameters were set for the excitation source to enhance the sensitivity to low-density filling media and high-density intact bedrock, respectively. Based on geological information, two amplitude parameters were set for the excitation source to detect the structures in the near-field strong attenuation region and the far-field weak attenuation region, respectively. According to a preset timing sequence, at least two different excitations consisting of two frequency parameters and two amplitude parameters are alternately applied to the excitation source to obtain a mechanical vibration excitation sequence.

4. The multi-method collaborative detection method for karst collapse hazard areas as described in claim 1, characterized in that, S1 further includes: The excitations in the electrical excitation sequence and the excitations in the mechanical vibration excitation sequence are arranged in time according to a preset alternation order; This ensures that both electrical and mechanical vibration excitations act on the same spatial region covered by the observation network. Through arrangement and interaction, a coordinated excitation sequence that is temporally interleaved and spatially unified is obtained.

5. The multi-method collaborative detection method for karst collapse hazard areas as described in claim 1, characterized in that, The process of inverting the potential difference data at each acquisition time to obtain a resistivity-time sequence includes: Differential processing is performed on the potential difference data before, during and after the same excitation in the electrical excitation sequence to obtain the potential difference change data caused by the excitation. Based on potential difference change data, power supply current parameters and electrode position geometric relationship, and according to the principle of electric field superposition, a set of linear constraint relationships on the conductivity distribution of underground medium is constructed. The linear constraint relationship is solved iteratively, and each iteration adjusts the equivalent characterization of the dielectric polarization effect based on the current solution until the solution converges, thereby outputting a sequence of resistivity changes over time.

6. The multi-method collaborative detection method for karst collapse hazard areas as described in claim 1, characterized in that, The process of inverting the vibration waveform data at each acquisition moment to obtain the sequence of elastic wave velocity changing with time includes: From the vibration waveform data corresponding to the same excitation in the mechanical vibration excitation sequence, extract the arrival time and waveform amplitude attenuation information of the direct wave and specific reflected wave; Based on the first arrival time and the spatial location of the observation point, and according to the integral relationship between travel time and slowness along the ray path, a set of travel time constraints on the slowness distribution of the medium is constructed. Based on waveform amplitude attenuation information and the exponential decay relationship between wave energy and quality factor, a set of attenuation constraints on the distribution of medium quality factor is constructed. The travel time constraint and the decay constraint are solved jointly. The distribution of slowness and quality factor is updated alternately and iteratively until synchronous convergence, thereby outputting the sequence of elastic wave velocity changing with time.

7. The multi-method collaborative detection method for karst collapse hazard areas as described in claim 1, characterized in that, The process of dividing regions with different response characteristics based on the spatial differences in the dynamic response modes of multiphysics includes: The synchronicity of the electrical response characteristics and elastic response characteristics of the multiphysics dynamic response mode at various locations in space is compared to obtain the inter-field synchronicity relationship at each location. The statistical distribution characteristics of inter-field synchronicity relationships are analyzed throughout the entire space. Based on these statistical distribution characteristics, the main distribution interval of inter-field synchronicity values ​​is defined, and positions where the synchronicity values ​​are outside the main distribution interval and are spatially continuous are selected to form preliminary anomaly zones. The spatial gradient changes of the field synchronicity relationship within the initial anomaly zone and between it and the surrounding area are analyzed. The boundary lines of the region are determined based on the extreme values ​​of the gradient changes, thus completing the division of the region.

8. A multi-method collaborative detection system for karst collapse hazard areas, used to implement the multi-method collaborative detection method for karst collapse hazard areas as described in any one of claims 1-7, characterized in that, The system includes: The network deployment and excitation setting module is used to deploy an observation network in the target hazard area and implement electrical excitation sequences and mechanical vibration excitation sequences respectively, wherein each excitation sequence contains at least two excitations with different parameter settings; The data acquisition module is used to collect potential difference data at a series of moments before, during, and after each excitation in the electrical excitation sequence; and to collect vibration waveform data at a series of moments before, during, and after each excitation in the mechanical vibration excitation sequence. The time-series parameter inversion module is used to invert the potential difference data at each acquisition time to obtain the resistivity change sequence over time; and to invert the vibration waveform data at each acquisition time to obtain the elastic wave velocity change sequence over time. The sequence analysis module is used to analyze the phase relationship and change pattern of resistivity sequence and elastic wave velocity sequence relative to their corresponding excitation time sequence, and extract the multi-physics dynamic response mode. The hazard zoning module is used to divide regions with different response characteristics based on the spatial differences in the dynamic response modes of multiphysics fields; and to perform geological genetic correlation and hazard type determination based on the characteristic combination of response modes in each region, generating a hazard zoning map.

Images