Method and system for intelligently reducing rockburst through advanced targeting non-explosion fracturing

Abstract

The invention discloses a method and a system for intelligently reducing rockburst through advanced targeting non-explosion fracturing. The method comprises the following steps of: arranging a monitoring system in a rockburst tendency area within a preset range from a tunnel face to obtain activity data of a rock mass; the installed monitoring system is connected to the micro-seismic monitoring equipment through a data acquisition line so as to collect the rock mass activity data acquired by the monitoring system to obtain micro-seismic monitoring data; a micro-seismic analysis device and an early warning device are in communication connection with a micro-seismic monitoring device, a time sequence model is established through the micro-seismic analysis device, micro-seismic monitoring data are analyzed, and the alarm device predicts a potential rockburst area of a rockburst tendency area near a tunnel construction face based on a micro-seismic analysis result. The grades of the potential rockburst areas are divided; and non-burst drill holes are formed in the potential rockburst area in a targeted mode according to the early warning signals, and a low-wave impedance propping agent is injected into the non-burst drill holes so as to eliminate potential rockburst disasters.

Description

Technical Field

[0001] The present invention relates to the technical field of rock engineering such as mines and tunnels, and in particular to a method and system for intelligent advanced targeted non-explosive fracturing to mitigate rockbursts. Background Art

[0002] During deep mining, mining and deep tunnel excavation, the hard surrounding rock is often in a high stress and strong disturbance environment, and dynamic disasters such as rock burst, rock burst, and mine earthquake occur frequently. When a rock burst occurs, a large amount of elastic energy will be released instantly, throwing the rock out, endangering the safety of personnel and equipment, and may cause serious accidents.

[0003] To this end, technicians in this field have designed and developed a variety of rock burst control methods, such as pre-blasting or rock fracturing to relieve pressure, and drilling pre-stress release holes to relieve stress concentration. By reducing the risk of rock bursts, the occurrence of impact accidents can be reduced, which has a certain effect in preventing and controlling disasters.

[0004] However, currently, the advanced stress release or pressure relief has become a complicated and one-size-fits-all approach. A large number of solutions are not targeted. They often work well in tests, but construction units are not very enthusiastic about applying them. The main reason is that taking measures in general on the working surface will increase costs dramatically, and in fact, it is unnecessary in some areas. In addition, due to problems such as complex processes, risk avoidance measures and ventilation procedures, most methods can only become "paper solutions."

[0005] At present, the advanced control technology of rock burst disasters has fallen into the inertial thinking of complexity, resulting in more and more new methods of pressure relief (stress release) and worse applicability. Therefore, under the premise of reducing the hidden dangers of rock burst, the accuracy of the scheme application should be improved as much as possible, and at the same time, the manufacturing cost and construction difficulty of the equipment should be greatly reduced, so that the implementation of the scheme can return to being simple and convenient. Summary of the invention

[0006] In response to the problems and needs raised above, this solution proposes a method and system for intelligent advanced targeted non-explosive fracturing to mitigate rockbursts. Due to the adoption of the following technical features, the above technical objectives can be achieved and many other technical effects can be brought about.

[0007] One object of the present invention is to provide a method for intelligent advanced targeted non-explosive fracturing to mitigate rockburst, comprising the following steps:

[0008] S10: arranging a monitoring system in a rockburst prone area within a predetermined range from the tunnel face to obtain activity data of the rock mass, wherein the monitoring system includes: a displacement sensor, a stress sensor, a temperature sensor, an acoustic wave sensor, and a vibration sensor;

[0009] S20: connecting the installed monitoring system to the microseismic monitoring equipment through a data acquisition line to collect the rock mass activity data acquired by the monitoring system to obtain microseismic monitoring data;

[0010] S30: connecting the microseismic analysis equipment and the early warning equipment to the microseismic monitoring equipment, establishing a time series model through the microseismic analysis equipment, analyzing the microseismic monitoring data, and using the alarm equipment to predict the potential rock burst area in the rock burst prone area near the tunnel construction face based on the microseismic analysis results, and classifying the potential rock burst area;

[0011] S40: opening non-explosion boreholes in the potential rock burst area in a targeted manner according to the early warning signal, and injecting low wave impedance proppants into the non-explosion fracture boreholes to eliminate potential rock burst disasters.

[0012] In an example of the present invention, in step S30, establishing a time series model includes the following steps:

[0013] S31: Use the historical and current data of the maximum energy of microseismic events in the monitoring area in a single day as samples to perform a stationary analysis. If the data is not stationary, perform a differential operation to obtain an approximately stationary sequence.

[0014] S32: Select an appropriate model to fit the maximum energy of the sample microseismic event in a single day, determine the order and estimate the parameters of the model, and perform an adaptability test on the model;

[0015] S33: Use the established model to predict the future value of the difference, and then obtain the future value of the maximum energy of the microseismic day through the reverse difference operation. During the mining process, use the established model to perform sample extrapolation prediction, and by comparing with the warning indicators, danger warning can be achieved.

[0016] In one example of the present invention, in step S31, the stationarity analysis includes the following steps:

[0017] First fit the adaptive model of the sequence, and then test whether there is a unit root in the characteristic equation composed of the autoregressive part parameters of the model. If there is no unit root, the sequence is stable, and if there is a unit root, the sequence is not stable.

[0018] The specific process is as follows: for the AR(p) process

[0019] X t =φ1X t-1 +φ2X t-2 +···+φ p X t-p +ε t

[0020] make The following hypothesis test is performed:

[0021]

[0022] Construct the ADF test statistic: in is the sample standard deviation of parameter ρ. ARMA(p,q) is the same as AR(p). If it is MA(q), the regression part has no parameters and is definitely stable.

[0023] In one example of the present invention, in step S32, determining the order of the model includes the following steps:

[0024] First, by judging whether formula (1) and formula (2) are valid, we can roughly judge the truncation of the autocorrelation coefficient of the stationary series and the truncation of the partial correlation coefficient of the stationary series and make a preliminary order determination; among them,

[0025]

[0026] Where: is the estimated autocorrelation coefficient of the i-th period; is the partial correlation coefficient estimate of the i-th period; N is the length of the sequence; M is the maximum lag period of the autocorrelation coefficient;

[0027] If formula (1) is established, it can be considered that the autocorrelation coefficient of the sequence is truncated at q steps, and the MA model can be established, and the initial order is q; if formula (2) is established, it can be considered that the partial correlation coefficient of the sequence is truncated at p steps, and the AR model can be established, and the initial order is p;

[0028] Then, the AIC criterion and the BIC criterion were used to determine the final order of the model.

[0029] In one example of the present invention, the mathematical expressions of the AIC criterion and the BIC criterion are respectively:

[0030] AIC(p,q)=nln(σ ε 2 )+2(p+q+1)

[0031] BIC(p,q)=nln(σ ε 2 )+ln(n)(p+q+1)

[0032] Where: n is the sequence length; σ ε 2 is the maximum likelihood estimate of the model residual variance.

[0033] In one example of the present invention, in step S32, performing adaptability test on the model includes the following steps:

[0034] Using χ 2 Test method, constructing statistics Q to test the adaptability of the model only requires white noise test of the model residuals. The expression of statistics Q is as follows:

[0035]

[0036] Where N is the number of samples, p ε (k) is the autocorrelation coefficient of the residual sequence, (p,q) is the model order, and M is the maximum lag period of the autocorrelation coefficient;

[0037] If the model is appropriate, the statistic Q approximately follows χ 2 (Mpq); at a given significance level a, if Q≤χ 2 (Mpq), the residual sequence is considered to be a white noise sequence, and the fitting model is significant; if Q>χ 2 (Mpq), then the residual sequence is considered not to be a white noise sequence and the fitting model is not significant.

[0038] In one example of the present invention, in step S32, when the ARMA model is selected to fit the maximum daily energy of sample microseismic events, the modeling steps are as follows:

[0039] Determine whether the original microseismic time series data is stable, that is, perform a stability test on the original data. If the series does not meet the stability condition, it can be made stable through differential transformation, logarithmic transformation or other transformations;

[0040] By calculating the statistics AC (autocorrelation) and PAC (partial autocorrelation) that describe the characteristics of the sequence, the order p and q values ​​of the ARMA model of the post-stationary microseismic time series are determined;

[0041] Estimate unknown parameters of the model, test the significance of the parameters and the rationality of the model;

[0042] Diagnostic analyses were performed to confirm that the resulting model matched the observed data characteristics.

[0043] In one example of the present invention, in the step S30, before establishing the time series model, it also includes: first, a white noise test is performed on the two groups of time series. If the test result is a non-white noise sequence, the stationarity of the microseismic data is judged. If the time series is non-stationary, it should be converted into a stationary sequence before establishing the model.

[0044] Another object of the present invention is to provide an intelligent advanced targeted non-explosive fracturing and rockburst mitigation system, comprising:

[0045] An activity data acquisition device is configured to arrange a monitoring system in a rockburst prone area within a predetermined range from a tunnel face to obtain activity data of the rock mass, wherein the monitoring system includes: a displacement sensor, a stress sensor, a temperature sensor, an acoustic wave sensor, and a vibration sensor;

[0046] A microseismic data acquisition device configured to connect the installed monitoring system to the microseismic monitoring equipment through a data acquisition line, so as to collect the rock mass activity data acquired by the monitoring system to obtain microseismic monitoring data;

[0047] A microseismic data analysis device is configured to connect the microseismic analysis device and the early warning device to the microseismic monitoring device for communication, establish a time series model through the microseismic analysis device, analyze the microseismic monitoring data, and predict the potential rock burst area in the rock burst prone area near the tunnel construction face based on the microseismic analysis results by the alarm device, and classify the potential rock burst area;

[0048] The rockburst elimination device is configured to open non-explosion boreholes in the potential rockburst area in a targeted manner according to the early warning signal, and inject low-wave impedance proppants into the non-explosion fracture boreholes to eliminate potential rockburst disasters.

[0049] In one example of the present invention, the microseismic data analysis device includes:

[0050] A stationary analysis unit is configured to use the historical and current data of the maximum energy of a single day of microseismic in the monitoring area as samples to perform stationary analysis, and if it is not stationary, perform a differential operation to obtain an approximately stationary sequence;

[0051] Adaptability test unit: configured to select an appropriate model to fit the maximum energy of a single day of sample microseismic events, determine the order and estimate the parameters of the model, and perform adaptability test on the model;

[0052] Prediction and warning unit: It is configured to use the established model to predict the future value of the difference, and then obtain the future value of the maximum energy of the microseismic day through the inverse difference operation. During the mining process, the established model is used to perform sample extrapolation prediction, and by comparing with the warning indicators, danger warning can be achieved.

[0053] The best embodiment for carrying out the present invention will be described in more detail below with reference to the accompanying drawings so that the features and advantages of the present invention can be easily understood. BRIEF DESCRIPTION OF THE DRAWINGS

[0054] In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the following briefly introduces the drawings of the embodiments of the present invention, wherein the drawings are only used to illustrate some embodiments of the present invention, but not to limit all embodiments of the present invention thereto.

[0055] Figure 1 Schematic diagram of the structure of intelligent advanced targeted non-explosive fracturing to mitigate rockburst according to an embodiment of the present invention;

[0056] Figure 2 This is a flow chart of a method for intelligent advanced targeted non-explosive fracturing to mitigate rockburst according to an embodiment of the present invention;

[0057] Figure 3 is a time series diagram of the maximum energy of a microseism in a single day according to an embodiment of the present invention;

[0058] Figure 4 This is a diagram of the unit root test results of the sequence Energy according to an embodiment of the present invention;

[0059] Figure 5 is a timing diagram after a primary difference according to an embodiment of the present invention;

[0060] Figure 6 A diagram showing the ADF test result of the first-order difference sequence according to an embodiment of the present invention;

[0061] Figure 7 is the fitting result of the ADF test model of the first-order difference sequence according to an embodiment of the present invention;

[0062] Figure 8 is a graph of an original sequence, a first-order difference sequence and a second-order difference sequence according to an embodiment of the present invention;

[0063] Fig. 9 The microseismic single-day maximum energy stationary series correlation diagram and partial correlation diagram according to an embodiment of the present invention;

[0064] Fig.10 A scatter plot of the correlation between the original sequence and the lagged sequence according to an embodiment of the present invention;

[0065] Fig.11 is the regression result of the MA(1) model according to an embodiment of the present invention;

[0066] Fig.12 The MA(1) model fitting effect according to an embodiment of the present invention;

[0067] Fig.13 The MA(1) residual curve and QQ graph according to an embodiment of the present invention;

[0068] Fig.14 1 is a time series diagram of the logarithmic difference transformation of the maximum energy of a microseismic day according to an embodiment of the present invention. DETAILED DESCRIPTION

[0069] In order to make the purpose, technical solution and advantages of the technical solution of the present invention clearer, the technical solution of the embodiment of the present invention will be clearly and completely described in conjunction with the drawings of specific embodiments of the present invention. The same figure marks in the drawings represent the same parts. It should be noted that the described embodiments are part of the embodiments of the present invention, not all of the embodiments. Based on the described embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0070] Unless otherwise defined, the technical terms or scientific terms used herein shall have the usual meanings understood by persons of ordinary skill in the field to which the present invention belongs. The words "first", "second" and similar words used in the patent application specification and claims of the present invention do not indicate any order, quantity or importance, but are only used to distinguish different components. Similarly, words such as "one" or "one" do not necessarily indicate a quantity limitation. Words such as "include" or "comprise" mean that the elements or objects appearing before the word include the elements or objects listed after the word and their equivalents, without excluding other elements or objects. Words such as "connect" or "connected" are not limited to physical or mechanical connections, but may include electrical connections, whether direct or indirect. "Up", "down", "left", "right" and the like are only used to indicate relative positional relationships. When the absolute position of the described object changes, the relative positional relationship may also change accordingly.

[0071] According to a first aspect of the present invention, a method for intelligent advanced targeted non-explosive fracturing to mitigate rockbursts comprises the following steps:

[0072] S10: arranging a monitoring system in a rockburst prone area within a predetermined range from the tunnel face to obtain activity data of the rock mass, wherein the monitoring system includes: a displacement sensor, a stress sensor, a temperature sensor, an acoustic wave sensor, and a vibration sensor;

[0073] S20: connecting the installed monitoring system to the microseismic monitoring equipment through a data acquisition line to collect the rock mass activity data acquired by the monitoring system to obtain microseismic monitoring data;

[0074] S30: connecting the microseismic analysis equipment and the early warning equipment to the microseismic monitoring equipment, establishing a time series model through the microseismic analysis equipment, analyzing the microseismic monitoring data, and using the alarm equipment to predict the potential rock burst area in the rock burst prone area near the tunnel construction face based on the microseismic analysis results, and classifying the potential rock burst area;

[0075] S40: opening non-explosion boreholes in the potential rock burst area in a targeted manner according to the early warning signal, and injecting low wave impedance proppants into the non-explosion fracture boreholes to eliminate potential rock burst disasters.

[0076] This method can not only accurately predict the potential rock burst area in the rock burst prone area near the tunnel construction face, but also actively reduce the stress concentration in deep hard rock tunnels and mining working faces, weaken stress wave disturbances, and has high accuracy, greatly reducing the manufacturing cost and construction difficulty of equipment. It plays an important guiding role in hard rock strata excavation and rock burst prediction and control, and is highly practical.

[0077] In an example of the present invention, in step S30, establishing a time series model includes the following steps:

[0078] S31: Use the historical and current data of the maximum energy of microseismic events in the monitoring area in a single day as samples to perform a stationary analysis. If the data is not stationary, perform a differential operation to obtain an approximately stationary sequence.

[0079] S32: Select an appropriate model to fit the maximum energy of the sample microseismic event in a single day, determine the order and estimate the parameters of the model, and perform an adaptability test on the model;

[0080] S33: Use the established model to predict the future value of the difference, and then obtain the future value of the maximum energy of the microseismic day through the reverse difference operation. During the mining process, use the established model to perform sample extrapolation prediction, and by comparing with the warning indicators, danger warning can be achieved.

[0081] In one example of the present invention, in step S31, the stationarity analysis includes the following steps:

[0082] First fit the adaptive model of the sequence, and then test whether there is a unit root in the characteristic equation composed of the autoregressive part parameters of the model. If there is no unit root, the sequence is stable, and if there is a unit root, the sequence is not stable.

[0083] The specific process is as follows: for the AR(p) process

[0084] X t =φ1X t-1 +φ2X t-2 +···+φ p X t-p +ε t

[0085] make The following hypothesis test is performed:

[0086]

[0087] Construct the ADF test statistic: in is the sample standard deviation of parameter ρ. ARMA(p,q) is the same as AR(p). If it is MA(q), the regression part has no parameters and is definitely stable.

[0088] In one example of the present invention, in step S32, determining the order of the model includes the following steps:

[0089] First, by judging whether formula (1) and formula (2) are valid, we can roughly judge the truncation of the autocorrelation coefficient of the stationary series and the truncation of the partial correlation coefficient of the stationary series and make a preliminary order determination; among them,

[0090]

[0091] Where: is the estimated autocorrelation coefficient of the i-th period; is the partial correlation coefficient estimate of the i-th period; N is the length of the sequence; M is the maximum lag period of the autocorrelation coefficient;

[0092] If formula (1) is established, it can be considered that the autocorrelation coefficient of the sequence is truncated at q steps, and the MA model can be established, and the initial order is q; if formula (2) is established, it can be considered that the partial correlation coefficient of the sequence is truncated at p steps, and the AR model can be established, and the initial order is p;

[0093] Then, the AIC criterion and the BIC criterion were used to determine the final order of the model.

[0094] In one example of the present invention, the mathematical expressions of the AIC criterion and the BIC criterion are respectively:

[0095] AIC(p,q)=nln(σ ε 2 )+2(p+q+1)

[0096] BIC(p,q)=nln(σ ε 2 )+ln(n)(p+q+1)

[0097] Where: n is the sequence length; σ ε 2 is the maximum likelihood estimate of the model residual variance.

[0098] In one example of the present invention, in step S32, performing adaptability test on the model includes the following steps:

[0099] Using χ 2 Test method, constructing statistics Q to test the adaptability of the model only requires white noise test of the model residuals. The expression of statistics Q is as follows:

[0100]

[0101] Where N is the number of samples, p ε (k) is the autocorrelation coefficient of the residual sequence, (p,q) is the model order, and M is the maximum lag period of the autocorrelation coefficient;

[0102] If the model is appropriate, the statistic Q approximately follows χ 2 (Mpq); at a given significance level a, if Q≤χ 2 (Mpq), the residual sequence is considered to be a white noise sequence, and the fitting model is significant; if Q>χ 2 (Mpq), then the residual sequence is considered not to be a white noise sequence and the fitting model is not significant.

[0103] In one example of the present invention, in step S32, when the ARMA model is selected to fit the maximum daily energy of sample microseismic events, the modeling steps are as follows:

[0104] Determine whether the original microseismic time series data is stable, that is, perform a stability test on the original data. If the series does not meet the stability condition, it can be made stable through differential transformation, logarithmic transformation or other transformations;

[0105] By calculating the statistics AC (autocorrelation) and PAC (partial autocorrelation) that describe the characteristics of the sequence, the order p and q values ​​of the ARMA model of the post-stationary microseismic time series are determined;

[0106] Estimate unknown parameters of the model, test the significance of the parameters and the rationality of the model;

[0107] Diagnostic analyses were performed to confirm that the resulting model matched the observed data characteristics.

[0108] In one example of the present invention, in the step S10, before establishing the time series model, it also includes: first, a white noise test is performed on the two groups of time series. If the test result is a non-white noise sequence, the stationarity of the microseismic data is judged. If the time series is non-stationary, it should be converted into a stationary sequence before establishing the model.

[0109] According to a second aspect of the present invention, an intelligent advanced targeted non-explosive fracturing and rockburst mitigation system comprises:

[0110] An activity data acquisition device is configured to arrange a monitoring system in a rockburst prone area within a predetermined range from a tunnel face to obtain activity data of the rock mass, wherein the monitoring system includes: a displacement sensor, a stress sensor, a temperature sensor, an acoustic wave sensor, and a vibration sensor;

[0111] A microseismic data acquisition device configured to connect the installed monitoring system to the microseismic monitoring equipment through a data acquisition line, so as to collect the rock mass activity data acquired by the monitoring system to obtain microseismic monitoring data;

[0112] A microseismic data analysis device is configured to connect the microseismic analysis device and the early warning device to the microseismic monitoring device for communication, establish a time series model through the microseismic analysis device, analyze the microseismic monitoring data, and predict the potential rock burst area in the rock burst prone area near the tunnel construction face based on the microseismic analysis results by the alarm device, and classify the potential rock burst area;

[0113] The rockburst elimination device is configured to open non-explosion boreholes in the potential rockburst area in a targeted manner according to the early warning signal, and inject low-wave impedance proppants into the non-explosion fracture boreholes to eliminate potential rockburst disasters.

[0114] The system can not only accurately predict potential rock burst areas in rock burst prone areas near the tunnel construction face, but also actively reduce stress concentration in deep hard rock tunnels and mining working faces, weaken stress wave disturbances, and has high accuracy, greatly reducing the manufacturing cost and construction difficulty of the equipment. It plays an important guiding role in hard rock strata excavation and rock burst prediction and control, and is highly practical.

[0115] In one example of the present invention, the microseismic data analysis device includes:

[0116] A stationary analysis unit is configured to use the historical and current data of the maximum energy of a single day of microseismic in the monitoring area as samples to perform stationary analysis, and if it is not stationary, perform a differential operation to obtain an approximately stationary sequence;

[0117] Adaptability test unit: configured to select an appropriate model to fit the maximum energy of a single day of sample microseismic events, determine the order and estimate the parameters of the model, and perform adaptability test on the model;

[0118] Prediction and warning unit: It is configured to use the established model to predict the future value of the difference, and then obtain the future value of the maximum energy of the microseismic day through the inverse difference operation. During the mining process, the established model is used to perform sample extrapolation prediction, and by comparing with the warning indicators, danger warning can be achieved.

[0119] In one example of the present invention, in the step S40, a low wave impedance proppant is injected into the non-explosive fracturing borehole by a hydraulic fracturing device;

[0120] Hydraulic fracturing equipment includes high-pressure pumps, nozzles, and piping systems. The precise arrangement of these devices ensures that the proppant can enter the advance borehole at the appropriate speed and direction to ensure that the cracks and pores in the rock are filled. The proppant usually adopts low wave impedance proppant, and these characteristics will be described in detail in the next section.

[0121] The high-pressure pump is one of the core components of the whole process. It is responsible for injecting the proppant into the advance borehole at high pressure to achieve rock fracturing. The high-pressure pump is equipped with a programmable pressure control system to allow the injection pressure of the proppant to be adjusted under different rock conditions. In addition, the design and arrangement of the nozzle ensure that the proppant enters the advance borehole at an appropriate speed and angle, minimizing the presence of unfilled areas.

[0122] In one example of the present invention, the proppant has the characteristic of low wave impedance, which helps to absorb the stress waves propagating in the rock mass and reduce the disturbance of the stress waves. The type of proppant selected and the adjustment of the parameters are crucial to ensure that the proppant can fill the cracks and pores in the rock, thereby enhancing the overall stability of the rock mass. In addition, the particle distribution and chemical properties of the proppant are precisely designed to ensure uniform distribution of the proppant during hydraulic fracturing and provide a lasting support effect. In addition, the physical properties and chemical composition of the proppant can be adjusted according to different rock mass types and engineering requirements.

[0123] Specific examples

[0124] When using the SOS microseismic monitoring system for ground pressure monitoring and early warning, in general, the following impact hazard indicators are adopted, and the monitoring data are continuously summarized and analyzed during the production process to formulate scientific and reasonable impact hazard indicators.

[0125] (1) During excavation, vibration energy exceeds 10 4 J, it can be predicted that there is a risk of impact in this area.

[0126] (2) During mining, vibration energy exceeds 10 5 J, it can be predicted that there is a risk of impact in this area.

[0127] (3) When the number of vibrations and the vibration energy increase at the same time, it is a sign of a dangerous increase in ground pressure.

[0128] (4) When the number of vibrations decreases but the vibration energy increases, it is a sign of dangerously rising ground pressure. The initial stages of several production cycles (mining working face advances at least 10m) when the number of vibrations begins to decrease may be particularly dangerous, indicating energy accumulation in the rock formation. The KJ551 high-precision microseismic monitoring system is used for ground pressure monitoring and early warning, and its early warning indicators are shown in Table 1:

[0129] Table 1 KJ551 dangerous microseismic evaluation index

[0130]

[0131] Note: The indicators in the table are reference indicators for microseismic early warning, which need to be improved in practical applications. For example, the relationship between the total number of daily microseismic events, the average energy of daily microseismic events and the ground pressure danger needs further study; E max is the maximum energy of daily microseismic events.

[0132] From the above analysis, it can be seen that when the energy of a single microseismic event on a single day exceeds the critical value or the cumulative energy increases within a cycle, it indicates that the possibility of ground pressure danger is higher. However, predicting the impact danger by the cumulative change of energy within a cycle is not intuitive and operational. A single microseismic energy event on a single day can meet this requirement, and high-energy microseismic events are also a necessary condition for the occurrence of ground pressure. The change in microseismic frequency (increase or decrease) plays a qualitative analysis role in the occurrence of ground pressure, and cannot be quantitatively identified and warned. Therefore, the maximum energy of a microseismic event on a single day is selected as the warning indicator for microseismic monitoring impact danger.

[0133] (1) ARIMA model building

[0134] The ARIMA (p, d, q) model is the abbreviation of the Auto Regressive Integrated Moving Average Model. It combines moving average, autoregressive analysis and difference. The model needs to determine three parameters, namely the autoregressive order (p), the number of differences (d) and the moving average order (q). Before building the model, it is necessary to determine whether the sequence has randomness, stationary and seasonality, and then to achieve the best fit before prediction, and finally to make predictions and evaluations. The specific modeling steps are as follows:

[0135] (1) Determine whether the original microseismic time series data is stable, that is, perform a stability test on the original data. If the series does not meet the stability condition, it can be made stable through differential transformation, logarithmic transformation or other transformations;

[0136] (2) By calculating the statistics AC (autocorrelation) and PAC (partial autocorrelation) diagrams that describe the characteristics of the sequence, the order p and q values ​​of the ARMA model of the post-stationary microseismic time series are determined;

[0137] (3) Estimate the unknown parameters of the model and test the significance of the parameters and the rationality of the model;

[0138] (4) Perform diagnostic analysis to confirm that the obtained model is consistent with the observed data characteristics;

[0139] (5) In steps (3) and (4), some statistics and tests are needed to analyze whether the model form chosen in step (2) is appropriate.

[0140] (2) Microseismic monitoring data analysis and processing

[0141] The microseismic monitoring data of the -290 middle working face of the mine from March 1 to July 10 were selected for sorting and analysis. The specific data are listed in Table 3-2. The data were imported into the database through the Pandas tool using Python software, and the results were shown in Table 2.

[0142] Table 2 Microseismic monitoring data of the 290 middle section working face from March 1 to July 10

[0143]

[0144]

[0145]

[0146] (3) Determination of stationarity of time series data

[0147] For time series data, the most important tests are whether the time series data is white noise data, whether the time series data is stable, and the analysis of the autocorrelation coefficient and partial autocorrelation coefficient of the time series data. If the time series data is white noise data, it means that it does not have any useful information. Many analysis methods for time series data require that the time series data under study is stable, so it is very important to determine whether the time series data is stable and how to convert non-stationary time series data into stationary series data. Modeling research on time series data is very important.

[0148] The premise of establishing the ARMA time series model is that the sequence to be analyzed must be stationary and non-white noise, because only non-white noise sequences can predict the future based on historical data. Therefore, the two sets of time series must first be tested for white noise. If the test result is a non-white noise sequence, the stationarity of the microseismic data is then judged. If the time series is non-stationary, it should be converted into a stationary sequence before building the model. This paper selects data from March 1 to July 7 as sample data for the training set for Python software model establishment, and data from July 8 to July 10 as the test set to test the accuracy of the test model.

[0149] (1) Time series white noise test

[0150] In the Python Jupyter Notebook work file window, you can use the Matplotlib tool to draw graphs. The data has been imported into Python before. You can use the plot() function in the Matplotlib tool to visualize the time series and obtain the time series diagram of the maximum energy of a microseismic day, as shown in the figure below: Figure 3 shown.

[0151] If a sequence is white noise (i.e., independent and identically distributed random data), then there is no need to establish a time series model for prediction, because predicting random numbers is meaningless.

[64] . Therefore, before establishing time series analysis, it is necessary to conduct a white noise test on it. The commonly used white noise test method is the Ljung-Box test (LB test for short), whose null hypothesis and alternative hypothesis are H0: the sequences with a delay period less than or equal to m are independent of each other (the sequence is white noise); H1: there is a correlation between the sequences with a delay period less than or equal to m (the sequence is not white noise). The Ljung-Box test can use the acorr_ljungbox() function that comes with the statsmodels tool in the Python software to perform a white noise test on the time series. The test results are shown in Table 3:

[0152] Table 3 White noise test results of sequence Energy

[0153]

[0154] From the above results, we can see that when the delay order is [4,8,16,32], the LB test p-value of the sequence Energy is less than 0.05, indicating that the null hypothesis that the sequence is white noise can be rejected, and it is believed that the data is not random data, that is, the data is not random, but has regular patterns and is worth analyzing.

[0155] (2) Time series stationarity test

[0156] Whether the time series is stationary is critical to the selection of a mathematical model for prediction. If a set of time series data is stationary, you can directly use the autoregressive moving average model (ARMA) for prediction. If the data is not stationary, you need to try to establish a differential autoregressive moving average model (ARIMA) for prediction.

[0157] There are two ways to determine whether a sequence is stationary: one is to make a judgment based on the characteristics shown in the time series diagram and autocorrelation diagram; the other is to construct a test statistic for hypothesis testing, such as a unit root test. The first method is more subjective, while the second is an objective method.

[0158] When a time series is a stationary series, its mean and variance are constants and the time series graph always fluctuates randomly around a constant. Figure 4 It can be seen that the mean of the maximum energy of a microseismic event in a single day decreases gradually with the increase of time. Therefore, it is preliminarily judged that the time series of the maximum energy of a microseismic event in a single day is non-stationary.

[0159] The adfuller() function of the statsmodels.tsa module in Python can be used to perform an ADF unit root test. The results of the unit root test for the Energy sequence are as follows: Figure 4 As shown:

[0160] From the output results of the unit root test above, we can find that the test P value of the Energy sequence is 0.172846, which is greater than 0.05, indicating that the Energy sequence is not stable. Combined with the previous white noise test results, it can be concluded that the Energy sequence is a non-white noise non-stationary time series. Therefore, it is necessary to process the non-stationary time series.

[0161] (3) Non-stationary time series processing

[0162] Difference is the most commonly used method for smoothing non-stationary time series data, especially in the ARIMA model. By using the diff() function in Python software, we can directly perform a difference operation on the Energy sequence to obtain the first-order difference sequence DEnergy. Then, by using the plot() function above, we can visualize the difference sequence and generate a time series diagram of the maximum energy of a single day of microseismic earthquakes after a difference, as shown below: Figure 5 shown.

[0163] from Figure 6 It can be seen that the time series diagram after difference basically fluctuates around the value of "0", and it is preliminarily determined that the time series after difference is a stationary series. The following also uses the quantitative method ADF test to further determine whether the sequence after difference is really a stationary series. The test results are as follows:

[0164] Combination Figure 7 The output results of the ADF test and the ADF test model fitting results show that the test P value of the first-order difference sequence DEnergy is 1.063001e-8, which is much smaller than 0.05. Figure 7Comparing the value of "Test Statistic" in the chart with the value of "Critical Value (5%)" in the chart, the result is -6.517377, which is less than -2.886151 (and also less than -3.486535 at the 1% significance level), indicating that the first-order difference sequence DEnergy under test is stable. Combining the previous test results, we can conclude that the original sequence Energy is not stable, and the first-order difference sequence DEnergy is stable. In order to further confirm the most suitable number of differences for the original sequence Energy, we perform another difference based on the first-order difference, and draw a curve of the second-order difference, as shown below: Figure 8 As shown in the figure, after the second-order difference is performed on the original sequence, the fluctuation range of the time series curve is obviously larger than that of the first-order, which will inevitably lead to a larger variance of the fitting result. Therefore, it is not appropriate to use the second-order difference for the original sequence Energy. On the whole, the microseismic single-day maximum energy sequence is a first-order difference stationary sequence.

[0165] (3) Model identification and verification

[0166] (1) Model identification

[0167] Autocorrelation analysis and partial autocorrelation analysis are methods used to determine the two parameters p and q in the ARMA (p, q) model. After determining that the sequence is a stationary non-white noise sequence, the truncation of the sequence can be analyzed by the values ​​of the autocorrelation coefficient and partial autocorrelation coefficient of the sequence.

[0168] For a time series, if the sample autocorrelation coefficient ACF is not equal to 0 until the lag period s = q, and the ACF is almost 0 when the lag period s> q, then it can be considered that the real data generation process is MA (q). If the sample partial autocorrelation coefficient PACF is not equal to 0 until the lag period s = p, and the PACF is almost 0 when the lag period s> p, then it can be considered that the real data generation process is AR (p). More generally, based on the performance of the sample ACF and PACF, a more appropriate ARMA (p, q) model can be fitted. Table 4 shows how to determine the parameters p and q in the model.

[0169] Table 4 Model selection discriminant table

[0170]

[0171]

[0172] For the autocorrelation coefficient and partial autocorrelation coefficient of the time series, you can use the plot_acf() function and plot_pacf() function to visualize them. Run the following program to obtain the autocorrelation coefficient and partial autocorrelation coefficient of the time series DEnergy. The results are as follows: Fig. 9 shown.

[0173] Depend on Fig. 9 It can be seen that the autocorrelation diagram is first-order truncated (the coefficient quickly approaches 0), and the partial autocorrelation diagram is tailed (the coefficient slowly approaches 0). The initial identification of the microseismic single-day maximum energy time series model is the ARIMA (0,1,1) model, that is, the MA (1) model. In addition, the scatter plot of the first-order difference sequence can be used to show the correlation of the microseismic single-day maximum energy, such as Fig.10 shown.

[0174] exist Fig.10 In the figure, the ordinate "Original" represents the original sequence of the first-order difference of the maximum energy of a microseismic day, the abscissa "tn" represents the time series of n-order lags based on the original sequence, and the title "corr" of each scatter plot is the correlation coefficient between the original sequence and the lagged sequence. The 9 scatter plot images as a whole represent the correlation between the original sequence and the time series of lags 1-9. It can be obtained that when n = 1, the maximum correlation coefficient between the original sequence and the lagged sequence is corr = 0.37, so it can be concluded that q = 1 in the corresponding ARIMA (p, d, q).

[0175] (2) Model verification

[0176] The MA(1) model regression results are as follows: Fig.11 As shown. Fig.11 It can be seen that the t-statistic of the constant term coefficient (column "z") is -1.160, with an absolute value less than 2, and the probability value corresponding to the "P>|z|" column is 0.246, which is greater than 0.05. The t-statistic of the explanatory variable coefficient of MA(1) has an absolute value greater than 2, which means that there is at least a 95% probability that the true value of the coefficient is not zero. The probability value corresponding to the "P>|z|" column is 0, which is less than 0.05, indicating that the parameter estimates are valid, that is, the explanatory variable estimated coefficients are significantly different from zero, and the MA(1) model parameters have passed the significance test.

[0177] from Fig.11 From the last line of the output, we can see that the characteristic root is 1.0896>1, which meets the requirement of stationarity. Using the fitted values ​​method, we can get the predicted value of the MA(1) model by fitting the ARIMA model, and then make a prediction curve. By putting the original time series curve and the prediction curve in the same coordinate system, we can intuitively show the fitting effect, such as Fig.12At the same time, the resid method is used to calculate the residual value of the fitting result, and the plot() function and qqplot() function are used in Python software to draw the corresponding residual curve graph and QQ graph, as shown in Fig.13 shown.

[0178] exist Fig.12 In the figure, the “diff” curve is the first-order difference sequence of the maximum energy value of a microseismic day, and the “forecast_diff” curve is the prediction curve of the MA(1) model. By comparison, it can be seen that the fitting effect of the MA(1) model is average, and the prediction curve and the sequence DEnergy curve partially overlap. Fig.13 It can be found that the residual curve of MA(1) fluctuates greatly, with the highest value even exceeding 20000. From the MA(1) residual QQ graph, it can be concluded that the model fitting quality is average, with some points located near the straight line and some points far away from the straight line. Therefore, the random error term of the model is a white noise sequence.

[0179] Although the MA(1) model passed the test, the goodness-of-fit value R calculated by Python software was 2 It is only 0.455, and the regression effect is poor and needs to be improved.

[0180] (4) Model modification

[0181] As can be seen from the previous text, the time series of the maximum energy of a microseismic day is a non-stationary series. After a simple difference operation, the series can be stabilized. Further observation of the original time series shows that with the increase of time, the maximum energy of a day shows a downward trend, which makes the variance larger and larger, that is, more and more unstable. In this case, the common analysis assumptions often do not meet the requirements of the normal distribution of errors obeying independent and identical distributions and the stability of the time series. This inevitably leads us to seek a way to make the data meet the assumptions as much as possible, keep the variance constant, that is, make the fluctuations relatively stable, and this purpose can be achieved through logarithmic transformation. Logarithmic transformation is a common way of data transformation. The purpose of data transformation is to make the presentation of data close to the premise assumptions we want, so as to better perform statistical inference.

[0182] The logarithmic transformation process of the microseismic daily maximum energy time series is similar to the difference process, so it will not be demonstrated here. DLnEnergy time series after logarithmic transformation and difference.

[0183] from Fig.14It can be seen that the fluctuation range of the microseismic daily maximum energy time series after taking the logarithmic first difference is between [-3,3]. Compared with the time series of direct difference, the fluctuation range is significantly reduced. This ensures that the variance is constant, that is, the fluctuation is relatively stable. The results of the model single root test are shown in Table 5. The Test Statistic value of ADF is significantly smaller than the value of "Critical Value (5%)", and is also smaller than the value at the 1% significance level. The p value is 3.022895e-11, which is much smaller than 0.05 and smaller than the p value of the DEnergy series. This shows that the tested series is stable and has better stability than the DEnergy series.

[0184] Table 5 Single root test results of DLnEnergy sequence model

[0185]

[0186] This method uses a microseismic monitoring system for rockburst early warning, and through comparative studies of microseismic event sequence models and various rockburst early warning models, the following conclusions are drawn:

[0187] (1) Use the historical and current data of the maximum daily energy of microseismic events in the monitoring area as samples to perform a stationary analysis. If the data is not stationary, perform a differential operation to obtain an approximately stationary sequence. Select an appropriate model to fit the maximum daily energy of the sample microseismic events, estimate the parameters of the model, and perform an adaptability test on the model. If the model is not adaptable, reselect the model to fit the maximum daily energy of the sample microseismic events.

[0188] (2) Finally, the established model is used to predict the future value of the difference, and then the future value of the maximum energy of the microseismic day is obtained through the reverse difference operation. During the mining process, the established model is used to perform sample extrapolation prediction, and by comparing it with the warning indicators, danger warning can be achieved.

[0189] The exemplary implementation scheme of the method for mitigating rockburst by intelligent advanced targeted non-explosive fracturing proposed in the present invention is described in detail above with reference to the preferred embodiments. However, those skilled in the art will appreciate that, without departing from the concept of the present invention, various modifications and variations can be made to the above-mentioned specific embodiments, and various technical features and structures proposed in the present invention can be combined in various ways without exceeding the protection scope of the present invention, which is determined by the attached claims.

Claims

1. A method for intelligent advanced targeted non-explosive fracturing to mitigate rockburst, characterized in that: The steps include: S10: arranging a monitoring system in a rockburst prone area within a predetermined range from the tunnel face to obtain activity data of the rock mass, wherein the monitoring system includes: a displacement sensor, a stress sensor, a temperature sensor, an acoustic wave sensor, and a vibration sensor; S20: connecting the installed monitoring system to the microseismic monitoring equipment through a data acquisition line to collect the rock mass activity data acquired by the monitoring system to obtain microseismic monitoring data; S30: connecting the microseismic analysis equipment and the early warning equipment to the microseismic monitoring equipment, establishing a time series model through the microseismic analysis equipment, analyzing the microseismic monitoring data, and using the alarm equipment to predict the potential rock burst area in the rock burst prone area near the tunnel construction face based on the microseismic analysis results, and classifying the potential rock burst area; S40: opening non-explosion boreholes in the potential rock burst area in a targeted manner according to the early warning signal, and injecting low wave impedance proppants into the non-explosion fracture boreholes to eliminate potential rock burst disasters.

2. The method for reducing rockburst by intelligent advanced targeted non-explosive fracturing according to claim 1 is characterized in that: In step S30, establishing a time series model includes the following steps: S31: Use the historical and current data of the maximum energy of microseismic events in the monitoring area in a single day as samples to perform a stationary analysis. If the data is not stationary, perform a differential operation to obtain an approximately stationary sequence. S32: Select an appropriate model to fit the maximum energy of the sample microseismic event in a single day, determine the order and estimate the parameters of the model, and perform an adaptability test on the model; S33: Use the established model to predict the future value of the difference, and then obtain the future value of the maximum energy of the microseismic day through the reverse difference operation. During the mining process, use the established model to perform sample extrapolation prediction, and by comparing with the warning indicators, danger warning can be achieved.

3. The method for reducing rockburst by intelligent advanced targeted non-explosive fracturing according to claim 2 is characterized in that: In step S31, the stability analysis includes the following steps: First fit the adaptive model of the sequence, and then test whether there is a unit root in the characteristic equation composed of the autoregressive part parameters of the model. If there is no unit root, the sequence is stable, and if there is a unit root, the sequence is not stable. The specific process is as follows: for the AR(p) process X t =φ1X t-1 +φ2X t-2 +···+φ p X t-p +e t make The following hypothesis test is performed: Construct the ADF test statistic: in is the sample standard deviation of parameter ρ. ARMA(p,q) is the same as AR(p). If it is MA(q), the regression part has no parameters and is definitely stable.

4. The method for reducing rockburst by intelligent advanced targeted non-explosive fracturing according to claim 2 is characterized in that: In step S32, determining the order of the model includes the following steps: First, by judging whether formula (1) and formula (2) are valid, we can roughly judge the truncation of the autocorrelation coefficient of the stationary series and the truncation of the partial correlation coefficient of the stationary series and make a preliminary order determination; among them, Where: is the estimated autocorrelation coefficient of the i-th period; is the partial correlation coefficient estimate of the i-th period; N is the length of the sequence; M is the maximum lag period of the autocorrelation coefficient; If formula (1) is established, it can be considered that the autocorrelation coefficient of the sequence is truncated at q steps, and the MA model can be established, and the initial order is q; if formula (2) is established, it can be considered that the partial correlation coefficient of the sequence is truncated at p steps, and the AR model can be established, and the initial order is p; Then, the AIC criterion and the BIC criterion were used to determine the final order of the model.

5. The method for reducing rockburst by intelligent advanced targeted non-explosive fracturing according to claim 4 is characterized in that: The mathematical expressions of the AIC criterion and the BIC criterion are: AIC(p,q)=nln(σ ε 2 )+2(p+q+1) BIC(p,q)=nln(σ ε 2 )+ln(n)(p+q+1) Where: n is the sequence length; σ ε 2 is the maximum likelihood estimate of the model residual variance.

6. The method for reducing rockburst by intelligent advanced targeted non-explosive fracturing according to claim 2 is characterized in that: In step S32, the adaptability test of the model includes the following steps: Using χ 2 Test method, constructing statistics Q to test the adaptability of the model only requires white noise test of the model residuals. The expression of statistics Q is as follows: Where N is the number of samples, p ε (k) is the autocorrelation coefficient of the residual sequence, (p,q) is the model order, and M is the maximum lag period of the autocorrelation coefficient; If the model is appropriate, the statistic Q approximately follows χ 2 (Mpq); at a given significance level a, if Q≤χ 2 (Mpq), the residual sequence is considered to be a white noise sequence, and the fitting model is significant; if Q>χ 2 (Mpq), then the residual sequence is considered not to be a white noise sequence and the fitting model is not significant.

7. The method for reducing rockburst by intelligent advanced targeted non-explosive fracturing according to claim 2 is characterized in that: In step S32, when the ARMA model is selected to fit the maximum daily energy of sample microseismic events, the modeling steps are as follows: Determine whether the original microseismic time series data is stable, that is, perform a stability test on the original data. If the series does not meet the stability condition, it can be made stable through differential transformation, logarithmic transformation or other transformations; By calculating the statistics AC and PAC describing the characteristics of the sequence, the order p and q values ​​of the ARMA model of the microseismic time series after stationarity are determined; Estimate unknown parameters of the model, test the significance of the parameters and the rationality of the model; Diagnostic analyses were performed to confirm that the resulting model matched the observed data characteristics.

8. The method for reducing rockburst by intelligent advanced targeted non-explosive fracturing according to claim 1 is characterized in that: In the step S30, before establishing the time series model, it also includes: first, a white noise test is performed on the two groups of time series. If the test result is a non-white noise sequence, the stationarity of the microseismic data is judged. If the time series is non-stationary, it should be converted into a stationary sequence before establishing the model.

9. An intelligent advanced targeted non-explosive fracturing and rockburst mitigation system, characterized in that: include: An activity data acquisition device is configured to arrange a monitoring system in a rockburst prone area within a predetermined range from a tunnel face to obtain activity data of the rock mass, wherein the monitoring system includes: a displacement sensor, a stress sensor, a temperature sensor, an acoustic wave sensor, and a vibration sensor; A microseismic data acquisition device configured to connect the installed monitoring system to the microseismic monitoring equipment through a data acquisition line, so as to collect the rock mass activity data acquired by the monitoring system to obtain microseismic monitoring data; A microseismic data analysis device is configured to connect the microseismic analysis device and the early warning device to the microseismic monitoring device for communication, establish a time series model through the microseismic analysis device, analyze the microseismic monitoring data, and predict the potential rock burst area in the rock burst prone area near the tunnel construction face based on the microseismic analysis results by the alarm device, and classify the potential rock burst area; The rockburst elimination device is configured to open non-explosion boreholes in the potential rockburst area in a targeted manner according to the early warning signal, and inject low-wave impedance proppants into the non-explosion fracture boreholes to eliminate potential rockburst disasters.

10. The intelligent advanced targeted non-explosive fracturing and rockburst mitigation system according to claim 9, characterized in that: The microseismic data analysis device comprises: A stationary analysis unit is configured to use the historical and current data of the maximum energy of a single day of microseismic in the monitoring area as samples to perform stationary analysis, and if it is not stationary, perform a differential operation to obtain an approximately stationary sequence; Adaptability test unit: configured to select an appropriate model to fit the maximum energy of a single day of sample microseismic events, determine the order and estimate the parameters of the model, and perform adaptability test on the model; Prediction and warning unit: It is configured to use the established model to predict the future value of the difference, and then obtain the future value of the maximum energy of the microseismic day through the inverse difference operation. During the mining process, the established model is used to perform sample extrapolation prediction, and by comparing with the warning indicators, danger warning can be achieved.

Images