An Adaptive Mesh Generation Method and System for Thin Reservoir Planar Mesh Based on Property Gradient

Abstract

This invention relates to the field of oil and gas reservoir geological modeling, and provides an adaptive meshing method and system for thin reservoir planar grids based on attribute gradients. The method includes: preprocessing seismic attribute surface data of the work area and constructing a planar attribute field; calculating the plane attribute gradient magnitude of the planar attribute field using the finite difference method to obtain the gradient field; mapping the gradient field to the center of corner grid cells using a bilinear interpolation algorithm to obtain the grid gradient field; calculating the statistical characteristic values ​​of the gradient magnitudes across the entire work area in the grid gradient field and generating an adaptive hierarchical threshold; automatically partitioning and labeling the grid cells of the entire work area according to the gradient magnitude to obtain the partitioning and labeling results; and performing hierarchical densification and smoothing optimization on the grid cells corresponding to the transition zones and abrupt change zones in the partitioning and labeling results to obtain an adaptively densified corner grid model. This invention improves the geometric quality and topological compliance of the corner grid model of thin reservoir planar grids.

Description

Technical Field

[0001] This invention relates to the field of oil and gas reservoir geological modeling and numerical simulation technology, and in particular to an adaptive meshing method and system for thin reservoir planar grids based on attribute gradients. Background Technology

[0002] In oil and gas reservoir exploration and development, thin reservoirs have become an important research target for increasing reserves and production due to their large reserves and wide distribution. Thin reservoirs generally have characteristics such as small single-layer thickness, rapid lateral facies changes, complex spatial distribution, strong heterogeneity, and well-developed faults and pinch-outs. Their geological modeling and numerical simulation place extremely high demands on grid partitioning: it is necessary to achieve detailed characterization in areas of abrupt property changes, sand body boundaries, near faults, and pinch-out zones, while controlling the number of grids in gentle facies zones and stable regions to ensure the efficiency of numerical simulation calculations.

[0003] The current mainstream grid for geological modeling is the corner grid. Although it can adapt to complex structures and faults, planar grid subdivision technology still has significant drawbacks. Existing technologies mainly rely on manual delineation of the densification area, depending on the experience of interpreters, resulting in low efficiency and poor consistency. Traditional uniform grid technology cannot achieve a balance between accuracy and efficiency. Globally uniform coarse grids cannot characterize subtle changes inside thin reservoirs, while globally uniform densification grids, although improving accuracy, cause the number of grids to increase exponentially, significantly increasing computational resource consumption. At the same time, existing adaptive methods are mostly oriented towards general grid design, incompatible with the unique characteristics of corner grids, and lack adaptive subdivision schemes specifically for thin reservoirs. Summary of the Invention

[0004] This invention provides an adaptive meshing method and system for thin reservoir planar grids based on attribute gradients, in order to overcome the shortcomings of existing technologies.

[0005] This invention provides an adaptive mesh generation method for thin reservoir planar grids based on attribute gradients, comprising: S1. Preprocess the seismic attribute surface data of the work area and construct a planar attribute field; S2. The gradient field is obtained by calculating the gradient magnitude of the planar attribute field using the finite difference method. S3. The gradient field is mapped to the center of the corner grid cell using a bilinear interpolation algorithm to obtain the grid gradient field; S4. Calculate the statistical characteristic value of the gradient magnitude of the entire working area in the grid gradient field, and generate an adaptive hierarchical threshold by weighting the statistical characteristic value; S5. Based on the adaptive hierarchical threshold, the grid cells of the entire work area are automatically partitioned and labeled according to the gradient magnitude to obtain the partitioning and labeling results. S6. Perform hierarchical densification and smoothing optimization on the grid cells corresponding to the transition zone and steep change zone in the partitioning marking results to obtain the adaptive densified corner grid model, so as to complete the adaptive partitioning of the thin reservoir planar grid.

[0006] According to the adaptive meshing method for thin reservoir planar grids based on attribute gradients provided by the present invention, step S2 further includes: S21. The partial derivatives of the planar attribute field are calculated in the first horizontal direction and the second horizontal direction using the finite difference method to obtain the partial derivative values ​​in the corresponding directions; wherein the finite difference method includes forward finite difference, backward finite difference, and central finite difference method; S22. Perform a square superposition operation and a square root operation on the partial derivative values ​​in the first horizontal direction and the partial derivative values ​​in the second horizontal direction to obtain a gradient field based on the gradient magnitude.

[0007] According to the adaptive meshing method for thin reservoir planar grids based on attribute gradients provided by the present invention, step S3 further includes: S31. Traverse each grid cell in the corner grid one by one to obtain the center coordinates of the corresponding grid cell; S32. Based on the position of the center coordinates of the mesh cell in the attribute surface rectangular mesh, find and determine the attribute surface mesh cell that contains the center point of the corresponding mesh cell; S33. Based on the coordinates of the four corner points of the attribute surface grid cell, the gradient magnitude at the four corner points of the attribute surface grid cell is mapped to the cell center of the corresponding grid cell through a bilinear interpolation function. After assigning gradient values ​​to all grid cells in the corner grid, the gradient field is mapped from the attribute surface grid to the corner grid, and the grid gradient field is obtained.

[0008] According to the adaptive meshing method for thin reservoir planar grids based on attribute gradients provided by the present invention, in step S33, the expression of the bilinear interpolation function is:

[0009] in, The x-coordinate value of the center coordinate of the grid cell. The ordinate value of the center coordinate of the grid cell. Center of corner grid cell The gradient magnitude is calculated at [location]. The minimum coordinate value of the attribute surface mesh cell in the first horizontal direction. The maximum coordinate value of the attribute surface mesh cell in the first horizontal direction. This represents the minimum coordinate value of the attribute surface mesh cell in the second horizontal direction. This represents the maximum coordinate value of the attribute surface mesh cell in the second horizontal direction. This represents the gradient magnitude at the bottom left corner of the attribute surface mesh cell. This represents the gradient magnitude at the bottom right corner of the attribute surface mesh cell. This represents the gradient magnitude at the top-left corner of the attribute surface mesh cell. This represents the gradient magnitude at the top right corner of the attribute surface mesh cell. The area weight of the bottom left corner point of the attribute surface mesh cell. The area weight of the bottom right corner point of the attribute surface mesh cell. The area weight of the top-left corner point of the attribute surface mesh cell. This represents the area weight of the top-right corner point of the attribute surface mesh cell.

[0010] According to the adaptive meshing method for thin reservoir planar grids based on attribute gradients provided by the present invention, step S4 further includes: S41. Perform statistical analysis on the gradient magnitude of all grid cells in the grid gradient field, and calculate the mean and standard deviation of the gradient magnitude of the entire work area. S42. Calculate a baseline threshold based on gradient statistical features, wherein the baseline threshold is used to filter out background fluctuation areas; S43. Introduce a first empirical coefficient and a second empirical coefficient, and scale the benchmark threshold proportionally to obtain an adaptive hierarchical threshold including a first-level threshold and a second-level threshold.

[0011] According to the adaptive meshing method for thin reservoir planar grids based on attribute gradient provided by the present invention, in step S43, the first empirical coefficient is 0.6 and the second empirical coefficient is 1.2.

[0012] According to the adaptive meshing method for thin reservoir planar grids based on attribute gradients provided by the present invention, step S5, the partitioning labeling result specifically includes: For grid cells whose gradient magnitude is less than or equal to the first threshold, they are marked as smooth regions; Mesh cells with gradient magnitudes greater than the first threshold and less than the second threshold are marked as transition regions; For grid cells with gradient magnitudes greater than or equal to the second threshold, they are marked as steep transition regions.

[0013] According to the adaptive meshing method for thin reservoir planar grids based on attribute gradients provided by the present invention, step S6 further includes: S61: Perform hierarchical densification and smoothing optimization on the grid cells corresponding to the transition and steep change regions in the partitioning marking results to obtain a densified grid; S62. For the encrypted mesh, the Laplacian smoothing algorithm is used to iteratively optimize the node coordinates within the encrypted mesh area to obtain the corner mesh model.

[0014] According to the present invention, a thin reservoir planar mesh adaptive partitioning method based on attribute gradient is provided, wherein step S61 specifically includes: For the marked grid cells in the transition zone, a first-level densification is performed in the first grid subdivision direction and the second grid subdivision direction respectively, and the first-level densification is subdivided by a densification factor of 2. For the marked grid cells in the steep transition zone, secondary densification is performed in the first grid subdivision direction and the second grid subdivision direction respectively, and the secondary densification is subdivided by a densification factor of 4.

[0015] The present invention also provides an adaptive mesh generation system for thin reservoir planar grids based on attribute gradients, for performing an adaptive mesh generation method for thin reservoir planar grids based on attribute gradients as described in any of the preceding claims, comprising: Preprocessing module: Used to preprocess the seismic attribute surface data of the work area and construct a planar attribute field; Calculation module: used to calculate the gradient magnitude of the planar attribute field using the finite difference method to obtain the gradient field; Mapping module: Used to map the gradient field to the center of the corner grid cell using a bilinear interpolation algorithm to obtain the grid gradient field; Generation module: used to calculate the statistical characteristic values ​​of the gradient magnitude of the entire work area in the grid gradient field, and to generate an adaptive hierarchical threshold by weighting the statistical characteristic values; The marking module is used to automatically partition and mark the grid cells of the entire work area according to the gradient magnitude based on the adaptive grading threshold, and obtain the partition marking results. Encryption module: used to perform hierarchical encryption and smoothing optimization on the grid cells corresponding to the transition zone and steep change zone in the partitioning marking results, to obtain the adaptively encrypted corner grid model, so as to complete the adaptive partitioning of thin reservoir planar grid. Attached Figure Description

[0016] To more clearly illustrate the technical solutions in this invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are some embodiments of this invention. For those skilled in the art, other drawings can be obtained from these drawings without creative effort.

[0017] Figure 1 This is a schematic diagram of a thin reservoir planar mesh adaptive partitioning method based on attribute gradient provided in an embodiment of the present invention; Figure 2A schematic diagram showing the comparison of meshes before and after adaptive densification, based on an adaptive thin reservoir planar mesh subdivision method based on attribute gradients, provided in an embodiment of the present invention. Figure 3 This is a schematic diagram of an adaptive subdivision system for thin reservoir planar grids based on attribute gradients, provided in an embodiment of the present invention. Detailed Implementation

[0018] To make the objectives, technical solutions, and advantages of this invention clearer, the technical solutions of this invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, embodiments of this invention, and should not be construed as limiting the invention. All other embodiments obtained by those skilled in the art based on the embodiments of this invention without creative effort are within the scope of protection of this invention. In the description of this invention, it should be understood that the terminology used is for descriptive purposes only and should not be construed as indicating or implying relative importance.

[0019] The embodiments of the present invention are described below with reference to the figures.

[0020] like Figure 1 As shown, this invention provides an adaptive meshing method for thin reservoir planar grids based on attribute gradients, comprising: S1. The seismic attribute surface data of the work area is preprocessed and a planar attribute field is constructed.

[0021] In step S1, this invention first imports multi-source data such as seismic attribute surfaces and sand thickness surfaces from the work area. Because the acquired raw data suffers from problems such as inconsistent coordinate systems, noise interference, and differences in numerical dimensions, this invention performs coordinate unification processing on these raw attribute data, unifying attribute surfaces from different sources to the same coordinate system. Subsequently, denoising processing is performed to eliminate outliers and high-frequency noise interference in the attribute surfaces. Finally, all attribute data are normalized, converting each attribute value to a unified numerical range, allowing attribute data with different dimensions to undergo unified mathematical operations. After preprocessing, a planar attribute field is formed, which is expressed as a functional relationship in two-dimensional space: ,in Indicates the position in planar coordinates. This indicates the attribute value at that location (such as earthquake amplitude, sand thickness, etc.).

[0022] S2. The gradient field is obtained by calculating the gradient magnitude of the planar attribute field using the finite difference method.

[0023] Step S2 further includes: S21. The partial derivatives of the planar attribute field are calculated in the first horizontal direction and the second horizontal direction using the finite difference method to obtain the partial derivative values ​​in the corresponding directions; wherein the finite difference method includes forward finite difference, backward finite difference and central finite difference method.

[0024] In step S21, the present invention applies the attribute field in the X and Y directions respectively. Finite difference operations are performed, and the finite difference method includes three forms: forward difference, backward difference, and central difference. These three forms can be used interchangeably. Taking central difference as an example, the partial derivative in the X direction is:

[0025] The partial derivative in the Y direction is:

[0026] in and Let be the step size of the attribute surface mesh in the corresponding direction. In this specific implementation, the basic mesh step size is taken as 60m. Furthermore, forward and backward differencing are used at node locations on the attribute surface boundary where two-sided data cannot be obtained. After processing in step S21, this invention obtains two independent partial derivative values ​​in the X and Y directions at each node of the attribute surface mesh, denoted as […]. and .

[0027] S22. Perform a square superposition operation and a square root operation on the partial derivative values ​​in the first horizontal direction and the partial derivative values ​​in the second horizontal direction to obtain a gradient field based on the gradient magnitude.

[0028] In step S22, after obtaining the partial derivative values ​​in the two directions, the present invention squares the partial derivative values ​​in the X direction and the Y direction respectively, then adds the two squared values, and then performs a square root operation on the summation to obtain the gradient magnitude at the node. The calculation formula is:

[0029] The gradient magnitude calculated using the above formula This reflects the degree of drastic change of the attribute at that spatial location; a larger value indicates a more drastic change, while a smaller value indicates a more gradual change. Furthermore, this invention performs the aforementioned square superposition and square root operations on each node in the attribute surface mesh, ultimately forming a complete gradient field across the entire attribute surface.

[0030] S3. The gradient field is mapped to the center of the corner grid cell using a bilinear interpolation algorithm to obtain the grid gradient field.

[0031] Step S3 further includes: S31. Traverse each grid cell in the corner grid one by one to obtain the center coordinates of the corresponding grid cell.

[0032] In step S31, the present invention iterates through each grid cell in the corner grid one by one, and obtains the planar coordinates of the center point of each grid cell, that is, the center position coordinates of the grid cell in the first horizontal direction and the second horizontal direction, as the target point coordinate input for subsequent interpolation calculation.

[0033] S32. Based on the position of the center coordinates of the grid cell in the attribute surface rectangular grid, find and determine the attribute surface grid cell that contains the center point of the corresponding grid cell.

[0034] For the center coordinates of each corner grid cell obtained in S31 This invention performs spatial location lookup within an attribute surface rectangular grid to determine which attribute surface grid cell spatially contains the center coordinate point. The determined attribute surface grid cell is defined by the coordinates of its four corner points. Let the coordinate range of this attribute surface grid cell in the X direction be... The coordinate range in the Y direction is Then when the target point coordinates satisfy and When this occurs, the attribute surface mesh cell is the target cell.

[0035] S33. Based on the coordinates of the four corner points of the attribute surface grid cell, the gradient magnitude at the four corner points of the attribute surface grid cell is mapped to the cell center of the corresponding grid cell through a bilinear interpolation function. After assigning gradient values ​​to all grid cells in the corner grid, the gradient field is mapped from the attribute surface grid to the corner grid, and the grid gradient field is obtained.

[0036] In step S33, the expression for the bilinear interpolation function is:

[0037] in, The x-coordinate value of the center coordinate of the grid cell. The ordinate value of the center coordinate of the grid cell. Center of corner grid cell The gradient magnitude is calculated at [location]. The minimum coordinate value of the attribute surface mesh cell in the first horizontal direction. The maximum coordinate value of the attribute surface mesh cell in the first horizontal direction. This represents the minimum coordinate value of the attribute surface mesh cell in the second horizontal direction. This represents the maximum coordinate value of the attribute surface mesh cell in the second horizontal direction. This represents the gradient magnitude at the bottom left corner of the attribute surface mesh cell. This represents the gradient magnitude at the bottom right corner of the attribute surface mesh cell. This represents the gradient magnitude at the top-left corner of the attribute surface mesh cell. This represents the gradient magnitude at the top right corner of the attribute surface mesh cell. The area weight of the bottom left corner point of the attribute surface mesh cell. The area weight of the bottom right corner point of the attribute surface mesh cell. The area weight of the top-left corner point of the attribute surface mesh cell. This represents the area weight of the top-right corner point of the attribute surface mesh cell.

[0038] In step S33, after determining the attribute surface mesh cell containing the target point, the present invention extracts the gradient magnitudes at the four corner points of the attribute surface mesh cell, namely the gradient magnitude of the lower left corner point. Gradient magnitude at the bottom right corner Gradient magnitude at the top left corner Gradient magnitude at the top right corner Subsequently, this invention calculates the area weights of the target point relative to the four corner points, with the area weights of the four corner points being as follows:

[0039]

[0040]

[0041]

[0042] Wherein, the numerator of each weight coefficient is the target point. The area of ​​the rectangle enclosed by the diagonal corners is the denominator of the total area of ​​the attribute surface mesh cells. The sum of the four weights is always equal to 1, and the closer the target point is to a corner point, the larger the corresponding weight. Subsequently, this invention multiplies the gradient magnitudes of the four corner points by their corresponding area weights and then sums them, i.e., it calculates the gradient magnitude at the target point using a bilinear interpolation function, as shown in the formula above. The calculated gradient magnitude is... This is the gradient magnitude at the center of the corner grid cell.

[0043] Furthermore, the present invention performs S31 to S33 on all grid cells in the corner grid, and after assigning values ​​one by one, finally obtains a grid gradient field that perfectly matches the corner grid.

[0044] S4. Calculate the statistical characteristic value of the gradient magnitude of the entire working area in the grid gradient field, and generate an adaptive hierarchical threshold by weighting the statistical characteristic value.

[0045] Step S4 further includes: S41. Perform statistical analysis on the gradient magnitude of all grid cells in the grid gradient field, and calculate the mean and standard deviation of the gradient magnitude of the entire work area.

[0046] In step S41, the present invention receives the grid gradient field output from S3 as input data, wherein each corner grid cell in the grid gradient field has been assigned a corresponding gradient magnitude. Subsequently, this invention applies to all areas within the entire work zone. Statistical analysis was performed on the gradient magnitudes of each grid cell. First, all... Gradient magnitude Add them together and divide by the total number of grid cells. The mean gradient magnitude of the entire work area was calculated. ; then with Using this as a baseline, calculate the gradient magnitude for each grid cell. with the mean The square of the difference will be all Sum of the squares of the differences and divide by Take the square root again to obtain the standard deviation of the gradient magnitude across the entire work area. Mean This reflects the overall level and standard deviation of the gradient amplitude across the entire work area. It reflects the degree of dispersion of the gradient magnitude of the entire work area relative to the mean, and the two together constitute the statistical feature input for subsequent threshold calculation.

[0047] S42. Calculate a baseline threshold based on gradient statistical features. The baseline threshold is used to filter out background fluctuation areas.

[0048] In obtaining the mean and standard deviation Then, the present invention directly adds the two together to calculate the benchmark threshold. The obtained baseline threshold In a physical sense, this represents the position one standard deviation above the mean, that is, the gradient magnitude in the entire work area is higher than... The grid cells belong to regions with relatively drastic changes, below The grid cells belong to the background fluctuation region. Because... It is calculated directly from the statistical characteristics of the gradient data of the entire work area, and its value adapts to the actual distribution of the work area data without relying on any manual input.

[0049] S43. Introduce a first empirical coefficient and a second empirical coefficient to scale the benchmark threshold proportionally to obtain an adaptive hierarchical threshold including a first-level threshold and a second-level threshold; wherein, the first empirical coefficient is 0.6 and the second empirical coefficient is 1.2.

[0050] Obtain the baseline threshold Subsequently, the present invention introduces a first empirical coefficient. Second empirical coefficient , respectively with Multiply the values ​​and scale them proportionally to obtain the first-level threshold. and secondary threshold It should be noted that, due to and It can be adjusted independently, and the ratio of the two threshold levels is... Fixed, with and The specific value is irrelevant, therefore the relative width of the gradation interval is not affected by the gradient distribution pattern.

[0051] S5. Based on the adaptive grading threshold, the grid cells of the entire work area are automatically partitioned and marked according to the gradient magnitude to obtain the partitioning and marking results.

[0052] In step S5, the partitioning and marking results specifically include: grid cells with gradient magnitude less than or equal to the first threshold are marked as smooth regions; grid cells with gradient magnitude greater than the first threshold and less than the second threshold are marked as transition regions; and grid cells with gradient magnitude greater than or equal to the second threshold are marked as steep transition regions.

[0053] The present invention receives the first-level threshold output by S43. and secondary threshold and the gradient magnitude of each grid cell in the grid gradient field. As input, all grid cells in the entire work area are individually partitioned and labeled. For gradient magnitudes that satisfy... The grid cells in this invention are marked as smooth regions, where the properties of the region where these cells are located change smoothly; for gradient magnitudes satisfying The grid cells in this invention are labeled as transition regions, and these cells are located in transition areas with moderate gradient changes; for gradient magnitudes that satisfy... The grid cells are labeled as abrupt change regions, which are located in areas with drastic gradient changes (such as sand body boundaries, pinch-out lines, and near faults). Ultimately, the three types of labeled regions do not overlap, covering all grid cells in the entire work area, and the output is a partitioned labeling result containing the region category labeling of each grid cell.

[0054] S6. Perform hierarchical densification and smoothing optimization on the grid cells corresponding to the transition zone and steep change zone in the partitioning marking results to obtain the adaptive densified corner grid model, so as to complete the adaptive partitioning of the thin reservoir planar grid.

[0055] Step S6 further includes: S61: Perform hierarchical encryption and smoothing optimization on the grid cells corresponding to the transition zone and abrupt change zone in the partitioning marking result to obtain a encrypted grid.

[0056] Specifically, step S61 includes: performing first-level encryption on the marked grid cells in the transition zone in the first grid subdivision direction and the second grid subdivision direction respectively, wherein the first-level encryption is subdivided by a 2x encryption factor; and performing second-level encryption on the marked grid cells in the abrupt change zone in the first grid subdivision direction and the second grid subdivision direction respectively, wherein the second-level encryption is subdivided by a 4x encryption factor.

[0057] In step S61, the present invention receives the partitioning marking result output in S5, and performs different multiples of encryption operations on the mesh cells marked as transition regions and abrupt change regions, while the mesh cells in the smooth region retain their original size. For the mesh cells in the transition region, the present invention subdivides them by a multiple of 2 in the first meshing direction (I direction) and the second meshing direction (J direction), that is, each original transition region mesh cell is evenly divided into 2 parts in the I direction and 2 parts in the J direction, and each original cell is subdivided into... Sub-cells; for grid cells in the abrupt change region, this invention subdivides them in both the I and J directions by a multiplier of 4, and each original cell is subdivided into... Subunits.

[0058] During the encryption process, this invention maintains the consistency of the mesh topology in the IJ direction, ensures the alignment of mesh nodes between adjacent encryption levels, and avoids the generation of dangling nodes, i.e., a mesh node is only connected to the midpoint of the edge of an adjacent cell instead of the corner point, resulting in abnormal nodes with discontinuous topology, or non-manifold structures, i.e., there are more than two faces in the mesh sharing the same edge, which are geometric topological errors.

[0059] S62. For the encrypted mesh, the Laplacian smoothing algorithm is used to iteratively optimize the node coordinates within the encrypted mesh area to obtain the corner mesh model.

[0060] After hierarchical densification, distortion problems such as uneven node distribution and excessively small element angles occur due to the influence of the densification operation on the grid node coordinates near the fault and in the densification boundary region. Therefore, this invention employs the Laplacian smoothing algorithm to iteratively optimize the node coordinates within the densified grid region. Specifically, the core operation of the Laplacian smoothing algorithm in this invention is to iteratively optimize the node coordinates within the densified grid region for each node... Obtain the coordinates of all directly connected neighboring nodes, calculate the arithmetic mean of these coordinates, and use this average as the node's coordinates. The new coordinates are used to complete one iteration update.

[0061] This invention performs the above coordinate update operation on all nodes in the region to complete one round of iteration. Then, the updated node coordinates are used as input to enter the next round of iteration. This process is repeated 3 to 5 times to make each node gradually move closer to the geometric center of its neighboring nodes, eliminate distorted cells and excessively small angles in the mesh, and finally output a corner mesh model with uniform node distribution and regular cell shape.

[0062] like Figure 2 The figure shown is a comparison of the mesh before and after adaptive densification based on the method of the present invention. It shows the comparison of the property surface distribution and mesh subdivision results before and after densification for two thin reservoir layers, NmIII-1 and NmIII-2.

[0063] Figure 2 The top two images on the left and right are two small-scale seismic attribute surface distribution maps. The colors from blue to red represent the changes in attribute values ​​from low to high. The areas with dense contour lines correspond to regions where the attribute values ​​change drastically in the lateral direction, such as sand body boundaries and pinch-out zones; the areas with sparse contour lines correspond to background areas where the attribute values ​​change gently. Figure 2 In the upper left figure (NmIII-1), the area marked by the red circle shows obvious blue-green low value anomalies, with drastic horizontal changes in attributes and extremely dense contour lines; Figure 2 The upper right image (NmIII-2) is predominantly warm in color, with relatively uniform property distribution, but some areas still exhibit abrupt band-like changes.

[0064] Figure 2 The two images at the bottom left and right are planar representations of the corresponding mesh partitioning results. Figure 2 The lower left image (NmIII-1) shows the corner mesh after adaptive encryption according to the present invention. The mesh is displayed in blue-green. The mesh is obviously finer in the area marked by the red circle, which corresponds to the low value anomaly area with steep attribute change in the upper left image. The mesh cell size in this area is significantly smaller than the surrounding gentle area, which intuitively reflects the hierarchical encryption result of 4 times encryption in the steep change area and 2 times encryption in the transition area. Figure 2 The lower right image (NmIII-2) shows the corresponding encrypted mesh, displayed in yellow. The spatial distribution of the encrypted mesh area closely matches the strip-shaped area with abrupt attribute changes in the upper right image, while the mesh size in the flat area remains unchanged from the original coarse mesh.

[0065] based on Figure 2 The comparison results of the two sub-layers show that the adaptive encrypted mesh output by this invention has a spatial location of the encrypted region that strictly corresponds to the spatial distribution of the gradient abrupt change region in the attribute surface. The encrypted range is automatically generated by the attribute gradient data without relying on manual intervention, and different sub-layers automatically generate different encrypted partition patterns due to different attribute distributions.

[0066] like Figure 3 As shown, the present invention also provides an adaptive meshing system for thin reservoir planar grids based on attribute gradients, comprising: Preprocessing module 100: Used to preprocess the seismic attribute surface data of the work area and construct a planar attribute field; Calculation module 200: used to calculate the gradient magnitude of the planar attribute field using the finite difference method to obtain the gradient field; Mapping module 300: used to map the gradient field to the center of the corner grid cell using a bilinear interpolation algorithm to obtain the grid gradient field; Generation module 400: used to calculate the statistical characteristic values ​​of the gradient magnitude of the entire work area in the grid gradient field, and to generate an adaptive hierarchical threshold by weighting the statistical characteristic values; Marking module 500: used to automatically partition and mark the grid cells of the entire work area according to the gradient magnitude based on the adaptive grading threshold, and obtain the partition marking result; Encryption module 600: used to perform hierarchical encryption and smoothing optimization on the grid cells corresponding to the transition zone and steep change zone in the partitioning marking result, to obtain the adaptively encrypted corner grid model, so as to complete the adaptive partitioning of thin reservoir planar grid.

[0067] The device embodiments described above are merely illustrative. The units described as separate components may or may not be physically separate. The components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the modules can be selected to achieve the purpose of this embodiment according to actual needs. Those skilled in the art can understand and implement this without any creative effort.

[0068] Through the above description of the embodiments, those skilled in the art can clearly understand that each embodiment can be implemented by means of software plus necessary general-purpose hardware platforms, and of course, it can also be implemented by hardware. Based on this understanding, the above technical solutions, in essence or the part that contributes to the prior art, can be embodied in the form of a software product. This computer software product can be stored in a computer-readable storage medium, such as ROM / RAM, magnetic disk, optical disk, etc., including several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute the attribute gradient-based thin reservoir planar grid adaptive partitioning method described in the various embodiments or some parts of the embodiments.

[0069] This invention achieves fully automated hierarchical densification of planar grids in thin reservoirs from a data-driven perspective by combining the seismic attribute surface gradient field with an adaptive subdivision mechanism for corner grids. In the attribute field construction stage, this invention performs coordinate unification, denoising, and normalization preprocessing on multi-source seismic attribute data, eliminating the impact of dimensional differences and noise interference in the original data on subsequent gradient calculations, ensuring the numerical consistency and spatial continuity of the planar attribute field. In the gradient calculation stage, this invention uses the finite difference method to calculate the partial derivatives of the attribute field in two horizontal directions, and obtains the gradient magnitude through square root extraction, ensuring clear numerical representation of abrupt attribute changes in areas such as sand body boundaries, pinch-out zones, and fault vicinity, providing accurate quantitative basis for subsequent zoning judgments. In the gradient mapping stage, this invention uses bilinear interpolation to accurately map the gradient magnitudes on the attribute surface grid to the center of the corner grid cells, achieving lossless transfer of gradient information between grid systems of different resolutions and ensuring the spatial consistency between the gradient field and the corner grid. In the threshold generation stage, this invention calculates a baseline threshold based on the mean and standard deviation of the gradient amplitude across the entire work area. A first empirical coefficient of 0.6 and a second empirical coefficient of 1.2 are introduced to proportionally scale the baseline threshold, generating primary and secondary thresholds. This allows the tiered thresholds to adaptively adjust according to the actual distribution of the work area data, overcoming the subjectivity and non-repeatability issues caused by manually setting thresholds in traditional methods. In the partitioning stage, this invention automatically divides the entire work area's grid cells into three categories: smooth zones, transition zones, and abrupt change zones based on the two-level thresholds. The partition boundaries are directly driven by gradient data, aligning with the geological change characteristics in the attribute surface, eliminating the need for manual delineation of the densification range. In the tiered densification stage, this invention performs double densification on the transition zone and quadruple densification on the abrupt change zone, while maintaining the original grid size in the smooth zone. This achieves refined characterization in critical areas with drastic attribute changes while effectively controlling the total number of grid cells across the entire work area, avoiding the exponential growth of grid size caused by globally uniform densification. During the smoothing optimization stage, this invention uses the Laplacian smoothing algorithm to iteratively update the node coordinates within the encrypted region 3 to 5 times, eliminating distorted cells and excessively small angles at the encrypted boundary, ensuring the geometric quality and topological compliance of the final corner mesh model, and enabling the output mesh to be directly connected to commercial numerical simulation software.

[0070] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. An adaptive mesh generation method for thin reservoir planar grids based on attribute gradients, characterized in that, include: S1. Preprocess the seismic attribute surface data of the work area and construct a planar attribute field; S2. The gradient field is obtained by calculating the gradient magnitude of the planar attribute field using the finite difference method. S3. The gradient field is mapped to the center of the corner grid cell using a bilinear interpolation algorithm to obtain the grid gradient field; S4. Calculate the statistical characteristic value of the gradient magnitude of the entire working area in the grid gradient field, and generate an adaptive hierarchical threshold by weighting the statistical characteristic value; S5. Based on the adaptive hierarchical threshold, the grid cells of the entire work area are automatically partitioned and labeled according to the gradient magnitude to obtain the partitioning and labeling results. S6. Perform hierarchical densification and smoothing optimization on the grid cells corresponding to the transition zone and steep change zone in the partitioning marking results to obtain the adaptive densified corner grid model, so as to complete the adaptive partitioning of the thin reservoir planar grid.

2. The adaptive meshing method for thin reservoir planar grids based on attribute gradients according to claim 1, characterized in that, Step S2 further includes: S21. The partial derivatives of the planar attribute field are calculated in the first horizontal direction and the second horizontal direction using the finite difference method to obtain the partial derivative values ​​in the corresponding directions; wherein the finite difference method includes forward finite difference, backward finite difference, and central finite difference method; S22. Perform a square superposition operation and a square root operation on the partial derivative values ​​in the first horizontal direction and the partial derivative values ​​in the second horizontal direction to obtain a gradient field based on the gradient magnitude.

3. The adaptive meshing method for thin reservoir planar grids based on attribute gradients according to claim 1, characterized in that, Step S3 further includes: S31. Traverse each grid cell in the corner grid one by one to obtain the center coordinates of the corresponding grid cell; S32. Based on the position of the center coordinates of the mesh cell in the attribute surface rectangular mesh, find and determine the attribute surface mesh cell that contains the center point of the corresponding mesh cell; S33. Based on the coordinates of the four corner points of the attribute surface grid cell, the gradient magnitude at the four corner points of the attribute surface grid cell is mapped to the cell center of the corresponding grid cell through a bilinear interpolation function. After assigning gradient values ​​to all grid cells in the corner grid, the gradient field is mapped from the attribute surface grid to the corner grid, and the grid gradient field is obtained.

4. The adaptive meshing method for thin reservoir planar grids based on attribute gradients according to claim 3, characterized in that, In step S33, the expression for the bilinear interpolation function is: in, The x-coordinate value of the center coordinate of the grid cell. The ordinate value of the center coordinate of the grid cell. Center of corner grid cell The gradient magnitude is calculated at [location]. The minimum coordinate value of the attribute surface mesh cell in the first horizontal direction. The maximum coordinate value of the attribute surface mesh cell in the first horizontal direction. This represents the minimum coordinate value of the attribute surface mesh cell in the second horizontal direction. This represents the maximum coordinate value of the attribute surface mesh cell in the second horizontal direction. This represents the gradient magnitude at the bottom left corner of the attribute surface mesh cell. This represents the gradient magnitude at the bottom right corner of the attribute surface mesh cell. This represents the gradient magnitude at the top-left corner of the attribute surface mesh cell. This represents the gradient magnitude at the top right corner of the attribute surface mesh cell. The area weight of the bottom left corner point of the attribute surface mesh cell. The area weight of the bottom right corner point of the attribute surface mesh cell. The area weight of the top-left corner point of the attribute surface mesh cell. This represents the area weight of the top-right corner point of the attribute surface mesh cell.

5. The adaptive meshing method for thin reservoir planar grids based on attribute gradients according to claim 1, characterized in that, Step S4 further includes: S41. Perform statistical analysis on the gradient magnitude of all grid cells in the grid gradient field, and calculate the mean and standard deviation of the gradient magnitude of the entire work area. S42. Calculate a baseline threshold based on gradient statistical features, wherein the baseline threshold is used to filter out background fluctuation areas; S43. Introduce a first empirical coefficient and a second empirical coefficient, and scale the benchmark threshold proportionally to obtain an adaptive hierarchical threshold including a first-level threshold and a second-level threshold.

6. The adaptive meshing method for thin reservoir planar grids based on attribute gradients according to claim 5, characterized in that, In step S43, the first empirical coefficient is set to 0.6, and the second empirical coefficient is set to 1.

2.

7. The adaptive meshing method for thin reservoir planar grids based on attribute gradients according to claim 5, characterized in that, Step S5, the partition marking result specifically includes: For grid cells whose gradient magnitude is less than or equal to the first threshold, they are marked as smooth regions; Mesh cells with gradient magnitudes greater than the first threshold and less than the second threshold are marked as transition regions; For grid cells with gradient magnitudes greater than or equal to the second threshold, they are marked as steep transition regions.

8. The adaptive meshing method for thin reservoir planar grids based on attribute gradients according to claim 1, characterized in that, Step S6 further includes: S61: Perform hierarchical densification and smoothing optimization on the grid cells corresponding to the transition and steep change regions in the partitioning marking results to obtain a densified grid; S62. For the encrypted mesh, the Laplacian smoothing algorithm is used to iteratively optimize the node coordinates within the encrypted mesh area to obtain the corner mesh model.

9. The adaptive meshing method for thin reservoir planar grids based on attribute gradients according to claim 8, characterized in that, Step S61 specifically includes: For the marked grid cells in the transition zone, a first-level densification is performed in the first grid subdivision direction and the second grid subdivision direction respectively, and the first-level densification is subdivided by a densification factor of 2. For the marked grid cells in the steep transition zone, secondary densification is performed in the first grid subdivision direction and the second grid subdivision direction respectively, and the secondary densification is subdivided by a densification factor of 4.

10. A thin reservoir planar mesh adaptive partitioning system based on attribute gradient, used to execute the thin reservoir planar mesh adaptive partitioning method based on attribute gradient as described in any one of claims 1 to 9, characterized in that, include: Preprocessing module: Used to preprocess the seismic attribute surface data of the work area and construct a planar attribute field; Calculation module: used to calculate the gradient magnitude of the planar attribute field using the finite difference method to obtain the gradient field; Mapping module: Used to map the gradient field to the center of the corner grid cell using a bilinear interpolation algorithm to obtain the grid gradient field; Generation module: used to calculate the statistical characteristic values ​​of the gradient magnitude of the entire work area in the grid gradient field, and to generate an adaptive hierarchical threshold by weighting the statistical characteristic values; The marking module is used to automatically partition and mark the grid cells of the entire work area according to the gradient magnitude based on the adaptive grading threshold, and obtain the partition marking results. Encryption module: used to perform hierarchical encryption and smoothing optimization on the grid cells corresponding to the transition zone and steep change zone in the partitioning marking results, to obtain the adaptively encrypted corner grid model, so as to complete the adaptive partitioning of thin reservoir planar grid.

Images