Open space well complex geology refined modeling method

Abstract

The invention relates to the technical field of geophysical exploration, in particular to an open space well complex geology refined modeling method. The invention provides a brand-new open space well complex geology refined modeling method in a tunnel exploration stage. A point set can be effectively optimized by interpolation according to the perimeter proportion based on the W configuration, the problem of sparsity and winding distribution discretization is solved, meanwhile, the winding difference of adjacent profile contour feature points is remarkably weakened, and a foundation is laid for enhancing the constraint of the feature point set on the profile shape and reconstructing the similarity between triangles in a triangulation network. The three-dimensional physical property model and the geological anomaly model of the research area can be efficiently and accurately established, the crushing and water-rich anomaly interpretation of the area is realized, and data reference is provided for unfavorable geology advanced prediction and tunnel safety construction.

Description

Technical Field

[0001] The invention relates to the field of geophysical detection technology, in particular to a method for fine modeling of complex geology of open-ground wells. Background Art

[0002] In recent years, railway construction in high-altitude areas of my country has developed rapidly. However, many lines are located in difficult and complex mountainous areas, and tunnel excavation is often required in areas with complex structures. The establishment of a three-dimensional model in the early stage of tunnel construction is particularly important. The interpretation of the broken and water-rich anomalies in this area provides reference data for the advanced prediction of adverse geology and safe tunnel construction.

[0003] Three-dimensional reconstruction and numerical simulation based on surface rendering are important research directions in the fields of geological information science, geophysical exploration, mining engineering, etc. In the context of geological big data, the model construction and visualization of scattered three-dimensional data fields such as aviation, ground, and tunnels are of great significance for massive geological data mining and information acquisition. At present, the isosurface extraction and reconstruction method based on distance function has defects such as empty cube retrieval, stepped surface, huge number of triangles and redundancy, which will lead to problems such as low efficiency of reconstruction process and simulation display and surface discontinuity; while surface rendering based on contour splicing has higher response speed and less GPU occupancy, and can simulate and describe dynamic processes in real time, which better meets the modeling needs in the context of geological big data. After comprehensive consideration of constraint sets and feature metrics, it is concluded that the objective function of a single or the same type of feature metrics and constraint sets cannot achieve both local and global optimality, and is insufficient to effectively constrain the surface reconstruction process and maintain the rationality of the optimal solution in situations such as complex contour shapes, discrete edge length distribution, sparse vertices, and large differences in projection length along the horizontal and vertical axes of the plane; while the objective function optimization algorithm based on multiple feature metrics and constraint sets can effectively enhance the constraints on contour shape, multi-dimensional Euclidean distance of the winding direction and normal, and inter-layer differences. Significantly improving the control of the reconstruction process has become the main research direction in the field of objective functions. Therefore, in response to the above problems, this method is based on the objective function optimal solution algorithm of multiple feature metrics and constraint sets, and uses local optimal constraints to constrain the unstudied object. A contour reconstruction algorithm based on normal approximation and multi-feature fusion is proposed to perform refined three-dimensional modeling of complex geology of open-ground wells. The interpretation of the broken and water-rich anomalies in this area is realized, providing data reference for the advanced prediction of adverse geology and safe tunnel construction.

[0004] A Chinese patent with publication number CN202010706373.6 and publication date of March 1, 2024 discloses a fast and refined three-dimensional geological modeling method based on multi-source data, which can accurately express the interior of the entity while taking into account the accurate expression of the surface of the geological body.

[0005] However, the patent uses the interpolation method to interpolate vector points at the grid points of the regular grid on each layer interface. This method lacks quantitative analysis based on statistics in the fusion of multiple features, resulting in a reduction in the correlation between the fusion of constraint balance and feature metrics of different rules and the distribution rationality of sparse feature point sets, making it difficult to achieve fine modeling of comprehensive geological information of open-ground wells. Therefore, a method for fine modeling of complex geology of open-ground wells is needed today that can overcome the problems of constraint defects, matching disorder caused by differences in inflection point distribution, and discrete distribution of circumferential directions. Summary of the invention

[0006] The purpose of the present invention is to overcome the problems of constraint defects, matching disorder caused by differences in inflection point distribution, and discrete directional distribution when performing three-dimensional modeling in the prior art, and to provide a method for fine modeling of complex geology of open-ground wells.

[0007] In order to achieve the above-mentioned object of the invention, the present invention provides the following technical solutions:

[0008] A complex geological refinement modeling method for an open-field well comprises the following steps:

[0009] S1: Acquire the vacant well data to be analyzed; the vacant well data to be analyzed includes a set of discrete contour points;

[0010] S2: define a winding vertex density, take the W-configuration triangulation as the target, determine the number of interpolation based on the winding vertex density, optimize the contour discrete point set through distance weighted interpolation, and unify the vertex density of different contour lines;

[0011] S3: defining a spatial quadrilateral as a basic unit for establishing a similarity function, and constructing a feature metric of the spatial quadrilateral;

[0012] S4: calculating the characteristic measure of the spatial quadrilateral according to the contour discrete point set, establishing a sample set matrix of the optimal spatial quadrilateral, and performing centralization processing on the sample set matrix;

[0013] S5: Generate a covariance matrix according to the sample set matrix after centralization, calculate the eigenvalues ​​and feature measures of the covariance matrix, calculate the PCA contribution, perform feature dimensionality reduction according to the PCA contribution, solve the optimal solution with the maximum similarity, and establish vertex matching with global optimal constraints;

[0014] S6: establishing a mesh composed of the spatial quadrilaterals according to the vertex matching result, performing normal calculation on two triangulation situations in a single spatial quadrilateral, calculating similarities of adjacent triangles according to the normals, and establishing a triangulated mesh with local optimal constraints based on the maximum similarity;

[0015] S7: After all grids have been calculated, the current triangulated network is output as a complex geological refinement model.

[0016] As a preferred embodiment of the present invention, the expression of the winding vertex density is:

[0017]

[0018] Constraints: ρ Contour1 =ρ Contour2 ≈C,

[0019] Where ρ is the vertex density, n is the number of feature points, side i is the distance between the adjacent points starting from i, C∈Z∩C>0, C is a constant, ρ Contour1 and ρ Contour2 They are the vertex densities of contours Contour1 and Contour2 respectively.

[0020] As a preferred solution of the present invention, the interpolation based on the winding vertex density in S2 includes the following steps:

[0021] S21: judging the sparseness and density of the geophysical exploration slices in the open-ground well data to be analyzed according to the vertex density, and setting a reference number N for interpolation;

[0022] S22: Fill the side length set of adjacent contour lines in the contour discrete point set, calculate the perimeters C1 and C2 of adjacent contour lines, and calculate the interpolation numbers of two contour lines based on the reference number N and the perimeter ratio:

[0023] N1 / N2=C1 / C2;

[0024] Among them, N1 and N2 are the interpolation numbers of adjacent contour lines;

[0025] S23: except for the first and last contour lines, copy the point sets of each contour line to the auxiliary space in turn, fill the edge length set of the current contour line, search for the maximum edge length and perform equal interpolation on the edge length set of the current contour line; update the edge length set of the current contour line;

[0026] Repeat this step until the number of points on the current contour line reaches N1;

[0027] S24: Traverse all contour lines until the number of points of all contour lines reaches N1.

[0028] As a preferred solution of the present invention, the spatial quadrilateral is composed of two groups of adjacent points of the same contour in adjacent contours, and its feature measurements include:

[0029] ① The line vector from the geometric center of the spatial quadrilateral outline to the feature point as the direction vector of the point;

[0030] ② The lines connecting the vertices of the starting span and the adjacent points of the contour in the clockwise direction are used as the two short sides of the spatial quadrilateral, and the lines connecting the starting span and the corresponding adjacent points are used as the long sides of the spatial quadrilateral;

[0031] ③ The direction vector of the line where the short side is located is used as the direction vector of the side, and the direction points to the adjacent point;

[0032] ④ The three-dimensional Euclidean distance d of the straight line on which the short side lies.

[0033] As a preferred solution of the present invention, the construction of the sample set matrix includes the following steps:

[0034] S41: Count the number of vertices of adjacent profiles as m and n respectively, and establish a feature point direction association degree set {A}=[A1, A2,…, A n ] and the distance association degree set {B}=[B1,B2,…,B n ], the edge direction association degree set {C} = [C1, C2, …, C n ] and the edge distance association degree set {D}=[D1,D2,…,D n ]; its expression is:

[0035]

[0036]

[0037] Among them, a ij is the directional correlation between feature point i and feature point j, b ij is the distance correlation between feature point i and feature point j; c ij is the correlation degree between feature point i and feature point j; d ij is the edge distance correlation between feature point i and feature point j, i∈[1,m]∩j∈[1,n];

[0038] S42: Quantize the set of each feature metric, and the expression is:

[0039] The two-dimensional coordinates of the feature point of the starting side of the spatial quadrilateral are (x i ,y i ),(x j ,y j ) Its direction vector is:

[0040]

[0041] The directional correlation between two points is:

[0042] The direction vector of the short side is:

[0043]

[0044] The direction association of the edge is:

[0045] Distance between two points on different contours ij for:

[0046]

[0047] Contour distance contour for:

[0048] distance contour =|Z contour1 -Z contour2 |

[0049] The correlation degree of the starting point distance is: c ij =distance ij -distance contour

[0050] The midpoint of the short side of the spatial quadrilateral (x' i ,y' i ) or (x' j ,y' j )for:

[0051]

[0052] The projection distance correlation is:

[0053] As a preferred embodiment of the present invention, the method for calculating the PCA contribution in S5 comprises the following steps:

[0054] S511: Establish a sample data matrix set G = {(A, B, C, D)} according to the sample set matrix; its mathematical expectation is the characteristic variable μ a ,μ b ,μ c ,μ d , whose expression is:

[0055]

[0056] Among them, a i ∈A,b i ∈B,c i ∈C,d i ∈D;

[0057] S512: Center the sample data and calculate the covariance matrix C of the characteristic variables, which is expressed as:

[0058] Centralization formula:

[0059]

[0060] Covariance matrix C:

[0061]

[0062] Among them, cov(X,Y) is the covariance of variables X and Y, X i and Y i are the variables of the X and Y matrices, respectively, with index i, and are the mathematical expectations of X and Y, respectively, a, b, c, d are four characteristic measures of the spatial quadrilateral;

[0063] S513: Calculate the eigenvalue of the covariance matrix and normalize it with the eigenvector. The expression is:

[0064] The eigenvalues ​​of the covariance matrix are:

[0065] eigenvalue=(λ1,λ2,λ3,λ4),

[0066] Its feature matrix is:

[0067]

[0068]

[0069] Among them, λ1, λ2, λ3, and λ4 are the eigenvalues ​​of the principal components PC1, PC2, PC3, and PC4 respectively;

[0070] S514: Determine the number of principal components based on the cumulative contribution rate and the variance contribution rate and construct the corresponding projection mapping matrix, and calculate the fused feature data set after PCA:

[0071]

[0072] Among them, PC1 m×n is the value of the m×nth maximum contribution variable**, PC2 m×n is the value of the m×nth second contribution variable**, A m×n is the variance contribution rate of the eigenvalue of the m×nth point, B m×n is the cumulative contribution rate of the eigenvalue of the m×nth point, are the eigenvectors corresponding to the two most important components;

[0073] S515: The fused new features can be used as variables of the objective function and measures of constraint rules in contour reconstruction.

[0074] As a preferred solution of the present invention, the vertex matching of establishing global optimal constraints in S5 includes the following steps:

[0075] ① Calculate the direction vector of each contour vertex according to the initial storage order of the linked list Traverse all points and get the direction vector set of the vertices of the two contour lines

[0076] ② According to the initial storage order of the linked list, calculate the point set Points of the i-th contour i [0] Points to the i+1th contour i+1 The distance d of all points on [] i , get the span matrix of the two contour lines

[0077] ③Traverse the values ​​obtained in step ① Calculate the standard direction difference Get the direction difference matrix ΔV mn

[0078]

[0079] ④Center the span matrix and vector difference matrix, perform PCA principal component analysis based on least squares constrained fitting, and obtain a fusion feature matrix that takes into account vertex directional distribution and profile normal distance

[0080] ⑤Search for fusion feature Q mn The minimum value of the corresponding index i, j, Points i [], Points i+1 [] Adjusted to Points i [i] and Points i+1 [j] is the first element, connecting two points to obtain the first control span;

[0081] ⑥Remove Points i [i] and Points i+1 [j] Repeat the above steps to traverse the contour line set with fewer vertices, and the algorithm ends.

[0082] As a preferred solution of the present invention, in S6, normal calculation is performed on two possible triangulation situations in a single spatial quadrilateral, and the specific calculation method is as follows:

[0083] The vertex coordinates of the triangle patch are:

[0084]

[0085] The feature metric is the normal vector of the triangle patch, expressed as:

[0086]

[0087]

[0088] The constraint set is:

[0089]

[0090] As a preferred solution of the present invention, the triangulated network splicing algorithm in S6 has the following specific steps:

[0091] ① From the starting span, generate the starting triangle along the detour direction according to the shortest diagonal rule;

[0092] ② Starting from the current span, calculate the current point p i Unit line vectors connecting the front and back points judge If it is established, execute step ④, otherwise execute step ③;

[0093] ③Calculate the weight of the circumference occupied by the knitted parts of the two contours Retrieve the contour with smaller weight and calculate its unit edge weight P = 1 / N, calculate the weight difference If P ≥ Δ p Then execute step ④, otherwise calculate the number of compensation points Q = round (Δ p ×N), and then i ,p i+1 …p i+Q} and p j Connect to generate new triangles and add them to the Mesh; side h is the perimeter of the hth contour, m and n are the number of vertices of the adjacent profile contours, i and j are the ordinal numbers of the points in the feature point set of the adjacent contour lines, and N is the number of vertices of the contour with a smaller weight;

[0094] ④ Starting from the current triangle, calculate the unit normal vectors of the current triangle and the two adjacent triangles along the direction of the detour and and calculate Select the adjacent triangle with the larger absolute value of the product as the new triangle and add it to the Mesh;

[0095] ⑤ Repeat steps ②③④. When the points on a contour line that are in the current control interval are traversed, the remaining points on the other contour line are connected to the last point in the current control interval of the contour line to generate a new triangle and add it to the Mesh;

[0096] ⑥ Repeat step ⑤ until all control intervals are traversed and the current Mesh is output as the final triangulated network.

[0097] A device for fine-grained modeling of complex geology of an open-air well comprises at least one processor and a memory in communication with the at least one processor; the memory stores instructions executable by the at least one processor, and the instructions are executed by the at least one processor so that the at least one processor can execute any of the methods described above.

[0098] Compared with the prior art, the present invention has the following beneficial effects:

[0099] This application proposes a new method for fine modeling of complex geology of open-ground wells during the tunnel survey stage. The interpolation based on the perimeter ratio of the W configuration can effectively optimize the point set, and while solving the sparse and discrete distribution of the detour, it significantly weakens the detour difference of adjacent profile feature points, laying the foundation for strengthening the constraint of the feature point set on the profile shape and reconstructing the similarity between triangles in the triangulation network. It can efficiently and accurately establish a three-dimensional physical property model and geological anomaly model of the study area, realize the interpretation of the broken and water-rich anomalies in the area, and provide reference data for the advanced prediction of adverse geology and safe tunnel construction. BRIEF DESCRIPTION OF THE DRAWINGS

[0100] Figure 1 This is a schematic diagram of a process of a complex geological refinement modeling method for an open-ground well according to Example 1 of the present invention;

[0101] Figure 2 A schematic diagram showing how to achieve consistency in corresponding side lengths after interpolation in a complex geological refinement modeling method for an open-ground well according to Example 2 of the present invention;

[0102] Figure 3 A schematic diagram illustrating the approximation problem when the Euclidean distances are equal in a complex geological refinement modeling method for an open-ground well according to Example 2 of the present invention;

[0103] Figure 4 A schematic diagram showing the spatial quadrilateral and characteristic measurement in a complex geological detailed modeling method for an open-ground well according to Example 2 of the present invention;

[0104] Figure 5 A schematic diagram illustrating the PCA process in a complex geological refinement modeling method for an open-ground well according to Example 2 of the present invention;

[0105] Figure 6 It is a schematic diagram for demonstrating the concave-convex contour model in the complex geological detailed modeling method for an open-ground well described in Example 2 of the present invention;

[0106] Figure 7 A schematic diagram of a triangle point index storage structure for rendering normals in a complex geological refined modeling method for an open-ground well described in Example 2 of the present invention;

[0107] Figure 8 A schematic diagram showing the comparison between the adaptive compensation control and the conventional algorithm in a complex geological refined modeling method for an open-ground well according to Example 2 of the present invention;

[0108] Fig. 9 A schematic diagram of a weight compensation mechanism process based on adaptive inflection point identification in a complex geological refined modeling method for an open-ground well according to Example 2 of the present invention;

[0109] Fig.10 A schematic diagram of mosaic roughness in a complex geological refinement modeling method for an open-ground well described in Example 3 of the present invention;

[0110] Fig.11 This is a schematic diagram of the distribution of adjacent contour feature points in a complex geological refined modeling method for an open-ground well described in Example 3 of the present invention;

[0111] Fig.12 This is a schematic diagram of the comparison of the distance distribution of adjacent points before and after interpolation in a complex geological refined modeling method for an open-ground well according to Example 3 of the present invention;

[0112] Fig.13 A schematic diagram of the distribution characteristics of PC1 in a complex geological refined modeling method for an open-ground well according to Example 3 of the present invention;

[0113] Fig.14 A complex geological refinement model diagram of a study area in a complex geological refinement modeling method for an open-ground well described in Example 3 of the present invention;

[0114] Fig.15 This is a schematic structural diagram of a device for fine-grained modeling of complex geology for open-ground wells described in Example 4 of the present invention, which utilizes the method for fine-grained modeling of complex geology for open-ground wells described in Example 1. DETAILED DESCRIPTION

[0115] The present invention is further described in detail below in conjunction with test examples and specific implementation methods. However, this should not be understood as the scope of the above subject matter of the present invention being limited to the following embodiments, and all technologies realized based on the content of the present invention belong to the scope of the present invention.

[0116] Example 1

[0117] like Figure 1 As shown, a complex geological refinement modeling method for an open-ground well includes the following steps:

[0118] S1: Acquire the vacant well data to be analyzed; the vacant well data to be analyzed includes a set of discrete contour points.

[0119] S2: Define a winding vertex density, take the W-configuration triangulation as the target, determine the interpolation quantity based on the winding vertex density, optimize the contour discrete point set through distance weighted interpolation, and unify the vertex density of different contour lines.

[0120] S3: Define a spatial quadrilateral as a basic unit for establishing a similarity function, and construct a feature metric of the spatial quadrilateral.

[0121] S4: Calculate the characteristic measure of the spatial quadrilateral according to the contour discrete point set, establish a sample set matrix of the optimal spatial quadrilateral, and perform centralization processing on the sample set matrix.

[0122] S5: Generate a covariance matrix based on the sample set matrix after centralization, calculate the eigenvalues ​​and feature measures of the covariance matrix, calculate the PCA contribution, perform feature dimensionality reduction based on the PCA contribution, solve the optimal solution with the maximum similarity, and establish vertex matching with global optimal constraints.

[0123] S6: Establishing a mesh composed of the spatial quadrilaterals according to the vertex matching results, performing normal calculations on two triangulation situations in a single spatial quadrilateral, calculating similarities of adjacent triangles based on the normals, and establishing a locally optimal constrained triangulated mesh based on the maximum similarity.

[0124] S7: After all grids have been calculated, the current triangulated network is output as a complex geological refinement model.

[0125] Example 2

[0126] This embodiment is a specific implementation method of the complex geological refinement modeling method of an open-ground well described in Example 1, comprising the following steps:

[0127] S1: Obtain the open well data to be analyzed.

[0128] Furthermore, the open-ground well data to be analyzed include tunnel aerial geophysical data, ground geophysical data, well logging data, tunnel advance prediction data and / or excavation geological data and contour discrete point sets.

[0129] S2: Define a winding vertex density, take the W-configuration triangulation as the target, determine the interpolation quantity based on the winding vertex density, optimize the contour discrete point set through distance weighted interpolation, and unify the vertex density of different contour lines.

[0130] The shape of the basic triangle facets and the shape and perimeter approximation of the triangle facets with adjacency and sub-adjacency relationships are the key factors in achieving local surface smoothness; restricted by the Fuchs rule, the solution set of the local optimal function of any triangle facet has only two points, so the distribution of the points will have an important impact on the rationality of the optimal solution.

[0131] The distance between adjacent points of the same contour feature point set is an important measure to characterize the point distribution characteristics and the change in the length of the triangle base. By introducing the adjacent distance distribution function, the differences between adjacent contours can be quantified, such as Figure 4 As shown in (a): Both are distributed discretely, with no local intersection and low overall similarity. This will reduce the similarity and congruence of adjacent triangles, significantly affecting the smooth transition of the triangulated network. This paper introduces a measure of the distance between adjacent points, the winding density ρ, as a feature measure to optimize the point set and construct local constraints. Its expression is:

[0132]

[0133] The optimization is to perform interpolation optimization on the feature point set based on the winding density ρ; at this time, the constraint condition of local control is:

[0134] ρ Contour1 =ρ Contour2 ≈C,

[0135] C∈Z∩C>0,

[0136] Where ρ is the vertex density, n is the number of feature points, side i is the distance between the neighboring points starting from i, ρ Contour1 and ρ Contour2 are the vertex densities of contours Contour1 and Contour2 respectively, C is a constant, and Z is a set of integers.

[0137] As attached Figure 2 As shown in (b), (c), and (d), the feature point set is interpolated and optimized based on the winding density. As the number of points increases, the balance and consistency of the point set are significantly improved.

[0138] Limited by the topological structure of the linked list of discrete feature points of the contour line, this interpolation algorithm has the problem of interpolation results being covered: except for the first and last contour lines, the remaining contour lines have the problem of double interpolation, that is, along the normal vector direction of the slice, the contour line will be interpolated with the two adjacent contour lines in succession, which will destroy the current topological structure of the contour line and cause the interpolation to be invalid. This problem is actually caused by the difference between the Mesh structure and the Tin structure. This embodiment introduces an auxiliary space to solve this problem, that is, all contour lines except the first and last contour lines are copied to the auxiliary space, and each interpolation is performed by sampling the auxiliary space of the original contour line. Although the two auxiliary spaces of the same contour line have different topological structures, the vertex coordinate information that controls the shape size and spatial position is the same, and the vertices of adjacent contour lines are not deduplicated, so that the output triangulated network can achieve effective fitting of the surface. Specifically, the interpolation based on the vertex density includes the following steps:

[0139] S21: judging the sparseness and density of the geophysical exploration slices in the open-ground well data to be analyzed according to the vertex density, and setting a reference number N for interpolation;

[0140] S22: Fill the side length set of adjacent contour lines in the contour discrete point set, calculate the perimeters C1 and C2 of adjacent contour lines, and calculate the interpolation numbers of two contour lines based on the reference number N and the perimeter ratio:

[0141] N1 / N2=C1 / C2;

[0142] Among them, N1 and N2 are the interpolation numbers of adjacent contour lines;

[0143] S23: except for the first and last contour lines, copy the point sets of each contour line to the auxiliary space in turn, fill the edge length set of the current contour line, search for the maximum edge length and perform equal interpolation on the edge length set of the current contour line; update the edge length set of the current contour line;

[0144] Repeat this step until the number of points on the current contour line reaches N1;

[0145] S24: Traverse all contour lines until the number of points of all contour lines reaches N1.

[0146] S3: Define a spatial quadrilateral as a basic unit for establishing a similarity function, and construct a feature metric of the spatial quadrilateral.

[0147] Searching for the starting span and generating the starting triangle from a set of discrete feature points are important contents of the global constraint of the local optimal objective function. The existing algorithm uses three-dimensional Euclidean distance as a measure, constructs the association matrix of the adjacent contour point set according to the absolute coordinates of the feature points, generates the starting span through the global nearest neighbor search and its evolution method, and then generates the starting triangle based on random winding and local constraints.

[0148] However, the existing algorithm has major defects: ① The information loss in the generation of the correlation matrix is ​​not considered: in the point set coordinate operation, the directional relationship represented by the three-dimensional coordinates is lost, and the Euclidean distance between two points only inherits the distance relationship represented by the three-dimensional coordinates, resulting in multiple solutions and lack of constraints on point distribution differences. Figure 3 As shown in (a), the rationality of the three solutions is not the same. Taking the direction into consideration can find a better rationality. ② The detour direction is not included in the constraint set: Due to the problems of offset accumulation and error propagation in the reconstruction process of the local optimal objective function, the rationality of the optimal solution caused by this accumulation and propagation is significantly reduced along the detour direction (the traversal order of the points).

[0149] Based on the above analysis, we can conclude that: ① The spatial quadrilateral where the starting triangle is located has a significant impact on the shape and tessellation of the starting triangle and its adjacent triangles. ② The spatial quadrilateral can be uniquely determined by the starting span and the winding direction, and both are indispensable.

[0150] Therefore, the present invention takes the spatial quadrilateral as the optimization object of the objective function to search for the optimal solution of the starting triangle, and selects the characteristic measurement that can completely, non-redundantly and quantitatively describe the properties of the spatial quadrilateral by deeply mining the distribution information represented by the coordinates of the feature points, the spatial parallel information of the adjacent sections, and the correlation information of each element of the triangular facet: ① The two-dimensional directional vector angle of the point describes the distribution difference; ② The three-dimensional Euclidean distance of the point describes the degree of far-near correlation; ③ The two-dimensional directional vector angle of the winding edge describes the distribution difference of the matching interval; ④ The three-dimensional Euclidean distance of the straight line where the winding edge is located describes the degree of far-near correlation of the matching interval; ⑤ The length difference of the winding edge represents the stability of various properties of the matching interval.

[0151] Among them, the spatial quadrilateral and characteristic measure are as shown in the attached Figure 4 As shown: (a) The line connecting the contour geometric center to the feature point The vector is used as the direction vector of the point. (b) The lines connecting the vertices of the starting span and the adjacent points of the contour along the clockwise direction are used as the two short sides of the spatial quadrilateral, and the lines connecting the starting span and the corresponding adjacent points are used as the long sides of the spatial quadrilateral. (c) The direction vector of the straight line where the short side is located is used as the direction vector of the side, and the direction points to the adjacent point. (d) The three-dimensional Euclidean distance d of the straight line where the short side is located.

[0152] S4: Calculate the characteristic measure of the spatial quadrilateral according to the contour discrete point set, establish a sample set matrix of the optimal spatial quadrilateral, and perform centralization processing on the sample set matrix.

[0153] The construction of the sample set matrix includes the following steps:

[0154] S41: Count the number of vertices of adjacent profiles as m and n respectively, and establish a feature point direction association degree set {A}=[A1, A2,…, A n ] and the distance association degree set {B}=[B1,B2,…,B n ], the edge direction association degree set {C} = [C1, C2, …, C n ] and the edge distance association degree set {D}=[D1,D2,…,D n ]; its expression is:

[0155]

[0156]

[0157] Among them, a ij is the directional correlation between feature point i and feature point j, b ij is the distance correlation between feature point i and feature point j; c ij is the correlation degree between feature point i and feature point j; d ij is the edge distance correlation between feature point i and feature point j, i∈[1,m]∩j∈[1,n];

[0158] S42: Quantize the set of each feature metric, and the expression is:

[0159] The two-dimensional coordinates of the feature point of the starting side of the spatial quadrilateral are (x i ,y i ),(x j ,y j ) Its direction vector is:

[0160]

[0161] The directional correlation between two points is:

[0162] The direction vector of the short side is:

[0163]

[0164] The direction association of the edge is:

[0165] Distance between two points on different contours ij for:

[0166]

[0167] Contour distance contour for:

[0168] distance contour =|Z contour1 -Z contour2 |

[0169] The correlation degree of the starting point distance is: c ij =distance ij -distance contour

[0170] The midpoint of the short side of the spatial quadrilateral (x' i ,y' i ) or (x' j ,y' j )for:

[0171]

[0172] The projection distance correlation is:

[0173] S5: Generate a covariance matrix based on the sample set matrix after centralization, calculate the eigenvalues ​​and feature measures of the covariance matrix, calculate the PCA contribution, perform feature dimensionality reduction based on the PCA contribution, solve the optimal solution with the maximum similarity, and establish vertex matching with global optimal constraints.

[0174] S511: Establish a sample data matrix set G = {({A}, {B}, {C}, {D})} according to the multi-feature constraint conditions, whose feature variable is μ a ,μ b ,μ c ,μ d .

[0175] The PCA solution process of the optimal spatial quadrilateral with multiple feature constraints is as follows: According to {A} = [A1, A2, …, A n ],{B}=[B1,B2,…,B n ],{C}=[C1,C2,…,C n ],{D}=[D1,D2,…,D n ] Establish a sample data matrix set G = {(A, B, C, D)}, a i ∈A,b i ∈B,c i ∈C,d i ∈D.

[0176] Furthermore, the multi-feature constraints are:

[0177]

[0178] a ij ∈{A},b ij∈ {B}, c ij ∈ {C}, d ij ∈ {D},

[0179] The set of sample data matrices is:

[0180] The characteristic variables are:

[0181]

[0182] Limited by the finiteness of the set and the finiteness of the search algorithm, the construction of the space-optimal quadrilateral should be transformed into the problem of finding the optimal solution of the constraint equation with a finite domain. There have been many studies on this type of mathematical model. PCA (principal component analysis) is a mathematical transformation model that can map a given set of n-dimensional features to a k-dimensional orthogonal space (k < n), and the quantity with the largest variance among these K new features is called the principal component. As shown in the appendix Figure 5 As shown, the four clusters are distributed in the left three-dimensional space according to three original features. By performing PCA dimensionality reduction, the two-dimensional space on the right can be obtained. If the vibration in the direction of the minor axis of the ellipse is ignored and the major axis of the ellipse is used as the fused principal component, feature fusion and dimensionality reduction can be achieved, which is beneficial to the separation of clustering. This method can analyze the weights of multiple features and avoid information loss and the decline in the rationality of the optimal solution caused by determining the feature priority based on prior theory and making selections and eliminations.

[0183] S312: Calculate the covariance matrix of the characteristic variables, and its expression is:

[0184]

[0185] Among them, cov(X, Y) is the covariance of variables X and Y, and X i and Y i are the variables of matrices X and Y respectively, and their index is i, and are the mathematical expectations of X and Y respectively, and a, b, c, and d are the four characteristic measures of the space quadrilateral mentioned above.

[0186] The sample data can be centralized through the characteristic variables to facilitate the calculation of the covariance matrix of the characteristic variables.

[0187] S513: Calculate the eigenvalues of the covariance matrix, and unitize them together with the eigenvectors. The expression is:

[0188] The eigenvalues of the covariance matrix are:

[0189] eigenvalue = (λ1, λ2, λ3, λ4),

[0190] Its feature matrix is:

[0191]

[0192]

[0193] Among them, λ1, λ2, λ3, and λ4 are the eigenvalues ​​of the principal components PC1, PC2, PC3, and PC4, respectively.

[0194] S514: Determine the number of principal components according to the cumulative contribution rate and the variance contribution rate and construct the corresponding projection mapping matrix, and calculate the fused feature data set after PCA, the expression of which is:

[0195]

[0196] Among them, PC1 m×n is the value of the m×nth maximum contribution variable, PC2 m×n is the value of the m×nth second contribution variable, A m×n is the variance contribution rate of the eigenvalue of the m×nth point, B m×n is the cumulative contribution rate of the eigenvalue of the m×nth point, are the eigenvectors corresponding to the two most important components respectively.

[0197] S515: The fused new features can be used as variables of the objective function and measures of constraint rules in contour reconstruction.

[0198] Furthermore, establishing the vertex matching of the global optimal constraint includes the following steps:

[0199] ① Calculate the direction vector of each contour vertex according to the initial storage order of the linked list Traverse all points and get the direction vector set of the vertices of the two contour lines

[0200] ② According to the initial storage order of the linked list, calculate the point set Points of the i-th contour i [0] Points to the i+1th contour i+1 The distance d of all points on [] i , get the span matrix of the two contour lines

[0201] ③Traverse the values ​​obtained in step ① Calculate the standard direction difference Get the direction difference matrix ΔV mn

[0202]

[0203] ④Center the span matrix and vector difference matrix, perform PCA principal component analysis based on least squares constrained fitting, and obtain a fusion feature matrix that takes into account vertex directional distribution and profile normal distance

[0204] ⑤Search for fusion feature Q mn The minimum value of the corresponding index i, j, Points i [], Points i+1 [] Adjusted to Points i [i] and Points i+1 [j] is the first element, connecting two points to obtain the first control span;

[0205] ⑥Remove Points i [i] and Points i+1 [j] Repeat the above steps to traverse the contour line set with fewer vertices, and the algorithm ends.

[0206] S6: establishing a mesh composed of the spatial quadrilaterals according to the vertex matching result, performing normal calculation on two triangulation situations in a single spatial quadrilateral, calculating similarities of adjacent triangles according to the normals, and establishing a triangulated mesh with local optimal constraints based on the maximum similarity;

[0207] This embodiment performs smooth tessellation based on the principle of approximation of the normal directions of adjacent faces. Based on the parallelism of the tangent planes of any adjacent points of the smooth surface, the reconstructed triangular facet set should meet the following constraints:

[0208]

[0209] i∈mesh, j∈mesh; lim{mesh}=∞

[0210] and is the unit normal vector of any adjacent triangle i, j of the triangulated mesh formed by adjacent contours. Limited by the accuracy of drilling instruments, the level of exploration technology and the influence of geological environment, the sparse profile obtained leads to significant differences in point distribution density along the winding direction and the normal direction. The non-smooth contour polyline causes the normal consistency of adjacent triangles to be discontinuous near the contour inflection point. Therefore, the above formula cannot be directly used as a constraint set. For non-continuous surface reconstruction, through the study of a large number of modeling examples, extensive literature analysis and algorithm simulation, it can be found that convex contours can obtain better reconstructed models. However, there are only two triangles in a solution set of the local objective function along the winding direction. Therefore, the algorithm should weaken the influence of concave contours while generating convex contours. As shown in the attached figure Figure 6As shown in the figure: The dihedral angle of adjacent faces with the common span as the axis under concave and convex contours is an important indicator to measure the fluctuation degree of contour concavity and transition. This algorithm uses the maximum dihedral angle of adjacent faces as the local constraint condition and has universal applicability. ① In the case of convex contours, the normal change of adjacent faces is small, which enhances the smoothness of surface transition. ② In the case of concave contours, the large dihedral angle weakens the amplitude of the concave and reduces the risk of local cross-section and misalignment. For concave contours, the minimum normal difference can weaken the influence of the concave area. For convex contours, the minimum normal difference can reduce the normal change between adjacent triangles and enhance the smoothness of the surface.

[0211] As attached Figure 6 As shown, the vertex coordinates of the triangle patch are:

[0212]

[0213] The feature metric is the normal vector of the triangle patch, expressed as:

[0214]

[0215]

[0216] The constraint set is:

[0217]

[0218] In order to unify the direction of the normal and avoid multiple solutions of the dihedral angle, this algorithm introduces a triangle point index storage structure that conforms to the rendering normal, as shown in the attached figure. Figure 7 As shown, the indices of the points are constructed in a counterclockwise manner to construct a triangle.

[0219] The present invention proposes a method for generating a control span based on the nearest neighbor search of the inflection point of the easier shape of adjacent contour lines to segment the matching interval, enhance the consistency of the matching interval in the winding trend and normal projection, weaken the influence of unreasonable spans, and eliminate the overall disorder of the matching correspondence and modeling errors caused by the uncertainty and randomness of the interpolation point and the non-starting span. The algorithm distinguishes the interpolation points and vertices of the contour line, as shown in the attached figure. Figure 8 and attached Fig. 9 As shown: the red points are the original vertices that represent the turning points and directions of the contours, and the blue points are the interpolation points. After unifying the direction of the contours, starting from the starting span, the turning points of the contours with low shape complexity are matched and searched one by one to the turning points of the contours with high complexity. The control span is generated according to the vertex index to control the triangulated interval, so as to effectively take into account the shape characteristics of the contours and enhance the robustness of the constraint rules.

[0220] Further, the following steps are included:

[0221] ① From the starting span, generate the starting triangle along the detour direction according to the shortest diagonal rule;

[0222] ② Starting from the current span, calculate the current point p i Unit line vectors connecting the front and back points judge If it is established, execute step ④, otherwise execute step ③;

[0223] ③Calculate the weight of the circumference occupied by the knitted parts of the two contours Retrieve the contour with smaller weight and calculate its unit edge weight P = 1 / N, calculate the weight difference If P ≥ Δ p Then execute step ④, otherwise calculate the number of compensation points Q = round (Δ p ×N), and then i ,p i+1 …p i+Q} and p j Connect to generate new triangles and add them to the Mesh; side h is the perimeter of the hth contour, m and n are the number of vertices of the adjacent profile contours, i and j are the ordinal numbers of the points in the feature point set of the adjacent contour lines, and N is the number of vertices of the contour with a smaller weight;

[0224] ④ Starting from the current triangle, calculate the unit normal vectors of the current triangle and the two adjacent triangles along the direction of the detour and and calculate Select the adjacent triangle with the larger absolute value of the product as the new triangle and add it to the Mesh;

[0225] ⑤ Repeat steps ②③④. When the points on a contour line that are in the current control interval are traversed, the remaining points on the other contour line are connected to the last point in the current control interval of the contour line to generate a new triangle and add it to the Mesh;

[0226] ⑥ Repeat step ⑤ until all control intervals are traversed and the current Mesh is output as the final triangulated network.

[0227] S7: After all grids have been calculated, the current triangulated network is output as a complex geological refinement model.

[0228] The difference between this embodiment and the above embodiment is that S7 also includes a model verification analysis process, which specifically includes the following steps:

[0229] S71: The present invention takes the triangular faces and their topological connections of the triangulated network that constitute the model as the research objects, combines the comparison of a large number of modeling examples with the analysis of the literature in the field of surface reconstruction, introduces a set of metrics that characterize the properties of the triangular faces in the triangulated network and the tessellation properties between the triangular faces, and uses the variable values ​​and the first and second order central moments to evaluate the local and global reconstructed surface respectively.

[0230] S72: The side length and angle of a triangle face represent its size and shape. Due to the sparseness of exploration data, the spacing between adjacent sections is much larger than the spacing between adjacent points on the section, making it difficult to form equilateral basic triangle faces. Isosceles triangle faces of similar shape and size are conducive to achieving smooth transition to the maximum extent, and are therefore used in triangulation. The distribution characteristics of span length are used as indicators to measure the proportion and size consistency of approximately isosceles triangles in a triangulation network; and the base and its corresponding angle are used as indicators to measure the shape of a single triangle face.

[0231] The topological relationship of the triangular facets is characterized by the tessellation of adjacent faces on the model, and the smoothness of the facets determines the fitting effect of the reconstructed surface. When the adjacent dihedral angle is an obtuse angle, the metric is negative, which will affect the effectiveness and intuitiveness of the overall data evaluation and increase the complexity of the analysis. The present invention proposes the concept of tessellation roughness as a measurement indicator. The roughness range is [0,2] and is negatively correlated with transition smoothness. When the roughness is greater than 1, unreasonable tessellation is prone to occur. Tessellation roughness is shown in the attached figure. Fig.10 As shown, the unit normal of the adjacent surface is The tessellation roughness formula is expressed as:

[0232]

[0233] S73: This section uses a set of feature points of a group of adjacent concave contour shape simulation slices as the data source, and the adjacent contour data distribution is as shown in the attached Fig.11 shown.

[0234] S74: Count the number and distribution of contour vertices, and its original neighbor point distance distribution function is shown in the attached figure. Fig.12 As shown in (a), the distance between the two points is discrete and the expected difference is large, and there is also the influence of extreme values.

[0235] For example, the interpolation algorithm based on directional density balance is used to increase the number of points of the two contour lines from 19 to 76 and from 17 to 69 respectively. The distribution function of the encrypted point spacing is shown in the attached figure. Fig.12 As shown in (b), the extreme values ​​are effectively eliminated, the distribution tends to be stable and the expectation is approximate, the oscillation amplitude along both ends is reduced from 1.75 to 0.5, and the covariance matrix of the balance and consistency feature measure of the point set distribution is:

[0236]

[0237] The eigenvalues ​​and eigenvectors of the covariance matrix are solved, and the four principal components to be selected are calculated, and their eigenvalues ​​and contribution rates are shown in Table 1. The quantities with eigenvalues ​​greater than 1 are selected as the principal components, which are significantly optimized.

[0238] Table 1 PCA contribution rates of the four principal components

[0239]

[0240] S75: Establish the correlation matrix of the starting direction, edge direction, starting distance, and end distance based on the optimized point set; calculate four PCAs. The fusion weights of the principal component and the original feature are shown in Table 2. Select the index of the variable corresponding to the PCA score closest to 0, which can be used as the starting span and starting quadrilateral. The starting triangle can be obtained by decomposition according to the shortest diagonal rule.

[0241] Table 2 Fusion weights of principal components and original features

[0242]

[0243] S76: Refined modeling after PCA analysis, Fig.13 is the distribution characteristic of PC1, Fig.14 For the built model.

[0244] Example 3

[0245] like Fig.15 As shown, a device for fine-grained modeling of complex geology of an open-ground well comprises at least one processor, a memory connected to the at least one processor in communication, and at least one input-output interface connected to the at least one processor in communication; the memory stores instructions executable by the at least one processor, and the instructions are executed by the at least one processor so that the at least one processor can execute the method for fine-grained modeling of complex geology of an open-ground well described in the above embodiment. The input-output interface may include a display, a keyboard, a mouse, and a USB interface for inputting and outputting data.

[0246] Those skilled in the art can understand that: all or part of the steps of implementing the above-mentioned method embodiment can be completed by hardware related to program instructions, and the aforementioned program can be stored in a computer-readable storage medium. When the program is executed, it executes the steps of the above-mentioned method embodiment; and the aforementioned storage medium includes: mobile storage devices, read-only memories (ROM), disks or optical disks, etc. Various media that can store program codes.

[0247] When the above-mentioned integrated unit of the present invention is implemented in the form of a software functional unit and sold or used as an independent product, it can also be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the embodiment of the present invention can be essentially or partly reflected in the form of a software product that contributes to the prior art. The computer software product is stored in a storage medium, including several instructions for a computer device (which can be a personal computer, a server, or a network device, etc.) to execute all or part of the methods described in each embodiment of the present invention. The aforementioned storage medium includes: various media that can store program codes, such as mobile storage devices, ROMs, magnetic disks, or optical disks.

[0248] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included in the protection scope of the present invention.

Claims

1. A complex geological refinement modeling method for open-air wells, characterized in that: The following steps are involved: S1: Acquire the vacant well data to be analyzed; the vacant well data to be analyzed includes a set of discrete contour points; S2: define a winding vertex density, take the W-configuration triangulation as the target, determine the number of interpolation based on the winding vertex density, optimize the contour discrete point set through distance weighted interpolation, and unify the vertex density of different contour lines; S3: defining a spatial quadrilateral as a basic unit for establishing a similarity function, and constructing a feature metric of the spatial quadrilateral; S4: calculating the characteristic measure of the spatial quadrilateral according to the contour discrete point set, establishing a sample set matrix of the optimal spatial quadrilateral, and performing centralization processing on the sample set matrix; S5: Generate a covariance matrix according to the sample set matrix after centralization, calculate the eigenvalues ​​and feature measures of the covariance matrix, calculate the PCA contribution, perform feature dimensionality reduction according to the PCA contribution, solve the optimal solution with the maximum similarity, and establish vertex matching with global optimal constraints; S6: establishing a mesh composed of the spatial quadrilaterals according to the vertex matching result, performing normal calculation on two triangulation situations in a single spatial quadrilateral, calculating similarities of adjacent triangles according to the normals, and establishing a triangulated mesh with local optimal constraints based on the maximum similarity; S7: After all grids have been calculated, the current triangulated network is output as a complex geological refinement model.

2. The complex geological detailed modeling method for open-air wells according to claim 1 is characterized in that: The expression of the winding vertex density is: Constraints: ρ Contour1 =ρ Contour2 ≈C, Where ρ is the vertex density, n is the number of feature points, side i is the distance between the neighboring points starting from i, C∈Z ∩C>0, C is a constant, ρ Contour1 and ρContou r2 They are the vertex densities of contours Contour1 and Contour2 respectively.

3. The complex geological detailed modeling method for open-air wells according to claim 1 is characterized in that: The interpolation based on the winding vertex density in S2 comprises the following steps: S21: judging the sparseness and density of the geophysical exploration slices in the open-ground well data to be analyzed according to the vertex density, and setting a reference number N for interpolation; S22: Fill the side length set of adjacent contour lines in the contour discrete point set, calculate the perimeters C1 and C2 of adjacent contour lines, and calculate the interpolation number of two contour lines according to the reference number N and the perimeter ratio: N1 / N2=C1 / C2; Among them, N1 and N2 are the interpolation numbers of adjacent contour lines; S23: except for the first and last contour lines, copy the point sets of each contour line to the auxiliary space in turn, fill the edge length set of the current contour line, search for the maximum edge length and perform equal interpolation on the edge length set of the current contour line; update the edge length set of the current contour line; Repeat this step until the number of points on the current contour line reaches N1; S24: Traverse all contour lines until the number of points of all contour lines reaches N1.

4. The complex geological detailed modeling method for open-air wells according to claim 1 is characterized in that: The spatial quadrilateral is composed of two groups of adjacent points of the same contour in adjacent contours, and its feature measures include: ① The line vector from the geometric center of the spatial quadrilateral outline to the feature point as the direction vector of the point; ② The lines connecting the vertices of the starting span and the adjacent points of the contour in the clockwise direction are used as the two short sides of the spatial quadrilateral, and the lines connecting the starting span and the corresponding adjacent points are used as the long sides of the spatial quadrilateral; ③ The direction vector of the line where the short side is located is used as the direction vector of the side, and the direction points to the adjacent point; ④ The three-dimensional Euclidean distance d of the straight line on which the short side lies.

5. The complex geological detailed modeling method for open-air wells according to claim 4 is characterized in that: The construction of the sample set matrix includes the following steps: S41: Count the number of vertices of adjacent profile contours as m and n respectively, and establish a feature point direction association degree set {A}=[A1, A2, ..., A n ] and the distance association degree set {B} = [B1, B2, ..., B n ], the edge direction association degree set {C} = [C1, C2, ..., C n ] and the edge distance association degree set {D} = [D1, D2, ..., D n ]; its expression is: Among them, a ij is the directional correlation between feature point i and feature point j, b ij is the distance correlation between feature point i and feature point j; c ij is the correlation degree between feature point i and feature point j; d ij is the edge distance correlation between feature point i and feature point j, i∈[1,m]∩j∈[1,n]; S42: Quantize the set of each feature metric, and the expression is: The two-dimensional coordinates of the feature point of the starting side of the spatial quadrilateral are (x i ,y i ), (x j ,y j ) Its direction vector is: The directional correlation between two points is: The direction vector of the short side is: The direction association of the edge is: Distance between two points on different contours ij for: Contour distance contour for: distance contour =|Z contour1 -FROM contour2 | The correlation degree of the starting point distance is: c ij =distance ij -distance contour The midpoint of the short side of the spatial quadrilateral (x′ i , y′ i ) or (x′ j , y′ j )for: The projection distance correlation is:

6. The complex geological detailed modeling method for open-air wells according to claim 5 is characterized in that: The method for calculating the PCA contribution in S5 comprises the following steps: S511: Establish a sample data matrix set G = {(A, B, C, D)} according to the sample set matrix; its mathematical expectation is the characteristic variable μ a , μ b , μ c , μ d , whose expression is: Among them, a i ∈A,b i ∈B,c i ∈C,d i ∈D; S512: Center the sample data and calculate the covariance matrix C of the characteristic variables, which is expressed as: Centralization formula: Covariance matrix C: Among them, cov(X,Y) is the covariance of variables X and Y, X i and Y i are the variables of the X and Y matrices, respectively, with index i, and are the mathematical expectations of X and Y, respectively, a, b, c, d are four characteristic measures of the spatial quadrilateral; S513: Calculate the eigenvalue of the covariance matrix and normalize it with the eigenvector. The expression is: The eigenvalues ​​of the covariance matrix are: eigenvalue=(λ1, λ2, λ3, λ4), Its feature matrix is: Among them, λ1, λ2, λ3, and λ4 are the eigenvalues ​​of the principal components PC1, PC2, PC3, and PC4 respectively; S514: Determine the number of principal components based on the cumulative contribution rate and the variance contribution rate and construct the corresponding projection mapping matrix, and calculate the fused feature data set after PCA: Among them, PC1 m×n is the value of the m×nth maximum contribution variable**, PC2 m×n is the value of the m×nth second contribution variable**, A m×n is the variance contribution rate of the eigenvalue of the m×nth point, B m×n is the cumulative contribution rate of the eigenvalue of the m×nth point, are the eigenvectors corresponding to the two most important components; S515: The fused new features can be used as variables of the objective function and measures of constraint rules in contour reconstruction.

7. The complex geological detailed modeling method for open-air wells according to claim 1 is characterized in that: The vertex matching of establishing the global optimal constraint in S5 comprises the following steps: ① Calculate the direction vector of each contour vertex according to the initial storage order of the linked list Traverse all points and get the direction vector set of the vertices of the two contour lines ② According to the initial storage order of the linked list, calculate the point set Points of the i-th contour i [0] Points to the i+1th contour i+1 The distance d of all points on [] i , get the span matrix of the two contour lines ③Traverse the values ​​obtained in step ① Calculate the standard direction difference Get the direction difference matrix ΔV mn ④Center the span matrix and vector difference matrix, perform PCA principal component analysis based on least squares constrained fitting, and obtain a fusion feature matrix that takes into account vertex directional distribution and profile normal distance ⑤Search fusion feature Q m The minimum value of the corresponding index i, j, Points i [], Points i+1 [] Adjusted to Points i [i] and Points i+1 [j] is the first element, connecting two points to obtain the first control span; ⑥Remove Points i [i] and Points i+1 [j] Repeat the above steps to traverse the contour line set with fewer vertices, and the algorithm ends.

8. The complex geological detailed modeling method for open-air wells according to claim 1 is characterized in that: In s6, normal calculation is performed for two possible triangulation situations in a single spatial quadrilateral, and the specific calculation method is as follows: The vertex coordinates of the triangle patch are: The feature metric is the normal vector of the triangle patch, expressed as: The constraint set is:

9. The complex geological detailed modeling method for open-ground wells according to claim 8 is characterized in that: The triangulated network splicing algorithm in S6 has the following specific steps: ① From the starting span, generate the starting triangle along the detour direction according to the shortest diagonal rule; ② Starting from the current span, calculate the unit line vector between the current point pi and the previous and next points judge If it is established, execute step ④, otherwise execute step ③; ③Calculate the weight of the circumference occupied by the knitted parts of the two contours Retrieve the contour with smaller weight and calculate its unit edge weight P = 1 / N, calculate the weight difference If P ≥ Δ p Then execute step ④, otherwise calculate the number of compensation points Q = round (Δ p ×N), and successively i , p i+1 …p i+Q } and p j Connect to generate new triangles and add them to the Mesh; side h is the perimeter of the hth contour, m and n are the number of vertices of the adjacent profile contours, i and j are the serial numbers of the points in the feature point set of the adjacent contour lines, and N is the number of vertices of the contour with a smaller weight; ④ Starting from the current triangle, calculate the unit normal vectors of the current triangle and the two adjacent triangles along the direction of the detour and and calculate Select the adjacent triangle with the larger absolute value of the product as the new triangle and add it to the Mesh; ⑤ Repeat steps ②③④. When the points on a contour line that are in the current control interval are traversed, the remaining points on the other contour line are connected to the last point in the current control interval of the contour line to generate a new triangle and add it to the Mesh; ⑥ Repeat step ⑤ until all control intervals are traversed and the current Mesh is output as the final triangulated network.

10. A complex geological detailed modeling device for open-ground wells, characterized in that: It includes at least one processor and a memory communicatively connected to the at least one processor; the memory stores instructions executable by the at least one processor, and the instructions are executed by the at least one processor so that the at least one processor can execute the method described in any one of claims 1 to 9.

Images