Underground environment data acquisition method and system for three-dimensional simulation environment modeling

Abstract

The invention relates to the technical field of geological exploration, in particular to an underground environment data acquisition method and system for three-dimensional simulation environment modeling, working condition parameters are acquired and input into a CFD model, and an underground dynamic risk map is generated based on shock wave signals; the method comprises the following steps: collecting water mist data to construct a humidity field, and embedding humidity field data into a dynamic risk map generation process, which comprises the following steps: in an injection compression period, solving an N-S equation set, detecting and encrypting a key area grid through a pressure gradient, and generating a burst dust cloud form and density distribution. The results of four stages of turbulence evolution are fused to form a five-dimensional tensor structure including space coordinates, duration time, humidity compensation risk levels, a water film distribution thermodynamic diagram and a water mist light attenuation coefficient, different risk areas can be clearly marked, switching of data acquisition strategies is guided, the problems of artifacts and distortion are effectively solved, and the accuracy of data acquisition is improved. Accurate point cloud data can be obtained in a complex environment.

Description

Technical Field

[0001] This invention relates to the field of geological exploration technology, specifically to a method and system for acquiring downhole environmental data for three-dimensional simulation environment modeling. Background Technology

[0002] With the development of computer graphics and virtual reality technology, 3D simulation environment modeling has been widely applied in many fields. Especially in the mining industry, 3D simulation environment modeling can achieve a virtual reproduction of the underground environment, which is of great significance for improving mining safety, optimizing mining plans, and training miners. Existing 3D simulation environment modeling technologies generally involve multi-source sensor fusion data acquisition, such as constructing a spatial point cloud model with centimeter-level accuracy using lidar and inertial navigation systems, and then utilizing distributed sensor networks to achieve dynamic monitoring of environmental parameters.

[0003] However, in coal mines, the high concentration of dust generated during coal mining operations forms a non-uniform suspended particle layer, causing multipath scattering in optical equipment. In particular, the collision between the cutting teeth of the coal mining machine and the coal seam generates an explosive dust cloud, forming a strong instantaneous attenuation barrier for the laser signal. This causes the laser beam to undergo total internal reflection or total absorption critical effects during penetration, resulting in pulsed tortuosity in the generated point cloud data.

[0004] Furthermore, the concentration gradient changes during the settling process of the dust cloud can lead to deceptive data acquisition by optical equipment. The high-pressure shock wave generated by the cutting teeth colliding with the coal seam instantly compresses the air, ejecting coal dust at high speed to form a cone-shaped dust cloud. As the dust cloud flows along the wall of the coal mining machine, it generates low-pressure vortices, forming the initial turbulent core. The turbulence evolves in four stages: the jet compression phase (dust cloud expansion compresses air), shear layer formation (a velocity difference arises between the dust cloud boundary and still air, Reynolds stress dominates turbulence development), vortex ring instability (Köhler instability causes vortex ring rupture), and full vortex development (forming multi-scale vortices, conforming to the turbulent energy spectrum law). Strong turbulence regions, through vortex entrainment, locally form turbulent void regions, i.e., low-concentration areas. When the laser beam passes through these turbulent void regions, a beam-like projection phenomenon occurs, leading to distortion of the spatial topology, including the generation of pseudo-geometric features. For example, a beam-like transmission channel may be misidentified as a physical channel, forming a non-existent tunnel-like pseudo-structure in the point cloud. Furthermore, the temporal continuity is disrupted, and the topological discontinuity between adjacent frames is caused by the instantaneous appearance of transparent windows, resulting in a discrepancy between the 3D simulation environment model and the actual physical environment. Summary of the Invention

[0005] To address the aforementioned shortcomings of existing technologies, this invention provides a method and system for acquiring downhole environmental data for three-dimensional simulation environment modeling. This method effectively solves the problems of point cloud distortion, artifact misjudgment, and topological fracture caused by explosive dust clouds and turbulence in existing technologies.

[0006] To achieve the above objectives, the present invention provides the following technical solution: This invention provides a downhole environment data acquisition method for three-dimensional simulation environment modeling. The technical solution adopted is as follows: acquiring working condition parameters and inputting them into a CFD model, and generating a downhole dynamic risk map based on shock wave signals; Collect water mist data to construct a humidity field, and embed the humidity field data into the generation process of the downhole dynamic risk map, including the injection compression stage, shear layer generation stage, vortex ring instability stage, and vortex full development stage; Finally, the simulation results of the above four stages are integrated into a reconstructed heat map, and the acquisition strategy of downhole environmental point cloud data is switched based on the risk area. By analyzing laser echo signals and detecting abnormal point cloud data in the downhole environment through distortion analysis, the vortex equations are reconstructed using CFD parameters to optimize the acquisition accuracy of downhole environmental point cloud data.

[0007] Furthermore, during the injection compression period, the NS equations are solved, the grid in the key area is refined by pressure gradient detection, the morphology and density distribution of the explosive dust cloud are generated, and the dust and steam mixing density field is output based on the volume expansion effect generated by the water vaporization phase change. When the shear layer is generated, the velocity gradient tensor of the boundary layer of the explosive dust cloud is calculated, and the vertical component of the velocity gradient tensor is corrected by the wet dust settling acceleration term to generate seed points for the turbulent development core containing water mist settling effect. When the vortex ring becomes unstable, the stability model of the vortex ring is reconstructed under the action of water mist viscous damping, the unstable region of explosive dust cloud flow is identified, and the location of vortex ring rupture is predicted. The inhibitory effect of water film on vortex ring rupture is quantified by adaptive vorticity threshold, and the turbulence morphology transformation feature map under humidity environment is output. When the eddies are fully developed, a water mist dissipation function is embedded in the K-ω turbulence model to simulate the evolution of multi-scale eddies. The real turbulent cavity region is identified by a dual threshold criterion of turbulent kinetic energy and humidity. Optical artifact regions and steam artifact regions are labeled simultaneously, and a distribution map of turbulent cavity regions with risk labels is output.

[0008] Furthermore, the method for encrypting the key area mesh is as follows: A virtual space for simulating explosive dust clouds is created. A Cartesian grid is created based on the geometry of the cutting teeth. The grid range of the virtual space is determined with the cutting teeth as the center, and the virtual space is divided into several grid units according to the preset grid size. Define the boundary and initial state of the explosive dust cloud flow, traverse all grid cells in the virtual space, calculate the pressure change rate and pressure gradient amplitude in three directions for each grid cell, compare the pressure gradient amplitude of the grid cell with the preset pressure gradient critical value, and divide the grid cell into normal grid and risk grid. The risk grid is decomposed into risk subgrids, and the data of the risk subgrids are integrated into the overall computation model according to the timestamps. The physical parameters of the affected grid cells are then iteratively updated based on the timestamps.

[0009] Furthermore, the method of decomposing into risk sub-grids is as follows: The risk grid is decomposed into several small cubes, denoted as risk sub-grids. Each risk sub-grid inherits the physical properties of its parent grid. The boundary data of the risk sub-grids is then smoothed using the following method: Locate and mark the junction between the risk subgrid and the normal grid, and symmetrically extend the transition zone on both sides of the junction; Key data points of risk subgrids and normal grids are collected. Using the width of the transition zone as a baseline, the asymptotic transition value between the normal grid and risk subgrids is calculated based on the bilinear interpolation method. Forced physical constraints are set for the transition zone, including mass conservation, momentum continuity and energy stability. Repeat the above steps for iterative optimization until the numerical mutation rate of all sampling points is lower than the preset threshold, then stop the optimization.

[0010] Furthermore, the method for constructing the humidity field is as follows: A global coordinate system with the cutter tooth as the origin is established in the virtual space; the pose transformation matrix is ​​obtained in real time through the UWB module, a grid index mapping table for different sensors is established, the sensor coordinates are mapped to the global coordinates, and then the global coordinates are used to project them to the nearest neighboring grid cell in the virtual space, and cross-cell assignment is performed on the boundary points. The water mist data projected onto the same cell are aligned according to the timestamp of the sampling time, and weights are assigned according to the sensor accuracy. The water mist data projected from multiple cells within the grid cell are weighted and fused to calculate the humidity gradient field within the grid area in the virtual space. Then, cross-grid correlation analysis is performed on the humidity gradient field to generate a water film distribution heatmap affected by water film diffusion; the humidity compensation risk level of each grid cell is marked in the water film distribution heatmap.

[0011] Furthermore, the method for performing cross-cell allocation of boundary points is as follows: Calculate the distance between the sensor projection point and the center points of all adjacent grid cells, denoted as . ; The minimum of all its distances Compared with the boundary influence threshold, the minimum value The point is determined as a boundary point; where P is the spatial coordinate of the sensor projection point; The spatial coordinates of the center point of the i-th adjacent grid cell; The size of the grid cell; Boundary influence factor; Select all mesh cells that satisfy the boundary point conditions and merge them into a candidate set; Calculate the assignment weight of the projection point to each adjacent grid cell in the candidate set, and then normalize all assignment weights to obtain the normalized assignment weights. The water mist data at sensor projection point P is allocated to each adjacent grid cell according to the normalized assignment weight; the component data of each boundary point within the same grid cell are accumulated, and boundary effects are eliminated according to intelligent gradient optimization.

[0012] Furthermore, the method for performing cross-grid correlation analysis on the humidity gradient field is as follows: A directed connectivity graph with 26 neighborhood indices is constructed for each grid cell, and dynamic association weights are calculated based on humidity gradient difference and wind flow vector angle; thus forming a weighted topology network that reflects spatial relationships and wind flow vectors. Humidity mass conservation is maintained based on gradient continuity constraint equations and flux equations across grid cells; Constructing the spatiotemporal correlation matrix Record the time delay of grid cell i. The threshold of the influence on grid cell j is used to predict the future state based on the current gradient field, and the actual measured value is used to correct the dynamic association weight. Construct specific propagation models for different types of risks within grid cells, and visualize the risk transmission paths and impact ranges.

[0013] Furthermore, the reconstructed heatmap adopts a five-dimensional time-space tensor structure, namely: in, To reconstruct the heatmap; The coordinates are those of the actual turbulent cavity region; The duration of the actual turbulent cavity region; Humidity compensation risk level; Thermographic diagram of water film distribution; Water mist light attenuation coefficient The criteria for identifying optical artifacts are: the water film thickness is greater than the critical value of the water film thickness and the light attenuation coefficient of the water mist is greater than the light attenuation threshold. The criteria for identifying steam artifact regions are: humidity H is greater than the humidity saturation threshold and turbulent kinetic energy is less than the critical value for turbulence disappearance.

[0014] Furthermore, the risk area switching-based data collection strategy includes: In the real turbulent cavity region, the laser sampling frequency is increased over a period of time to establish a spatiotemporal interpolation model to fill the topological fracture region, and a dynamic confidence factor C is used to perform weighted fusion of the filled point cloud. In the optical artifact region, the optical path offset is calculated based on the water film thickness and refractive index to perform beam path correction; the laser power is enhanced based on the optical attenuation coefficient; when the water film thickness exceeds the critical value, a moving median filter is used to eliminate signal flicker noise caused by water film reflection. In the vapor artifact region, the laser is switched to the 1570nm infrared band and pulse code modulation is used.

[0015] A downhole environment data acquisition system for three-dimensional simulation environment modeling includes a data acquisition module for shock wave signals, operating parameters, and water mist data; The CFD simulation module performs four-stage turbulence evolution calculations based on operating parameters and water mist data. The dynamic risk assessment module generates a five-dimensional reconstructed heatmap with risk labels based on the calculation results of the CFD simulation module. The acquisition and control module is based on the region type switching acquisition strategy in the five-dimensional reconstructed heat map; The data update module is used to analyze the waveform characteristics of the laser echo signal, detect anomalies, and optimize the turbulence model by combining CFD parameters.

[0016] In summary, due to the adoption of the above technical solution, the beneficial effects of the present invention are: 1. This invention embeds humidity field data into CFD simulations, considering the influence of water mist in each of the four stages of turbulent evolution, including volume expansion, wet dust settling, viscous damping, and water mist dissipation. It achieves end-to-end optimization of explosive dust clouds from diffusion source to settling endpoint. Through dynamic coupling of gas-liquid-solid three-phase flow, it improves the prediction of dust settling trajectory, solves the interference caused by water mist dust settling (such as optical transparency, misjudgment caused by vapor filling, and changes in trajectory due to wet dust settling), improves the accuracy of turbulent cavity region identification, and provides a reliable basis for subsequent data acquisition.

[0017] 2. This invention integrates the results of four stages of turbulent evolution to form a five-dimensional tensor structure that includes spatial coordinates, duration, humidity compensation risk level, water film distribution thermogram, and water mist light attenuation coefficient. This structure can clearly mark different risk areas (real void areas and artifact areas), guide the switching of data acquisition strategies, effectively overcome artifact and distortion problems, and ensure accurate point cloud data in complex environments.

[0018] 3. When constructing the humidity field, this invention uses 26-neighborhood index, dynamic weight calculation, and spatiotemporal correlation matrix to achieve accurate simulation of the humidity field and prediction of risk propagation. By considering factors such as airflow, humidity gradient, and time delay, an accurate humidity field model is constructed to predict water film diffusion and risk propagation paths, providing a basis for humidity compensation risk levels and enhancing the system's adaptability to dynamic environments. Attached Figure Description

[0019] Figure 1 This is a simplified flowchart illustrating the downhole environmental data acquisition method in Embodiment 1 of the present invention. Detailed Implementation

[0020] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention.

[0021] Example 1 Reference Figure 1 This embodiment proposes an underground environmental data acquisition method for 3D simulation environment modeling. This method can acquire accurate 3D point cloud data in environments with high concentrations of dust generated during coal mining machine operations, providing reliable data input for subsequent 3D simulation environment modeling in coal mines. It avoids problems such as multipath scattering, total reflection, or total absorption critical effects caused by dust clouds, as well as artifacts caused by voids due to turbulence, which distort the acquired point cloud data and thus affect the accuracy of the 3D simulation environment model.

[0022] As a further preferred implementation of this embodiment, the downhole environment data acquisition method includes: A sensor array is installed at the cutting edge of the coal mining machine to monitor the shock wave signal generated at the moment of collision between the cutting edge and the coal seam in real time, thereby determining whether turbulence has occurred. Based on dust generation dynamics, the shock wave pressure and dust ejection rate of the cutting edge colliding with the coal seam are exponentially related. When the shock wave pressure exceeds a certain threshold, it indicates the formation of an explosive dust cloud. Furthermore, the evolution of turbulence is predicted. Because the time required from collision to the evolution of an explosive dust cloud is extremely short, a simulation strategy of triggering the fluid dynamics CFD model in advance through the shock wave signal is adopted. This provides the acquisition system with additional response time, ensuring that the system can react more quickly and thus improving its ability to cope with the impact of potential explosive dust clouds.

[0023] Simultaneously, operating parameters such as the rotational speed of the cutting teeth, feed rate, and coal seam hardness are acquired. By inputting the measured operating parameters into the CFD model as boundary conditions, simulation strategy errors caused by purely theoretical assumptions can be avoided, thereby improving the accuracy of turbulence risk level assessment.

[0024] The CFD model generates a dynamic risk map of underground coal mines based on operating parameters, marking high-concentration dust areas and turbulent cavitation areas, and indicating the risk level of each area. The dynamic risk map is generated as follows: during the four-stage evolution of turbulence: During the jet compression phase, the compressible Navier-Stokes equations are solved based on the large eddy model to dynamically simulate the initial expansion process of the explosive dust cloud; the grid in key areas is automatically refined by pressure gradient detection, and the dust cloud morphology and density distribution are output.

[0025] During the generation of the shear layer, the velocity gradient tensor of the boundary layer of the explosive dust cloud is calculated, and seed points for the core of turbulence development are generated based on the Reynolds stress model, providing initial conditions for subsequent turbulence evolution.

[0026] When the vortex ring becomes unstable, the Richardson number is used to identify the unstable region of the explosive dust cloud flow, and the Kelvin-Helmholtz instability model is applied to predict the location of the vortex ring rupture, thus capturing the key nodes of turbulent morphology transformation.

[0027] When the eddies are fully developed, the K-ω turbulence model is used to simulate the evolution of multi-scale eddies. The turbulence cavitation region is automatically identified and marked by the turbulence kinetic energy threshold, and the prediction of the entire life cycle of turbulence is completed.

[0028] Finally, the simulation results of the above four stages are integrated into a spatiotemporal four-dimensional heat map, i.e., a dynamic risk map; the three-dimensional coordinates and duration of the turbulent cavity region in the dynamic risk map are marked, and the risk level of the pseudo-structure is calculated based on the dust concentration gradient and duration.

[0029] It is worth noting that water spraying is commonly used to suppress dust during coal mining operations, a fundamental and crucial measure for safe coal mine production. While water spraying significantly reduces dust concentration, it also significantly alters the physical properties and movement patterns of the dust, leading to the misidentification of pseudo-turbulent cavitation zones. Specifically: water mist forms a protective water film on the dust surface, causing localized areas to exhibit optical transparency (high actual dust concentration, but encased in the water film); simultaneously, steam fills the actual turbulent cavitation zones, maintaining turbulence but reducing light attenuation, causing low-concentration areas to be incorrectly labeled as high-concentration areas; furthermore, the accelerated settling of moist dust alters the boundary layer evolution trajectory, interfering with the spatiotemporal location of cavitation zones. To address these issues, the method for correcting and constructing dynamic risk maps based on water mist interference is as follows: When acquiring operating parameters, water mist data is collected to construct a humidity field, and humidity field data is embedded in the four-stage evolution of turbulence. During the jet compression phase, water mist data is incorporated into the initial expansion stage of the explosive dust cloud. A steam expansion equation is introduced to calculate the volumetric expansion effect caused by the water vaporization phase change, simulating the expansion process of the steam-dust mixture in real time. The calculation formula is as follows: In the formula, This is the steam expansion source term, reflecting the mass of steam generated per unit volume per unit time, and characterizing the water mist vaporization intensity; The phase transition coefficient reflects the water mist vaporization efficiency. Water mist mass flow rate represents the mass of water mist passing through a unit area per unit time; The temperature difference between dust and water mist is represented. An improved dust-vapor mixing density field is output to capture the initial turbulence generation characteristics under humid conditions; the output formula is: In the formula, The density of the mixture of dust and steam; The velocity field of the mixed medium; During shear layer formation, a term for the settling acceleration of moist dust is added to quantify the effect of water mist-induced coal powder agglomeration and weight gain on the settling velocity; Moist dust settling acceleration In the formula, g is the acceleration due to gravity. The particle density of the dry dust; The density of the water mist; The equivalent particle size of the wetted dust represents the increase in equivalent diameter after dust agglomeration. By correcting the vertical component of the velocity gradient tensor, a corrected vertical component is obtained that accurately describes the evolution trajectory of the wetted dust boundary layer, outputting the core parameters of turbulent development with water mist settling effects. Corrected vertical component In the formula, This is the vertical component of the original velocity gradient tensor; The coefficient of adhesion of the water film; For time step.

[0030] When the vortex ring becomes unstable, the stability model of the vortex ring is reconstructed under the viscous damping effect of water mist, and the Richardson number after humidity correction is used. Recalibrate the flow instability region; where, represents the water film damping coefficient, reflecting the water film's efficiency in suppressing turbulence; RI is the original Charison number. Using an adaptive vorticity thresholding technique, the inhibitory effect of the water film on the vortex ring rupture process is quantified, and a turbulence morphology transformation characteristic map under humid conditions is output.

[0031] When the eddies are fully developed, a three-phase flow turbulence model is established, and the water mist dissipation function is embedded in the turbulence kinetic energy equation. The expression of the three-phase flow turbulence model is: In the formula, K is the turbulent kinetic energy; w is the specific dissipation rate, which reflects the rate of turbulent kinetic energy dissipation, i.e., the vortex breaking efficiency. This is the turbulent kinetic energy generation term; These are empirical constants, i.e., the coefficients of the standard K-ω model; The water mist dissipation function is used. A dual threshold criterion of turbulent kinetic energy and humidity is used to rigorously identify the real turbulent cavity region, and simultaneously determine the optical artifact region (water film influence region) and vapor artifact region (high humidity and low pressure region), and output a turbulent cavity region distribution map with risk labels.

[0032] Finally, the simulation results from the above four stages are merged into a reconstructed heatmap, using a five-dimensional time-space tensor structure: in, The coordinates are those of the actual turbulent cavity region; The duration of the actual turbulent cavity region; Humidity compensation risk level; Thermographic diagram of water film distribution; is the light attenuation coefficient of water mist.

[0033] The criteria for identifying optical artifact regions are as follows: and In the formula, This is the critical value for water film thickness. This is the light attenuation threshold.

[0034] The criteria for identifying steam artifact regions are: and In the formula, The humidity saturation threshold. This is the critical value for the disappearance of turbulence.

[0035] A data acquisition strategy based on switching between real turbulent cavitation regions, optical artifact regions, and steam artifact regions in the reconstructed heatmap; specifically: For real turbulent cavity regions, the laser sampling frequency is increased over a duration to establish a spatiotemporal interpolation model to fill the topological break regions. A dynamic confidence factor C is then used to perform weighted fusion of the filled point cloud to ensure the geometric topological integrity of the cavity region.

[0036] For optical artifact areas, the optical path offset is calculated based on the water film thickness and refractive index to perform beam path correction; secondly, the laser power is enhanced based on the optical attenuation coefficient to overcome signal attenuation; when the water film thickness exceeds the critical value, multi-angle data is acquired through triple redundancy acquisition, and moving median filtering is used to eliminate signal flicker noise caused by water film reflection.

[0037] For the vapor artifact region, when the humidity H is greater than or equal to the humidity saturation threshold, the system automatically switches to the 1570nm infrared band laser and uses pulse code modulation (PCM) technology to enhance signal penetration. At the same time, a phase change compensation model is applied: the physical parameters are determined by a vapor density and temperature lookup table, and a vapor compensation factor is introduced to eliminate the ranging error caused by vapor phase change, ensuring the distance measurement accuracy in high temperature and high humidity environments.

[0038] When multiple regions are superimposed, they are processed in the following priority order: steam artifact region (which causes complete attenuation of laser signal and irreversible data loss), optical artifact region (which causes spatial coordinate shift, which can be corrected but affects accuracy), and turbulent cavity region (which only has local data loss and can be recovered by interpolation).

[0039] Finally, based on quantum noise suppression technology, the waveform characteristics of the laser echo signal are analyzed in real time, and the distortion of the laser echo signal is calculated. When the distortion of the laser echo signal is greater than the distortion threshold, an anomaly marker is triggered. The spatiotemporal coordinates of the point cloud data of the anomaly marker are collected, and the CFD turbulence field parameters of the corresponding region are extracted to generate a sample set with environmental labels. The sample set is injected into the CFD model to reconstruct the vortex equation, and the compensation parameters are adaptively adjusted by the gradient descent method to minimize the prediction error.

[0040] It should be noted that using the LES model to solve the Navier-Stokes equations can accurately describe the expansion process of dust clouds under high-pressure impact, providing a physical basis for subsequent calculations. Its expression is: In the formula, The mass conservation equation states that the rate of change of mass of a dust cloud per unit volume plus the net outflow mass caused by the flow of the dust cloud equals the mass of the dust cloud generated at that point; where... The density of the dust cloud; For velocity vector field, Where u is the velocity component in the x-direction, v is the velocity component in the y-direction, and w is the velocity component in the z-direction; The density-time change rate of a dust cloud reflects the increasing or decreasing trend of the dust cloud density over time. As a velocity vector, it reflects the flow velocity of the dust cloud, including direction; Mass flux divergence represents the net outflow of dust cloud mass per unit volume per unit time. This is the divergence operator, representing the degree of divergence of a vector field; This is a mass source term, representing the mass of dust clouds generated per unit volume per unit time.

[0041] The momentum conservation equation states that the rate of change of momentum per unit volume plus the net outflow rate of momentum due to fluid motion equals the net force acting on that point, consisting of the pressure gradient and viscous force; where... The momentum density of a dust cloud reflects the magnitude of the momentum of the dust cloud's motion. This is the tensor product operator; The momentum convection term describes the migration of momentum as the dust cloud moves; This is the pressure gradient force, the driving force from the high-pressure area to the low-pressure area; For viscous stress tensor, it represents the stress generated due to the viscous effect of dust cloud.

[0042] The method for automatically refining the mesh in key areas using pressure gradient detection is as follows: A virtual space for simulating explosive dust clouds is created. A Cartesian mesh is generated based on the geometry of the cutting edge, and the mesh range of the virtual space is determined with the cutting edge as the center. The virtual space is divided into several mesh cells according to a preset mesh size (usually 0.5-5 mm). The boundaries and initial states of the explosive dust cloud flow are defined; the boundaries include the cutting edge surface (where dust and solid surfaces are in relative motion), the coal seam interface (the source of dust generation), and the outer boundary (where dust can flow out freely). The mathematical expression for the initial state of the explosive dust cloud flow is: In the formula, The position vector represents the three-dimensional coordinates (x, y, z) of any position in the virtual space, that is, the grid cell at the three-dimensional coordinates (x, y, z); Location in virtual space Dust concentration at the initial moment; The three-dimensional coordinates of the center point of the explosive dust cloud eruption; This represents the peak concentration at the dust source. The diffusion radius of the explosive dust cloud is determined based on experimental data.

[0043] The algorithm iterates through all grid cells in the virtual space, calculating the pressure change rate and pressure gradient magnitude in three directions for each cell. It compares the pressure gradient magnitude of each cell with a preset pressure gradient threshold. Cells with pressure gradient magnitudes less than the threshold are marked as normal cells, while those with magnitudes greater than or equal to the threshold are marked as risk cells. Risk cells are then decomposed into several small cubes, denoted as risk sub-cells. These sub-cells inherit the physical properties of their parent cell (the original risk cell), including dust concentration (the average value of the parent cell) and flow velocity direction (maintaining consistency with the parent cell). Next, the boundary data of the risk subgrid is smoothed. The smoothing method is as follows: locate and mark the junction between the risk subgrid and the normal grid, and symmetrically extend the transition zone on both sides of the junction; collect key data points (such as corner data) of the risk subgrid and the normal grid, use the width of the transition zone as a baseline, calculate the asymptotic transition value between the normal grid and the risk subgrid based on bilinear interpolation, and set mandatory physical constraints for the transition zone, including mass conservation (ensuring that the total amount of dust in the transition zone does not change abruptly), momentum continuity (eliminating abrupt changes in airflow direction), and energy stability (preventing the generation of false hot spots / cold spots); repeat the above steps for iterative optimization until the numerical mutation rate of all sampling points is lower than the preset threshold, at which point the optimization stops. Then, the physical field of the virtual space is updated, the data of the risk subgrid is integrated into the overall calculation model according to the timestamp, and the physical parameters of the affected grid cells are iteratively updated based on the timestamp.

[0044] Finally, the transient explosive dust cloud morphology and density distribution are output. t represents the time scale.

[0045] It should be noted that the formula for the velocity gradient tensor is: In the formula, It represents the rate of change of the velocity component u in the x-direction along the x-axis (0 for lateral stretching or compression, which determines the expansion rate of the explosive dust cloud). This represents the rate of change of the y-direction velocity component u along the y-axis (lateral shear), which drives the generation of shear layer vortices. It represents the rate of change of the z-direction velocity component u along the z-axis (vertical mixing), which affects the upward diffusion of explosive dust clouds; It represents the rate of change of the velocity component v in the x-direction along the x-axis (lateral momentum transfer), which controls the energy exchange at the boundary of explosive dust clouds; It represents the rate of change of the y-direction velocity component v along the y-axis (normal deformation), which dominates the anisotropy of turbulence; This represents the rate of change of the z-direction velocity component v along the z-axis (rotational contribution term), which is directly related to the formation of the vortex ring. The rate of change of the velocity component w along the x-axis (vertical momentum transfer) determines the strength of the collateral vortex; This represents the rate of change of the y-direction velocity component w along the y-axis (curl contribution term), which affects the morphology of the turbulent cavity region; The rate of change of the z-direction velocity component w along the z-axis (vertical compression or expansion) determines the dust settling rate.

[0046] More specifically, the humidity field is constructed as follows: Water mist data was collected, including installing a millimeter-wave sensor at the leading edge of the cutting teeth to measure the water film thickness using the dielectric constant; deploying a laser projection rollout on the top of the coal mine roadway to measure the light attenuation characteristics of the water mist; and embedding temperature and humidity sensors at the shock wave detection point to collect basic temperature and humidity data. The captured water mist data was then analyzed according to timestamps.

[0047] In the virtual space simulating explosive dust clouds, positioning beacons are placed at the boundaries of grid cells to establish a global coordinate system for the roadway with the cutter tooth as the origin. The roadway direction is defined as the X-axis (extension direction), and the vertical direction as the Z-axis. The pose transformation matrix is ​​acquired in real time through the UWB module, and a grid index mapping table for different sensors is established to map the sensor coordinates to the global coordinate system. Then, the global coordinate system is used to project the sensor coordinates onto the nearest neighboring grid cell in the virtual space, and cross-cell assignment is performed on the boundary points. Water mist data projected onto the same cell are aligned according to the timestamp of the sampling time, and weights are assigned according to the sensor accuracy. The water mist data projected from multiple grid cells are weighted and fused to calculate the humidity gradient field within the grid range in the virtual space. Cross-grid correlation analysis is then performed to solve the data continuity and correlation problems between different grid cells, generating a risk propagation path map of water film diffusion effects, i.e., a water film distribution heat map. The humidity compensation risk level of each grid cell is marked in the water film distribution heat map. The risk level is calculated based on the weighted calculation of different risk values ​​in the grid cell and compared with the preset risk threshold range, and is divided into low-risk, medium-risk, and high-risk levels.

[0048] The cross-cell allocation method is as follows: calculate the distance between the sensor projection point and the center points of all adjacent grid cells, denoted as . The minimum value among all its distances And compare it with the boundary influence threshold, It is determined to be a boundary point. Here, P represents the spatial coordinates of the sensor projection point; The spatial coordinates of the center point of the i-th adjacent grid cell; The size of the grid cell; The boundary influence factor determines the scaling factor of the boundary's effective range. All mesh elements satisfying the boundary point conditions are selected and merged into a candidate set, forming a set with a radius of [radius value missing] centered at the projection point. The spherical influence region. Then calculate the assignment weight of the projection point to each adjacent grid cell in the candidate set, using the following formula: In the formula, The weight assigned to the i-th adjacent grid cell; The spatial attenuation coefficient controls the rate attenuation of weights with distance; a larger value indicates a more concentrated allocation range, and vice versa. All assigned weights are then normalized to obtain the normalized assigned weights. Water mist data from sensor projection point P is allocated to adjacent grid cells according to the normalized assigned weights. Component data from boundary points within the same grid cell are accumulated, and intelligent gradient optimization is used to eliminate boundary effects and ensure the continuity of the physical field.

[0049] The cross-grid correlation analysis process is as follows: An adjacency matrix is ​​constructed for each grid cell, forming a basic network of connections between grid cells in physical space, providing a structural foundation for subsequent field variable transfer. A 26-neighborhood index is created for each grid cell; a directed connection graph structure is then constructed, with the center point of the grid cell as the vertices and adjacency relationships as edges; connection weights are determined based on the wind flow vector, with a higher weight for the downwind direction than the upwind direction; dynamic association weights are calculated based on the humidity gradient difference and the angle between the wind flow vectors; a weighted topology network reflecting spatial relationships and wind flow vectors is formed, prioritizing the strengthening of downwind associations. The formula for calculating the dynamic weight factor is: In the formula, The dynamic correlation weight between grid cell i and grid cell j reflects the spatial transmission characteristics of the humidity field. The larger the value, the stronger the physical correlation between the two grid cells. The norm of the humidity gradient vector interpolation for grid cell i and grid cell j quantifies the degree of abrupt change in humidity between the two grid cells. The larger the value, the more drastic the change in humidity field in that direction. Let i be the humidity gradient vector of grid cell i; The distance between the center points of grid cell i and grid cell j is the spatial distance between them; the ratio of the two measures the difference in humidity gradient per unit distance. Represents the wind vector in virtual space Direction of the line connecting the two grid cells The included angle; This is the wind direction factor; the higher the value, the stronger the dominant role of the wind vector, and vice versa. It is a nonlinear decay function based on the airflow angle.

[0050] Construct gradient continuity constraint equations at the boundaries of adjacent grid cells: In the formula, Let i be the humidity gradient of grid cell i. Let be the humidity gradient of grid cell j. Ideally, the humidity gradients of adjacent grid cells at the boundary should be equal to conform to the basic principle of smooth transmission of physical fields. However, in coal mines, wind turbulence can disrupt the continuity under ideal conditions, so a wind turbulence disturbance correction term is introduced. ,when This indicates that the airflow causes the humidity gradient on the i-th side of the grid cell to be greater than the humidity gradient on the j-th side. Conversely, by solving for the optimal smoothness of the global gradient field using the least squares method, the humidity gradient is ensured to transition smoothly between adjacent grids, avoiding abrupt changes. Simultaneously, based on the phase transition coupling mechanism, the liquid water grid transfers mass to the high-temperature steam grid, and the steam grid transfers condensed mass back to the low-temperature grid. This achieves cross-grid transfer of water molecules in different phases, maintaining the physical consistency of the phase transition process. By ensuring the physical continuity of the humidity gradient at grid cell boundaries, the mass transfer relationship during the phase transition process can be accurately quantified, thus establishing a precise coupling mechanism for establishing a dynamic water-vapor equilibrium.

[0051] Then construct the flux equation across grid cells: In the formula, This indicates that the calculation is performed on all directly adjacent grid cells j of grid cell i; The cross-boundary flux, i.e. the mass of humidity flowing from grid cell i to grid cell j per unit time, is calculated using the finite volume method. Let be the humidity change rate, and represent the rate of change of absolute humidity within grid cell i over time. Let i be the volume of grid cell i; maintain the total humidity mass conservation of the system and ensure the cross-grid conservation of water vapor mass under the action of airflow, so as to maintain physical rationality when data is missing.

[0052] Constructing the spatiotemporal correlation matrix Record the time delay of grid cell i. The influence threshold on grid cell j is set, and the future state is predicted based on the current gradient field. Upon actual arrival, the measured value is used to correct the dynamically associated weights, forming a prediction-correction closed loop for the evolution of the humidity field over time, thus improving the accuracy of dynamic process simulation. The formula for predicting the future state is: In the formula, Predict humidity values ​​for future moments; The actual humidity value measured at the current moment; This represents a three-dimensional airflow velocity field. For humidity gradient tensor; Where is the predicted time step. The formula for correcting the dynamic association weights is: In the formula, The corrected dynamic correlation weights; The sensitivity coefficient was calibrated based on experimental data to correct it. Here is the predicted humidity value for grid cell j; The actual measured humidity value for grid cell j; This is an exponential decay function used to convert prediction errors into weight decay ratios.

[0053] Specific propagation models for different types of risks within grid cells are constructed to quantify the propagation patterns of risks, visualize the risk transmission paths and impact ranges, and generate comprehensive early warning signals for the synergistic effects of multiple risks. The construction method is as follows: Risk types are defined, including three categories: optical deception zone, vapor artifact zone, and sensor blind zone. A propagation equation and termination threshold are established for each risk. Specifically, the risk propagation in the optical deception zone depends on the direction and magnitude of the moist dust gradient, and a diffusion equation is used: In the formula, The risk value of the optical deception zone quantifies the risk of optical measurements being interfered with by water film. The diffusion coefficient reflects the water film's ability to transmit risk. The thickness of the water film; The risk gradient for the water film is defined; and once the water film thickness falls below a set safety threshold (0.1-0.5 mm), the risk value is forcibly reset to zero. Risk propagation in the steam artifact zone depends on the airflow vector and temperature gradient, and is calculated using the advection equation: In the formula, The risk value of the steam artifact region is quantified, which quantifies the risk value of measurement error caused by high-temperature steam. For steam risk gradient; The vapor decay coefficient reflects the rate of vapor condensation; and vapor condenses below the dew point, forcing the risk value to zero. Risk propagation in sensor blind zones depends on the device's movement trajectory, employing a path tracing model. In the formula, For position Risk value of sensor blind zone; path is the device movement path; For Dirac functions; For mobile device trajectories; when location Risk is activated when a path point coincides with a location, and the risk value is accumulated based on the dwell time. The risk propagation weight of each connection direction in the adjacency matrix of the grid cells is calculated based on the airflow vector and humidity. Then, the risk propagation intensity of each propagation path is calculated based on the risk propagation weight. Primary and secondary propagation paths are identified using a graph algorithm and merged to construct a risk propagation path graph. In the formula, The final risk propagation path map represents a comprehensive risk field distribution map that integrates all risk types and their propagation paths. It identifies the risk level and propagation path of each grid cell; N is the number of risk types; k is the index of the risk type. Let M be the propagation intensity of risk k from grid cell i to grid cell j; M is the number of grid cells along the propagation path. Risk field data for risk k; This is a risk type superposition operator used to apply the risk propagation intensity in the risk propagation path to the corresponding generic type.

[0054] Example 2 Building upon Example 1, this paper proposes a downhole environment data acquisition system for 3D simulation environment modeling, comprising modules for data acquisition, CFD simulation, dynamic risk assessment, acquisition control, and data updating. The data acquisition module acquires the shock wave signal generated at the moment of collision between the cutting tooth and the coal seam, operating parameters such as the cutting tooth's rotational speed, feed rate, and coal seam hardness, as well as water mist data. The CFD simulation module performs four-stage turbulence evolution calculations based on the operating parameters and water mist data. The dynamic risk assessment module generates a five-dimensional reconstructed heatmap with risk labels based on the calculation results from the CFD simulation module. The acquisition control module switches acquisition strategies based on the region type in the five-dimensional reconstructed heatmap. The data updating module analyzes the waveform characteristics of the laser echo signal, detects anomalies, and optimizes the turbulence model using CFD parameters to improve the accuracy of fluid prediction in complex environments.

[0055] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A method for acquiring downhole environmental data for three-dimensional simulation environment modeling, characterized in that, Includes the following steps: The operating parameters are collected and input into the CFD model, and a downhole dynamic risk map is generated based on the shock wave signal; Collect water mist data to construct a humidity field, and embed the humidity field data into the generation process of the downhole dynamic risk map, including the injection compression stage, shear layer generation stage, vortex ring instability stage, and vortex full development stage; Finally, the simulation results of the above four stages are integrated into a reconstructed heat map, and the acquisition strategy of downhole environmental point cloud data is switched based on the risk area. By analyzing laser echo signals and detecting abnormal point cloud data in the downhole environment through distortion, and by reconstructing vortex equations using CFD parameters, the acquisition accuracy of downhole environmental point cloud data can be optimized.

2. The downhole environment data acquisition method as described in claim 1, characterized in that: During the jet compression period, the Navier-Stokes equations are solved, the grid in the key area is refined by pressure gradient detection, the morphology and density distribution of explosive dust cloud are generated, and the density field of dust and steam mixture is output based on the volume expansion effect generated by water vaporization phase change. When the shear layer is generated, the velocity gradient tensor of the boundary layer of the explosive dust cloud is calculated, and the vertical component of the velocity gradient tensor is corrected by the wet dust settling acceleration term to generate seed points for the turbulent development core containing water mist settling effect. When the vortex ring becomes unstable, the stability model of the vortex ring is reconstructed under the action of water mist viscous damping, the unstable region of explosive dust cloud flow is identified, and the location of vortex ring rupture is predicted. The inhibitory effect of water film on vortex ring rupture is quantified by adaptive vorticity threshold, and the turbulence morphology transformation feature map under humidity environment is output. When the eddies are fully developed, a water mist dissipation function is embedded in the K-ω turbulence model to simulate the evolution of multi-scale eddies. The real turbulent cavity region is identified by a dual threshold criterion of turbulent kinetic energy and humidity. Optical artifact regions and steam artifact regions are labeled simultaneously, and a distribution map of turbulent cavity regions with risk labels is output.

3. The downhole environment data acquisition method as described in claim 1, characterized in that: The method for encrypting the key area mesh is as follows: A virtual space for simulating explosive dust clouds is created. A Cartesian grid is created based on the geometry of the cutting teeth. The grid range of the virtual space is determined with the cutting teeth as the center, and the virtual space is divided into several grid units according to the preset grid size. Define the boundary and initial state of the explosive dust cloud flow, traverse all grid cells in the virtual space, calculate the pressure change rate and pressure gradient amplitude in three directions for each grid cell, compare the pressure gradient amplitude of the grid cell with the preset pressure gradient critical value, and divide the grid cell into normal grid and risk grid. The risk grid is decomposed into risk subgrids, and the data of the risk subgrids are integrated into the overall computation model according to the timestamps. The physical parameters of the affected grid cells are then iteratively updated based on the timestamps.

4. The downhole environment data acquisition method as described in claim 2, characterized in that: The method of decomposing into risk subgrids is as follows: The risk grid is decomposed into several small cubes, denoted as risk sub-grids. Each risk sub-grid inherits the physical properties of its parent grid. The boundary data of the risk sub-grids is then smoothed using the following method: Locate and mark the junction between the risk subgrid and the normal grid, and symmetrically extend the transition zone on both sides of the junction; Key data points of risk subgrids and normal grids are collected. Using the width of the transition zone as a baseline, the asymptotic transition value between the normal grid and risk subgrids is calculated based on the bilinear interpolation method. Forced physical constraints are set for the transition zone, including mass conservation, momentum continuity and energy stability. Repeat the above steps for iterative optimization until the numerical mutation rate of all sampling points is lower than the preset threshold, then stop the optimization.

5. The downhole environment data acquisition method as described in claim 3, characterized in that: The method for constructing the humidity field is as follows: A global coordinate system with the cutter tooth as the origin is established in the virtual space; the pose transformation matrix is ​​obtained in real time through the UWB module, a grid index mapping table for different sensors is established, the sensor coordinates are mapped to the global coordinates, and then the global coordinates are used to project them to the nearest neighboring grid cell in the virtual space, and cross-cell assignment is performed on the boundary points. The water mist data projected onto the same cell are aligned according to the timestamp of the sampling time, and weights are assigned according to the sensor accuracy. The water mist data projected from multiple cells within the grid cell are weighted and fused to calculate the humidity gradient field within the grid area in the virtual space. Then, cross-grid correlation analysis is performed on the humidity gradient field to generate a water film distribution heatmap affected by water film diffusion; the humidity compensation risk level of each grid cell is marked in the water film distribution heatmap.

6. The downhole environment data acquisition method as described in claim 4, characterized in that: The method for performing cross-cell allocation of boundary points is as follows: Calculate the distance between the sensor projection point and the center points of all adjacent grid cells, denoted as . ; Minimum of all its distances Compared with the boundary influence threshold, the minimum value The point is determined as a boundary point; where P is the spatial coordinate of the sensor projection point; The spatial coordinates of the center point of the i-th adjacent grid cell; The size of the grid cell; Boundary influence factor; Select all mesh cells that satisfy the boundary point conditions and merge them into a candidate set; Calculate the assignment weight of the projection point to each adjacent grid cell in the candidate set, and then normalize all assignment weights to obtain the normalized assignment weights. The water mist data at sensor projection point P is allocated to each adjacent grid cell according to the normalized assignment weight; the component data of each boundary point within the same grid cell are accumulated, and boundary effects are eliminated according to intelligent gradient optimization.

7. The downhole environment data acquisition method as described in claim 5, characterized in that: The method for performing cross-grid correlation analysis on the humidity gradient field is as follows: A directed connectivity graph with 26 neighborhood indices is constructed for each grid cell, and dynamic association weights are calculated based on humidity gradient difference and wind flow vector angle; thus forming a weighted topology network that reflects spatial relationships and wind flow vectors. Humidity mass conservation is maintained based on gradient continuity constraint equations and flux equations across grid cells; Constructing the spatiotemporal correlation matrix Record the time delay of grid cell i. The threshold of the influence on grid cell j is used to predict the future state based on the current gradient field, and the actual measured value is used to correct the dynamic association weight. Construct specific propagation models for different types of risks within grid cells, and visualize the risk transmission paths and impact ranges.

8. The downhole environment data acquisition method as described in claim 7, characterized in that: The reconstructed heatmap adopts a five-dimensional time-space tensor structure, namely: in, To reconstruct the heatmap; The coordinates are those of the actual turbulent cavity region; The duration of the actual turbulent cavity region; Humidity compensation risk level; Thermographic diagram of water film distribution; The light attenuation coefficient is the water mist light attenuation coefficient. The criteria for identifying optical artifacts are: the water film thickness is greater than the critical value of the water film thickness and the light attenuation coefficient of the water mist is greater than the light attenuation threshold. The criteria for identifying steam artifact regions are: humidity H is greater than the humidity saturation threshold and turbulent kinetic energy is less than the critical value for turbulence disappearance.

9. The downhole environment data acquisition method as described in claim 8, characterized in that: The risk zone switching-based data collection strategy includes: In the real turbulent cavity region, the laser sampling frequency is increased over a period of time to establish a spatiotemporal interpolation model to fill the topological fracture region, and a dynamic confidence factor C is used to perform weighted fusion of the filled point cloud. In the optical artifact region, the optical path offset is calculated based on the water film thickness and refractive index to perform beam path correction; the laser power is enhanced based on the optical attenuation coefficient; when the water film thickness exceeds the critical value, a moving median filter is used to eliminate signal flicker noise caused by water film reflection. In the vapor artifact region, the laser is switched to the 1570nm infrared band and pulse code modulation is used.

10. A downhole environment data acquisition system for three-dimensional simulation environment modeling, characterized in that, include: The data acquisition module is used for shock wave signals, operating parameters, and water mist data. The CFD simulation module performs four-stage turbulence evolution calculations based on operating parameters and water mist data. The dynamic risk assessment module generates a five-dimensional reconstructed heatmap with risk labels based on the calculation results of the CFD simulation module. The acquisition and control module is based on the region type switching acquisition strategy in the five-dimensional reconstructed heat map; The data update module is used to analyze the waveform characteristics of the laser echo signal, detect anomalies, and optimize the turbulence model by combining CFD parameters.

Images