A method, apparatus and equipment for rapid modeling of shield tunneling areas in complex geological formations

By acquiring point set data to construct geological domain entities and performing adaptive meshing and layer assignment, the problem of insufficient accuracy in stratum modeling was solved, enabling rapid modeling of complex strata and automated model updates, thereby improving the accuracy of tunnel boring machine (TBM) construction.

CN121505208BActive Publication Date: 2026-04-03SHENZHEN UNIV
View PDF 2 Cites 0 Cited by

Patent Information

Authority / Receiving Office
CN · China
Patent Type
Patents(China)
Current Assignee / Owner
Filing Date
2026-01-12
Publication Date
2026-04-03

AI Technical Summary

Technical Problem

Existing geological modeling methods cannot effectively reflect the nonlinear fluctuations and local abrupt changes in strata, resulting in insufficient model accuracy and the inability to achieve automated construction from point set data to "geometry-mesh-layer assignment", which seriously affects the real-time performance and consistency of geological models in shield tunneling areas.

Method used

By acquiring point set data containing three-dimensional coordinates and stratigraphic labels, a three-dimensional entity of the geological domain is constructed, Boolean cutting and mesh generation are performed, a voxel-level stratigraphic prior field is established, stratigraphic categories are adaptively assigned, and a finite element model file is generated.

Benefits of technology

It improved the accuracy and realism of the geological model in the shield tunneling area, realized the automated construction of the geological domain and the rapid updating of the model, and improved the real-time performance and consistency of the model.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN121505208B_ABST
    Figure CN121505208B_ABST
Patent Text Reader

Abstract

This invention discloses a rapid modeling method, apparatus, and equipment for shield tunneling areas in complex geological formations, relating to the fields of 3D geological modeling, tunnel and underground engineering technology. Driven by 3D point set data containing stratigraphic categories, the invention constructs a modeling domain through data verification and spatial bounding box calculation, generating geological entities and shield tunnel models. The tunnel space is then subtracted to form a geological modeling domain containing chambers. A continuous grid size field is constructed based on the tunnel distance field and stratigraphic complexity, automatically determining the partition width and representative size to achieve adaptive zoning and densification. Simultaneously, an adaptive stratigraphic determination mechanism is constructed, automatically assigning stratigraphic values ​​to grid units for stable and continuous operation through joint decision-making. This invention avoids the problems of boundary ambiguity and stratigraphic misjudgment, improves the accuracy and consistency of complex geological modeling, supports dynamic model updates, and can output 3D geological models that can be directly used for construction analysis and parameter optimization, providing efficient and reliable technical support for all aspects of shield tunneling construction.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This invention relates to the fields of three-dimensional geological modeling, tunnel and underground engineering technology, and in particular to a method, apparatus and equipment for rapid modeling of shield tunnel areas in complex strata. Background Technology

[0002] In tunnel and underground engineering construction, to perform in-situ stress analysis, stratum stability assessment, and shield tunneling parameter optimization, it is usually necessary to establish a three-dimensional geological model containing cavern structure and layer information in numerical analysis platforms such as finite element methods. However, existing stratum modeling often only extracts stratum depth based on limited geological borehole data and delineates strata boundaries through linear or planar extrapolation. This simplified approach fails to reflect the nonlinear fluctuations and local abrupt changes in the strata in space, resulting in insufficient model accuracy and significant deviations from reality. In shield tunneling scenarios, geological conditions change frequently and dynamically. Existing solutions cannot achieve an automated construction process from point set data to "geometry-mesh-layer assignment," severely restricting the real-time performance and consistency of model updates.

[0003] Therefore, there is an urgent need for an automated method that can efficiently complete geological domain construction, mesh generation, and layer assignment based directly on 3D point set data with stratum labels, so as to improve the accuracy, realism, and update efficiency of geological models in shield tunneling areas. Summary of the Invention

[0004] To address the aforementioned technical problems in related technologies, this invention proposes a method, apparatus, and equipment for rapid modeling of shield tunneling areas in complex geological formations.

[0005] In a first aspect, the present invention provides a method for rapid modeling of shield tunneling areas in complex geological formations, comprising the following steps:

[0006] S1. Obtain point set data and lithological attribute data containing three-dimensional coordinates and stratigraphic labels, and determine the three-dimensional envelope boundary of the modeling domain based on the coordinate range of the point set.

[0007] S2. Construct a three-dimensional entity of the geological domain based on the three-dimensional envelope boundary, and construct a tunnel cylinder based on the shield tunnel parameters and spatially align it with the geological domain;

[0008] S3. Perform Boolean cutting on the three-dimensional entity of the geological domain and the tunnel cylinder to obtain a geological domain modeling entity containing the tunnel chamber space;

[0009] S4. Construct a continuous grid size field based on the tunnel distance field and stratigraphic complexity of the geological domain modeling entity, and then generate grid cells by adaptive grid division and zoning densification control of the geological domain modeling entity based on the continuous grid size field.

[0010] S5. In the geological domain modeling entity, construct a voxel grid with a step size not less than the global grid control size. Perform stratigraphic clustering statistics on the point set in each voxel grid to obtain stratigraphic prior information. Establish a voxel-level stratigraphic prior field that maps the voxel grid to stratigraphic prior information.

[0011] S6. For each grid cell in the geological domain modeling entity, determine its voxel cell by using its centroid coordinates as the query point. Then, retrieve the stratigraphic prior information corresponding to the voxel cell based on the voxel-level stratigraphic prior field. If the confidence level of the stratigraphic prior information is greater than or equal to a preset confidence threshold, then use the stratigraphic category of the stratigraphic prior information as the stratigraphic category of the grid cell; otherwise, determine the stratigraphic category of the grid cell based on the nearest neighbor retrieval and distance-weighted voting strategy.

[0012] S7. Automatically generate the material and cross-section of the geological domain modeling entity according to the layer category and assign values ​​to each layer;

[0013] S8. Output a finite element model file containing mesh topology and layered assignment information based on the geological domain model entity, and use the finite element model file for shield tunneling construction analysis.

[0014] Specifically, step S2 includes: the shield tunnel parameters include the tunnel radius R, centerline position parameters, and axis; the construction of the tunnel cylinder based on the shield tunnel parameters and spatial alignment with the geological domain specifically includes: generating a tunnel cylinder with the same spatial orientation as the tunnel radius R and centerline position parameters, and realizing the spatial positioning of the tunnel configuration in the geological domain through axis alignment and position registration.

[0015] Specifically, step S4, which involves constructing a continuous grid size field based on the tunnel distance field and stratigraphic complexity of the geological domain modeling entity, includes the following steps:

[0016] S41. First, apply global mesh control dimensions within the geological domain modeling entity to determine the baseline dimensions of the overall mesh of the geological domain modeling entity;

[0017] S42. Then, using the tunnel wall as the zero level set, calculate the minimum Euclidean distance from each point in the geological modeling entity to the tunnel wall to construct the tunnel distance field.

[0018] S43. Construct a stratigraphic complexity index C(x) based on the degree of stratigraphic variation and local lithological fluctuations of the point set, as shown in the following formula:

[0019]

[0020] Among them, C layer (x) represents the dispersion of the stratigraphic variation within the local neighborhood of point x, used to characterize the degree of stratigraphic variation of the point set; C litho(x) represents the degree of fluctuation of lithological properties in the local neighborhood of point x, and is used to characterize the local lithological fluctuation of the point set; λ is the stratigraphic variation weighting coefficient;

[0021] S44. Introduce the tunnel distance field and the geological complexity index into the grid size continuous function to form a continuous grid size field S(x);

[0022] The grid size continuity function is shown in the following equation:

[0023] ;

[0024] Among them, h min h is the minimum allowed mesh size for the modeling domain. max The maximum allowed mesh size for the modeling domain; The encrypted weight of point x is based on the tunnel distance field; This represents the maximum value in the tunnel distance field. This represents the value corresponding to the distance x from the midpoint of the field to the tunnel.

[0025] The encryption weight of point x is based on the stratigraphic complexity index; This represents the maximum value of the stratigraphic complexity index. This represents the minimum value of the stratigraphic complexity index. Let x be the stratigraphic complexity index at point x.

[0026] Specifically, the process of adaptively dividing and zoning the geological domain modeling entity based on the continuous grid size field to generate grid cells includes the following steps:

[0027] S45. Select M representative points on the tunnel wall. Starting from each representative point, sample along its unit outward normal direction at preset fixed intervals to obtain M target size curves L(t). Then, statistically process the M target size curves L(t) according to their distance and position to obtain representative size curves. The target size curves L(t) are shown in the following formula:

[0028] ;

[0029] Where x0 is a representative point selected on the tunnel wall, n(x0) is the unit outward normal of the representative point x0, t is the preset sampling distance along the outward normal; S(·) represents the continuous grid size field;

[0030] S46. Based on the representative dimension curve, determine the boundary points of the three zones according to the trend of dimension change with normal distance, divide the three zones according to the boundary points of the three zones and obtain the distance intervals corresponding to the three zones; the three zones include the wall refinement zone, the inner zone and the outer zone;

[0031] S47. On the representative size curve, calculate the average or median of the target size values ​​in each distance interval of the three zones to obtain the representative size of each distance interval, and use them as the target mesh size of the three zones respectively.

[0032] S48. Construct zone boundaries using isosurfaces of the distance field based on the distance intervals corresponding to the three zones. Apply the corresponding target grid size to the distance intervals corresponding to each zone during the grid generation stage. Set preset size transition parameters between zones to complete the grid cell division of the geological domain modeling entity.

[0033] Specifically, the prior information of the strata includes the stratum category and confidence level.

[0034] Specifically, the nearest neighbor retrieval and distance-weighted voting strategy described in step S6 includes:

[0035] Calculate the average neighbor distance of the m nearest neighbor points around the query point with the centroid coordinates of the grid cell as the query point, set an initial search radius, and if the number of sampling points within the search radius meets the standard or the search radius reaches the maximum limit radius, otherwise, use a linear expansion formula to adaptively expand the search radius until the number of sampling points within the search radius meets the standard or the search radius reaches the maximum limit radius. Finally, the layer category of the grid cell is obtained through layer consistency determination. The layer consistency determination is to accumulate the weights of the same candidate layers, take the layer with the largest sum of weights as the layer of the grid cell, and map it to the corresponding grid cell.

[0036] Specifically, the degree of stratum change C mentioned in step S43 layer (x) passes through the local neighborhood B of point x. r The occurrence ratio of different layers or the number of layer jumps within (x) is calculated.

[0037] Specifically, the local lithological fluctuations C mentioned in step S43 litho (x) passes through the local neighborhood B of point x. r (x) is obtained by weighted variance or degree of variation of lithological attribute values; specifically, for point x, its neighborhood B is taken. r The surrounding spatial point x within (x) i lithological properties a(x) i The fluctuations are statistically analyzed using distance weighting, and the local lithological fluctuation C is defined. litho The formula for calculating (x) is shown below:

[0038]

[0039] in This is the weighted average of the lithological properties surrounding point x; Let x be the surrounding space point iThe weight of is determined by its distance from point x, and is calculated using the distance decay function;

[0040] The distance attenuation function is shown in the following equation:

[0041] ,

[0042] Where, d i Let x be a point and its surrounding space x i The distance; ε is the weighted stabilizing term; ρ>0 is the distance decay exponent.

[0043] Secondly, the present invention provides a rapid modeling device for shield tunneling areas in complex geological formations, and a rapid modeling method for shield tunneling areas in complex geological formations as described in any one of the first aspects above, comprising the following units:

[0044] The point cloud and attribute acquisition unit is used to acquire point set data and lithological attribute data containing three-dimensional coordinates and stratigraphic labels, and to determine the three-dimensional envelope boundary of the modeling domain based on the range of point set coordinates.

[0045] The geological domain and tunnel construction unit is used to construct a three-dimensional entity of the geological domain according to the three-dimensional envelope boundary, and to construct a tunnel cylinder according to the shield tunnel parameters and spatially align it with the geological domain.

[0046] The geological domain modeling entity construction unit is used to perform Boolean cutting on the three-dimensional entity of the geological domain and the tunnel cylinder to obtain a geological domain modeling entity containing the tunnel chamber space.

[0047] Mesh generation and zoning densification units are used to construct a continuous mesh size field based on the tunnel distance field and stratigraphic complexity of the geological domain modeling entity. Then, based on the continuous mesh size field, adaptive mesh generation and zoning densification control are applied to the geological domain modeling entity to generate mesh units.

[0048] The voxel-level stratigraphic prior field construction unit is used to construct voxel grids in the geological domain modeling entity with a step size not less than the global grid control size. Stratigraphic prior information is obtained by performing stratigraphic clustering statistics on the point set in each voxel grid, and a voxel-level stratigraphic prior field is established to establish the mapping relationship between voxel grid and stratigraphic prior information.

[0049] The grid layer assignment unit is used to determine the voxel cell in which each grid cell in the geological domain modeling entity belongs based on its centroid coordinates as the query point. Then, it retrieves the layer prior information corresponding to the voxel cell based on the voxel-level layer prior field. If the confidence of the layer prior information is greater than or equal to a preset confidence threshold, the layer category of the layer prior information is taken as the layer category of the grid cell; otherwise, the layer category of the grid cell is determined based on the nearest neighbor retrieval and distance-weighted voting strategy.

[0050] The material and solid section assignment unit is used to automatically generate the material and solid section of the geological domain modeling entity according to the layer category and assign values ​​to each layer.

[0051] The model export and analysis unit is used to output a finite element model file containing mesh topology and layered assignment information based on the geological domain modeling entity, and to use the finite element model file for shield tunneling construction analysis.

[0052] Thirdly, the present invention also provides an electronic device including a processor, a memory, a communication interface, and one or more programs stored in the memory and configured to be executed by the processor, the programs including instructions for performing the steps of the method described in any one of the first aspects.

[0053] This invention provides a method, apparatus, and equipment for rapid modeling of shield tunneling areas in complex geological formations. Attached Figure Description

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

[0055] Figure 1 A schematic diagram of a rapid modeling method for shield tunneling areas in complex geological formations provided by an embodiment of the present invention;

[0056] Figure 2 This is a schematic representation of stratigraphic information provided in the embodiments of the present invention;

[0057] Figure 3 This is a schematic diagram of geological grid points in a point set region provided in an embodiment of the present invention;

[0058] Figure 4 This is a schematic diagram of the geological domain modeling entity after adaptive mesh division and zoning refinement provided in an embodiment of the present invention;

[0059] Figure 5 This is a schematic diagram of the geological domain modeling entity for automatic stratigraphic construction and assignment provided in an embodiment of the present invention;

[0060] Figure 6 A schematic diagram of a rapid modeling device for shield tunneling areas in complex geological formations provided in an embodiment of the present invention;

[0061] Figure 7 This is a schematic diagram of a rapid modeling device for shield tunneling areas in complex geological formations, provided as an embodiment of the present invention. Detailed Implementation

[0062] The present invention will be explained in detail through the following embodiments. The purpose of this invention is to protect all technical improvements within its scope. In the description of this invention, it should be understood that the terms "first" and "second" are used for descriptive purposes only and should not be construed as indicating or implying relative importance or implicitly specifying the number of indicated technical features. Therefore, a feature defined as "first" or "second" may explicitly or implicitly include one or more of that feature. In the description of this invention, "a plurality of" means two or more, unless otherwise explicitly specified.

[0063] To make the technical problems, technical solutions and advantages of the present invention clearer, a detailed description will be given below in conjunction with the accompanying drawings and specific embodiments.

[0064] Example 1

[0065] refer to Figure 1 This embodiment provides a rapid modeling method for shield tunneling areas in complex geological formations, including the following steps:

[0066] S1. Obtain point set data and lithological attribute data containing three-dimensional coordinates and stratigraphic labels, and determine the three-dimensional envelope boundary of the modeling domain based on the coordinate range of the point set.

[0067] Before determining the three-dimensional envelope boundary of the modeling domain based on the coordinate range of the point set, the method further includes: performing format consistency verification and outlier removal on the point set;

[0068] refer to Figure 2 First, point set data is read from a 3D point set data file (such as CSV / Excel / Parquet, etc.) containing stratigraphic labels (layer categories). The point set data contains the 3D coordinates (x, y, z) and layer category of each sampling point, which is used to describe the spatial geological distribution of the modeling area. Then, the data is format-validated and cleaned, including field integrity verification, data type matching check, and outlier removal. When the file does not exist, fields are missing, or the format is incorrect, an error message is output to ensure that the input data for subsequent geometric construction is reliable and consistent.

[0069] After successfully reading the point set, the 3D bounding box of the modeling domain is calculated based on the coordinate range of all points to obtain the minimum and maximum coordinate positions, and the spatial boundary conditions of the geological modeling domain are determined accordingly. This bounding box is used to limit the construction range of subsequent geometric entities and also provides a spatial basis for voxel lattice generation, distance field generation, and mesh size field calculation.

[0070] The lithological attribute data can be obtained from the preliminary geological survey report. When importing the point set data file at the initial stage of modeling, the lithological attribute data should be imported together. The lithological attributes can include many rock properties, such as density, uniaxial compressive strength, modulus, etc. The lithological attribute data to be imported can be determined according to actual needs. There is a correspondence between the lithological attribute data and the point set data. The corresponding lithological attribute can be obtained through the three-dimensional coordinates or serial number of the point.

[0071] S2. Construct a three-dimensional entity of the geological domain based on the three-dimensional envelope boundary, and construct a tunnel cylinder based on the shield tunnel parameters and spatially align it with the geological domain; the shield tunnel parameters include the tunnel radius R, centerline position parameters, and axis.

[0072] Constructing a tunnel cylinder based on shield tunnel parameters and spatially aligning it with the geological domain specifically includes: generating a tunnel cylinder with a spatial orientation consistent with the tunnel radius R and centerline position parameters, and achieving spatial positioning of the tunnel configuration in the geological domain through axis alignment and position registration;

[0073] refer to Figure 3 The spatial geological field constructed based on point set data forms the basic distribution characteristics of the modeling domain. First, in the finite element modeling environment, a cuboid 3D geological domain entity is automatically generated according to the 3D bounding box size of the modeling domain to cover the spatial distribution range of all geological points. On this basis, a cylindrical geometric model consistent with the tunnel axis, i.e., the tunnel cylinder, is generated according to the preset shield tunnel radius R and centerline position parameters. The tunnel cylinder is initially constructed along the Z-axis, and its axis is aligned with the Y-axis direction of the global coordinate system through rotation transformation. Spatial translation is also performed according to the coordinates of the tunnel center point to ensure that the tunnel configuration is accurately and consistently positioned in the 3D geological domain.

[0074] S3. Perform Boolean cutting on the three-dimensional entity of the geological domain and the tunnel cylinder to obtain a geological domain modeling entity containing the tunnel chamber space;

[0075] Based on the spatial intersection relationship between the three-dimensional entity of the geological domain and the tunnel cylinder, geometric set operations are performed in the modeling environment to subtract the overlapping entity region of the tunnel cylinder in the geological domain using Boolean subtraction, forming a geological modeling domain entity containing the tunnel chamber space. To maintain the consistency and integrity of the geometric topology, the original entity is automatically cleaned up and new shaped parts are generated after Boolean cutting, which serve as the basic geometric input for subsequent distance field generation, size field construction and mesh discretization.

[0076] S4. Construct a continuous grid size field based on the tunnel distance field and stratigraphic complexity of the geological domain modeling entity, and then generate grid cells by adaptive grid division and zoning densification control of the geological domain modeling entity based on the continuous grid size field.

[0077] refer to Figure 4 First, global mesh control parameters are applied within the modeling domain to determine the baseline size of the overall mesh. Then, a tunnel distance field is constructed with the tunnel wall as the zero level set. Combined with the stratigraphic complexity index formed by the point set layer variation and local lithological fluctuations, a continuous mesh size field is constructed to control the degree of mesh refinement. This allows smaller target mesh sizes to be automatically obtained near the tunnel wall and in areas with drastic stratigraphic changes, while larger target mesh sizes are obtained far from the tunnel or in geologically homogeneous areas.

[0078] Subsequently, based on the change of the mesh size field along the normal direction, the effective band width of the wall refinement area, the inner band and the outer band and the size of the representative target unit are automatically obtained, and they are transformed into edge lines, surfaces and global mesh control parameters in the modeling environment, so as to realize the adaptive determination of the band boundary and mesh size, replacing the traditional fixed ratio or manually specified densification method.

[0079] Regarding the selection of element types, based on features such as the local curvature of the modeling domain, the ratio of feature size to mesh size, a hybrid mesh is automatically assigned, mainly composed of reduced integral hexahedral elements (C3D8R), supplemented by prism elements (C3D6) and tetrahedral elements (C3D4), to ensure the generativeability of the mesh in complex geometric regions and the overall numerical convergence stability.

[0080] The construction of a continuous grid size field based on the tunnel distance field and stratigraphic complexity of the geological domain modeled entity specifically includes the following steps:

[0081] S41. First, apply global mesh control dimensions within the geological domain modeling entity to determine the baseline dimensions of the overall mesh of the geological domain modeling entity;

[0082] The global mesh control size is a unified control quantity that determines the overall mesh reference unit size of the modeling domain. It serves as a global scale reference for mesh generation and is used to provide a basis for mesh size before local refinement. In this embodiment, the global mesh control size is 1m.

[0083] S42. Then, using the tunnel wall as the zero level set, calculate the minimum Euclidean distance from each point in the geological modeling entity to the tunnel wall to construct the tunnel distance field.

[0084] It is understandable that the tunnel distance field is the set of minimum Euclidean distances from each point within the geological modeling entity to the tunnel wall. Then, each point within the geological modeling entity can obtain its minimum Euclidean distance from the tunnel wall through the tunnel distance field.

[0085] The tunnel wall is the wall of the tunnel chamber space;

[0086] S43. Construct a stratigraphic complexity index based on the degree of stratigraphic variation and local lithological fluctuations in the point set;

[0087] For any point x in the point set of a geological modeling entity, statistically analyze the changes in stratigraphic category and lithological properties within its local neighborhood, and define the stratigraphic complexity index C(x) based on these two factors:

[0088]

[0089] Among them, C layer (x) represents the dispersion of the stratigraphic variation within the local neighborhood of point x, used to characterize the degree of stratigraphic variation of the point set; C litho (x) represents the degree of fluctuation of lithological properties (such as density, uniaxial compressive strength, etc.) in the local neighborhood of point x, and is used to characterize the local lithological fluctuation of the point set; λ∈[0,1] is the stratigraphic variation weight coefficient, which is used to control the contribution of the two parts to the complexity. It can be selected according to actual needs. In this embodiment, the value of λ is 0.5.

[0090] Because in actual geological conditions, the distribution of different strata is complex, not simply a neatly stacked layer as ideally envisioned. In reality, strata are often interlayered, and there are abrupt geological changes. Therefore, this embodiment describes the complexity of the strata distribution, i.e., the strata fluctuation, by the degree of stratigraphic variation. Simultaneously, the physical properties of different strata vary, sometimes being similar and sometimes significantly different. Large differences can have a significant impact on engineering and computational analysis; therefore, lithological fluctuations are used to describe this complexity. When modeling and meshing, if the mesh is coarse in key areas with dense stratigraphic interfaces and dramatic lithological changes, it will affect the accuracy and precision of the calculations. Therefore, a stratigraphic complexity index is constructed using the degree of stratigraphic variation of point sets and local lithological fluctuations to regulate the mesh division.

[0091] Specifically, for any point x in the point set of geological modeling entities, define its neighborhood sphere B. r (x) as a local neighborhood of x:

[0092] Taking any point x in the geological modeling entity point set as the center, points whose distance from x does not exceed the local search radius are included in its neighborhood sphere to obtain the local neighborhood. The local search radius is selected according to actual needs.

[0093] In one possible implementation, the degree of stratigraphic variation C layer (x) passes through the local neighborhood B of point x. r The proportion of different layers or the number of layer jumps within (x) are calculated;

[0094] Reference: Altieri L, Cocchi D, Roli G. Efficient Computation of Spatial Entropy Measures. Entropy. 2023; 25(12):1634. https: / / doi.org / 10.3390 / e25121634. This paper provides a technical approach, mainly using a spatial entropy calculation framework based on neighborhood / distance (co-occurrence). This embodiment innovatively applies the idea of ​​characterizing heterogeneity with the dispersion of category distribution to a three-dimensional point set with layer categories, through the local neighborhood B of point x. r The proportion of different strata within (x) or the number of strata jumps measures the degree of local strata variation and serves as part of the stratigraphic complexity index.

[0095] In one possible implementation, local lithological fluctuations C litho (x) passes through the local neighborhood B of point x. r (x) is obtained by weighted variance or degree of variation of lithological attribute values;

[0096] Specifically, for point x, take its neighborhood B. r The surrounding spatial point x within (x) i lithological properties a(x) i (One of the lithological properties can be selected, such as density or uniaxial compressive strength), and its fluctuation is statistically analyzed using distance weighting, defining the local lithological fluctuation C. litho The formula for calculating (x) is shown below:

[0097]

[0098] in This is the weighted average of the lithological properties surrounding point x; Let x be the surrounding space point i The weight of is determined by its distance from point x, and is calculated using the distance decay function;

[0099] The distance attenuation function is shown in the following equation:

[0100] ,

[0101] Where, d i Let x be a point and its surrounding space x i The distance, ε is the weighted stabilization term, taken as a small constant to prevent division by zero; ρ>0 is the distance decay exponent, taken as 2 in this embodiment; local lithological fluctuations C litho The larger the (x) value, the more obvious the fluctuation of lithological properties around point x (the stronger the local lithological fluctuation).

[0102] For example: in the neighborhood B of point x r (x) contains 3 surrounding spatial points: x1, x2, x3;

[0103] The distances between point x and three surrounding spatial points x1, x2, and x3 are d = [1, 2, 3] m, i.e., d1 = 1, d2 = 2, d3 = 3. The lithological properties (uniaxial compressive strength) are: a = [100, 120, 80] MPa, i.e., a1 = 100, a2 = 120, a3 = 80. The parameters are: ε = 0.1, ρ = 2.

[0104] Calculate the weights of the three surrounding spatial points: w1 = 1 / (1 + 0.1) 2 =0.82, w2=1 / (2+0.1) 2 =0.23, w3=1 / (3+0.1) 2 =0.10, therefore ∑w i =w1+w1+w1=1.15;

[0105] ∑w i (a i =(0.826446×100)+(0.226757×120)+(0.104058×80)=118.18;

[0106] =118.18 / 1.15=102.12 MPa

[0107] Substitute into the formula =

[0108] Therefore, in this embodiment, the uniaxial compressive strength property in the neighborhood of point x exhibits a local fluctuation on the order of approximately 10.48 MPa.

[0109] S44. Introduce the tunnel distance field and the geological complexity index into the grid size continuous function to form a continuous grid size field S(x);

[0110] The grid size continuity function is shown in the following equation:

[0111] ;

[0112] Among them, h min h is the minimum allowed mesh size for the modeling domain. max The maximum allowed mesh size for the modeling domain; The encrypted weight of point x is based on the tunnel distance field; This represents the maximum value in the tunnel distance field. This represents the value corresponding to the distance x from the midpoint of the field to the tunnel.

[0113] The encryption weight of point x is based on the stratigraphic complexity index; This represents the maximum value of the stratigraphic complexity index. This represents the minimum value of the stratigraphic complexity index. Let x be the stratigraphic complexity index at point x;

[0114] The continuous grid size field S(x) automatically presents a smaller target size in areas close to the cavern and areas with significant geological changes, while a larger size is presented in locations far from the cavern or in areas with uniform geological formations. The target size value can be obtained based on the point coordinate x through the continuous grid size field S(x).

[0115] The cavern boundary is a crucial control interface for structural analysis, often exhibiting stress concentration or significant deformation gradients. High boundary curvature can lead to geometric approximation errors and numerical stiffness drift if the mesh is too large. In finite element simulations, elements closer to the cavern require higher precision to accurately capture mechanical behavior. Therefore, in areas with significant geological variations (high complexity), mesh refinement is necessary to prevent excessive computational errors, while in areas with homogeneous strata, mesh enlargement can be appropriately used to reduce computational load.

[0116] The process of adaptively dividing and zoning the geological domain modeling entity based on the continuous grid size field to generate grid cells specifically includes the following steps:

[0117] S45. Select M representative points on the tunnel wall. Starting from each representative point, sample along its unit outward normal direction at a preset fixed interval to obtain M target size curves L(t). Then, perform statistical processing on the M target size curves L(t) according to their distance and position to obtain representative size curves.

[0118] To determine the actual variation characteristics of the continuous grid size field along the normal direction, M points need to be randomly selected on the tunnel wall as representative points. Starting from these representative points, samples are taken point by point along the unit outward normal direction of the representative points at a fixed preset sampling distance. The target size values ​​of each sampling point are read from the continuous grid size field to form M target size curves L(t) that show the variation of target size values ​​with normal distance.

[0119] ;

[0120] Where x0 is a representative point selected on the tunnel wall, n(x0) is the unit outward normal of the representative point x0, t is the preset sampling distance along the outward normal; S(·) represents the continuous grid size field;

[0121] It is understandable that the unit outward normal of each point on the tunnel wall is the direction "outward from the wall" where the point is located. For circular tunnels, the direction from the center of the circle to the point is used as the unit outward normal of that point; for tunnels with a centerline, the radial direction within the cross-section is used as the unit outward normal of the point; if it is a mesh model, the patch normal is used as the unit outward normal of the point. This data can be obtained in the software through code, which is existing technology and will not be elaborated here.

[0122] The value of M is a positive integer greater than or equal to 8, generally determined based on the perimeter of the cavern wall and the minimum target grid size in the continuous grid size field, ensuring that the distance between adjacent representative points is no greater than the preset sampling resolution. In this embodiment, M=32. The value of the preset sampling distance t is determined based on the target encryption band width to be covered along the unit outward normal direction and the desired sampling resolution. In this embodiment, t is taken as Δt=h. min Maximum sampling distance t max =10h min That is, t takes the value 0 and h min 2h min ...until 10h min .

[0123] Then, the M target size curves L(t) are statistically processed according to their distance and position to obtain a typical "representative size curve of size change with normal distance", which is used to judge the overall trend of size increase. This curve reflects how the mesh should gradually become thicker from the wall to the outside, and is the basis for automatic three-zone division.

[0124] Statistical processing of M target size curves L(t) according to their distance and position to obtain representative size curves is a common technique for data fusion, which will not be elaborated here.

[0125] S46. Based on the representative dimension curve, determine the boundary points of the three zones according to the trend of dimension change with normal distance, divide the three zones according to the boundary points of the three zones and obtain the distance intervals corresponding to the three zones; the three zones include the wall refinement zone, the inner zone and the outer zone;

[0126] After obtaining representative dimension curves, the three-zone boundary positions can be automatically determined based on the trend of dimension variation with normal distance.

[0127] In this embodiment, the three zones include a wall refinement zone, an inner zone, and an outer zone. Typically, the dimensional changes near the tunnel wall are relatively gradual, but a significant initial increase occurs at a certain distance. The first inflection point of the representative dimensional curve (i.e., the curve transitions from a gradual to a growth phase) is recorded as the boundary t1 of the inner zone of the wall refinement zone, serving as its boundary. Continuing outwards, the dimensional growth gradually stabilizes. The second inflection point of the representative dimensional curve (i.e., the point where the curve transitions from accelerated to slow growth) is recorded as the boundary point t2 between the inner and outer zones, serving as the dividing line between them. This yields the distance intervals corresponding to the three zones, including:

[0128] ;

[0129] Where T1 is the distance interval of the wall refinement zone; T2 is the distance interval of the inner zone; and T3 is the interval of the outer zone. The maximum sampling distance along the unit outward normal direction of the representative point;

[0130] The flat sections of the curve are characterized by slow changes in size, smooth curves, and approximately constant slopes; while the increasing sections of the curve are intervals where the slope changes significantly and the curve becomes steeper. This is common knowledge in the field of mathematics. Detecting the inflection point of the curve is a current technology and will not be elaborated here.

[0131] S47. On the representative size curve, calculate the average or median of the target size values ​​in each distance interval of the three zones to obtain the representative size of each distance interval, and use them as the target mesh size of the three zones respectively.

[0132] It is understandable that the distance intervals corresponding to the three zones divide the representative size curve into three segments, each segment containing multiple sampling points. Each sampling point corresponds to a target size value. That is, the target size value of each sampling point can be read from the continuous grid size field, and all target size values ​​contained in each distance interval can be counted. The average value of all target size values ​​is taken as the representative size of each distance interval, or the median value of all target size values ​​is taken as the representative size of each distance interval.

[0133] Then, the representative dimensions of each distance interval of the three zones are used as the target mesh dimensions of the three zones respectively. That is, the representative dimensions of the wall refinement zone are used as the target mesh dimensions of the wall refinement zone, the representative dimensions of the inner zone are used as the target mesh dimensions of the inner zone, and the representative dimensions of the outer zone are used as the target mesh dimensions of the outer zone, which are used for the specific control of subsequent mesh generation.

[0134] S48. Construct zone boundaries using isosurfaces of the distance field based on the distance intervals corresponding to the three zones, and apply the corresponding target grid size to the distance intervals corresponding to each zone during the grid generation stage, and set preset size transition parameters between zones to complete the grid cell division of the geological domain modeling entity.

[0135] Step S48 further includes: preferentially generating reduced integral hexahedral elements in geometrically smooth regions, using prism elements for transition near boundaries with normal features, and using tetrahedral elements for filling areas with large curvature or dense geometric details.

[0136] Finally, after determining the distance intervals corresponding to the three zones and the target grid size of the three zones, the zone boundaries can be constructed using the isosurface of the distance field based on the distance intervals corresponding to the three zones.

[0137] Constructing zone boundaries based on the distance intervals corresponding to the three zones using the isosurface of the distance field can be achieved through isosurface extraction and triangulation surface reconstruction techniques, which are existing technologies and will not be elaborated here.

[0138] During the mesh generation stage, corresponding target mesh sizes are applied to regions with different distance intervals in the three zones. The preset size transition parameter is used to constrain the smooth transition of target mesh sizes between adjacent zones at the zone boundaries. Near the zone boundaries, the value of a continuous mesh size field S(x) is used as the local target size to avoid abrupt changes in mesh sizes between different zones, thereby ensuring mesh quality.

[0139] In the actual mesh generation process, it is also necessary to select the appropriate element type based on the local geometry: in geometrically flat areas, reduce integral hexahedral elements (C3D8R) are generated first; prism elements (C3D6) are used for transition near the boundary with normal features; and tetrahedral elements (C3D4) are used to fill areas with large curvature or dense geometric details, so as to balance mesh generation capability and overall numerical stability.

[0140] Local geometry (geometrically flat regions, regions with normal characteristics, and areas with dense curvature and geometric details) can be obtained through techniques such as normal calculation, curvature estimation, and local mesh morphology analysis.

[0141] In geometry, the normal vector is the direction of a straight line perpendicular to the tangent plane at a point on a surface. In three-dimensional space, it is defined as the direction of the normal vector of the surface at that point. Normal vector calculation can be achieved through plane equation coefficient extraction, surface partial derivative cross product calculation, or eigenvalue decomposition techniques of the covariance matrix of neighborhood points in point cloud data. It has wide applications in mathematics, computer graphics, and engineering.

[0142] Curvature describes the degree of bending of a local surface and is an important indicator of surface geometry. In point cloud processing, curvature estimation is usually based on principal component analysis (PCA) to decompose the covariance matrix into eigenvalues ​​and eigenvectors, and then to calculate the maximum and minimum curvature.

[0143] Local mesh morphology analysis involves parameters such as mesh density, element type, and shape quality, and is used to evaluate the mesh's adaptability to geometric features. These techniques are widely used in fields such as finite element analysis, computer-aided design (CAD), and computer-aided engineering (CAE), and will not be elaborated here.

[0144] This step first constructs a tunnel distance field with the tunnel wall as the zero level set, and then combines the stratigraphic complexity index formed by the changes in the point set stratigraphic position and local lithological fluctuations to construct a continuous grid size field based on the tunnel distance field and stratigraphic complexity. This results in smaller target unit sizes near the tunnel wall and in areas with drastic stratigraphic changes, while gradually transitioning to larger sizes in areas far from the tunnel or in areas with uniform stratigraphy.

[0145] This invention automatically calculates the wall refinement zone, the width of the inner and outer zones, and the target mesh size corresponding to each zone based on the spatial gradient of the size field along the tunnel normal direction. It then converts these into local mesh control parameters for the edge lines (constructing the zone boundaries using the isosurface of the distance field based on the distance intervals corresponding to the three zones) and the surface (applying the corresponding target mesh size to the distance intervals corresponding to each zone during the mesh generation stage). This enables adaptive definition of the zone range and the degree of densification without the need for manually setting a fixed ratio or absolute size, thus achieving automatic division of the three zones and automatic densification construction of the mesh.

[0146] Regarding the selection of mesh element types, the system automatically selects a hybrid mesh format based on local geometry and different regions, with reduced integral hexahedral elements (C3D8R) as the main type, supplemented by prism elements (C3D6) and tetrahedral elements (C3D4) to ensure the mesh generation capability of complex geometric regions, while maintaining the numerical accuracy of the tunnel near-field region and the stability of the overall problem-solving process.

[0147] Through the aforementioned adaptive zonal densification strategy, this embodiment effectively controls the number of global cells while ensuring the near-field grid resolution of the cavern, achieving a comprehensive balance between grid quality, computational efficiency, and geological representation accuracy, and providing a structurally unified and reliable grid foundation for subsequent spatial stratification assignment.

[0148] S5. In the geological domain modeling entity, construct a voxel grid with a step size not less than the global grid control size. Perform stratigraphic clustering statistics on the point set in each voxel grid to obtain stratigraphic prior information. Establish a voxel-level stratigraphic prior field that maps the voxel grid to stratigraphic prior information.

[0149] After mesh generation is completed, the system enters the voxel construction and layer mapping stage, such as... Figure 5 As shown, a regular voxel grid is constructed by discretizing the modeling domain with a step size no smaller than the global grid size. Stratigraphic clustering is then performed on the point sets within each voxel grid to obtain stratigraphic prior information (including stratigraphic category and confidence level), forming a voxel-level stratigraphic prior field reflecting the stability of stratigraphic distribution. This voxel-level stratigraphic prior field is dominated by the majority class and automatically exhibits null values ​​or low confidence levels in sparsely distributed point sets, thus constructing a coarse-scale preliminary stratigraphic distribution structure in space, forming a priori mapping table of "voxel grid-stratigraphic prior information". This voxel-level stratigraphic prior field is not only used for direct determination of the stratigraphic layer by the centroid of the grid cell, but also serves as the initial constraint and candidate set basis for subsequent local nearest neighbor retrieval and adaptive stratigraphic refinement, thereby improving spatial consistency and robustness of stratigraphic determination under complex stratigraphic conditions.

[0150] Example 1 exists, where the modeling domain in three-dimensional space ranges from 0 to 1000 meters in the x-direction, 0 to 800 meters in the y-direction, and -500 to 0 meters in the z-direction (the negative sign indicates the underground depth).

[0151] Set the global mesh control size to 1 meter, and discretize the modeling domain with a step size no smaller than this size (here, we take 5 meters). Then, 1000 / 5=200 intervals will be divided in the x direction, 800 / 5=160 intervals in the y direction, and 500 / 5=100 intervals in the z direction. In this way, a total of 200×160×100=3,200,000 regular voxel grids will be constructed, and each voxel grid has a fixed position and size in space (5 meters × 5 meters × 5 meters).

[0152] Each sampling point in the point set of each voxel grid contains three-dimensional coordinates (x, y, z) and a layer category. Layer clustering statistics are performed on the point set within each voxel grid to obtain layer prior information, that is, to count the layer category of all points falling within it. For example, for a voxel grid with 10 points, statistics show that 7 points have a layer category of 1, 2 points have a layer category of 2, and 1 point has a layer category of 3.

[0153] Based on the majority class principle, the stratigraphic category of the prior information for this voxel is determined to be stratum 1. Simultaneously, the confidence level of this prior information is calculated based on the distribution of the point set, such as the proportion of points from different stratigraphic layers. In the example above, since points from stratum 1 account for 70%, the confidence level is 0.7, which is relatively high. If the number of points in a voxel is small, or if the distribution of points from different stratigraphic layers is relatively uniform and there is no clear majority class, then the confidence level of the stratigraphic prior information for that voxel will be low, or even null.

[0154] By integrating the stratigraphic prior information from all voxel lattices, a voxel-level stratigraphic prior field is formed, reflecting the stability of stratigraphic distribution. This prior field spatially presents a coarse-scale preliminary stratigraphic distribution structure. In areas with dense point set distribution and relatively stable stratigraphy, the prior field can clearly reflect the stratigraphic distribution; while in areas with sparse point set distribution and drastic stratigraphic changes, the prior field may exhibit null values ​​or low confidence.

[0155] Simultaneously, a priori mapping table of voxel grid-stratum category-confidence is established. The table records the number of each voxel grid (which can be determined by its index in the x, y, z directions) and the corresponding prior information of the stratum (stratum category and confidence). The priori mapping table is shown in Table 1.

[0156] Table 1 Prior Mapping Table

[0157]

[0158] S6. For each grid cell in the geological domain modeling entity, determine its voxel cell by using its centroid coordinates as the query point, and then retrieve the stratigraphic prior information corresponding to the voxel cell based on the voxel-level stratigraphic prior field. If the confidence of the stratigraphic prior information is greater than or equal to the preset confidence threshold, then the stratigraphic category of the stratigraphic prior information is taken as the stratigraphic category of the grid cell.

[0159] When determining the stratigraphic level of a grid cell, the spatial coordinates of the cell's centroid are first used as the query point to determine the voxel cell in which it belongs. Then, the stratigraphic prior information corresponding to that voxel cell is retrieved. If the prior information has sufficient confidence (confidence greater than or equal to a preset confidence threshold), the stratigraphic category in the prior information is directly used as the stratigraphic level determination result for the grid cell. For example, if the centroid of a grid cell falls within voxel cell numbered (1,4,5), and the stratigraphic prior information for that voxel cell has a stratigraphic category of 3 and a confidence level of 0.8, then the stratigraphic category of that grid cell is directly determined to be 3.

[0160] The preset reliability threshold is selected according to actual needs, and is generally selected in the range of [0.6,1). In this embodiment, the preset reliability threshold is 0.6.

[0161] If the confidence level of the prior information of the layer is less than the preset confidence threshold, the layer category of the grid cell is determined based on the nearest neighbor retrieval and distance-weighted voting strategy.

[0162] The nearest neighbor retrieval and distance-weighted voting strategy specifically includes:

[0163] Calculate the average neighbor distance of the m nearest points around the query point with the centroid coordinates of the grid cell as the query point, set the initial search radius, and if the number of sampling points within the search radius meets the standard or the search radius reaches the maximum limit radius, otherwise use the linear expansion formula to adaptively expand the search radius until the number of sampling points within the search radius meets the standard or the search radius reaches the maximum limit radius. Finally, the layer category of the grid cell is obtained through layer consistency determination.

[0164] If the confidence level of the prior information on the layer position of the voxel grid is insufficient (the confidence level is less than the preset confidence threshold), or if the point set within the voxel grid is sparse and there is no majority class (no direct layer position information), a local refinement process needs to be initiated. At the same time, an initial search radius is set, and then neighboring sampling points are collected as candidate points in the voxel grid within the initial search radius around the centroid. The layer position of the candidate points is used as the candidate layer position. The weights of the candidate points are calculated according to the distance decay function (such as the reciprocal or power function) to form a locally weighted layer position information set. If the number of sampling points within the initial search radius reaches the preset minimum number of samples, the layer position consistency judgment is performed on the weighted candidate layers to obtain the layer position result of the grid cell.

[0165] The layer consistency determination is to accumulate the weights of the same candidate layers, take the layer with the largest sum of weights as the layer of the grid cell, and map it to the corresponding grid cell.

[0166] The distance decay function described in this step is the same as the distance decay function described in step S4:

[0167] If the point distribution within the initial search radius is still sparse or insufficient to achieve a stable judgment, that is, if the local point density within the initial search radius is insufficient to satisfy k samples, the search radius is dynamically expanded based on the average point distance of the neighborhood, and the above process is repeated until the number of sampling points within the search radius reaches the preset minimum number of samples or the search radius reaches the preset maximum limit radius. Then, a layer consistency judgment is performed, and finally, a final layer result with spatial continuity and stability is generated within the adaptive neighborhood (search radius) and mapped to the corresponding grid cell.

[0168] The nearest neighbor retrieval and distance-weighted voting strategy can be implemented through the following steps:

[0169] a) First, calculate the local point density. Calculate the average neighbor distance ρ(x) of the m nearest points around each query point x. This reflects the local density of the point. The value of m is selected according to actual needs, and is generally m=10.

[0170] b) Set the initial search radius r = c1 × ρ(x) based on the average neighbor distance ρ(x); where the first expansion factor c1 is a constant, which is set autonomously according to the distance of the model;

[0171] c) If the number of sampling points found within the search radius r is n ≥ k, or the search radius reaches the preset maximum limit radius r max Proceed to step f); otherwise proceed to step d).

[0172] d) Adaptively expand the current search radius using the linear expansion formula, then proceed to step c); the linear expansion formula is: r = r + c2 × ρ(x), where the second expansion factor c2 is a constant, set autonomously based on the model's distance;

[0173] e) Finally, weighted voting is performed, and the weight w is calculated based on the distance to each point within the search radius. i For each possible level, calculate its weighted sum, and finally select the level with the largest weight as the level of the query point.

[0174] The stratum of the query point is returned based on the weighted majority vote result.

[0175] Example 2 exists. Assume the centroid coordinates of the query point are (x0, y0, z0). Within the voxel lattice containing this point, there is no direct stratigraphic information or the confidence level of the prior stratigraphic information is insufficient (the confidence level is less than a preset confidence threshold). Therefore, it is necessary to obtain the corresponding information from neighboring sampling points among several nearby voxels. Assume the three found neighboring sampling points are A, B, and C, and their distance is d. A =1.0 (distance), d B =2.0, d C =3.0;

[0176] The weight of each query point is calculated using a distance decay function. In this embodiment, w is set to 0.01, i.e.: w A =1 / (1+0.01) 2 ≈0.9803, w B =1 / (2+0.01) 2 ≈0.2475, w C =1 / (3+0.01) 2 =0.1104;

[0177] If the layer of point A is "Layer1", the layer of point B is "Layer2", and the layer of point C is "1", then according to the weights, the weight of "Layer1" is 0.9803 + 0.1104 = 1.0907, while the weight of "Layer2" is 0.2475. Finally, the layer of the grid cell is determined to be "Layer1".

[0178] Example 3 exists. Assume that the average neighbor distance ρ(x) of the query point x is 2.0 m, k=8, c1=1.25, the initial search radius r=c1×ρ(x)=2.5 m, and the expansion step of the search radius is a linear expansion c2=0.5.

[0179] The initial search begins with a search radius of r = 2.5 m. If 5 points are found, the number of points is insufficient for k = 8, so the radius needs to be expanded. The new radius is r = 2.5 + 0.5 × 2.5 = 3.75 m. If 7 points are found, the number of points is still insufficient for k = 8, so the search continues. The expanded radius is r = 3.75 + 0.5 × 2.5 = 5.0 m. If 8 points are found, the requirement is met, and weighted voting begins. For the 8 found points, refer to Example 2 to calculate the distance d from each point to the query point. i And calculate the weight w of each query point according to the distance decay function. i A weighted vote is performed on each level, and the level with the most votes is selected as the query point.

[0180] S7. Automatically generate the material and cross-section of the geological domain modeling entity according to the layer category and assign values ​​to each layer;

[0181] A material parameter dictionary is pre-established with layer category as the key. Each material category includes preset density, elastic modulus, and Poisson's ratio parameters, meaning that each layer corresponds to a set of density, elastic modulus, and Poisson's ratio parameters. After the layer determination of the mesh unit is completed, material objects and homogeneous solid sections are generated according to the layer. Based on the final layer mapping result, the layer-unit set correspondence is constructed, and the corresponding solid sections are assigned to each layer unit set in batches. This realizes the automatic and consistent association between material properties and spatial layering structure, effectively avoiding omissions or mismatches in manual operation and improving the efficiency and reliability of layer assignment.

[0182] S8. Output a finite element model file containing mesh topology and layered assignment information based on the geological domain modeling entity, and use the finite element model file for shield tunneling construction analysis.

[0183] Understandably, the geological domain modeling entity is the "mother body" for geometry and mesh generation, while the finite element model element is the "analysis unit" generated from this entity and inheriting its geometric and stratigraphic information. Simply put, the geological domain modeling entity is a 3D geometric model (a virtual soil mass with stratigraphic attributes) used to build the mesh, and the finite element model element is a specific mesh unit (with numbering, nodes, topology, and material) generated from this entity.

[0184] Export finite element model files containing mesh and layer assignments (compatible with solvers such as ABAQUS and ANSYS); generate 3D visualization data (such as VTK / PLY) when needed to display the geological distribution and element layering results around the tunnel. The system also records log information (bounding boxes, mesh statistics, number of elements in each layer, miss rate, etc.) for quality verification and experiment reproduction. The generated 3D geological model accurately presents the geometric features of the tunnel chamber and achieves layered mesh discretization consistent with the geological structure throughout the entire domain, which can be directly used for shield tunneling construction simulation, geological response analysis, and engineering design optimization.

[0185] This embodiment provides a rapid modeling method for shield tunneling areas in complex geological formations. Driven by 3D point set data containing stratigraphic categories, the method constructs a modeling domain through data verification and spatial bounding box calculations, generating geological entities and shield tunnel models. It then identifies and subtracts tunnel space through geometric set relationships to form a geological modeling domain containing caverns. Regarding mesh generation, this invention constructs a continuous mesh size field based on the tunnel distance field and stratigraphic complexity, automatically determining the partition width and representative mesh size of the wall refinement zone, inner zone, and outer zone, achieving adaptive zoning and densification for shield tunneling areas. Simultaneously, it constructs an adaptive stratigraphic determination mechanism for the modeling space. Through a joint decision-making process involving coarse-scale priors, local refinement, and dynamic neighborhood expansion, it achieves stable and continuous stratigraphic automatic assignment to mesh cells. Compared to traditional modeling methods that rely on simplified stratigraphic interfaces or linear partitioning, this invention effectively avoids boundary ambiguity and stratigraphic misjudgment, significantly improving modeling accuracy and consistency under complex geological conditions, and supporting dynamic model updates based on incremental point sets. This method can output a three-dimensional geological model that can be directly used for construction analysis and parameter optimization, providing efficient and reliable technical support for shield tunneling construction design, risk assessment and operation management.

[0186] Example 2

[0187] refer to Figure 6 This embodiment provides a rapid modeling device for shield tunneling areas with complex geological formations, comprising the following units:

[0188] The point cloud and attribute acquisition unit is used to acquire point set data and lithological attribute data containing three-dimensional coordinates and stratigraphic labels, and to determine the three-dimensional envelope boundary of the modeling domain based on the range of point set coordinates.

[0189] The geological domain and tunnel construction unit is used to construct a three-dimensional entity of the geological domain according to the three-dimensional envelope boundary, and to construct a tunnel cylinder according to the shield tunnel parameters and spatially align it with the geological domain.

[0190] The geological domain modeling entity construction unit is used to perform Boolean cutting on the three-dimensional entity of the geological domain and the tunnel cylinder to obtain a geological domain modeling entity containing the tunnel chamber space.

[0191] Mesh generation and zoning densification units are used to construct a continuous mesh size field based on the tunnel distance field and stratigraphic complexity of the geological domain modeling entity. Then, based on the continuous mesh size field, adaptive mesh generation and zoning densification control are applied to the geological domain modeling entity to generate mesh units.

[0192] The voxel-level stratigraphic prior field construction unit is used to construct voxel grids in the geological domain modeling entity with a step size not less than the global grid control size. Stratigraphic prior information is obtained by performing stratigraphic clustering statistics on the point set in each voxel grid, and a voxel-level stratigraphic prior field is established to establish the mapping relationship between voxel grid and stratigraphic prior information.

[0193] The grid layer assignment unit is used to determine the voxel cell in which each grid cell in the geological domain modeling entity belongs based on its centroid coordinates as the query point. Then, it retrieves the layer prior information corresponding to the voxel cell based on the voxel-level layer prior field. If the confidence of the layer prior information is greater than or equal to a preset confidence threshold, the layer category of the layer prior information is taken as the layer category of the grid cell; otherwise, the layer category of the grid cell is determined based on the nearest neighbor retrieval and distance-weighted voting strategy.

[0194] The material and solid section assignment unit is used to automatically generate the material and solid section of the geological domain modeling entity according to the layer category and assign values ​​to each layer.

[0195] The model export and analysis unit is used to output a finite element model file containing mesh topology and layered assignment information based on the geological domain modeling entity, and to use the finite element model file for shield tunneling construction analysis.

[0196] Example 3

[0197] refer to Figure 7 , Figure 7 This is a schematic diagram of the structure of a rapid shield tunneling area modeling device for complex geological formations according to this embodiment. The rapid shield tunneling area modeling device 20 of this embodiment includes a processor 21, a memory 22, and a computer program stored in the memory 22 and executable on the processor 21. When the processor 21 executes the computer program, it implements the steps in the above method embodiments. Alternatively, when the processor 21 executes the computer program, it implements the functions of each module / unit in the above device embodiments.

[0198] For example, the computer program can be divided into one or more modules / units, which are stored in the memory 22 and executed by the processor 21 to complete the present invention. The one or more modules / units can be a series of computer program instruction segments capable of performing specific functions, which describe the execution process of the computer program in the shield tunneling area rapid modeling device 20 for complex geological formations. For example, the computer program can be divided into the modules described in Embodiment 2. The specific functions of each module are described in the working process of the device described in the above embodiments, and will not be repeated here.

[0199] The rapid modeling device 20 for shield tunneling areas with complex geological formations may include, but is not limited to, a processor 21 and a memory 22. Those skilled in the art will understand that the schematic diagram is merely an example of the rapid modeling device 20 for shield tunneling areas with complex geological formations and does not constitute a limitation on the device. It may include more or fewer components than illustrated, or combine certain components, or use different components. For example, the rapid modeling device 20 for shield tunneling areas with complex geological formations may also include input / output devices, network access devices, buses, etc.

[0200] The processor 21 can be a central processing unit (CPU), or other general-purpose processors, digital signal processors (DSPs), application-specific integrated circuits (ASICs), field-programmable gate arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc. The general-purpose processor can be a microprocessor or any conventional processor. The processor 21 is the control center of the rapid modeling equipment 20 for shield tunneling areas facing complex geological formations, connecting all parts of the equipment 20 via various interfaces and lines.

[0201] The memory 22 can be used to store the computer programs and / or modules. The processor 21 implements various functions of the shield tunneling area rapid modeling device 20 for complex strata by running or executing the computer programs and / or modules stored in the memory 22 and calling the data stored in the memory 22. The memory 22 may mainly include a program storage area and a data storage area. The program storage area may store the operating system, at least one application program required for a function (such as sound playback function, image playback function, etc.), etc.; the data storage area may store data created according to the use of the mobile phone (such as audio data, phonebook, etc.). In addition, the memory 22 may include high-speed random access memory, and may also include non-volatile memory, such as hard disk, memory, plug-in hard disk, smart media card (SMC), secure digital (SD) card, flash card, at least one disk storage device, flash memory device, or other volatile solid-state storage device.

[0202] The modules / units integrated in the shield tunneling area rapid modeling equipment 20 for complex geological formations, if implemented as software functional units and sold or used as independent products, can be stored in a computer-readable storage medium. Based on this understanding, all or part of the processes in the above embodiments of the present invention can also be implemented by a computer program instructing related hardware. The computer program can be stored in a computer-readable storage medium, and when executed by the processor 21, it can implement the steps of the various method embodiments described above. The computer program includes computer program code, which can be in the form of source code, object code, executable files, or certain intermediate forms. The computer-readable medium can include: any entity or device capable of carrying the computer program code, recording media, USB flash drives, portable hard drives, magnetic disks, optical disks, computer memory, read-only memory (ROM), random access memory (RAM), electrical carrier signals, telecommunication signals, and software distribution media, etc. It should be noted that the content contained in the computer-readable medium may be appropriately added to or subtracted from the content as required by the legislation and patent practice in the jurisdiction. For example, in some jurisdictions, according to legislation and patent practice, the computer-readable medium may not include electrical carrier signals and telecommunication signals.

[0203] It should be noted that the device embodiments described above are merely illustrative. The units described as separate components may or may not be physically separate, and the components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the modules can be selected to achieve the purpose of this embodiment according to actual needs. Furthermore, in the accompanying drawings of the device embodiments provided by this invention, the connection relationships between modules indicate that they have communication connections, which can be specifically implemented as one or more communication buses or signal lines. Those skilled in the art can understand and implement this without any creative effort.

[0204] This specification is described with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of this specification. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, create a machine for implementing the flowchart illustrations and / or block diagrams. Figure 1 A process, multiple processes, and / or boxes Figure 1 Devices that specify the functions in one or more boxes.

[0205] These computer program instructions may also be stored in a computer-readable storage medium that can direct a computer or other programmable data processing device to function in a particular manner, such that the instructions stored in the computer-readable storage medium produce an article of manufacture including an instruction device, which is implemented in a process Figure 1 One or more processes and / or boxes Figure 1 The function specified in one or more boxes.

[0206] These computer program instructions may also be loaded onto a computer or other programmable data processing equipment to cause a series of operational steps to be performed on the computer or other programmable equipment to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable equipment for implementing the process. Figure 1 One or more processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.

[0207] The parts of this invention not described in detail are prior art. It will be apparent to those skilled in the art that this invention is not limited to the details of the above exemplary embodiments, and that the invention can be implemented in other specific forms without departing from the spirit or essential characteristics of the invention. Therefore, the embodiments should be regarded as exemplary and non-limiting in all respects, and are intended to encompass all changes falling within the meaning and scope of equivalents within this invention.

Claims

1. A rapid modeling method for shield tunneling areas in complex geological formations, characterized in that, Includes the following steps: S1. Obtain point set data and lithological attribute data containing three-dimensional coordinates and stratigraphic labels, and determine the three-dimensional envelope boundary of the modeling domain based on the coordinate range of the point set. S2. Construct a three-dimensional entity of the geological domain based on the three-dimensional envelope boundary, and construct a tunnel cylinder based on the shield tunnel parameters and spatially align it with the geological domain; S3. Perform Boolean cutting on the three-dimensional entity of the geological domain and the tunnel cylinder to obtain a geological domain modeling entity containing the tunnel chamber space; S4. Construct a continuous grid size field based on the tunnel distance field and stratigraphic complexity of the geological domain modeling entity, and then generate grid cells by adaptive grid division and zoning densification control of the geological domain modeling entity based on the continuous grid size field. The construction of a continuous grid size field based on the tunnel distance field and stratigraphic complexity of the geological domain modeled entity specifically includes the following steps: S41. First, apply global mesh control dimensions within the geological domain modeling entity to determine the baseline dimensions of the overall mesh of the geological domain modeling entity; S42. Then, using the tunnel wall as the zero level set, calculate the minimum Euclidean distance from each point in the geological modeling entity to the tunnel wall to construct the tunnel distance field. S43. Construct a stratigraphic complexity index C(x) based on the degree of stratigraphic variation and local lithological fluctuations of the point set, as shown in the following formula: ; Among them, C layer (x) represents the dispersion of the stratigraphic variation within the local neighborhood of point x, used to characterize the degree of stratigraphic variation of the point set; C litho (x) represents the degree of fluctuation of lithological properties in the local neighborhood of point x, and is used to characterize the local lithological fluctuation of the point set; λ is the stratigraphic variation weighting coefficient; S44. Introduce the tunnel distance field and the geological complexity index into the grid size continuous function to form a continuous grid size field S(x); The grid size continuity function is shown in the following equation: ; Among them, h min h is the minimum allowed mesh size for the modeling domain. max The maximum allowed mesh size for the modeling domain; The encrypted weight of point x is based on the tunnel distance field; This represents the maximum value in the tunnel distance field. This represents the value corresponding to the distance x from the midpoint of the field to the tunnel. The encryption weight of point x is based on the stratigraphic complexity index; This represents the maximum value of the stratigraphic complexity index. This represents the minimum value of the stratigraphic complexity index. Let x be the stratigraphic complexity index at point x; S5. In the geological domain modeling entity, construct a voxel grid with a step size not less than the global grid control size. Perform stratigraphic clustering statistics on the point set in each voxel grid to obtain stratigraphic prior information. Establish a voxel-level stratigraphic prior field that maps the voxel grid to stratigraphic prior information. S6. For each grid cell in the geological domain modeling entity, determine its voxel cell by using its centroid coordinates as the query point. Then, retrieve the stratigraphic prior information corresponding to the voxel cell based on the voxel-level stratigraphic prior field. If the confidence level of the stratigraphic prior information is greater than or equal to a preset confidence threshold, then use the stratigraphic category of the stratigraphic prior information as the stratigraphic category of the grid cell; otherwise, determine the stratigraphic category of the grid cell based on the nearest neighbor retrieval and distance-weighted voting strategy. S7. Automatically generate the material and cross-section of the geological domain modeling entity according to the layer category and assign values ​​to each layer; S8. Output a finite element model file containing mesh topology and layered assignment information based on the geological domain model entity, and use the finite element model file for shield tunneling construction analysis.

2. The rapid modeling method for shield tunneling areas in complex geological formations according to claim 1, characterized in that, Step S2 specifically includes: the shield tunnel parameters include the tunnel radius R, centerline position parameters, and axis; the construction of the tunnel cylinder based on the shield tunnel parameters and spatial alignment with the geological domain specifically includes: generating a tunnel cylinder consistent with its spatial orientation based on the tunnel radius R and centerline position parameters, and realizing the spatial positioning of the tunnel configuration in the geological domain through axis alignment and position registration.

3. The rapid modeling method for shield tunneling areas in complex geological formations according to claim 1, characterized in that, The process of adaptively dividing and zoning the geological domain modeling entity based on the continuous grid size field to generate grid cells specifically includes the following steps: S45. Select M representative points on the tunnel wall. Starting from each representative point, sample along its unit outward normal direction at preset fixed intervals to obtain M target size curves L(t). Then, statistically process the M target size curves L(t) according to their distance and position to obtain representative size curves. The target size curves L(t) are shown in the following formula: ; Where x0 is a representative point selected on the tunnel wall, n(x0) is the unit outward normal of the representative point x0, t is the preset sampling distance along the outward normal; S(·) represents the continuous grid size field; S46. Based on the representative dimension curve, determine the boundary points of the three zones according to the trend of dimension change with normal distance, divide the three zones according to the boundary points of the three zones and obtain the distance intervals corresponding to the three zones; the three zones include the wall refinement zone, the inner zone and the outer zone; S47. On the representative size curve, calculate the average or median of the target size values ​​in each distance interval of the three zones to obtain the representative size of each distance interval, and use them as the target mesh size of the three zones respectively. S48. Construct zone boundaries using isosurfaces of the distance field based on the distance intervals corresponding to the three zones. Apply the corresponding target grid size to the distance intervals corresponding to each zone during the grid generation stage. Set preset size transition parameters between zones to complete the grid cell division of the geological domain modeling entity.

4. The rapid modeling method for shield tunneling areas in complex geological formations according to claim 1, characterized in that, The prior information about the strata includes the stratum category and confidence level.

5. The rapid modeling method for shield tunneling areas in complex geological formations according to claim 4, characterized in that, The nearest neighbor retrieval and distance-weighted voting strategy described in step S6 specifically includes: Calculate the average neighbor distance of the m nearest neighbor points around the query point with the centroid coordinates of the grid cell as the query point, set an initial search radius, and if the number of sampling points within the search radius meets the standard or the search radius reaches the maximum limit radius, otherwise, use a linear expansion formula to adaptively expand the search radius until the number of sampling points within the search radius meets the standard or the search radius reaches the maximum limit radius. Finally, the layer category of the grid cell is obtained through layer consistency determination. The layer consistency determination is to accumulate the weights of the same candidate layers, take the layer with the largest sum of weights as the layer of the grid cell, and map it to the corresponding grid cell.

6. The rapid modeling method for shield tunneling areas in complex geological formations according to claim 1, characterized in that, The degree of stratigraphic change C mentioned in step S43 layer (x) is calculated by the proportion of different layers or the number of layer jumps within the local neighborhood Bᵣ(x) of point x.

7. The rapid modeling method for shield tunneling areas in complex geological formations according to claim 1, characterized in that, The local lithological fluctuations C mentioned in step S43 litho (x) is obtained by weighting the variance or variability of lithological attribute values ​​within the local neighborhood Bᵣ(x) of point x; specifically, it includes taking the neighborhood B of point x. r The surrounding spatial point x within (x) i lithological properties a(x) i The fluctuations are statistically analyzed using distance weighting, and the local lithological fluctuation C is defined. litho The formula for calculating (x) is shown below: , in This is the weighted average of the lithological properties surrounding point x; Let x be the surrounding space point i The weight of is determined by its distance from point x, and is calculated using the distance decay function; The distance attenuation function is shown in the following equation: , Where, d i Let x be a point and its surrounding space x i The distance is ε, where ε is the weight stabilization term and ρ is the distance decay exponent.

8. A rapid modeling device for shield tunneling areas in complex geological formations, characterized in that, The rapid modeling method for shield tunneling areas in complex geological formations based on any one of claims 1-7 includes the following components: The point cloud and attribute acquisition unit is used to acquire point set data and lithological attribute data containing three-dimensional coordinates and stratigraphic labels, and to determine the three-dimensional envelope boundary of the modeling domain based on the range of point set coordinates. The geological domain and tunnel construction unit is used to construct a three-dimensional entity of the geological domain according to the three-dimensional envelope boundary, and to construct a tunnel cylinder according to the shield tunnel parameters and spatially align it with the geological domain. The geological domain modeling entity construction unit is used to perform Boolean cutting on the three-dimensional entity of the geological domain and the tunnel cylinder to obtain a geological domain modeling entity containing the tunnel chamber space. Mesh generation and zoning densification units are used to construct a continuous mesh size field based on the tunnel distance field and stratigraphic complexity of the geological domain modeling entity. Then, based on the continuous mesh size field, adaptive mesh generation and zoning densification control are applied to the geological domain modeling entity to generate mesh units. The voxel-level stratigraphic prior field construction unit is used to construct voxel grids in the geological domain modeling entity with a step size not less than the global grid control size. Stratigraphic prior information is obtained by performing stratigraphic clustering statistics on the point set in each voxel grid, and a voxel-level stratigraphic prior field is established to establish the mapping relationship between voxel grid and stratigraphic prior information. The grid layer assignment unit is used to determine the voxel cell in which each grid cell in the geological domain modeling entity is located based on its centroid coordinates as the query point. Then, the layer prior information corresponding to the voxel cell is retrieved based on the voxel-level layer prior field. If the confidence of the layer prior information is greater than or equal to the preset confidence threshold, the layer category of the layer prior information is taken as the layer category of the grid cell. Otherwise, the layer category of the grid cell is determined based on the nearest neighbor retrieval and distance-weighted voting strategy; The material and solid section assignment unit is used to automatically generate the material and solid section of the geological domain modeling entity according to the layer category and assign values ​​to each layer. The model export and analysis unit is used to output a finite element model file containing mesh topology and layered assignment information based on the geological domain modeling entity, and to use the finite element model file for shield tunneling construction analysis.

9. An electronic device, comprising a storage medium, a processor, and a computer program stored on the storage medium and executable on the processor, characterized in that, When the processor executes the computer program, it implements the rapid modeling method for shield tunneling areas oriented towards complex strata as described in any one of claims 1-7.

Citation Information

Patent Citations

  • Risk assessment method based on combination of fuzzy analytic hierarchy process and TOPSIS

    CN120598376A

  • Seismic attribute inversion and ancient karst form modeling fused reservoir identification method

    CN120742423A