Layered simulation geologic modeling method under anisotropic tensor and uncertainty constraint

Abstract

The application discloses a hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraint, relates to the technical field of geological modeling, and comprises the following steps: inputting selected profile sample grade data, fitting a variation function and defining simulation parameters; performing conditional simulation according to the simulation parameters to generate a plurality of simulation implementations reflecting spatial anisotropy; controlling a hierarchical decomposition process through checking a termination condition, performing threshold smoothing processing on a sub-region after termination, saving an optimal threshold value and a class probability field of the sub-region, updating a global probability field by using a back transmission mechanism, outputting a grid of a final class probability field, and forming a rasterized ore body shape. The application obtains local threshold optimization through continuous hierarchical simulation and iteration, replaces unoptimized results, and forms optimal simulation results of the ore body, so that the uncertainty transmission of a traditional method is effectively reduced, and the accuracy and rationality of geological domain division are improved.

Description

Technical Field

[0001] This invention relates to the field of geological modeling technology, and in particular to a layered simulation geological modeling method under anisotropic tensor and uncertainty constraints. Background Technology

[0002] In geological modeling, accurate determination of orebody boundaries is crucial for resource estimation and mining design. Existing modeling methods typically rely on techniques such as geometric modeling, distance function modeling, geostatistical methods, and data-driven classification, using stochastic simulation, interpolation, and classification algorithms to estimate orebody boundaries.

[0003] Geometric modeling and distance function modeling methods determine boundary morphology by extracting isosurfaces or isodomains, but they lack the ability to handle uncertainties and struggle to quantify the reliability of boundaries. Geostatistical methods, especially conditional simulation, can generate multiple realizations and calculate the probability of orebody occurrence; however, existing methods often suffer from insufficient uncertainty control, boundary instability, and loss of local details. These shortcomings limit the reliability of orebody boundaries and their practical application in engineering. Summary of the Invention

[0004] In view of the aforementioned existing problems, the present invention is proposed.

[0005] Therefore, this invention provides a layered simulation geological modeling method under anisotropic tensor and uncertainty constraints, which solves the problems of uncertainty transmission and boundary instability in the determination of ore body boundaries in traditional methods.

[0006] To solve the above-mentioned technical problems, the present invention provides the following technical solution:

[0007] This invention provides a hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints. The method includes: inputting selected profile sample grade data, fitting a variogram function, and defining simulation parameters; performing conditional simulation based on the simulation parameters to generate multiple simulation implementations reflecting spatial anisotropy; rasterizing the simulation implementations into a unified grid to form a set of simulated grades, setting an initial threshold to calculate the class probability field for each pixel, and calculating the Shannon entropy field based on the class probability field; using a quadtree decomposition method to hierarchically divide the study domain into sub-regions, constructing an anisotropic positive definite matrix tensor based on the variogram function parameters, and calculating the norm of the gradient matrix of the class probability field under the positive definite tensor; for each sub-region, calculating the average classification uncertainty and norm-average gradient continuity based on the Shannon entropy field and norm, constructing a local energy function coupling uncertainty and continuity, and solving for the optimal threshold of the sub-region by minimizing the local energy function; controlling the hierarchical decomposition process by checking termination conditions, performing threshold smoothing on the terminated sub-regions, saving the optimal threshold and class probability field of the sub-regions, updating the global probability field using a backpropagation mechanism, and outputting the final class probability field grid to form a rasterized orebody morphology.

[0008] As a preferred embodiment of the hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints described in this invention, wherein: when the norm average gradient continuity of a sub-region is lower than a preset gradient continuity judgment threshold, a first penalty weight coefficient is assigned to the classification uncertainty term in the local energy function;

[0009] When the norm-average gradient continuity is higher than the preset gradient continuity judgment threshold, a second penalty weight coefficient is assigned to the classification uncertainty term; the first penalty weight coefficient is greater than the second penalty weight coefficient.

[0010] As a preferred embodiment of the hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints described in this invention, the norm-averaged gradient continuity is obtained by calculating the arithmetic mean of the spatial gradient of the probability field class within the subregion under the anisotropic tensor-induced norm.

[0011] The anisotropic tensor adjusts the gradient contribution weights based on the main continuous direction of the variogram function.

[0012] The gradient contribution weights indirectly trigger the dynamic allocation of the first and second penalty weight coefficients in the local energy function by determining the value of the norm-averaged gradient continuity.

[0013] As a preferred embodiment of the hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints described in this invention, the principal axis direction parameters of the variogram are converted into rotation matrices, and the major and minor axis range parameters in the variogram are converted into eigenvalue matrices through scaling mapping, and positive definite tensors are generated through matrix operations.

[0014] As a preferred embodiment of the hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints described in this invention, wherein: the optimal threshold of the sub-region is determined by minimizing the local energy function that couples classification uncertainty and gradient continuity;

[0015] Based on the optimal threshold of the sub-region, the class probability value of all pixels in the corresponding sub-region is recalculated, and the updated class probability value is sent back to the previous level.

[0016] As a preferred embodiment of the hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints described in this invention, the study domain hierarchy refers to the set of sub-regions formed by dividing the entire study area into a pyramid-like hierarchical structure using a quadtree decomposition method.

[0017] In the quadtree decomposition method, the current level of decomposition is terminated if any of the following conditions are met: the local energy function value is lower than the preset energy convergence threshold, the energy improvement is less than the preset tolerance, the number of effective samples in the sub-region is less than the preset sample number threshold, and the decomposition level reaches the maximum pyramid setting level.

[0018] As a preferred embodiment of the hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints described in this invention, the threshold smoothing process refers to interpolating the overlapping region by using the weighted average of the optimal thresholds of adjacent sub-regions.

[0019] As a preferred embodiment of the hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints described in this invention, the backpropagation mechanism replaces the optimized sub-region class probability field with the parent class probability field at the pixel level for updating and transferring probability values ​​between layers.

[0020] As a preferred embodiment of the hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints described in this invention, the determination of the rasterized ore body morphology is based on the comparison result of the class probability field and the preset deterministic threshold. If the class probability value of the class probability field is greater than the preset deterministic threshold for a pixel, it is identified as an ore body region.

[0021] As a preferred embodiment of the hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints described in this invention, the calculation of the average classification uncertainty and the norm average gradient continuity are respectively performed by the arithmetic mean method, and the entropy value and gradient norm value of all pixels in the sub-region are statistically averaged.

[0022] The beneficial effects of this invention are as follows: it utilizes original data conditions to simulate multiple solutions to construct a probability field of probability class, an entropy field of uncertainty, and an anisotropic tensor; it uses pyramid layering for local partitioning, balancing the continuity of the uncertainty entropy field and the gradient of the anisotropic tensor to determine the local orebody boundary threshold; through continuous layered simulation and iteration, it obtains local threshold optimization, replaces unoptimized results, and forms the optimal simulation result of the orebody, thereby effectively reducing the uncertainty transmission of traditional methods and improving the accuracy and rationality of geological domain partitioning. Attached Figure Description

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

[0024] Figure 1 The flowchart shows a hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints.

[0025] Figure 2 This is a diagram of the simulation test process.

[0026] Figure 3 The diagram shows the local energy function optimization process during simulation.

[0027] Figure 4 This is a comparison chart of simulation results from simulation methods and other methods with real production data.

[0028] Figure 5 This is a comparison chart of the statistical characteristics of simulation results from the simulation method and other methods. Detailed Implementation

[0029] To make the above-mentioned objects, features and advantages of the present invention more apparent and understandable, the specific embodiments of the present invention will be described in detail below with reference to the accompanying drawings.

[0030] Many specific details are set forth in the following description in order to provide a full understanding of the invention. However, the invention may also be practiced in other ways different from those described herein, and those skilled in the art can make similar extensions without departing from the spirit of the invention. Therefore, the invention is not limited to the specific embodiments disclosed below.

[0031] Secondly, the term "one embodiment" or "embodiment" as used herein refers to a specific feature, structure, or characteristic that may be included in at least one implementation of the present invention. The phrase "in one embodiment" appearing in different places in this specification does not necessarily refer to the same embodiment, nor is it a single or selective embodiment that is mutually exclusive with other embodiments.

[0032] Reference Figures 1-4 As one embodiment of the present invention, this embodiment provides a layered simulation geological modeling method under anisotropic tensor and uncertainty constraints, comprising the following steps:

[0033] The sample has a strong anisotropic structure and indistinct boundary grade transitions, and simultaneously exhibits fracture-controlled abrupt boundary and alteration-controlled gradual boundary.

[0034] like Figure 1 As shown, this invention proposes a hierarchical simulation geological modeling method based on anisotropic tensor constraints and uncertainty optimization, comprising the following steps:

[0035] Input the grade data of the selected profile sample, fit the original variogram function, and set the required parameters for simulation based on the variogram, including simulation type, search neighborhood, realization number, etc.

[0036] Multiple implementations are generated using conditional randomization simulation (50 implementations are used in this embodiment). Each implementation maintains local statistical characteristics and continuous anisotropy direction under the constraint of the fitted variogram function.

[0037] The profile area is standardized to form a rectangle S. Each grid is rasterized to a uniform resolution, and these grids are then stacked to form a simulated grade set. Using the boundary grade as the initial threshold, the frequency of each pixel in the simulation set that is higher than the threshold is statistically analyzed to obtain the pixel-level ore probability.

[0038] ;

[0039] In the formula, Let k be a specific point in the rasterization, and k be a specific implementation in N sets of simulation implementations. For this realization The taste of the location The threshold for classifying pixels as either inside or outside the ore body. This indicates the probability that the location is identified as an ore body across multiple implementations.

[0040] Based on pixel-class probability Calculate Shannon entropy The expression is:

[0041] ;

[0042] In the formula, for Point threshold Entropy under certain conditions.

[0043] like Figure 1As shown, to better accommodate the uncertainty and continuity of pixel calculations in profile rasterization, the non-perpendicular principal and secondary axis vectors of the variogram are mapped from the natural coordinate system of the geological body to the standard Cartesian coordinate system, thus constructing an anisotropic tensor. The expression is:

[0044] ;

[0045] In the formula, 'a' represents the range of the variogram in the principal direction (the direction of strongest change), 'b' represents the range of the secondary direction (the direction of weakest change), and the rotation matrix R is used to rotate and align the principal and secondary axes of the variogram to the Cartesian coordinate system, making their directions consistent without changing their length. The diagonal matrix... Extract the standardized length weights of the anisotropic vectors along the Cartesian coordinate axes to describe the relative intensity of change along each coordinate axis. Tensor A is the positive definite matrix of R, representing the weights of the variogram function along different directions in space.

[0046] The spatial gradient matrix is ​​calculated from the probability field of ore, and the local direction and magnitude of the probability change are obtained.

[0047] To reflect the property connectivity between a single point and its neighboring points within a sub-region, maintain the continuity of the simulation image, and control the occurrence of noise and abrupt changes, a property gradient function is constructed. The expression is:

[0048] ;

[0049] In the formula, , Indicates two directions, Represents the geological properties at a specific pixel or node in space. (In this case, grade g) the rate of change along both directions;

[0050] Then use tensors We obtain the anisotropy norms of each term of the variogram by weighting the gradient matrix, i.e.:

[0051] ;

[0052] In the formula, the gradient matrix is ​​adjusted using a positive definite matrix A. This makes the geological properties The rate of change is consistent with the anisotropic characteristics of the geological body, that is, the change is smaller along the main direction (the direction of the maximum range) and larger along the secondary direction (the direction of the small range);

[0053] For the input raw data grid The study domain is divided into pyramid-shaped layers using quadtree decomposition, with the layers progressing from coarse to fine: firstly, it is divided into four overlapping sub-regions. , where i represents the number of layers in the iteration, and j represents the number of the four rectangular sub-regions contained in each layer; adjacent sub-regions are kept to overlap by 20% of their side length to avoid threshold mutation.

[0054] For the region Shannon entropy of all points within the range Summation is performed to obtain sub-regions. Total Entropy The expression is:

[0055] ;

[0056] Calculate the average class uncertainty for each sub-region The expression is:

[0057] ;

[0058] In the formula, sub-region The area of ​​the region, sub-region Total entropy;

[0059] Gradient continuity is obtained by summing and averaging the anisotropic norms of all points within each sub-region. The expression is:

[0060] ;

[0061] In the formula, sub-region The area of ​​the region, Let be the anisotropic norm.

[0062] For the current level sub-region Constructing an energy function for coupled uncertainties and gradient continuity The expression is:

[0063] ;

[0064] In the formula, For the class of uncertainty calculated based on Shannon entropy for sub-regions, To ensure gradient continuity based on anisotropic norm calculations, Assigning penalty weights to the uncertainties in the classification;

[0065] The optimization process of the energy function involves finding the optimal threshold for the sub-region by minimizing the local energy function, such as... Figure 3 As shown.

[0066] ;

[0067] In the formula, argmin represents the expression that makes the subregion The independent variable that minimizes the energy function The value of is the optimal grade threshold for classifying pixel x as either inside or outside the ore body, taking into account both the smoothness and uncertainty of the sub-region.

[0068] Another optimization involves assigning a first penalty weight coefficient to the classification uncertainty term in the local energy function when the average gradient continuity of the sub-region is lower than a preset gradient continuity threshold. When the norm-average gradient continuity is higher than the preset gradient continuity threshold, a second penalty weight coefficient is assigned to the classification uncertainty term; the first penalty weight coefficient is greater than the second penalty weight coefficient.

[0069] ;

[0070] In the formula, The preset threshold for judging gradient continuity; The first penalty weighting coefficient, This is the second penalty weighting coefficient.

[0071] It should be noted that the preset gradient continuity judgment threshold is determined by taking the statistical quantile value of the norm average gradient continuity of all sub-regions in the current study domain. It is preferred to take the median or a value in the normalized interval of 0.3 to 0.7 as an example value.

[0072] Both methods use gradient continuity to set thresholds. It was determined that when the norm-average gradient continuity of a sub-region is lower than a preset gradient continuity threshold, It is assigned a large value (e.g., ranging from 0.5 to 1.0), and when the norm-averaged gradient continuity is higher than the gradient continuity judgment threshold, Smaller (e.g., ranging from 0.1 to 0.3); specific weight coefficient values ​​can be adjusted based on the main continuity direction of the anisotropic tensor of the gradient continuity and uncertainty level in the current subregion, according to the variogram function.

[0073] By using weight coefficients to determine the value of the norm-averaged gradient continuity, the dynamic allocation of the first and second penalty weight coefficients in the local energy function is indirectly triggered.

[0074] The two-level optimization process terminates if any of the following stopping conditions are met: a. Energy value is below a threshold; b. Energy improvement is negligible; c. Insufficient sample size in the sub-region; d. Reaching the maximum pyramid level; e. Threshold deviates too much from economic feasibility level. If any termination condition is met, the process proceeds to the threshold smoothing step. If the termination condition is not met, the process returns to the quadtree decomposition step, continuing to decompose the current sub-region and repeating the subsequent optimization process.

[0075] Specifically, "energy value below the threshold" means that the local energy function value of the current sub-region is lower than the energy convergence threshold set proportionally to the initial value of the global energy function during the initialization phase, preferably 1% to 5% of the initial energy value as the example value range. "Insufficient number of samples in the sub-region" means that the number of valid samples participating in the statistics in the sub-region is less than the sample number threshold preset according to the overall sample size at the beginning of modeling, preferably no less than 10 to 30 valid samples as the example value range.

[0076] After the decomposition terminates, the overlapping regions of the last sub-region are smoothed using the average of the optimal thresholds of adjacent sub-regions to ensure the spatial continuity of the threshold field. Figure 3 As shown.

[0077] Save the optimal threshold and the corresponding class probability for each sub-region. Replace the previous level's class probability field with the final optimized sub-region class probability field to achieve backpropagation and transmission control of uncertainty between levels. Output a grid containing the final class probabilities, and form a raster-based geological domain model of the ore body boundary with the set of pixels whose determinism is greater than a set lower limit.

[0078] In the process of forming the geological domain model of the ore body boundary, the class probability value of each pixel in the final class probability field is compared with the deterministic threshold (i.e., the lower limit) preset in the modeling initialization stage based on the global class probability field statistical distribution. If the class probability value is greater than the deterministic threshold, the corresponding pixel is identified as the ore body region. The deterministic threshold is preferably taken as an example value range of 0.5 to 0.8.

[0079] like Figure 4 The figure shows a comparison between the simulation results of the present invention and the implicit modeling method of the distance function, the manual ore delineation ordinary kriging method, and the a posteriori boundary of supplementary drilling. The geological boundary modeling is more realistic and does not exhibit excessive smoothing; in areas with faults and low-ore interlayers, the simulation reveals boundary changes through high entropy. In stable and gradually changing regions of the ore body, the high continuity of the boundary simulation is consistent with the ore body structure.

[0080] like Figure 5As shown, the boundary of this simulation, the boundary of the distance function implicit modeling method, and the boundary of the manual mineral delineation ordinary Kriging method are selected and compared with the posterior boundary of the supplementary drilling. The intersection-union ratio, overlap rate, and Hausdorf distance index are calculated. The present invention shows a higher overlap rate and a smaller boundary deviation.

[0081] In summary, this invention utilizes original data conditions to simulate multiple solutions to construct a probability field of probability class, an entropy field of uncertainty, and an anisotropic tensor; it uses pyramid layering for local partitioning, balancing the continuity of the uncertainty entropy field and the gradient of the anisotropic tensor to determine the local orebody boundary threshold; through continuous layered simulation and iteration, it obtains local threshold optimization, replaces unoptimized results, and forms the optimal simulation result of the orebody, thereby effectively reducing the uncertainty transmission of traditional methods and improving the accuracy and rationality of geological domain partitioning.

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

Claims

1. A hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints, characterized in that: include, Input the grade data of the selected profile sample, fit the variogram function and define the simulation parameters; Based on the simulation parameters, conditional simulations are performed to generate multiple simulation implementations that reflect spatial anisotropy. The simulation is rasterized into a unified grid to form a set of simulated grades. An initial threshold is set to calculate the class probability field of each pixel, and the Shannon entropy field is calculated based on the class probability field. The study domain is hierarchically divided into sub-regions using the quadtree decomposition method. An anisotropic positive definite matrix tensor is constructed based on the variogram parameters. The norm of the gradient matrix of the probability field class under the positive definite tensor is calculated. For each sub-region, the average classification uncertainty and norm-average gradient continuity are calculated based on the Shannon entropy field and norm. A local energy function that couples uncertainty and continuity is constructed, and the optimal threshold of the sub-region is solved by minimizing the local energy function. By checking the termination condition to control the hierarchical decomposition process, threshold smoothing is performed on the terminated sub-regions, the optimal threshold and class probability field of the sub-regions are saved, and the global probability field is updated using a backpropagation mechanism to output the grid of the final class probability field, forming a rasterized ore body morphology.

2. The hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints as described in claim 1, characterized in that: When the norm-average gradient continuity of a sub-region is lower than the preset gradient continuity judgment threshold, a first penalty weight coefficient is assigned to the classification uncertainty term in the local energy function. When the norm-average gradient continuity is higher than the preset gradient continuity judgment threshold, a second penalty weight coefficient is assigned to the classification uncertainty term; the first penalty weight coefficient is greater than the second penalty weight coefficient.

3. The hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints as described in claim 1, characterized in that: The continuity of the norm-averaged gradient is obtained by calculating the arithmetic mean of the spatial gradients of the probability field class within the subregion under the anisotropic tensor-induced norm; The anisotropic tensor adjusts the gradient contribution weights based on the main continuous direction of the variogram function. The gradient contribution weights indirectly trigger the dynamic allocation of the first and second penalty weight coefficients in the local energy function by determining the value of the norm-averaged gradient continuity.

4. The hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints as described in claim 1, characterized in that: The principal axis direction parameter of the variogram is converted into a rotation matrix, and the major and minor axis range parameters in the variogram are converted into eigenvalue matrices through scaling mapping, and positive definite tensors are generated through matrix operations.

5. The hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints as described in claim 1, characterized in that: The optimal threshold for the sub-region is determined by minimizing the local energy function that minimizes the coupling classification uncertainty and gradient continuity. Based on the optimal threshold of the sub-region, the class probability value of all pixels in the corresponding sub-region is recalculated, and the updated class probability value is sent back to the previous level.

6. The hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints as described in claim 1, characterized in that: The research domain hierarchy refers to the set of sub-regions formed by dividing the entire research area into a pyramid-shaped hierarchical structure using the quadtree decomposition method. In the quadtree decomposition method, the current level of decomposition is terminated if any of the following conditions are met: the local energy function value is lower than the preset energy convergence threshold, the energy improvement is less than the preset tolerance, the number of effective samples in the sub-region is less than the preset sample number threshold, and the decomposition level reaches the maximum pyramid setting level.

7. The hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints as described in claim 1, characterized in that: The threshold smoothing process refers to interpolating the overlapping region using a weighted average of the optimal thresholds of adjacent sub-regions.

8. The hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints as described in claim 1, characterized in that: The backhaul mechanism replaces the optimized sub-region class probability field with the parent class probability field at the pixel level, which is used for updating and transferring probability values ​​between layers.

9. The hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints as described in claim 1, characterized in that: The determination of the rasterized ore body morphology is based on the comparison result between the class probability field and the preset deterministic threshold. If the class probability value of the class probability field is greater than the preset deterministic threshold for a pixel, it is identified as an ore body region.

10. The hierarchical simulation geological modeling method under anisotropic tensor and uncertainty constraints as described in claim 3, characterized in that: The calculation of the average classification uncertainty and the norm average gradient continuity are respectively performed using the arithmetic mean method, which statistically averages the entropy value and gradient norm value of all pixels in the sub-region.

Images