A method and system for managing the recovery of caking ore in a stope

Abstract

This application discloses a method and system for managing the recovery of bonded ore in a stope, belonging to the field of computer technology. The technical solution mainly involves digital signal processing and general control or regulation systems. By fusing multi-source sensing data, a heterogeneous coupled digital twin model is constructed, achieving a complete understanding of the internal structure, bonding causes, and mechanical properties of the bonded ore pile. Through intelligent simulation optimization, a target scheme balancing ore recovery efficiency and operational safety can be generated, avoiding the lack of basis and potential damage associated with traditional empirical blasting. Simultaneously, by combining real-time construction feedback for model assimilation and dynamic adjustment of the scheme, the efficiency, safety, and resource recovery rate of bonded ore recovery operations are effectively improved.

Description

Technical Field

[0001] This application relates to the field of computer technology, and in particular to a method and system for managing the recovery of bonded ore in a mining area. Background Technology

[0002] In mines employing mining methods such as large-diameter deep-hole backfilling, large-scale bonded ore heaps often form at the bottom of the stope due to the high sulfur content and abundance of fine-grained minerals in the ore. These ore heaps are strong and tightly bonded to the surrounding rock, making it impossible for conventional shoveling and transport equipment to operate directly, resulting in severe resource losses and affecting the safety of subsequent operations in the stope.

[0003] Currently, the industry generally relies on empirical blasting methods, such as vibratory blasting with fan-shaped holes in adjacent roadways. These methods have several drawbacks: First, they treat the bonded ore body as a "black box," lacking knowledge of its internal structure, bonding origin, and strength distribution, leading to blindly setting blasting parameters and highly unstable results. Second, in pursuit of loosening, excessive explosive charges are often used, which can easily cause uncontrollable damage to critical structures such as the bottom structure of the stope and ore extraction roadways, posing significant safety hazards. Finally, the operational plan is isolated and static, failing to assess its disturbance impact within the global stress field and unable to be dynamically adjusted according to real-time geological conditions.

[0004] Therefore, existing technologies lack a systematic solution that can integrate multi-source sensing, quantitative diagnosis of the cause of adhesion, and intelligent simulation and dynamic control within the framework of a global stress field. This has led to the long-term severe challenges of low efficiency, poor safety, and low resource recovery rate in the recovery of bonded ore. Summary of the Invention

[0005] This application provides a method and system for managing the recovery of bonded ore in a stope, the technical solution of which is as follows:

[0006] On the one hand, a method for managing the recovery of caking ore in a stope is provided, the method comprising:

[0007] The three-dimensional laser scanning data, internal physical detection data, and core genetic diagnosis data of the bonded ore pile are fused to obtain the spatial morphology and internal structural defect distribution of the bonded ore pile. Based on the differences in chemical and mineral composition indicated by the core genetic diagnosis data, multiple bonding type zones are defined, and the quantitative rock mass mechanical parameters corresponding to each bonding type zone are determined.

[0008] Based on the spatial morphology, the distribution of internal structural defects, the multiple bonding type partitions, and the quantified rock mass mechanical parameters, the stress background model of the mining site is configured to obtain a heterogeneous coupled digital twin model. The heterogeneous coupled digital twin model is used to characterize the spatial variation mechanical properties of the bonded ore pile.

[0009] When simulating in the heterogeneous coupled digital twin model, based on the multiple bonding type partitions and the quantified rock mass mechanical parameters, multi-field coupled simulation is performed on at least one candidate recovery operation scheme. With ore recovery index and safety risk index as optimization objectives, the operation parameters of the candidate recovery operation scheme are iteratively optimized to generate the target recovery operation scheme.

[0010] During the execution of the target recovery operation plan, the real-time collected construction feedback data is dynamically compared and assimilated with the predicted data at the corresponding location in the heterogeneous coupled digital twin model, and the operation parameters that have not yet been executed in the target recovery operation plan are dynamically adjusted based on the model assimilation results. Attached Figure Description

[0011] To more clearly illustrate the technical solutions in the embodiments of this application, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments of this application. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0012] Figure 1 This is a flowchart of a method for managing the recovery of bonded ore in a mining area, provided in an embodiment of this application.

[0013] Figure 2 This is a flowchart of another method for managing the recovery of bonded ore in a mining area, provided in an embodiment of this application.

[0014] Figure 3 This is a schematic diagram of a mining site bonded ore recovery management system provided in an embodiment of this application. Detailed Implementation

[0015] To make the objectives, technical solutions, and advantages of this application clearer, the embodiments of this application will be described in further detail below with reference to the accompanying drawings.

[0016] In this application, the terms "first," "second," etc., are used to distinguish identical or similar items with essentially the same function. It should be understood that there is no logical or temporal dependency between "first," "second," and "nth," nor are there any restrictions on quantity or execution order.

[0017] It should be noted that the information (including but not limited to user device information, user personal information, etc.), data (including but not limited to data used for analysis, data stored, data displayed, etc.) and signals involved in this application are all authorized by the user or fully authorized by all parties, and the collection, use and processing of related data must comply with the relevant laws, regulations and standards of the relevant countries and regions.

[0018] Traditional operations for recovering bonded ore in stopes generally rely on experience-based blasting, treating the bonded ore body as a "black box" with a lack of understanding of its internal structure, bonding mechanisms, and strength distribution. This results in inconsistent blasting parameters and fluctuating performance. Furthermore, overloading can cause uncontrollable damage to the stope's bottom structure, posing significant safety hazards. In addition, the operational plans are independent and fixed, unable to be dynamically adjusted according to real-time geological conditions, leading to low efficiency, poor safety, and low resource recovery rates.

[0019] In response, this application proposes a method for managing the recovery of caking ore in a stope, see [link to relevant documentation]. Figure 1 Taking the server as the executing entity as an example, the following steps are included.

[0020] 101. By fusing the three-dimensional laser scanning data, internal physical detection data and core genetic diagnosis data of the bonded ore pile, the spatial morphology and internal structural defect distribution of the bonded ore pile are obtained.

[0021] 102. Based on the differences in chemical and mineral composition indicated by the genetic diagnostic data of the rock core, define multiple bonding type zones and determine the quantitative rock mass mechanical parameters corresponding to each bonding type zone.

[0022] 103. Based on spatial morphology, distribution of internal structural defects, multiple bonding type partitions, and quantified rock mass mechanical parameters, the stress background model of the mining site is configured to obtain a heterogeneous coupled digital twin model. This heterogeneous coupled digital twin model is used to characterize the spatial variation mechanical properties of bonded ore piles.

[0023] 104. When performing simulations in a heterogeneous coupled digital twin model, based on multiple bonding type partitions and quantified rock mass mechanical parameters, multi-field coupled simulations are performed on at least one candidate recovery operation scheme. The operation parameters of the candidate recovery operation scheme are iteratively optimized with ore recovery index and safety risk index as optimization objectives to generate the target recovery operation scheme.

[0024] 105. During the execution of the target recovery operation plan, the real-time collected construction feedback data is dynamically compared and assimilated with the predicted data at the corresponding location in the heterogeneous coupled digital twin model, and the operation parameters that have not yet been executed in the target recovery operation plan are dynamically adjusted based on the model assimilation results.

[0025] For ease of understanding, the following explains some key terms in this embodiment:

[0026] Three-dimensional laser scanning data refers to point cloud data of the surface of a bonded ore pile obtained through laser scanning equipment, which is used to construct the external geometry of the ore pile.

[0027] Internal physical exploration data refers to the internal structural data of a bonded ore pile obtained through physical exploration methods such as ground-penetrating radar and microseismic monitoring, which is used to identify defects such as cracks and cavities.

[0028] Genetic diagnostic data from rock cores refers to the data obtained by analyzing the chemical composition, mineral composition, and microstructure of rock cores drilled from bonded ore deposits. This data is used to reveal the formation mechanism and mechanical properties of bonded ore.

[0029] Bonding type zoning refers to the division of bonded ore deposits into regions with different bonding characteristics based on the differences in chemical and mineral composition indicated by the genetic diagnostic data of the core.

[0030] Quantitative rock mechanics parameters refer to rock mechanical property parameters that are determined experimentally or empirically and correspond to each type of bonding zone, such as compressive strength, elastic modulus, Poisson's ratio, etc.

[0031] The mining area stress background model refers to a numerical model that characterizes the initial stress distribution in the mining area and provides boundary conditions for the mechanical behavior analysis of bonded ore piles.

[0032] Heterogeneous coupled digital twin model refers to a digital model that integrates the spatial morphology, internal structural defects, mechanical properties of heterogeneous materials, and stress background of the mining site of a bonded ore pile, and can characterize the mechanical response of the ore pile under different working conditions in real time.

[0033] Ore recovery indicators refer to quantitative indicators used to evaluate the effectiveness of ore recovery operations, such as the average size of fragments and the percentage of qualified fragments.

[0034] Safety risk indicators refer to quantitative indicators used to evaluate the safety of ore recovery operations, such as the maximum plastic zone depth and the peak displacement of key points.

[0035] Construction feedback data refers to construction parameters and environmental monitoring data collected in real time during the recovery operation, such as drilling parameters, vibration data, and displacement data.

[0036] Model assimilation refers to the process of comparing real-time collected construction feedback data with the prediction data of a digital twin model, and adjusting and updating the model parameters based on the comparison results to improve the model's accuracy.

[0037] This embodiment provides a method for the management of the recovery of bonded ore in a mining area.

[0038] The spatial morphology and internal structural defect distribution of the bonded ore pile are obtained by fusing 3D laser scanning data, internal physical exploration data, and core genetic diagnostic data. Specifically, the 3D laser scanning data can be used to obtain the surface morphology of the ore pile by directly converting point cloud data into a triangular mesh model, or by using a voxel-based method to discretize the point cloud data into a regular mesh before surface reconstruction. Internal physical exploration data can be used to identify macroscopic fractures, cavities, or stress concentration areas by manually interpreting ground-penetrating radar images or microseismic event distribution maps, or by using a simple threshold segmentation method to extract abnormal areas as defects from the exploration data. Core genetic diagnostic data can be used to identify the dominant bonding material type by manually analyzing the core's chemical composition and mineral composition reports, or by using simple statistical methods to classify the content of different components. These data can be fused by simply overlaying feature layers from different data sources, for example, by visualizing the surface model and internal defect point clouds in the same coordinate system to intuitively present the spatial morphology and internal structural defects.

[0039] Based on the differences in chemical and mineral composition indicated by the genetic diagnostic data of the core samples, multiple bonding type zones are defined, and quantitative rock mass mechanical parameters corresponding to each bonding type zone are determined. Specifically, bonding type zones can be roughly divided into regional divisions based on human experience and significant differences in chemical or mineral composition reported in the core analysis report, such as "siliceous bonding zone" and "carbonate bonding zone." Quantitative rock mass mechanical parameters can be obtained by consulting a general rock mechanics handbook and selecting representative mechanical parameters approximating the bonding type of rock as the quantitative rock mass mechanical parameters. For example, for the "siliceous bonding zone," the average mechanical parameters of quartz sandstone can be selected.

[0040] Furthermore, based on spatial morphology, distribution of internal structural defects, multiple bonding type zones, and quantified rock mass mechanics parameters, a heterogeneous coupled digital twin model is configured for the mining area stress background model, resulting in a heterogeneous coupled digital twin model. This model is used to characterize the spatially variable mechanical properties of the bonded ore pile. Specifically, the mining area stress background model can be obtained by directly embedding the geometric model of the bonded ore pile into a pre-established mining area geostress model and manually adjusting the boundary conditions to match the external contour of the ore pile. The heterogeneous material mechanical properties can be obtained by directly assigning the quantified rock mass mechanics parameters corresponding to each bonding type zone to the material properties of the corresponding area in the digital twin model, without considering the influence of internal defects. Coupling calculation and initialization can be performed in general finite element software by simply splicing the bonded ore pile model and the surrounding rock model, applying initial geostress, and performing a static equilibrium calculation to obtain the initial stress field and displacement field.

[0041] Based on this, during simulation in a heterogeneous coupled digital twin model, multi-field coupled simulations are performed on at least one candidate recovery operation scheme based on multiple bonding type partitions and quantified rock mass mechanical parameters. The operational parameters of the candidate recovery operation schemes are iteratively optimized with ore recovery indicators and safety risk indicators as optimization objectives to generate the target recovery operation scheme. Specifically, dynamic disturbance sources can be manually set in the digital twin model by manually setting one or more point or line force sources to simulate the instantaneous effects of blasting or hydraulic fracturing, and setting their action time and intensity. Multi-field coupled simulations can be performed independently in the simulation software, including mechanical and fluid simulations, and then the correlation between the results can be manually analyzed, such as observing the correspondence between the mechanical damage area and the fluid seepage path. Predicted indicator values ​​can be extracted by manually visually inspecting the simulation results (animations or images) to roughly determine the size distribution of fragments and the deformation of key areas, thereby estimating the ore recovery indicators and safety risk indicators. Iterative optimization involves manually adjusting the operational parameters, rerunning the simulation, and repeatedly trying different schemes based on the manually evaluated indicator values ​​until a satisfactory solution is found.

[0042] During the execution of the target recovery operation plan, real-time collected construction feedback data is dynamically compared and assimilated with the predicted data at corresponding locations in the heterogeneous coupled digital twin model. Based on the model assimilation results, the parameters of the unexecuted operations in the target recovery operation plan are dynamically adjusted. Specifically, dynamic comparison and model assimilation can be achieved by manually comparing the data collected by sensors at the construction site with the predicted values ​​of the digital twin model at the same location and time point. If significant differences are found, the material parameters of local areas in the model are manually modified to better match the actual observation values. Dynamic adjustment of operation parameters can be achieved by manually re-evaluating the parameters of the unexecuted operation units based on the updated model state after model assimilation and manually modifying the operation plan.

[0043] This method integrates multi-source sensing data to construct a heterogeneous coupled digital twin model, achieving a comprehensive understanding of the internal structure, bonding causes, and mechanical properties of bonded ore piles. Through intelligent simulation optimization, it can generate target schemes that balance ore recovery efficiency and operational safety, avoiding the lack of basis and potential damage associated with traditional experience-based blasting. Simultaneously, by incorporating real-time construction feedback for model assimilation and dynamic scheme adjustment, it effectively improves the efficiency, safety, and resource recovery rate of bonded ore recovery operations.

[0044] This application further proposes a method for fusing three-dimensional laser scanning data, internal physical exploration data, and core genetic diagnostic data of a bonded ore pile to obtain the spatial morphology and internal structural defect distribution of the bonded ore pile. See [link to relevant documentation]. Figure 2 Taking the server as the executing entity as an example, the method includes:

[0045] 201. Based on the three-dimensional laser scanning data, a three-dimensional surface model is generated, which is used to characterize the external contour of the bonded ore pile.

[0046] 202. By fusing and analyzing the ground-penetrating radar data and the microseismic monitoring event sequence in the internal physical exploration data, the spatial distribution data of the internal defect elements of the bonded ore pile were obtained.

[0047] 203. Spatial registration and fusion of the three-dimensional surface model and the spatial distribution data of the internal defect elements are performed to obtain the spatial morphology and internal structural defect distribution of the bonded ore pile.

[0048] The three-dimensional laser scanning data refers to point cloud data acquired through laser scanning equipment to describe the geometry of an object's surface. This data typically contains a large number of discrete three-dimensional coordinate points, accurately reflecting the external morphology of the object being measured. Acquiring this three-dimensional laser scanning data can be achieved through methods including, but not limited to, scanning the mining area using a Terrestrial Laser Scanner (TLS) or using a drone equipped with a LiDAR (Light Detection and Ranging) system for aerial scanning. Generating the three-dimensional surface model involves converting the discrete three-dimensional laser scanning data into a continuous surface representation that can be used for geometric analysis and visualization. This three-dimensional surface model characterizes the external contour of the bonded ore pile, providing accurate geometric boundaries for subsequent modeling. The internal physical detection data refers to data acquired through non-invasive or minimally invasive detection techniques to reveal the internal structural characteristics of the bonded ore pile. This data provides information on defects such as fractures, cavities, and weak surfaces within the ore pile. This internal physical detection data can include ground-penetrating radar data, microseismic monitoring event sequences, acoustic emission data, resistivity tomography data, etc.

[0049] The ground-penetrating radar (GPR) data is obtained by transmitting and receiving electromagnetic waves using GPR, and by analyzing the reflection and attenuation characteristics of electromagnetic waves in different media to detect underground structures. This GPR data can effectively identify macroscopic fissures, cavities, and other media discontinuities within the ore pile. The microseismic monitoring event sequence refers to the real-time monitoring and recording of micro-fracture events (micro-seismic events) within the rock mass by deploying sensor arrays in or around the ore pile. This microseismic monitoring event sequence can reflect stress concentration areas, active weak points, and potential instability trends within the rock mass.

[0050] The fusion analysis of the ground-penetrating radar (GPR) data and the microseismic monitoring event sequence refers to the comprehensive processing and interpretation of internal detection data from these two different sources to obtain more comprehensive and accurate information on internal defects in the bonded ore pile. Fusion analysis methods may include: spatiotemporally correlating structural anomalies identified by GPR with the epicenter locations of microseismic events to determine active structural defects; or using techniques such as data overlay, feature extraction, and pattern recognition to comprehensively determine the type, scale, and activity of internal defects. Through fusion analysis, spatial distribution data of internal defect elements in the bonded ore pile can be obtained, which details the location, geometry, and mechanical characteristics of defects such as cracks, cavities, and weak surfaces within the ore pile.

[0051] Spatial registration and fusion of the 3D surface model and the spatial distribution data of the internal defect elements refers to aligning and integrating the 3D surface model representing the external contour of the ore pile with the spatial distribution data representing the internal defects of the ore pile in a unified coordinate system. Spatial registration ensures an accurate spatial correspondence between the external geometry and the internal defects. Fusion integrates this information into a unified data structure, forming a complete digital model of the ore pile that includes both the external morphology and internal structural defects. After fusion, the spatial morphology and internal structural defect distribution of the bonded ore pile can be obtained, and this data forms the basis for subsequent construction of a heterogeneous coupled digital twin model.

[0052] Through the above technical solutions, this application can integrate multi-source heterogeneous data, solving the problem of how to comprehensively characterize the external contour and internal defects of bonded ore piles. Generating a three-dimensional surface model based on high-precision three-dimensional laser scanning data enables accurate capture of the external geometry of the ore pile, avoiding morphological deviations that may exist in traditional methods. By fusing and analyzing ground-penetrating radar detection data with microseismic monitoring event sequences, it is possible to complementaryly identify defect elements such as macroscopic cracks, cavities, active weak surfaces, and stress concentration zones within the ore pile, overcoming the limitations of single detection methods and obtaining more comprehensive and reliable spatial distribution data of internal defects. Spatial registration and fusion of the external three-dimensional surface model and the spatial distribution data of internal defect elements achieves the integration of external morphology and internal structural defects within a unified spatial framework, ensuring a high degree of consistency and integrity in the foundational data of the subsequent digital twin model. This comprehensive data fusion method provides a solid foundation for the configuration of subsequent heterogeneous coupled digital twin models, enhances the model's ability to represent the real state of bonded ore piles, and thus provides more accurate input for subsequent simulation and operation plan optimization. It avoids model prediction bias caused by insufficient or distorted data, thereby improving the efficiency and safety of ore recovery operations.

[0053] This application proposes a method for generating a three-dimensional surface model based on three-dimensional laser scanning data. The method includes: extracting geometric feature parameters from the three-dimensional laser scanning data, whereby these geometric feature parameters characterize the surface morphology of the bonded ore pile, including local curvature and surface roughness; dividing the surface of the bonded ore pile into different geometric feature partitions based on the distribution of these geometric feature parameters; and performing partition-adaptive surface reconstruction on the three-dimensional laser scanning data based on these different geometric feature partitions to obtain the three-dimensional surface model.

[0054] Geometric feature parameters, including local curvature and surface roughness, are extracted from the 3D laser scanning data to characterize the surface morphology of the bonded ore pile. These parameters are quantitative indicators describing the geometric properties of the bonded ore pile surface, reflecting the complexity and detail variations in its morphology. Local curvature describes the degree of bending of the surface at a given point, while surface roughness reflects the severity of microscopic undulations. The purpose of extracting these parameters is to provide a refined description of the bonded ore pile surface, offering a data foundation for subsequent zoning and adaptive reconstruction. Specifically, local curvature can be estimated by calculating the normal vector of each point in the point cloud and the rate of change of the normal vectors of its neighboring points. For example, principal component analysis (PCA) can be used to calculate the eigenvalues ​​of the local covariance matrix of the point cloud; the ratio between eigenvalues ​​reflects the local curvature. Alternatively, curvature parameters can be extracted by fitting a local surface (such as a quadratic surface) and then extracting the curvature parameters from the fitted surface. Surface roughness can be quantified by calculating the standard deviation of the distance from each point in the point cloud to its local fitted plane (such as a least-squares plane). A larger standard deviation indicates greater local undulation and higher roughness. Alternatively, a Fourier transform-based method can be used to analyze the frequency components of the surface height data; the more high-frequency components, the higher the roughness.

[0055] Based on the distribution of these geometric feature parameters, the surface of the bonded ore pile is divided into different geometric feature zones. Geometric feature partitioning refers to dividing the entire surface into several sub-regions with similar geometric characteristics based on the numerical range or distribution pattern of geometric feature parameters (such as local curvature and surface roughness) in different regions of the bonded ore pile surface. The purpose of this is to identify the heterogeneity of surface morphology, distinguishing between different types of regions such as high roughness, low roughness, flatness, and sharpness, providing a basis for subsequent targeted surface reconstruction. For example, clustering algorithms such as K-means or DBSCAN can be used to group points with similar local curvature and surface roughness values ​​into one class, forming different geometric feature zones. Before clustering, these parameters can be normalized. Alternatively, thresholds can be set; for example, regions with surface roughness above a certain threshold can be classified as high-roughness zones, regions below another threshold as low-roughness zones, and regions in between as medium-roughness zones. Furthermore, a more refined division can be achieved by combining a threshold for local curvature; for example, regions with both high curvature and high roughness may represent sharp edges.

[0056] Based on these different geometric feature partitions, partition-adaptive surface reconstruction is performed on the 3D laser scanning data to obtain the 3D surface model. Partition-adaptive surface reconstruction refers to using different surface reconstruction algorithms or parameter configurations for different geometric feature partitions to maximize the preservation of geometric details or achieve a smoothing effect in each partition. Its purpose is to overcome the limitations of traditional uniform reconstruction strategies, enabling the capture of sharp edges and complex textures in high-roughness regions, while effectively suppressing noise and generating smooth surfaces in low-roughness regions, thereby improving the accuracy and realism of the overall 3D surface model. For example, for high-roughness partitions, a detail-preserving reconstruction algorithm can be used, such as implicit function or voxel-based methods, with a fine parameter set (e.g., smaller voxel size, lower smoothness weight) to faithfully reconstruct sharp edges and complex textures. For low-roughness partitions, a smoothing and noise-suppression-focused reconstruction algorithm can be used, such as local surface fitting or mesh optimization methods, with a relatively relaxed parameter set (e.g., larger local fitting radius, higher smoothness weight) to effectively remove noise and generate a smooth surface. In addition, different mesh generation strategies can be dynamically selected based on the geometric characteristics of the partitions. For example, non-uniform meshing can be used for regions with complex structures, while uniform meshing can be used for flat regions to optimize computational efficiency and model accuracy.

[0057] By extracting geometric feature parameters such as local curvature and surface roughness from 3D laser scanning data, the surface morphology of bonded ore piles can be precisely and quantitatively described, providing an objective basis for subsequent processing. Based on the distribution of these geometric feature parameters, the surface of the bonded ore pile is divided into different geometric feature zones, thereby identifying the heterogeneity of the surface morphology and distinguishing between different regions such as high roughness and low roughness. On this basis, zone-adaptive surface reconstruction is performed for different geometric feature zones, employing the most suitable reconstruction algorithm and parameter configuration according to the specific geometric characteristics of each zone. For example, in high-roughness regions, a reconstruction strategy that preserves sharp edge features can be used to avoid loss of detail. In low-roughness regions, a reconstruction strategy that effectively suppresses noise and achieves a smoothing effect can be used to improve the smoothness of the model. This zone-adaptive reconstruction strategy enables the generated 3D surface model to reflect the complex surface morphology of the bonded ore pile with high fidelity, improving the accuracy and realism of the model. This is crucial for the subsequent spatial registration and fusion of the three-dimensional surface model and the spatial distribution data of internal defect elements, as well as for obtaining the accuracy of the spatial morphology and internal structural defect distribution of the bonded ore pile. This lays a solid foundation for the configuration of the mining site stress background model and the construction of the heterogeneous coupled digital twin model, thereby improving the reliability and effectiveness of the entire mining site bonded ore recovery management method.

[0058] This application further proposes a method for performing partition-adaptive surface reconstruction on the 3D laser scanning data based on the different geometric feature partitions to obtain the 3D surface model. Specifically, it includes: configuring a first reconstruction strategy for the high-roughness partitions within the geometric feature partitions. This first reconstruction strategy employs a Poisson-based surface reconstruction algorithm and a first parameter set, which includes a first voxel size and a first smoothing coefficient. The Poisson-based surface reconstruction algorithm is configured to generate a surface model that retains sharp edge features under the first parameter set. For the low-roughness partitions within the geometric feature partitions, a second reconstruction strategy is configured. This second reconstruction strategy employs a rolling sphere-based surface reconstruction algorithm and a second parameter set, which includes a second voxel size and a second smoothing coefficient. The second voxel size is larger than the first voxel size, and the second smoothing coefficient is larger than the first smoothing coefficient. Based on the first and second reconstruction strategies, point cloud data within the corresponding geometric feature partitions of the 3D laser scanning data are calculated to generate and fuse surface models for each partition to obtain the 3D surface model.

[0059] Specifically, configuring a first reconstruction strategy for the high-roughness partition within the geometric feature partition refers to identifying areas with significant surface undulations and rich details based on pre-set geometric feature parameter thresholds, and assigning a specific reconstruction method to these areas. This configuration process can be achieved by associating a specific processing function or algorithm identifier with the data blocks of the high-roughness partition in the data processing flow, or by allowing users to manually or automatically select a reconstruction mode suitable for high-roughness regions in the reconstruction software interface.

[0060] The first reconstruction strategy employs a Poisson-based surface reconstruction algorithm and a first parameter set. The Poisson-based surface reconstruction algorithm is a method for reconstructing 3D surfaces from point cloud data. Its core idea is to infer the implicit function that best represents the input point cloud by solving the Poisson equation, thereby generating a closed and smooth surface. This algorithm is particularly adept at handling noisy and incomplete point cloud data and effectively preserves surface details. The first parameter set is a set of configuration parameters tailored for this Poisson reconstruction algorithm, designed to optimize its reconstruction performance in high-roughness regions.

[0061] The first parameter set includes a first voxel size and a first smoothing coefficient. The first voxel size defines the minimum spatial unit size used by the Poisson reconstruction algorithm when constructing the octree; a smaller voxel size can capture finer geometric details. The first smoothing coefficient controls the degree of surface smoothing during reconstruction; a smaller smoothing coefficient means less smoothing, helping to preserve sharp features and rough textures in the original point cloud. For example, the first voxel size can be set to 0.01 meters, and the first smoothing coefficient can be set to 0.1 to enable fine characterization of high-roughness regions.

[0062] This Poisson-based surface reconstruction algorithm is configured to generate surface models that preserve sharp edge features under this first parameter set. This means that when applying the Poisson reconstruction algorithm, its internal parameters and optimization objectives are adjusted to prioritize the accuracy of edges and details. In addition to the first voxel size and the first smoothing coefficient, parameters such as the maximum depth of the octree and the minimum number of sample points can be adjusted to enable the algorithm to more sensitively identify and reconstruct sharp features such as rock fissures and edges when dealing with high-roughness regions, avoiding the loss of details due to over-smoothing.

[0063] Configuring a second reconstruction strategy for the low-roughness partition within the geometric feature partition involves identifying areas with relatively flat surfaces and minimal variations, and assigning a different, specialized reconstruction method to these areas. This configuration process can be achieved by associating a different processing function or algorithm identifier with the data blocks of the low-roughness partition during the data processing flow, or by allowing the user to manually or automatically select a reconstruction mode suitable for the low-roughness region within the reconstruction software interface.

[0064] The second reconstruction strategy employs a rolling sphere-based surface reconstruction algorithm and a second parameter set. The rolling sphere-based surface reconstruction algorithm connects adjacent points and constructs a triangular mesh by simulating a virtual sphere rolling across a point cloud surface. This algorithm performs well on relatively smooth surfaces, generating continuous and topologically correct meshes. The second parameter set consists of configuration parameters customized for this rolling sphere algorithm, designed to optimize its reconstruction performance in low-roughness regions. For example, the mesh density and smoothness are controlled by adjusting the rolling sphere radius and the point cloud density threshold.

[0065] The second parameter set includes a second voxel size and a second smoothing coefficient, where the second voxel size is larger than the first voxel size, and the second smoothing coefficient is larger than the first smoothing coefficient. The second voxel size (e.g., the rolling ball radius) defines the maximum distance considered by the rolling ball algorithm when connecting points. A larger voxel size helps generate a smoother, more continuous surface and effectively filters out local noise. The second smoothing coefficient controls the smoothness of the reconstructed surface; a larger smoothing coefficient means stronger smoothing, eliminating minor undulations in low-roughness areas and making them appear smoother. For example, the second voxel size can be set to 0.15 meters, and the second smoothing coefficient can be set to 0.7 to achieve smooth reconstruction of low-roughness areas.

[0066] Based on the first and second reconstruction strategies, point cloud data within corresponding geometric feature partitions of the 3D laser scanning data are calculated to generate and fuse surface models of each partition to obtain the 3D surface model. This step divides the original 3D laser scanning data according to its geometric feature partitions, and then independently inputs the data of each partition into its corresponding reconstruction strategy for calculation, generating its own local surface model. Through spatial registration and topological stitching techniques, these local surface models are seamlessly stitched together to form a complete and unified 3D surface model of the bonded ore pile.

[0067] Through the above technical solutions, for high-roughness partitions, an algorithm based on Poisson reconstruction combined with a first parameter set of small voxel size and small smoothing coefficient can accurately capture and preserve sharp edge features, avoiding the loss of details. For low-roughness partitions, an algorithm based on the rolling sphere method combined with a second parameter set of larger voxel size and larger smoothing coefficient is used to optimize the smooth surface, effectively filtering noise and avoiding the generation of excessive details. Through this differentiated strategy, this application can more accurately characterize the external contour of the bonded ore pile, improving the accuracy and realism of the 3D surface model, thus laying a solid foundation for the accurate acquisition of the spatial morphology and internal structural defect distribution of the bonded ore pile, and further improving the reliability of the heterogeneous coupled digital twin model.

[0068] This application further proposes a method for fusing and analyzing ground-penetrating radar (GPR) data and microseismic monitoring event sequences from internal physical exploration data to obtain spatial distribution data of internal defect elements within the bonded ore pile. This process includes: extracting reflected wave amplitude and electromagnetic wave attenuation characteristics from the GPR data, and identifying a first type of defect element based on these characteristics. This first type of defect element includes macroscopic cracks and cavities. From the microseismic monitoring event sequence, extracting the source location and energy release characteristics of the microseismic events, and identifying a second type of defect element based on these characteristics. This second type of defect element includes active weak surfaces and stress concentration zones. A spatial coupling relationship model is established between the first and second types of defect elements. Based on this spatial coupling relationship model, the first and second types of defect elements are fused and interpolated to generate spatial distribution data of internal defect elements within the bonded ore pile.

[0069] Specifically, ground-penetrating radar (GPR) data refers to data acquired by transmitting high-frequency electromagnetic waves underground using a ground-penetrating radar (GPR) and receiving the reflected echoes. This data records information about the propagation, reflection, and attenuation of electromagnetic waves in the underground medium, reflecting electrical differences in the underground medium and thus revealing the spatial distribution of underground structures, defects, and anomalies. Acquisition methods typically involve placing a GPR antenna on the surface of a bonded ore pile or inside a borehole, scanning along a pre-defined survey line, and constructing two-dimensional or three-dimensional detection profiles by recording echo signals at different depths and locations. The reflected wave amplitude refers to the intensity of the reflected electromagnetic wave signal received by the ground-penetrating radar; its magnitude is closely related to the differences in dielectric constant and conductivity of the underground medium, as well as the geometry of the reflecting interface. Electromagnetic wave attenuation characteristics refer to the degree of energy loss during electromagnetic wave propagation in the medium, mainly affected by factors such as the medium's conductivity, dielectric constant, and frequency. For example, when an electromagnetic wave encounters an interface with a large difference in dielectric constant or conductivity, a strong reflected wave amplitude is generated. In high-conductivity media, electromagnetic wave attenuation increases significantly. These characteristics are important physical bases for identifying underground defects. The first type of defect element, namely macroscopic fissures and cavities, refers to large-scale, typically static, structural discontinuities or spaces within the bonded ore pile. Identifying these defects based on reflected wave amplitude and electromagnetic wave attenuation characteristics is usually achieved by analyzing anomalous strong reflection signals, reflected waveform distortion, signal interruptions, or regions of rapid energy attenuation in the GPR profile. For example, macroscopic fissures may manifest as highly reflective interfaces with poor continuity, while cavities may appear as highly reflective anomalies with specific morphologies or as signal attenuation regions behind them. Identification methods can employ threshold-based image processing techniques, such as setting amplitude thresholds to delineate highly reflective areas, or analyzing signal energy attenuation rates to identify high-absorption regions. Alternatively, pattern recognition algorithms can be used to automatically identify specific reflection patterns associated with fissures and cavities in GPR images through model training.

[0070] A microseismic monitoring event sequence refers to a collection of information on the time, space (source location), and energy (magnitude or energy release) of microseismic events obtained through real-time recording of elastic wave signals caused by minor fractures or stress adjustments within the rock mass, obtained by a network of seismic sensors deployed around the mining area or bonded ore pile. These microseismic events are typically direct manifestations of local instability, crack propagation, or frictional slippage within the rock mass under stress. Acquisition methods usually involve deploying high-sensitivity accelerometers or velocity sensor arrays, picking up, locating, and calculating the magnitude of the waveform data received by the sensors to form a series of discrete microseismic event records. The source location refers to the spatial coordinates of the microseismic event, calculated by locating the waveforms received by multiple sensors based on time differences. Energy release characteristics are usually characterized by the magnitude of the microseismic event or directly calculated energy, reflecting the intensity and scale of the microseismic event. The source location directly indicates the specific area within the rock mass where fractures or stress concentrations occur, while the energy release characteristics reflect the activity level and potential destructive power of these areas. For example, high-energy microseismic events often indicate significant rock mass fracturing or substantial stress adjustment. The second type of defect element, namely active weak surfaces and stress concentration zones, refers to structural surfaces within a bonded ore pile that are currently under active stress and may undergo further deformation or failure, as well as potentially unstable areas with significantly higher stress than surrounding regions. Identifying these defects based on source location and energy release characteristics is typically achieved through spatiotemporal clustering analysis and energy density mapping of microseismic event sequences. For example, high-density clusters of microseismic events within a specific region may indicate active weak surfaces, while concentrated areas of high-energy microseismic events may correspond to stress concentration zones. Identification methods can include density-based clustering algorithms (such as DBSCAN) to identify clusters of microseismic events, or methods such as kernel density estimation to transform discrete source location and energy information into a continuous stress concentration probability map.

[0071] A spatial coupling model is a mathematical or conceptual model used to describe the spatial interrelationships, influences, or symbiotic relationships between first-type and second-type defect elements. The model aims to understand how static structural defects affect dynamic stress response and how dynamic stress activity indicates or activates static defects. For example, the presence of macroscopic cracks may lead to stress concentration, thereby inducing microseismic activity. Conversely, microseismic activity may occur along existing weak surfaces, making them active weak surfaces. This model can be built using various methods, such as constructing an empirical model by statistically analyzing the spatial co-occurrence, distance correlation, or geometric intersection of the two types of defect elements. Alternatively, based on rock mechanics principles, numerical simulations (such as finite element or discrete element methods) can be used to simulate the influence of cracks on the stress field, thereby establishing a physical mechanism model. Fusion refers to integrating information from different data sources into a unified data structure or spatial model to form a more comprehensive and accurate description of defects within a bonded ore pile. Interpolation refers to estimating data values ​​at unknown locations using mathematical methods based on known discrete data points, thereby generating continuous, high-resolution spatial distribution data. For example, geostatistical methods such as Kriging interpolation, inverse distance weighted interpolation, or radial basis function interpolation can be used to transform discrete defect element information into a continuous defect intensity field or defect probability field. The fusion process may include unifying the spatial location and attribute information of the two types of defect elements into a three-dimensional mesh model, and then using interpolation algorithms to fill in areas of the mesh that were not directly detected, thereby generating continuous spatial distribution data of internal defect elements. Spatial distribution data of internal defect elements in bonded ore piles refers to a comprehensive dataset based on three-dimensional spatial coordinates, describing the location, geometry, size, type, activity, or intensity of various defects within the bonded ore pile. This data is typically presented in the form of a three-dimensional mesh model, voxel model, or point cloud attribute model, and can intuitively and quantitatively reflect the heterogeneity and complexity of the internal structure of bonded ore piles.

[0072] By extracting the amplitude of reflected waves and the attenuation characteristics of electromagnetic waves from ground-penetrating radar data, the above-mentioned technical solutions can accurately identify and quantify static structural defects such as macroscopic cracks and cavities within the bonded ore pile, providing a solid geometric foundation for internal structural analysis. Simultaneously, extracting the source location and energy release characteristics from microseismic monitoring event sequences allows for the dynamic capture of active weak surfaces and stress concentration zones, revealing the stress state and potential instability areas within the bonded ore pile. Crucially, by establishing a spatial coupling relationship model between the first and second types of defect elements, this application can gain a deeper understanding of the intrinsic connection between static defects and dynamic stress activity, avoiding information conflicts and inconsistencies that may result from simple superposition. Based on this, the fusion and interpolation of the two types of defect elements using this spatial coupling relationship model not only integrates static and dynamic information but also generates continuous and complete spatial distribution data of defect elements within the bonded ore pile through interpolation techniques. This overcomes the limitations of single detection methods and provides a more comprehensive and coordinated characterization of structural defects within the bonded ore pile. This internal defect distribution data provides high-quality input for the subsequent construction of heterogeneous coupled digital twin models, improving the digital twin model's ability to characterize the spatial variation mechanical properties of bonded ore piles and the accuracy of simulation deduction, thereby providing a more reliable basis for optimizing the mining site bonded ore recovery operation scheme.

[0073] This application further proposes a method for spatially registering and fusing the three-dimensional surface model with the spatial distribution data of the internal defect elements to obtain the spatial morphology and internal structural defect distribution of the bonded ore pile. Specifically, the method includes: projecting the geometric center of each defect element in the spatial distribution data of the internal defect elements onto the three-dimensional surface model to obtain the surface influence area corresponding to each defect element on the three-dimensional surface model; determining the mechanical influence intensity value of the defect element on the corresponding surface influence area based on the type and size of each defect element, and assigning the mechanical influence intensity value as an attribute to the surface influence area; and generating and drawing a surface defect influence intensity distribution map on the three-dimensional surface model based on the surface influence areas corresponding to all defect elements and the mechanical influence intensity values ​​of the surface influence areas. This surface defect influence intensity distribution map is used to represent the spatial morphology and internal structural defect distribution of the bonded ore pile.

[0074] This process involves projecting the geometric center of each defect element in the spatial distribution data of the internal defect elements onto the 3D surface model to obtain the surface influence region corresponding to each defect element on the 3D surface model. The aim is to establish the spatial relationship between internal defects and the external surface. This step is crucial for understanding how internal structural defects potentially affect the external behavior of bonded ore piles. Specifically, a ray casting method can be used, in which rays are emitted from the geometric center of each defect element along a predetermined direction (e.g., perpendicular to the surface or along the principal stress direction), and the intersection of the ray and the 3D surface model is the center of the surface influence region. Alternatively, a nearest-point search method can be used, in which the Euclidean distance from the geometric center of each defect element to all points on the 3D surface model is calculated, and the nearest point or group of points is determined as the center of the surface influence region.

[0075] Based on the type and size of each defect element, the mechanical influence intensity value of that defect element on the corresponding surface influence area is determined, and this mechanical influence intensity value is assigned as an attribute to the surface influence area, aiming to quantitatively assess the potential mechanical influence of internal defects. This step combines the geometric information and mechanical properties of the defect, providing more accurate input for subsequent simulation analysis. In practice, an empirical mapping relationship can be established by consulting a pre-defined defect type-size-intensity lookup table or empirical formula to map different types (such as macroscopic cracks, cavities, and active weak surfaces) and sizes (such as length, area, and volume) of defects to corresponding mechanical influence intensity values. Alternatively, a simplified mechanical model can be used for estimation. For example, for crack-type defects, the stress concentration factor can be calculated based on fracture mechanics theory and normalized to a mechanical influence intensity value. The determined mechanical influence intensity value is then used as an additional attribute and bound to the corresponding surface influence area.

[0076] Based on the surface influence areas corresponding to all defect elements and the mechanical influence intensity values ​​of these areas, a surface defect influence intensity distribution map is generated and plotted on the 3D surface model. This map represents the spatial morphology and internal structural defect distribution of the bonded ore pile, aiming to provide an intuitive and comprehensive visualization clearly showing the combined influence of internal defects on the surface mechanical behavior of the bonded ore pile. Specifically, interpolation algorithms can be used, such as spatial interpolation (e.g., inverse distance weighted interpolation or kriging interpolation) of the mechanical influence intensity values ​​of all surface influence areas to generate a continuous intensity field, which is then plotted on the 3D surface model as a color gradient map. Alternatively, a discrete region coloring method can be used, assigning a specific color or texture to each surface influence area based on its mechanical influence intensity value and rendering it directly on the 3D surface model. Overlapping areas can be handled by intensity superposition or taking the maximum value.

[0077] By projecting the geometric center of internal defect elements onto a 3D surface model and obtaining the corresponding surface influence area, this application establishes a direct and spatial relationship between internal defects and the external surface, avoiding the limitations of treating the internal structure as a "black box" in traditional methods. Furthermore, by determining the mechanical influence intensity value based on the type and size of the defect and assigning it as an attribute to the surface influence area, a quantitative assessment of the mechanical influence of internal defects is achieved, enabling the model to more realistically reflect the actual mechanical behavior of the bonded ore pile. By generating and plotting the surface defect influence intensity distribution map, this application provides an intuitive and comprehensive representation of the spatial morphology and internal structural defect distribution of the bonded ore pile, providing high-precision and high-reliability basic data for the subsequent configuration of the mining site stress background model and the construction of the heterogeneous coupled digital twin model. This improves the accuracy of subsequent recovery operation simulation and optimization, thereby effectively reducing potential safety risks and improving resource recovery efficiency.

[0078] This application further proposes defining multiple bonding type zones based on the differences in chemical and mineral composition indicated by the genetic diagnostic data of the core, and determining the quantitative rock mass mechanical parameters corresponding to each bonding type zone. The scheme includes the following steps:

[0079] The genetic diagnostic data from this core was interpreted to distinguish at least two dominant chemical or mineralogical modes. This genetic diagnostic data refers to data obtained through laboratory analysis of core samples drilled from bonded ore deposits. It includes information on the ore's chemical composition, mineral composition, and microstructure, reflecting the formation mechanism and cementation characteristics of the bonded ore. The interpretation of this data aims to identify the main chemical reactions or mineralogical cementation types leading to ore bonding. For example, mineral phase composition can be determined through X-ray diffraction (XRD), and the microstructure and elemental distribution of the cement can be observed through scanning electron microscopy (SEM) combined with energy dispersive spectroscopy (EDS). Distinguishing the dominant chemical or mineralogical modes can be achieved using statistical analysis methods, such as principal component analysis (PCA) or cluster analysis, to identify representative chemical or mineralogical assemblages from multidimensional data, thereby summarizing different bonding genetic types. For example, one mode might be sulfate cementation formed by the oxidation of high sulfides, while another mode might be carbonate cementation formed by the precipitation of carbonate minerals. In addition, experts can use pre-defined mineralogical and geochemical discrimination rules to logically deduce different dominant genetic patterns from the diagnostic data.

[0080] Based on the characteristic differences reflected in the genetic diagnostic data by at least two dominant chemical or mineralogical genetic models, the bonded ore deposit is divided into multiple bond type partitions corresponding to each of the at least two dominant genetic models. These characteristic differences refer to the distinguishable quantitative indicators exhibited by different dominant genetic models in diagnostic data such as chemical composition, mineral composition, and microstructure. For example, a sulfate cementation model may be characterized by high contents of sulfur, calcium, and gypsum minerals, while a carbonate cementation model may be characterized by high contents of carbon, calcium, magnesium, and calcite and dolomite minerals. Based on these differences, geostatistical methods, such as kriging interpolation or co-kriging interpolation, can be used to extrapolate the diagnostic data of core points to the entire bonded ore deposit space, forming a continuous characteristic difference distribution map. By setting thresholds or applying classification algorithms (such as support vector machines or random forests), the bonded ore deposit space is divided into several discrete bond type partitions according to these characteristic differences, with each partition corresponding to a dominant chemical or mineralogical genetic model.

[0081] For each bonding type zone, based on the dominant chemical or mineralogical genetic model corresponding to that bonding type zone, the corresponding quantitative rock mass mechanical parameters are retrieved from a pre-stored mapping table of genetic mechanical parameters. This pre-stored mapping table is a database or knowledge base that stores the correspondence between different chemical or mineralogical genetic models and their corresponding quantitative rock mass mechanical parameters. These relationships are obtained through the accumulation of extensive laboratory test data (such as uniaxial compressive strength, Brazilian splitting strength, shear strength, etc.) and field test data, and are formed through statistical analysis and empirical summarization. For example, this mapping table might contain a set of mechanical parameter ranges corresponding to "high sulfide oxidative cementation" (e.g., compressive strength 20-50 MPa, elastic modulus 10-30 GPa), while "carbonate cementation" corresponds to another set of mechanical parameter ranges (e.g., compressive strength 50-100 MPa, elastic modulus 30-60 GPa). Once the dominant genetic pattern of a certain bonding type partition is determined, the mapping relationship table can be queried according to the pattern to obtain the corresponding quantified rock mechanics parameters, such as compressive strength, elastic modulus, Poisson's ratio, cohesion, and internal friction angle, thereby providing accurate material property input for subsequent mechanical models.

[0082] By employing the aforementioned technical solutions and deeply interpreting the genetic diagnostic data, at least two dominant chemical or mineralogical genetic modes can be accurately distinguished. This provides a comprehensive understanding of the formation mechanism of bonded ore, laying the foundation for subsequent zoning. Zoning based on the characteristic differences of these dominant genetic modes allows for the precise capture and characterization of the spatial variability of the mechanical properties of the ore pile, avoiding model errors caused by ambiguous zoning. Furthermore, by querying quantified rock mass mechanical parameters from a pre-stored mapping table of genetic mechanical parameters, the objectivity and consistency of parameter acquisition are ensured, eliminating biases from human estimation and providing reliable and matching mechanical input for the heterogeneous coupled digital twin model. This characterization of mechanical properties improves the accuracy of subsequent configuration of the mining site stress background model, enabling the heterogeneous coupled digital twin model to more realistically reflect the actual mechanical behavior of bonded ore piles. This provides a more solid foundation for multi-field coupled simulation of candidate recovery operation schemes, ultimately optimizing ore recovery indicators and safety risk indicators, and improving the efficiency and safety of recovery operations.

[0083] This application further proposes a method for managing the recovery of bonded ore in a stope. Based on spatial morphology, distribution of internal structural defects, multiple bond type zones, and quantified rock mechanics parameters, a stope stress background model is configured to obtain a heterogeneous coupled digital twin model. Specifically, the method includes the following steps: coupling spatial morphology with the distribution of internal structural defects to construct the geometric and structural entity boundaries of the bonded ore pile within the stope stress background model; spatially assigning values ​​to multiple bond type zones and quantified rock mechanics parameters to define heterogeneous material mechanical properties for different regions within the geometric and structural entity boundaries; and performing coupled calculations and initialization on the stope stress background model based on the geometric and structural entity boundaries and the heterogeneous material mechanical properties to obtain a heterogeneous coupled digital twin model for simulation.

[0084] Specifically, this method couples spatial morphology with the distribution of internal structural defects, constructing the geometric and structural entity boundaries of the bonded ore pile within the stress background model of the mining site. The aim is to integrate the external geometry and internal structural weaknesses of the bonded ore pile into the model. One implementation involves integrating a three-dimensional surface model obtained from processing three-dimensional laser scanning data with spatial distribution data of internal defect elements obtained from analyzing internal physical detection data.

[0085] This method spatially assigns values ​​to multiple bonding type partitions and quantified rock mass mechanical parameters, defining heterogeneous material mechanical properties for different regions within the geometric and structural entity boundaries. The aim is to ensure that the mechanical parameters of each region in the model can truly reflect the variability of mechanical properties caused by differences in bonding causes.

[0086] Based on the geometric and structural entity boundaries and the material mechanical properties of heterogeneous materials, this method performs coupled calculations and initialization on the stress background model of the mining site to obtain a heterogeneous coupled digital twin model for simulation. The aim is to integrate the constructed bonded ore pile model with the stress background model of the mining site and bring it to a stable initial mechanical equilibrium state, providing a reliable starting point for subsequent dynamic simulation.

[0087] Through the above technical solution, this application can finely couple the external geometry (spatial morphology) of the bonded ore pile with its internal structural defects such as cracks and cavities (distribution of internal structural defects), constructing a highly realistic geometric and structural entity boundary in the mining site stress background model. This boundary not only reflects the macroscopic outline of the ore pile, but more importantly, it directly embeds the internal structural weaknesses and discontinuities into the model geometry, thus avoiding the simplification of treating the ore pile as a homogeneous continuum in traditional modeling. This improves the model's accuracy in representing the true structural characteristics of the ore pile, providing a more reliable geometric basis for subsequent mechanical simulations and avoiding simulation distortion and errors caused by model simplification. Furthermore, by spatially assigning multiple bonding type partitions defined based on core genetic diagnostic data to the corresponding quantified rock mass mechanical parameters of each partition, different regions within the aforementioned geometric and structural entity boundary can be endowed with heterogeneous material mechanical properties consistent with their true bonding genesis and mechanical characteristics. This refined property definition overcomes the limitations of homogeneous material properties in traditional models, enabling the model to realistically reflect the spatial variability of mechanical properties within the bonded ore pile due to differences in its formation. This allows simulation results to more accurately predict the actual response of the ore pile under disturbance, further reducing simulation errors. Based on the constructed geometric and structural entity boundaries and the finely defined heterogeneous material mechanical properties, coupled calculations and initialization are performed on the mining site stress background model. This process fully couples the bonded ore pile model with the initial geostress field of the surrounding rock, achieving mechanical equilibrium through iterative calculations, resulting in a stable and physically consistent heterogeneous coupled digital twin model. This heterogeneous coupled digital twin model not only integrates the complex geometry, internal defects, and heterogeneous mechanical properties of the ore pile, but its initial state also fully considers the constraints of the geostress background, providing highly reliable and accurate initial conditions for subsequent multi-field coupled simulations. This significantly improves the predictive power of the simulation results and the effectiveness of guiding recovery operation plans.

[0088] This application further proposes coupling spatial morphology with the distribution of internal structural defects to construct the geometric and structural entity boundaries of the bonded ore pile in the mining site stress background model. Specifically, this includes: converting the spatial morphology into an initial geometric entity in the mining site stress background model; defining key defect elements in the distribution of internal structural defects as internal boundary conditions and attribute constraints of the initial geometric entity; and reconstructing and generating the geometric and structural entity boundaries in the mining site stress background model based on the internal boundary conditions and attribute constraints, which embed the distribution of internal structural defects.

[0089] The process of converting the spatial morphology into initial geometric entities in the mining site stress background model aims to transform the external geometry of the bonded ore pile, obtained through methods such as 3D laser scanning, into initial 3D entities that the mining site stress background model can recognize and process. This is the foundation for constructing a digital twin model, enabling the model to accurately reflect the actual external contours of the bonded ore pile. Specifically, this can be achieved in various ways. For example, the initial geometric entity can be generated by directly importing the 3D surface model (such as STL, OBJ formats) into the finite element or discrete element modeling software used for the mining site stress background model and performing mesh generation. Alternatively, the 3D surface model can be parametrically fitted, for example, using NURBS surfaces or implicit surfaces to describe the external morphology of the bonded ore pile, and then these parametric descriptions can be converted into geometric entities in the mining site stress background model.

[0090] Furthermore, key defect elements in the distribution of internal structural defects are defined as internal boundary conditions and property constraints of the initial geometric entity. This step aims to transform the key defect information such as cracks, cavities, and weak surfaces inside the bonded ore pile obtained through internal physical exploration data analysis into internal boundary conditions and material property constraints that can be understood and applied by the mining site stress background model. This is crucial for accurately simulating the mechanical behavior of bonded ore piles, as internal defects are key factors affecting their overall stability and failure mode. In practice, key defect elements (such as crack surfaces and cavity boundaries) can be directly defined as the internal geometric boundaries of the model, and specific contact properties (such as friction coefficient and adhesion force) or fracture criteria can be assigned to these boundaries. At the same time, the rock mass region surrounding the defects can be assigned specific material properties (such as reduced elastic modulus and tensile strength). Another approach is to embed key defect elements into the initial geometric entity through implicit functions or level set methods, and generate discontinuous surfaces or weakened regions inside the model according to the type and location of the defects. The material properties of these regions can be set to anisotropic or have damage evolution characteristics, thus serving as property constraints.

[0091] Building upon this foundation, based on internal boundary conditions and property constraints, geometric and structural entity boundaries are reconstructed and generated in the mining site stress background model. These geometric and structural entity boundaries embed the distribution of internal structural defects. This step aims to refine the initial geometric entity using the internal boundary conditions and property constraints defined in the previous step, thereby constructing geometric and structural entity boundaries containing internal defect information in the mining site stress background model. This allows the digital twin model to not only possess an external morphology but also accurately reflect internal heterogeneity and structural weaknesses, providing a high-fidelity foundation for subsequent mechanical simulations. In specific implementation, Boolean operations or cutting operations can be used to "cut out" or "embed" the geometric boundaries defined by internal defect elements from the initial geometric entity, forming new geometric entities with complex internal structures. Simultaneously, the corresponding material property constraints are mapped onto these newly generated internal structures. Alternatively, adaptive mesh generation techniques can be used to refine the mesh near the internal defect elements, and based on internal boundary conditions and property constraints, the heterogeneity of material properties or the mechanical behavior of contact interfaces can be directly defined on the mesh cells, thereby achieving defect embedding without explicitly changing the overall geometric topology.

[0092] The above technical solution transforms the spatial morphology of the bonded ore pile into an initial geometric entity within the mining site stress background model, laying the foundation for subsequent refined modeling. Furthermore, key defect elements in the internal structural defect distribution are defined as internal boundary conditions and attribute constraints of the initial geometric entity. This allows the model to accurately capture key structural features such as cracks, voids, and weak surfaces within the bonded ore pile, avoiding the problem of treating the bonded ore pile as a homogeneous "black box." Based on these internal boundary conditions and attribute constraints, geometric and structural entity boundaries embedded with the internal structural defect distribution are reconstructed and generated within the mining site stress background model, thus constructing a high-fidelity heterogeneous coupled digital twin model. This model not only reflects the external contour of the bonded ore pile but, more importantly, characterizes the influence of internal structural defects on mechanical behavior, providing a more realistic and reliable physical basis for subsequent multi-field coupled simulations, improving the accuracy of simulation results, and thereby optimizing the prediction and optimization effects of ore recovery indicators and safety risk indicators.

[0093] This application further proposes spatial assignment of multiple bond type partitions and quantified rock mass mechanical parameters, defining heterogeneous material mechanical properties for different regions within the geometric and structural entity boundaries, specifically including:

[0094] The spatial extent of multiple bonding type partitions is mapped to the boundaries of geometric and structural entities, and corresponding material property partitions are defined within the boundaries of geometric and structural entities.

[0095] For each material property partition, based on the bond type partition corresponding to the material property partition, the parameter value matching the bond type partition is called from the quantified rock mechanics parameters, and the parameter value is used as the material property value of the material property partition.

[0096] The material property values ​​of each material property partition are assigned to the corresponding spatial units in the mining site stress background model, thereby generating material mechanical properties within the geometric and structural entity boundaries.

[0097] This process involves mapping the spatial extent of multiple bond-type partitions to the boundaries of geometric and structural entities. Within these boundaries, corresponding material property partitions are then defined, aiming to establish a spatial correspondence between the bond-type partitions and the geometric representation of the ore pile in the digital twin model. This step ensures the accuracy of subsequent material property assignments and avoids model distortion caused by spatial misalignment. Specifically, a geometric topology mapping method can be used, employing Boolean or intersection operations to superimpose the geometric model of the bond-type partitions onto the entity boundaries, thereby defining the material property partitions. Alternatively, a spatial mesh discretization method can be used, dividing the space within the entity boundaries into mesh cells and labeling each cell according to its corresponding bond-type partition.

[0098] For each material property partition, based on the corresponding bond type partition, parameter values ​​matching the bond type partition are retrieved from the quantitative rock mechanics parameters. These parameter values ​​are then used as the material property values ​​for the material property partition. The purpose is to assign quantitative mechanical parameters consistent with the bond type to each spatially defined material property partition. This is crucial for accurately simulating the heterogeneous mechanical behavior of bonded ore piles. For example, a database or lookup table containing different bond types and their corresponding quantitative rock mechanics parameters can be pre-established. When processing each material property partition, the corresponding parameter values ​​can be directly retrieved from this table based on its bond type. Alternatively, a rule engine can be designed to dynamically infer and generate the corresponding quantitative rock mechanics parameter values ​​based on the bond type characteristics of the material property partition, combined with expert knowledge or machine learning models.

[0099] The material property values ​​for each material property partition are assigned to the corresponding spatial elements in the mining site stress background model, thereby generating material mechanical properties within the geometric and structural entity boundaries. This step integrates discrete material property values ​​into a continuous numerical simulation model, forming a material mechanical property field with spatial variability. This allows the digital twin model to realistically reflect the heterogeneous characteristics of bonded ore piles. Specifically, if the mining site stress background model uses the finite element method, the property values ​​of the material property partitions can be directly assigned to the finite element mesh elements they cover. If the discrete element method is used, the property values ​​can be assigned to the discrete element particles or particle clusters they contain. Alternatively, the model can be discretized into voxels, and the property values ​​can be assigned to the corresponding voxel elements.

[0100] The above technical solution maps the spatial range of the bonding type partition to the geometric and structural entity boundaries, and divides the corresponding material property partitions within them. This ensures that the partitioning strictly corresponds to the actual structure of the bonded ore pile, avoiding regional misalignment or omissions caused by empirical division. Furthermore, for each material property partition, matching parameter values ​​are retrieved from the quantified rock mechanics parameters based on its corresponding bonding type partition and used as material property values. This fully utilizes the correlation between bonding origin diagnostic data and mechanical parameters, effectively preventing parameter selection bias or type mismatch, ensuring that the attribute values ​​of each partition accurately reflect the mechanical behavior of the bonded ore. Assigning the attribute values ​​of each material property partition to the corresponding spatial units in the mining site stress background model achieves spatial positioning of the attribute values, ensuring a one-to-one match between model units and the actual ore pile area, thereby improving the realism of the heterogeneous coupled digital twin model and the reliability of simulation results. This definition of heterogeneous material properties provides a solid foundation for subsequent multi-field coupled simulations in the heterogeneous coupled digital twin model, enabling simulation results to more accurately predict ore recovery indicators and safety risk indicators, thus guiding the generation of more optimized target recovery operation schemes.

[0101] This application further proposes a method for performing coupled calculations and initialization on a mining site stress background model based on geometric and structural entity boundaries and heterogeneous material mechanical properties to obtain a heterogeneous coupled digital twin model for simulation. The method includes the following steps: setting the geometric and structural entity boundaries as new embedded geometric boundaries in the mining site stress background model; applying heterogeneous material mechanical properties to the model region defined by the new embedded geometric boundaries; establishing a contact relationship between the new embedded geometric boundaries and the surrounding rock region of the mining site stress background model, and performing iterative calculations based on the initial geostress conditions provided by the mining site stress background model and the contact relationship to solve the coupled stress balance within and outside the new embedded geometric boundaries, obtaining a converged mechanical equilibrium solution; setting the mechanical equilibrium solution as the initial simulation condition for the heterogeneous coupled digital twin model to obtain the heterogeneous coupled digital twin model.

[0102] Specifically, setting the geometric and structural entity boundaries as new embedded geometric boundaries in the mining site stress background model aims to integrate the geometric shape and internal structural defect information of the bonded ore pile into the model, thereby providing accurate physical boundaries for subsequent mechanical analysis. This can be achieved in several ways. For example, using the geometric modeling function of professional CAD / CAE software, a pre-generated 3D geometric model (including internal defects) can be imported into the mining site stress background model and marked as a new embedded geometry, typically involving Boolean operations or geometric merging operations. Alternatively, in the mining site stress background model, based on the spatial morphology and internal structural defect distribution of the bonded ore pile, the mesh of the corresponding area can be locally reconstructed and refined, so that the reconstructed mesh boundary fits the geometric and structural entity boundaries of the bonded ore pile.

[0103] Applying heterogeneous material mechanical properties to the model region defined by the newly embedded geometric boundary aims to accurately assign spatially variable mechanical parameters, determined based on the actual conditions of the bonded ore pile, to the corresponding region in the model, thus reflecting the true heterogeneous characteristics of the bonded ore pile. This can be achieved through attribute partitioning and database mapping. Based on multiple predefined bond type partitions, these partitions are mapped to corresponding spatial regions within the newly embedded geometric boundary. Then, from a pre-stored material database or parameter mapping table, the corresponding quantized rock mechanics parameters are automatically retrieved and assigned to the finite element elements or mesh nodes of that region according to the bond type of each partition. Alternatively, a script can be written to read the spatial coordinate information and corresponding mechanical parameters of the bond type partitions, and then programmatically traverse all elements or nodes within the newly embedded geometric boundary, dynamically assigning them the corresponding heterogeneous material mechanical properties according to their bond type partitions.

[0104] Based on this, a contact relationship is established between the newly embedded geometric boundary and the surrounding rock region in the mining site stress background model. This serves to simulate the actual mechanical interaction between the bonded ore pile and the surrounding rock mass, ensuring that stress transfer and deformation behavior are accurately reflected at the interface. This can be achieved by defining a surface contact algorithm. For example, in finite element analysis software, a contact pair can be defined between the outer surface of the newly embedded geometric boundary and the inner surface of the adjacent surrounding rock region in the mining site stress background model, and contact properties such as friction coefficient, normal contact stiffness, and tangential contact behavior can be configured. Alternatively, the interface element method can be used, inserting a layer of special interface elements between the newly embedded geometric boundary and the surrounding rock region. These interface elements can simulate the mechanical behavior of the contact surface, including normal compression, tension, shear, and friction, thereby achieving mechanical coupling between the two.

[0105] Furthermore, iterative calculations are performed based on the initial geostress conditions and contact relationships provided by the mining site stress background model to solve the coupled stress balance in the regions inside and outside the newly embedded geometric boundary, obtaining a convergent mechanical equilibrium solution. This step is crucial for ensuring the model reaches a stable and balanced state before simulation. It iteratively solves a complex set of mechanical equations, enabling the entire model to reach mechanical equilibrium under the initial geostress field. This can be achieved using a stepwise loading and equilibration method: an initial geostress field is applied to the mining site stress background model, then the properties of the newly embedded geometric boundary and heterogeneous materials are introduced into the model, and the boundary conditions are gradually loaded or adjusted while performing nonlinear iterative calculations until the entire model reaches mechanical equilibrium. Alternatively, an implicit solver and convergence criteria can be used. Iterative calculations are performed using an implicit finite element solver, and convergence criteria such as residual force, energy change rate, or displacement increment are monitored to determine whether coupled stress equilibrium has been reached.

[0106] The mechanical equilibrium solution is set as the initial simulation condition for the heterogeneous coupled digital twin model. This step uses the rigorously calculated and verified equilibrium state as the starting point for all subsequent dynamic simulation analyses, ensuring the stability of the simulation process and the reliability of the results. This can be achieved by storing the convergent mechanical equilibrium solution (including displacement, stress, and strain fields) obtained through iterative calculations as the initial state variables in the digital twin model's data structure. Alternatively, the mechanical equilibrium solution can be exported in a specific data format via file export and import, and then imported into the simulation model as input when starting subsequent dynamic simulation tasks.

[0107] By employing the aforementioned technical solution, the geometric and structural entity boundaries are set as newly embedded geometric boundaries, defining the physical boundary framework of the model. This ensures that the geometric shape and internal structural defects of the bonded ore pile are accurately incorporated, providing a solid geometric constraint for coupled calculations. The application of heterogeneous material mechanical properties assigns mechanical parameters reflecting the actual heterogeneous characteristics of the bonded ore pile, avoiding property distortion caused by model simplification and providing a realistic material basis for subsequent simulations. Establishing a contact relationship between the newly embedded geometric boundary and the surrounding rock region simulates the actual mechanical interaction between the ore pile and the surrounding rock, preventing interaction omission errors caused by boundary simplification. Iterative calculations based on initial geostress conditions and contact relationships fully consider the original stress state of the stope. By iteratively adjusting and solving the coupled stress balance, the stress in the internal and external regions of the model is coordinated, eliminating the divergence risk caused by unbalanced forces. Setting the convergent mechanical equilibrium solution as the initial condition for simulation provides a stable and realistic starting point for dynamic deduction, avoiding the failure of the simulation process due to the instability of the initial state. This improves the accuracy and reliability of the heterogeneous coupled digital twin model, provides a solid foundation for the optimization of subsequent recovery operation schemes, and thus enhances the safety and efficiency of recovery operations.

[0108] This application further proposes an iterative calculation based on the aforementioned initial geostress conditions and contact relationships to solve the coupled stress balance within and outside the newly embedded geometric boundary, thereby obtaining a convergent mechanical equilibrium solution. This process includes:

[0109] Based on the established contact relationship, the contact force and relative displacement between the newly embedded geometric boundary and the surrounding rock region are determined. This step aims to quantify the interaction between the newly embedded geometric boundary and the surrounding rock region. In numerical simulations, contact relationships are typically handled using specific algorithms to avoid unreasonable penetration or separation between model elements. Determining the contact force and relative displacement is a key input for establishing the mechanical equilibrium equations. One implementation method is the penalty function method. In this method, when penetration occurs between the newly embedded geometric boundary and the surrounding rock region, a virtual "penalty" force is introduced, which is proportional to the penetration amount, simulating the repulsive effect of the contact surface. The relative displacement can be obtained by calculating the penetration amount, and the penalty force is the contact force. Another implementation method is the Lagrange multiplier method. This method introduces Lagrange multipliers as unknowns, directly incorporating the contact constraints into the equilibrium equations for solution. The Lagrange multipliers physically represent the contact force, while the difference in nodal displacements on the contact surface directly reflects the relative displacement. In addition, the augmented Lagrange method can be used, which combines the advantages of the penalty function method and the Lagrange multiplier method. By iteratively updating the Lagrange multipliers and penalty terms, it can improve convergence efficiency while ensuring computational accuracy.

[0110] Based on the aforementioned contact forces and relative displacements, the displacement and stress fields of the regions inside and outside the newly embedded geometric boundary are iteratively adjusted to gradually reduce the unbalanced forces on the newly embedded geometric boundary. This step is the core iterative process for solving the coupled stress equilibrium. Its purpose is to achieve equilibrium of all forces on the newly embedded geometric boundary by repeatedly correcting the displacement and stress states at each point in the model, i.e., the unbalanced forces approach zero. One implementation method is the Newton-Raphson iterative method. In each iteration, the tangent stiffness matrix is ​​calculated based on the current unbalanced forces, and the displacement increment is solved, then the displacement and stress fields of the entire model are updated. This method has a fast convergence speed, but requires the calculation and assembly of the tangent stiffness matrix. Another implementation method is the explicit dynamic relaxation method. This method treats the unbalanced forces as virtual inertial forces and, by introducing damping terms, gradually approaches static equilibrium within the virtual time step. This method does not require assembling the global stiffness matrix and is suitable for large-scale problems, but the convergence speed may be relatively slow. Alternatively, the conjugate gradient method can be used. This method minimizes the energy functional by iteratively searching for conjugate directions, thereby gradually adjusting the displacement field and stress field until mechanical equilibrium is reached.

[0111] When both the non-equilibrium force norm and the rate of change of strain energy on the newly embedded geometric boundary are below a preset convergence threshold, the coupled stress equilibrium is determined to be reached, and the displacement and stress fields at this point are output as the converged mechanical equilibrium solution. This step defines the termination condition for the iterative calculation, ensuring that the obtained mechanical equilibrium solution is stable and sufficient. Employing a dual convergence criterion, i.e., simultaneously monitoring the non-equilibrium force norm and the rate of change of strain energy, can avoid false convergence or premature termination of iteration, thereby improving the reliability of the results. For the non-equilibrium force norm, the preset convergence threshold can be set as a percentage (e.g., 0.01% to 1%) of the sum of the non-equilibrium force vector magnitudes of all nodes in the model relative to the total external or internal forces, or an absolute small value (e.g., 10^-6 N). For the rate of change of strain energy, the preset convergence threshold can be set as the relative value of the change in total strain energy between two adjacent iterations (e.g., less than 0.001%), or an absolute small value. One determination method is that both conditions are met simultaneously. After each iteration, the non-equilibrium force norm and the rate of change of strain energy are calculated and checked. Only when both of these indicators are simultaneously below their respective preset convergence thresholds is it determined that coupled stress equilibrium has been reached. Another method is a weighted comprehensive determination. The non-equilibrium force norm and the rate of change of strain energy can be weighted and combined to form a comprehensive convergence index. When this comprehensive index is below a preset threshold, coupled stress equilibrium is determined to have been reached.

[0112] Through the above technical solution, the contact force and relative displacement are determined based on the established contact relationship. This step, by quantifying physical interaction parameters, provides a basic input for the iterative process and avoids deviations caused by ambiguous boundary conditions. The displacement and stress fields are iteratively adjusted based on the contact force and relative displacement. By dynamically optimizing field variables, a gradual optimization of stress distribution is achieved, improving computational efficiency and accuracy. When both the non-equilibrium force norm and the rate of change of strain energy are below the preset convergence threshold, equilibrium is determined, and the displacement and stress fields are output. This dual convergence condition, by simultaneously monitoring force and energy indices, ensures robust convergence of the calculation results, preventing premature or invalid convergence, thus laying a mechanical foundation for subsequent simulations. Overall, this solution significantly enhances the reliability of model initialization through structured iteration and intelligent convergence determination, thereby improving the accuracy and safety of the heterogeneous coupled digital twin model in subsequent simulations, and ultimately optimizing the overall performance of the bonded ore recovery operation scheme.

[0113] This application further proposes a method for managing the recovery of caking ore in a stope, which includes:

[0114] When simulating in this heterogeneous coupled digital twin model, based on the multiple bonding type partitions and the quantified rock mass mechanical parameters, multi-field coupled simulations are performed on at least one candidate recovery operation scheme. The operation parameters of the candidate recovery operation scheme are iteratively optimized with ore recovery indicators and safety risk indicators as optimization objectives to generate the target recovery operation scheme. Specifically, the following steps are included:

[0115] In this heterogeneous coupled digital twin model, a dynamic disturbance source corresponding to the candidate recovery operation scheme is defined, and its spatiotemporal evolution parameters and energy release parameters are configured. The dynamic disturbance source refers to the external force or energy input applied to the bonded ore pile by simulating actual mining operations (such as blasting, hydraulic fracturing, static dilatation, etc.) in the digital twin model. Its purpose is to reproduce the destructive process of the ore body caused by real operations in a virtual environment. The dynamic disturbance source can be defined by applying instantaneous or continuous forces, pressures, or energy density fields within a pre-defined geometric region in the model. The spatiotemporal evolution parameters describe how the dynamic disturbance source changes and acts in the spatial and temporal dimensions; these parameters determine the range, duration, and propagation speed of the disturbance. The energy release parameters describe the magnitude and distribution characteristics of the energy released by the dynamic disturbance source; these parameters directly affect the degree of damage caused to the ore body by the disturbance.

[0116] The heterogeneous coupled digital twin model was run to perform multi-field coupled simulations of the stress wave propagation, rock fracture, and fragment transport processes triggered by the dynamic disturbance source. Simulation results were obtained, and the estimated values ​​of ore recovery indicators and safety risk indicators were extracted from these results. This multi-field coupled simulation refers to simultaneously considering and simulating the interactions and influences between different physical fields (such as mechanical fields, fluid fields, and thermal fields) during the simulation process. For mining operations, stress fields, displacement fields, damage fields, and fluid fields (such as explosive gases or water flows) are typically involved. The simulation results can include mechanical field data such as stress, strain, displacement, velocity, and acceleration at various points in the model, as well as the degree of rock mass damage, crack propagation paths, fragment geometry, size distribution, and movement trajectories. The estimated ore recovery indicators are used to quantitatively evaluate the performance of the recovery operation scheme in terms of ore recovery efficiency and quality. For example, these can be the average fragment size, the proportion of qualified fragments (e.g., the mass or volume percentage of fragments within the particle size range that meets the requirements of subsequent beneficiation processes), ore loss rate, and dilution rate. This safety risk indicator estimate is used to quantitatively assess the impact of the recovery operation plan on the stability of the mining area structure, personnel and equipment safety, etc. For example, it can be the maximum displacement of the surrounding rock of the mining area, the depth of the plastic zone, the stress concentration factor, the degree of damage to key structures (such as ore extraction roadways and support structures), the frequency or energy of microseismic events, etc.

[0117] Using the estimated values ​​of the ore recovery index and the safety risk index as the objective function values, a multi-objective optimization algorithm is driven to adjust the spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source, generating an updated set of parameters. The dynamic disturbance source is reconfigured using these updated parameters, and the steps of running the heterogeneous coupled digital twin model, extracting the index estimates from the simulation results, and driving the multi-objective optimization algorithm for adjustment are repeated until a preset optimization convergence criterion is met. The operation plan corresponding to the parameter set determined in the final iteration is output as the target ore recovery operation plan. The objective function value refers to the numerical value used to quantitatively evaluate the merits of candidate solutions during the optimization process. Multi-objective optimization typically involves multiple conflicting objectives, requiring these objectives to be transformed into computable function values. For example, these could be the estimated values ​​of the ore recovery index and the safety risk index themselves, or their weighted combination, points on the Pareto front, etc. This multi-objective optimization algorithm can simultaneously handle two or more conflicting optimization objectives and find a set of Pareto optimal solutions. This target recovery operation plan was determined through multi-objective optimization to achieve the optimal balance between ore recovery and safety risks. The plan includes the spatiotemporal evolution parameters and energy release parameters of dynamic disturbance sources, such as the borehole parameters, charge structure, detonation sequence, and delay time for blasting operations, or the nozzle type, water pressure, flow rate, and scanning path for hydraulic fracturing operations.

[0118] By defining dynamic disturbance sources corresponding to candidate recovery operation schemes and configuring their spatiotemporal evolution parameters and energy release parameters, the disturbance sources can accurately simulate the dynamic changes in actual operations, avoiding the blind setting of parameters and thus providing a realistic basis for simulation. Running a heterogeneous coupled digital twin model, multi-field coupled simulations are performed on the stress wave propagation, rock mass fracturing, and fragment transport processes triggered by the dynamic disturbance sources. This comprehensively simulates the physical interaction process and provides detailed simulation results, laying a data foundation for subsequent evaluation. The estimated values ​​of ore recovery indicators and safety risk indicators from the simulation results are extracted, quantifying the scheme's performance in resource recovery and safety, providing an objective basis for optimization, and avoiding biases caused by subjective experience. Using these indicators as the optimization objective function values, a multi-objective optimization algorithm is driven to adjust the parameters of the dynamic disturbance sources, automatically balancing multiple optimization objectives, resulting in systematic and efficient parameter adjustment. The simulation, extraction, and optimization steps are repeated until the preset convergence criteria are met. Through an iterative mechanism, the parameters gradually approach the optimal solution, avoiding the risks of overloading or inefficient schemes. The final output is a target recovery operation plan, which generates a reliable operation plan based on simulation optimization, improving the overall feasibility and safety of the operation.

[0119] This application further proposes a step in defining a dynamic disturbance source corresponding to the candidate recovery operation scheme in the heterogeneous coupled digital twin model, and configuring the spatiotemporal evolution parameters and energy release parameters of the disturbance source. This includes: selecting matching disturbance generation mechanisms for different types of bonded partitions based on the mechanical properties of the multiple bonded partition types, including blasting, high-pressure water jetting, or static fracturing. For each selected disturbance generation mechanism, initial spatiotemporal evolution parameters and energy release parameters are set based on the quantified rock mass mechanical parameters of the bonded partition corresponding to the disturbance generation mechanism. The set initial spatiotemporal evolution parameters and energy release parameters are mapped to the corresponding spatial location in the heterogeneous coupled digital twin model, completing the definition and configuration of the dynamic disturbance source.

[0120] Specifically, based on the mechanical properties of these multiple bonded ore zones, matching disturbance mechanisms are selected for different types of bonded ore zones. This step aims to choose the most suitable ore crushing or loosening method according to the actual mechanical characteristics of different areas within the bonded ore pile. For example, for bonded ores with extremely high strength and poor toughness, blasting operations may be preferred to achieve rapid and large-scale crushing. For areas with moderate strength and many micro-cracks, high-pressure water jetting may be more efficient and cause less disturbance to the surrounding rock. For areas requiring precise control and avoiding severe vibration, static fracturing technology can be considered. By comprehensively evaluating the mechanical properties of each bonded ore zone, such as compressive strength, tensile strength, elastic modulus, and Poisson's ratio, a decision-making model or expert can be established to automatically recommend or assist manual selection of the most suitable disturbance mechanism. Another approach is to construct a machine learning-based classifier by analyzing historical operational data and the effects of different disturbance mechanisms. The classifier takes the mechanical properties of the bonded ore zone as input and outputs the optimal disturbance mechanism. The disturbance mechanisms include blasting, high-pressure water jetting, and static fracturing, which are specific techniques used to break or loosen bonded ores. Blasting typically involves drilling holes inside the ore body and loading explosives. The shock wave and high-temperature, high-pressure gas generated by the explosion break the rock. Its advantages include concentrated energy and high breaking efficiency, but it causes significant disturbance to the surrounding environment. High-pressure water jetting technology uses a high-pressure pump to pressurize water to extremely high speeds, forming a high-speed water jet. The impact, shearing, and cavitation effects of the water jet cut or break the ore. Its advantages include no sparks, no vibration, and environmental friendliness, making it suitable for precision operations. Static fracturing technology involves injecting an expanding agent or hydraulic medium into the borehole, causing it to slowly expand and generate enormous pressure, thus creating cracks inside the ore and gradually fracturing them. Its advantages include safety, no vibration, and high controllability, but its operational efficiency is relatively low.

[0121] For each selected disturbance mechanism, initial spatiotemporal evolution parameters and energy release parameters are set based on the quantified rock mass mechanical parameters of the corresponding bonding type zone. This step aims to provide an initial parameter configuration based on actual geological conditions for the selected disturbance mechanism. For example, if a blasting mechanism is selected, initial spatiotemporal evolution parameters such as explosive consumption, hole pattern parameters (hole spacing, row spacing), detonation sequence, and delay time, as well as energy release parameters such as explosive type and charge amount, can be set according to the quantified rock mass mechanical parameters of the bonding type zone (such as rock density, wave velocity, and blast resistance). If a high-pressure water jet mechanism is selected, initial spatiotemporal evolution parameters such as water pressure, flow rate, nozzle diameter, jet velocity, and jet path, as well as energy release parameters such as water jet power, can be set according to parameters such as ore hardness and toughness. These initial parameter settings can be based on pre-established empirical formulas, numerical simulation result databases, or expert knowledge bases to ensure that the parameters match the mechanical properties of the ore. The initial spatiotemporal evolution parameters and energy release parameters are mapped to their corresponding spatial locations in the heterogeneous coupled digital twin model, thus completing the definition and configuration of the dynamic disturbance source. This step accurately imports the specific operational parameters determined above into the digital twin model, making them part of the model simulation. For example, in the digital twin model, parameters such as the location, depth, charge amount, and initiation time of the blast hole can be defined in the form of load boundary conditions, initial stress fields, or material property changes in the Discrete Element Method (DEM) or Finite Element Method (FEM). For high-pressure water jets, they can be simulated as moving loads or erosion boundary conditions in a specific region of the model. Static fracturing can be simulated as volume expansion forces or material property changes in a specific region of the model. Through this mapping, the digital twin model can simulate the mode, range, and intensity of the disturbance source in the mining space, providing accurate input for subsequent multi-field coupled simulations.

[0122] By employing the aforementioned technical solution and selecting matching disturbance mechanisms for different types of bonded ore zones based on their mechanical properties, this application can choose the most suitable crushing or loosening method according to the actual heterogeneity within the bonded ore pile. This avoids the blind "one-size-fits-all" approach of traditional empirical methods, thereby improving the applicability and efficiency of the operational plan. For example, for hard bonded ore, blasting can be selected to achieve efficient crushing. For ore susceptible to hydraulic forces, high-pressure water jets can be used to reduce vibration and damage to the surrounding rock. Furthermore, for each selected disturbance mechanism, initial spatiotemporal evolution parameters and energy release parameters are set based on the quantified rock mechanics parameters of the bonded ore zone corresponding to that disturbance mechanism. This ensures that the setting of disturbance source parameters has a solid scientific basis and data support. This makes the initial parameter configuration more accurate, reduces simulation errors caused by improper parameters, and provides a more reliable starting point for subsequent iterative optimization. By mapping the initial spatiotemporal evolution parameters and energy release parameters to their corresponding spatial locations within the heterogeneous coupled digital twin model, the definition and configuration of the dynamic disturbance source were completed. This ensures that the simulation model accurately reflects the actual operation process, making the simulation results more reliable. In summary, this application, through a refined, data-driven method for defining and configuring disturbance sources, enables simulations in the heterogeneous coupled digital twin model to more accurately predict ore recovery effects and potential safety risks. This lays a solid foundation for generating efficient and safe mining site ore recovery operation schemes, thereby effectively improving ore recovery efficiency and operational safety.

[0123] This application further proposes to run the heterogeneous coupled digital twin model to perform multi-field coupled simulation of the stress wave propagation, rock mass fracturing, and fragment migration processes triggered by the dynamic disturbance source, and obtain simulation results. Specifically, the method includes: in the heterogeneous coupled digital twin model, converting the spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source into dynamic load boundary conditions acting on the corresponding region of the model; executing an explicit dynamic analysis algorithm to solve the transient dynamic response of the heterogeneous coupled digital twin model under the dynamic load boundary conditions, simulating the continuous dynamic process, which includes stress wave propagation, rock mass fracturing, and fragment generation and migration; recording and outputting the spatiotemporal evolution data of key physical fields and the geometric information of fragments after rock mass failure during the continuous dynamic process, and using the spatiotemporal evolution data and the fragment geometric information as the simulation results.

[0124] Specifically, the spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source are converted into dynamic load boundary conditions acting on the corresponding region of the model. This aims to transform abstract operational parameters into physical inputs that the model can recognize. Dynamic load boundary conditions refer to forces, pressures, or displacements that vary over time and are applied to a specific region of the model, used to simulate the actual physical effects of the disturbance source. For example, in blasting operations, the spatiotemporal evolution parameters and energy release parameters such as the explosive charge, blast hole arrangement, and detonation delay can be converted into time-varying pressure pulses acting on the blast hole wall. The peak value, duration, and attenuation law of these pressure pulses can be set according to the characteristics and charge of the explosive. Similarly, in high-pressure water jet operations, parameters such as the water jet pressure, flow rate, and nozzle trajectory can be converted into local pressure loads or erosion boundary conditions moving on the rock surface.

[0125] An explicit dynamic analysis algorithm is executed to solve the transient dynamic response of the heterogeneous coupled digital twin model under dynamic load boundary conditions, simulating the continuous dynamic process, including stress wave propagation, rock mass fracturing, and fragment generation and migration. The aim is to capture the rapid, nonlinear response of the rock mass under dynamic disturbances. Explicit dynamic analysis is a numerical computation method particularly suitable for simulating transient dynamic problems involving large deformations, high strain rates, and complex contact conditions, such as impact, explosion, and fracture processes. This algorithm iteratively solves the equations of motion with small time steps, effectively tracking the propagation path and attenuation characteristics of stress waves in heterogeneous rock masses, simulating crack initiation, propagation, and eventual rock mass fracturing, and further simulating the movement, collision, and accumulation of fragments after fracturing. For example, explicit dynamic software based on the finite element method (such as LS-DYNA or Abaqus / Explicit) can be used in conjunction with constitutive models of the rock mass (such as elastoplastic models, damage models, or fracture mechanics models) to simulate the rock mass fracturing process. In addition, discrete element method (such as PFC or UDEC) can be used to discretize the rock mass into a series of interacting particles or blocks, directly simulating their contact, separation and movement, thus more intuitively showing the generation and migration of fragments.

[0126] This simulation records and outputs the spatiotemporal evolution data of key physical fields and the geometric information of fragments after rock mass failure during the continuous dynamic process. This spatiotemporal evolution data and fragment geometric information are used as the simulation results, aiming to provide a comprehensive and reliable data foundation for subsequent index extraction and optimization decisions. Key physical fields include, but are not limited to, stress, strain, displacement, velocity, and damage fields. Their spatiotemporal evolution data records the changes of these physical quantities at different spatial locations and time points. Fragment geometric information includes the size, shape, spatial distribution, and average size of the fractured rock fragments. For example, during the simulation, specific monitoring points or areas can be set to record the time history curves of stress, strain, and displacement in real time. Simultaneously, after the simulation, post-processing tools are used to identify and statistically analyze the fracture elements or discrete elements in the model, extracting geometric features such as fragment volume, surface area, and aspect ratio, and calculating the fragment particle size distribution. This data can be output in the form of structured files (such as CSV, HDF5) or visualization files (such as VTK).

[0127] Through the above technical solution, the abstract parameters of dynamic disturbance sources are transformed into dynamic load boundary conditions that the model can recognize, thereby improving the physical realism and operability of the simulation input and avoiding errors caused by parameter abstraction. Simultaneously, the execution of explicit dynamic analysis algorithms can capture highly nonlinear and transient continuous dynamic processes such as stress wave propagation, rock mass fracturing, and fragment migration, overcoming the problems of process fragmentation or simplification in traditional methods, thus improving the accuracy of simulation results. Furthermore, comprehensive recording of the spatiotemporal evolution data of key physical fields and the geometric information of fragments after rock mass failure provides multi-dimensional simulation output, offering a solid data foundation for the subsequent extraction and optimization decisions of ore recovery indicators and safety risk indicators, and solving the problem of insufficient prediction reliability caused by incomplete data recording. Overall, this application achieves refined control of the simulation process through the synergy of parameter transformation, algorithm execution, and data recording, laying a reliable foundation for optimizing recovery operation schemes and solving the technical problem of how to simulate continuous dynamic processes triggered by dynamic disturbance sources and effectively record key information to ensure the accuracy and comprehensiveness of simulation results.

[0128] This application further proposes a method for extracting the estimated values ​​of ore recovery indicators and safety risk indicators from the simulation results. Specifically, this includes: statistically analyzing the geometric information of the fragments in the simulation results to obtain the estimated value of the ore recovery indicator, which includes the average size of the fragments and the proportion of qualified fragments; and processing the spatiotemporal evolution data of the stress field and plastic strain field in key areas of the simulation results to obtain the estimated value of the safety risk indicator, which includes the maximum plastic zone depth and the peak displacement of key points.

[0129] The fragment geometry information in this simulation result refers to the spatial attribute data, such as shape, size, and location, of the broken ore blocks obtained during multi-field coupled simulation by simulating rock mass fracturing and fragment transport. This information can be directly output by simulation software. For example, in methods such as Discrete Element Method (DEM) or Smoothed Particle Hydrodynamics (SPH), the geometric information of fragments is usually represented in the form of discrete particles or particle assemblies, including the diameter, shape, and centroid coordinates of each particle. Another approach is to post-process the damage or crack propagation region in the continuous medium model (such as the finite element method) to identify and extract independent fragmentation elements and calculate the geometric characteristics of these elements.

[0130] Statistical analysis of the fragment geometry information refers to the systematic quantification of a large amount of fragment geometry data to reveal its overall distribution characteristics and key parameters. One approach is to use methods such as histogram analysis and cumulative distribution functions to statistically analyze the size distribution of the fragments and calculate statistical quantities such as the mean, median, and variance. Another approach is to utilize image processing or 3D point cloud processing techniques to identify and classify the shapes of the fragments and calculate the proportion of fragments with different shapes.

[0131] The average size of the fragments is a macroscopic indicator of the degree of ore crushing, reflecting the average level of ore block size after crushing. One calculation method is to take the average of the cube roots of the equivalent diameters or volumes of all fragments as the average size. Another method is to integrate the fragment size distribution curve or use a weighted average method to consider the contribution of fragments of different sizes to the total mass or volume, thus obtaining a more representative average size.

[0132] The percentage of acceptable block size refers to the proportion of all crushed ore that meets preset recycling standards (e.g., blocks within a specific size range). One calculation method is to set one or more block size thresholds, count the number or volume of all blocks falling within these threshold ranges, and then calculate the ratio of this to the total number or volume of blocks. Another calculation method is to define a continuous acceptable block size function based on the feed requirements of the ore processing equipment, evaluate the acceptable level of each block, and then perform a weighted sum to obtain the overall acceptable block size percentage.

[0133] The simulation results include spatiotemporal evolution data of stress and plastic strain fields in key areas. The stress field refers to the distribution of internal forces at various points within the rock mass during the simulation, while the plastic strain field represents the degree and distribution of irreversible deformation of the rock mass under stress. Spatiotemporal evolution data refers to the records of these physical quantities changing over time and spatially. Key areas typically refer to areas with significant impacts on operational safety, such as the bottom of the mining area, the ore outlet, and the contact surface with the surrounding rock. These data can be calculated and recorded at each time step and each grid point using numerical simulation methods such as the finite element method and the finite difference method. For example, in simulation software, specific monitoring points or surfaces can be specified, and their stress tensor, plastic strain tensor, and other data can be output throughout the simulation process.

[0134] Processing the spatiotemporal evolution data of the stress field and plastic strain field involves screening, calculating, and analyzing the simulated data to extract quantitative indicators directly related to safety risks. One approach is to perform time-series analysis on stress or strain data in specific regions to identify peak values, rates of change, and other characteristics. Another approach is to use data visualization technology to display the stress field and plastic strain field as cloud maps, and combine this with threshold judgments to automatically identify high stress concentration areas or plastic deformation regions.

[0135] The maximum plastic zone depth refers to the maximum extension depth of the rock mass undergoing plastic deformation along a specific direction (e.g., perpendicular to the stope floor or parallel to the ore extraction roadway) during the simulation. One method to determine this is to identify all elements in the simulation results whose plastic strain values ​​exceed a preset threshold, and then calculate the maximum projected distance of these elements in the specified direction. Another method is to find the plastic zone boundary through contour analysis of the plastic strain field and measure the maximum distance from this boundary to a reference plane.

[0136] The peak displacement of this key point refers to the maximum displacement of pre-defined key monitoring points (e.g., roof of ore extraction roadways, surrounding rock slopes, support structures, etc.) in the mining area structure during the simulation process. One method is to pre-set multiple monitoring points in the simulation model, record the three-dimensional displacement time series of these points throughout the simulation process, and then extract the maximum absolute value of the displacement at each point. Another method is to post-process the entire displacement field to identify several points with the largest displacements and use their peak displacements as safety risk indicators.

[0137] The above technical solution utilizes statistical analysis of the geometric information of fragments in the simulation results. By leveraging the size and shape data of the fragments generated in the simulation, the physical state of the ore after crushing is directly reflected, thus generating a predicted ore recovery index. This predicted index includes the average fragment size and the proportion of qualified fragments. The average fragment size quantifies the overall degree of ore crushing, and the proportion of qualified fragments assesses the percentage that meets the recovery standards, providing a concrete basis for optimizing recovery efficiency. By processing the spatiotemporal evolution data of stress field and plastic strain field in key areas of the simulation results, focusing on the changes in safety-related physical fields in the model's dynamic response, a predicted safety risk index is generated. This predicted index includes the maximum plastic zone depth and the peak displacement of key points. The maximum plastic zone depth reveals the depth of rock mass failure, and the peak displacement of key points quantifies potential deformation risks, providing measurable indicators for optimizing safety performance. The entire solution, by combining fragment geometry and physical field data, achieves accurate extraction of recovery and safety indicators, supporting the effective operation of the optimization algorithm. This allows subsequent optimization processes to be based on quantitative and objective data, improving the accuracy and reliability of the recovery operation optimization plan.

[0138] This application further proposes using the estimated values ​​of ore recovery index and safety risk index as the optimization objective function values ​​to drive a multi-objective optimization algorithm to adjust the spatiotemporal evolution parameters and energy release parameters of dynamic disturbance sources, generating an updated set of parameters. The specific steps include: constructing a multi-objective fitness function to evaluate the comprehensive performance of the scheme, using the estimated values ​​of ore recovery index and safety risk index as inputs. During the initial execution, the spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source are used as the initial parameter set. Based on the evaluation results of the initial parameter set using the multi-objective fitness function, the multi-objective optimization algorithm is driven to search within a preset parameter constraint space, generating and outputting the first set of candidate parameters. It is then determined whether the search process of the multi-objective optimization algorithm satisfies a preset convergence criterion. If satisfied, the newly generated candidate parameter set is output as an updated set of parameters. If not satisfied, the newly generated candidate parameter set is used as new input, and the steps of evaluation based on the fitness function, algorithm search, and convergence determination are repeated.

[0139] Specifically, constructing a multi-objective fitness function is a mathematical expression used to quantitatively evaluate the comprehensive performance of candidate recovery operation schemes. It takes the predicted values ​​of ore recovery indicators and safety risk indicators as inputs, aiming to balance these two often conflicting objectives. For example, the two indicators can be linearly combined using a weighted summation method, where the weights can be adjusted according to actual production needs and safety priorities. Alternatively, a Pareto dominance-based method can be used to identify a set of non-dominated solutions, i.e., the Pareto front, where neither objective can be improved without sacrificing the other. The construction of this function allows the optimization process to comprehensively consider recovery efficiency and operational safety, providing a clear evaluation basis for subsequent optimization algorithms.

[0140] The spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source are used as the initial parameter set, which serves as the starting point for the multi-objective optimization algorithm's search. Its selection significantly impacts the efficiency of the optimization process and the quality of the final result. One approach is to establish a set of validated parameters as initial values ​​based on historical operational data or expert experience, combined with the characteristics of the bonded ore pile in the stope. For example, this could be based on the charge quantity, detonation sequence, and interval time from previous successful blasting cases. Another approach is to conduct small-scale pre-simulation or sensitivity analysis to quickly determine a relatively high-performing parameter combination as the initial parameter set. This provides a reasonable starting point for subsequent iterative optimization, avoiding a completely random parameter search and improving optimization efficiency.

[0141] Driving a multi-objective optimization algorithm to search within a predefined parameter constraint space, generating and outputting the first set of candidate parameters, refers to using the algorithm's mechanism to find a better solution within a defined parameter range. Multi-objective optimization algorithms can take various forms, such as evolutionary computation-based algorithms like Non-Dominated Sorting Genetic Algorithm II (NSGA-II) or Multi-Objective Particle Swarm Optimization (MOPSO), which explore the solution space by simulating natural selection and genetic mechanisms. Alternatively, decomposition-based multi-objective evolutionary algorithms (MOEA / D) can be used, breaking down the multi-objective problem into a series of single-objective sub-problems for solution. The predefined parameter constraint space refers to the reasonable range of values ​​set for the spatiotemporal evolution parameters (such as blast hole spacing and detonation delay) and energy release parameters (such as charge amount and explosive type) of the dynamic disturbance source. These constraints are typically based on engineering experience, equipment capabilities, and safety regulations, ensuring the algorithm searches within a practically feasible range.

[0142] To determine whether a multi-objective optimization algorithm's search process satisfies a pre-defined convergence criterion, this criterion is used to determine when the algorithm stops iterating, indicating that a sufficiently good solution has been found. A common convergence criterion is that if the improvement in the objective function value is less than a pre-defined small threshold over N consecutive iterations, the algorithm is considered close to a local or global optimum. Another convergence criterion could be reaching a pre-defined maximum number of iterations to control computational resources and time costs. Furthermore, the uniformity or diversity of the Pareto front can be evaluated; when it no longer improves significantly, the algorithm is considered convergent. These criteria allow the optimization process to achieve the desired results while avoiding unnecessary computational overhead.

[0143] If the conditions are met, the newly generated candidate parameter set is output as an updated set of parameters. If not, the newly generated candidate parameter set is used as new input, and the steps of evaluation based on the fitness function, driving the algorithm search, and determining convergence are repeated. This step describes the iterative working mechanism of the multi-objective optimization algorithm. When the convergence criterion is met, the algorithm stops and outputs the current optimal or non-dominated parameter set as the final updated parameters to guide the actual operation. If the convergence criterion is not yet met, the algorithm will use the newly generated candidate parameter set as input for the next iteration to continue fitness evaluation, algorithm search, and convergence determination. This iterative feedback mechanism enables the optimization algorithm to gradually approach the optimal solution, continuously improve the operation plan, until the expected optimization objective is reached or the preset stopping condition is met, thereby ensuring the reliability and effectiveness of the optimization results.

[0144] The above technical solution constructs a multi-objective fitness function to evaluate the overall performance of a solution by using the estimated values ​​of ore recovery indicators and safety risk indicators as inputs. This allows the optimization process to comprehensively and objectively balance recovery efficiency and operational safety, avoiding the bias that may result from optimizing a single indicator and improving the overall benefits of the optimization results. During the initial execution, the spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source are used as the initial parameter set. This allows the optimization algorithm to start its search from a practically meaningful starting point, improving search efficiency and the rationality of the initial solution, and avoiding blind exploration. Based on the evaluation results of the initial parameter set using the multi-objective fitness function, the multi-objective optimization algorithm is driven to search within the preset parameter constraint space, generating and outputting the first set of candidate parameters. This ensures that the optimization process is conducted within a feasible range and can efficiently explore the solution space, quickly discovering potential superior solutions. By determining whether the search process of the multi-objective optimization algorithm meets the preset convergence criteria, and deciding whether to repeat the evaluation, search, and judgment steps based on the judgment results, this application establishes an adaptive iterative optimization mechanism. This mechanism controls the stability and convergence of the optimization process, avoiding premature stopping or infinite loops, thus ensuring that the final output set of updated parameters is fully optimized and reliable. Overall, these steps work synergistically, making the simulation optimization performed in the heterogeneous coupled digital twin model more efficient, stable, and reliable. The resulting target recovery operation scheme can better balance ore recovery rate and operational safety risks, improving the overall efficiency and safety of the mining area's bonded ore recovery operation.

[0145] This application further proposes a dynamic comparison and model assimilation process between real-time collected construction feedback data and corresponding predicted data in a heterogeneous coupled digital twin model. Specifically, this includes: during construction, real-time collection of construction feedback data via sensors mounted on the work equipment, including at least drilling process parameters and vibration monitoring data. This construction feedback data is then synchronously compared with the predicted mechanical response data of the same spatial location at the same construction stage in the heterogeneous coupled digital twin model to determine the deviation value between the two. Based on this deviation value, it is determined whether a model assimilation condition is triggered. If triggered, the mechanical property parameters of a local area in the heterogeneous coupled digital twin model are reverse-calibrated and updated based on the construction feedback data to complete model assimilation and obtain the model assimilation result.

[0146] During construction, sensors mounted on the operating equipment collect construction feedback data in real time. This feedback data includes at least drilling process parameters and vibration monitoring data. "Real-time sensor acquisition of construction feedback data" refers to the continuous acquisition of performance or environmental data directly related to the mining operation (e.g., drilling, blasting, hydraulic fracturing) using physical or electronic devices. For example, accelerometers, strain gauges, and pressure sensors can be installed on the drilling equipment to monitor parameters such as drill bit stress, drilling speed, torque, and borehole pressure in real time. Alternatively, microseismic sensor arrays and acoustic emission sensors can be deployed near the blasting or hydraulic fracturing area to capture vibration signals such as rock fracture and stress wave propagation in real time. "Drilling process parameters" reflect the real-time status of the drilling operation and the rock mass response, serving as crucial evidence for evaluating drilling efficiency and rock mass properties. These parameters may include drilling speed, drilling pressure, torque, impact frequency, and footage, which can be directly output by the drilling rig's built-in controls. Alternatively, they may include borehole mud pressure, flow rate, and temperature, acquired through sensors installed near the drill pipe or borehole opening. Vibration monitoring data reflects the dynamic response of rock mass under operational disturbances and is a key indicator for assessing rock mass stability, fracturing degree, and safety risks. For example, this data can include vibration velocity, vibration acceleration, and vibration frequency, acquired through vibration sensors (such as three-component accelerometers) installed on the operating equipment, surrounding rock, or surrounding structures. Alternatively, it can include the occurrence time, source location, and energy magnitude of microseismic events, enabling real-time location and analysis through microseismic monitoring.

[0147] The construction feedback data is synchronously compared with the predicted mechanical response data of the same spatial location and construction stage in the heterogeneous coupled digital twin model to determine the deviation between the two. "Synchronous comparison" aims to ensure consistency between actual observations and model predictions in time and space, thereby accurately assessing the model's prediction accuracy and the actual changes in rock mass condition. For example, real-time collected construction feedback data (such as drilling parameters at a borehole location) can be matched with the predicted data (such as predicted rock mass strength and stress response) of the digital twin model at that borehole location, corresponding drilling depth, and time step through timestamp alignment and spatial coordinate mapping. Alternatively, data fusion algorithms can be used to interpolate or resample data with different sampling frequencies and spatial resolutions to achieve synchronous comparison. "Predicted mechanical response data" refers to the physical quantities such as stress, strain, displacement, and failure mode of the rock mass at corresponding time and space points, calculated by the heterogeneous coupled digital twin model based on its internal rock mass mechanical properties and boundary conditions when simulating specific construction operations (such as drilling and blasting). For example, the model predicts the stress distribution, plastic zone range, and displacement field of the rock mass surrounding the borehole. Alternatively, it predicts the initial velocity, block size distribution, and vibration intensity of rock fragments under blasting. The "deviation value" quantifies the difference between actual observations and model predictions, providing a quantitative basis for subsequent model assimilation and operational adjustments. For example, it can be an absolute difference, relative difference, root mean square error (RMSE), or specific statistical indicators (such as correlation coefficients or difference coefficients). Alternatively, it can compare the actual measured drilling speed with the model-predicted drilling speed and calculate the difference, or compare the actual monitored vibration peak value with the model-predicted vibration peak value and calculate the difference.

[0148] Based on this deviation value, it is determined whether to trigger the model assimilation condition. The "model assimilation condition" is a preset rule or threshold used to determine when to initiate or execute a digital twin model update or calibration. For example, model assimilation is triggered when the deviation value exceeds a preset tolerance threshold (such as a drilling speed deviation greater than 10% or a vibration peak deviation greater than 20%). Alternatively, model assimilation is triggered when the deviation values ​​at multiple consecutive time steps show a certain trend (such as continuous increase), or when the deviation value reaches a certain critical value.

[0149] If triggered, the mechanical property parameters of a local area in the heterogeneous coupled digital twin model are reverse-calibrated and updated based on the construction feedback data to complete model assimilation and obtain the model assimilation result. The role of "reverse calibration and update" is to correct the parameters inside the model based on actual observation data, so that the model more accurately reflects the physical state of the real world. For example, data assimilation algorithms such as Kalman filtering, particle filtering, and ensemble Kalman filtering can be used to integrate construction feedback data into the model and inversely update the rock mechanical parameters (such as elastic modulus, Poisson's ratio, compressive strength, cohesion, internal friction angle, etc.) of the local area in the model. Alternatively, optimization algorithms (such as genetic algorithms and neural networks) or Bayesian inversion methods can be used to iteratively adjust the model parameters with the goal of minimizing the deviation value until the model prediction highly matches the actual observation. "Local area" refers to the finite spatial range in the digital twin model that is directly related to or affected by the current construction activity. For example, the rock mass area around the current borehole, the rock mass unit near the blasting hole, and the ore area within the range of hydraulic fracturing. Alternatively, based on the spatial distribution and impact range of construction feedback data, the model mesh or elements that require parameter updates can be dynamically determined.

[0150] The above technical solution synchronously compares construction feedback data with predicted mechanical response data at the same spatial location and construction stage in a heterogeneous coupled digital twin model. Strict correspondence between specified locations and stages ensures spatial and temporal consistency in the comparison, reduces error accumulation, and improves the accuracy of deviation detection. Then, the deviation value between the two is determined, quantifying the difference between actual and predicted data and providing an objective basis for judging update needs. Based on the deviation value, the model assimilation condition is triggered, intelligently determining the update timing and avoiding unnecessary waste of computational resources, initiating updates only when there are critical differences. Subsequently, the mechanical property parameters of local areas in the heterogeneous coupled digital twin model are reverse-calibrated and updated based on the construction feedback data. Targeted adjustments to local areas rather than the global model save computational overhead and enable the model to quickly respond to actual changes. This completes model assimilation, obtaining the model assimilation result and achieving real-time model calibration. This lays an accurate foundation for dynamically adjusting operational parameters, thus solving problems such as incomplete data collection, asynchronous comparison, or delayed model updates, and improving the safety and efficiency of mining site bonded ore recovery operations.

[0151] This application further proposes a dynamic adjustment of unexecuted operational parameters in the target recovery operation plan based on model assimilation results. This adjustment process includes: identifying local regions in the model assimilation results where parameter updates have been completed and their corresponding updated mechanical property parameters; matching the spatial location of the local regions with unexecuted operational units in the target recovery operation plan to determine the affected operational units and their corresponding operational parameters to be adjusted; recalculating the operational parameters to be adjusted for the affected operational units based on the updated mechanical property parameters; and replacing the corresponding parameters in the target recovery operation plan with the recalculated operational parameters to complete the dynamic adjustment.

[0152] Specifically, identifying the local regions where parameter updates were completed in the model assimilation results and the corresponding updated mechanical property parameters aims to pinpoint the specific spatial range of the digital twin model where mechanical property parameter updates occurred due to model assimilation. Model assimilation is a process of comparing real-time construction feedback data with model prediction data and calibrating model parameters. By identifying these local regions, it can be determined which parts of the model parameters have been corrected, thus providing an accurate basis for subsequent operational parameter adjustments. For example, one implementation is to record the IDs of all modified mesh elements or finite element elements and their new mechanical property parameter values ​​during the model assimilation process. These modified element sets constitute the local regions where parameter updates were completed. Another implementation is to maintain a parameter update log, recording the bounding box or geometry of the region involved in each model assimilation operation, as well as the type and old and new values ​​of the updated mechanical property parameters within that region. By querying this log, the latest updated region can be identified.

[0153] The spatial location of the local area is matched with the unexecuted work units in the target recovery operation plan to determine the affected work units and their corresponding operational parameters to be adjusted. The purpose of this step is to establish a spatial correlation between the model parameter update area and the actual operation plan, thereby identifying which unexecuted work units will be affected by these model updates. Work units can be blasting holes, hydraulic cutting paths, or static fracturing points, etc., and each work unit contains a series of operational parameters, such as location, depth, charge amount, and water pressure. For example, one implementation is to spatially intersect or include the geometric location of each unexecuted work unit in three-dimensional space (e.g., the center point or trajectory of a blasting hole) with the identified local area. If the geometric range of a work unit overlaps with the local area, the work unit is considered affected. Another implementation is to predefine the influence range of each work unit and then determine whether the local area intersects with the influence range of any work unit. Once the affected work units are identified, the operational parameters related to the updated mechanical property parameters within that work unit are further identified; for example, if the updated rock mass strength parameter may require adjustment of the blasting charge amount or water pressure.

[0154] Based on this, the adjusted operational parameters for the affected work units are recalculated using the updated mechanical property parameters. This step involves re-optimizing the operational parameters of the affected work units based on the latest and more accurate local mechanical property parameters obtained after model assimilation. This allows subsequent operational parameters to better adapt to the current geological conditions, improving the accuracy and effectiveness of operations. For example, one approach is to use a pre-established operational parameter optimization model or empirical formula, taking the updated mechanical property parameters as input, to recalculate the optimal operational parameters. For blasting operations, the optimal charge amount, hole spacing, and detonation sequence can be re-determined using a blasting engineering calculation model based on new rock mass strength and toughness parameters. Another approach is to again call upon the heterogeneous coupled digital twin model for local simulation, but only for the affected work units and their surrounding areas. The updated mechanical property parameters are then substituted into the model, and the parameters of the work unit are iteratively optimized on a small scale, using ore recovery indicators and safety risk indicators as optimization targets, to obtain new operational parameters to be adjusted.

[0155] The corresponding parameters in the target recovery operation plan are replaced with the recalculated operational parameters to achieve dynamic adjustment. This step updates the original target recovery operation plan with the recalculated and optimized operational parameters, thereby enabling dynamic adjustment of the operation plan. This replacement operation ensures that subsequent construction is executed based on the latest and most accurate parameters, thus improving operational efficiency and safety. For example, one implementation is to directly write the recalculated parameters into the database or configuration file of the target recovery operation plan in the operation plan management and mark it as updated. Another implementation is to generate a patch file containing all updated parameters, which is automatically loaded and applied when the operation equipment executes the corresponding operation unit, ensuring the real-time performance and accuracy of the operation.

[0156] By identifying the local areas where parameters have been updated in the model assimilation results and their corresponding updated mechanical property parameters, the key areas after model calibration can be accurately located. This ensures that subsequent adjustments are based on the latest reliable data, avoiding blind modifications. Matching the spatial location of the local areas with the unexecuted work units in the target recovery operation plan achieves spatial correlation mapping, allowing adjustments to be made only to the affected areas, reducing resource consumption from ineffective adjustments and improving the targeting of adjustments. Based on this, the adjusted operation parameters of the affected work units are recalculated using the updated mechanical property parameters, and the parameters are re-optimized using the latest rock mass properties, making the operation plan more closely match actual mechanical conditions. Replacing the corresponding parameters in the target recovery operation plan with the recalculated operation parameters completes the dynamic adjustment of the operation plan, thereby ensuring the efficiency and safety of subsequent operations, improving the accuracy of ore recovery, and enhancing the overall adaptability of the operation. This dynamic adjustment mechanism enables the entire recovery management method to respond in real time to changes in the geological environment, continuously optimizing operation strategies, and achieving optimal recovery results and safety levels in complex and ever-changing environments.

[0157] Specifically, in one embodiment of this application, the bonding ore recovery project of the 1301 mining area in the 500-600m section of a copper mine in Xinjiang is used as an example for detailed explanation.

[0158] A terrestrial 3D laser scanner (RIEGL VZ-400i) was used to perform a 3D scan of the bonded ore pile in stope 1301, acquiring 3D point cloud data with a resolution of 5cm. After processing, a 3D surface model of the ore pile was obtained, accurately representing the external contour of the ore pile, which is 20 meters long and 18 meters high at its highest point. Ground-penetrating radar (MALA ProEx, equipped with a 100MHz shielded antenna) was used to continuously probe the bottom of the ore pile along the walls of the ore extraction roadway, acquiring internal electromagnetic reflection profiles and identifying the distribution of macroscopic cracks and cavities in the contact zone between the ore pile and the floor and sidewalls. Simultaneously, microseismic events from the microseismic monitoring system (IMS) deployed around the stope over the past week were retrieved, and the location of the seismic source (accuracy ±5 meters) and energy release (10²~10⁻⁶ m²) were analyzed. 4 J) revealed an active stress concentration zone in the central part of the ore pile. Five sets of core samples were drilled at different heights of the bonded ore pile and subjected to X-ray diffraction (XRD) and scanning electron microscopy-energy dispersive spectroscopy (SEM-EDS). Diagnostic data showed that the ore pile was mainly bonded by two genesis: sulfate cementation (mainly gypsum) at the bottom and carbonate cementation (mainly calcite and dolomite) in the middle and upper parts. Using a spatial registration algorithm, the 3D surface model, internal defect elements (fractures, cavities, stress zones), and genetic data of the core sampling points were fused to obtain a complete digital model that includes both external morphology and internal structural defects.

[0159] Data quality control and preprocessing: To ensure the accuracy of the fused data, strict quality control and preprocessing were carried out on various types of raw data: (1) Three-dimensional laser scanning data: Three different stations were set during scanning to ensure that the overlap of the point cloud was ≥30%. After collection, multi-station registration was performed using the ICP algorithm, and the registration error was controlled within 2cm. Statistical filtering was performed on the point cloud data to remove outlier noise points, and finally a complete and continuous three-dimensional surface point cloud was obtained. (2) Ground radar data: The common center point method was used for collection during detection, and each survey line was measured twice to ensure data repeatability. Zero-time correction, background removal, gain compensation and bandpass filtering were performed on the raw radar data to highlight the effective reflection signal and eliminate system noise. (3) Microseismic data: Events with a signal-to-noise ratio >3 and a positioning residual <5m were used for analysis, and non-rock mass rupture events such as blasting and mechanical interference were removed. Clustering algorithm was used to identify the spatiotemporal clustering characteristics of microseismic events to ensure that the identified stress concentration areas were statistically significant. (4) Core data: Three sets of parallel tests were performed on each core sample, and the average value was taken as the final data. Core segments with obvious cracks or fractures were not included in the mechanical parameter statistics to ensure the reliability of the parameters. Through the above quality control measures, the authenticity and accuracy of the input data were ensured, laying a reliable foundation for the subsequent construction of the digital twin model.

[0160] Based on the XRD and SEM-EDS analyses of the five core samples, the region with sulfur and calcium as the main chemical elements and gypsum mineral content of 15% was defined as the "sulfate-cemented zone" (mainly located at the bottom of the ore pile); the region with carbon, calcium, and magnesium as the main chemical elements and calcite and dolomite content of 20% was defined as the "carbonate-cemented zone" (mainly located in the upper middle and sides of the ore pile). Using Kriging interpolation, these two genetic models were extrapolated in three-dimensional space to form continuous bonded zone types. Uniaxial compressive strength tests and Brazilian splitting tests were conducted on the core samples from both zones. The results showed: sulfate-cemented zone: average compressive strength 25 MPa, elastic modulus 15 GPa, Poisson's ratio 0.28; carbonate-cemented zone: average compressive strength 65 MPa, elastic modulus 40 GPa, Poisson's ratio 0.22. These quantitative parameters were entered into the "Genetic-Mechanical Parameter Mapping Table" for subsequent model use.

[0161] Construction of the heterogeneous coupled digital twin model: The three-dimensional surface model and internal structural defect distribution obtained in the above steps are imported into the mining area stress background model (which already includes the initial geostress field of the mining area). Boolean operations are used to embed the geometric entity of the bonded ore pile into the model. Then, the two previously divided bonding type zones and their corresponding quantified mechanical parameter spaces are assigned to the corresponding regions in the model: the sulfate cemented zone is assigned a compressive strength of 25 MPa, and the carbonate cemented zone is assigned a compressive strength of 65 MPa. A contact relationship is established between the bonded ore pile and the surrounding rock (friction coefficient set to 0.6), and iterative calculations are performed based on the initial geostress conditions (vertical stress approximately 15 MPa, horizontal stress approximately 8 MPa) to solve for the coupled stress balance. After 200 iterations, the non-equilibrium force norm drops below 0.01%, obtaining a convergent mechanical equilibrium solution, which is used as the initial condition for subsequent dynamic simulations, thus obtaining the heterogeneous coupled digital twin model.

[0162] Optimization of Blasting Scheme Based on Digital Twin Model: In this embodiment, the candidate recovery operation scheme is fan-shaped medium-deep hole blasting. In the digital twin model, the dynamic disturbance source is defined as a fan-shaped medium-deep hole, with initial parameters referencing empirical values: hole diameter 80mm, row spacing 3.5m, hole bottom distance 3.0m, inclination angle 40°~75°, hole depth 6.5~16.5m, and the charge structure is continuous coupled loading of powdered ammonium nitrate explosive at the bottom 2.0m of the hole. Explicit dynamic analysis software (such as LSDYNA) is used for multi-field coupled simulation to simulate the stress wave propagation, rock mass fracturing, and fragment migration process after the explosive detonation. After each simulation, the following indicators are extracted: Ore recovery indicator: the proportion of ore volume with a block size <500mm after blasting (the proportion of qualified block size); Safety risk indicator: peak vibration velocity (PPV) and maximum plastic zone depth at key points on the roof of the mining roadway. With the goal of maximizing the percentage of qualified blocks and minimizing PPV, a multi-objective optimization algorithm (NSGA-II) was used to iteratively optimize the blasting parameters. The parameter constraints were: row spacing 2.5–4.0 m, hole bottom distance 2.5–3.5 m, single-hole charge 15–25 kg, and row-to-row delay 25–50 ms. After 50 iterations, the algorithm converged, and a set of parameters was selected from the Pareto front as the target scheme: row spacing 3.0 m, hole bottom distance 2.8 m, single-hole charge 18 kg, and row-to-row delay 30 ms. Row-by-row detonation was adopted, with lower-level detonators (e.g., 3-level) used for rows near the ore outlet, and higher-level detonators (e.g., 7-level) used for rows further away.

[0163] Verification and Sensitivity Analysis of Digital Twin Model: To ensure the prediction accuracy of the digital twin model, model verification was carried out before formal optimization: (1) Historical Data Inversion Verification: The actual data of a small-scale blasting operation in the early stage of the mining area were selected, and the actual blasting parameters were input into the model for simulation. The vibration waveform and fracture range obtained by simulation were compared with the actual monitoring data. The results showed that the average error between the peak vibration velocity obtained by simulation and the actual monitoring value was 12.3%, and the fracture range was more than 85% consistent, which verified the effectiveness of the model. (2) Sensitivity Analysis: Sensitivity analysis was carried out on key parameters (rock mass compressive strength, elastic modulus, charge amount, etc.). The results showed that when the rock mass compressive strength changed by ±20%, the proportion of qualified block size changed by about ±15%, which was highly sensitive; when the charge amount changed by ±20%, the PPV changed by about ±25%, which was highly sensitive; when the inter-row delay changed by ±20%, the fracture effect changed by about ±8%, which was moderately sensitive. The sensitivity analysis results provide a basis for adjusting the weights of parameters in the subsequent optimization process, namely, prioritizing the accuracy of rock mass parameters, optimizing the charge amount, and then adjusting the delay parameters.

[0164] Real-time feedback and dynamic adjustment: During actual construction in the 1301 mining area, sensors were installed on the YGZ-90 drilling rig to collect parameters such as drilling speed, drilling pressure, and torque in real time. When drilling reached the second row of blast holes, the actual drilling speed (0.8 m / min) deviated from the predicted value (1.2 m / min) of the digital twin model by -33%, exceeding the preset threshold (±15%), triggering model assimilation. The reverse calibration algorithm was activated, and the rock mass compressive strength of the local area (the location of the second row of holes) was corrected from 65 MPa to 78 MPa based on the actual drilling data. After model assimilation, it was identified that the subsequent third row and subsequent blast holes were located within the influence range of the corrected area (distance <5m). Therefore, the charge amount was recalculated based on the updated compressive strength of 78 MPa, increasing the single-hole charge amount from 18 kg to 22 kg, and automatically updating the parameters of the corresponding blast holes in the target recovery operation plan.

[0165] Effectiveness Verification: Following the dynamically adjusted plan, 30 fan-shaped medium-deep holes were drilled, blasted in two stages, consuming a total of 540 kg of explosives. After blasting, the bonded ore pile was effectively loosened, and a large amount of ore fell into the extraction roadway. A 4m³ loader was used for continuous ore extraction, recovering a total of 3900 tons of bonded ore. Sampling and analysis showed that the average grade of the recovered ore was 2.59% copper and 1.72% zinc, demonstrating significant economic benefits. Post-blasting inspection of the bottom structure of the stope revealed that the extraction roadway and the surrounding rock of the cross-cutting veins were stable, with only minor cracks and no structural damage, verifying the safety of the plan. The reduced height of the remaining ore pile did not affect the subsequent construction of the backfill retaining wall or the flow of backfill grout, creating safe conditions for the second-stage mining preparation project.

[0166] Before recovering the bonded ore in this mining area, two traditional methods were tried, but the results were unsatisfactory: (1) Exposed blasting method: The explosive charge was placed directly on the surface of the ore pile for blasting. A total of 300 kg of explosives and 150 detonators were consumed. Only more than 400 tons of ore were shaken off at a height of 3-4 meters above the bottom of the mining roadway. The bonded ore pile at higher levels was not loosened, and there was a significant safety risk in the subsequent placement of explosive charges. (2) Robotic hole-cutting blasting method: A robot was used to construct shallow holes with a diameter of 50-60 cm at the bottom of the ore pile. Explosive rolls were tied to bamboo poles and placed in the holes for blasting. Due to the difficulty of hole cutting and insufficient charge depth, only a small amount of ore was shaken off, and the top was suspended, posing a risk of collapse.

[0167] Universality and application prospects of the method: This embodiment takes the high-sulfur ore bonding recovery of a copper mine in Xinjiang as an example, but the scope of application of this method is not limited to this: (1) Different bonding causes: Whether it is sulfate cementation caused by sulfide oxidation, carbonate cementation caused by carbonate precipitation, or physical bonding caused by clay mineral hydration, this method can identify its causes through core diagnostic data and match the corresponding mechanical parameters, which has wide applicability. (2) Different mining methods: This method is not only applicable to the large-diameter deep hole subsequent backfilling mining method, but also to the bottom ore bonding problem in other mining methods such as segmented open field method and room and pillar method. (3) Different recovery operation mode: In addition to blasting vibration, the dynamic disturbance source in this method can also be replaced by different crushing modes such as high-pressure water jet, static expansion and cracking, and mechanical impact. Only by defining the corresponding disturbance mechanism and parameters in the model can the simulation optimization of different operation modes be realized. (4) Different geological conditions: The geostress background model in this method can be configured according to the measured geostress data of different mines, and is applicable to mining conditions with different burial depths and different tectonic stress fields. Therefore, the method provided in this application has high versatility and portability, and can be widely used in the recovery operation of bonded ore in various mines, and has good prospects for promotion and application.

[0168] All of the above-mentioned optional technical solutions can be combined in any way to form the optional embodiments of this application, and will not be described in detail here.

[0169] Figure 3 This is a schematic diagram of a mining area bonded ore recovery management system provided in an embodiment of this application. See also... Figure 3 The system includes one or more processors and one or more memories, in which at least one computer program is stored. The computer program is loaded and executed by the one or more processors to implement the above-mentioned method for managing the recovery of bonded ore in the mining area. The specific implementation process is detailed in the method embodiment and will not be repeated here.

[0170] The above are merely optional embodiments of this application and are not intended to limit this application. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of this application should be included within the protection scope of this application.

Claims

1. A method for managing the recovery of caking ore in a stope, characterized in that, The method includes: The three-dimensional laser scanning data, internal physical exploration data, and core genetic diagnosis data of the bonded ore pile are fused to obtain the spatial morphology and internal structural defect distribution of the bonded ore pile. Based on the differences in chemical and mineral composition indicated by the core genetic diagnosis data, multiple bonding type zones are defined, and the quantitative rock mass mechanical parameters corresponding to each bonding type zone are determined. Based on the spatial morphology, the distribution of internal structural defects, the multiple bonding type partitions, and the quantified rock mass mechanical parameters, the stress background model of the mining site is configured to obtain a heterogeneous coupled digital twin model. The heterogeneous coupled digital twin model is used to characterize the spatial variation mechanical properties of the bonded ore pile. Based on the mechanical properties of the multiple bond-type partitions, matching disturbance mechanisms are selected for different types of bond-type partitions. These disturbance mechanisms include blasting, high-pressure water jetting, or static fracturing. For each selected disturbance mechanism, initial spatiotemporal evolution parameters and energy release parameters are set based on the quantified rock mass mechanical parameters of the bond-type partition corresponding to the disturbance mechanism. The set initial spatiotemporal evolution parameters and energy release parameters are mapped to the corresponding spatial locations in the heterogeneous coupled digital twin model to complete the definition and configuration of the dynamic disturbance source. The heterogeneous coupled digital twin model is run to perform multi-field coupled simulation of the stress wave propagation, rock mass fracturing, and fragment transport processes triggered by the dynamic disturbance source. The simulation results are obtained, and the estimated values ​​of ore recovery indicators and safety risk indicators are extracted from the simulation results. Using the estimated values ​​of ore recovery indicators and safety risk indicators as the optimization objective function values, a multi-objective optimization algorithm is driven to adjust the spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source, generating an updated set of parameters. The dynamic disturbance source is reconfigured using the updated parameters, and the steps of running the heterogeneous coupled digital twin model, extracting the estimated values ​​of indicators from the simulation results, and driving the multi-objective optimization algorithm for adjustment are repeated until the preset optimization convergence criterion is met. The operation plan corresponding to the parameter set finally determined by the iteration is output as the target recovery operation plan. During the execution of the target recovery operation plan, the real-time collected construction feedback data is dynamically compared and assimilated with the predicted data at the corresponding location in the heterogeneous coupled digital twin model, and the operation parameters that have not yet been executed in the target recovery operation plan are dynamically adjusted based on the model assimilation results.

2. The method according to claim 1, characterized in that, The fusion processing of three-dimensional laser scanning data, internal physical detection data, and core genetic diagnostic data of the bonded ore pile yields the spatial morphology and internal structural defect distribution of the bonded ore pile, including: Based on the three-dimensional laser scanning data, a three-dimensional surface model is generated, which is used to characterize the external contour of the bonded ore pile. By fusing and analyzing the ground-penetrating radar data and the microseismic monitoring event sequence in the internal physical exploration data, the spatial distribution data of the internal defect elements of the bonded ore pile are obtained. Spatial registration and fusion are performed between the three-dimensional surface model and the spatial distribution data of the internal defect elements to obtain the spatial morphology and internal structural defect distribution of the bonded ore pile.

3. The method according to claim 2, characterized in that, The process of generating a three-dimensional surface model based on the three-dimensional laser scanning data includes: Geometric feature parameters are extracted from the three-dimensional laser scanning data. These geometric feature parameters are used to characterize the surface morphology of the bonded ore pile, including local curvature and surface roughness. Based on the distribution of the geometric feature parameters, the surface of the bonded ore pile is divided into different geometric feature zones; Based on the different geometric feature partitions, partition-adaptive surface reconstruction is performed on the three-dimensional laser scanning data to obtain the three-dimensional surface model.

4. The method according to claim 1, characterized in that, Based on the differences in chemical and mineral composition indicated by the genetic diagnostic data of the rock core, multiple bonding type zones are defined, and quantitative rock mass mechanical parameters corresponding to each bonding type zone are determined, including: The genetic diagnostic data of the core were interpreted to distinguish at least two dominant chemical or mineralogical modes. Based on the characteristic differences reflected in the genetic diagnostic data by the at least two dominant chemical or mineralogical genetic modes, the bonded ore pile is divided into multiple bonded type partitions corresponding to the at least two dominant genetic modes respectively. For each bonding type partition, based on the dominant chemical or mineralogical genetic mode corresponding to the bonding type partition, the corresponding quantitative rock mass mechanical parameters are retrieved from a pre-stored genetic mechanical parameter mapping table.

5. The method according to claim 1, characterized in that, Based on the spatial morphology, the distribution of internal structural defects, the multiple bonding type partitions, and the quantified rock mass mechanical parameters, the stress background model of the mining area is configured to obtain a heterogeneous coupled digital twin model, including: The spatial morphology is coupled with the distribution of internal structural defects to construct the geometric and structural entity boundary of the bonded ore pile in the stress background model of the mining site; Spatially assign values ​​to the multiple bonding type partitions and the quantified rock mass mechanical parameters, and define heterogeneous material mechanical properties for different regions within the geometric and structural entity boundaries; Based on the geometric and structural entity boundaries and the heterogeneous material mechanical properties, coupled calculations and initialization are performed on the stress background model of the mining site to obtain the heterogeneous coupled digital twin model for simulation.

6. The method according to claim 5, characterized in that, Based on the geometric and structural entity boundaries and the heterogeneous material mechanical properties, the stress background model of the mining site is coupled and initialized to obtain the heterogeneous coupled digital twin model for simulation derivation, including: Set the geometric and structural entity boundary as a new embedded geometric boundary in the mining site stress background model; The heterogeneous material mechanical properties are applied to the model region defined by the newly embedded geometric boundary; A contact relationship is established between the newly embedded geometric boundary and the surrounding rock region of the mining site stress background model. Based on the initial geostress conditions provided by the mining site stress background model and the contact relationship, iterative calculations are performed to solve the coupled stress balance between the inner and outer regions of the newly embedded geometric boundary, and a converged mechanical equilibrium solution is obtained. The mechanical equilibrium solution is set as the initial simulation condition for the heterogeneous coupled digital twin model, thus obtaining the heterogeneous coupled digital twin model.

7. The method according to claim 1, characterized in that, The heterogeneous coupled digital twin model is used to perform multi-field coupled simulations of the stress wave propagation, rock mass fracturing, and fragment transport processes triggered by the dynamic disturbance source, and the simulation results are obtained, including: In the heterogeneous coupled digital twin model, the spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source are converted into dynamic load boundary conditions acting on the corresponding region of the model. An explicit dynamic analysis algorithm is executed to solve the transient dynamic response of the heterogeneous coupled digital twin model under the dynamic load boundary conditions, and the continuous dynamic process is simulated, which includes stress wave propagation, rock mass fracturing, and fragment generation and migration. Record and output the spatiotemporal evolution data of key physical fields and the geometric information of fragments after rock mass failure during the continuous dynamic process, and use the spatiotemporal evolution data and the fragment geometric information as the simulation results.

8. The method according to claim 1, characterized in that, The method uses the estimated values ​​of the ore recovery index and the estimated values ​​of the safety risk index as the optimization objective function values ​​to drive a multi-objective optimization algorithm to adjust the spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source, generating a set of updated parameters, including: Using the estimated values ​​of the ore recovery index and the estimated values ​​of the safety risk index as inputs, a multi-objective fitness function is constructed to evaluate the overall performance of the scheme. During the initial execution, the spatiotemporal evolution parameters and energy release parameters of the dynamic disturbance source are used as the initial parameter set; based on the evaluation results of the initial parameter set by the multi-objective fitness function, the multi-objective optimization algorithm is driven to search within the preset parameter constraint space to generate and output the first set of candidate parameters. Determine whether the search process of the multi-objective optimization algorithm meets the preset convergence criterion; if it does, output the newly generated candidate parameter set as a set of updated parameters; if it does not, use the newly generated candidate parameter set as new input and repeat the steps of evaluation based on fitness function, driving algorithm search and determining convergence.

9. The method according to claim 1, characterized in that, The step of dynamically comparing and assimilating the real-time collected construction feedback data with the predicted data at the corresponding location in the heterogeneous coupled digital twin model includes: During construction, construction feedback data is collected in real time by sensors mounted on the working equipment. The construction feedback data includes at least drilling process parameters and vibration monitoring data. The construction feedback data is synchronously compared with the predicted mechanical response data of the same spatial location and the same construction stage in the heterogeneous coupled digital twin model to determine the deviation value between the two. Based on the deviation value, it is determined whether the model assimilation condition is triggered; if triggered, the mechanical property parameters of the local area in the heterogeneous coupled digital twin model are reverse-calibrated and updated according to the construction feedback data to complete the model assimilation and obtain the model assimilation result.

10. A management system for recovering ore bonded in a mining area, characterized in that, The system includes one or more processors and one or more memories, wherein at least one computer program is stored in the one or more memories, and the computer program is loaded and executed by the one or more processors to implement the mining area bonded ore recovery management method according to any one of claims 1-9.

Images