A dynamic simulation method for vine cubic slope repair effect

Abstract

This invention relates to the field of 3D modeling technology, specifically to a dynamic simulation method for slope restoration using a multi-dimensional model. The method includes: reconstructing a slope mesh based on remote sensing imagery and point clouds, according to abrupt slope changes and consistent slope aspect; extracting the 3D coordinates and contact layout of anchor points; analyzing stress and density to mark disturbance levels; selecting continuous deformation nodes to form deformation paths; and drawing tension offset trajectories and intersection sequences. This invention integrates the grayscale differences between remote sensing imagery and point cloud elevations, extracts topographic polylines based on abrupt slope changes and consistent slope aspect for mesh reconstruction, matches topographic features with engineering layout, extracts the 3D coordinates of anchor points and groups them according to their placement direction to generate a spatial topology layout, identifies disturbance units based on principal tensile stress and node density to construct a multi-dimensional disturbance level system, selects continuous deformation paths, draws tension offset trajectories and indexes intersection sections, and achieves visualization of multi-segment restoration evolution.

Description

Technical Field

[0001] This invention relates to the field of 3D modeling technology, and in particular to a dynamic simulation method for slope repair effects using a granular material. Background Technology

[0002] The field of 3D modeling technology encompasses various techniques for digitally representing and reconstructing real-world objects or scenes. Its core content involves generating 3D virtual models with spatial depth, geometric structure, and visual realism through data acquisition, spatial reconstruction, graphics rendering, and dynamic updates. It widely utilizes input information such as point cloud data, image measurement data, and terrain scan data, combined with graphical geometric algorithms for modeling and simulation processing. 3D modeling is widely applied in various fields such as engineering construction, environmental monitoring, virtual reality, disaster assessment, and geographic information systems. It enables the digital restoration, visualization analysis, and dynamic demonstration of real-world environmental characteristics, exhibiting high technical systematization and application complexity.

[0003] The dynamic simulation method for slope restoration using the Tengcube structure refers to constructing a dynamic evolution model of the restoration effect under different temporal and physical effects during slope support and restoration. This mainly involves extracting original slope structure data through remote sensing image interpretation and slope cross-section measurement, spatially reconstructing the data by combining the geometric parameters and layout of Tengcube components, and using a 3D modeling engine to set different natural environmental parameters such as rainfall, wind erosion, and gravity to simulate the temporal response changes. Furthermore, a grid-based method based on terrain data is used to segment and define the overall deformation trend of the slope, and dynamic evolution rules are established by setting fixed boundary conditions and material mechanics parameters. Finally, the simulation and reconstruction of the geometric evolution, deformation response path, and structural stability trend during the Tengcube slope restoration process are completed.

[0004] Existing technologies rely on single terrain scan data and fixed boundary conditions for mesh generation, resulting in insufficient matching between slope structure zoning and terrain features. This makes it difficult to reflect the impact of complex slope aspect changes on engineering layout. Traditional component mechanics analysis uses static contact area calculation and homogenized material parameter settings, which cannot capture the local stress concentration phenomenon caused by changes in the component layout direction and slope aspect angle. Deformation trend prediction is based on the elevation change analysis of discrete nodes, lacking the coordinated observation of horizontal displacement trajectory and tension direction offset. This results in lag and fragmentation in deformation path identification, and the dynamic simulation results fail to present the continuous evolution law of structural response in the spatiotemporal dimension. Summary of the Invention

[0005] To address the technical problems existing in the prior art, this invention provides a dynamic simulation method for slope restoration using the Tenglifang method. The technical solution is as follows:

[0006] A dynamic simulation method for slope restoration using a granular material includes the following steps:

[0007] S1: Obtain remote sensing images and point cloud data of the area where the vine cube structure is laid, sort the grid nodes according to the elevation gradient, extract the terrain polyline to reconstruct the slope grid points, divide the area with consistent slope aspect, and generate a slope structure partition unit matrix.

[0008] S2: Based on the slope structure partition unit matrix, extract the three-dimensional coordinate data of the anchor points of the vine cube structure, calculate the slope contact area and overlap length, group and mark the components according to the angle between the component layout direction and the slope direction vector, and generate a component slope coupling contact layout map.

[0009] S3: Call the component slope coupling contact layout map, obtain the principal tensile stress direction and contact normal stress distribution of the node in the contact area of ​​the component, calculate the angle between the component direction vector and the principal stress direction of the node, determine the consistency distribution, count the contact rate and distribution frequency in the unit grid, mark the disturbance level of the area where the angle deviates from the set limit, and output the slope response interference level classification layer.

[0010] S4: Call the disturbance level marker area in the slope response disturbance level classification layer, obtain the elevation change rate of the slope grid nodes and the horizontal displacement trajectory within a continuous time period, extract the nodes according to the consistency of deformation direction, and form a slope structure deformation trend path group.

[0011] As a further embodiment of the present invention, the slope structure partitioning unit matrix includes reconstructed grid boundaries, slope aspect consistency identifiers, and unit topological relationships; the component slope coupling contact layout map specifically includes a three-dimensional coordinate set of anchor points, a contact area parameter matrix, and a component grouping coding table; the slope response interference level classification layer includes a principal stress direction vector field, a contact rate density distribution cloud map, and a disturbance unit boundary coordinate set; and the slope structure deformation trend path group specifically includes an elevation change rate gradient map, a displacement direction sequence matrix, and a deformation node topological network.

[0012] As a further aspect of the present invention, the step of obtaining the slope structure partitioning unit matrix is ​​as follows:

[0013] S101: Obtain the gray value matrix of remote sensing image and the three-dimensional point cloud elevation dataset. Compare the image gray value and point cloud elevation value at the corresponding position of each grid node pixel by pixel, calculate the absolute value of the difference between the two values, establish a heat map of the distribution of grid node difference values, filter the node positions where the difference value exceeds the terrain change threshold, and generate a terrain difference coefficient matrix.

[0014] S102: Based on the terrain difference coefficient matrix, the slope change rate of adjacent nodes is calculated along the longitude direction. The sliding window method is used to detect abrupt regions where the slope change rate exceeds a set number of consecutive nodes. The coordinate sequence of abrupt points is recorded. Adjacent abrupt points are connected to form a broken line trajectory. The curvature continuity of the broken line is optimized by cubic spline interpolation to obtain the terrain broken line feature set.

[0015] S103: Call the terrain polyline feature set as the grid division boundary, perform topological reconstruction on the original grid, calculate the elevation gradient vector of all nodes in the new grid cell, count the deviation value of the gradient vector azimuth angle in the cell, retain the grid cells with deviation values ​​less than the set angle threshold, merge adjacent cells to form a continuous region, and output the slope structure partition cell matrix.

[0016] As a further aspect of the present invention, the step of obtaining the component slope coupling contact layout diagram is as follows:

[0017] S201: Call the slope structure partition unit matrix, extract the three-dimensional coordinate dataset of the anchor points of the vine cube structure, spatially match the coordinates of each anchor point with the corresponding grid node, calculate the elevation surface fitting plane of all nodes in the grid unit where the anchor point is located, count the number of contact points projected onto the fitting plane of the anchor point, and generate a contact area vector set.

[0018] S202: Based on the contact area vector set, extract the spacing data of adjacent anchor points along the component boundary, measure the length difference of the overlapping area at the end of the component, record the overlapping length variation curve according to the component number order, establish a component connection sequence index table, and output the connection sequence.

[0019] S203: Based on the connection sequence, calculate the cosine value of the spatial angle between the layout direction vector of each component and the corresponding slope aspect vector, classify the components with a cosine value greater than a set threshold into the same group, perform continuous area numbering mapping on the components in the same group, integrate the contact area and connection sequence data, and generate a component slope coupling contact layout map.

[0020] As a further aspect of the present invention, the step of obtaining the slope response disturbance level classification layer is as follows:

[0021] S301: Call the component slope coupling contact layout map, extract the principal tensile stress direction vector and normal stress scalar value of the node in the contact area, calculate the cosine value of the spatial angle between the principal tensile stress direction of each node and the corresponding component direction vector, count the frequency of cosine value distribution, and generate a set of force direction consistency coefficients.

[0022] S302: Based on the set of force consistency coefficients, the node density values ​​are statistically analyzed according to the grid cells, the number of contact nodes and the frequency of normal stress distribution within a unit area are calculated, the product of the number of contact nodes and the frequency of normal stress is used as the contact rate index, the grid contact rate distribution matrix is ​​established, and the contact rate gradient field is output.

[0023] S303: Based on the contact rate gradient field, calculate the absolute value of the difference in force direction consistency coefficient between adjacent grid cells, filter out regions where the difference exceeds the stress disturbance threshold, perform topological merging on continuous exceeding regions, assign values ​​according to disturbance intensity, and generate a slope response disturbance level classification layer.

[0024] As a further aspect of the present invention, the step of obtaining the slope structure deformation trend path group is as follows:

[0025] S401: Call the slope response disturbance level classification layer, extract the time series data of elevation change rate of grid nodes within the disturbance boundary, calculate the ratio of node elevation increment in adjacent time periods, establish node elevation change trend vector, and generate elevation dynamic trend set.

[0026] S402: Based on the elevation dynamic trend set, obtain the horizontal displacement trajectory coordinates of the same node in a continuous time period, calculate the cosine value of the angle between the displacement direction vectors of adjacent time periods, count the number of consecutive time periods where the cosine value is greater than the displacement same direction threshold, and output the horizontal displacement same direction sequence.

[0027] S403: Based on the horizontal displacement homogeneity sequence, select nodes whose elevation increment ratio changes synchronously with the cosine value of the horizontal displacement direction, perform spatial clustering analysis on the nodes that meet the synchronization condition, connect adjacent nodes to form a continuous path, and generate a slope structure deformation trend path group.

[0028] As a further aspect of the present invention, the method further includes:

[0029] S5: Based on the slope structure deformation trend path group, record the change data of structural tension direction of continuous deformation nodes in the time series, draw tension offset trajectory diagram according to the tension direction vector transformation order in the path, mark the node coordinates and component connection order of each offset point in the path, mark the tension turning node and the intersection of the path, organize the path sequence of the whole area and output the dynamic sequence diagram of multi-segment structural response of the repair area.

[0030] The dynamic sequence diagram of the multi-segment structural response in the repair area includes a tension direction vector trajectory diagram, a set of spatial coordinates of turning nodes, and a topological index of path intersection segments.

[0031] As a further aspect of the present invention, the step of obtaining the dynamic sequence diagram of the multi-segment structural response in the repair region is as follows:

[0032] S501: Call the slope structure deformation trend path group, extract the time series tension direction vector of continuous deformation node, calculate the angle change between vectors in adjacent time periods in time order, draw a line graph of tension direction offset trajectory in each path, and generate tension offset trajectory set.

[0033] S502: Based on the tension offset trajectory set, associate the node coordinates of the path where each offset point is located with the component connection sequence number, establish a mapping relationship table between the spatial position of the offset point and the component sequence, mark the node position where the tension direction changes abruptly by more than a set angle, and output the tension turning node index.

[0034] S503: Combining the tension turning node index, the density of turning nodes in the path intersection area is statistically analyzed, the intersection segment level is divided according to the density threshold, the spatiotemporal evolution data of the entire path is integrated, and a dynamic sequence diagram of the multi-segment structural response of the repair area is generated.

[0035] The beneficial effects of the technical solutions provided in the embodiments of the present invention include at least the following:

[0036] This invention integrates remote sensing image grayscale values ​​and 3D point cloud elevation difference analysis to establish a terrain polyline extraction mechanism based on slope abrupt change point sorting and slope aspect consistency constraints. This enables grid reconstruction of slope structure zoning, improving the matching accuracy between terrain features and engineering layout. By combining a dynamic grouping strategy of anchor point 3D coordinates and component layout direction vector angle, a contact layout map of spatial topological association is constructed, enhancing the pertinence of structural mechanical response analysis. A disturbance unit identification method based on principal tensile stress direction and node density dual criteria is adopted to form a multi-dimensional disturbance level classification system, improving the spatiotemporal resolution of slope deformation early warning. Based on the analysis of the elevation change rate and displacement direction symmetry of continuous deformation nodes, a dynamic association model of tension offset trajectory and path intersection segment is established, realizing a multi-segmented visualization of the repair effect evolution. Attached Figure Description

[0037] Figure 1 This is a flowchart of the method of the present invention;

[0038] Figure 2 This is a flowchart illustrating the process of obtaining the slope structure partitioning unit matrix of the present invention.

[0039] Figure 3 This is a flowchart illustrating the process of obtaining the inclined plane coupling contact layout diagram of the components of this invention.

[0040] Figure 4 This is a flowchart illustrating the process of obtaining the slope response interference level classification layer of the present invention.

[0041] Figure 5 This is a flowchart illustrating the process of obtaining the deformation trend path group of the slope structure according to the present invention.

[0042] Figure 6 This is a flowchart of the process for obtaining the dynamic sequence diagram of the multi-segment structural response in the repair region of the present invention. Detailed Implementation

[0043] The technical solution of the present invention will now be described with reference to the accompanying drawings.

[0044] In embodiments of the present invention, words such as "exemplarily," "for example," etc., are used to indicate that something is an example, illustration, or description. Any embodiment or design described as "exemplary" in the present invention should not be construed as being more preferred or advantageous than other embodiments or designs. Specifically, the use of the word "exemplary" is intended to present the concept in a concrete manner. Furthermore, in embodiments of the present invention, the meaning expressed by "and / or" can be both, or either one.

[0045] In the embodiments of this invention, the terms "image" and "picture" may sometimes be used interchangeably. It should be noted that, without emphasizing the distinction between them, they convey the same meaning. Similarly, the terms "of," "corresponding (relevant)," and "corresponding" may sometimes be used interchangeably. It should be noted that, without emphasizing the distinction between them, they convey the same meaning.

[0046] In this embodiment of the invention, sometimes a subscript such as W1 may be written in a non-subscript form such as W1. When the difference is not emphasized, the meaning they express is the same.

[0047] 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.

[0048] Please see Figure 1 This invention provides a technical solution: a dynamic simulation method for slope restoration using a granular material, comprising the following steps:

[0049] S1: Obtain remote sensing images and terrain scanning data of the area where the vine cube structure is laid. By comparing the difference between the gray value of the image and the elevation of the three-dimensional point cloud, sort the surface grid nodes according to the location of the slope change point. Extract the terrain polyline based on the elevation gradient and slope aspect consistency between adjacent nodes. Use the polyline as the grid boundary to reconstruct the grid points, divide the area into slope aspect consistent grid groups, and generate a slope structure partition unit matrix.

[0050] S2: Based on the slope structure partition unit matrix, extract the three-dimensional coordinate data of the anchor points of the vine cube structure, match the anchor points with the grid nodes and calculate the slope contact area, measure the overlap length on the component boundary and record the connection sequence according to the component number, group the components according to the angle between the component layout direction and the slope aspect vector (the slope aspect vector is a standard topographic term that represents the horizontal projection direction of a unit normal of a slope, often used for slope orientation calculation, and has a standard definition in GIS and 3D modeling), number and mark the areas where the component layout direction is consistent with the slope aspect, and generate a component slope coupling contact layout map.

[0051] S3: Call the component slope coupling contact layout map to obtain the principal tensile stress direction and contact normal stress distribution of the nodes in the contact area of ​​the component (the principal tensile stress direction and contact normal stress are commonly used parameters in geotechnical mechanics, which can be obtained by decomposing the eigenvalues ​​of the stress tensor). Analyze the uniformity of force distribution at each node based on the angle between the component direction vector and the principal stress direction of the node. Calculate the contact rate and distribution frequency in the unit grid according to the node density. Mark the area where the force difference in each cell exceeds the set limit as a disturbance unit and output the slope response disturbance level classification layer.

[0052] S4: Call the disturbance level marker area in the slope response disturbance level classification layer, obtain the elevation change rate and horizontal displacement trajectory of the slope grid nodes within the disturbance boundary (elevation change rate and horizontal displacement trajectory are standard observation indicators in slope deformation monitoring, widely used in GNSS, InSAR, tilt sensors and other monitoring methods), calculate the elevation increment ratio and horizontal displacement direction sequence of the nodes in adjacent time periods respectively, and select the continuous deformation node set based on the same direction of the two changes to form the slope structure deformation trend path group;

[0053] S5: Based on the deformation trend path group of the slope structure, record the change data of the structural tension direction of the continuous deformation nodes in the time series (structural tension direction is a commonly used term in the field of structural engineering, used to describe the direction of tension transmitted by tension members (such as cables and rods) in the structure, which can be calculated based on the node force data), draw the tension offset trajectory diagram according to the order of tension direction vector transformation in the path, and mark the node coordinates and component connection sequence of each offset point in the path, mark the tension turning point and the intersection of the path, organize the path sequence of the whole area and output the dynamic sequence diagram of the multi-segment structural response of the repair area.

[0054] The slope structure partitioning unit matrix includes reconstructed grid boundaries, slope aspect consistency identifiers, and unit topological relationships. The component slope coupling contact layout map specifically includes a three-dimensional coordinate set of anchor points, a contact area parameter matrix, and a component grouping coding table. The slope response disturbance level classification layer includes a principal stress direction vector field, a contact rate density distribution cloud map, and a disturbance unit boundary coordinate set. The slope structure deformation trend path group specifically includes an elevation change rate gradient map, a displacement direction sequence matrix, and a deformation node topological network. The multi-segment structural response dynamic sequence map of the repair area includes a tension direction vector trajectory map, a turning node spatial coordinate set, and a path intersection segment topological index.

[0055] Please see Figure 2 The steps for obtaining the slope structure partitioning element matrix are as follows:

[0056] S101: Obtain the gray value matrix of remote sensing image and the three-dimensional point cloud elevation dataset. Compare the image gray value and point cloud elevation value at the corresponding position of each grid node pixel by pixel, calculate the absolute value of the difference between the two values, establish a heat map of the distribution of grid node difference values, filter the node positions where the difference value exceeds the terrain change threshold, and generate a terrain difference coefficient matrix.

[0057] The call retrieves the coverage coordinate range (east longitude). ,north latitude ) to (East Longitude) ,north latitude The remote sensing image grayscale matrix and 3D point cloud elevation dataset for the region were first preprocessed to a resolution of [missing information]. The grayscale values ​​of the remote sensing images range from The 8-bit integer range is linearly normalized to The floating-point range is used, and the 3D laser point cloud data is processed into a regular shape using Kriging interpolation. Grid, generate elevation value matrix, for specified grid nodes, such as node ( The corresponding image grayscale value is After normalization, it becomes The elevation value after interpolation of its point cloud data is Next, the two values ​​are compared, and the absolute value of the difference is calculated. ,in This is the elevation scaling factor. Set as the lowest elevation in the region The regional elevation range is ,but Calculations yielded Traverse all grid nodes and calculate The values ​​are organized into a matrix corresponding to the grid, and a heatmap of the distribution of grid node differences is established. Based on this, a terrain abrupt change threshold is set. This threshold was determined through statistical analysis of all the differences. Values ​​are all mean of values Plus Double standard deviation For example, statistical calculations of the entire region and ,but Finally, all were filtered out. The node position, for example, node of If the node is marked as a terrain abrupt change point, then the logical values ​​of all marked nodes ( Indicates mutation, (Indicates non-mutation) is stored in a binary matrix to generate a terrain difference coefficient matrix.

[0058] S102: Based on the terrain difference coefficient matrix, the slope change rate of adjacent nodes is calculated along the longitude direction. The sliding window method is used to detect abrupt regions where the slope change rate exceeds a set number of consecutive nodes. The coordinate sequence of abrupt points is recorded. Adjacent abrupt points are connected to form a polyline trajectory. The curvature continuity of the polyline is optimized by cubic spline interpolation to obtain the terrain polyline feature set.

[0059] Based on the topographic difference coefficient matrix, along the longitude direction of the geographic coordinates ( The slope change rate is calculated for adjacent nodes marked as terrain abrupt change points in the matrix (axis). This process uses a width of... Nodes ( A sliding window is used to scan column by column from west to east. Within the window, the center node is calculated first. Its two adjacent nodes before and after and slope and ,in For node elevation, Node spacing ( Then calculate the slope change rate. For example, the center node of a window elevation The elevations of its preceding and following nodes are respectively and ,but ,and The corresponding slope change rate Set a threshold for the rate of change of slope. This threshold is based on samples of stable slopes and historical landslide areas in the region. Comparative analysis of value distributions was conducted, and values ​​that could distinguish between them were selected. The critical value of the sample is determined when the sliding window detects a continuous... Each node The values ​​are all greater than At that time, the first and last nodes of the region, for example and Mark the endpoints of the mutation region and record their coordinate sequence. Then, record the coordinates of all adjacent mutation points, such as... The initial polyline trajectory is formed by direct connection. To make the trajectory curve smoother at the nodes, the polyline trajectory is processed by cubic spline interpolation. The coefficients of each cubic polynomial are determined by solving a tridiagonal linear equation system to ensure the continuity of the second derivative at the connection points, and the final optimized terrain polyline feature set is obtained.

[0060] S103: Call the terrain polyline feature set as the grid division boundary, perform topological reconstruction on the original grid, calculate the elevation gradient vector of all nodes in the new grid cell, count the deviation of the gradient vector azimuth angle in the cell, retain the grid cells with deviation values ​​less than the set angle threshold, merge adjacent cells to form a continuous region, and output the slope structure partition cell matrix.

[0061] The obtained terrain polyline feature set is used as the mandatory grid boundary for the original [grid]. For square meshes, topology reconstruction is performed. The original polyline-traversed mesh cells are divided into two or more new irregular polygonal mesh cells. After reconstruction, the elevation gradient vector of all nodes within each new mesh cell is calculated. ,in and The nodes are respectively and The rate of change of elevation in a direction, for example, for a reconstructed triangular element. It contains three nodes. Calculate the azimuth angle of its gradient vector respectively. Statistics for this unit Standard deviation of the azimuth angle of the gradient vector of all nodes Calculated Set an azimuth deviation threshold. This threshold is based on the definition of "slope aspect consistency" in geotechnical engineering, which means that the slope aspect swing range is less than... (Corresponding standard deviation approximately) The area within (within) can be considered the same slope, because Therefore, the grid cells are retained. Perform this check on all grid cells, and then process all cells that pass the check and are spatially adjacent, such as... and (That ), perform topological merging to form a larger, continuous region with consistent internal slope, and assign unique numbers to all merged continuous regions, for example The final output is a two-dimensional array that records the numbers of the slope-oriented regions to which each grid cell belongs, i.e., the slope structure partitioning cell matrix.

[0062] Please see Figure 3 The steps for obtaining the component slope coupling contact layout diagram are as follows:

[0063] S201: Call the slope structure partition unit matrix, extract the three-dimensional coordinate dataset of the anchor points of the vine cube structure, spatially match the coordinates of each anchor point with the corresponding grid node, calculate the elevation surface fitting plane of all nodes in the grid unit where the anchor point is located, count the number of contact points projected onto the fitting plane of the anchor point, and generate a set of contact area vectors.

[0064] The slope structure partitioning element matrix is ​​invoked, and a pre-defined 3D coordinate dataset of the anchor points of the vine cube structure is extracted. This dataset contains the ID of each anchor point and its... Coordinates, such as anchor points The coordinates are The coordinates of the anchor point are spatially matched with the grid nodes in the slope structure zoning element matrix to determine... It falls on the number Within the partitioned unit, then, calculation All nodes within the given grid cell (e.g., with) The elevation surface fitting plane (with nodes) is given by the equation of the plane. Its coefficient Through this The three-dimensional coordinates of each node are fitted using the least squares method. This aims to obtain a reference plane that can represent the local slope morphology of this small area, and to statistically analyze the anchor points (assuming their foundation is...). The number of contact points (of a square) projected onto the fitting plane is determined by using the area within the projection region. The spacing generates a sampling point array. If the distance between the sampling point and the fitting plane in the Z direction is less than... If a contact point is not found, it is considered a valid contact point, such as an anchor point. Co-generated within its projection area There are 10 sampling points, of which 10 are 100. If any one is counted as a valid contact point, then the total number of contact points is 1. Organize the IDs of all anchor points and the corresponding number of contact points into a vector, for example... Generate a vector set of contact areas.

[0065] S202: Based on the contact area vector set, extract the spacing data of adjacent anchor points along the component boundary, measure the length difference of the overlap area at the end of the component, record the overlap length variation curve in the order of component number, establish a component connection sequence index table, and output the connection sequence sequence.

[0066] Based on the contact area vector set, the spacing data between adjacent anchor points is extracted along the layout path of the rattan cube component. For example, the component... Connecting anchor points and Its standard design spacing is The actual three-dimensional spatial distance is obtained through coordinate calculation. Meanwhile, the end of the measuring component is at and The difference in the overlap length of the anchorage zone is obtained by analyzing images of the relative positions of the component and the anchorage foundation. For example, the component... exist The overlap length at the end is ,exist The overlap length at the end is Then the difference in their lengths is According to the actual connection sequence of the components (e.g., from the top of the slope to the bottom of the slope). (), sequentially record the change in overlap length of each component to form an overlap length change curve, for example, the sequence data is Based on this connection order, a component connection order index table is created. This table clearly indicates the preceding and following component numbers of each component, as shown in Table 1. The index table and the corresponding overlap length and other data are output, which together constitute the connection order sequence.

[0067] Table 1: Example of the connection sequence index for rattan cube components

[0068]

[0069] As shown in Table 1, this table records the connection relationship of components along a path and the slight changes in their overlap length, providing basic data for subsequent structural continuity analysis and mechanical transfer path identification.

[0070] S203: Combine the connection sequence to calculate the cosine value of the spatial angle between the layout direction vector of each component and the corresponding slope aspect vector. Components with a cosine value greater than a set threshold are classified into the same group. Components in the same group are mapped with continuous area numbers. The contact area and connection sequence data are integrated to generate a component slope coupling contact layout map.

[0071] Based on the generated connection sequence, calculate the layout direction vector for each component. For example, components The direction vector is determined by the anchor point. point to The normalized vector is determined, and at the same time, the region where the component is located is extracted from the slope structure partitioning unit matrix. Average slope aspect vector This vector is obtained by averaging the gradient vectors of all nodes within the region, and its direction represents the macroscopic flow direction or the direction of maximum gradient in that partition. Next, the cosine of the spatial angle between these two vectors is calculated. For example, calculating the components of ,member of And another component of Set a cosine threshold for the included angle. This threshold corresponds to approximately The included angle is set based on on-site construction specifications, which require that the angle between the main direction of the component and the main slope direction of the slope should not be too large. Components, such as and Components are categorized into the same group (e.g., "slope group"), and spatially continuous components within the same group are assigned a unified area numbering system. This grouping and numbering information is then integrated with the contact area vector set and connection sequence data obtained in the previous steps. Specifically, a "group ID" field and a "slope consistency" field are added to each component record. Finally, a multi-dimensional data set containing the component geometry, connection, contact and coupling relationship with the terrain is generated, namely the component slope coupling contact layout map.

[0072] Please see Figure 4 The steps for obtaining the slope response disturbance level classification layer are as follows:

[0073] S301: Call the component slope coupling contact layout map, extract the principal tensile stress direction vector and normal stress scalar value of the contact area node, calculate the spatial cosine value of the angle between the principal tensile stress direction of each node and the corresponding component direction vector, count the frequency of cosine value distribution, and generate a set of force direction consistency coefficients.

[0074] The component slope coupling contact layout map is called. For each component and its contact area with the slope recorded in the map, the principal tensile stress direction vectors of all mesh nodes in the area are extracted through the finite element analysis software interface. and contact normal stress scalar value For example, for components A contact node below Its principal tensile stress direction vector is Its contact normal stress value is Next, calculate the principal tensile stress direction vector of this node. With the component to which this node belongs Direction vector Cosine value of the spatial angle between (Assuming the vectors have been normalized), the calculation yields... , for components All under the contact area Repeat this calculation for each node to obtain a set of cosine values. Calculate the frequency distribution of values ​​in this set, for example, by... The interval is divided into Each sub-interval is counted, and the number of occurrences within each interval is recorded. The number of values ​​is used to obtain a frequency distribution histogram, which includes all nodes. The values ​​are collected to form an array that corresponds one-to-one with the node position, generating a set of force direction consistency coefficients.

[0075] S302: Based on the force uniformity coefficient set, the node density value is statistically calculated according to the grid cell, the number of contact nodes and the frequency of normal stress distribution within a unit area are calculated, the product of the number of contact nodes and the frequency of normal stress is used as the contact rate index, the grid contact rate distribution matrix is ​​established, and the contact rate gradient field is output.

[0076] Based on the force orientation consistency coefficient set, The analysis grid is used as a unit for statistical analysis. First, the node density value within each analysis grid unit is calculated. That is, the number of contact nodes falling into the unit, for example, unit Inside If there are 1 contact node, then Simultaneously, calculate the normal stress values ​​of all contact nodes within the unit. The distribution frequency, specifically statistical The value is greater than a certain reference value (e.g.) The proportion of nodes that have made effective contact to the total number of contact nodes. For example, unit Inside The normal stress at each node is greater than ,but The contact ratio index of the element is obtained by multiplying the node density with the frequency of normal stress distribution. Perform the same calculation on all analyzed mesh elements to build a two-dimensional mesh contact rate distribution matrix. Each element of this matrix represents the contact rate distribution of the corresponding element. Finally, gradient calculation is performed on the matrix to calculate the contact rate of each unit. and The rate of change in direction is output as a visualization layer, namely the contact rate gradient field.

[0077] S303: Combine the contact rate gradient field to calculate the absolute value of the difference in force direction consistency coefficient between adjacent grid cells, screen regions where the difference exceeds the stress disturbance threshold, perform topological merging on continuous exceeding regions, assign values ​​according to disturbance intensity, and generate a slope response disturbance level classification layer.

[0078] Based on the output contact rate gradient field, the absolute value of the difference in force direction consistency coefficients between adjacent analysis grid elements is calculated. Specifically, for two adjacent elements... and Calculate the force uniformity coefficient of all nodes within it. average and Then calculate the difference. ,For example, , ,but Set a stress disturbance threshold The threshold value was set with reference to the gradient variation characteristics of stress concentration regions in materials mechanics. It was derived through ROC analysis of simulation data from numerous stable and unstable slope repair cases. This value can achieve [the desired effect]. The accuracy of the calculation distinguishes stress anomaly transition zones. The boundary between these two units is then marked as a perturbation region. This filtering is performed on all adjacent units, and all continuous regions marked as perturbations are topologically merged to form independent perturbation patches. The perturbation patches are then analyzed based on their internal structure. The average value is used to classify and assign values ​​to the disturbance intensity. For example, the average Defined as a Level 1 disturbance, average Defined as a second-order disturbance, average The disturbances are classified as Level 3, and this classification information is rendered onto the geographic layer with different colors or attribute values ​​to generate a slope response disturbance level classification layer.

[0079] Please see Figure 5 The steps for obtaining the slope structure deformation trend path group are as follows:

[0080] S401: Call the slope response disturbance level classification layer, extract the time series data of elevation change rate of grid nodes within the disturbance boundary, calculate the ratio of node elevation increment in adjacent time periods, establish node elevation change trend vector, and generate elevation dynamic trend set.

[0081] The generated slope response disturbance level classification layer is invoked, and all areas marked as level two or higher disturbances are identified. All grid nodes within the boundaries of these areas are then extracted in continuous... One monitoring cycle (e.g., once a month), to Elevation change rate time series data, which are derived from simulation or actual monitoring, for example, node... exist The elevation at that moment is ,exist Time for ,exist Time for Their elevation change rates are respectively , Next, the ratio of the node elevation increments within adjacent time periods is calculated. ,for time, This ratio reflects the trend of deformation rate. If the ratio is greater than 1, it indicates accelerated deformation, and if it is less than 1, it indicates decelerated deformation. The ratio is calculated for all consecutive time periods of all nodes in the disturbance zone to form a set containing node ID, timestamp, and incremental ratio. This is used to establish the node elevation change trend vector, and the vectors of all nodes are collected to generate the elevation dynamic trend set.

[0082] S402: Based on the dynamic trend set of elevation, obtain the horizontal displacement trajectory coordinates of the same node in a continuous time period, calculate the cosine value of the angle between the displacement direction vectors of adjacent time periods, count the number of consecutive time periods with the cosine value greater than the displacement in the same direction threshold, and output the horizontal displacement in the same direction sequence.

[0083] Based on the dynamic trend set of elevation, we focus on those exhibiting a continuously accelerating deformation trend (i.e.) For nodes in two consecutive periods, obtain their values ​​within the same consecutive time period ( to The coordinates of the horizontal displacement trajectory within the node, for example. The horizontal coordinates from of Become of , and then become of Calculate its in displacement direction vector during the time period and displacement direction vector during the time period Calculate the cosine of the angle between these two vectors. Set a displacement in the same direction threshold. This value corresponds to approximately The included angle is determined based on statistics of the deformation characteristics of creep landslides, whose early displacement directions are highly consistent, for nodes. The entire time series ( arrive Perform this calculation and statistics. Values ​​continuously greater than The number of time periods, in terms of nodes and For example, the calculation results of the horizontal displacement unidirectionality over a continuous period are shown in Table 2. The statistical quantity of all target nodes is used as their unidirectionality index, and the horizontal displacement unidirectionality sequence is output.

[0084] Table 2: Examples of Nodal Horizontal Displacement Same Direction Sequence

[0085]

[0086] As shown in Table 2, nodes exist The cosine value of the included angle is greater than the threshold for three consecutive time periods. Therefore, the number of consecutive same-direction time periods is counted as And nodes exist cosine of the angle between time periods When the value is less than the threshold, continuity is interrupted, and the number of consecutive same-direction time periods is [number missing]. .

[0087] S403: Combining the horizontal displacement homogeneity sequence, nodes whose elevation increment ratio changes synchronously with the cosine value of the horizontal displacement direction are selected. Spatial clustering analysis is performed on nodes that meet the synchronization condition, and adjacent nodes are connected to form a continuous path to generate a slope structure deformation trend path group.

[0088] Based on the output horizontal displacement unidirectional sequence, the nodes are finally selected. The selection criteria are: in continuous... The ratio of elevation increments of nodes within one or more time periods logarithm The direction of change and the cosine value of the horizontal displacement direction The direction of change increases synchronously; the specific judgment process is to compare... and At two moments, if and If so, then synchronization is considered to occur within a certain time period. The logarithmic transformation here is to handle... Nonlinear growth, for example, a node at Moment They are respectively , , To ensure synchronous growth across two consecutive time periods, all nodes meeting this condition are marked as continuously deforming nodes. Density-based spatial clustering analysis (DBSCAN) is then performed on all marked nodes, with a search radius set. (equivalent to) (multiple component length) and minimum number of neighbors The core nodes within the same cluster are connected according to their spatial proximity to form one or more continuous paths. These paths represent the potential, physically related overall deformation areas of the slope, generating a group of slope structure deformation trend paths.

[0089] Please see Figure 6 The steps for obtaining the dynamic sequence diagram of the multi-segment structural response in the repair region are as follows:

[0090] S501: Call the slope structure deformation trend path group, extract the time series tension direction vector of continuous deformation node, calculate the angle change between vectors in adjacent time periods in time order, draw the tension direction offset trajectory polyline graph in each path, and generate tension offset trajectory set.

[0091] Call the slope structure deformation trend path group, and for each path in the path group, for example path Extract all continuously deformed nodes on the path in the time series ( to The structural tension direction vector calculated by the structural analysis module in () This vector represents the main direction of tension transmission of the rattan cube component at this node. The change in the angle between this vectors in adjacent time periods is calculated sequentially. For example, nodes In the path superior, The tension vector at time t is , Time becomes The change in the included angle This will display the data of all nodes along a path at different times. The values ​​are connected and plotted as a line graph of the tension direction offset trajectory, where the X-axis represents the cumulative length along the path, the Y-axis represents time, and the Z-axis (or represented by color) represents... The size is used to generate such charts for all paths, which are then aggregated into a set of tension offset trajectories.

[0092] S502: Based on the tension offset trajectory set, associate the node coordinates of the path where each offset point is located with the component connection sequence number, establish a mapping relationship table between the spatial position of the offset point and the component sequence, mark the node position where the tension direction changes abruptly by more than a set angle, and output the tension turning node index.

[0093] Based on the tension offset trajectory set, for each offset point in the figure, i.e., each Combine and associate their three-dimensional spatial coordinates within the path. And the component connection sequence number to which the node belongs (from S202), for example, node The offset point, whose coordinates are... The component number is First, establish a mapping table from the offset point to its physical location and structural affiliation. Then, set a threshold for the angle of sudden change in tension direction. This threshold is set by analyzing the nonlinear mechanical behavior of the component material to determine the minimum rotation angle that may cause local stress redistribution, and all offset points are checked. Value, if If the node is marked as a "tension inflection node" at that moment, and after filtering all the nodes on the path, the partial content of the generated tension inflection node index is shown in Table 3. The index file is then output.

[0094] Table 3: Example of Tension Inflection Node Index

[0095]

[0096] See Table 3. This index records the nodes identified as tension inflection points and their key spatiotemporal and structural information, such as nodes. exist The change in the angle between the tension direction and the moment Exceed They are marked based on thresholds, and this information will be used for subsequent high-risk segment identification.

[0097] S503: Combining the tension turning point index, the density of turning points in the path intersection area is statistically analyzed, the intersection segment level is divided according to the density threshold, the spatiotemporal evolution data of the entire path is integrated, and a dynamic sequence diagram of the multi-segment structural response in the repair area is generated.

[0098] By combining the tension inflection node index, we analyze and statistically analyze the areas where paths intersect in space, such as paths and In a certain There is overlap within the region. The number density of all tension inflection points within this intersection region is statistically analyzed. For example, a total of [number missing] were found within this region. If there are several tension inflection points, then its density is... Based on this density value, a threshold is set to classify the response level of the intersection segment. For example, if the density is less than... In the low response zone, arrive The area between is the medium response zone, greater than The high-response zone is then divided into a medium-response zone. Finally, the spatiotemporal evolution data of all paths are integrated, including the geometry of the path, the deformation rate of nodes on the path (from S401), the tension offset trajectory (from S501), and the level labels of turning points and intersection sections (generated in this step). Through a 3D visualization engine, this multi-dimensional information is superimposed on the slope model and dynamically played using the time axis as the controller, thereby generating a dynamic sequence diagram of the multi-segment structural response in the repair area.

[0099] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in the present invention should be included within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be determined by the scope of the claims.

Claims

1. A dynamic simulation method for slope restoration effect using a granular material, characterized in that, Includes the following steps: S1: Obtain remote sensing images and point cloud data of the area where the vine cube structure is laid, sort the grid nodes according to the elevation gradient, extract the terrain polyline to reconstruct the slope grid points, divide the area with consistent slope aspect, and generate a slope structure partition unit matrix. The steps for obtaining the slope structure partitioning unit matrix are as follows: S101: Obtain the gray value matrix of remote sensing image and the three-dimensional point cloud elevation dataset. Compare the image gray value and point cloud elevation value at the corresponding position of each grid node pixel by pixel, calculate the absolute value of the difference between the two values, establish a heat map of the distribution of grid node difference values, filter the node positions where the difference value exceeds the terrain change threshold, and generate a terrain difference coefficient matrix. S102: Based on the terrain difference coefficient matrix, the slope change rate of adjacent nodes is calculated along the longitude direction. The sliding window method is used to detect abrupt regions where the slope change rate exceeds a set number of consecutive nodes. The coordinate sequence of abrupt points is recorded. Adjacent abrupt points are connected to form a broken line trajectory. The curvature continuity of the broken line is optimized by cubic spline interpolation to obtain the terrain broken line feature set. S103: Call the terrain polyline feature set as the grid division boundary, perform topological reconstruction on the original grid, calculate the elevation gradient vector of all nodes in the new grid cell, count the deviation of the gradient vector azimuth angle in the cell, retain the grid cells with deviation values ​​less than the set angle threshold, merge adjacent cells to form a continuous region, and output the slope structure partition cell matrix. S2: Based on the slope structure partition unit matrix, extract the three-dimensional coordinate data of the anchor points of the rattan cube structure, calculate the slope contact area and overlap length, group and mark the rattan cube structure components according to the angle between the layout direction of the rattan cube structure components and the slope direction vector, and generate a component slope coupling contact layout map. S3: Call the component slope coupling contact layout map, obtain the principal tensile stress direction and contact normal stress distribution of the nodes in the contact area of ​​the vine cube structure component, calculate the angle between the direction vector of the vine cube structure component and the principal stress direction of the node, determine the consistency distribution, count the contact rate and distribution frequency in the unit grid, mark the disturbance level of the area where the angle deviates from the set limit, and output the slope response interference level classification layer. S4: Call the disturbance level marker area in the slope response disturbance level classification layer, obtain the elevation change rate of the slope grid nodes and the horizontal displacement trajectory within a continuous time period, extract the nodes according to the consistency of deformation direction, and form a slope structure deformation trend path group. S5: Based on the slope structure deformation trend path group, record the change data of the structural tension direction of the continuous deformation nodes in the time series, draw the tension offset trajectory diagram according to the tension direction vector transformation order in the path, and mark the node coordinates of each offset point in the path and the connection order of the Cube structural components in parallel, mark the tension turning nodes and the intersection of the path, organize the path sequence of the whole area and output the dynamic sequence diagram of the multi-segment structural response of the repair area. The dynamic sequence diagram of the multi-segment structural response in the repair area includes a tension direction vector trajectory diagram, a set of spatial coordinates of turning nodes, and a topological index of path intersection segments.

2. The dynamic simulation method for slope restoration effect of Tenglifang according to claim 1, characterized in that: The slope structure partitioning unit matrix includes reconstructed grid boundaries, slope aspect consistency identifiers, and unit topological relationships. The component slope coupling contact layout map specifically includes a three-dimensional coordinate set of anchor points, a contact area parameter matrix, and a component grouping coding table. The slope response interference level classification layer includes a principal stress direction vector field, a contact rate density distribution cloud map, and a disturbance unit boundary coordinate set. The slope structure deformation trend path group specifically includes an elevation change rate gradient map, a displacement direction sequence matrix, and a deformation node topological network.

3. The dynamic simulation method for slope restoration effect of Tenglifang according to claim 1, characterized in that, The steps for obtaining the component ramp coupling contact layout diagram are as follows: S201: Call the slope structure partition unit matrix, extract the three-dimensional coordinate dataset of the anchor points of the vine cube structure, spatially match the coordinates of each anchor point with the corresponding grid node, calculate the elevation surface fitting plane of all nodes in the grid unit where the anchor point is located, count the number of contact points projected onto the fitting plane of the anchor point, and generate a contact area vector set. S202: Based on the contact area vector set, extract the spacing data between adjacent anchor points along the boundary of the rattan cube structure component, measure the length difference of the overlapping area at the end of the rattan cube structure component, record the overlap length variation curve according to the rattan cube structure component numbering order, establish a rattan cube structure component connection sequence index table, and output the connection sequence sequence. S203: Based on the connection sequence, calculate the cosine value of the spatial angle between the layout direction vector of each rattan cube structure component and the corresponding slope aspect vector. Classify rattan cube structure components with a cosine value greater than a set threshold into the same group. Perform continuous area numbering mapping on rattan cube structure components in the same group, integrate contact area and connection sequence data, and generate a component slope coupling contact layout map.

4. The dynamic simulation method for slope restoration effect of Tenglifang according to claim 1, characterized in that, The steps for obtaining the slope response disturbance level classification layer are as follows: S301: Call the component slope coupling contact layout map, extract the principal tensile stress direction vector and normal stress scalar value of the contact area node, calculate the spatial cosine value of the angle between the principal tensile stress direction of each node and the direction vector of the corresponding rattan cube structure component, count the frequency of cosine value distribution, and generate a set of force direction consistency coefficients. S302: Based on the set of force consistency coefficients, the node density values ​​are statistically analyzed according to the grid cells, the number of contact nodes and the frequency of normal stress distribution within a unit area are calculated, the product of the number of contact nodes and the frequency of normal stress is used as the contact rate index, the grid contact rate distribution matrix is ​​established, and the contact rate gradient field is output. S303: Based on the contact rate gradient field, calculate the absolute value of the difference in force direction consistency coefficient between adjacent grid cells, filter out regions where the difference exceeds the stress disturbance threshold, perform topological merging on continuous exceeding regions, assign values ​​according to disturbance intensity, and generate a slope response disturbance level classification layer.

5. The dynamic simulation method for slope restoration effect of Tenglifang according to claim 1, characterized in that, The steps for obtaining the deformation trend path group of the slope structure are as follows: S401: Call the slope response disturbance level classification layer, extract the time series data of elevation change rate of grid nodes within the disturbance boundary, calculate the ratio of node elevation increment in adjacent time periods, establish node elevation change trend vector, and generate elevation dynamic trend set. S402: Based on the elevation dynamic trend set, obtain the horizontal displacement trajectory coordinates of the same node in a continuous time period, calculate the cosine value of the angle between the displacement direction vectors of adjacent time periods, count the number of consecutive time periods where the cosine value is greater than the displacement same direction threshold, and output the horizontal displacement same direction sequence. S403: Based on the horizontal displacement homogeneity sequence, select nodes whose elevation increment ratio changes synchronously with the cosine value of the horizontal displacement direction, perform spatial clustering analysis on the nodes that meet the synchronization condition, connect adjacent nodes to form a continuous path, and generate a slope structure deformation trend path group.

6. The dynamic simulation method for slope restoration effect of Tenglifang according to claim 1, characterized in that, The steps for obtaining the dynamic sequence diagram of the multi-segment structural response in the repair region are as follows: S501: Call the slope structure deformation trend path group, extract the time series tension direction vector of continuous deformation node, calculate the angle change between vectors in adjacent time periods in time order, draw a line graph of tension direction offset trajectory in each path, and generate tension offset trajectory set. S502: Based on the tension offset trajectory set, associate the node coordinates of the path where each offset point is located with the component connection sequence number, establish a mapping relationship table between the spatial position of the offset point and the component sequence, mark the node position where the tension direction changes abruptly by more than a set angle, and output the tension turning node index. S503: Combining the tension turning node index, the density of turning nodes in the path intersection area is statistically analyzed, the intersection segment level is divided according to the density threshold, the spatiotemporal evolution data of the entire path is integrated, and a dynamic sequence diagram of the multi-segment structural response of the repair area is generated.

Images