A regional scale landslide hydraulic triggering modeling method considering deep heterogeneity

Abstract

The application discloses a regional scale landslide hydraulic triggering modeling method considering depth heterogeneity, comprising the following steps: acquiring meteorological, terrain, soil and vegetation characteristics and the like as input layer data of a subsequent model; obtaining soil water content and pore water pressure of each time step of a target region based on a GEOtop 3.0 model; constructing a landslide hydraulic triggering process model VH-LHT considering soil property along depth heterogeneity, and based on the input layer data and output results of the GEOtop 3.0 model at a current time step, the safety factor and damage depth of a landslide unit at each layer are synchronously obtained. The application adopts the above regional scale landslide hydraulic triggering modeling method considering depth heterogeneity, which not only considers the change of parameters along the depth, but also realizes the force transmission between soil columns, realizes seamless coupling with the GEOtop 3.0 model, improves the modeling ability of depth heterogeneity, and improves the spatial prediction precision of group landslide.

Description

Technical Field

[0001] This invention relates to the field of landslide hydraulic triggering process monitoring technology, and in particular to a regional-scale landslide hydraulic triggering modeling method that considers depth heterogeneity. Background Technology

[0002] Physics-based process modeling is achieved by coupling a distributed hydrological model with a slope stability model. The distributed hydrological model calculates the soil water potential and moisture content of the target area. After cell mapping, the slope stability model calculates the slope stability state of each cell, and finally evaluates the impact of soil moisture dynamics on the timing, rate, volume, and location of shallow landslides.

[0003] Due to limitations in computational efficiency and disciplinary constraints, early hydrological models focused on calculating the sources and sinks of the entire watershed, treating groundwater flow and saturated-unsaturated seepage in the soil layer as black boxes. However, with the development of computer performance and the increasing interdisciplinarity of disciplines, hydrological models have gradually transformed into process models with clear physical concepts and have become increasingly coupled with models in fields such as climate, geology, and ecology.

[0004] Based on these hydrological models, landslide researchers have conducted physical modeling of slope instability at the watershed scale. To explain the hydraulic triggering process of rainfall-induced landslides, many studies have employed more refined coupled modeling, such as coupled models using the 3DRichards equation solver and infinite slope stability analysis. Subsequently, optimization of slope stability models has also been gradually implemented. A landslide hydraulic-mechanical coupled triggering model (LHT model) was established based on the theory of sand pile element automata, discretizing the soil mass above the bedrock into hexagonal prisms interconnected by mechanical "bonds." This model improves upon the traditional infinite slope model, reproducing the interaction process between local soil strength failure and forces between soil columns, and has been widely applied in areas with thin weathered layers. To expand the model's application scope, a slip surface search method was coupled into the LHT model, improving the problem of overestimating landslide volume and thickness in areas with thick weathered layers.

[0005] However, the coupled modeling method of the deep layered landslide hydraulic-mechanical coupled triggering model (D-LHT model) based on LHT and combined with the SWAT model does not spatially discretize along the slope normal direction (geocentric direction), and therefore cannot store heterogeneous information of soil at different depths (e.g., soil dry density, hydraulic conductivity, root density, etc.). This makes the model appear to achieve dynamic search of the slip surface, but in reality, the soil column stability calculation results of the model at one time step are only functions of depth and water content, which is not applicable to areas with thick weathered layers and high vertical heterogeneity.

[0006] For example, the high vegetation cover in the low mountain and hilly areas of southeastern coastal provinces of China results in root density decreasing with depth, significantly increasing the gradient of parameters such as hydraulic conductivity and cohesion between root soil and underlying rootless soil. If all subsurface heterogeneity is attributed to average or "effective" parameters, it may fail to capture the "fault" phenomenon in soil hydraulic parameters caused by vegetation and the root-controlled sliding surface location. Secondly, the SWAT+D-LHT modeling strategy inputs the soil moisture content and surface water head of the hydrological response units calculated by the SWAT model into the D-LHT model, and then divides saturated and unsaturated zones by calculating the wetting front depth. This method may lead to greater deviations between the inversion results and the true values ​​due to overestimation of soil unit weight and simplification of the unsaturated seepage process.

[0007] Based on the problems of the aforementioned models, this invention aims to construct a dynamically graded landslide hydraulic-mechanical coupled triggering model (VH-LHT) that considers depth heterogeneity, achieving three-dimensional hierarchical data mapping with a comprehensive hydrological model. Then, the VH-LHT model is coupled with the GEOtop 3.0 model to establish a coupled process model for landslide hydraulic triggering processes with high spatiotemporal resolution and considering depth heterogeneity, specifically for landslides caused by frequent rainfall. Summary of the Invention

[0008] The purpose of this invention is to provide a regional-scale hydraulic triggering modeling method for landslides that considers depth heterogeneity, solving the transient seepage error caused by the reliance on the "wetting front approximation method" in traditional models and the slope stability calculation error caused by the homogenization of input parameters along the depth; supporting automatic switching between saturated and unsaturated states to accurately reflect the hysteresis effect of rainfall infiltration; improving the ability to model depth heterogeneity and increasing the accuracy of spatial prediction of cluster landslides.

[0009] To achieve the above objectives, this invention provides a regional-scale hydraulic triggering modeling method for landslides that considers depth heterogeneity, comprising the following steps:

[0010] Step S1: Obtain meteorological data, topographic data, soil data, and vegetation feature data as input layer data for subsequent models;

[0011] Step S2: Based on the GEOtop 3.0 model, obtain the soil moisture content and pore water pressure of the target area at each time step, and use them as input for subsequent models;

[0012] Step S3: Construct a landslide hydraulic triggering process model VH-LHT that considers the heterogeneity of soil properties along depth. Based on the input layer data and the output results of the GEOtop 3.0 model at the current time step, obtain the safety factor and failure depth of each landslide unit layer by layer.

[0013] Preferably, in step S1, meteorological data, topographic data, soil data, and vegetation characteristics are acquired as input layer data for subsequent models;

[0014] Meteorological data includes, but is not limited to, temperature and rainfall; topographic data includes, but is not limited to, digital elevation models, aspect, and slope; soil data includes, but is not limited to, hydrological characteristics and soil and rock characteristics; vegetation characteristics include, but are not limited to, leaf area index, evapotranspiration, and root area ratio; in addition to the above data types, it also includes, but is not limited to, land use and soil type.

[0015] Preferably, in step S2, the input layer data is preprocessed based on the GEOtop 3.0 model, and the specific process is as follows:

[0016] First, land use types are mapped to a two-dimensional grid, and data is assigned to each data type;

[0017] Secondly, the soil types are mapped to the initial values ​​of the three-dimensional grid, and data is assigned to each data type;

[0018] Finally, numerical settings are made to simulate control conditions, including processing time step, maximum number of iterations, convergence residual, and matching simulation start and end dates.

[0019] The time step can be adjusted while adhering to the CFL condition. Each grid cell at a certain time step uses the D8 topology algorithm to calculate surface runoff as a surface computing unit, while the three-dimensional grid below it is a subsurface computing unit. By calculating the water and energy flux of the computing units, the inversion and prediction of the real world can be achieved.

[0020] Preferably, the nonlinear partial differential equations are solved by a three-dimensional Richards equation solver, and the runoff on the grid surface in the target area and the soil moisture content and pore water pressure at different depths below it are output along the time series.

[0021] Numerical solutions are obtained based on the Newton-Raphson iterative method, and the Jacobian matrix is ​​calculated iteratively using the BI-CGSTAB method.

[0022] Preferably, in step S3, based on the sand pile unit automata model, the implementation method of the LHT model is optimized, and the construction of the landslide hydraulic triggering process model VH-LHT, which considers the heterogeneity of soil properties along the depth, is realized by using object-oriented programming method; the VH-LHT model will be called after each time step of the GEOtop 3.0 model calculation.

[0023] Preferably, based on the constructed VH-LHT model, the input layer data and the output results of the GEOtop 3.0 model at the current time step are preprocessed, and the specific process is as follows:

[0024] First, instantiate landslide elements, map element static parameters, construct element matrix, and assign static and dynamic parameters to each layer of landslide elements;

[0025] Then, numerical settings are made, including processing time step, spatial resolution, maximum number of iterations, and convergence residual;

[0026] Finally, load the dynamic parameters for the current time step and map them to the landslide element;

[0027] To balance the information abundance of the integrated hydrological model output and the mechanical transmission mechanism of the LHT model, the VH-LHT model is based on the concept of spatiotemporal cubes. It discretizes the plane space with the grid resolution of the digital elevation model and regards the entire soil column under the digital elevation model as a landslide unit. Each landslide unit can store Z-dimensional information in layers, achieving seamless coupling with the GEOtop 3.0 model.

[0028] Preferably, based on the established landslide unit object, each soil layer of the landslide unit can be mapped one-to-one with unsaturated water content and dry soil unit weight. When calculating the total stress borne by each layer of the landslide unit, the stress contributed by the dry soil is used. Stress contributed by water Additional stress from vegetation get;

[0029] Among them, the stress contributed by dry soil Stress contributed by water The value is obtained through hierarchical cumulative summation, as shown below:

[0030] ;

[0031] in, The density is measured in dry soil. The specific gravity of water; The depth of the current layer is obtained by adding half the thickness of the current layer to the thickness of each layer above it, and is aligned with the direction of gravity. For the first Layer thickness; This is the index of the layer above the current layer; The current layer thickness; This is the current layer number; For the first Layer moisture content; This represents the volumetric water content of the current layer.

[0032] Preferably, the safety factor FS of each landslide unit is obtained by the limit equilibrium method, assuming the sliding surface is parallel to the slope surface, and the safety factor is the resistance force. With sliding force The ratio;

[0033] In the VH-LHT model, anti-slip force The expression is as follows:

[0034] ;

[0035] in, Based on soil shear strength Root system reinforcement effect Additional anti-sliding force provided by surrounding soil columns It consists of three parts; Effective normal stress; The effective friction angle; For effective cohesion; This is a correction factor for the overestimation in the model; It is the average tensile strength of the root system per unit area of ​​soil; Root orientation factor; This represents the root area ratio;

[0036] Sliding force The expression is as follows:

[0037] ;

[0038] in, From soil self-weight stress to component The slope-direction component of additional stress from vegetation Additional sliding force provided by surrounding soil columns composition;

[0039] Therefore, the safety factors of the VH-LHT model under the three states of unsaturation, transient saturation, and stable saturation are as follows:

[0040] ;

[0041] in, This represents the pore water pressure in the current layer, which is negative under unsaturated conditions. Slope; This represents the residual volumetric moisture content. This represents the saturated volumetric water content. This refers to the depth of the groundwater level.

[0042] Preferably, based on the initial safety factor distribution of a single time step, the landslide element is subjected to force transmission iteration, and the safety factor FS and failure depth of the landslide element in each time step are finally obtained;

[0043] The landslide element transmits forces in eight directions and has one potential instability direction; the eight landslide elements adjacent to the current landslide element are considered as neighboring elements, corresponding one-to-one with the direction; the potential instability direction of the current landslide element and A 45° deviation is considered downhill, while the other five directions are uphill.

[0044] If the current landslide unit becomes unstable, the stabilizing unit on the upslope provides the current landslide unit with additional anti-sliding force not exceeding the root tensile strength, and the stabilizing unit on the downslope provides the current landslide unit with additional anti-sliding force other than that used for its own stability.

[0045] If the current landslide unit remains unstable after reinforcement, the residual sliding force will be transferred to the landslide unit in the downslope.

[0046] Therefore, if the current landslide element exists in the downslope of an unstable element in the upslope, the unstable element will not only fail to provide additional anti-sliding force, but will also exert compressive stress on it.

[0047] Therefore, the present invention employs the above-mentioned regional-scale landslide hydraulic triggering modeling method that considers depth heterogeneity, and the beneficial effects are as follows:

[0048] (1) The soil parameters are stored by three-dimensional hierarchical discretization, which overcomes the homogeneity assumption of traditional models and improves the ability to model deep heterogeneity; it accurately captures the parameter gradient changes caused by vegetation roots and avoids ignoring the phenomenon of "hydraulic parameter faults".

[0049] (2) The unsaturated zone water transport is dynamically calculated by the three-dimensional Richards equation solver, which solves the transient seepage error caused by the traditional model relying on the "wetting front approximation method"; it supports automatic switching between saturated and unsaturated states and accurately reflects the hysteresis effect of rainfall infiltration.

[0050] (4) Based on the infinite slope method, the safety factor is calculated in layers and the mechanical formula is dynamically selected according to the soil condition, which significantly reduces the overestimation of landslide volume and thickness; through the eight-directional force transmission of hexagonal landslide units, the progressive failure process of soil is reproduced, and the spatial prediction accuracy of group landslides is improved.

[0051] (5) Using the netCDF data interface of the GEOtop 3.0 model, the hydrological parameters and the mechanical model are mapped in layers to avoid data redundancy; the Newton-Raphson iterative method and the BI-CGSTAB sparse matrix solver are used to optimize the computational efficiency while ensuring accuracy.

[0052] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. Attached Figure Description

[0053] Figure 1This is a flowchart of the GEOtop 3.0+VH-LHT coupling model of the present invention;

[0054] Figure 2 This is an overview of the framework of the GEOtop 3.0 model of this invention; wherein, (a) is a three-dimensional spatial discretization map based on DEM; (b) is the computational unit structure and data mapping mechanism; and (c) is the code architecture of the hydrological process module.

[0055] Figure 3 This is a schematic diagram of the modeling process of the VH-LHT model of the present invention; wherein, (a) is the output data of the GEOtop 3.0 model; (b) is the construction of landslide elements; (c) is the distribution diagram of the initial safety factor FS; (d) is the force transmission mechanism between landslide elements; and (e) is the output of the VH-LHT model. Detailed Implementation

[0056] The technical solution of the present invention will be further described below with reference to the accompanying drawings and embodiments.

[0057] like Figure 1 As shown, a regional-scale hydraulic triggering modeling method for landslides considering depth heterogeneity includes the following steps:

[0058] Step S1: Acquire meteorological data, topographic data, soil data, and vegetation characteristic data, etc., as input layer data for subsequent models;

[0059] Step S2: Based on the GEOtop 3.0 model, obtain the soil moisture content and pore water pressure of the target area at each time step, and use them as input for subsequent models;

[0060] Step S3: Construct a regional-scale landslide hydraulic triggering process model (VH-LHT model) that considers the heterogeneity of soil properties along depth. Based on the input layer data and the output results of the GEOtop 3.0 model at the current time step, obtain the safety factor and failure depth of each landslide unit.

[0061] Example

[0062] Step S1: Acquire meteorological data, topographic data, soil data, and vegetation characteristic data, etc., as input layer data for subsequent models.

[0063] Meteorological data includes temperature and rainfall; topographic data includes digital elevation models, aspect, and slope; soil data includes hydrological characteristics and soil and rock properties; vegetation characteristics include leaf area index, evapotranspiration, and root area ratio; in addition to the above data types, it also includes land use and soil type.

[0064] Step S2: Based on the GEOtop 3.0 model, obtain the soil moisture content and pore water pressure of the target area at each time step, and use them as input for subsequent models.

[0065] GEOtop is built on the C / C++ language. Its improved version 3.0 uses object-oriented programming to re-encapsulate the code of modules such as hydrological processes and energy exchange, making it easier for users to deploy and perform secondary development on the Linux operating system.

[0066] The GEOtop 3.0 model integrates the advantages of hydrological models and land surface models (LSMs), realizing the coupled inversion of energy transfer and hydrological processes from the surface to the subsurface under complex terrain conditions at the watershed scale. Based on the raster resolution of the digital elevation model (DEM), the real space is discretized in planar space, and then the subsurface space is discretized layer by layer along the geocentric direction, such as... Figure 2 As shown in (a) in the figure.

[0067] Step S21: Based on the GEOtop 3.0 model, preprocess the input layer data.

[0068] First, land use types are mapped to a two-dimensional grid, and data is assigned to each data type.

[0069] Secondly, the soil types are mapped to the initial values ​​of the three-dimensional grid, and data is assigned to each data type.

[0070] Finally, numerical settings are performed to simulate control conditions, including processing time step, maximum number of iterations, convergence residual, etc., and to match the simulation start and end dates.

[0071] The time step can be adjusted while adhering to the CFL conditions. At a given time step, each grid cell uses the D8 topology algorithm to calculate surface runoff as surface computational units, while the underlying 3D grid serves as the subsurface computational unit. For example... Figure 2 As shown in (b) of the diagram. , Corresponding to the spatial resolution of the DEM, The discrete thickness corresponding to each soil layer can be adjusted according to data accuracy and computational efficiency. The model achieves real-world inversion and prediction by calculating the water and energy fluxes of the calculation units.

[0072] Step S22: Perform energy and mass balance using the meteorological event-driven model at each time step.

[0073] This invention focuses on changes in pore water pressure (matrix suction in the unsaturated zone) and soil moisture content caused by groundwater hydrological processes. The GEOtop 3.0 model uses finite difference equations to represent the general partial differential equations of the Richards equations and energy conservation equations, employing an implicit backward Euler method to ensure numerical stability. Time discretization is performed on a specific three-dimensional mesh element within the spatially discretized computational domain. The calculation point is located at the geometric center of the three-dimensional mesh.

[0074] Therefore, the constructed difference equation expression is as follows:

[0075] ;

[0076] in, and They represent time steps, yes and At some point in the middle, It is the time step; It is a set of adjacent three-dimensional mesh cells ( (index of) , and These represent the current adjacent 3D mesh cell sets along... , and The amount of directional deflection; This represents the distance between two calculation points in a grid cell, calculated using dx, dy, and dz. For source and sink residuals; For mass / energy conservation residuals; Representing the field variables that vary with time and space, for the Richards equation corresponding to soil water potential For the energy conservation equation corresponding to temperature ; Indicates from grid ( ) points to adjacent grids ( The flux in the direction of the Richards equation corresponds to the permeability coefficient. For the energy conservation equation, the corresponding thermal conductivity coefficient ; It is a nonlinear function of the state variable, corresponding to the soil volumetric water content in the Richards equation. For the energy conservation equation, the corresponding internal energy density .

[0077] Step S23: Solve the nonlinear partial differential equations using a three-dimensional Richards equation solver, and output the runoff on the grid surface in the target area and the soil moisture content and pore water pressure at different depths below it along the time series.

[0078] (1) Numerical solution is performed based on the Newton-Raphson iteration method. Since the model considers the lateral seepage process, the value of N is 6, including a three-dimensional mesh in two vertical directions and four lateral directions, such as Figure 2 As shown in (b) of the diagram.

[0079] (2) The Jacobian matrix is ​​iteratively calculated using the BI-CGSTAB method to reduce the waste of computational resources caused by sparse matrices. For example... Figure 2 As shown in (c), the key code logic architecture of the reconstructed hydrological process module is presented.

[0080] Because the model uses a three-dimensional Richards equation solver for continuous transient solutions from the saturated to the unsaturated state, it sacrifices computational efficiency in the steady-state computational domain (saturated zone). However, the dynamically solved hydraulic conductivity, soil moisture content, and matrix suction allow the model to better simulate the actual seepage conditions in the unsaturated zone, both horizontally and vertically, solving the problem that traditional hydrological models can only allow the wetting front to advance further after the soil enters the saturated state. Therefore, the GEOtop 3.0 model is very suitable for uncovering the "black box" of soil moisture content changes in the unsaturated zone, which is precisely the prerequisite for landslide researchers to model the hydraulic triggering process of landslides.

[0081] The GEOtop 3.0 model employs a netCDF (Network Common Data Format) external library for efficient storage of grid data. Watershed-scale GIS data, including land use, slope, and soil properties, can be directly mapped to a computational domain constructed based on a DEM discrete grid after preprocessing. Point-scale field or experimental data can be mapped to different land use types and layer-by-layer soil types, ultimately achieving the mapping of all data and parameters. For example... Figure 2 As shown on the left side of (b) in the image.

[0082] This storage and mapping method differs from models that classify hydrological units based on posterior knowledge (such as SWAT). It stores more raw and heterogeneous subsurface information in each 3D grid, making it possible for landslide researchers to develop progressive destructive hydraulic triggering models with soil columns as computational units.

[0083] Step S3: Construct a landslide hydraulic triggering process model (VH-LHT model) that considers the heterogeneity of soil properties along depth. Based on the input layer data and the output results of the GEOtop 3.0 model at the current time step, obtain the safety factor and failure depth of each landslide unit layer by layer.

[0084] Step S31: Construct a landslide hydraulic triggering process model that considers the heterogeneity of soil properties along depth (VH-LHT model).

[0085] Since traditional hydrological models are not conducive to revealing small-scale spatial heterogeneity and surface-to-ground hydrological processes, this invention develops a dynamic graded landslide hydraulic-mechanical coupling triggering model (VH-LHT model) to meet the coupling requirements of fine hydromechanical models and achieve high spatial resolution mapping of water content and pore water pressure calculation units that cover vegetation disturbance.

[0086] Based on the sand pile automata model, the implementation method of the LHT model was optimized, and object-oriented programming was used to construct the model. The VH-LHT model will be called after each time step of the GEOtop 3.0 model computation.

[0087] Step S32: Based on the constructed VH-LHT model, preprocess the input layer data and the output results of the GEOtop 3.0 model.

[0088] First, instantiate landslide elements, map element static parameters, construct element matrix, and assign static and dynamic parameters to each layer of the landslide elements.

[0089] Then, numerical settings are made, including processing time step, spatial resolution, maximum number of iterations, convergence residual, etc.

[0090] Finally, load the dynamic parameters for the current time step and map them to the landslide element.

[0091] To balance the information abundance of the integrated hydrological model output with the mechanical transmission mechanism of the LHT model, the VH-LHT model is based on the concept of a spatiotemporal cube. It discretizes the plane space with the grid resolution of the DEM and treats the entire soil column under the DEM as a landslide unit. Each landslide unit can store Z-dimensional information in layers, which can not only consider the changes of parameters along the depth, but also realize the transmission of forces between soil columns, thus achieving more accurate prediction.

[0092] like Figure 3 As shown in (a) and (b), the GEOtop 3.0 model output contains multi-dimensional information, which is mapped layer by layer to the landslide element objects. It is precisely because the VH-LHT model considers the heterogeneity of soil properties along depth that it can be seamlessly coupled with the GEOtop 3.0 model, allowing the computational results of the GEOtop 3.0 model to be fully utilized, and ensuring one-to-one mapping of computational units without simplification. , Maintaining the same output resolution as the GEOtop 3.0 model. It can be consistent with the GEOtop 3.0 model, or it can be customized by the user.

[0093] Because the landslide element object is established in this way, each soil layer of the landslide element can be mapped one-to-one with unsaturated water content and dry soil unit weight. When calculating the total stress borne by each layer of the landslide element, the stress contributed by the dry soil can be used. Stress contributed by water Additional stress from vegetation This method yields stresses that better avoid errors caused by homogenization.

[0094] Among them, the stress contributed by dry soil Stress contributed by water The value is obtained through hierarchical cumulative summation, as shown below:

[0095] ;

[0096] in, The density is measured in dry soil. The specific gravity of water; The depth of the current layer is obtained by adding half the thickness of the current layer to the thickness of each layer above it, and is aligned with the direction of gravity. For the first Layer thickness; This is the index of the layer above the current layer; The current layer thickness; This is the current layer number; For the first Layer moisture content; This represents the volumetric water content of the current layer. Further parameter details are shown in the table.

[0097] Table 1 Parameter Description

[0098] ;

[0099] Step S33: Calculate the safety factor FS of each layer of the landslide unit using the limit equilibrium method.

[0100] Calculating the safety factor of landslide elements using the Infinity Slope Method .when At that time, the landslide unit is considered unstable; when The landslide unit is considered stable. It is assumed that the sliding surface is parallel to the slope surface, and the safety factor is equal to the resistance force. With sliding force The ratio.

[0101] In the VH-LHT model, anti-slip force The expression is as follows:

[0102] ;

[0103] in, Based on soil shear strength Root system reinforcement effect Additional anti-sliding force provided by surrounding soil columns It consists of three parts; The effective normal stress is a function of soil moisture content and pore water pressure or matrix suction, and has different functional relationships under unsaturated, transiently saturated and steady-state saturated conditions. The effective friction angle; For effective cohesion, root reinforcement follows the Wu / Waldron model. This is a correction factor for the overestimation in the model; It is the average tensile strength of the root system per unit area of ​​soil; Root orientation factor; This represents the root area ratio.

[0104] Sliding force The expression is as follows:

[0105] ;

[0106] in, From soil self-weight stress to component The slope-direction component of additional stress from vegetation Additional sliding force provided by surrounding soil columns composition.

[0107] Therefore, the safety factors of the VH-LHT model under the three states of unsaturation, transient saturation, and stable saturation are as follows:

[0108] ;

[0109] in, This represents the pore water pressure in the current layer, which is negative under unsaturated conditions. Slope; This represents the residual volumetric moisture content. This represents the saturated volumetric water content. This represents the groundwater level depth. Further parameter details are shown in Table 1. Due to the layered calculation method, the buoyant unit weight can be used when calculating the self-weight stress of soil layers below the groundwater level, avoiding overestimation of self-weight stress caused by homogenization calculations.

[0110] Based on the above parameter mapping and landslide element stability calculation ( =0, =0), thus obtaining the safety factor distribution map of the target calculation region along the depth and time series variation, such as Figure 3 As shown in (c) in the figure.

[0111] Step S34: Based on the initial safety factor distribution map of a single time step, the landslide element is subjected to force transmission iteration to finally obtain the safety factor FS and failure depth of the landslide element in each time step.

[0112] The VH-LHT model's landslide elements can transmit forces in eight directions, with one potential instability direction. The eight adjacent landslide elements are considered as neighborhood elements, corresponding one-to-one with the direction, such as... Figure 3 As shown in (d) in the diagram. The potential instability direction of the current landslide unit and... The direction with a 45° deviation is considered downhill, while the other five directions are uphill.

[0113] If the current landslide unit becomes unstable, the stabilizing units on the upslope can provide the current landslide unit with additional anti-sliding force not exceeding the tensile strength of the root system, and the stabilizing units on the downslope can provide the current landslide unit with additional anti-sliding force other than that used for its own stability.

[0114] If the current landslide unit remains unstable after reinforcement, the residual sliding force will be transferred to landslide units in the downslope. It can be inferred that if the current landslide unit exists in the downslope of an unstable unit in the upslope, the unstable unit will not only fail to provide additional anti-sliding force but will also exert compressive stress on it, making the current landslide unit even more difficult to stabilize. The model follows the natural law of force transmission and employs a multi-round, simultaneous allocation algorithm for additional anti-sliding force, more realistically simulating the stress change process of the soil.

[0115] Therefore, this invention adopts the above-mentioned regional-scale landslide hydraulic triggering modeling method that considers depth heterogeneity. It not only considers the changes in parameters along the depth, but also realizes the force transmission between soil columns. Furthermore, it achieves seamless coupling with the GEOtop 3.0 model, so that the calculation results of the GEOtop 3.0 model can be fully utilized. The calculation units are mapped one-to-one without simplification, achieving more accurate prediction, improving the deep heterogeneity modeling capability, improving the spatial prediction accuracy of cluster landslides, supporting automatic switching between saturated and unsaturated states, and accurately reflecting the hysteresis effect of rainfall infiltration.

[0116] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit them. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the technical solutions of the present invention, and these modifications or equivalent substitutions cannot cause the modified technical solutions to deviate from the spirit and scope of the technical solutions of the present invention.

Claims

1. A regional-scale hydraulic triggering modeling method for landslides considering depth heterogeneity, characterized in that, Includes the following steps: Step S1: Obtain meteorological data, topographic data, soil data, and vegetation feature data as input layer data for subsequent models; Step S2: Based on the GEOtop 3.0 model, obtain the soil moisture content and pore water pressure of the target area at each time step, and use them as input for subsequent models; The nonlinear partial differential equations are solved by a three-dimensional Richards equation solver, and the runoff on the grid surface in the target area and the soil moisture content and pore water pressure at different depths below it are output along the time series. Numerical solution is obtained based on the Newton-Raphson iterative method, and the Jacobian matrix is ​​calculated iteratively using the BI-CGSTAB method. Step S3: Construct a landslide hydraulic triggering process model VH-LHT that considers the heterogeneity of soil properties along depth. Based on the input layer data and the output results of the GEOtop 3.0 model at the current time step, simultaneously obtain the safety factor and failure depth of each landslide unit. In order to balance the output information abundance of the integrated hydrological model and the mechanical transmission mechanism of the LHT model, the VH-LHT model is based on the concept of spatiotemporal cubes. It discretizes the plane space with the grid resolution of the digital elevation model (DEM), and regards the entire soil column under the digital elevation model as a landslide unit. Each landslide unit stores Z-dimensional information in layers to achieve seamless coupling with the GEOtop 3.0 model. The safety factor of the landslide unit is calculated using the infinite slope method.

2. The regional-scale hydraulic triggering modeling method for landslides considering depth heterogeneity according to claim 1, characterized in that, In step S1, meteorological data, topographic data, soil data, and vegetation characteristics are acquired as input layer data for subsequent models. Meteorological data includes, but is not limited to, temperature and rainfall; topographic data includes, but is not limited to, digital elevation models, aspect, and slope; soil data includes, but is not limited to, hydrological characteristics and soil and rock characteristics; vegetation characteristics include, but are not limited to, leaf area index, evapotranspiration, and root area ratio; in addition to the above data types, it also includes, but is not limited to, land use and soil type.

3. The regional-scale hydraulic triggering modeling method for landslides considering depth heterogeneity according to claim 1, characterized in that, In step S2, the input layer data is preprocessed based on the GEOtop 3.0 model. The specific process is as follows: First, land use types are mapped to a two-dimensional grid, and data is assigned to each data type; Secondly, the soil types are mapped to the initial values ​​of the three-dimensional grid, and data is assigned to each data type; Finally, numerical settings are made to simulate control conditions, including processing time step, maximum number of iterations, convergence residual, and matching simulation start and end dates. The time step can be adjusted while adhering to the CFL condition. Each grid cell at a certain time step uses the D8 topology algorithm to calculate surface runoff as a surface computing unit, while the three-dimensional grid below it is a subsurface computing unit. By calculating the water and energy flux of the computing units, the inversion and prediction of the real world can be achieved.

4. The regional-scale hydraulic triggering modeling method for landslides considering depth heterogeneity according to claim 1, characterized in that, In step S3, based on the sand pile unit automata model, the implementation method of the LHT model is optimized. The VH-LHT model, which considers the heterogeneity of soil properties along depth, is constructed using an object-oriented programming method. The VH-LHT model will be called after each time step of the GEOtop 3.0 model calculation.

5. The regional-scale hydraulic triggering modeling method for landslides considering depth heterogeneity according to claim 4, characterized in that, Based on the constructed VH-LHT model, the input layer data and the output results of the GEOtop 3.0 model at the current time step are preprocessed. The specific process is as follows: First, instantiate landslide elements, map element static parameters, construct element matrix, and assign static and dynamic parameters to each layer of landslide elements; Then, numerical settings are made, including processing time step, spatial resolution, maximum number of iterations, and convergence residual; Finally, load the dynamic parameters for the current time step and map them to the landslide element.

6. The regional-scale hydraulic triggering modeling method for landslides considering depth heterogeneity according to claim 5, characterized in that, Based on the established landslide unit object, each soil layer of the landslide unit can be mapped one-to-one with unsaturated water content and dry soil unit weight. When calculating the total stress borne by each layer of the landslide unit, the stress contributed by the dry soil is used. Stress contributed by water Additional stress from vegetation get; Among them, the stress contributed by dry soil Stress contributed by water The value is obtained through hierarchical cumulative summation, as shown below: ； in, The density is measured in dry soil. The specific gravity of water; The depth of the current layer is obtained by adding half the thickness of the current layer to the thickness of each layer above it, and is aligned with the direction of gravity. For the first Layer thickness; This is the index of the layer above the current layer; The current layer thickness; This is the current layer number; For the first Layer moisture content; This represents the volumetric water content of the current layer.

7. A regional-scale landslide hydraulic triggering modeling method considering depth heterogeneity according to claim 6, characterized in that, The safety factor FS for each layer of the landslide unit is obtained using the limit equilibrium method. Assuming the sliding surface is parallel to the slope surface, the safety factor is the resistance force. With sliding force The ratio; In the VH-LHT model, anti-slip force The expression is as follows: ； in, Based on soil shear strength Root system reinforcement effect Additional anti-sliding force provided by surrounding soil columns It consists of three parts; Effective normal stress; The effective friction angle; For effective cohesion; This is a correction factor for the overestimation in the model; It is the average tensile strength of the root system per unit area of ​​soil; Root orientation factor; This represents the root area ratio; Sliding force The expression is as follows: ； in, From soil self-weight stress to component The slope-direction component of additional stress from vegetation Additional sliding force provided by surrounding soil columns composition; Therefore, the safety factors of the VH-LHT model under the three states of unsaturation, transient saturation, and stable saturation are as follows: ； in, This represents the pore water pressure in the current layer, which is negative under unsaturated conditions. Slope; This represents the residual volumetric moisture content. This represents the saturated volumetric water content. This refers to the depth of the groundwater level. This represents the volumetric water content of the current layer. The stress contributed to dry soil; The stress contributed by the water; Adding stress to vegetation; The depth of the current layer is obtained by adding half the thickness of the current layer to the thickness of each layer above it, and is aligned with the direction of gravity.

8. The regional-scale hydraulic triggering modeling method for landslides considering depth heterogeneity according to claim 7, characterized in that, Based on the initial safety factor distribution of a single time step, the landslide element is subjected to force transmission iteration, and finally the safety factor FS and failure depth of the landslide element in each time step are obtained. The landslide element transmits forces in eight directions and has one potential instability direction; the eight landslide elements adjacent to the current landslide element are considered as neighboring elements, corresponding one-to-one with the direction; the potential instability direction of the current landslide element and A 45° deviation is considered downhill, while the other five directions are uphill. If the current landslide unit becomes unstable, the stabilizing unit on the upslope provides the current landslide unit with additional anti-sliding force not exceeding the root tensile strength, and the stabilizing unit on the downslope provides the current landslide unit with additional anti-sliding force other than that used for its own stability. If the current landslide unit remains unstable after reinforcement, the residual sliding force will be transferred to the landslide unit in the downslope. Therefore, if the current landslide element exists in the downslope of an unstable element in the upslope, the unstable element will not only fail to provide additional anti-sliding force, but will also exert compressive stress on it.

Images