Rock burst early warning method for surrounding rock of cavern / tunnel under water inrush disaster
By constructing a fully coupled simulation model of hydrodynamics and a multi-index fusion hierarchical early warning threshold system, the problem of the inability of existing technologies to accurately capture the precursors of surrounding rock disasters under water inrush disasters has been solved, achieving high-precision and rapid disaster early warning and ensuring the safety of deep cavern/tunnel engineering.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- JINAN UNIVERSITY
- Filing Date
- 2026-05-20
- Publication Date
- 2026-06-19
AI Technical Summary
Existing early warning methods cannot accurately capture the catastrophic precursors of surrounding rock under water inrush disasters, and cannot meet the high-precision, quantitative, and efficient requirements of thermal-water-mechanical coupled disasters in deep cavern/tunnel engineering. Traditional methods cannot be adapted to the catastrophic patterns of heterogeneous surrounding rock.
A fully coupled hydrodynamic simulation model was constructed. Through multi-physics field energy interpolation calculation and crack network reconstruction, a multi-index fusion hierarchical early warning threshold system was established to achieve refined characterization and rapid simulation of surrounding rock damage evolution.
It improves the accuracy and reliability of disaster early warning, takes into account computing efficiency, and provides safety assurance for deep underground engineering construction.
Smart Images
Figure SMS_4 
Figure SMS_12 
Figure SMS_24
Abstract
Description
Technical Field
[0001] This invention relates to the field of disaster early warning technology for underground engineering, and in particular to a method for early warning of water inrush disasters caused by rock impact in the surrounding rock of caverns / tunnels. Background Technology
[0002] In the construction of deep caverns / tunnels, the surrounding rock is often in a complex environment characterized by high ground stress, high pore water pressure, and high ground temperature. When a water inrush disaster occurs, the sudden increase in pore pressure leads to a sharp decrease in the effective stress of the surrounding rock. Combined with dynamic impact loads from drilling and blasting, and tunneling, this can easily trigger the rapid initiation and propagation of cracks in the surrounding rock, ultimately leading to impact disasters and seriously threatening the safety of the project. Therefore, achieving accurate early warning of impact disasters in the surrounding rock under water inrush conditions is a core requirement for ensuring the safety of underground engineering construction.
[0003] Existing early warning methods are mostly based on monitoring data from a single mechanical field, without considering the thermo-hydraulic-mechanical coupling effect caused by water inrush. They cannot accurately capture the precursors of disasters under the coupling effect. The fracturing of surrounding rock under water inrush disasters is essentially a strongly coupled process of mechanical deformation, fluid seepage, and heat transfer. Surrounding rock materials have natural heterogeneity, and traditional parameter characterization methods are prone to destroying their geometric symmetry and physical correlation, resulting in early warning models that cannot truly reflect the disaster laws of heterogeneous surrounding rocks. The simulation of water inrush-thermal-mechanical coupling involves strong coupling and non-local interaction of multiple physics fields. Traditional solution methods have low computational efficiency and cannot quickly simulate the disaster process under different working conditions to support the calibration of early warning thresholds.
[0004] In summary, existing early warning methods cannot meet the high-precision, quantitative, and efficient requirements for early warning of sudden thermal-hydraulic-mechanical coupled disasters in deep cavern / tunnel engineering. There is an urgent need to develop an early warning method based on multi-field coupled simulation, adapted to heterogeneous surrounding rock, and capable of threshold quantification. Summary of the Invention
[0005] The purpose of this invention is to overcome the shortcomings of the prior art and provide an early warning method for water inrush disasters in caverns / tunnels, solving the technical defects of existing methods and achieving accurate and timely early warning of water inrush and impact coupled disasters.
[0006] To achieve the above objectives, the present invention adopts the following technical solution: Early warning method for water inrush disaster in caverns / tunnels due to rock impact, including the following steps: Obtain surrounding rock parameters; A fully coupled hydrodynamic simulation model was constructed. The surrounding rock parameters were simulated using the fully coupled hydrodynamic simulation model to obtain discrete damage point clouds, and the discrete damage point clouds were reconstructed into a continuous crack network composed of triangular elements. The weighted quantization method is used to perform multi-physics energy interpolation calculations on the continuous crack network. The multi-physics field includes mechanical field, thermal field, and flow field. The mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element are obtained. Based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element, the energy density weight and scale transition factor are calculated. The crack fractal evolution parameters are then calculated based on the energy density weight and scale transition factor. Based on the energy release rate criterion, the energy density weight and scale transition factor are further solved to obtain the crack propagation parameters; The energy distribution ratio of the multiphysics field is calculated based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element. Based on crack fractal evolution parameters, crack propagation parameters, and energy distribution ratio, a multi-index fusion hierarchical early warning threshold system is constructed and used for early warning.
[0007] Furthermore, a fully coupled hydrodynamic simulation model is constructed. This model is used to simulate the surrounding rock parameters and obtain a discrete damage point cloud. The discrete damage point cloud is then reconstructed into a continuous crack network composed of triangular elements. The specific steps include: A fully coupled simulation model of hydrodynamics was constructed based on conventional peri-field dynamics. A fully coupled hydrodynamic simulation model was used to simulate the surrounding rock parameters and obtain a discrete damage point cloud. The discrete damage point cloud is reconstructed into a continuous crack network using the Delaunay triangulation method. The expression for the continuous crack network is as follows: ,in It is a triangular unit. k This is the index for the triangular unit number; Calculate the geometric parameters of a continuous crack network, which includes two-dimensional and three-dimensional continuous crack networks. The three-dimensional continuous crack mesh is constructed based on the triangular elements of the two-dimensional continuous crack mesh. The geometric parameters include at least the area of the triangular elements of the two-dimensional continuous crack network and the crack volume of the closed surface elements of the three-dimensional continuous crack network.
[0008] Further, the geometric parameters of the continuous crack network are calculated, specifically including the following steps: The area of triangular elements in a two-dimensional continuous crack network is calculated using the vector cross product.
[0009] Area of triangular unit The calculation formula is:
[0010] in, , , These are the three vertices of a triangular unit. , , All are material points; Calculate the crack volume of the closed surface element in a three-dimensional continuous crack network based on the area of the triangular element in the two-dimensional continuous crack network. Crack volume of closed surface element The calculation formula is:
[0011] in, Let be the normal vector of the closed surface element. This represents the centroid displacement of the triangular element. It is a triangular unit.
[0012] Furthermore, the weighted quantization method is used to perform multiphysics energy interpolation calculations on the continuous crack network to obtain the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element. Specifically, this includes the following steps: Calculate any vertex of a triangular unit Mechanical strain energy density at the point Thermal strain energy density Fluid potential energy density ; According to the vertex Calculate the mechanical strain energy density, thermal strain energy density, and fluid potential energy density of the triangular element. Thermal strain energy Fluid potential energy The calculation formulas are as follows:
[0013]
[0014]
[0015] in, h The thickness of the crack surface; As vertices The volume of the spatial region at the initial moment. .
[0016] Furthermore, the energy density weight and scale transition factor are calculated based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element. This process includes the following steps: Based on the mechanical strain energy of the triangular element Thermal strain energy Fluid potential energy Calculate the energy density weight. The calculation formula is:
[0017] in, , , These are the coupling coefficients for the mechanical field, thermal field, and flow field, respectively, used to balance the contributions of each physical field; The feature length of the unit. ; The scale transition factor is calculated based on the energy density weight. The calculation formula is:
[0018] in, The characteristic scale of the crack, The total number of triangular units. The characteristic length of the material.
[0019] Furthermore, the crack fractal evolution parameters are calculated based on the energy density weight and the scale transition factor, specifically including the following steps: The fractal dimension is determined based on energy density weights and scale transition factors. The crack fractal evolution parameters are obtained by differentiating the fractal dimension with respect to time. These parameters represent the crack fractal dimension evolution rate. fractal dimension The calculation formula is:
[0020] in, Energy-driven weights; Crack fractal dimension evolution rate The calculation formula is:
[0021] in, , , All are regression coefficients. For spatial gradient norm, Second-order coupling term Spatial variance.
[0022] Furthermore, based on the energy release rate criterion, the energy density weight and scale transition factor are further solved to obtain the crack propagation parameters, specifically including the following steps: The energy release rate of the triangular unit is calculated based on the energy density weight. The calculation formula is:
[0023] in, For crack propagation increment, For time variables, For time step; The energy release rate and a preset energy release rate threshold are used to determine whether the triangular element participates in the crack propagation process; if the energy release rate is greater than or equal to the preset energy release rate threshold, the triangular element is determined to be an active crack element and participates in the crack propagation process. The crack propagation rate is calculated based on the energy release rate and a preset energy release rate threshold. The formula for calculating the crack propagation rate is as follows:
[0024] in, The preset energy release rate threshold; , Rayleigh wave speed; This is a sensitivity index for crack propagation rate.
[0025] Furthermore, the energy distribution ratio of the multiphysics field is calculated based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element. This process includes the following steps: The sensitivity of the multiphysics field is calculated based on the mechanical strain energy and / or thermal strain energy and / or fluid potential energy of the triangular elements. The sensitivity of each physics field includes thermal field sensitivity and flow field sensitivity. Thermal field sensitivity The expression is:
[0026] Flow field sensitivity The expression is:
[0027] Determine whether the thermal field or the flow field is dominant based on the sensitivity of each physical field. like When the temperature is high, it is determined that the thermal field dominates; otherwise, it is determined that the flow field dominates. Using thermal field sensitivity Flow field sensitivity The second-order coupling term is obtained by calculating the cross derivative, which is then used to analyze the inter-field interactions. The formula for calculation is:
[0028] in, For temperature, Pore pressure, For pore pressure variables; like A positive value indicates synergistic coupling between the thermal field and the flow field, meaning that the fields reinforce each other. like A negative value indicates antagonistic coupling between the thermal field and the flow field; The energy distribution ratio is calculated based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element. According to the energy distribution ratio Perform energy path analysis and energy distribution ratio. The expression is:
[0029] like This indicates that the mechanical field dominates. like This indicates that the thermal field and the flow field are dominant.
[0030] Furthermore, based on crack fractal evolution parameters, crack propagation parameters, and energy distribution ratio, a multi-index fusion hierarchical early warning threshold system is constructed, specifically including the following steps: The preset value ranges for fractal evolution parameters, crack propagation parameters, and energy distribution ratio of the crack are defined. The warning level is set according to whether the values of one or more of the following parameters fall within the corresponding preset parameter range: crack fractal evolution parameters, crack propagation parameters, and energy distribution ratio. The warning levels are at least divided into normal, attention, preliminary warning, Level II warning, and Level I alarm.
[0031] Compared with the prior art, the beneficial effects of the present invention are as follows: By constructing a fully coupled thermal-hydraulic-mechanical simulation model, the multi-field coupled evolution law of surrounding rock impact disaster under water inrush disaster is restored; through crack network reconstruction and accurate interpolation of multi-field energy, the refined characterization of surrounding rock damage evolution is realized; a multi-dimensional integrated hierarchical early warning system is established, which greatly improves the accuracy and reliability of disaster early warning, takes into account both computational efficiency and engineering adaptability, and can provide effective support for the construction safety of deep underground engineering. Detailed Implementation The present invention will be further described in detail below with reference to specific embodiments. The scope of protection of the present invention is not limited to the following embodiments. All equivalent transformations made based on the content of the present invention are within the scope of protection of the present invention.
[0032] It should be noted that the following detailed descriptions are illustrative and intended to provide further explanation of this application. Unless otherwise specified, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this application pertains.
[0033] This invention provides an early warning method for rock impact disasters in lower chambers / tunnels during water inrush, comprising the following steps: Obtain surrounding rock parameters; A fully coupled hydrodynamic simulation model was constructed. The surrounding rock parameters were simulated using the fully coupled hydrodynamic simulation model to obtain discrete damage point clouds, and the discrete damage point clouds were reconstructed into a continuous crack network composed of triangular elements. First, a hydrodynamic fully coupled (THM-PD) simulation model is constructed. The construction process is as follows: A fully coupled simulation model of THM-PD was constructed based on conventional near-field dynamics (OSB-PD).
[0034] The THM-PD simulation model includes a nonlocal integral THM coupled control equation set, which includes mechanical field control equations, seepage field control equations, and thermal field control equations to characterize the coupling relationship between the three fields of heat, water, and force.
[0035] Mechanical field governing equations seepage field governing equations and thermal field governing equations The expressions are as follows:
[0036] in, Mass density; For displacement field; It is a physical field.
[0037]
[0038] For mechanical volumetric strain rate, The rate of temperature change; It is the osmotic pressure; The fluid diffusion coefficient, The mechanical coupling coefficient is... Let be the thermo-coupling coefficient, where Dynamic viscosity; M is the permeability tensor; M is the Biot modulus. For Biot coefficients; It is porosity; The coefficient of thermal expansion of fluid volume; The coefficient of thermal expansion of the surrounding rock skeleton is given. This refers to the water pressure storage capacity of porous media.
[0039]
[0040] in, For mechanical energy consumption, It is the Taylor-Quinni coefficient, which represents the proportion of inelastic work converted into heat; For convective heat transfer energy, It is fluid density. For specific heat capacity, It is the fluid velocity divergence; is the thermal conductivity.
[0041] Secondly, based on the simulation of surrounding rock parameters using a fully coupled hydrodynamic simulation model, discrete damage point clouds are obtained, and these discrete damage point clouds are reconstructed into a continuous crack network composed of triangular elements. Specifically, this mainly includes the following: The THM-PD simulation model is used to simulate the surrounding rock parameters and obtain the geometry of the material points, i.e., the discrete damage point cloud. The expression for the discrete damage point cloud is as follows: ;in, For the discrete damage point cloud, the first i Material point; For three-dimensional Euclidean space, N is the total number of matter points; The discrete damage point cloud is reconstructed into a continuous crack network using Delaunay triangulation. The expression for the continuous crack network is as follows: ,in Triangular elements are used to realize the geometric solidification of the crack network and to calculate the geometric parameters of the continuous crack mesh. The geometric parameters of the continuous crack mesh include the area of the triangular elements of the two-dimensional continuous crack mesh and the crack volume of the three-dimensional continuous crack mesh. The three-dimensional continuous crack mesh is constructed based on the triangular elements of the two-dimensional continuous crack mesh. The geometric parameters of the continuous crack network are calculated as follows: Each triangular unit Composed of three vertices , , constitute, , , All of them are material points.
[0042] Determine the area of the triangular elements in a two-dimensional continuous crack network, where the area of the triangular elements is... The area is calculated using the cross product of vectors. The calculation formula is:
[0043] In a three-dimensional continuous crack network, multiple triangular elements constitute a closed surface element. Closed surface element Crack volume The crack volume Vc can be approximated using the Gaussian divergence theorem, and the formula for calculating it is:
[0044] in, Let be the normal vector of the closed surface element. This represents the centroid displacement of the triangular element. As vertices The displacement vector.
[0045] A weighted quantization method is used to perform multiphysics energy interpolation calculations on a continuous crack network. The multiphysics fields include mechanical, thermal, and fluid fields, yielding the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular elements. Based on these values, energy density weights and scale transition factors are calculated, and crack fractal evolution parameters are then obtained. The specific steps include: The weighted quantization method maps the mechanical strain energy density, thermal strain energy density, and fluid potential energy density at a material point to a triangular element by interpolating the centroid coordinates, thereby obtaining the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element. First, determine the centroid coordinates at the vertices of the triangular unit, expressed as:
[0046]
[0047]
[0048] in, , , Given the three barycentric coordinates of a vertex of a triangular unit, satisfying... and ; Let be any point within the triangular unit. This represents the area of the sub-triangular unit.
[0049] The physical quantities at a point in the material are discrete, while those within the triangular mesh are continuous. For any point within a triangular element... Its field variable value Quantities such as displacement, temperature, and pore pressure need to be obtained by interpolating from the known values at the three vertices using the centroid coordinates.
[0050] Calculate the vertices of the triangular unit mechanical strain energy density Thermal strain energy density Fluid potential energy density and the mechanical strain energy of the triangular unit Thermal strain energy Fluid potential energy Specifically, it includes: For vertex Its mechanical strain energy density With force state vector Related to the deformation state, considering linear elastic materials, vertex mechanical strain energy density Represented as:
[0051] in, for The near field; To be relative to position The force state vector is a variable; To be relative to position The strain state vector is a variable; For the vertex The associated representative volume element; t is the integral variable of the representative volume element; t is the time variable; for The race point, of which the near field region is excluded All points outside the family are family points; The relative positions of vertices and family points. = - ; Represents vertices Influenced by its ethnic points The force; Represents vertices Influenced by its ethnic points The stress; The mechanical strain energy of the triangular element is obtained by weighted averaging of the vertex values:
[0052] in As vertices The associated volume, The vertices are determined using the Voronoi diagram partitioning method. The volume of the spatial region at time t=0, which is the initial configuration. At t=0, the material exhibits a pristine, stress-free, and deformation-free geometric shape and spatial configuration unaffected by any external loads, temperature changes, or fluid pressure. Its value is equal to the Voronoi cell corresponding to that vertex. The volume, that is , In the case of uniform discreteness, It can be taken as the intrinsic volume of the background discretization unit; this volume is used to convert the density of the node-based field variables into the corresponding physical quantity contribution.
[0053] Mechanical strain energy of triangular element Based on the integration of the centroid coordinates, we obtain:
[0054] Using the properties of integrals, we obtain the analytical form:
[0055] Taking the crack surface thickness in two dimensions .
[0056] vertex thermal strain energy density Derived from the constitutive relation of thermal expansion, it can be expressed as:
[0057] The coefficient of thermal expansion; These are thermoelastic parameters; For material points At any moment Temperature increment; coefficients in the above formula These are the conventional coefficients in the elastic strain energy density formula, derived from linear elasticity theory; Thermal strain energy of triangular element Based on the integration of the centroid coordinates, we obtain:
[0058] Using the properties of integrals, we obtain the analytical form:
[0059] vertex Fluid potential energy density Represented as:
[0060] in For Biot coefficient, Pore pressure, For volumetric strain.
[0061] Based on the centroid coordinate integral, the fluid potential energy of the triangular element Similarly, it can be expressed as:
[0062] The crack fractal evolution parameters are calculated based on the energy density weight and the scale transition factor, specifically including the following steps: Based on the mechanical strain energy of the triangular element Thermal strain energy Fluid potential energy Calculate energy density weights Among them, energy density weight Represented as:
[0063] in, , , These are the coupling coefficients for the mechanical field, thermal field, and flow field, used to balance the contributions of each physical field. For cases lacking experimental calibration or where the contributions of the three physical fields are roughly equal under the desired operating conditions, the coupling coefficients can be further simplified to an equal-weighted approach, i.e., directly taken as follows: . The feature length of the unit. ,in, , , It represents the length between any two vertices of a triangular unit.
[0064] Energy density weights are used to quantify the activity of crack elements, with high values corresponding to areas of concentrated energy, such as crack tips.
[0065] Further calculations are made of crack size parameters, i.e., scale transition factors. Define the scale transition factor It is used to bridge micro and macro scales.
[0066] Scale transition factor The expression is:
[0067] in, The characteristic scale of the crack is used as the energy-weighted average crack size. The total number of triangular units. The characteristic length of the material can be taken as the near-field radius. .
[0068] By energy density weight and scale transition factor Determine the fractal dimension Through fractal dimension Achieve topological representation and address fractal dimension Differentiating over time yields the crack fractal evolution parameters, i.e., the crack fractal dimension evolution rate. ; Among them, fractal dimension The expression is:
[0069] in, Energy-driven weights ensure that at the microscale, i.e. To suppress noise, at macroscale... , To strengthen the dominant crack. Dynamically updated. and This allows for tracing crack propagation paths to reveal the evolution of crack networks under THM coupling.
[0070] And for fractal dimension Determining the crack fractal dimension evolution rate by taking the derivative over time. Crack fractal dimension evolution rate The formula for calculation is:
[0071] in , , All are regression coefficients. For spatial gradient norm, Let be the average state variable of the hydrothermal coupled field. Second-order coupling term Spatial variance This is a second-order coupling term between the thermal field and the flow field.
[0072] when The values are all negative, indicating that the spatial gradient norm has an inhibitory effect on the crack evolution rate; The value is significantly positive, reflecting the positive driving effect of the thermal-fluid coupling variance; A stable positive value indicates a shift in the base axis. If the temperature increases... An increase in absolute value indicates that the thermal softening effect strengthens the gradient inhibition; an increase in pressure causes... The decrease indicates that the confining pressure nonlinearly modulates the thermal-fluid coupling. Furthermore, the occurrence of extreme responses under these conditions suggests that the surrounding rock material undergoes microstructural reorganization within this temperature-pressure window, enhancing the coupling mechanism.
[0073] Based on the energy release rate criterion, the energy density weight and scale transition factor are further solved to obtain the crack propagation parameters; specifically, the following steps are included: Define the energy release rate of the triangular unit Through energy release rate Determine whether the triangular elements participate in the crack propagation process.
[0074] The dynamic behavior of crack propagation can be described by the energy release rate criterion. Within the near-field dynamics framework, crack propagation is essentially a bond breaking process, and its critical condition is determined by the energy state. In this embodiment, the triangular unit... Energy release rate Defined as:
[0075] in, The spatial gradient of the energy weights; For crack propagation increment, For time variables, For time steps.
[0076] Determining the energy release rate: using the energy release rate With the preset energy release rate threshold Determining whether a triangular unit participates in the crack propagation process involves an energy release rate threshold. Determined by the fracture toughness of the material; Crack propagation conditions are If the energy release rate Greater than or equal to the preset energy release rate threshold If the condition is met, the triangular element is determined to be an active crack element and participates in the crack propagation process. And based on the energy release rate With the preset energy release rate threshold The crack propagation rate was calculated. Among them, crack propagation rate The expression is:
[0077] in, , Rayleigh wave speed, The index is to be determined and is obtained by performing double log-linear regression on simulated data points.
[0078] By obtaining the crack propagation rate The probability density distribution can be used to determine the strong driving region and thus identify the tip of the main crack.
[0079] Crack propagation rate In addition, crack propagation parameters include crack propagation direction, crack propagation rate, crack path, and crack propagation stability index, which are used for crack propagation analysis.
[0080] Crack propagation direction The crack propagation direction is determined by the spatial gradient of the energy weight. The expression is:
[0081] in The spatial gradient of the energy weights.
[0082] The spatial gradient of the energy weight is defined by the local coordinate system at the crack front, using the normal vector of the triangular element. For reference, the gradient of energy weights for:
[0083] In the formula, As vertices The energy weight gradient at a given location is obtained by fitting the energy weight values of neighboring cells using the least squares method.
[0084] Actual expansion direction Defined by the displacement direction of the centroid of the crack front in adjacent time steps, the expression is:
[0085] Calculate the crack propagation direction With actual expansion direction deviation angle The expression is:
[0086] By analyzing the direction of crack propagation With actual expansion direction deviation angle It can analyze areas where cracks intersect or branch.
[0087] Crack path tracing is achieved by updating the triangular mesh. At each time step, active cells are identified based on the energy release rate condition and tracing is performed along the crack propagation direction. Extend the crack front. Increase the extension increment. The expression is: It is determined by the speed and the time step.
[0088] The generation of new crack surfaces is accomplished through local re-division of the Delaunay triangulation, inserting new points in the propagation direction and re-triangulating to maintain the open-circle characteristic. Crack length The dynamic update is based on energy-weighted accumulation:
[0089] Crack length With scale transition factor This creates a feedback loop to ensure consistency across multiple scales.
[0090] The expression for the crack propagation stability index is:
[0091] in, For energy curvature.
[0092] when When the energy release rate decreases with crack propagation, crack propagation becomes stable; when At this time, the expansion may become unstable, leading to branching or acceleration.
[0093] The energy distribution ratio of the multiphysics field is calculated based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element; specifically including... The sensitivity of field variables in a multiphysics hydrodynamic field is calculated, the inter-field coupling effect is determined based on the sensitivity of the field variables, and the energy distribution ratio of the multiphysics field is calculated. Spatial differentiation and temporal evolution characteristics are determined through inter-field interactions, parameter sensitivity, and energy path analysis. The specific steps include: Based on the mechanical strain energy of the triangular element and / or thermal strain energy and / or fluid potential energy Calculate the sensitivity of the field variables of the hydrodynamic multiphysics field, where the sensitivity of each physics field includes the sensitivity of the thermal field. Flow field sensitivity And based on the sensitivity of each physical field, determine whether the thermal field or the flow field is dominant, where: Thermal sensitivity The expression is:
[0094] in, Directly derived from thermal strain energy. It is necessary to perform thermal stress constitutive calculations based on the thermal components of the PD force state. Force state caused by thermal stress The expression is ; This originates from the effect of temperature on pore pressure and volumetric strain. Scale transition factor. right The derivative reflects the modulation of the characteristic scale of the crack by temperature.
[0095] Flow field sensitivity The expression is:
[0096] in, Defined by fluid potential energy; This is based on the principle of effective stress, which states that pressure affects mechanical strain energy.
[0097] Sensitivity analysis revealed the dominant coupling mechanism: when When the thermal field is dominant, the thermal field dominates; conversely, when the thermal field is dominant, the flow field dominates; the coupling strength is determined by the ratio. Characterization. In addition, sensitivity parameters. and The spatial distributions of these differ significantly. The peak value is strictly limited to the vicinity of the crack tip, which is due to local thermal softening caused by the plastic work-thermal conversion under high strain rate (Taylor-Quinni effect), constituting the core driving force for crack propagation; The high-value areas are mainly distributed in the main channels of the opened cracks and the connected damage zones, reflecting that the fluid pressure maintains the crack opening and promotes the flow of the fluid network by reducing the effective confining pressure and splitting action.
[0098] And adopt thermal field sensitivity Flow field sensitivity The second-order coupling term is obtained by calculating the cross derivative. To analyze the inter-field interactions, the second-order coupling term The expression is:
[0099] In the formula, Pore pressure, For pressure variables; like A positive value indicates cooperative coupling between physical fields, meaning that the fields mutually reinforce each other; if... A negative value indicates antagonistic coupling between physical fields; Based on the mechanical strain energy of the triangular element and / or thermal strain energy and / or fluid potential energy Determine the energy distribution ratio To achieve energy path analysis and energy distribution ratio The expression is:
[0100] Among them, when Indicates that the mechanical field dominates, when This indicates that the thermal-fluid field is dominant. Dynamic tracking. The transition of coupling mechanisms can be identified, such as the transition from thermal dominance to flow dominance under impact loading.
[0101] Based on crack fractal evolution parameters, crack propagation parameters, and energy distribution ratio, a multi-index fusion hierarchical early warning threshold system is constructed and used for early warning. Specifically, this includes the following: The construction process of the multi-indicator fusion hierarchical early warning threshold system includes: Preset parameter levels for crack fractal evolution parameters, crack propagation parameters, and energy distribution ratio are defined separately; crack fractal evolution parameters, crack propagation parameters, and energy distribution ratio serve as early warning indicators.
[0102] The warning level is set based on whether the values of one or more of the following parameters—crack fractal evolution parameters, crack propagation parameters, and energy distribution ratio—fall within the corresponding preset parameter range. The warning levels are at least divided into Normal, Attention, Preliminary Warning, Level II Warning, and Level I Alarm. In this embodiment, the preset value parameter level range of the early warning indicators is divided into four threshold levels: Level I, Level II, Level III, and Level IV. The preset value parameter level range of each early warning indicator in this embodiment is shown in Table 1.
[0103] Table 1. Preset value parameter level ranges for each early warning indicator
[0104] Based on the aforementioned preset parameter range, the specific criteria for determining the warning level are as follows: If any of the warning indicators reaches Level II or above, a preliminary warning is triggered. If two or more of the warning indicators reach Level II at the same time, the warning will be upgraded to Level II. If any one of the warning indicators reaches Level I, or two or more indicators simultaneously reach Level II and include... This triggers a Level I alarm.
[0105] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention in any way. Any simple modifications, equivalent changes and alterations made to the above embodiments based on the technical essence of the present invention shall still fall within the scope of the technical solution of the present invention.
Claims
1. A method for early warning of rock impact disasters in caverns / tunnels during water inrush disasters, characterized in that, Includes the following steps: Obtain surrounding rock parameters; A fully coupled hydrodynamic simulation model is constructed. The surrounding rock parameters are simulated using the fully coupled hydrodynamic simulation model to obtain a discrete damage point cloud. The discrete damage point cloud is then reconstructed into a continuous crack network composed of triangular elements. The continuous crack network is subjected to multi-physics energy interpolation calculation using a weighted quantization method. The multi-physics field includes mechanical field, thermal field, and flow field. The mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular unit are obtained. Based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular unit, the energy density weight and scale transition factor are calculated. The crack fractal evolution parameters are then calculated based on the energy density weight and the scale transition factor. Based on the energy release rate criterion, the energy density weight and the scale transition factor are further solved to obtain the crack propagation parameters; The energy distribution ratio of the multiphysics field is calculated based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular unit. Based on the crack fractal evolution parameters, the crack propagation parameters, and the energy distribution ratio, a multi-index fusion hierarchical early warning threshold system is constructed, and early warning is issued.
2. The method for early warning of rock impact disasters in caverns / tunnels caused by sudden water inrush as described in claim 1, characterized in that, A fully coupled hydrodynamic simulation model is constructed. This model is used to simulate the surrounding rock parameters to obtain a discrete damage point cloud. The discrete damage point cloud is then reconstructed into a continuous crack network composed of triangular elements. The specific steps include: A fully coupled simulation model of hydrodynamics was constructed based on conventional peri-field dynamics. The surrounding rock parameters are simulated using the aforementioned fully coupled hydrodynamic simulation model to obtain a discrete damage point cloud; The discrete damage point cloud is reconstructed into a continuous crack network using the Delaunay triangulation method. The expression for the continuous crack network is as follows: ,in It is a triangular unit. k This is the index for the triangular unit number; Calculate the geometric parameters of the continuous crack network, which includes a two-dimensional continuous crack network and a three-dimensional continuous crack network. The three-dimensional continuous crack mesh is constructed based on the triangular elements of the two-dimensional continuous crack mesh. The geometric parameters include at least the area of the triangular elements of the two-dimensional continuous crack network and the crack volume of the closed surface elements of the three-dimensional continuous crack network.
3. The method for early warning of rock impact disasters in caverns / tunnels caused by sudden water inrush as described in claim 2, characterized in that, The calculation of the geometric parameters of the continuous crack network specifically includes the following steps: The area of the triangular elements in the two-dimensional continuous crack network is calculated using the vector cross product. The area of the triangular unit The calculation formula is: in, , , These are the three vertices of a triangular unit. , , All are material points; The crack volume of the closed surface element of the three-dimensional continuous crack network is calculated based on the area of the triangular element in the two-dimensional continuous crack network. Crack volume of the closed curved surface unit The calculation formula is: in, Let be the normal vector of the closed surface element. This represents the centroid displacement of the triangular element. It is a triangular unit.
4. The method for early warning of rock impact disasters in caverns / tunnels caused by sudden water inrush as described in claim 3, characterized in that, The energy interpolation of the continuous crack network using the weighted quantization method is performed to obtain the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element. Specifically, the following steps are included: Calculate any vertex of a triangular unit Mechanical strain energy density at the point Thermal strain energy density Fluid potential energy density ; According to the vertex Calculate the mechanical strain energy density, thermal strain energy density, and fluid potential energy density of the triangular element. Thermal strain energy Fluid potential energy The calculation formulas are as follows: in h The thickness of the crack surface; As vertices The volume of the spatial region at the initial moment. .
5. The method for early warning of rock impact disasters in caverns / tunnels caused by sudden water inrush as described in claim 4, characterized in that, The energy density weight and scale transition factor are calculated based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular unit, specifically including the following steps: Based on the mechanical strain energy of the triangular element Thermal strain energy Fluid potential energy Calculate the energy density weight, the energy density weight The calculation formula is: in, , , These are the coupling coefficients for the mechanical field, thermal field, and flow field, respectively, used to balance the contributions of each physical field; The feature length of the unit. ; The scale transition factor is calculated based on the energy density weight. The calculation formula is: in, The characteristic scale of the crack, The total number of triangular units. The characteristic length of the material.
6. The method for early warning of rock impact disasters in caverns / tunnels under water inrush disasters according to claim 5, characterized in that, The crack fractal evolution parameters are calculated based on the energy density weight and the scale transition factor, specifically including the following steps: The fractal dimension is determined based on the energy density weight and the scale transition factor. The crack fractal evolution parameter is obtained by differentiating the fractal dimension with respect to time; this parameter is the crack fractal dimension evolution rate. The fractal dimension The calculation formula is: in, Energy-driven weights; The crack fractal dimension evolution rate The calculation formula is: in, , , All are regression coefficients. For spatial gradient norm, Second-order coupling term Spatial variance.
7. The method for early warning of rock impact disasters in caverns / tunnels under the category of water inrush disasters according to claim 6, characterized in that, Based on the energy release rate criterion, the energy density weight and the scale transition factor are further solved to obtain the crack propagation parameters, specifically including the following steps: The energy release rate of the triangular unit is calculated based on the energy density weight. The calculation formula is: in, For crack propagation increment, For time variables, For time step; The energy release rate and a preset energy release rate threshold are used to determine whether the triangular unit participates in the crack propagation process. If the energy release rate is greater than or equal to the preset energy release rate threshold, then the triangular unit is determined to be an active crack unit and participates in the crack propagation process; The crack propagation rate is calculated based on the energy release rate and the preset energy release rate threshold. The calculation formula is: in, The preset energy release rate threshold; , Rayleigh wave speed; This is a sensitivity index for crack propagation rate.
8. The method for early warning of water inrush disaster in caverns / tunnels according to claim 7, characterized in that, The energy distribution ratio of the multiphysics field is calculated based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular unit, specifically including the following steps: The sensitivity of the multiphysics field is calculated based on the mechanical strain energy and / or thermal strain energy and / or fluid potential energy of the triangular elements, wherein the sensitivity of each physics field includes the thermal field sensitivity. Flow field sensitivity The thermal field sensitivity The expression is: The flow field sensitivity The expression is: Determine whether the thermal field or the flow field is dominant based on the sensitivity of each physical field. like When the temperature is high, it is determined that the thermal field dominates; otherwise, it is determined that the flow field dominates. Using thermal field sensitivity Flow field sensitivity The second-order coupling term is obtained by calculating the cross derivative to analyze the inter-field interaction. The formula for calculation is: in, For temperature, Pore pressure, For pore pressure variables; like A positive value indicates synergistic coupling between the thermal field and the flow field, meaning that the fields reinforce each other. like A negative value indicates antagonistic coupling between the thermal field and the flow field; The energy distribution ratio is calculated based on the mechanical strain energy, thermal strain energy, and fluid potential energy of the triangular element. According to the energy distribution ratio Perform energy path analysis and energy distribution ratio. The expression is: like This indicates that the mechanical field dominates; if This indicates that the thermal field and the flow field are dominant.
9. The method for early warning of rock impact disasters in caverns / tunnels under the category of water inrush disasters according to claim 1, characterized in that, Based on the crack fractal evolution parameters, the crack propagation parameters, and the energy distribution ratio, a multi-index fusion hierarchical early warning threshold system is constructed, specifically including the following steps: The preset value parameter levels of the crack fractal evolution parameters, the crack propagation parameters, and the energy distribution ratio are defined. The warning level is set according to whether the values of one or more of the crack fractal evolution parameters, crack propagation parameters and energy distribution ratio fall within the corresponding preset parameter range. The warning level is at least divided into normal, attention, preliminary warning, level II warning and level I alarm.