Intelligent centralized control system and method for paste filling mining in coal mine

Abstract

The application discloses a coal mine paste filling mining intelligent centralized control system and method, relates to the technical field of coal mining control, and comprises the following steps: arranging a sensor array on the surface of a coal mining working face, a roof stratum and a filled area, and synchronously collecting surface subsidence, roof pressure and paste compactness data. A three-dimensional geomechanics digital model is constructed, the monitoring data is coupled with the model in real time for analysis, the dynamic deformation trend of the roof of a goaf is predicted through an improved deformation collaborative algorithm, and the paste filling opportunity and initial filling strength are calculated accordingly. The paste preparation and conveying subsystem is dispatched to run, the material ratio is configured according to the initial filling strength, the compactness data is compared with a preset threshold in real time during the filling process, the material ratio and the conveying working pressure are dynamically adjusted, and a closed-loop intelligent control is formed. The application is suitable for intelligent control of coal mine paste filling mining by predicting and driving filling parameter configuration and dynamically adjusting operation parameters in combination with compactness feedback.

Description

Technical Field

[0001] This invention belongs to the field of coal mine mining control technology, specifically based on an intelligent centralized control system and method for coal mine paste backfilling mining. Background Technology

[0002] Currently, coal mine paste backfilling operations mostly employ decentralized monitoring methods, collecting data only on roof pressure or the state of the backfill body. This results in a single data source and a lack of regional coordination. Geomechanical models are often statically constructed, failing to integrate with real-time on-site monitoring data for comprehensive analysis. Backfilling operation parameters rely heavily on empirical presets, with paste material ratios and conveying pressures remaining relatively fixed during operations, making it difficult to adapt to changes in the mining environment.

[0003] Existing technologies struggle to predict roof deformation trends in goaf areas in advance. There is a lack of corresponding control logic linking the roof condition to backfilling operations, and the determination of backfilling timing and intensity exhibits significant lag, failing to match the dynamic changes in the roof. Furthermore, the monitoring of the compaction of the backfill paste during the backfilling process lacks continuity, and the test results cannot be applied in real-time to adjust operational parameters. The backfill preparation and transportation processes also lack a closed-loop feedback mechanism.

[0004] This invention requires the simultaneous acquisition of multiple types of monitoring data from the surface, roof, and filling area. The monitoring data is coupled in real-time with a three-dimensional geomechanical digital model. An improved deformation coordination algorithm is used to predict the dynamic deformation trend of the roof, and the filling timing and initial filling intensity are determined based on the prediction results. During the filling process, the material ratio and conveying pressure are dynamically adjusted based on real-time density data, forming a complete intelligent closed-loop control process. This overcomes the shortcomings of existing filling control methods in data coordination, dynamic prediction, and real-time regulation. Summary of the Invention

[0005] This invention aims to solve at least one of the technical problems existing in the prior art; Therefore, this invention proposes an intelligent centralized control method for coal mine paste backfilling mining, comprising: Multiple sensor arrays are deployed in the surface area, roof strata, and filled area of ​​the target coal mining face. The sensor arrays simultaneously collect surface settlement data, roof pressure data, and compaction data of the filled paste. A three-dimensional geomechanical digital model of the target coal mining face is constructed, and the surface subsidence data and roof pressure data are coupled and analyzed in real time with the three-dimensional geomechanical digital model. The dynamic deformation trend of the roof in the goaf is predicted by an improved deformation coordination algorithm. Based on the dynamic deformation trend, the timing and initial filling strength of the paste filling used to maintain the stability of the top plate are calculated. Based on the timing of paste filling, the paste preparation and delivery subsystem is scheduled to start working, and the material ratio parameters in the preparation subsystem are set according to the initial filling intensity. During the filling process of the paste, the density data of the filled paste is acquired in real time and compared with the preset density threshold. Based on the comparison results, the material proportioning parameters and the working pressure of the conveying subsystem are dynamically adjusted to form a closed-loop control.

[0006] Furthermore, a three-dimensional geomechanical digital model of the target coal mining face is constructed, including: Collect geological exploration data of the target coal mining face, including coal seam thickness distribution, roof and floor lithological sequence, geological structure location and rock mechanical parameters; Based on the coal seam thickness distribution and the lithological sequence of the roof and floor, a three-dimensional stratigraphic structure grid model is established; The geological structure location is used as a constraint boundary and embedded into the three-dimensional stratigraphic mesh model; Using the rock strata mechanical parameters, each element in the three-dimensional stratigraphic structure mesh model is assigned elastic modulus, Poisson's ratio, and tensile strength properties; By integrating the three-dimensional stratigraphic structure mesh model and its unit cell properties, a three-dimensional geomechanical digital model capable of performing numerical simulations of stress and deformation is formed.

[0007] Furthermore, the surface subsidence data, roof pressure data, and the three-dimensional geomechanical digital model are coupled and analyzed in real time. An improved deformation coordination algorithm is used to predict the dynamic deformation trend of the goaf roof, including: The surface subsidence data is used as a displacement boundary condition and applied to the corresponding surface nodes of the three-dimensional geomechanical digital model. The roof pressure data is used as a stress load condition and applied to the corresponding unit of the three-dimensional geomechanical digital model simulating the roof of the goaf. Initiate transient mechanical calculations of the three-dimensional geomechanical digital model to simulate the stress redistribution process of the roof under mining disturbance; During the transient mechanical calculation process, the improved deformation coordination algorithm is called to process the deformation coordination relationship of the overlying rock strata in the goaf. The improved deformation coordination algorithm corrects the displacement transmission between rock strata by searching for contact between rock strata structures and judging slippage. Based on the corrected displacement transmission results, the displacement vector and subsidence velocity of each point on the roof of the goaf are output within a specified future time period. The displacement vector and subsidence velocity together constitute the dynamic deformation trend.

[0008] Furthermore, the improved deformation coordination algorithm is invoked to process the deformation coordination relationship of the overlying strata in the goaf, including: In the three-dimensional geomechanical digital model, lithological interfaces and potential weak structural surfaces are identified and marked as internal contact pairs; Define the normal contact stiffness, tangential contact stiffness, and coefficient of friction for each of the internal contact pairs; In each mechanical calculation iteration step, the internal contact pairs are searched for contact state to determine whether they are in an adhesive, slipping or separating state. For the internal contact pair in a slipping state, the maximum tangential resistance is calculated based on the friction coefficient and normal pressure, and the tangential displacement constraint is adjusted accordingly. For the internal contact pairs that are in a separated state, release their displacement constraints in the current iteration step, allowing the rock strata elements to separate freely; The updated contact state and displacement constraints are substituted into the overall stiffness matrix, and the mechanical equilibrium calculation is performed again to simulate the deformation discontinuity and energy dissipation process between rock layers, thereby correcting the displacement transmission between rock layers.

[0009] Furthermore, based on the dynamic deformation trend, the timing and initial filling strength of the paste filling used to maintain the stability of the roof slab are calculated, including: Extract the sinking speed time series of key points of the top plate from the dynamic deformation trend. When the rate of change of the sinking speed time series exceeds the preset warning threshold for the first time, trigger the filling timing signal. The filling timing signal is the filling timing of the paste. Based on the instantaneous subsidence velocity and cumulative subsidence of the key point of the roof when the filling timing signal is triggered, combined with the self-weight stress and overlying load of the rock strata at this location in the three-dimensional geomechanical digital model, the required support force is calculated in reverse. Divide the required support force by the designed filling area to convert it into the required strength of the paste. Based on the strength growth characteristic curve of the paste material, the initial solid mass concentration and cement content required to meet the strength requirements are determined. This solid mass concentration and cement content together define the initial filling strength.

[0010] Furthermore, based on the timing of the paste filling, the paste preparation and delivery subsystem is scheduled to start working, and the material proportioning parameters in the preparation subsystem are set according to the initial filling intensity, including: When the central controller receives the ointment filling timing signal, it sends a start command to the ointment preparation and delivery subsystem; Retrieve the formula table corresponding to the initial filling strength from the database. The formula table includes the set values ​​for solid mass concentration, cement mass ratio, fly ash mass ratio, and aggregate gradation parameters. The set values ​​for solid mass concentration, cement mass ratio, fly ash mass ratio, and aggregate gradation parameters are sent to the automatic batching controller of the preparation subsystem. The automatic batching controller controls the opening degree of the feeders in the aggregate bin, cement bin, and fly ash bin according to the received parameters, and controls the valve opening degree of the water supply pipeline to complete the weighing and mixing of materials and form a paste mixture that meets the initial filling strength.

[0011] Furthermore, during the paste filling process, the density data of the filled paste is acquired in real time and compared with a preset density threshold, including: In the filled area of ​​the target coal mining face, density sensors in the pre-embedded sensor array are used to periodically collect paste density values ​​at different depths. The density values ​​of the paste at different depths within the same vertical profile are weighted and averaged to obtain an average density value representing the overall compactness of the filled area. The preset compaction threshold corresponding to the current filling stage is read from the control strategy table. The preset compaction threshold is positively correlated with the mining depth and the roof pressure. The calculated average density value is compared with the preset density threshold to determine whether the average density value is lower than the preset density threshold, and the difference is calculated.

[0012] Furthermore, the step of dynamically adjusting the material proportioning parameters and the working pressure of the conveying subsystem based on the comparison results to form a closed-loop control includes: When the average density value is lower than the preset density threshold, the adjustment level is determined according to the magnitude of the difference; Based on the adjustment level, query the corresponding material ratio adjustment scheme in the control strategy table. The material ratio adjustment scheme includes the increase in solid mass concentration, the increase in cement ratio, and the adjustment value of additive dosage. The adjusted values ​​for the increase in solid mass concentration, the increase in cement ratio, and the adjustment of additive dosage are sent to the preparation subsystem to change the fluidity and early strength growth characteristics of the subsequently prepared paste. At the same time, based on the difference and the current length and height difference of the conveying pipeline, the working pressure setting value of the conveying subsystem is calculated and increased to ensure that the paste with increased viscosity after formula adjustment is smoothly conveyed to the goaf area. After the adjustment is executed, the density data of the filled paste in the next cycle is collected, and the comparison and adjustment steps are repeated until the average density value reaches or exceeds the preset density threshold.

[0013] Furthermore, the method also includes a global optimization step based on the filling process: In a complete coal mining cycle, record the time of each triggering of backfilling, the backfilling intensity used, the actual density of the backfill body achieved, and the stabilization time of the roof deformation after backfilling. The filling timing, initial filling strength, actual filling density, and stabilization time of the top plate deformation are used as input samples to train a filling efficiency evaluation model. Using the trained backfilling efficiency evaluation model, the predicted dynamic deformation trend in subsequent mining cycles is forward extrapolated to find the optimal combination of backfilling parameters that balances backfilling material cost and roof control effect. The optimal filling parameter combination is updated in the control strategy table to guide the filling control of subsequent mining cycles; The step of using the trained backfilling efficiency evaluation model to perform forward extrapolation of the predicted dynamic deformation trend in subsequent mining cycles includes: Key features of the dynamic deformation trend predicted by the improved deformation coordination algorithm in subsequent mining cycles are extracted. These key features include the peak value of the roof subsidence acceleration, the time to reach the peak value, and the deformation influence range. Input the key features along with the current material cost parameters into the filling performance evaluation model; The filling efficiency evaluation model outputs the expected evaluation score of multiple candidate filling parameter combinations, and the expected evaluation score is the weighted sum of the filling cost score and the roof stability prediction score. The candidate filling parameter combination with the highest expected evaluation score is selected as the optimal filling parameter combination that balances filling material cost and roof control effect.

[0014] Furthermore, the present invention also includes an intelligent centralized control system for coal mine paste backfilling mining, the system including a memory, a processor, and a computer program stored in the memory and running on the processor, wherein when the processor executes the computer program, it implements the steps of the intelligent centralized control method for coal mine paste backfilling mining as described above.

[0015] Compared with the prior art, the beneficial effects of the present invention are: Multiple sensor arrays are simultaneously deployed in the surface area, roof strata, and filled areas to simultaneously collect various types of data, including surface subsidence, roof pressure, and the density of the filled paste. A three-dimensional geomechanical digital model is continuously coupled with real-time monitoring data. An improved deformation coordination algorithm directly applies to the calculation of the dynamic deformation trend of the goaf roof, allowing for continuous tracking and quantification of roof state changes. The filling timing and initial filling intensity are directly calculated based on the roof's dynamic deformation trend. The paste preparation and conveying subsystems start operation according to the calculation results, establishing a correspondence between material proportioning parameters and roof stability requirements. The filling operation start point and intensity configuration are aligned with the real-time mechanical state of the stope.

[0016] During the filling process, the density data of the filled paste is continuously acquired and compared with a preset threshold. The material ratio parameters are changed in real time based on the comparison results, and the working pressure of the conveying subsystem is adjusted synchronously. The monitoring data and the actuator form a continuous interactive response. Changes in the density state directly drive the adjustment of parameters in the preparation and conveying stages. Deviations between the performance indicators of the filled body and the on-site working conditions are continuously corrected. The filling operation is in a dynamic adjustment state throughout the entire process, and the parameter settings are consistent with the actual filling effect. Attached Figure Description

[0017] Figure 1 This is a state diagram of the intelligent centralized control method for coal mine paste backfilling mining described in this invention; Figure 2 Flowchart for the construction of a three-dimensional geomechanical digital model; Figure 3 A flowchart for predicting the dynamic deformation trend of the roof in the goaf. Detailed Implementation

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

[0019] See Figure 1 This invention provides an intelligent centralized control method for coal mine paste backfilling mining, the specific method including: Multiple sensor arrays are deployed in the surface area, roof strata, and filled areas of the target coal mining face. These sensor arrays simultaneously collect surface subsidence data, roof pressure data, and the density data of the filled paste. A three-dimensional geomechanical digital model of the target coal mining face is constructed. The surface subsidence data and roof pressure data are coupled and analyzed in real time with the three-dimensional geomechanical digital model. An improved deformation coordination algorithm is used to predict the dynamic deformation trend of the roof in the goaf. Based on the dynamic deformation trend, the timing and initial filling intensity of the paste filling to maintain roof stability are calculated. Based on the paste filling timing, the paste preparation and conveying subsystem is scheduled to start working, and the material ratio parameters in the preparation subsystem are set according to the initial filling intensity. During the paste filling process, the density data of the filled paste is acquired in real time and compared with a preset density threshold. Based on the comparison results, the material ratio parameters and the working pressure of the conveying subsystem are dynamically adjusted to form a closed-loop control.

[0020] In one embodiment of the present invention, see [reference] Figure 2 The geological exploration data of the target coal mining face is collected, including coal seam thickness distribution, roof and floor lithological sequence, geological structure location, and rock mechanics parameters. Based on the coal seam thickness distribution and roof and floor lithological sequence, a three-dimensional stratigraphic structure mesh model is established. The geological structure location is used as a constraint boundary and embedded into the three-dimensional stratigraphic structure mesh model. Using the rock mechanics parameters, each unit cell in the three-dimensional stratigraphic structure mesh model is assigned elastic modulus, Poisson's ratio, and tensile strength properties. The three-dimensional stratigraphic structure mesh model and its unit cell properties are integrated to form the three-dimensional geomechanical digital model capable of performing numerical simulations of stress and deformation.

[0021] In practice, geological exploration data of the target coal mining face is collected. This data includes coal seam thickness distribution, roof and floor lithological sequence, geological structural location, and rock strata mechanical parameters. The coal seam thickness distribution is obtained through a combination of borehole logging and roadway exposure measurements, recording the roof and floor elevation data of the coal seam within the entire working face area with a grid precision of 10m × 10m. The roof and floor lithological sequence is determined based on the core logging of exploration boreholes and the interpretation results of geophysical logging curves, dividing the vertical stratification interfaces and thicknesses of different lithologies such as sandstone, mudstone, and limestone. The geological structural location includes the location coordinates and occurrence elements of fault fracture zones and fold axes, which are comprehensively located through three-dimensional seismic exploration interpretation and downhole advanced detection data. The rock strata mechanical parameters are obtained from uniaxial compression tests, Brazilian splitting tests, and triaxial shear tests conducted on various core specimens in the laboratory, obtaining the elastic modulus, Poisson's ratio, tensile strength, cohesion, and internal friction angle values ​​of each rock stratum. Based on the coal seam thickness distribution and the lithological sequence of the roof and floor, a three-dimensional stratigraphic structure mesh model was established. The modeling process adopted a tetrahedral and hexahedral hybrid element discretization technique, using the curved surfaces of the coal seam roof and floor as the main control surfaces. The target coal mining face and the rock strata space within a range of 100 meters above and 50 meters below were divided into a finite number of elements. The element size gradually increased with the distance from the coal seam being mined. The smallest element side length in the near-coal seam region was 0.5 meters, and the largest element side length in the far-field region was 5 meters. The geometric topology of the three-dimensional stratigraphic structure mesh model completely recorded the vertex coordinates, adjacency relationships, and lithological identifiers of each stratigraphic element. Geological structures are used as boundary constraints and embedded into the 3D stratigraphic mesh model. For major faults crossing the working face, zero-thickness contact unit layers are set at the contact surfaces of the hanging wall and footwall, and the spatial orientation of the contact units is adjusted according to the dip and strike of the fault. For fold structures, the attitude changes of the strata on both flanks are transformed into displacement constraint equations for the unit node, forcing the model's deformation pattern in this region to conform to the fold morphology. The embedding of geological structures ensures that the boundary conditions of the 3D stratigraphic mesh model can realistically reflect the influence of geological discontinuities on the stress transmission path. Using rock mechanics parameters, each unit in the 3D stratigraphic mesh model is assigned elastic modulus, Poisson's ratio, and tensile strength properties; the mathematical expressions for these property assignments are: ; in: This represents the equivalent elastic modulus, Poisson's ratio, or tensile strength of the element, where n represents the number of different lithological types contained within the element. This represents the volume occupied by the i-th lithology within a unit cell. This represents the numerical value of the mechanical parameter corresponding to the i-th lithology. This formula represents the total volume of a unit cell; it is applicable to handling mixed units spanning two or more lithological interfaces, reflecting the macroscopic mechanical behavior of composite rock strata through volume-weighted averaging. Integrating a three-dimensional stratigraphic structure mesh model and its unit cell properties, a three-dimensional geomechanical digital model capable of stress and deformation numerical simulation is formed. The data structure of the three-dimensional geomechanical digital model includes a set of node coordinates, a set of unit topologies, a set of material properties, and a set of constraint boundaries, which can be directly imported into finite element analysis software for gravity initialization, excavation simulation, and infill response calculation. The three-dimensional geomechanical digital model realizes the transformation from static geological description to dynamic mechanical simulation, becoming the basic platform for subsequent deformation collaborative analysis and infill decision-making.

[0022] In one embodiment of the present invention, see [reference] Figure 3 The surface subsidence data is loaded as displacement boundary conditions onto the corresponding surface nodes of the three-dimensional geomechanical digital model; the roof pressure data is loaded as stress load conditions onto the corresponding units simulating the roof of the goaf in the three-dimensional geomechanical digital model; the transient mechanical calculation of the three-dimensional geomechanical digital model is initiated to simulate the stress redistribution process of the roof under mining disturbance; during the transient mechanical calculation, the improved deformation coordination algorithm is called to handle the deformation coordination relationship of the overlying rock strata in the goaf. The improved deformation coordination algorithm corrects the displacement transmission between rock strata through rock stratum structural surface contact search and slip judgment; based on the corrected displacement transmission results, the displacement vector and subsidence velocity of each point of the roof of the goaf are output within a specified future time period. The displacement vector and subsidence velocity together constitute the dynamic deformation trend.

[0023] In the three-dimensional geomechanical digital model, lithological interfaces and potentially weak structural surfaces are identified and marked as internal contact pairs. Normal contact stiffness, tangential contact stiffness, and friction coefficient are defined for each internal contact pair. In each mechanical calculation iteration step, the contact state of the internal contact pairs is searched to determine whether they are in a bonded, slipped, or separated state. For internal contact pairs in a slipped state, the maximum tangential resistance is calculated based on the friction coefficient and normal pressure, and the tangential displacement constraint is adjusted accordingly. For internal contact pairs in a separated state, the displacement constraint in the current iteration step is released, allowing the rock strata elements to separate freely. The updated contact state and displacement constraint are substituted into the overall stiffness matrix, and a new mechanical equilibrium calculation is performed to simulate the deformation discontinuity and energy dissipation process between rock strata, thereby correcting the displacement transmission between rock strata.

[0024] In practice, surface subsidence data is used as displacement boundary conditions and applied to the corresponding surface nodes of the three-dimensional geomechanical digital model. The surface subsidence data comes from GNSS monitoring stations deployed above the working face and InSAR remote sensing interpretation results, recording the displacement of surface points in the east-west, north-south, and vertical directions with millimeter-level precision. In the three-dimensional geomechanical digital model, model nodes that match the plane coordinates of the measured surface points are selected. The vertical subsidence value obtained from GNSS monitoring is directly assigned to the Z-direction displacement degree of freedom of the corresponding node, and the horizontal displacement component is used as the displacement constraint input of the node in the X and Y directions, so that the surface deformation morphology of the model is consistent with the actual monitoring data. The roof pressure data is used as a stress load condition and applied to the corresponding unit of the roof in the three-dimensional geomechanical digital model simulating the goaf. The roof pressure data is collected by pressure sensors and roof delamination instruments installed on the hydraulic support columns, reflecting the vertical force and delamination deformation of the roof rock strata on the support system during mining. In the three-dimensional geomechanical digital model, the upper surface of the roof unit corresponding to the sensor installation position is located, the pressure sensor readings are converted into a uniformly distributed load applied to the surface of the unit, and the delamination deformation is converted into the relative displacement difference between the upper and lower parts of the unit, which serves as the basis for adjusting the initial stress field.

[0025] Transient mechanical calculations were initiated using a three-dimensional geomechanical digital model to simulate the stress redistribution process of the roof under mining disturbance. The transient mechanical calculations employed an explicit dynamic solver, with the mining advance time as the independent variable and a single-step time increment of 0.1 hours. In each step of the calculation, the life and death status of the excavation units was updated according to the position of the coal mining machine, the material contribution of the mined area units was removed, and release loads were applied to the boundaries of the newly formed goaf. During the calculation process, the magnitude and direction of the principal stresses of the roof units, the unit strain tensors, and the nodal displacement vectors were recorded for each step, and the migration trajectory of the stress concentration zone and the expansion range of the plastic zone were tracked. During transient mechanics calculations, an improved deformation coordination algorithm is invoked to handle the deformation coordination relationship of the overlying rock strata in the goaf. The improved deformation coordination algorithm corrects the displacement transmission between rock strata through contact search and slip judgment of rock strata structure surfaces. Based on the corrected displacement transmission results, the displacement vector and subsidence velocity of each point on the roof of the goaf are output within a specified time period in the future. The displacement vector and subsidence velocity together constitute the dynamic deformation trend. The output results are stored in matrix form, with the rows of the matrix corresponding to the grid node numbers on the roof surface, the columns corresponding to the hourly times in the next 24 hours, and the matrix elements representing the predicted subsidence amount and subsidence rate of the node at that time.

[0026] In the three-dimensional geomechanical digital model, lithological interfaces and potentially weak structural surfaces are identified and marked as internal contact pairs. Lithological interfaces refer to the junctions where there are significant differences in lithology between upper and lower strata, such as the interface between sandstone and mudstone. Potentially weak structural surfaces include known densely jointed zones and well-developed bedding planes. The identification process involves traversing the material property identifiers of adjacent units. When adjacent units belong to different lithologies or have weak interlayer labels, their common surface is registered as an internal contact pair, and a unique ID number is assigned to this contact pair. Normal contact stiffness, tangential contact stiffness, and friction coefficient are defined for each internal contact pair. Normal contact stiffness characterizes the contact surface's resistance to normal compressive deformation, and its value is proportional to the combined elastic modulus of the upper and lower strata. The calculation formula is: ; in: Indicates the normal contact stiffness. and These represent the elastic moduli of the upper and lower rock strata, respectively. and This represents the equivalent thickness of the upper and lower rock strata near the contact surface; the tangential contact stiffness is taken as 0.3 times the normal stiffness, and the friction coefficient is selected between 0.35 and 0.65 based on the roughness and water content of the upper and lower rock strata.

[0027] In each iteration of the mechanical calculation, a contact state search is performed on the internal contact pairs to determine whether they are in an adhesive, slip, or separation state. The contact state search is based on the relative displacement of the nodes on both sides of the contact surface and the magnitude of the contact force. When the normal pressure of the contact surface is positive and the tangential force is less than the frictional force, it is determined to be in an adhesive state; when the normal pressure is positive and the tangential force is greater than or equal to the frictional force, it is determined to be in a slip state; and when the normal pressure is negative, it is determined to be in a separation state. For internal contact pairs in the slip state, the maximum tangential resistance is calculated based on the friction coefficient and the normal pressure, and the tangential displacement constraint is adjusted accordingly. The calculation of the maximum tangential resistance adopts the Coulomb friction criterion, i.e. ,in This represents the maximum tangential resistance, and μ represents the coefficient of friction. This represents the normal contact pressure in the current iteration step; during adjustment, the limit value of the tangential relative displacement of the contact surface is set to the same as... The corresponding critical displacement is defined, and any displacement exceeding this value is considered irreversible plastic slip. For internal contact pairs in a separated state, their displacement constraints in the current iteration step are removed, allowing the rock strata elements to separate freely. In the separated state, the contact surface no longer transmits compressive and shear stresses, and the elements on both sides respond independently to external loads and inertial forces. The updated contact state and displacement constraints are substituted into the global stiffness matrix, and a new mechanical equilibrium calculation is performed to simulate the deformation discontinuity and energy dissipation process between rock strata, thereby correcting the displacement transmission between rock strata. Updating the global stiffness matrix involves increasing or decreasing the degrees of freedom of the contact elements and modifying the constraint equations. After the equilibrium calculation converges, the corrected nodal displacement field is output as the initial condition for the next calculation.

[0028] In one embodiment of the present invention, the subsidence velocity time series of key points of the roof is extracted from the dynamic deformation trend. When the rate of change of the subsidence velocity time series first exceeds a preset warning threshold, a filling timing signal is triggered, which is the filling timing of the paste. Based on the instantaneous subsidence velocity and cumulative subsidence of the key points of the roof when the filling timing signal is triggered, combined with the self-weight stress and overlying load of the rock strata at this location in the three-dimensional geomechanical digital model, the required support force is calculated. The required support force is divided by the designed filling area to convert it into the strength requirement of the paste. Based on the strength growth characteristic curve of the paste material, the initial solid mass concentration and cement content required to meet the strength requirement are determined. The solid mass concentration and cement content together define the initial filling strength.

[0029] When the central controller receives the paste filling timing signal, it sends a start command to the paste preparation and conveying subsystem; it retrieves the formula table corresponding to the initial filling strength from the database, the formula table including the solid mass concentration setting value, cement mass ratio setting value, fly ash mass ratio setting value, and aggregate gradation parameters; it sends the solid mass concentration setting value, cement mass ratio setting value, fly ash mass ratio setting value, and aggregate gradation parameters to the automatic batching controller of the preparation subsystem; the automatic batching controller controls the opening degree of the feeders of the aggregate bin, cement bin, and fly ash bin according to the received parameters, and controls the valve opening degree of the water supply pipeline to complete the weighing and mixing of materials to form a paste mixture that meets the initial filling strength.

[0030] In practice, the subsidence velocity time series of key points on the roof is extracted from the dynamic deformation trend. When the rate of change of the subsidence velocity time series exceeds the preset warning threshold for the first time, a filling timing signal is triggered, which is the timing of paste filling. The key points on the roof are located in the middle of the goaf, 120 meters away from the opening and 15 meters away from the return airway. The subsidence velocity time series is sampled at minute intervals, and the warning threshold is set at a subsidence acceleration of 5 mm per hour. When the real-time monitored subsidence velocity increases from 1.2 mm / hour to 6.5 mm / hour, the corresponding rate of change reaches 5.3 mm / hour², which exceeds the warning threshold. The central controller immediately generates a filling timing signal. Based on the instantaneous subsidence velocity and cumulative subsidence of key points in the roof at the time of triggering the filling signal, combined with the self-weight stress and overlying load of the rock strata at this location in the three-dimensional geomechanical digital model, the required support force is calculated. The instantaneous subsidence velocity is 6.5 mm / h, the cumulative subsidence is 38 mm, and the three-dimensional geomechanical digital model shows that the self-weight stress of the roof rock strata at this location is 12.5 MPa, and the additional load of the 300-meter overlying rock strata is 7.2 MPa. The calculation formula is: ; in: Indicates the required support force. This represents the stress due to the self-weight of the top strata. Indicates the additional load on top. Indicates the exposed area of ​​the roof slab. This indicates that the stress transfer efficiency coefficient is taken as 0.85. Indicates rock density, Represents gravitational acceleration. Indicates instantaneous sinking speed. This indicates the filling lag time; the required support force is calculated to be 28,500 kN. Dividing the required support force by the designed filling area converts it into the required strength of the paste; the designed filling area is 180 square meters, and the required paste strength is 158 kPa; based on the strength growth characteristic curve of the paste material, the initial solid mass concentration and cement content required to meet the strength requirements are determined, and the solid mass concentration and cement content together define the initial filling strength; the strength growth characteristic curve shows that at an ambient temperature of 20℃, a solid mass concentration of 76% corresponds to a 28-day paste strength of 165 kPa, and the cement content accounts for 11% of the solid mass (see Table 1).

[0031] Table 1: Comparison of Formulation Parameters Corresponding to Initial Filling Strength ; In some embodiments, when the central controller receives the paste filling timing signal, it sends a start command to the paste preparation and conveying subsystem. Within 500 milliseconds of generating the filling timing signal, the central controller completes signal verification and command encapsulation, and transmits the start command to the PLC controller of the preparation subsystem via industrial Ethernet. The central controller retrieves the formula table corresponding to the initial filling strength from the database. The formula table includes the set values ​​for solid mass concentration, cement mass percentage, fly ash mass percentage, and aggregate gradation parameters. The database stores formula records arranged according to strength gradients, indexed at 10 kPa intervals. The central controller queries the formula at the 160 kPa level with the closest strength requirement and retrieves the corresponding parameter set. The set values ​​for solid mass concentration, cement mass percentage, fly ash mass percentage, and aggregate gradation parameters are sent to the automatic batching controller of the preparation subsystem. Data transmission uses the Modbus TCP protocol, and parameter precision is retained to two decimal places. The automatic batching controller receives a solid mass concentration set value of 77%, a cement mass percentage set value of 12%, a fly ash mass percentage set value of 21%, and an aggregate gradation parameter of 2.6. The automatic batching controller, based on received parameters, controls the opening of the feeders in the aggregate bins, cement bins, and fly ash bins, as well as the valve opening of the water supply pipeline, to complete the weighing and mixing of materials, forming a paste mixture that meets the initial filling strength. The automatic batching controller sets the feeder opening to 78%, adjusts the cement bin screw conveyor speed to 45 rpm, controls the fly ash pneumatic conveying flow rate at 1.2 cubic meters per minute, adjusts the water supply valve opening to 42%, maintains the mixer operating power at 55 kW, and controls the slump of the produced paste within the range of 210 ± 10 mm. The paste formulations selected in this embodiment have all been verified by indoor mixing tests. The relationship between the short-term support strength and the final strength at 28 days during the initial setting stage of the paste is recorded simultaneously for each combination of solid mass concentration and cement dosage. When the system retrieves the formulation table, it simultaneously reads the early strength reference value corresponding to the ratio. If the instantaneous deformation rate of the current roof is too high, the system will automatically increase the proportion of cement and early strength additives to prioritize improving the early bearing capacity of the paste during initial setting, taking into account both short-term roof support and long-term surrounding rock stability requirements, and eliminating the time scale contradiction between the final strength calculation and the early stress condition of the fluid paste on site.

[0032] It is understandable that, in the paste filling timing judgment stage, the calculation of the rate of change of the sinking velocity time series adopts the five-point central difference method to reduce the risk of false triggering caused by noise interference. The warning threshold is dynamically adjusted based on the working face mining depth and roof lithology classification, with a higher threshold for shallow hard rock working faces and a lower threshold for deep soft rock working faces. It is also understandable that, in the formula parameter execution stage, the automatic batching controller has an online calibration function. When the weighing sensor detects that the actual feed rate deviates from the set value by more than 2%, it automatically corrects the feeder opening compensation coefficient to ensure that the final paste ratio error is controlled within ±1.5%.

[0033] In one embodiment of the present invention, in the filled area of ​​the target coal mining face, density sensors in the pre-embedded sensor array are used to periodically collect the paste density values ​​at different depths. The specific deployment method is as follows: before the filling operation, shallow hole anchors are constructed in advance on the bottom plate of the filling section. Multiple sets of sensor mounting brackets are vertically and layered upwards using rigid metal support rods. The brackets are rigidly locked to the support rods by locking buckles. Each density sensor is covered with a hollow steel protective cage. The protective cage is welded and fixed to the bracket to form an integrated rigid monitoring unit. The entire monitoring unit relies on the bottom plate anchors and vertical support rods to form a stable support system, completely isolating the sensor body from lateral and vertical impact loads of the flowing paste. The filling operation proceeds slowly in layers along one side of the tunnel. The paste is poured in layers and the monitoring units are buried in layers. The thickness of each layer of paste is controlled within the spacing between the layers of the bracket. The impact load of the paste flow is transmitted to the ground anchor on the bottom plate through the outer steel protective cage and rigid support rod. The sensor body does not directly bear the impact force of the slurry. At the same time, the sensor signal line is laid along the inner side of the support rod and fixed in sections with pipe clamps to prevent the slurry from washing away and pulling the cable, causing the sensor to shift.

[0034] The density values ​​of the paste at different depths within the same vertical profile are weighted and averaged to calculate an average density value representing the overall compactness of the filled area. The preset density threshold corresponding to the current filling stage is read from the control strategy table. The preset density threshold is positively correlated with the mining depth and roof pressure. The calculated average density value is compared with the preset density threshold to determine whether the average density value is lower than the preset density threshold, and the difference is calculated.

[0035] When the average density value is lower than the preset density threshold, the adjustment level is determined based on the difference. Based on the adjustment level, the corresponding material ratio adjustment scheme in the control strategy table is queried. This scheme includes adjustments for the solid mass concentration increase, cement ratio increase, and additive dosage. The adjusted solid mass concentration increase, cement ratio increase, and additive dosage adjustments are sent to the preparation subsystem to alter the fluidity and early strength growth characteristics of the subsequently prepared paste. Simultaneously, based on the difference and the current length and elevation difference of the conveying pipeline, the working pressure setting of the conveying subsystem is calculated and increased to ensure that the paste with increased viscosity after formula adjustment is smoothly conveyed to the goaf. After the adjustment is executed, the density data of the filled paste for the next cycle is collected, and the comparison and adjustment steps are repeated until the average density value reaches or exceeds the preset density threshold.

[0036] In practical implementation, density sensors in a pre-embedded sensor array are used to periodically collect paste density values ​​at different depths in the filled area of ​​the target coal mining face. A monitoring section is arranged every 30 meters along the filling area, with sensors installed at three height levels (bottom, middle, and top) of the filling body at each section. The data collection frequency is set to one set every 10 minutes. The collected paste density values ​​range from 1.85 tons / cubic meter to 2.05 tons / cubic meter, reflecting the compaction degree of the filling body at different solidification stages. A weighted average is calculated on the paste density values ​​at different depths within the same vertical profile to obtain the average compaction value representing the overall compaction degree of the filling body in the filled area. The weighting coefficient is determined based on the paste volume ratio corresponding to the sensor depth: 0.5 for the bottom sensor, 0.3 for the middle sensor, and 0.2 for the top sensor. The calculation formula is as follows: ; in: This represents the average density value. This represents the weighting coefficient of the k-th layer sensor. This represents the paste density value collected by the sensor at layer k. The preset compaction threshold corresponding to the current filling stage is read from the control strategy table. The preset compaction threshold is positively correlated with the mining depth and roof pressure. The control strategy table stores threshold parameters for different combinations of mining depth ranges (200-300 meters, 300-400 meters, 400-500 meters) and roof pressure levels (low, medium, high). For example, the threshold is 1.92 tons / cubic meter when the mining depth is 350 meters and the roof pressure is medium. The calculated average compaction value is compared with the preset compaction threshold to determine if the average compaction value is lower than the preset threshold, and the difference is calculated. If the calculated average compaction value is 1.88 tons / cubic meter, which is lower than the threshold of 1.92 tons / cubic meter, the difference is -0.04 tons / cubic meter.

[0037] In some embodiments, when the average density value is lower than a preset density threshold, the adjustment level is determined based on the magnitude of the difference. The adjustment level is divided into three levels: Level 1 adjustment is 0-0.02 tons / cubic meter, Level 2 adjustment is 0.02-0.04 tons / cubic meter, Level 3 adjustment is greater than 0.04 tons / cubic meter, and Level 3 adjustment is a difference of -0.04 tons / cubic meter. Based on the adjustment level, the corresponding material proportion adjustment scheme is queried from the control strategy table. The material proportion adjustment scheme includes the increase in solid mass concentration, the increase in cement proportion, and the adjustment value of additive dosage. The parameters for the Level 3 adjustment scheme in the control strategy table are shown in Table 2. Table 2: Material Proportion Adjustment Scheme Corresponding to the Three Adjustment Levels ; The adjusted increases in solid mass concentration, cement ratio, and additive dosage are sent to the preparation subsystem to alter the fluidity and early strength growth characteristics of the prepared paste. Upon receiving the adjustment instructions, the preparation subsystem increases the aggregate feed rate by 12%, the cement silo rotation speed by 18%, the early strength agent metering pump frequency by 8 Hz, the water-reducing agent dosage by 240 ml / m³, and the mixing time by 15 seconds, reducing the paste slump from 220 mm to 190 mm and shortening the initial setting time by 35 minutes. Simultaneously, based on the differences and the current length and elevation difference of the delivery pipeline, the working pressure setting of the delivery subsystem is calculated and increased to ensure the smooth delivery of the increased viscosity paste to the goaf. The total length of the delivery pipeline is 850 meters, with a vertical elevation difference of 120 meters. The original working pressure was 8.5 MPa, which is increased to 11.2 MPa after adjustment, while the pumping flow rate remains at 90 m³ / h. After the adjustment is executed, the density data of the filled paste in the next cycle is collected, and the comparison and adjustment steps are repeated until the average density value reaches or exceeds the preset density threshold. The average density value collected in the next cycle is 1.94 tons / cubic meter, which is higher than the threshold of 1.92 tons / cubic meter. The system stops adjusting and maintains the current ratio parameters.

[0038] Optionally, the installation location of the density sensor can be dynamically selected according to the shape of the filling area. For irregular filling areas, the sensor density in the inclined boreholes can be increased to ensure that the monitoring section covers the key load-bearing area of ​​the filling body. For rapid filling of large-flow paste, additional lateral limiting steel cables can be added. The two ends of the steel cables are anchored to the sides of the roadway on both sides of the filling area, laterally restraining the horizontal displacement of the sensor brackets at each layer, further improving the overall impact resistance stability of the pre-embedded monitoring unit, and ensuring that the density sensor maintains the preset installation coordinates throughout the entire process of the paste from the flow state to the initial setting. The monitoring data can stably reflect the true density of the filling body at different depths. Optionally, the weighted average calculation can use moving window filtering technology, taking the density values ​​of the most recent three periods for moving average, eliminating the interference of random measurement errors on the density determination. It can be understood that the preset density threshold in the control strategy table will be dynamically updated as mining progresses. When the working face enters the fault influence zone, the system automatically lowers the threshold by 0.03-0.05 tons / cubic meter to adapt to roof pressure fluctuations. It is understandable that the working pressure adjustment of the delivery subsystem adopts a segmented and gradual strategy, with each pressure increase not exceeding 3 MPa, to avoid a sudden increase in pipeline pressure that could lead to a pipe burst.

[0039] In one embodiment of the present invention, during a complete coal mining cycle, the time of triggering backfilling, the backfilling intensity used, the actual density of the backfill body achieved, and the stabilization time of roof deformation after backfilling are recorded for each instance. The timing of paste backfilling, the initial backfilling intensity, the actual density of the backfill body achieved, and the stabilization time of roof deformation are used as input samples to train a backfilling efficiency evaluation model. Using the trained backfilling efficiency evaluation model, the predicted dynamic deformation trend in subsequent mining cycles is forward extrapolated to find the optimal combination of backfilling parameters that balances backfilling material cost and roof control effect. The optimal combination of backfilling parameters is updated in the control strategy table to guide backfilling control in subsequent mining cycles.

[0040] Key features of the dynamic deformation trend predicted by the improved deformation coordination algorithm in subsequent mining cycles are extracted. These key features include the peak value of the roof subsidence acceleration, the time to reach the peak value, and the deformation influence range. The key features are input into the backfilling efficiency evaluation model along with the current material cost parameters. The backfilling efficiency evaluation model outputs expected evaluation scores for multiple candidate backfilling parameter combinations. The expected evaluation scores are the weighted sum of the backfilling cost score and the roof stability prediction score. The candidate backfilling parameter combination with the highest expected evaluation score is selected as the optimal backfilling parameter combination that balances backfilling material cost and roof control effect.

[0041] In practice, during a complete coal mining cycle, the time of each backfilling trigger, the backfilling intensity used, the actual density of the backfill, and the stabilization time of the roof deformation after backfilling are recorded. A coal mining cycle is defined as the entire process from the opening cut to the stop line, which includes 52 backfilling operations. Recorded fields include the backfilling trigger time (accurate to the second), the initial solid mass concentration of the paste and the cement content, the measured average density 24 hours after the end of backfilling, and the duration for which the roof monitoring point subsidence rate recovers to below 0.1 mm / hour. An example of a single record is as follows: the backfilling time is 14:26:18 on the 3rd day, the backfilling intensity corresponds to a solid mass concentration of 77% and a cement content of 12%, the measured density is 1.95 tons / cubic meter, and the roof stabilization time is 48 hours. A backfilling efficiency evaluation model was trained using the backfilling timing, initial backfilling strength, actual backfill density, and roof deformation stabilization time as input samples. A total of 52 input samples were used, each containing four feature dimensions and one target variable. After feature dimension normalization, the samples were input into a three-layer fully connected neural network. The target variable was the comprehensive efficiency score. The network output layer used a linear activation function, and the training iterations were set to 1000 times with a fixed learning rate of 0.001. After training, the weight parameters of the backfilling efficiency evaluation model were fixed and deployed to the real-time inference module of the central controller. Using the trained backfilling efficiency evaluation model, the predicted dynamic deformation trend in subsequent mining cycles was forward-engineered to find the optimal combination of backfilling parameters that balances backfilling material cost and roof control effect. The forward-engineering cycle was set to the estimated time span of the next coal mining cycle. The input data came from historical records of the previous cycle and the predicted dynamic deformation trend characteristics of the current cycle. The output was the ranking results of multiple candidate backfilling parameter combinations. The optimal filling parameter combination is updated to the control strategy table to guide the filling control of subsequent mining cycles. The control strategy table is version-managed, with new parameter combinations overwriting existing entries and being synchronously distributed to the local configuration cache of each subsystem to ensure that it takes effect immediately when the next filling is triggered.

[0042] Key features of the dynamic deformation trend predicted by the improved deformation coordination algorithm in subsequent mining cycles are extracted. These key features include the peak value of the roof subsidence acceleration, the time of reaching the peak value, and the deformation influence range. The prediction data comes from the roof deformation simulation results after the current working face advances 120 meters eastward. The peak subsidence acceleration is 8.2 mm / h², occurring 16 hours after mining. The deformation influence range is an elliptical area with a dip width of 45 meters and a strike length of 60 meters. The key features, along with current material cost parameters, are input into the backfilling efficiency evaluation model. These material cost parameters include cement cost of 520 yuan / ton, fly ash cost of 180 yuan / ton, aggregate cost of 65 yuan / ton, and electricity cost of 0.68 yuan / kWh. These parameters are concatenated as a vector after the key features to form an extended input vector. The backfilling efficiency evaluation model outputs the expected evaluation score for multiple combinations of candidate backfilling parameters. The expected evaluation score is a weighted sum of the backfilling cost score and the roof stability prediction score. The scoring formula is as follows: ; in: Indicates the expected assessment score. This indicates the filling cost score (the lower the cost, the higher the score). This indicates the roof stability prediction score (the longer the stability time, the higher the score). and The weighting coefficients for cost and stability are 0.35 and 0.65, respectively. The candidate filling parameter combination with the highest expected evaluation score is selected as the optimal filling parameter combination that balances the cost of filling materials and the control effect on the roof. The candidate combinations include permutations of solid mass concentration of 76%-82% (step size 2%) and cement content of 11%-15% (step size 1%). After calculation, the combination with solid mass concentration of 78% and cement content of 13% scored the highest and was determined as the optimal parameter.

[0043] In some embodiments, a regularization term is added to the training process of the filling efficiency evaluation model to prevent overfitting. The L2 regularization coefficient is set to 0.01, the validation set accounts for 20% of the total samples, and the training stops when the validation set loss no longer decreases after 50 consecutive rounds. The model inference stage adopts a batch processing mode, calculating the scores of all candidate combinations at once, with sorting time controlled within 200 milliseconds. Optionally, key feature extraction can employ principal component analysis to reduce the dimensionality of the displacement time-series data of each point on the roof to three principal component factors, retaining 97% of the variance information of the original data and reducing model input redundancy. It is understood that the weight coefficients of the expected evaluation scores can be dynamically adjusted according to the mine safety policy. When the priority of roof control is increased, the stability weight β is increased to 0.8, and the cost weight α is correspondingly decreased.

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

Claims

1. A method for intelligent centralized control of coal mine paste backfilling mining, characterized in that, include: Multiple sensor arrays are deployed in the surface area, roof strata, and filled area of ​​the target coal mining face. The sensor arrays simultaneously collect surface settlement data, roof pressure data, and compaction data of the filled paste. A three-dimensional geomechanical digital model of the target coal mining face is constructed, and the surface subsidence data and roof pressure data are coupled and analyzed in real time with the three-dimensional geomechanical digital model. The dynamic deformation trend of the roof in the goaf is predicted by an improved deformation coordination algorithm. Based on the dynamic deformation trend, the timing and initial filling strength of the paste filling used to maintain the stability of the top plate are calculated. Based on the timing of paste filling, the paste preparation and delivery subsystem is scheduled to start working, and the material ratio parameters in the preparation subsystem are set according to the initial filling intensity. During the filling process of the paste, the density data of the filled paste is acquired in real time and compared with the preset density threshold. Based on the comparison results, the material proportioning parameters and the working pressure of the conveying subsystem are dynamically adjusted to form a closed-loop control.

2. The intelligent centralized control method for coal mine paste backfilling mining according to claim 1, characterized in that, Constructing a three-dimensional geomechanical digital model of the target coal mining face includes: Collect geological exploration data of the target coal mining face, including coal seam thickness distribution, roof and floor lithological sequence, geological structure location and rock mechanical parameters; Based on the coal seam thickness distribution and the lithological sequence of the roof and floor, a three-dimensional stratigraphic structure grid model is established; The geological structure location is used as a constraint boundary and embedded into the three-dimensional stratigraphic mesh model; Using the rock strata mechanical parameters, each element in the three-dimensional stratigraphic structure mesh model is assigned elastic modulus, Poisson's ratio, and tensile strength properties; By integrating the three-dimensional stratigraphic structure mesh model and its unit cell properties, a three-dimensional geomechanical digital model capable of performing numerical simulations of stress and deformation is formed.

3. The intelligent centralized control method for coal mine paste backfilling mining according to claim 1, characterized in that, The surface subsidence data, roof pressure data, and the three-dimensional geomechanical digital model are coupled and analyzed in real time. An improved deformation synergy algorithm is used to predict the dynamic deformation trend of the goaf roof, including: The surface subsidence data is used as a displacement boundary condition and applied to the corresponding surface nodes of the three-dimensional geomechanical digital model. The roof pressure data is used as a stress load condition and applied to the corresponding unit of the three-dimensional geomechanical digital model simulating the roof of the goaf. Initiate transient mechanical calculations of the three-dimensional geomechanical digital model to simulate the stress redistribution process of the roof under mining disturbance; During the transient mechanical calculation process, the improved deformation coordination algorithm is called to process the deformation coordination relationship of the overlying rock strata in the goaf. The improved deformation coordination algorithm corrects the displacement transmission between rock strata by searching for contact between rock strata structures and judging slippage. Based on the corrected displacement transmission results, the displacement vector and subsidence velocity of each point on the roof of the goaf are output within a specified future time period. The displacement vector and subsidence velocity together constitute the dynamic deformation trend.

4. The intelligent centralized control method for coal mine paste backfilling mining according to claim 3, characterized in that, The improved deformation coordination algorithm is invoked to process the deformation coordination relationship of the overlying strata in the goaf, including: In the three-dimensional geomechanical digital model, lithological interfaces and potential weak structural surfaces are identified and marked as internal contact pairs; Define the normal contact stiffness, tangential contact stiffness, and coefficient of friction for each of the internal contact pairs; In each mechanical calculation iteration step, the internal contact pairs are searched for contact state to determine whether they are in an adhesive, slipping or separating state. For the internal contact pair in a slipping state, the maximum tangential resistance is calculated based on the friction coefficient and normal pressure, and the tangential displacement constraint is adjusted accordingly. For the internal contact pairs that are in a separated state, release their displacement constraints in the current iteration step, allowing the rock strata elements to separate freely; The updated contact state and displacement constraints are substituted into the overall stiffness matrix, and the mechanical equilibrium calculation is performed again to simulate the deformation discontinuity and energy dissipation process between rock layers, thereby correcting the displacement transmission between rock layers.

5. The intelligent centralized control method for coal mine paste backfilling mining according to claim 1, characterized in that, Based on the dynamic deformation trend, the timing and initial filling strength of the paste filling used to maintain the stability of the roof slab are calculated, including: Extract the sinking speed time series of key points of the top plate from the dynamic deformation trend. When the rate of change of the sinking speed time series exceeds the preset warning threshold for the first time, trigger the filling timing signal. The filling timing signal is the filling timing of the paste. Based on the instantaneous subsidence velocity and cumulative subsidence of the key point of the roof when the filling timing signal is triggered, combined with the self-weight stress and overlying load of the rock strata at this location in the three-dimensional geomechanical digital model, the required support force is calculated in reverse. Divide the required support force by the designed filling area to convert it into the required strength of the paste. Based on the strength growth characteristic curve of the paste material, the initial solid mass concentration and cement content required to meet the strength requirements are determined. This solid mass concentration and cement content together define the initial filling strength.

6. The intelligent centralized control method for coal mine paste backfilling mining according to claim 1, characterized in that, Based on the timing of paste filling, the paste preparation and delivery subsystem is scheduled to start working, and the material proportioning parameters in the preparation subsystem are set according to the initial filling intensity, including: When the central controller receives the ointment filling timing signal, it sends a start command to the ointment preparation and delivery subsystem; Retrieve the formula table corresponding to the initial filling strength from the database. The formula table includes the set values ​​for solid mass concentration, cement mass ratio, fly ash mass ratio, and aggregate gradation parameters. The set values ​​for solid mass concentration, cement mass ratio, fly ash mass ratio, and aggregate gradation parameters are sent to the automatic batching controller of the preparation subsystem. The automatic batching controller controls the opening degree of the feeders in the aggregate bin, cement bin, and fly ash bin according to the received parameters, and controls the valve opening degree of the water supply pipeline to complete the weighing and mixing of materials and form a paste mixture that meets the initial filling strength.

7. The intelligent centralized control method for coal mine paste backfilling mining according to claim 1, characterized in that, During the paste filling process, the density data of the filled paste is acquired in real time and compared with a preset density threshold, including: In the filled area of ​​the target coal mining face, density sensors in the pre-embedded sensor array are used to periodically collect paste density values ​​at different depths. The density values ​​of the paste at different depths within the same vertical profile are weighted and averaged to obtain an average density value representing the overall compactness of the filled area. The preset compaction threshold corresponding to the current filling stage is read from the control strategy table. The preset compaction threshold is positively correlated with the mining depth and the roof pressure. The calculated average density value is compared with the preset density threshold to determine whether the average density value is lower than the preset density threshold, and the difference is calculated.

8. The intelligent centralized control method for coal mine paste backfilling mining according to claim 7, characterized in that, The process of dynamically adjusting the material proportioning parameters and the working pressure of the conveying subsystem based on the comparison results to form a closed-loop control includes: When the average density value is lower than the preset density threshold, the adjustment level is determined according to the magnitude of the difference; Based on the adjustment level, query the corresponding material ratio adjustment scheme in the control strategy table. The material ratio adjustment scheme includes the increase in solid mass concentration, the increase in cement ratio, and the adjustment value of additive dosage. The adjusted values ​​for the increase in solid mass concentration, the increase in cement ratio, and the adjustment of additive dosage are sent to the preparation subsystem to change the fluidity and early strength growth characteristics of the subsequently prepared paste. At the same time, based on the difference and the current length and height difference of the conveying pipeline, the working pressure setting value of the conveying subsystem is calculated and increased to ensure that the paste with increased viscosity after formula adjustment is smoothly conveyed to the goaf area. After the adjustment is executed, the density data of the filled paste in the next cycle is collected, and the comparison and adjustment steps are repeated until the average density value reaches or exceeds the preset density threshold.

9. The intelligent centralized control method for coal mine paste backfilling mining according to claim 8, characterized in that, The method also includes a global optimization step based on the filling process: In a complete coal mining cycle, record the time of each triggering of backfilling, the backfilling intensity used, the actual density of the backfill body achieved, and the stabilization time of the roof deformation after backfilling. The filling timing, initial filling strength, actual filling density, and stabilization time of the top plate deformation are used as input samples to train a filling efficiency evaluation model. Using the trained backfilling efficiency evaluation model, the predicted dynamic deformation trend in subsequent mining cycles is forward extrapolated to find the optimal combination of backfilling parameters that balances backfilling material cost and roof control effect. The optimal filling parameter combination is updated in the control strategy table to guide the filling control of subsequent mining cycles; The step of using the trained backfilling efficiency evaluation model to perform forward extrapolation of the predicted dynamic deformation trend in subsequent mining cycles includes: Extract key features of the dynamic deformation trend predicted by the improved deformation coordination algorithm in subsequent mining cycles. The key features include the peak value of the roof subsidence acceleration, the time to reach the peak value, and the deformation influence range. Input the key features along with the current material cost parameters into the filling performance evaluation model; The filling efficiency evaluation model outputs the expected evaluation score of multiple candidate filling parameter combinations, and the expected evaluation score is the weighted sum of the filling cost score and the roof stability prediction score. The candidate filling parameter combination with the highest expected evaluation score is selected as the optimal filling parameter combination that balances filling material cost and roof control effect.

10. A centralized intelligent control system for coal mine paste backfilling mining, comprising a memory, a processor, and a computer program stored in the memory and running on the processor, characterized in that, When the processor executes the computer program, it implements the steps of the intelligent centralized control method for coal mine paste backfilling mining as described in any one of claims 1 to 9.

Images