A flood routing method for steep slope-short flow path basin of tropical island

Abstract

The present application relates to the technical field of hydrological simulation, and particularly relates to a flood routing simulation method for a steep slope-short process basin of a tropical island. The method comprises the following steps: calculating a steep slope index and completing partitioning based on a basin river network and sub-basin boundaries; constructing a variable dimension hydrodynamic model of a motion wave, a one-dimensional Saint-Venant equation and a two-dimensional shallow water equation for different units; introducing distributed rainfall to drive runoff calculation of each unit; applying Courant number constraint to steep slope units for local sub-time step iteration correction; and finally solving to obtain flood peak arrival time, inundation range and water depth distribution results. The present application realizes high-precision and stable simulation of a flood rapid concentration process of a tropical island steep slope-short process basin by constructing a variable dimension hydrodynamic coupling model based on steep slope partitioning, and combining high-precision automatic extraction of a basin river network and local sub-time step iteration calculation under Courant number constraint, thereby significantly improving the accuracy and adaptability of flood peak evolution and inundation range prediction.

Description

Technical Field

[0001] This invention relates to the field of hydrological simulation technology, and in particular to a flood evolution simulation method for steep slope-short flow basins in tropical islands. Background Technology

[0002] Tropical island watersheds are generally characterized by steep terrain, short confluence paths, concentrated rainfall, and extremely rapid flood response times. This allows floods to generate and converge rapidly and evolve quickly to downstream areas, exhibiting significant nonlinearity and spatial variability. Existing flood evolution simulation methods mostly employ uniform-dimensional hydrodynamic models and fixed-time-step numerical solutions, which struggle to simultaneously adapt to the differentiated hydrodynamic characteristics of rapid slope confluence and concentrated river transport. This results in insufficient characterization of flood propagation processes under complex terrain conditions, limiting simulation accuracy.

[0003] Meanwhile, traditional methods lack adaptive control mechanisms for local velocity changes in steep slope areas during the time progression, which can easily lead to decreased numerical stability or accumulation of calculation errors in high-slope units, thus affecting the reliable prediction of flood peak arrival time and inundation range. Summary of the Invention

[0004] Therefore, it is necessary to provide a flood evolution simulation method for steep slope-short flow basins in tropical islands to solve at least one of the above-mentioned technical problems.

[0005] To achieve the above objectives, a flood evolution simulation method for steep slope-short flow basins in tropical islands is proposed, the method comprising the following steps: Step S1: Calculate the slope and flow length based on the obtained watershed river network and sub-watershed boundaries to obtain the steep slope index; use the steep slope index to divide the watershed into steep slope units, fast-flowing river units, and downstream slow-flowing or floodplain units. Step S2: Establish a variable-dimensional hydrodynamic model based on different units. The steep slope unit adopts the kinematic wave model, the rapid river unit adopts the one-dimensional Saint-Venant equation, and the downstream slow-flow or floodplain unit adopts the two-dimensional local inertial shallow water equation. Step S3: Construct distributed rainfall input and calculate the runoff generation process of each unit; Step S4: Use the slope information of steep slope units to perform Coulomb number constraints, and perform local sub-time step iterative calculations on steep slope units based on time step constraints to obtain the runoff generation process of steep slope units after iteration. Step S5: Solve the runoff generation process of each unit to obtain the peak arrival time, inundation range and water depth distribution results.

[0006] The present invention has the following beneficial effects: I. This invention introduces a watershed zoning based on steep slope index and a variable-dimensional hydrodynamic coupling modeling mechanism to achieve differentiated characterization of the "steep slope-short flow-rapid confluence" watershed of tropical islands. This enables slopes, channels and floodplains to be simulated using appropriate physical equations within a unified framework, thereby significantly improving the ability to express the non-uniform response characteristics of complex terrain and the overall consistency of calculations in the flood evolution process.

[0007] Second, this invention constructs an automatic extraction method for river networks and sub-basins based on high-resolution lidar point clouds, and combines it with a parallel accumulation mechanism of flow direction matrix and confluence intensity. This achieves high-precision automatic identification of river network structure and stable tracking of water flow path, effectively enhancing the topological connectivity expression ability of water flow convergence process under complex slope conditions, and improving the analytical accuracy of the model for rapid confluence process of short-duration floods.

[0008] Third, this invention introduces a local sub-time step iteration mechanism for slope units based on Courant number constraints, which enables the time step to be adaptively adjusted according to the local slope and propagation speed. This improves the computational resolution of steep slope areas while ensuring numerical stability, effectively avoids the numerical dissipation and computational distortion problems that are prone to occur in traditional unified time step methods, and significantly improves the accuracy and reliability of flood peak arrival time and inundation range prediction. Attached Figure Description

[0009] Figure 1 A schematic diagram of the steps for a flood evolution simulation method for steep slope-short flow basins in tropical islands; Figure 2 for Figure 1 A detailed flowchart illustrating the implementation steps of step S4. Figure 3 This is a schematic diagram of a round-trip flight path scan in one embodiment; Figure 4 This is a schematic diagram illustrating the generation of a flow direction matrix in one embodiment; The realization of the objective, functional features and advantages of the present invention will be further explained in conjunction with the embodiments and with reference to the accompanying drawings. Detailed Implementation

[0010] The technical method of the present invention will now be clearly and completely described with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of the present invention. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without inventive effort are within the scope of protection of the present invention.

[0011] Furthermore, the accompanying drawings are merely illustrative of the invention and are not necessarily drawn to scale. The same reference numerals in the drawings denote the same or similar parts, and therefore repeated descriptions of them will be omitted. Some block diagrams shown in the drawings are functional entities and do not necessarily correspond to physically or logically independent entities. These functional entities can be implemented in software, in one or more hardware modules or integrated circuits, or in different network and / or processor methods and / or microcontroller methods.

[0012] It should be understood that although the terms "first," "second," etc., may be used herein to describe various units, these units should not be limited by these terms. These terms are used merely to distinguish one unit from another. For example, without departing from the scope of the exemplary embodiments, a first unit may be referred to as a second unit, and similarly, a second unit may be referred to as a first unit. The term "and / or" as used herein includes any and all combinations of one or more of the associated listed items.

[0013] To achieve the above objectives, please refer to Figures 1 to 4 A flood evolution simulation method for steep slope-short flow basins in tropical islands, the method comprising the following steps: Step S1: Calculate the slope and flow length based on the obtained watershed river network and sub-watershed boundaries to obtain the steep slope index; use the steep slope index to divide the watershed into steep slope units, fast-flowing river units, and downstream slow-flowing or floodplain units. Step S2: Establish a variable-dimensional hydrodynamic model based on different units. The steep slope unit adopts the kinematic wave model, the rapid river unit adopts the one-dimensional Saint-Venant equation, and the downstream slow-flow or floodplain unit adopts the two-dimensional local inertial shallow water equation. Step S3: Construct distributed rainfall input and calculate the runoff generation process of each unit; Step S4: Use the slope information of steep slope units to perform Coulomb number constraints, and perform local sub-time step iterative calculations on steep slope units based on time step constraints to obtain the runoff generation process of steep slope units after iteration. Step S5: Solve the runoff generation process of each unit to obtain the peak arrival time, inundation range and water depth distribution results.

[0014] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. This watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of 38°, characterized by steep slopes, short course, and dramatic changes in river cross-section. First, the digital elevation model (DEM), river network data, and sub-watershed boundary information of the watershed are acquired. Based on the DEM data, the slope and course length of each grid cell are calculated, and the steepness index is calculated according to the formula: ,in For unit slope, This represents the length of the flow from the unit to the confluence point of the river channel. Based on the steepness index, the watershed is divided into three types of units: steep slope units, fast-flowing river units, and downstream slow-flowing or floodplain units. Steep slope units are mainly distributed in the upstream mountainous areas, accounting for approximately 32% of the watershed area; fast-flowing river units extend along the main river channels, accounting for 20%; and downstream slow-flowing or floodplain units account for 48%, mainly consisting of flat land and low-lying coastal areas.

[0015] Corresponding variable-dimensional hydrodynamic models were established for different units: the steep slope unit used a kinematic wave model to describe the rapid propagation of floodwater along the steep slope, considering gravity-driven and frictional effects; the fast-flowing river unit used the one-dimensional Saint-Venant equation, combined with unsteady flow boundary conditions to solve for changes in river level and discharge; the downstream slow-flowing and floodplain units used two-dimensional local inertial shallow water equations, with a gridded water depth field to simulate the floodplain process. The models were coupled with boundary discharge and water level to ensure the continuity of hydrodynamics between the upstream and downstream areas of the basin.

[0016] Rainfall input uses high-resolution radar rainfall data or historical heavy rainfall event data for distributed rainfall simulation, with values ​​assigned according to grid cells. In steep slope cells, the runoff generation coefficient is calculated based on local slope and hydrodynamic characteristics, and the initial runoff is generated by combining rainfall intensity. To ensure the computational stability of steep slope cells, the Coulomb number (CFL) constraint is introduced, and the local time step of the calculation cell is: In the formula, For time step, As an empirical safety factor, For local slope, For the speed of water propagation, For spatial grid scale, It is the acceleration due to gravity. The water is deep.

[0017] For elements that do not meet the time step constraints, local sub-time step iterative calculations are used until the global time step is completed, in order to ensure the stability and accuracy of flood propagation.

[0018] In the fast-flowing zone and downstream floodplain, a set of hydrodynamic equations was established using gridded river sections and floodplains, and an explicit-implicit coupling method was employed to solve for water level evolution. The river nodes and the floodplain grid boundaries were coupled through water volume continuity to ensure consistency between runoff generation and floodplain processes. The entire flood evolution process of the watershed was recorded using time series data, including water depth, flow velocity, and peak arrival time for each grid.

[0019] To verify the simulation accuracy, the model was calibrated using historical rainfall and tidal observation data. Simulation results show that, under the design rainfall intensity of 120 mm / h, the average time for the flood peak to reach the downstream estuary in the upstream steep slope unit is approximately 18 minutes, the peak water depth in the rapid current zone is approximately 1.4 m, and the maximum inundation depth in the downstream floodplain is approximately 0.85 m. The flood inundation area covers low-lying residential areas and parts of commercial areas. These simulation results can be used for further flood control scheduling, dam design, and emergency response plan development.

[0020] Preferred methods for obtaining the river network and sub-basin boundaries of a watershed include: The control drone platform is equipped with a lidar scanning component to perform round-trip route scanning of the watershed, continuously acquire surface point cloud data, and cache it through the onboard storage unit; A continuous digital elevation model is generated by registration and fusion of surface point cloud data. The elevation difference operation is performed on the digital elevation model grid by grid to output the water flow direction information of each unit and form a flow direction matrix. The flow accumulation operation is performed based on the flow direction matrix. The flow intensity distribution is generated by superimposing the contributions of upstream units in parallel, and the river network structure of the watershed is extracted. Based on the river network structure of the basin, the location of the flow direction boundary is determined and the region is divided to obtain the sub-basin boundary.

[0021] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. The watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river cross-sections within the watershed change dramatically, and the course is short and steep, exhibiting typical characteristics of steep-slope flood evolution. To accurately obtain the watershed's river network and sub-watershed boundaries, this embodiment uses a UAV equipped with a high-precision lidar scanning component to perform round-trip scanning of the entire watershed, collecting high-resolution point cloud data covering the watershed surface at a sampling interval of approximately 0.5 m. The data is cached in real-time using an onboard storage unit, and the flight path timestamps and location information are recorded simultaneously.

[0022] The collected point cloud data, after feature point matching, iterative nearest point (ICP) fine registration, and anomaly removal, generates a continuous digital elevation model (DEM) covering the entire watershed. For missing points or outliers in the DEM, spatial interpolation methods are used to complete them, ensuring the complete representation of steep slopes, river cross-sections, and tidal flat topographic features. Subsequently, the DEM is rasterized with each grid cell having a side length of 5m. The elevation difference between each grid cell and its eight neighboring grid points is calculated to determine the direction of the fastest descent of water flow, generating a flow direction matrix to accurately depict the spatial distribution path of water flow in steep slopes, rapids, and slow-flowing areas.

[0023] Based on the flow direction matrix, a flow accumulation calculation is performed, and the contribution of upstream grid cells is superimposed in parallel to form a flow intensity distribution map. Within this map, thresholds are set to extract the main channel and tributary network. Simultaneously, considering channel slope, width, and cross-sectional variation characteristics, rapid flow zones and slow flow zones are identified, providing foundational river network data for subsequent hydrodynamic simulations. Based on the extracted river network structure, the flow direction boundaries are further analyzed and regions are segmented to generate sub-basin boundaries. Each sub-basin corresponds to a set of catchment units, which can be used for subsequent division of steep slope units, rapid flow units, and floodplain units, thus providing high-precision foundational data for flood evolution simulation and flood control planning.

[0024] Preferably, controlling the unmanned aerial vehicle platform equipped with a lidar scanning component to perform round-trip route scanning of the watershed includes: Import watershed boundary data into the ground control terminal, generate non-equidistant flight strips along the main confluence direction based on the slope aspect grid calculation results, and write the start and end coordinates of the flight strips into the flight control system; After receiving the flight path coordinates, the UAV flight control system drives the power unit to take off and uses the inertial measurement unit and satellite positioning module for combined navigation to enable the UAV to fly stably along the first flight path; During flight, the lidar scanning component is controlled to perform continuous lateral sweeps and simultaneously acquire attitude angle data to compensate for the scanning angle in real time. When the flight control system identifies the coordinates of the current flight path's end point, it controls the UAV to perform a turnaround maneuver based on slope constraints, switching the flight direction to the adjacent flight path; During the flight of adjacent flight strips, the flight altitude and lateral offset are adjusted to make the adjacent scan strips form a preset overlapping area; Repeatedly execute the flight strip and turnaround switching until the watershed is fully covered and continuous point cloud data is output.

[0025] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. This watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river channel cross-section varies dramatically, and the course is short and steep. To obtain high-precision surface topographic information, one can refer to... Figure 3 First, the UAV platform, equipped with a lidar scanning component, conducts aerial surveys of the entire watershed. Watershed boundary information is imported into the ground control terminal, and non-equidistant flight strips are generated along the main confluence direction based on slope aspect grid calculations. The spacing of the flight strips is dynamically adjusted according to the slope and a preset overlap ratio, and the start and end coordinates are written into the flight control system. The UAV's flight altitude is set to an average altitude of 120m relative to the ground, the flight speed is 5~7m / s, the lidar scanning frequency is set to 200kHz, and the scanning angle is ±30° to ensure a point cloud density of at least 30 points per square meter.

[0026] While the UAV flies stably along the first flight strip, the lidar continuously performs lateral scans, simultaneously collecting attitude angle data (pitch, yaw, roll), and corrects for laser point projection errors in real time. Upon reaching the end of the flight strip, the flight control system controls the UAV to perform a slope-turn maneuver, adjusting its flight altitude and lateral offset to ensure at least 30% overlap between adjacent scan strips. This process is repeated until full coverage is achieved. Flight data is cached in the onboard storage unit, recording GPS position, laser reflection intensity, and flight attitude information in real time.

[0027] After acquiring continuous point cloud data, point cloud registration is performed based on lidar echo intensity and GPS / IMU information. A frame-by-frame and global optimization algorithm is used to fuse the data, eliminating inter-line errors and generating a continuous digital elevation model (DEM) with a spatial resolution of 0.5m × 0.5m. Subsequently, elevation differences are calculated grid-by-grid on the DEM to obtain the flow direction information for each unit and construct a flow direction matrix. Then, the contribution of upstream units is accumulated in parallel to generate the confluence intensity distribution. Based on the confluence intensity threshold, the main river network structure is extracted, confirming river inflow points and boundary locations. Combined with the flow direction matrix, regional division is performed to obtain complete sub-basin boundaries.

[0028] Preferably, the elevation difference operation is performed grid-by-grid on the digital elevation model to output the water flow direction information of each unit and form a flow direction matrix, including: Load the digital elevation model; The terrain is divided into regular grid cells with a fixed spatial resolution, and a grid neighborhood index relationship is established. For each grid cell and its eight neighboring cells, the elevation difference is calculated to obtain the direction of the minimum elevation neighboring cell, and the direction of the minimum elevation neighboring cell is used as the main water flow direction of the grid cell. The remaining grid cells are arranged in a matrix based on the direction of the main water flow and the spatial position of the grid cell to form a flow direction matrix.

[0029] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. The watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river cross-section within the watershed varies dramatically, the course is short and steep, and the terrain is complex, making it suitable for high-precision three-dimensional watershed topographic analysis. The watershed's digital elevation model (DEM) was obtained through the aforementioned UAV lidar scanning and point cloud registration fusion, with a spatial resolution of 0.5 m × 0.5 m.

[0030] For detailed procedures, please refer to [link / reference]. Figure 4 First, the DEM is divided into regular raster cells with a fixed spatial resolution, and an eight-neighbor cell index is created for each raster cell, with the neighbor order being Northeast, East, Southeast, South, Southwest, West, Northwest, and North. Then, the elevation drop is calculated by performing a difference operation on the elevation values ​​of each raster cell and its eight neighbor cells. ,in( ) represents the current grid cell, ( () represents an eight-neighborhood unit, with elevation values and All data are taken in meters from the DEM data. Among the differences between eight neighboring units, the direction with the largest elevation drop is selected as the main flow direction for the current unit. This direction number is mapped to a value in the flow direction matrix (e.g., D8 coding, where 1-8 correspond to the eight neighboring directions). If the maximum elevation drop is equal in multiple directions, the direction consistent with the overall main flow direction of the watershed is prioritized to ensure the continuity of the flow.

[0031] To ensure the accuracy of boundary cell processing, external sentinel cells are set for the watershed boundary and river outlet grids to prevent flow direction calculations from exceeding the boundaries. Simultaneously, a slope weighting factor is introduced for steep slope areas with a gradient exceeding 30°, and flow direction is selected through weighted elevation differences to reflect the rapid runoff characteristics of steep slopes.

[0032] After determining the main flow direction grid by grid, all grid cells are arranged in the order of the DEM row and column indices to form a two-dimensional flow direction matrix. Each cell in the flow direction matrix records the main flow direction and the corresponding grid number, which can be directly used for subsequent runoff accumulation calculations and river network extraction. This matrix can not only accurately represent the spatial distribution of water flow in short-course steep-slope watersheds, but also provide accurate grid-level hydrodynamic input for flood runoff simulation.

[0033] Preferably, calculating the elevation difference between each grid cell and its eight neighboring cells further includes: Identify the center grid based on the current grid cell; A 3×3 neighborhood calculation window is constructed based on the central grid. The grid elevation values ​​of the eight neighboring cells are read sequentially and differential operations are performed with the central grid cell in each direction to obtain the elevation difference sequence in eight directions. The elevation difference sequences in eight directions are sorted, the direction with the largest elevation decrease is extracted as the candidate mainstream direction, and the corresponding neighborhood grid position is marked. When there are equal elevation differences in the eight neighboring directions, or when the elevation differences between the center grid and the neighboring grids are all close to zero, the current grid cell is marked as a flat cell. For flat cells, perform neighborhood expansion and slope aspect updates.

[0034] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. The watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river channel cross-section within the watershed varies dramatically, and the course is short and steep. The digital elevation model (DEM) is generated through UAV lidar scanning and point cloud registration fusion, with a spatial resolution of 0.5 m × 0.5 m.

[0035] For each raster cell in the DEM, first identify it as the center raster, and then construct a 3×3 neighborhood calculation window centered on it. The center of the window is the current raster, and the eight cells surrounding the window are the eight-neighborhood cells. Then, sequentially read the elevation values ​​of the eight-neighborhood cells. , and the center grid elevation value Performing difference operations yields an elevation difference sequence in eight directions: , where the direction number The D8 encoding rule is adopted, corresponding to the northeast, east, southeast, south, southwest, west, northwest, and north directions respectively.

[0036] The elevation difference sequences in eight directions are sorted, and the direction with the largest elevation decrease is selected as the candidate mainstream direction. The corresponding neighboring grid positions are marked as the initial flow direction for the subsequent flow direction matrix. If multiple directions have equal elevation differences or the elevation difference between the central grid and its neighboring grids is close to zero (less than a preset threshold of 0.01m), the central grid is marked as a flat cell to distinguish it from cells with a clear slope aspect. For flat cells, a neighborhood expansion slope aspect update is implemented: the nearest non-flat cell is searched from the eight neighbors of the flat cell, and its flow direction information is obtained along the shortest spatial distance direction; this flow direction information is assigned to the flat cell, and the distance weight is recorded simultaneously for flow allocation in subsequent flow accumulation calculations; if all neighbors around the flat cell are flat cells, the process is expanded to a larger neighborhood until a slope aspect information that can be assigned is found, ensuring that each grid cell has a definite flow direction and avoiding flow direction breaks in flat areas of the DEM.

[0037] Ultimately, each grid cell obtains information on the main flow direction, marking whether it is a flat cell and the slope direction source of its neighborhood, forming a complete grid flow direction dataset.

[0038] Preferably, the neighborhood expansion slope aspect update for flat cells specifically involves: For flat cells, the calculation window is expanded step by step with the current grid as the center to form a 5×5 or 7×7 extended neighborhood, and the elevation difference is recalculated to obtain the set of slope gradient vectors. The slope gradient vectors corresponding to the eight neighboring directions are classified and aggregated to form directional gradient grouped data. The slope gradients within each directional group are accumulated item by item, and a distance attenuation coefficient is introduced to weight the gradients of distant grid cells in a decreasing manner. By comparing and analyzing the cumulative slope gradients in each direction, the direction with the largest cumulative slope value is selected as the dominant water flow direction. The determined dominant water flow direction is converted into a standard direction code and written into the corresponding grid position of the flow direction matrix to complete the flow direction update process.

[0039] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. The watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river cross-section within the watershed varies dramatically, and the course is short and steep. For the flat cells identified in the digital elevation model (DEM), a neighborhood expansion slope aspect update method is employed to ensure that each grid cell obtains a definite main flow direction.

[0040] Specifically, for each central grid cell marked as a flat unit, the calculation window range is first expanded step by step with it as the center to form a 5×5 or 7×7 extended neighborhood. The elevation values ​​of the grid cells in the extended neighborhood are obtained in turn, and the elevation difference of each grid cell relative to the central grid cell is calculated to form a set of slope gradient vectors.

[0041] Subsequently, the slope gradient vectors corresponding to the eight directions in the extended neighborhood are categorized and aggregated to obtain directional gradient grouping data for eight directions. The slope gradients within each directional group are then summed item by item, while a distance decay coefficient is introduced. We apply a weighted decrease to the gradient of distant grid cells to reflect the greater influence of the neighboring slope on the mainstream direction of the central grid cell.

[0042] After accumulation, the slope gradient values ​​in eight directions are compared and analyzed, and the direction with the largest accumulated slope value is selected as the dominant flow direction of the central grid. This direction is then converted into a standard direction code (such as D8 code) and written into the corresponding grid position in the flow direction matrix, thereby completing the flow direction update process for the flat cell.

[0043] Preferably, the confluence accumulation operation is performed based on the flow direction matrix, and the confluence intensity distribution is generated by superimposing the contributions of upstream units in parallel, and the watershed river network structure is extracted, including: Read the flow direction matrix and construct the directed connection relationship of the grid cells. By traversing the flow direction of each grid cell, count the number of upstream cells pointing to the same grid cell and generate the corresponding upstream association count. The grid cells are classified according to the upstream correlation count, and the grid cells with an upstream correlation count of zero are selected as the starting cells and assigned an initial flow contribution. The process proceeds step by step from upstream to downstream, transferring the flow contribution of each grid cell to the downstream cell along the flow path, and superimposing the contributions from different upstream directions. During the transmission process, the upstream association count of the corresponding downstream unit is reduced synchronously. When the upstream association count decreases to zero, the grid unit is included in the next round of transmission sequence. After completing the global transfer, the cumulative results of the flow of each grid cell are obtained to form the flow intensity distribution, and the watershed river network structure is extracted based on the continuity of the flow intensity.

[0044] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. This watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river cross-section within the watershed varies dramatically, and the course is short and steep. For each grid cell obtained from the flow direction matrix, a confluence accumulation calculation method is used to generate the confluence intensity distribution and extract the watershed's river network structure.

[0045] Specifically, the previously generated flow direction matrix is ​​first read, and a directed connection relationship is constructed for each grid cell based on the matrix, meaning each grid cell points to its downstream grid cell. By traversing the flow direction matrix, the number of upstream grid cells pointing to the same grid cell is counted, generating a corresponding upstream association count. Grid cells with an upstream association count of zero are identified as watershed initiating units. Initial confluence contributions are assigned to these initiating units to initiate confluence propagation.

[0046] During the flow confluence process, the flow is advanced step-by-step from upstream to downstream. Each grid cell transmits its own flow contribution along its flow path to the downstream cell. Simultaneously, flow from different upstream directions is superimposed to achieve multi-source flow accumulation. During the flow, the upstream association count of downstream cells is reduced in real time. When the upstream association count of a downstream cell decreases to zero, it is included in the next round of flow confluence, ensuring that no cell is missed in the step-by-step flow confluence process from upstream to downstream.

[0047] After the entire region is transferred, the cumulative flow results of each grid cell form a continuous flow intensity distribution. Based on the spatial continuity of flow intensity and a threshold filtering method, the main and tributary channels and confluence channels can be identified, thereby extracting the complete watershed river network structure, including the river initiation point, river branches, and the connectivity of the confluence network.

[0048] Preferably, the extraction of the watershed river network structure based on the continuity of confluence intensity includes: The current intensities of each grid cell are sorted based on the current intensities, and regions with continuous high values ​​are identified. Based on the connectivity of each grid cell, continuous high-flow-intensity cells are aggregated to form candidate river network segments; Candidate river network segments are extended and connected along the direction of increasing confluence intensity to construct a continuous river network path; The connectivity of the river network paths is checked, isolated units are eliminated, and the watershed river network structure is generated.

[0049] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. This watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river channel cross-section varies dramatically, and the course is short and steep. Based on the aforementioned distribution of flow intensity in the grid cells obtained through flow accumulation, this embodiment performs river network extraction.

[0050] First, the runoff intensities of all grid cells are sorted, and continuous high-value areas are identified from high to low. The specific method is as follows: the cumulative runoff value of each grid cell is read, and grid cells with runoff intensities exceeding the threshold are selected as candidate river cells by setting a threshold. At the same time, the threshold can be adaptively adjusted according to changes in sub-basin area or slope to take into account both main channels and secondary tributaries.

[0051] Subsequently, based on the spatial connectivity of the grid cells, adjacent continuous high-convergence-intensity cells are aggregated to form candidate river network segments. The aggregation process employs an eight-neighborhood connection rule to ensure that diagonally connected cells can also be correctly assigned to the same river segment; discontinuous and isolated high-intensity cells are temporarily retained, with subsequent connectivity checks determining whether to remove them.

[0052] Next, the candidate river network segments are extended and connected along the direction of increasing confluence intensity. Specifically, starting from the upstream high-convergence unit, adjacent downstream high-intensity units are added sequentially in the direction of the flow direction matrix, connecting the candidate segments into a continuous river network path. During the connection process, considering slope changes and flow direction consistency, units that deviate too much from the flow direction are removed or redistributed to ensure the continuity of the river network path and the true flow direction.

[0053] Finally, the connectivity of the constructed river network path is checked. This check includes removing isolated points, short tributary segments, and non-mainstream units, while merging tributary segments with the main channel to form a complete and continuous watershed river network structure. The generated river network path includes both the main channel and retains information on important tributaries.

[0054] As an example of the present invention, reference is made to... Figure 2 As shown, step S4 in this example includes: Step S41: Extract the slope information and corresponding water depth and flow velocity information of steep slope units to calculate the local propagation velocity of each unit; Step S42: Determine the time step constraint based on the local propagation speed and grid scale; Step S43: Based on the time step constraint, identify steep slope elements and mark elements that do not meet the time step constraint as elements that need to be refined. Step S44: Divide the calculation unit to be refined into multiple sub-time steps, and perform multiple recursive update calculations along the water flow direction within the same global time step to generate the flow generation process of the steep slope unit after iteration.

[0055] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. This watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river channel cross-section varies dramatically, and the course is short and steep. For the flood runoff process in steep-slope units within the watershed, this embodiment employs a time-step constraint calculation method combining local propagation velocity and grid scale.

[0056] First, the slope information, along with the corresponding water depth and flow velocity information, is extracted from steep slope units to calculate the local propagation velocity of each unit. Specifically, the local slope of each steep slope unit is obtained based on the digital elevation model and watershed division results; water depth and flow velocity information are obtained using the previous calculation step or earlier simulation results; and the wave velocity within the unit is calculated using the Saint-Venant kinematic wave theory or local shallow water dynamic equations to obtain the local propagation velocity.

[0057] Subsequently, based on the local propagation velocity of each element and the spatial size of the grid elements, time step constraints were determined. Specifically, according to the CFL conditions (Courant-Friedrich-Levi conditions), the time step was restricted to no more than the ratio of the element length to the local propagation velocity, thereby ensuring the stability and accuracy of the hydrodynamic numerical calculations.

[0058] Next, steep slope elements are identified based on the time step constraint, and elements that do not meet the global time step condition are marked as elements requiring further refinement. These elements are then divided into multiple sub-time steps, with each global time step further divided into several sub-time steps, to ensure that the local water flow propagation process remains numerically stable and accurate.

[0059] Finally, multiple recursive update calculations are performed on the computationally refined cells along the flow direction. Within each sub-time step, the water depth and flow velocity of the computational cell are updated, and the calculation results are passed to the downstream cells. After the iteration is completed, the flow generation process of the steep slope cell throughout the entire global time step is generated.

[0060] Preferably, the formula for the time step constraint is as follows: In the formula, For time step, As an empirical safety factor, For local slope, For the speed of water propagation, For spatial grid scale, It is the acceleration due to gravity. The water is deep.

[0061] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. This watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river cross-section varies dramatically, and the course is short and steep. For flood runoff calculations of steep-slope units, this embodiment employs a method combining local propagation velocity with grid-scale time step constraints to ensure the stability and accuracy of the numerical calculations.

[0062] First, the slope information and corresponding water depth and flow velocity information of steep slope units are extracted to calculate the local propagation velocity. Based on the digital elevation model and hydrodynamic calculation results, the local slope of each unit is obtained. , water depth and water flow velocity and combined with grid space scale The formula for establishing the time step constraint is as follows: In the formula, For time step, As an empirical safety factor, For local slope, For the speed of water propagation, For spatial grid scale, It is the acceleration due to gravity. The water depth is given. This formula comprehensively considers the influence of hydrodynamic propagation speed and topographic slope on numerical stability, and ensures that the CFL condition is satisfied by selecting a smaller time step to limit the update rate of local cells.

[0063] Next, the calculated Compared with the global time step, steep slope units that do not meet the conditions are divided into units requiring further calculation. Within the same global time step, these units are further divided into multiple sub-time steps, and recursive updates are performed along the flow direction to gradually generate the local runoff process. This ensures the numerical accuracy of rapid flood propagation while also considering computational efficiency.

[0064] Of particular importance is the generation of non-equidistant flight strips along the main confluence direction based on slope grid calculation results, and the writing of the start and end coordinates of the flight strips into the flight control system, including: The main confluence direction field of the watershed is constructed based on the slope aspect grid calculation results; The starting point of the flight strip is divided along the main confluence direction, and the flight strip spacing is dynamically adjusted according to the slope of each grid and the local confluence potential, so that the flight strips are dense in steep slopes and areas with concentrated confluence, while the flight strip spacing is increased in gentle slopes and areas with sparse confluence, thereby generating dynamic flight strip spacing data. The coordinate sequence of each flight path is generated by using dynamic flight path spacing data, and then converted into flight control commands and written into the UAV flight control system.

[0065] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. The watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river cross-section varies dramatically, and the course is short and steep. The watershed contains areas of concentrated steep-slope confluence, gentle-slope areas, and low-lying riverbeds. To obtain high-precision topographic point cloud data, an unmanned aerial vehicle (UAV) platform equipped with a lidar scanning component deploys non-equidistant flight paths along the main confluence direction of the watershed for full-coverage scanning.

[0066] First, slope aspect raster data is calculated based on the watershed digital elevation model to generate the main confluence direction field of the watershed. By performing slope aspect analysis on the elevation value of each raster cell and the elevation difference with neighboring raster cells, the local dominant confluence direction is determined. This is then combined with 3×3 or 5×5 neighborhood weighted average smoothing to form a continuous main confluence direction field. This direction field not only reflects the overall confluence direction of the water flow but also shows the flow direction changes on steep slopes and in areas of concentrated confluence, providing a basis for navigation strip layout.

[0067] Secondly, the starting point of the flight strip is generated along the main confluence direction, and the spacing of the flight strips is dynamically adjusted according to the local slope and confluence potential. Specifically, a fixed spacing is initially defined for each flight strip, and the local slope θ and confluence potential P (e.g., a function based on slope and confluence accumulation) are calculated at each grid cell. ), reducing the spacing of flight strips in areas with steep slopes and concentrated runoff to The spacing between gentle slopes and sparsely populated runoff areas can be expanded to By progressively scanning the entire watershed, a dynamic spacing sequence for each flight strip is generated, enabling non-equidistant deployment.

[0068] Subsequently, the path coordinate sequence of each flight strip is converted into UAV flight control commands. First, the flight strip path coordinates are mapped to flight altitude, flight speed, and attitude angle control commands, taking into account UAV dynamics constraints and lidar scanning coverage requirements. Second, the flight commands for the first flight strip are set and uploaded to the flight control system, which then uses integrated navigation (inertial measurement unit + GNSS positioning) to maintain stable flight along the flight strip. At the end of the flight strip, the flight control system executes a slope-constrained turnaround maneuver according to a preset turnaround strategy, causing the UAV to switch to an adjacent flight strip while adjusting altitude and lateral offset to ensure that a preset overlapping area is formed between adjacent flight strips, avoiding blind spots caused by terrain obstruction.

[0069] During flight, the lidar scanning component continuously performs lateral scanning and acquires the UAV's attitude angle information in real time for scanning angle compensation, ensuring the accuracy and spatial consistency of the point cloud data. It repeatedly executes flight stripe flight, turnaround switching, and overlap control until the entire watershed is covered. After the flight mission is completed, the system automatically outputs continuous point cloud data, including the elevation information, reflection intensity, and scan timestamp of each grid cell, providing foundational data for subsequent digital elevation model generation and watershed hydrological analysis.

[0070] Most importantly, the dynamic adjustment of the flight strip spacing based on the slope of each grid and the local confluence potential is specifically as follows: The cumulative upstream flow rate of each grid cell is calculated based on the flow direction matrix to obtain the preliminary flow intensity distribution; By using the slope data of each grid cell, the cumulative runoff intensity is weighted to obtain the local runoff potential value of each grid cell. The grid cells with larger slopes or higher upstream cumulative flow have greater potential. The watershed is divided into high potential areas and low potential areas, and the spacing between the flight strips is automatically adjusted according to the local confluence potential, so that the flight strips are dense in the high potential areas and sparse in the low potential areas, thus obtaining dynamic flight strip spacing data.

[0071] In one embodiment, a steep-slope, short-course watershed on a tropical island is used as the target study area. The watershed has a total area of ​​approximately 23 km², an average slope of about 16°, and a maximum slope of up to 38°. The river channel cross-section varies dramatically, and the course is short and steep. To achieve high-precision terrain scanning, an unmanned aerial vehicle (UAV) platform equipped with a lidar scanning module deploys non-equidistant flight strips along the main confluence direction, and dynamically adjusts the flight strip spacing to cover the entire watershed.

[0072] First, slope aspect raster data is calculated based on the digital elevation model to generate the main runoff direction field of the watershed. The cumulative upstream runoff volume is calculated for each raster cell using the runoff direction matrix to obtain a preliminary runoff intensity distribution. This distribution is then weighted using slope information to generate the local runoff potential value for each raster cell. Raster cells with steeper slopes or higher upstream cumulative flow have greater potential. Subsequently, the watershed is divided into high-potential and low-potential zones: high-potential zones are typically located on steep slopes and in areas of concentrated runoff, while low-potential zones are typically located on gentle slopes or in areas of sparse runoff.

[0073] During the deployment of flight strips, the spacing between the flight strips is dynamically adjusted based on the potential for convergence of currents within the grid pattern. Specifically, the spacing between flight strips in high-potential areas is reduced to a preset minimum value. To ensure scanning accuracy and point cloud density in steep slopes and areas of concentrated runoff; the spacing between flight strips in low potential areas is increased to the preset maximum value. This improves scanning efficiency and reduces redundant data. By progressively scanning the entire watershed, a dynamic spacing sequence for each flight strip is generated, forming complete dynamic flight strip spacing data.

[0074] After generating the flight path coordinate sequence, it is converted into UAV flight control commands. The flight control system generates attitude and speed control commands based on the flight path coordinates, flight altitude, and scan parameters to ensure stable flight of the UAV along the first flight path. When the UAV reaches the end of the flight path, the flight control system performs a slope-constrained turnaround maneuver, switching the flight direction to an adjacent flight path and adjusting the altitude and lateral offset to form a preset scan overlap area. This process of flight path switching and turnaround is repeated until the entire watershed scan is completed.

[0075] During flight, the lidar scanning component continuously performs lateral scans and acquires the UAV's attitude angle in real time for scanning angle compensation, ensuring the accuracy and spatial consistency of the point cloud data. The final output continuous point cloud data contains the elevation information and reflection intensity of each grid cell, which can be used for digital elevation model generation, watershed slope analysis, and flood evolution simulation.

[0076] Therefore, the embodiments should be considered as exemplary and non-limiting in all respects, and the scope of the invention is defined by the appended claims rather than the foregoing description. Thus, all variations falling within the meaning and scope of the equivalents of the application are intended to be included within the invention.

[0077] The above description is merely a specific embodiment of the present invention, enabling those skilled in the art to understand or implement the invention. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of the invention. Therefore, the present invention is not to be limited to the embodiments shown herein, but is to be accorded the widest scope consistent with the principles and novel features of the invention herein.

Claims

1. A flood evolution simulation method for steep slope-short flow basins in tropical islands, characterized in that, Includes the following steps: Step S1: Calculate the slope and flow length based on the obtained watershed river network and sub-watershed boundaries to obtain the steep slope index; use the steep slope index to divide the watershed into steep slope units, fast-flowing river units, and downstream slow-flowing or floodplain units. Step S2: Establish a variable-dimensional hydrodynamic model based on different units. The steep slope unit adopts the kinematic wave model, the rapid river unit adopts the one-dimensional Saint-Venant equation, and the downstream slow-flow or floodplain unit adopts the two-dimensional local inertial shallow water equation. Step S3: Construct distributed rainfall input and calculate the runoff generation process of each unit; Step S4: Use the slope information of steep slope units to perform Coulomb number constraints, and perform local sub-time step iterative calculations on steep slope units based on time step constraints to obtain the runoff generation process of steep slope units after iteration. Step S5: Solve the runoff generation process of each unit to obtain the peak arrival time, inundation range and water depth distribution results.

2. The flood evolution simulation method for steep slope-short flow basins in tropical islands according to claim 1, characterized in that, Methods for obtaining the river network and sub-basin boundaries of a watershed include: The control drone platform is equipped with a lidar scanning component to perform round-trip route scanning of the watershed, continuously acquire surface point cloud data, and cache it through the onboard storage unit; A continuous digital elevation model is generated by registration and fusion of surface point cloud data. The elevation difference operation is performed on the digital elevation model grid by grid to output the water flow direction information of each unit and form a flow direction matrix. The flow accumulation operation is performed based on the flow direction matrix. The flow intensity distribution is generated by superimposing the contributions of upstream units in parallel, and the river network structure of the watershed is extracted. Based on the river network structure of the basin, the location of the flow direction boundary is determined and the region is divided to obtain the sub-basin boundary.

3. The flood evolution simulation method for steep slope-short flow basins in tropical islands according to claim 2, characterized in that, The control of the unmanned aerial vehicle platform, equipped with a lidar scanning component, to perform round-trip route scanning of the watershed includes: Import watershed boundary data into the ground control terminal, generate non-equidistant flight strips along the main confluence direction based on the slope aspect grid calculation results, and write the start and end coordinates of the flight strips into the flight control system; After receiving the flight path coordinates, the UAV flight control system drives the power unit to take off and uses the inertial measurement unit and satellite positioning module for combined navigation to enable the UAV to fly stably along the first flight path; During flight, the lidar scanning component is controlled to perform continuous lateral sweeps and simultaneously acquire attitude angle data to compensate for the scanning angle in real time. When the flight control system identifies the coordinates of the current flight path's end point, it controls the UAV to perform a turnaround maneuver based on slope constraints, switching the flight direction to the adjacent flight path; During the flight of adjacent flight strips, the flight altitude and lateral offset are adjusted to make the adjacent scan strips form a preset overlapping area; Repeatedly execute the flight strip and turnaround switching until the watershed is fully covered and continuous point cloud data is output.

4. The flood evolution simulation method for steep slope-short flow basins in tropical islands according to claim 2, characterized in that, The digital elevation model is subjected to grid-by-grid elevation difference calculations, and the flow direction information of each cell is output to form a flow direction matrix, including: Load the digital elevation model; The terrain is divided into regular grid cells with a fixed spatial resolution, and a grid neighborhood index relationship is established. For each grid cell and its eight neighboring cells, the elevation difference is calculated to obtain the direction of the minimum elevation neighboring cell, and the direction of the minimum elevation neighboring cell is used as the main water flow direction of the grid cell. The remaining grid cells are arranged in a matrix based on the direction of the main water flow and the spatial position of the grid cell to form a flow direction matrix.

5. The flood evolution simulation method for steep slope-short flow basins in tropical islands according to claim 4, characterized in that, The calculation of the elevation difference between each grid cell and its eight neighboring cells also includes: Identify the center grid based on the current grid cell; A 3×3 neighborhood calculation window is constructed based on the central grid. The grid elevation values ​​of the eight neighboring cells are read sequentially and differential operations are performed with the central grid cell in each direction to obtain the elevation difference sequence in eight directions. The elevation difference sequences in eight directions are sorted, the direction with the largest elevation decrease is extracted as the candidate mainstream direction, and the corresponding neighborhood grid position is marked. When there are equal elevation differences in the eight neighboring directions, or when the elevation differences between the center grid and the neighboring grids are all close to zero, the current grid cell is marked as a flat cell. For flat cells, perform neighborhood expansion and slope aspect updates.

6. The flood evolution simulation method for steep slope-short flow basins in tropical islands according to claim 5, characterized in that, The specific steps for neighborhood expansion and slope aspect update for flat cells are as follows: For flat cells, the calculation window is expanded step by step with the current grid as the center to form a 5×5 or 7×7 extended neighborhood, and the elevation difference is recalculated to obtain the set of slope gradient vectors. The slope gradient vectors corresponding to the eight neighboring directions are classified and aggregated to form directional gradient grouped data. The slope gradients within each directional group are accumulated item by item, and a distance attenuation coefficient is introduced to weight the gradients of distant grid cells in a decreasing manner. By comparing and analyzing the cumulative slope gradients in each direction, the direction with the largest cumulative slope value is selected as the dominant water flow direction. The determined dominant water flow direction is converted into a standard direction code and written into the corresponding grid position of the flow direction matrix to complete the flow direction update process.

7. The flood evolution simulation method for steep slope-short flow basins in tropical islands according to claim 2, characterized in that, The confluence accumulation operation is performed based on the flow direction matrix. The confluence intensity distribution is generated by superimposing the contributions of upstream units in parallel, and the watershed river network structure is extracted, including: Read the flow direction matrix and construct the directed connection relationship of the grid cells. By traversing the flow direction of each grid cell, count the number of upstream cells pointing to the same grid cell and generate the corresponding upstream association count. The grid cells are classified according to the upstream correlation count, and the grid cells with an upstream correlation count of zero are selected as the starting cells and assigned an initial flow contribution. The process proceeds step by step from upstream to downstream, transferring the flow contribution of each grid cell to the downstream cell along the flow path, and superimposing the contributions from different upstream directions. During the transmission process, the upstream association count of the corresponding downstream unit is reduced synchronously. When the upstream association count decreases to zero, the grid unit is included in the next round of transmission sequence. After completing the global transfer, the cumulative results of the flow of each grid cell are obtained to form the flow intensity distribution, and the watershed river network structure is extracted based on the continuity of the flow intensity.

8. The flood evolution simulation method for steep slope-short flow basins in tropical islands according to claim 7, characterized in that, The watershed river network structure extracted based on the continuity of confluence intensity includes: The current intensities of each grid cell are sorted based on the current intensities, and regions with continuous high values ​​are identified. Based on the connectivity of each grid cell, continuous high-flow-intensity cells are aggregated to form candidate river network segments; Candidate river network segments are extended and connected along the direction of increasing confluence intensity to construct a continuous river network path; The connectivity of the river network paths is checked, isolated units are eliminated, and the watershed river network structure is generated.

9. The flood evolution simulation method for steep slope-short flow basins in tropical islands according to claim 1, characterized in that, Step S4 includes the following steps: Step S41: Extract the slope information and corresponding water depth and flow velocity information of steep slope units to calculate the local propagation velocity of each unit; Step S42: Determine the time step constraint based on the local propagation speed and grid scale; Step S43: Based on the time step constraint, identify steep slope elements and mark elements that do not meet the time step constraint as elements that need to be refined. Step S44: Divide the calculation unit to be refined into multiple sub-time steps, and perform multiple recursive update calculations along the water flow direction within the same global time step to generate the flow generation process of the steep slope unit after iteration.

10. The flood evolution simulation method for steep slope-short flow basins in tropical islands according to claim 1, characterized in that, The formula for the time step constraint is shown below: In the formula, For time step, As an empirical safety factor, For local slope, For the speed of water propagation, For spatial grid scale, It is the acceleration due to gravity. The water is deep.

Images