Natural resource three-dimensional visualization modeling method based on digital twinning

Abstract

This invention relates to the field of digital simulation technology, specifically to a method for three-dimensional visualization modeling of natural resources based on digital twins. The method includes: constructing a discrete computational map based on a real-world three-dimensional model of the target area and geographic information system (GIS) data; projecting the design base of the building information model (BIM) onto the discrete computational map to determine design boundary nodes; performing constraint propagation calculations on the discrete computational map using the design boundary nodes as propagation starting points to determine the cumulative compliance consumption and cumulative horizontal projection distance of each node in the discrete computational map; determining the topology type of each node based on its cumulative compliance consumption, cumulative horizontal projection distance, original terrain elevation, and preset engineering slope parameters; and geometrically reconstructing the discrete computational map based on the topology type of each node to generate a digital twin model that integrates the BIM and the real-world three-dimensional model. This invention simultaneously achieves automatic compliance avoidance and structural closure generation of the slope model.

Description

Technical Field

[0001] This invention relates to the field of digital simulation technology, specifically to a three-dimensional visualization modeling method for natural resources based on digital twins. Background Technology

[0002] With the deep application of digital twin technology in natural resource management and land spatial planning engineering, the demand for the integration of heterogeneous models in earthwork engineering simulation and compliance review is becoming increasingly urgent. It is necessary to achieve effective combination of Building Information Modeling (BIM) and real-world 3D models that reflect the current surface to support the accurate implementation of core processes such as engineering design and earthwork volume estimation.

[0003] In existing technologies, the fusion of heterogeneous models largely relies on traditional geometric calculation methods, such as Boolean operations or mesh smoothing algorithms based on spatial distance. These methods aim to eliminate intersections and conflicts between models at the geometric level, attempting to construct 3D models that meet engineering requirements and satisfy the basic need for combining terrain and design intent during earthwork scheme simulation. However, these traditional methods fail to fully consider the coupling relationship between surface semantic constraints and engineering physical constraints, resulting in fused 3D models often exhibiting problems such as insufficient compliance, structural incompleteness, and geometric accuracy that does not meet engineering specifications. Consequently, they are unable to reliably support compliance reviews and accurate simulations in actual engineering projects. Summary of the Invention

[0004] To address the technical problem of poor performance of existing heterogeneous model fusion methods in 3D visualization modeling of natural resources, the present invention aims to provide a 3D visualization modeling method for natural resources based on digital twins. The specific technical solution adopted is as follows:

[0005] Firstly, a method for three-dimensional visualization modeling of natural resources based on digital twins is provided. This method includes: for a target area to be simulated in earthwork engineering, constructing a discrete computational map based on the real-world three-dimensional model of the target area and geographic information system data. Nodes in the discrete computational map record the original terrain elevation, and edges in the discrete computational map are configured with traffic resistance weights, which characterize the compliance costs of traversing different types of surface areas; projecting the design base plate of the building information model onto the discrete computational map, determining design boundary nodes and assigning them design reference elevations corresponding to the design base plate, with the design boundary nodes being nodes within the projection area of ​​the design base plate; and performing constraint propagation calculations on the discrete computational map using the design boundary nodes as the propagation starting point. The cumulative compliance consumption and cumulative horizontal projection distance of each node in the discrete computational graph are determined. The cumulative compliance consumption is obtained by accumulating the traffic resistance weights of each side on the propagation path, and the cumulative horizontal projection distance is obtained by accumulating the projection length of the propagation path on the horizontal plane. For each node in the discrete computational graph, the theoretical design elevation of the node is determined according to the node's cumulative horizontal projection distance and preset engineering slope parameters. The topology type of the node is determined according to the first comparison result between the node's theoretical design elevation and the original terrain elevation, and the second comparison result between the node's cumulative compliance consumption and the preset compliance consumption threshold. The discrete computational graph is geometrically reconstructed according to the topology type of each node to generate a digital twin model that integrates the building information model and the real-world 3D model.

[0006] In one possible design, constructing a discrete computational atlas includes: discretizing the triangular mesh of the real-world 3D model, extracting and removing the vertices of all triangular faces to form a node set; defining connecting edges between nodes in the node set according to the edge connection relationships of the triangular faces to form an edge set; and recording the elevation value of each node in 3D space as the original terrain elevation.

[0007] In one possible design, determining the traffic resistance weight corresponding to each edge in the discrete computational map includes: for each edge in the discrete computational map, determining the land cover category corresponding to the edge in the geographic information system data; obtaining the resistance coefficient corresponding to the land cover category according to a preset resistance coefficient mapping relationship reflecting the engineering constraint strength of different land cover categories; and determining the traffic resistance weight of the edge based on the three-dimensional geometric length of the edge and the resistance coefficient.

[0008] In one possible design, constraint propagation calculations are performed on the discrete computation graph, including: initializing the cumulative compliance consumption and cumulative horizontal projection distance of the design boundary nodes to zero, and initializing the cumulative compliance consumption and cumulative horizontal projection distance of other nodes in the discrete computation graph (excluding the design boundary nodes) to preset maximum values, while recording the self-identity of each design boundary node as its design affiliation; establishing a minimum priority queue based on the cumulative compliance consumption of nodes, and adding all the design boundary nodes to the minimum priority queue; if the minimum priority queue is not empty, retrieving the node with the minimum cumulative compliance consumption as the current node; for each adjacent node of the current node, updating the cumulative compliance consumption and cumulative horizontal projection distance of the adjacent node based on the traffic resistance weight and horizontal projection distance between the current node and the adjacent node, and if the update is successful, recording the design affiliation of the current node as the design affiliation of the adjacent node, and adding or reordering the adjacent node in the minimum priority queue based on the updated cumulative compliance consumption of the adjacent node.

[0009] In one possible design, updating the cumulative compliance cost and cumulative horizontal projection distance of the adjacent node includes: determining the attempted compliance cost of the adjacent node based on the cumulative compliance cost of the current node and the passage resistance weight of the edge connecting the current node and the adjacent node; determining the attempted horizontal projection distance of the adjacent node based on the cumulative horizontal projection distance of the current node and the two-dimensional Euclidean distance between the current node and the adjacent node on the horizontal projection plane; and updating the cumulative compliance cost of the adjacent node based on the attempted compliance cost and the cumulative horizontal projection distance of the adjacent node based on the attempted horizontal projection distance if the attempted compliance cost of the adjacent node is less than the current cumulative compliance cost of the adjacent node.

[0010] In one possible design, the theoretical design elevation of the node is determined based on the cumulative horizontal projection distance of the node and preset engineering slope parameters, including: obtaining the design elevation associated with the node, wherein the design elevation associated with the node is the design reference elevation of the design boundary node corresponding to the design affiliation of the node record; and determining the theoretical design elevation of the node based on the design elevation associated with the node, the cumulative horizontal projection distance of the node, and the preset engineering slope parameters.

[0011] In one possible design, determining the topology type of the node includes: determining the cut / fill type of the node based on the design elevation associated with the node and the original terrain elevation of the node; determining the node as a landing node if the cut / fill type is fill and the first comparison result is that the theoretical design elevation is less than or equal to the original terrain elevation; determining the node as the landing node if the cut / fill type is cut and the theoretical design elevation is greater than or equal to the original terrain elevation, and determining the topology type of the landing node as a first topology type; for non-landing nodes, determining the topology type of the node as a second topology type if the second comparison result is that the cumulative compliance consumption is less than or equal to the compliance consumption threshold; and determining the topology type of the node as a third topology type if the second comparison result is that the cumulative compliance consumption is greater than the compliance consumption threshold.

[0012] In one possible design, the discrete computational graph is geometrically reconstructed according to the topology type of each node, including: for nodes of the second topology type, updating the spatial elevation of the node according to the theoretical design elevation of the node; for nodes of the third topology type, performing a vertex splitting operation to generate vertices representing the top of the designed slope and vertices representing the foundation of the retaining structure, adjusting the vertex indices of the triangular facets connected to the original node to reconstruct the connection relationship, and constructing new triangular facets between the newly generated vertices of adjacent third topology type nodes to form a solidified retaining structure surface; for nodes of the first topology type, calculating the normal vectors of all triangular facets associated with the node, and performing a weighted average of the normal vectors to generate and assign a smoothed vertex normal vector to the node.

[0013] In one possible design, adjusting the vertex indices of the triangular facets connected to the original node to reconstruct the connection relationship includes: for each triangular facet connected to the original node, obtaining the topology type of the other vertices in the triangular facet besides the original node; if the topology type of the other vertices contains the second topology type, modifying the vertex index of the triangular facet that points to the original node to point to the vertex representing the top of the designed slope; if the topology type of the other vertices does not contain the second topology type, modifying the vertex index of the triangular facet that points to the original node to point to the vertex representing the retaining structure foundation.

[0014] In one possible design, a pre-defined resistance coefficient mapping relationship reflecting the engineering constraint strength of different land cover categories is configured as follows: for the area category identified as prohibited from engineering encroachment in the geographic information system data, the corresponding resistance coefficient is a first value; for the area category identified as permitted for earthwork engineering, the corresponding resistance coefficient is a second value, and the first value is greater than the second value.

[0015] The present invention has the following beneficial effects:

[0016] In the digital twin-based 3D visualization modeling method for natural resources provided by this invention, by integrating the real-world model, GIS semantics, and BIM design into a discrete computational map, and utilizing a dual-channel constraint propagation mechanism, the generated slope model can automatically avoid intrusive areas such as roads and waterways based on preset traffic resistance weights, naturally conforming to land use boundary requirements and eliminating the need for extensive manual interaction and adjustments required by traditional methods. Simultaneously, the theoretical design elevation is determined based on the cumulative horizontal projection distance and preset slope ratio, ensuring that the geometry of the earthwork slope strictly adheres to engineering specifications and guaranteeing the accuracy of earthwork volume calculations. Crucially, the introduction of priority-based topology type adjudication logic automatically triggers the physical generation of retaining structures when the slope cannot naturally land (such as touching the boundary line or reaching the slope limit), thereby completely eliminating geometric defects such as model suspension and cracks generated by traditional algorithms in complex terrain, outputting a high-fidelity digital twin model with a closed structure that can be directly used for analysis and simulation. Attached Figure Description

[0017] To more clearly illustrate the technical solutions and advantages in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are 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.

[0018] Figure 1 This is a flowchart illustrating a method for three-dimensional visualization modeling of natural resources based on digital twins, provided in one embodiment of the present invention.

[0019] Figure 2 This is a schematic diagram of the structure of a three-dimensional visualization modeling system for natural resources based on digital twins, provided as an embodiment of the present invention. Detailed Implementation

[0020] To further illustrate the technical means and effects adopted by the present invention to achieve its intended purpose, the following, in conjunction with the accompanying drawings and preferred embodiments, details the specific implementation, structure, features, and effects of a digital twin-based three-dimensional visualization modeling method for natural resources proposed according to the present invention. In the following description, different "one embodiment" or "another embodiment" do not necessarily refer to the same embodiment. Furthermore, specific features, structures, or characteristics in one or more embodiments can be combined in any suitable form.

[0021] In embodiments of the present invention, the terms "exemplary" or "for example" are used to indicate that something is an example, illustration, or description. Any embodiment or design described as "exemplary" or "for example" in embodiments of the present invention should not be construed as being more preferred or advantageous than other embodiments or designs. Specifically, the use of the terms "exemplary" or "for example" is intended to present the relevant concepts in a specific manner.

[0022] In the description of this invention, unless otherwise stated, " / " means "or". For example, A / B can mean A or B. The term "and / or" in this document is merely a description of the relationship between related objects, indicating that three relationships can exist. For example, A and / or B can represent: A alone, A and B simultaneously, and B alone. Furthermore, "at least one" and "more than one" refer to two or more. The terms "first," "second," etc., do not limit the quantity or order of execution, and "first," "second," etc., do not necessarily imply differences.

[0023] Unless otherwise defined, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains.

[0024] The following description, in conjunction with the accompanying drawings, details a specific scheme for a digital twin-based three-dimensional visualization modeling method for natural resources provided by this invention.

[0025] Please see Figure 1 The diagram illustrates a flowchart of a method for three-dimensional visualization modeling of natural resources based on digital twins, according to an embodiment of the present invention, including the following steps S101-S106.

[0026] S101. For the target area to be simulated in earthwork engineering, a discrete computation map is constructed based on the real-world 3D model of the target area and geographic information system data.

[0027] In the discrete computational graph, the nodes record the original terrain elevation, and the edges are configured with traffic resistance weights. The traffic resistance weights are used to characterize the compliance costs of traversing different types of surface areas.

[0028] As one possible approach, firstly, for the target area where the earthwork engineering simulation is to be conducted, acquire the collected real-world 3D model data within the target area. This real-world 3D model data is typically a triangular mesh model (such as OSGB or OBJ format) generated based on oblique photogrammetry technology, consisting of continuous triangular facets. Simultaneously, acquire Geographic Information System (GIS) data for the same area. This GIS data contains multiple vectorized layers used to characterize the types and spatial distribution of land cover, such as roads, water systems, vegetation, construction land, and ecological protection red lines.

[0029] Furthermore, the triangular mesh of the real-world 3D model is discretized, and the vertices of all triangular faces are extracted and deduplicated to form a node set.

[0030] In some embodiments, all triangular faces in the real-world 3D model are traversed, and the coordinates of all vertices of each triangular face are extracted. Subsequently, spatial deduplication is performed on all extracted vertices, that is, vertices with identical 3D coordinates are merged to form a set of nodes in the discrete computational graph.

[0031] Furthermore, the connections between nodes are established by inheriting the edge information of the original triangular facets. For each triangular facet, connecting edges between nodes are defined in the node set based on the edge connections of the triangular facets, thus forming an edge set. In this way, the continuous triangular mesh is transformed into a graph structure composed of discrete nodes and connecting edges, that is, the topological skeleton of the discrete computational graph.

[0032] After forming the node set, each node is assigned its geometric attributes. Optionally, the Z-axis coordinate value of each node in its original realistic 3D model world coordinate system is extracted and recorded as the node's original terrain elevation (denoted as ). This attribute objectively reflects the actual ground elevation at the location of the node. In a preferred embodiment, if the real-world 3D model contains non-ground features such as vegetation, point cloud filtering algorithms such as Cloth Simulation Filter (CSF) can be used to preprocess the original data before extracting the original terrain elevation to remove vegetation noise and ensure that the original terrain elevation represents the true ground elevation.

[0033] Subsequently, semantic attributes, namely traffic resistance weights (denoted as ), are assigned to the edges in the discrete computational map based on geographic information system data. For each edge in the discrete computational map, perform the following operations: determine the spatial location of the midpoint of the edge; query the geographic information system data to obtain the corresponding land cover category (e.g., "paved road", "natural water system", "bare soil", etc.). Based on the preset drag coefficient mapping relationship reflecting the engineering constraint strength of different land cover categories, assign corresponding drag coefficients (denoted as ) to different land cover categories. The mapping relationship is configured as follows: for rigidly constrained areas that are prohibited from encroachment in engineering projects (such as roads, waterways, and red-line protection zones), the allocated size is the first value (e.g., ...). The resistance coefficient is assigned as follows: for areas where earthwork is permitted (such as grasslands or wastelands), a second resistance coefficient (e.g., 1.0) is assigned. Then, based on the three-dimensional geometric length of the edge and the resistance coefficient, the passage resistance weight of the edge is determined.

[0034] In some embodiments, the formula for calculating the passage resistance weight of an edge is as follows:

[0035]

[0036] In the formula, For nodes in the discrete computation graph With nodes The passage resistance weight corresponding to the connecting edge. For nodes in the discrete computation graph With nodes The three-dimensional geometric length of the connecting edge. For nodes in the discrete computation graph With nodes The resistance coefficient corresponding to the connecting edge.

[0037] Through the above process, the semantic constraints (compliance requirements) in the geographic information system data are quantified into weight values ​​on the discrete computational graph edges, thereby embedding the engineering compliance cost into the topological relationship of the graph structure.

[0038] In some embodiments, global engineering parameters, such as the standard slope ratio, can also be initialized when constructing the discrete computational map. (defined as the ratio of horizontal distance to vertical elevation difference), and the maximum slope limit. The maximum slope limit can be based on the maximum horizontal slope distance allowed by the engineering design. Derivation, for example, through formulas Calculate, where, This is the baseline drag coefficient, which is the second-value drag coefficient corresponding to landforms such as grassland and wasteland. The safety redundancy factor can be empirically set to 1.2. The maximum slope limit is used to prevent the unlimited spread of the slope range and defines the maximum total compliance cost that a natural slope can be extended under unobstructed and ideal terrain conditions.

[0039] S102. Project the design base plate of the building information model onto the discrete calculation map, determine the design boundary nodes, and assign the design boundary nodes the design reference elevation corresponding to the design base plate.

[0040] Among them, the design boundary nodes are the nodes within the projection area of ​​the design base plate.

[0041] As one possible approach, firstly, the data of the design slab in the Building Information Model (BIM) is acquired and coordinate system one is established. Then, the closed two-dimensional contour data of the design slab in the proposed project's BIM is extracted; this contour data accurately represents the planar extent of the design slab. Simultaneously, the design elevation value corresponding to this design slab is read; this value is the pre-set earthwork slope reference height for the project.

[0042] In some embodiments, since BIM design data is usually defined in an independent local coordinate system, while discrete calculation maps are constructed based on a world coordinate system, it is necessary to first use a coordinate transformation algorithm (such as the seven-parameter transformation method) to uniformly transform the two-dimensional outline and design elevation values ​​of the BIM design base plate to the world coordinate system corresponding to the discrete calculation map, so as to ensure the spatial consistency of the projection operation.

[0043] Secondly, using the XY horizontal plane of the discrete calculation map as the projection reference plane, the closed two-dimensional contour of the BIM design base plate after coordinate system one is vertically projected along the Z-axis direction, forming a projection contour area within the plane range of the discrete calculation map that is consistent with the plane range of the design base plate. This projection contour area is the initial mapping range of the engineering design intent in geographic space.

[0044] Next, a geometric determination algorithm for "points within polygons" (such as the ray method or rotation number method) is selected, and all nodes in the discrete computation graph are traversed. Each node's XY plane coordinates are determined to be within the above-mentioned projected contour area or on the contour boundary. All nodes that meet the determination conditions are selected to form a set of design boundary nodes. The nodes in this set are the core source points that directly carry the engineering design intent and are the starting basis for subsequent constraint propagation calculations.

[0045] Finally, for each design boundary node in the design boundary node set, based on its corresponding position in the BIM design base plate, the previously read design elevation value of the design base plate is assigned one-to-one to that node as the design reference elevation (denoted as ). The design benchmark elevation is used to characterize the starting elevation of earthwork slope protection at the corresponding node and is the core benchmark for subsequent theoretical design elevation calculations. For other nodes in the discrete calculation map that are not included in the design boundary node set, their design benchmark elevations are initialized to null values. During subsequent constraint propagation calculations, the corresponding design benchmark basis will be dynamically obtained through design attribution tracing to ensure the integrity of the design elevation association logic for all nodes in the domain.

[0046] S103. Using the design boundary node as the propagation starting point, perform constraint propagation calculations on the discrete computation graph to determine the cumulative compliance consumption and cumulative horizontal projection distance of each node in the discrete computation graph.

[0047] Among them, the cumulative compliance consumption is obtained by accumulating the passage resistance weights on each side of the propagation path, and the cumulative horizontal projection distance is obtained by accumulating the projection length of the propagation path on the horizontal plane.

[0048] As one possible implementation, state initialization is performed first: for each design boundary node Accumulated compliance costs With cumulative horizontal projection distance Initialize to zero, that is , This indicates that no compliance costs or horizontal extensions have been incurred at this design origin. Simultaneously, at each design boundary node... Record a design attribution This is used to trace the source of the design elevation later and can be initialized as a unique identifier for the node itself. For nodes in the discrete computational graph other than the design boundary nodes... Accumulated compliance costs With cumulative horizontal projection distance Initialize it to a maximum value (e.g., the maximum floating-point number) to indicate that it has not yet been reached by any design source point propagation, i.e. , Its design belongs to the logo Initialize to empty.

[0049] Furthermore, after initialization, a system is established based on the cumulative compliance consumption of nodes. Minimum priority queue as the sorting criterion and all design boundary nodes Insert into this minimum priority queue .

[0050] Subsequently, the iterative propagation and update process is initiated in the least priority queue. If the value is not empty, repeat the following steps S1031-S1034.

[0051] S1031, From the lowest priority queue Extract (i.e., pop up) accumulated compliance consumption The smallest node is selected as the current node. .

[0052] S1032, Traversal and Current Node All adjacent nodes directly connected by an edge .

[0053] S1033, For each adjacent node The "trial value" is calculated in two dimensions (cumulative compliance consumption and cumulative horizontal projection distance).

[0054] Among them, according to the current node The cumulative compliance costs, and the connection to the current node. and adjacent nodes Determine the travel resistance weight of the edges to identify adjacent nodes. The formula for calculating the compliance costs of the attempt is as follows:

[0055]

[0056] In the formula, To connect the current node and adjacent nodes The weight of the passage resistance on the edge, For the current node The cumulative compliance costs, Adjacent nodes The compliance consumption of the attempt represents the cost if the current node is passed. The path to the adjacent node The total accumulated compliance costs.

[0057] Based on the current node The cumulative horizontal projection distance, and the current node With adjacent nodes Determine adjacent nodes using two-dimensional Euclidean distance on the horizontal projection plane. The formula for calculating the horizontal projection distance is as follows:

[0058]

[0059] In the formula, For the current node and adjacent nodes The formula for calculating the two-dimensional Euclidean distance on a horizontal plane is as follows: , For the current node The cumulative horizontal projection distance, Adjacent nodes The attempted horizontal projection distance represents the cumulative unfolded length along the same path on the horizontal plane.

[0060] S1034, Status Update.

[0061] In adjacent nodes The attempt to comply with the consumption smaller than adjacent nodes Current cumulative compliance costs In this case, it indicates that a path has been found from the design boundary node to the node. The optimal path with lower compliance costs. At this point, execute the state update, based on the compliance costs of attempting the update. Update adjacent nodes Cumulative compliance costs (and soon) The value is updated to (value), and based on the attempted horizontal projection distance Update adjacent nodes Cumulative horizontal projection distance (and soon) The value is updated to (Value). Simultaneously, the current node... Design ownership Record as adjacent nodes Design ownership This ensures that it can be traced back to the correct design boundary node elevation.

[0062] Finally, based on adjacent nodes Updated cumulative compliance costs If node Not in the lowest priority queue In the middle, it is in the lowest priority queue. Add adjacent nodes If node Already in the lowest priority queue In the middle, adjust the nodes. In the lowest priority queue Maintain the position in the minimum priority queue. Cumulative compliance consumption of China-Israel nodes The order of sorting from smallest to largest.

[0063] The above iterative process continues until the minimum priority queue is reached. Empty. When the algorithm terminates, for all nodes in the discrete computational graph reachable from the design boundary nodes... All achieved a definite and optimal cumulative compliance cost. With cumulative horizontal projection distance The value was recorded, along with its design ownership identifier. These calculation results constitute a comprehensive propagation field containing both compliance and geometric information, providing precise numerical basis for subsequent engineering decisions and geometric generation.

[0064] S104. For each node in the discrete calculation map, determine the theoretical design elevation of the node based on the cumulative horizontal projection distance of the node and the preset engineering slope parameters.

[0065] As one possible implementation, iterate through all nodes in the discrete computation graph that have been affected by the propagated computation (i.e., accumulated compliance consumption). (Nodes that are not infinite). For each such node... According to the node Record design ownership This is traced back to the source of its propagation path, namely the corresponding design boundary node. The design baseline elevation assigned during the initialization phase is read from this design boundary node and denoted as... This elevation value is the node. Related design elevation.

[0066] Furthermore, based on the nodes Related design elevation and nodes Original topographic elevation Determine the node The type of cut and fill at the location. Among them, if Then determine the node The corresponding project type is fill, indicating that the design slope needs to be achieved by raising the ground level; if Then determine the node. The corresponding engineering type is cut, which indicates that the design slope needs to be achieved by lowering the ground level; clarifying the cut-fill type provides a basis for determining the direction of elevation change in subsequent formulas.

[0067] Next, the preset engineering slope parameters are the standard slope ratio. It is defined as the ratio of the horizontal distance to the vertical height difference required by the engineering design, and is a non-zero positive number (e.g., 1:1.5 corresponds to...). This parameter is used to limit the slope inclination to ensure that the generated slope surface complies with civil engineering specifications. This parameter can be customized according to the specific project's safety level, soil type, and other requirements. After configuration, it is stored as a global calculation parameter for all nodes to share and call.

[0068] Furthermore, by combining the cut / fill type, the design elevation associated with the node, the cumulative horizontal projection distance of the node, and the preset engineering slope parameters, the theoretical design elevation of the node is determined. The calculation formula is as follows:

[0069]

[0070] In the formula, For nodes The theoretical design elevation, that is, the target elevation that conforms to the slope protection specifications. For nodes Associated design elevation, For nodes The cumulative horizontal projection distance, determined in the constraint propagation calculation, characterizes the horizontal unfolding length originating from the design boundary node. The preset standard slope ratio can be taken as an empirical value of 1.5, which is a non-zero positive number with no denominator and zero risk. In the formula... The selection is determined by the node. The type of fill and cut is determined. If it is a fill type, a minus sign is used, which means that the theoretical design elevation decreases as the cumulative horizontal projection distance increases, which is consistent with the law that the fill slope extends from the reference elevation to the ground. If it is a cut type, a plus sign is used, which means that the theoretical design elevation increases as the cumulative horizontal projection distance increases, which is consistent with the logic that the cut slope is excavated from the reference elevation to the ground.

[0071] Based on this, the theoretical design elevation of all nodes in the discrete computational graph that have been affected by the propagation computation can be obtained.

[0072] S105. Based on the first comparison result between the node's theoretical design elevation and the original terrain elevation, and the second comparison result between the node's cumulative compliance consumption and the preset compliance consumption threshold, determine the node's topology type.

[0073] As one possible implementation, the core attributes of each node are first obtained, such as: theoretical design elevation. Original topographic elevation Accumulated compliance costs And a preset global compliance consumption threshold. Compliance consumption threshold The maximum allowable compliance cost budget for the project is defined, which is usually calculated based on the maximum allowable slope distance and the benchmark resistance coefficient. Refer to step S101 above, which will not be repeated here.

[0074] Furthermore, when determining the topology type of each node, for each node... The following three levels of priority logic are executed sequentially, and the corresponding topology type is assigned a flag once the condition is met. And terminate the subsequent determination of that node.

[0075] The first priority is to determine "landing" and mark it as the first topology type (zero fill-dig line node).

[0076] The first priority reflects the priority of the natural termination of engineering physics. The criterion is: whether the designed slope surface geometrically contacts or cuts into the actual ground. This is determined by comparing the theoretical design elevation. Compared with the original terrain elevation accomplish.

[0077] If the node is of the fill type, then in In the case of determining the node "Landing".

[0078] If the node is of the cut type, then in In the case of determining the node "Landing".

[0079] As long as the above "implementation" conditions are met, regardless of the accumulated compliance costs... Does it exceed the compliance consumption threshold? Immediately set the node Topology type This is identified as the first topology type. These nodes form the geometric boundary line (i.e., the toe line or excavation line) between the designed slope and the natural terrain, marking the physical boundary of the earthwork.

[0080] Second priority: Determine "compliant extension" and mark it as the second topology type (slope extension node).

[0081] The second priority is activated only if the node is determined to be "not yet logged in". The determination condition is: node Whether the extended path is still within the compliance budget's allowance. Specifically, examine its cumulative compliance costs. Does it meet the requirements? If this condition is met, it means that the optimal path from the design boundary to this point does not touch an insurmountable high-resistance area (such as the road right-of-way) and does not exceed the maximum slope range set by the project. In this case, the node... Topology type This is identified as the second topology type. These nodes represent the main area where the design slope is compliant and naturally extends.

[0082] Third priority: Determine "forced truncation" and mark it as the third topology type (retaining structure node).

[0083] The third priority is the default exit, which handles nodes that neither meet the "landing" conditions nor exceed the compliance budget. That is, when a node is "not landed" and its accumulated compliance costs... satisfy At that time, the node Topology type This is identified as the third topology type. This situation typically occurs when the slope path encounters rigid constraints such as roads or water systems, or when the design limit is reached due to excessive horizontal distance. In engineering practice, the natural slope must be terminated at this point, and a retaining structure must be constructed.

[0084] Understandably, through the above-mentioned order of decision based on strict priority, each node in the discrete computation graph is assigned to one of the first, second, or third topological types without omission or conflict.

[0085] S106. Based on the topology type of each node, the discrete computation graph is geometrically reconstructed to generate a digital twin model that integrates the building information model and the real-world 3D model.

[0086] As one possible implementation, for nodes with a topology type of the second topology type, the spatial elevation of the node is updated based on the node's theoretical design elevation.

[0087] In some embodiments, the second topology type node (slope extension node) represents the main area of ​​the designed slope compliant extension, and the goal of reconstruction is to adjust its spatial position to a theoretical elevation that strictly conforms to the design slope. Read each second topology type node. (Slope extension node) The determined cut / fill type (fill or cut) and the calculated theoretical design elevation. In this process, the theoretical design elevation for fill operations decreases as the cumulative horizontal projection distance increases, while the theoretical design elevation for cut operations increases as the cumulative horizontal projection distance increases, ensuring that the elevation adjustment conforms to the engineering slope logic. The spatial elevation of the nodes is further updated, keeping the XY plane coordinates of each slope extension node unchanged (consistent with the original spatial position in the discrete calculation map), only replacing the Z-axis coordinate of the nodes with the theoretical design elevation. This operation forces the original undulating terrain grid to be stretched to meet a preset standard slope ratio. The regular geometric surface allows the slope to form a coherent engineering structure with the BIM-designed base plate.

[0088] In some embodiments, slope geometric continuity is verified by traversing all nodes of the second topology type (slope extension nodes), retrieving triangular patches formed by adjacent nodes, verifying the vertex connection relationship and normal vector direction of the patches, and eliminating overlapping patches or reverse normal vector patches caused by elevation adjustment to ensure that the triangular patches on the slope surface are free of cracks and gaps, forming a continuous and closed slope geometry that conforms to the physical form of the slope in actual engineering.

[0089] For nodes of the third topology type, a vertex splitting operation is performed to generate vertices representing the top of the designed slope and vertices representing the foundation of the retaining structure. The vertex indices of the triangular facets connected to the original nodes are adjusted to reconstruct the connection relationship. New triangular facets are constructed between the new vertices generated by adjacent third topology type nodes to form a solidified retaining structure surface.

[0090] In some embodiments, for each third topology type node (For the retaining structure node) A vertex splitting operation is performed, generating a completely overlapping copy node at the XY plane coordinate position of the original node, forming a vertical node pair of "original node - copy node". The Z-axis coordinate of the original node is updated to its theoretical design elevation. (i.e., the top elevation of the retaining structure), the Z-axis coordinate of the replica node remains the original terrain elevation of that node. (i.e., the elevation of the base of the retaining structure wall), the vertical difference between the two is equal to the elevation separation. This ensures that the height of the retaining structure is precisely matched with the actual needs of the project.

[0091] Furthermore, all associated triangular faces of each retaining structure node are retrieved and categorized into inner and outer faces based on their orientation. Inner faces are those facing the design boundary node and connected to the slope extension node, while outer faces are those facing away from the design boundary node and connected to the undisturbed original terrain node. Vertex index redirection is performed on both types of faces: in the vertex index of the inner face, the entries originally pointing to the retaining structure node remain unchanged (still pointing to the original node, stretched to the top of the wall along with the original node's elevation); in the vertex index of the outer face, the entries originally pointing to the retaining structure node are modified to point to the copy node (maintaining the connection with the original terrain). This operation achieves the geometric separation and connection of the inner and outer sides of the retaining structure.

[0092] Optionally, the indices of other vertices can be determined by the topology type of the other vertices of the associated triangular facet. If a second topology type exists for the other vertices, the vertex index in the triangular facet that points to the original node is modified to point to the vertex representing the top of the designed slope (i.e., still pointing to the original node, stretched to the top of the wall along with the elevation of the original node). If a second topology type does not exist for the other vertices, the vertex index in the triangular facet that points to the original node is modified to point to the vertex representing the foundation of the retaining structure (i.e., maintaining the connection with the original terrain).

[0093] Furthermore, based on the distribution order of the third topology type nodes (retaining structure nodes) in the discrete computational graph, adjacent retaining structure node pairs are extracted sequentially. , ) and their corresponding replica node pairs ( , ), construct two complementary triangular patches between these four nodes (such as and This forms a closed vertical wall surface. The stitching surfaces of all adjacent node pairs are connected in sequence to form a continuous and solid geometric model of the retaining structure. This model accurately represents the gravity retaining wall or support structure in the project in the digital twin, completely solving the problem of model holes caused by the suspended design slope.

[0094] For a node with a topology type of 1, calculate the normal vectors of all triangles associated with the node, and perform a weighted average of the normal vectors to generate and assign a smoothed vertex normal vector to the node.

[0095] In some embodiments, for nodes of the first topology type (zero-cut / fill line nodes), the core of the reconstruction is to eliminate the geometric hard edges at the boundary between the designed slope and the original terrain, achieving a natural transition between the two. First, for each node of the first topology type (zero-cut / fill line node), all associated triangular faces within its first-order neighborhood are retrieved. The first-order neighborhood is the region formed by nodes directly connected to the zero-cut / fill line node via edges. Associated faces include faces partially located on the designed slope and faces partially located in the original terrain; the difference in normal vectors among these faces is the core reason for the geometric hard edges. Second, for each associated triangular face, its face normal vector is calculated using the coordinates of its three vertices. The face normal vector is used to characterize the orientation of the facet in 3D space. During calculation, it is ensured that the normal vector direction is consistent (all pointing outwards from the model) to provide a consistent vector basis for subsequent weighted averaging. Furthermore, the area of ​​each associated triangular facet is used... Using the weights as values, calculate the weighted average of all face normal vectors to obtain the average normal vector of the zero-fill / cut line node. The calculation logic involves balancing the influence of the design slope and the original terrain's normal vectors through patch area weighting, making the transition more closely resemble natural landform features. The introduction of area weighting aims to ensure that larger patches have a stronger impact on the node normal vectors, guaranteeing a smoothing effect consistent with the actual terrain shape. Finally, the calculated average normal vector... The value is assigned to the corresponding zero-fill / cut line node, replacing the node's original normal vector. During the lighting and rendering stage of the digital twin model, this average normal vector will make the light and shadow transition in the zero-fill / cut line area natural, completely eliminating the geometric hard edges caused by elevation differences, achieving seamless integration of the designed slope and the original terrain, and meeting the dual requirements of ecological restoration and engineering aesthetics.

[0096] After completing the targeted reconstruction of three types of topological nodes, the discrete computational map is fully geometrically integrated. The 3D model of the BIM design base plate is aligned with the reconstructed slopes, retaining structures, and natural terrain by coordinate alignment and topological stitching to ensure that the BIM design intent is fully coupled with the real terrain features. Global mesh optimization is performed on the integrated model to simplify redundant triangular faces and retain key geometric details, thereby improving the rendering efficiency and interactive smoothness of the model. Finally, a digital twin model is generated that integrates BIM design information, real terrain data, engineering compliance constraints, and physical structural features. This model can accurately map the earthwork engineering status of the target area and provide a high-fidelity 3D visualization carrier for subsequent engineering simulation, compliance review, and operation and maintenance management.

[0097] Understandably, in the digital twin-based three-dimensional visualization modeling method for natural resources provided in this invention, by integrating the real-world model, GIS semantics, and BIM design into a discrete computational map, and utilizing a dual-channel constraint propagation mechanism, the generated slope model can automatically avoid intrusive areas such as roads and waterways based on preset traffic resistance weights, naturally conforming to land use boundary requirements and eliminating the need for extensive manual interaction and adjustments required by traditional methods. Simultaneously, the theoretical design elevation is determined based on the cumulative horizontal projection distance and preset slope ratio, ensuring that the geometric shape of the earthwork slope strictly adheres to engineering specifications and guaranteeing the accuracy of earthwork volume calculation. Crucially, the introduction of priority-based topology type adjudication logic automatically triggers the physical generation of retaining structures when the slope cannot naturally fall (such as touching the boundary or reaching the slope limit), thereby completely eliminating geometric defects such as model suspension and cracks generated by traditional algorithms in complex terrain, outputting a high-fidelity digital twin model with a closed structure that can be directly used for analysis and simulation.

[0098] Please see Figure 2 The diagram illustrates a structural schematic of a three-dimensional visualization modeling system for natural resources based on digital twins, according to an embodiment of the present invention. Figure 2 As shown, the digital twin-based three-dimensional visualization modeling system for natural resources 20 includes a data processing unit 21, a constraint propagation calculation unit 22, a topology type determination unit 23, and a model generation unit 24.

[0099] The data processing unit 21 is used to construct a discrete calculation map for the target area to be simulated in earthwork engineering, based on the real-world 3D model of the target area and geographic information system data, and to project the design base plate of the building information model onto the discrete calculation map, determine the design boundary nodes and assign the design boundary nodes with the design reference elevation corresponding to the design base plate; its specific implementation can be referred to the description of steps S101 and S102 in the above embodiment of the present invention, and will not be repeated here.

[0100] The constraint propagation calculation unit 22 is used to perform constraint propagation calculation on the discrete calculation graph with the design boundary node as the propagation starting point, and to determine the cumulative compliance consumption and cumulative horizontal projection distance of each node in the discrete calculation graph. The cumulative compliance consumption is obtained based on the cumulative passage resistance weight of each side on the propagation path, and the cumulative horizontal projection distance is obtained based on the cumulative projection length of the propagation path on the horizontal plane. Its specific implementation can be referred to the description of step S103 in the above embodiment of the present invention, and will not be repeated here.

[0101] The topology type determination unit 23 is used to determine the theoretical design elevation of each node in the discrete computational map based on the cumulative horizontal projection distance of the node and the preset engineering slope parameters; and to determine the topology type of the node based on the first comparison result between the theoretical design elevation of the node and the original terrain elevation, and the second comparison result between the cumulative compliance consumption of the node and the preset compliance consumption threshold. Its specific implementation can be referred to the description of steps S104 and S105 in the above embodiment of the present invention, and will not be repeated here.

[0102] The model generation unit 24 is used to perform geometric reconstruction of the discrete computational graph according to the topology type of each node, and generate a digital twin model that integrates the building information model and the real-world 3D model; its specific implementation can be referred to the description of step S106 in the above embodiment of the present invention, and will not be repeated here.

[0103] It should be noted that the order of the above embodiments of the present invention is merely for descriptive purposes and does not represent the superiority or inferiority of the embodiments. The processes depicted in the accompanying drawings do not necessarily require a specific or sequential order to achieve the desired result. In some embodiments, multitasking and parallel processing are also possible or may be advantageous.

[0104] The various embodiments in this specification are described in a progressive manner. The same or similar parts between the various embodiments can be referred to each other. Each embodiment focuses on describing the differences from other embodiments.

Claims

1. A natural resource three-dimensional visualization modeling method based on digital twinning, characterized in that, The method includes: For the target area to be simulated for earthwork engineering, a discrete computation map is constructed based on the real-world 3D model of the target area and geographic information system data. The nodes in the discrete computation map record the original terrain elevation, and the edges in the discrete computation map are configured with traffic resistance weights. The traffic resistance weights are used to characterize the compliance costs of traversing different types of surface areas. The design base plate of the building information model is projected onto the discrete computation map, the design boundary nodes are determined, and the design boundary nodes are assigned design reference elevations corresponding to the design base plate. The design boundary nodes are nodes within the projection area of ​​the design base plate. Using the design boundary node as the propagation starting point, constraint propagation calculation is performed on the discrete computation graph to determine the cumulative compliance consumption and cumulative horizontal projection distance of each node in the discrete computation graph. For each node in the discrete computational graph, the theoretical design elevation of the node is determined based on the cumulative horizontal projection distance of the node and the preset engineering slope parameters. Based on the first comparison result between the theoretical design elevation and the original terrain elevation of the node, and the second comparison result between the cumulative compliance consumption of the node and the preset compliance consumption threshold, the topology type of the node is determined. The discrete computational graph is geometrically reconstructed based on the topology type of each node to generate a digital twin model that integrates the building information model and the real-world 3D model.

2. The natural resource three-dimensional visualization modeling method based on digital twinning according to claim 1, characterized in that, Constructing a discrete computational graph includes: The triangular mesh of the real-world 3D model is discretized, and all vertices of the triangular facets are extracted and deduplicated to form a node set. Based on the edge connection relationship of the triangular facets, define the connecting edges between nodes in the node set to form an edge set; The elevation value of each node in three-dimensional space is recorded as the original terrain elevation.

3. The natural resource three-dimensional visualization modeling method based on digital twinning according to claim 1, wherein, Determining the traffic resistance weight corresponding to each edge in the discrete computational graph includes: For each edge in the discrete computational map, determine the land cover category corresponding to the edge in the geographic information system data; Based on a preset drag coefficient mapping relationship that reflects the engineering constraint strength of different land cover categories, the drag coefficient corresponding to the land cover category is obtained; The passage resistance weight of the edge is determined based on the three-dimensional geometric length of the edge and the resistance coefficient.

4. The natural resource three-dimensional visualization modeling method based on digital twinning according to claim 1, wherein, Performing constraint propagation calculations on the discrete computational graph includes: The cumulative compliance consumption and cumulative horizontal projection distance of the design boundary node are initialized to zero, and the cumulative compliance consumption and cumulative horizontal projection distance of other nodes in the discrete computation graph other than the design boundary node are initialized to a preset maximum value. At the same time, each design boundary node is recorded with its own identifier as the design affiliation. Establish a minimum priority queue based on the cumulative compliance consumption of nodes, and add all the design boundary nodes to the minimum priority queue; If the minimum priority queue is not empty, the node with the lowest cumulative compliance consumption is taken from the minimum priority queue as the current node. For each neighboring node of the current node, the cumulative compliance consumption and cumulative horizontal projection distance of the neighboring node are updated based on the traffic resistance weight and horizontal projection distance between the current node and the neighboring node. If the update is successful, the design ownership of the current node is recorded as the design ownership of the neighboring node. Based on the updated cumulative compliance consumption of the neighboring node, the neighboring node is added to the minimum priority queue or reordered.

5. The natural resource three-dimensional visualization modeling method based on digital twinning according to claim 4, characterized in that, Updating the cumulative compliance consumption and cumulative horizontal projection distance of the adjacent nodes includes: The attempt compliance cost of the adjacent node is determined based on the cumulative compliance cost of the current node and the passage resistance weight of the edge connecting the current node and the adjacent node. The attempted horizontal projection distance of the adjacent node is determined based on the cumulative horizontal projection distance of the current node and the two-dimensional Euclidean distance between the current node and the adjacent node on the horizontal projection plane. If the attempted compliance consumption of the adjacent node is less than the current cumulative compliance consumption of the adjacent node, the cumulative compliance consumption of the adjacent node is updated based on the attempted compliance consumption, and the cumulative horizontal projection distance of the adjacent node is updated based on the attempted horizontal projection distance.

6. The natural resource three-dimensional visualization modeling method based on digital twinning according to claim 4, characterized in that, Based on the cumulative horizontal projection distance of the node and the preset engineering slope parameters, the theoretical design elevation of the node is determined, including: Obtain the design elevation associated with the node, which is the design reference elevation of the design boundary node corresponding to the design ownership of the node record; The theoretical design elevation of the node is determined based on the design elevation associated with the node, the cumulative horizontal projection distance of the node, and the preset engineering slope parameters.

7. The method for three-dimensional visualization modeling of natural resources based on digital twins according to claim 6, characterized in that, Determining the topology type of the node includes: The cut / fill type of the node is determined based on the design elevation associated with the node and the original terrain elevation of the node. If the cut / fill type is fill and the first comparison result is that the theoretical design elevation is less than or equal to the original terrain elevation, the node is determined to be a landing node; if the cut / fill type is cut and the theoretical design elevation is greater than or equal to the original terrain elevation, the node is determined to be the landing node, and the topology type of the landing node is determined to be the first topology type. For non-landing nodes, if the second comparison result is that the cumulative compliance consumption is less than or equal to the compliance consumption threshold, the topology type of the node is determined to be the second topology type. If the second comparison result indicates that the cumulative compliance consumption is greater than the compliance consumption threshold, the topology type of the node is determined to be the third topology type.

8. The natural resource three-dimensional visualization modeling method based on digital twinning according to claim 7, characterized in that, The discrete computation graph is geometrically reconstructed based on the topology type of each node, including: For nodes with a topology type of 2, update the spatial elevation of the node according to the theoretical design elevation of the node; For nodes of the third topology type, a vertex splitting operation is performed to generate vertices representing the top of the designed slope and vertices representing the foundation of the retaining structure. The vertex indices of the triangular facets connected to the original nodes are adjusted to reconstruct the connection relationship. New triangular facets are constructed between the new vertices generated by adjacent third topology type nodes to form a solidified retaining structure surface. For a node with a topology type of the first topology type, calculate the normal vectors of all triangles associated with the node, and perform a weighted average of the normal vectors to generate and assign a smoothed vertex normal vector to the node.

9. The natural resource three-dimensional visualization modeling method based on digital twinning according to claim 8, characterized in that, Adjust the vertex indices of the triangles connected to the original node to reconstruct the connectivity, including: For each triangular facet connected to the original node, obtain the topological type of the other vertices in the triangular facet besides the original node; If the other vertices have the second topology type, modify the vertex index in the triangle patch that points to the original node to point to the vertex representing the top of the designed slope. If the second topology type does not exist in the topology type of the other vertices, the vertex index in the vertex index of the triangle that points to the original node is modified to point to the vertex representing the retaining structure foundation.

10. The natural resource three-dimensional visualization modeling method based on digital twinning according to claim 3, wherein, A preset resistance coefficient mapping relationship reflecting the engineering constraint strength of different land cover categories is configured as follows: for the area category marked as prohibited from engineering encroachment in the geographic information system data, the corresponding resistance coefficient is a first value; for the area category marked as permitted for earthwork engineering, the corresponding resistance coefficient is a second value, and the first value is greater than the second value.

Images