Method for constructing three-dimensional gridding stratigraphic interface with geologic contact geometric consistency

Abstract

The invention discloses a geological contact geometry consistent three-dimensional gridding stratigraphic interface construction method, which comprises the following steps of: constructing a three-dimensional base grid, and performing three-dimensional gridding on each fault plane to obtain a fault plane polygon; forming a geological relation graph, and performing node marking on the fault plane polygon; constructing a split grid unit binary tree, and performing virtual and real marking; constructing a grid unit edge splitting binary tree, and adding and numbering virtual grid points; establishing a grid point neighborhood relation table, and solving to obtain the value of each numbered grid point; and extracting stratigraphic interface polygons, and cutting the stratigraphic interface polygons and the fault plane polygons to obtain stratigraphic interfaces with geologic contact geometry consistent. According to the method, by considering the fault cutting relation of multiple fault plane polygons, the geological contact geometric consistency is ensured, and the accuracy of formation interface construction is improved.

Description

Technical Field

[0001] The present application relates to the technical field of three-dimensional geological modeling, and in particular to a method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry. Background Art

[0002] The method of constructing the stratigraphic interface is an important part of geological modeling. The geometric consistency of geological contact is the relationship between the geometric form and position of the stratigraphic interface at the contact point, especially whether they can be in close contact, without gaps or intersections. In geology, the geometric consistency of geological contact is of great significance for understanding geological processes such as crustal movement, sedimentation, and tectonic deformation. Therefore, the geometric consistency of geological contact needs to be considered when studying the method of constructing the stratigraphic interface.

[0003] In the prior art, Chinese patent CN103500468A discloses a spatial surface reconstruction method based on boundary constraints, which collects original stratigraphic data points and fault data points, draws an upper plate polygon and a lower plate polygon; obtains the height attributes, relationship attributes and projection attributes of the original stratigraphic data points and the points on the fault polygon; connects the original stratigraphic data points and the points on the fault polygon using constrained triangular meshing; completes the splicing of the inner and outer areas of the upper plate polygon and the lower plate polygon, and thus completes the reconstruction of the three-dimensional spatial surface.

[0004] However, the above-mentioned prior art only considers the tangency between the upper plate and the lower plate of a fault, but does not consider the fault cutting relationship of multiple fault plane polygons, resulting in poor geometric consistency of geological contact and poor accuracy of stratigraphic interface construction. Summary of the invention

[0005] The present application provides a three-dimensional gridded stratum interface construction method with consistent geological contact geometry, which is used to solve the problems in the prior art that the fault cutting relationship of multiple fault plane polygons is not considered, the geological contact geometry consistency is poor, and the accuracy of stratum interface construction is poor.

[0006] On the one hand, the present application provides a method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry, comprising the following steps: Step 1: Acquire geological modeling source data, construct a three-dimensional base grid according to the geological modeling source data, and use the three-dimensional base grid to three-dimensionally grid each fault plane to obtain a fault plane polygon.

[0007] Step 2: forming a geological relationship map according to the fault cutting relationship and the fault plane normal relationship of the fault plane polygon, and marking nodes of the fault plane polygon according to the geological relationship map.

[0008] Step three: construct a split grid unit binary tree according to the geological relationship map, and mark the grid points of each grid unit in the split grid unit binary tree as real or virtual.

[0009] Step 4: Check the validity of the split grid unit according to the node mark and the virtual and real mark. If there is no invalid situation, go to step 5.

[0010] Step 5: construct a grid unit edge splitting binary tree, and add and number virtual grid points to the split grid units according to the edge splitting results of the grid unit edge splitting binary tree.

[0011] Step six, establishing different types of grid point neighborhood relationship tables, and solving the value of each numbered grid point according to the grid point number and the grid point neighborhood relationship table.

[0012] Step seven, extracting the stratigraphic interface polygon according to the values ​​of the numbered grid points, cutting the stratigraphic interface polygon and the fault plane polygon to obtain a stratigraphic interface with consistent geological contact geometry.

[0013] In a possible implementation, in step 1, the three-dimensional gridding of each fault plane using the three-dimensional base grid to obtain the fault plane polygon includes: Each fault plane is intersected with the three-dimensional base grid, and the grid unit edge of the three-dimensional base grid forms a grid unit edge intersection with each fault plane. It is stipulated that each grid unit edge has only one intersection point, and the intersection points of each grid unit edge with more than one intersection point are merged by algorithm weighted average to form a closed polygon.

[0014] The broken line edges of the closed polygon on the grid unit surface are transformed into straight lines to obtain the fault plane polygon.

[0015] In a possible implementation, in step 2, forming a geological relationship map according to the fault cutting relationship and the fault plane normal relationship of the fault plane polygon includes: The positive and negative directions of the fault plane are set according to the fault cutting relationship and the normal direction relationship of the fault plane polygon to form a geological relationship map.

[0016] In a possible implementation, in step 2, marking nodes of the fault plane polygon according to the geological relationship map includes: The nodes of each fault plane polygon are traversed according to the geological relationship map, and the nodes of other fault plane polygons that are cut by the current fault plane polygon and are not located in the cutting direction of the current fault plane polygon are marked as invalid nodes.

[0017] In a possible implementation, step three includes: The nodes of each fault plane polygon are traversed according to the geological relationship diagram, and a split grid unit binary tree is constructed according to the cutting relationship between the positive part and the negative part of each fault plane polygon.

[0018] The split grid unit binary tree is traversed, and the grid points of each grid unit are marked as virtual grid points and real grid points according to the direction in which the grid points of the corresponding fault plane polygon are located.

[0019] In a possible implementation, step five includes: The split grid unit binary tree is traversed, and a grid unit edge split binary tree is constructed according to the intersection relationship between each fault plane polygon and the grid unit edge.

[0020] According to the edge splitting result of the grid unit edge splitting binary tree, virtual grid points are added to the split grid unit, wherein the virtual grid points include edge splitting grid points, fault termination splitting grid points and completion effect grid points, and the added virtual grid points are numbered.

[0021] In a possible implementation, in step 6, establishing different types of grid point neighborhood relationship tables includes: A six-neighborhood relationship table between the original grid points and the added virtual grid points is established.

[0022] In a possible implementation, in step 6, solving the value of each numbered grid point according to the grid point number and the grid point neighborhood relationship table includes: According to the grid point number and the grid point neighborhood relationship table, a finite difference format for solving the implicit geological field interpolation function and a linear function relationship between the control point and the grid node are constructed, a solution equation is established, and a multi-grid method is used for iterative solution to obtain the value of each numbered grid point.

[0023] A three-dimensional gridded stratum interface construction method with consistent geological contact geometry in this application has the following advantages: By combining the three-dimensional base grid, geological relationship map, split grid unit binary tree, grid unit edge split binary tree, grid point neighborhood relationship table and polygon cutting processing, the fault cutting relationship of multiple fault plane polygons is considered to ensure the geometric consistency of geological contact and improve the accuracy of stratigraphic interface construction. BRIEF DESCRIPTION OF THE DRAWINGS

[0024] In order to more clearly illustrate the embodiments of the present application or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present application. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0025] Figure 1 A schematic flow chart of a method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry provided in an embodiment of the present application; Figure 2 A schematic diagram of constructing a three-dimensional base grid based on geological modeling source data provided in an embodiment of the present application; Figure 3 A schematic diagram of a grid unit provided in an embodiment of the present application; Figure 4 A schematic diagram of a fragmentary graphic of a fault plane in each grid unit provided in an embodiment of the present application; Figure 5 A schematic diagram of the intersection of a fault plane and a three-dimensional base grid provided in an embodiment of the present application; Figure 6 A schematic diagram of a process for forming a closed polygon provided in an embodiment of the present application; Figure 7 A schematic diagram of a process of obtaining a fault plane polygon by converting a straight line provided in an embodiment of the present application; Figure 8 A schematic diagram of cutting three faults F1, F2 and F3 provided in an embodiment of the present application; Fig. 9 A geological relationship diagram of the three faults F1, F2 and F3 provided in the embodiment of the present application; Fig.10 A schematic diagram of polygonal node markings on a fault plane provided in an embodiment of the present application; Fig.11 A schematic diagram of a split grid unit binary tree provided in an embodiment of the present application; Fig.12 A grid unit edge splitting binary tree provided in an embodiment of the present application; Fig.13 A schematic diagram of extracting stratum interface polygons provided in an embodiment of the present application; Fig.14 A geological relationship diagram of the faults F1, F2 and the formation interface H provided in the embodiment of the present application; Fig.15 A schematic diagram of a stratum interface with consistent geological contact geometry provided in an embodiment of the present application; Fig.16 A schematic diagram of fault plane data of 5 fault planes provided in an embodiment of the present application; Fig.17 The embodiment of the present application provides Fig.16 Schematic diagram of the five fault planes generated; Fig.18 The embodiment of the present application provides Fig.17 Schematic diagram of the cutting of the five generated fault planes; Fig.19 The embodiment of the present application provides Fig.18 A schematic diagram of the stratigraphic interface with consistent geological contact geometry generated through steps 2 to 7; Fig. 20 The embodiment of the present application provides Fig.19 Schematic diagram of the stratigraphic interfaces without showing the fault planes. DETAILED DESCRIPTION

[0026] The following will be combined with the drawings in the embodiments of the present application to clearly and completely describe the technical solutions in the embodiments of the present application. Obviously, the described embodiments are only part of the embodiments of the present application, not all of the embodiments. Based on the embodiments in the present application, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of this application.

[0027] like Figure 1 As shown, the embodiment of the present application provides a method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry, comprising the following steps: Step 1: Acquire geological modeling source data, construct a three-dimensional base grid according to the geological modeling source data, and use the three-dimensional base grid to three-dimensionally grid each fault plane to obtain a fault plane polygon.

[0028] Step 2: forming a geological relationship map according to the fault cutting relationship and the fault plane normal relationship of the fault plane polygon, and marking nodes of the fault plane polygon according to the geological relationship map.

[0029] Step three: construct a split grid unit binary tree according to the geological relationship map, and mark the grid points of each grid unit in the split grid unit binary tree as real or virtual.

[0030] Step 4: Check the validity of the split grid unit according to the node mark and the virtual and real mark. If there is no invalid situation, go to step 5.

[0031] Step 5: construct a grid unit edge splitting binary tree, and add and number virtual grid points to the split grid units according to the edge splitting results of the grid unit edge splitting binary tree.

[0032] Step six, establishing different types of grid point neighborhood relationship tables, and solving to obtain the value of each numbered grid point according to the grid point number and the grid point neighborhood relationship table.

[0033] Step seven, extracting the stratigraphic interface polygon according to the values ​​of the numbered grid points, cutting the stratigraphic interface polygon and the fault plane polygon to obtain a stratigraphic interface with consistent geological contact geometry.

[0034] like Figure 2 As shown, specifically, in this embodiment, constructing a three-dimensional base grid according to the geological modeling source data includes: calculating the range of the minimum hexahedron according to the geological modeling source data, and then expanding the minimum hexahedron along the XYZ direction to 1 / 10 of the length of each side of the hexahedron to form a modeling range. The three-dimensional modeling area is divided into a three-dimensional regular grid of M*N*T along the XYZ direction according to the modeling range, that is, a three-dimensional base grid.

[0035] In this embodiment, the geological modeling source data srcData is described as a binary tuple srcData=(F,H), where F represents the set of all fault surface data, H represents the set of all horizon data, and each horizon includes the horizon control point of the horizon. The range of the minimum hexahedron envelope is calculated based on all fault surface data and horizon data in the geological modeling source data srcData, and then the minimum hexahedron envelope is expanded along the XYZ direction to 1 / 10 of the length of each side of the hexahedron to form a modeling range. The three-dimensional modeling area is divided into a three-dimensional regular grid of M*N*T (in this embodiment, M=N=T=5) along the XYZ direction according to the modeling range, that is, a three-dimensional base grid.

[0036] like Figure 3 As shown in the figure, each hexahedron that makes up the three-dimensional base grid is called a grid unit, which has 12 edges and 8 grid points. Except for the boundary edges, each edge is shared by 4 grid units. The 12 edges of the grid unit are numbered in counterclockwise and top-to-bottom order. The default order is as follows: Figure 3 As shown, it is implicitly stored. According to this order, the neighborhood relationship of the edge can be queried by constructing a lookup table. Each edge has an edge id number, each grid cell has a cell id number, and each grid cell records the edge id numbers of 12 edges. Each edge is composed of two grid points. The grid points have grid point id numbers, and the edges reference the grid point id numbers of the grid points. Each edge also requires an original edge id number, and the original edge id number can find the edge splitting binary tree. During initialization, each grid cell of the three-dimensional base grid is assigned an original cell id number, each edge is assigned an original edge id number, and each grid point is assigned an original grid point id number.

[0037] Exemplarily, in step 1, the three-dimensional gridding of each fault plane using the three-dimensional base grid to obtain the fault plane polygons includes: Each fault plane is intersected with the three-dimensional base grid, and the grid unit edge of the three-dimensional base grid forms a grid unit edge intersection with each fault plane. It is stipulated that each grid unit edge has only one intersection point, and the intersection points of each grid unit edge with more than one intersection point are merged by algorithm weighted average to form a closed polygon.

[0038] The broken line edges of the closed polygon on the grid unit surface are transformed into straight lines to obtain the fault plane polygon.

[0039] like Figure 4 As shown, specifically, within each grid unit, there are four possibilities for the fragment graphics of the fault plane: triangle, quadrilateral, pentagon and hexagon. Each fragment graphic is formed by a corresponding node, and the node state is marked as a valid state. A hash lookup table is used for the edge where the node is located to establish node information. The node information includes the fault plane ID number, node state, normal information and whether it is a boundary node. This embodiment stipulates that each grid unit edge has only one intersection, and the line of each grid unit face is a straight line segment. The purpose of three-dimensional gridding of the fault plane is that each grid unit is split into two independent areas by a fault plane.

[0040] like Figure 5 As shown, each fault plane is intersected with the three-dimensional base grid. The grid unit edge of the three-dimensional base grid forms a grid unit edge intersection with each fault plane. In this embodiment, each grid unit edge has only one intersection point, such as Figure 6 As shown in , the intersection points of each grid unit with more than one edge are merged by algorithm weighted average to form a closed polygon. Figure 7 As shown, the polyline edge of the closed polygon on the grid unit surface is transformed into a straight line to obtain a fault plane polygon. If a node of a fault plane polygon has no adjacent fault plane polygon, the node is set as a boundary node.

[0041] Exemplarily, in step 2, forming a geological relationship map according to the fault cutting relationship and the fault plane normal relationship of the fault plane polygon includes: The positive and negative directions of the fault plane are set according to the fault cutting relationship and the normal direction relationship of the fault plane polygon to form a geological relationship map.

[0042] Specifically, according to the semantic characteristics of geology, the positive direction of the stratum or rock interface points to the newer direction of the geological era, and the negative direction points to the older direction of the geological era. The positive direction of the fault plane points to the direction of the upper wall of the fault, and the negative direction points to the direction of the lower wall of the fault. Any two intersecting fault planes are cut in the grid unit according to the relationship between the old and new cutting of the fault and the positive and negative normal directions.

[0043] like Figure 8 As shown, there are three faults, F1, F2 and F3, where F1 cuts F2 and F3, F2 cuts F3, F2 is in the positive direction of F1, and F3 is in the negative direction of F2. Fig. 9 The geological relationship diagram shown. Fig. 9In the geological relationship diagram, the cutting is transitive, that is, F1 cuts F2, F2 cuts F3, and F1 cuts F3. The normal relationship is inclusive, F2 is in the positive direction of F1 and F2 cuts F3, then F3 is in the positive direction of F1. Figure 8 The three faults in the figure divide the hexahedral grid cell into four sub-regions, each of which is a polyhedron.

[0044] Exemplarily, in step 2, marking nodes of the fault plane polygon according to the geological relationship map includes: The nodes of each fault plane polygon are traversed according to the geological relationship map, and the nodes of other fault plane polygons that are cut by the current fault plane polygon and are not located in the cutting direction of the current fault plane polygon are marked as invalid nodes.

[0045] Specifically, each fault plane is a complete polygon within the grid unit, the nodes of the polygon are located on the edges of the grid unit, and its initial state is valid. Fig. 9 The geological relationship diagram is constructed by traversing each node Fi of the geological relationship diagram one by one, and the section nodes that are cut by Fi and are not located in the cutting direction of Fi are marked as invalid. Fig. 9 and Fig.10 As shown, traverse F1, F1 cuts F2 and F3, F2 and F3 are located in the positive direction of F1, then the nodes of F2 and F3 located in the negative direction of F1 are set as invalid nodes; traverse F2, F2 cuts F3, F3 is located in the negative direction of F2, then the node of F3 located in the positive direction of F2 is set as invalid node; traverse F3, F3 does not cut other faults, and the traversal terminates. Fig.10 In the figure, black nodes are valid nodes and white nodes are invalid nodes.

[0046] Exemplarily, step three includes: The nodes of each fault plane polygon are traversed according to the geological relationship diagram, and a split grid unit binary tree is constructed according to the cutting relationship between the positive part and the negative part of each fault plane polygon.

[0047] The split grid unit binary tree is traversed, and the grid points of each grid unit are marked as virtual grid points and real grid points according to the direction in which the grid points of the corresponding fault plane polygon are located.

[0048] Specifically, each fault plane polygon divides the grid unit into two parts, the part located in the positive direction of the polygon is called the positive part, and the part located in the negative direction of the polygon is called the negative part. By traversing each fault node in the geological relationship map in turn, and cutting according to the positive and negative parts of the fault plane polygon, a binary tree of split grid units can be constructed. The root node is the original grid unit. Except for the leaf nodes, each binary tree node has two subtree nodes, namely the negative part tree node and the positive part tree node. The grid unit corresponding to the negative part tree node is located in the negative direction of the currently traversed fault, and the grid unit corresponding to the positive part tree node is located in the positive direction of the currently traversed fault. The leaf node represents the real part of the current grid unit that is split, and the number of leaf nodes is consistent with the number of polyhedrons after the grid unit is split. For example Fig.11 As shown, Fig.11 The four leaf nodes in Figure 8 The three faults in the figure divide the hexahedral grid cell into four sub-regions, each of which is a polyhedron.

[0049] Specifically, traverse the split grid unit binary tree, first mark the grid point of the cell corresponding to the root node as a real grid point, that is, the Boolean value v is 1; then traverse the left and right subtrees respectively, and the default value of the grid point of the cell corresponding to the subtree node inherits the Boolean value v of the parent node. Then, according to the currently traversed fault Fi, determine whether the grid point is in the positive direction or the negative direction of Fi. The specific method is: for the positive part of the tree node, the grid point located in the positive direction of Fi is marked as a real grid point, and the Boolean value is 1∩v, otherwise it is marked as a virtual grid point, and the Boolean value is 0; for the negative part of the tree node, the grid point located in the positive direction of Fi is marked as a virtual grid point, and the Boolean value is 0, otherwise it is marked as a real grid point, and the Boolean value is 1∩v.

[0050] Specifically, in step 4, for the grid cell corresponding to the leaf node in the binary tree of split grid cells, if the grid points of the grid cell are all virtual grid points (i.e., the Boolean value is 0), then the grid cell is a small grid cell. If the two endpoints of a straight line segment of a fault on a certain surface of the grid cell are invalid nodes, then the grid cell divided this time is an invalid grid cell, that is, the grid cell is not associated with all the surrounding grid cells, and it is necessary to increase the grid resolution or modify the fault surface morphology and geological relationship, otherwise the modeling process cannot continue. If there is no invalid situation, go to step 5.

[0051] Exemplarily, step five includes: The split grid unit binary tree is traversed, and a grid unit edge split binary tree is constructed according to the intersection relationship between each fault plane polygon and the grid unit edge.

[0052] According to the edge splitting result of the grid unit edge splitting binary tree, virtual grid points are added to the split grid unit, wherein the virtual grid points include edge splitting grid points, fault termination splitting grid points and completion effect grid points, and the added virtual grid points are numbered.

[0053] Specifically, edge splitting is to break the edge according to the old and new relationship of the fault, generate new edges and new grid points to form split grid units, and the neighborhood relationship between edges and grid points needs to be considered. The root of the grid unit edge splitting binary tree is the current edge. Except for the leaf node, each tree node has a left and right subtree, corresponding to the left and right parts of the edge split. The definition of left and right is based on the positive direction of the rectangular coordinate, with the left smaller and the right larger. The leaf node represents the edge splitting result.

[0054] like Fig.12 As shown, the specific process is: traverse the split grid unit binary tree, if the fault plane polygon being traversed has an intersection with the edge (that is, the intersection of the edge and the fault), the treatment of the intersection as a valid intersection and an invalid intersection is the same, the edge is broken to form two new edges, and the endpoints of the two new edges at the intersection are added with virtual grid points (that is, the edge splitting grid points), the left endpoint of the left subtree inherits the left endpoint of the parent node, and the right endpoint is the newly added virtual grid point, and the right endpoint of the right subtree inherits the right endpoint of the parent node, and the left endpoint is the newly added virtual grid point.

[0055] Specifically, virtual grid points are added to the split grid units according to the edge splitting result of the binary tree of grid unit edge splitting, and the virtual grid points include edge splitting grid points, fault termination splitting grid points and completion grid points. The specific process is: according to the edge splitting result of the binary tree of grid unit edge splitting, the original grid points are inherited from the three-dimensional base grid. If the grid unit is split, the intersection of the edge and the fault at the split is processed according to whether the fault is terminated. The fault termination splitting grid point is added at the fault termination, otherwise the edge splitting grid point is added. If the grid points of the formed hexahedron are processed, the completion grid points are added. The added virtual grid points are numbered, and the adjacent grids at the fault termination splitting grid points have the same number. Adjacent grid units share edges, and their edge splitting grid points have the same number.

[0056] Exemplarily, in step six, establishing different types of grid point neighborhood relationship tables includes: A six-neighborhood relationship table between the original grid points and the added virtual grid points is established.

[0057] Specifically, in this embodiment, the grid point neighborhood relationship table is a six-neighborhood relationship table between the original grid point and the added virtual grid point, and each record in the table includes: grid point number, upper grid point number, lower grid point number, front grid point number, rear grid point number, left grid point number, and right grid point number. The specific construction process is: in the current split grid, for any vertex of the eight vertices, initialize its six-neighborhood relationship to -1; then find its six-neighborhood vertex number and set it in the corresponding six-neighborhood record.

[0058] In order to find the neighborhood relationship of the grid cells, a six-neighborhood relationship table of all original grid cells and split grid cells is also established. Each record in the table includes: grid cell number, upper grid cell number, lower grid cell number, front grid cell number, back grid cell number, left grid cell number, and right grid cell number.

[0059] Exemplarily, in step 6, solving the value of each numbered grid point according to the grid point number and the grid point neighborhood relationship table includes: According to the grid point number and the grid point neighborhood relationship table, a finite difference format for solving the implicit geological field interpolation function and a linear function relationship between the control point and the grid node are constructed, a solution equation is established, and a multi-grid method is used for iterative solution to obtain the value of each numbered grid point.

[0060] Specifically, solving the implicit geological field interpolation function is an optimization problem, as shown in the following formula: J

[0061] st:

[0062] Among them, J represents the energy functional, represents the roughness factor, For spatial location , , The implicit function of Indicates a known point The function value of .

[0063] The construction process of the finite difference format considers the differences in the three directions of x, y, and z. In any direction, assuming that the function value at a certain point is , then the function value at this point will be affected by the six grid points, and the second-order partial derivatives and mixed partial derivatives can be expressed by the following equations:

[0064]

[0065]

[0066]

[0067]

[0068]

[0069] in, , , They represent the step size of the grid in each direction, and i, j, and k represent the grid coordinates.

[0070] The difference matrix of the second-order partial derivative is constructed as shown below:

[0071]

[0072]

[0073] in, represents the second-order partial differential, N represents the number of grids, , , represents the difference matrix of different grid directions, u is the matrix of grid point function values, and b represents the constant matrix.

[0074] Similarly, the construction method of the difference matrix of the mixed partial derivative is the same as that of the second-order partial derivative, which is represented by H, I, and J. Then the finite difference format for solving the implicit geological field interpolation function is summarized as:

[0075] The linear function relationship between the control points and the grid nodes uses a three-dimensional linear interpolation formula as shown below:

[0076]

[0077]

[0078]

[0079] in, represents the function value at the corner points of the cube, and , , are the normalized coordinates of P relative to the cube boundary.

[0080] The Lagrange multiplier method is used to find local extreme values ​​in the optimization problem with constraints. The specific method is to convert the constrained problem into an unconstrained problem, so that the optimization problem can be solved by finding the critical point (the point where the gradient is zero) of the unconstrained problem. The above formula after differential discretization is combined into The linear equation system format is used to solve the corresponding attribute value of each numbered grid point by an iterative solution method. In this embodiment, an algebraic multigrid method is used for solving.

[0081] Specifically, step seven includes: according to the values ​​of the numbered grid points, all original grid cells and split grid cells are subjected to a moving cube algorithm to extract isosurfaces to obtain a stratum interface polygon, and according to the geological relationship between the fault and the stratum interface, the stratum interface polygon and the fault plane polygon are cut to obtain a stratum interface with consistent geological contact geometry. Among them, the moving cube algorithm is a commonly used isosurface extraction algorithm, which will not be described in detail here.

[0082] like Fig.13 As shown, the polygon of the stratum interface H is extracted.

[0083] like Fig.14 As shown in Figure 1, it is the geological relationship between faults F1, F2 and the formation interface H. Fig.14 The geological relationship is cut into Fig.15 As shown, the stratigraphic interface H after cutting is a stratigraphic interface with consistent geological contact geometry.

[0084] In a possible embodiment, actual geological modeling source data is collected and experimented with the method of the present application to verify the reliability of the method of the present application. Figures 16 to 20 As shown, Fig.16 The fault plane data of the five fault planes collected are: Fig.17 These are the five fault planes generated based on the fault plane data of the five fault planes. Fig.18 To cut the five fault planes, Fig.19 To perform steps 2 to 7 of the method of the present application on the five fault planes after cutting, to generate stratum interfaces with consistent geological contact geometry, Fig. 20 is to Fig.19 The five fault planes in the Fig. 20 It can be seen intuitively that the stratum interface is disconnected at the 5 fault planes according to the morphology of the fault planes, and the disconnection points are perfectly matched with the 5 fault planes, that is, the geological contact is consistent. Among them, F1`, F2`, F3`, F4`, and F5` represent 5 fault planes respectively, and H` represents the stratum interface with consistent geological contact geometry in this embodiment.

[0085] The embodiment of the present application combines a three-dimensional base grid, a geological relationship map, a split grid unit binary tree, a grid unit edge split binary tree, a grid point neighborhood relationship table and polygon cutting processing, considers the fault cutting relationship of multiple fault plane polygons, ensures the geometric consistency of geological contact, and improves the accuracy of stratigraphic interface construction.

[0086] Although the preferred embodiments of the present application have been described, those skilled in the art may make other changes and modifications to these embodiments once they have learned the basic creative concept. Therefore, the appended claims are intended to be interpreted as including the preferred embodiments and all changes and modifications falling within the scope of the present application.

[0087] Obviously, those skilled in the art can make various changes and modifications to the present application without departing from the spirit and scope of the present application. Thus, if these modifications and variations of the present application fall within the scope of the claims of the present application and their equivalents, the present application is also intended to include these modifications and variations.

Claims

1. A method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry, characterized in that: The following steps are involved: Step 1, obtaining geological modeling source data, constructing a three-dimensional base grid according to the geological modeling source data, and using the three-dimensional base grid to three-dimensionally grid each fault plane to obtain a fault plane polygon; Step 2: forming a geological relationship map according to the fault cutting relationship and the fault plane normal relationship of the fault plane polygon, and marking nodes of the fault plane polygon according to the geological relationship map; Step 3, constructing a split grid unit binary tree according to the geological relationship map, and marking the grid points of each grid unit in the split grid unit binary tree as virtual or real; Step 4: Check the validity of the split grid unit according to the node mark and the virtual and real mark. If there is no invalid situation, proceed to step 5; Step 5: construct a grid unit edge splitting binary tree, and add and number virtual grid points to the split grid units according to the edge splitting result of the grid unit edge splitting binary tree; Step 6, establishing different types of grid point neighborhood relationship tables, and solving the value of each numbered grid point according to the grid point number and the grid point neighborhood relationship table; Step seven, extracting the stratigraphic interface polygon according to the values ​​of the numbered grid points, cutting the stratigraphic interface polygon and the fault plane polygon to obtain a stratigraphic interface with consistent geological contact geometry.

2. A method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry according to claim 1, characterized in that: In step 1, the three-dimensional gridding of each fault plane using the three-dimensional base grid to obtain the fault plane polygons includes: Intersecting each fault plane with the three-dimensional base grid, the grid unit edge of the three-dimensional base grid and each fault plane form a grid unit edge intersection, each grid unit edge is required to have only one intersection, and the intersections of each grid unit edge with more than one intersection are combined by algorithm weighted average to form a closed polygon; The broken line edges of the closed polygon on the grid unit surface are transformed into straight lines to obtain the fault plane polygon.

3. The method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry according to claim 1, characterized in that: In step 2, forming a geological relationship map according to the fault cutting relationship and the fault plane normal relationship of the fault plane polygon includes: The positive and negative directions of the fault plane are set according to the fault cutting relationship and the normal direction relationship of the fault plane polygon to form a geological relationship map.

4. A method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry according to claim 1, characterized in that: In step 2, the node marking of the fault plane polygon according to the geological relationship map includes: The nodes of each fault plane polygon are traversed according to the geological relationship map, and the nodes of other fault plane polygons that are cut by the current fault plane polygon and are not located in the cutting direction of the current fault plane polygon are marked as invalid nodes.

5. The method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry according to claim 1, characterized in that: Step three includes: Traversing the nodes of each fault plane polygon according to the geological relationship diagram, and constructing a split grid unit binary tree according to the cutting relationship between the positive part and the negative part of each fault plane polygon; The split grid unit binary tree is traversed, and the grid points of each grid unit are marked as virtual grid points and real grid points according to the direction in which the grid points of the corresponding fault plane polygon are located.

6. The method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry according to claim 1, characterized in that: Step five includes: Traversing the split grid unit binary tree, and constructing a grid unit edge split binary tree according to the intersection relationship between each fault plane polygon and the grid unit edge; According to the edge splitting result of the grid unit edge splitting binary tree, virtual grid points are added to the split grid unit, wherein the virtual grid points include edge splitting grid points, fault termination splitting grid points and completion effect grid points, and the added virtual grid points are numbered.

7. The method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry according to claim 1, characterized in that: In step six, the establishment of different types of grid point neighborhood relationship tables includes: A six-neighborhood relationship table between the original grid points and the added virtual grid points is established.

8. The method for constructing a three-dimensional gridded stratum interface with consistent geological contact geometry according to claim 1, characterized in that: In step 6, the value of each numbered grid point is obtained by solving the grid point number and the grid point neighborhood relationship table, including: According to the grid point number and the grid point neighborhood relationship table, a finite difference format for solving the implicit geological field interpolation function and a linear function relationship between the control point and the grid node are constructed, a solution equation is established, and a multi-grid method is used for iterative solution to obtain the value of each numbered grid point.

Images