Hydrodynamic simulation pavement drainage structure optimization method and system
By constructing a refined three-dimensional model and combining two-dimensional shallow water wave equation simulation with a neural network proxy model, and by combining field test calibration and adaptive genetic algorithm optimization, the problems of refinement and low efficiency in road drainage design in existing technologies have been solved, and the scientific design and reliable prediction of drainage structures have been realized.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- ZHEJIANG JIAOTOU EXPRESSWAY CONSTR MANAGEMENT CO LTD
- Filing Date
- 2026-03-19
- Publication Date
- 2026-05-19
AI Technical Summary
Existing road drainage designs lack detailed 3D models, resulting in low simulation efficiency and insufficient calibration of simulation model parameters. This makes it difficult to accurately describe the complex hydrodynamic processes of road surface runoff, leading to a lack of targeted drainage structure design and an inability to effectively cope with varying rainfall conditions.
A refined three-dimensional digital model of the road surface was constructed using UAV oblique photography and point cloud data processing. Combined with two-dimensional shallow water wave equation simulation and neural network proxy model, and calibrated by field test, multi-objective optimization was carried out through adaptive constraint penalty genetic algorithm to achieve scientific design of drainage structure parameters.
It enables rapid and accurate assessment of drainage performance, improves the reliability of simulation predictions and the effectiveness of engineering implementation, and obtains drainage structure parameter schemes that balance performance and cost through online monitoring system verification.
Smart Images

Figure CN122065411A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of intelligent road drainage technology, and more specifically, to a method and system for optimizing road drainage structure using fluid dynamics simulation. Background Technology
[0002] With the increasing frequency of extreme rainfall events, waterlogging and flooding problems in road engineering are becoming more and more prominent, seriously affecting traffic safety and road operation efficiency. Traditional road drainage design mainly relies on empirical formulas and standard parameters, which are difficult to accurately reflect the complex hydrodynamic processes of road surface runoff. This results in drainage structure design lacking specificity and failing to effectively cope with variable rainfall conditions.
[0003] Existing road drainage optimization methods have the following shortcomings: they lack detailed three-dimensional road surface models, making it difficult to accurately describe the impact of road surface micro-topography on runoff collection; simulation calculation efficiency is low, making it difficult to support large-scale parameter optimization; traditional CFD simulation is time-consuming, making it difficult to meet engineering optimization needs; simulation model parameters rely on empirical values, lack field test calibration, making it difficult to guarantee prediction accuracy, and limiting the reliability and applicability of the model.
[0004] Therefore, there is an urgent need for a road drainage structure optimization method that integrates refined 3D modeling, efficient fluid dynamics simulation, intelligent optimization algorithms, and on-site test calibration, in order to achieve the scientific design and precise modification of drainage structure parameters and improve the drainage capacity and flood control level of road engineering. Summary of the Invention
[0005] This invention provides a method and system for optimizing road drainage structures using fluid dynamics simulation, which solves the technical problems in related technologies such as lack of refined simulation support for road drainage structure design, low optimization efficiency, and insufficient calibration of simulation model parameters.
[0006] This invention provides a method for optimizing road drainage structures using fluid dynamics simulation, comprising the following steps:
[0007] S1. Collect geographic information of the target road section and use UAV oblique photography and point cloud data processing methods to obtain a refined three-dimensional digital model of the road surface.
[0008] S2, acquire historical rainfall data, and use frequency analysis and design rainfall pattern generation methods to obtain a set of rainfall boundary conditions;
[0009] S3. Based on the refined three-dimensional digital model of the road surface and the set of rainfall boundary conditions, the drainage performance index is obtained by using the two-dimensional shallow water wave equation simulation method.
[0010] S4. Based on the drainage performance indicators, a drainage structure parameter scheme is obtained by using a neural network surrogate model and an improved genetic algorithm optimization method.
[0011] S5. Based on the drainage structure parameter scheme, the parameters of the calibrated simulation model are obtained by using field tests and parameter calibration methods.
[0012] S6. Based on the calibrated simulation model parameters and drainage structure parameter scheme, engineering implementation and online monitoring methods are used to obtain the engineering effect verification results.
[0013] In a preferred embodiment, S1 includes:
[0014] Based on a high-resolution oblique image dataset, an incremental motion recovery structure algorithm is used to recover camera extrinsic parameters and sparse 3D point clouds. Combined with ground control points, a bundle adjustment algorithm is used for overall optimization.
[0015] A semi-global matching algorithm is used to generate a high-density point cloud through stereo matching. Statistical filtering is then applied to the point cloud to remove noise points and outliers, resulting in a high-density three-dimensional point cloud model of the road surface.
[0016] A threshold segmentation method based on elevation and intensity features is used to identify different road elements in the point cloud. By setting elevation and density thresholds, the point cloud of the road surface is automatically extracted, and the location of rainwater inlets and road boundary lines are identified.
[0017] The road is divided into multiple segments along its longitudinal direction. For each segment, the least squares method is used to fit the plane to calculate the cross slope and longitudinal slope, resulting in a parametric pavement geometric model that includes slope parameters and roughness coefficients.
[0018] In a preferred embodiment, S2 includes:
[0019] The maximum annual rainfall was obtained by extracting the maximum rainfall over different durations using the sliding time window method. The coefficient of variation and skewness coefficient were calculated by fitting the data using the Pearson type III distribution. The formula for rainstorm intensity was fitted using the least squares method.
[0020] The Chicago rainfall pattern method was used to design the rainfall process line. The peak rainfall location coefficient was set, and the rainfall intensity sequence was calculated step by step according to the set time step. The sequence was normalized and adjusted so that the total rainfall was equal to the design value, and the time history rainfall intensity process line was obtained.
[0021] Runoff volume is calculated using a reasoning formula. Runoff volume equals the comprehensive runoff coefficient multiplied by the rainfall intensity and then multiplied by the catchment area. The inflow volume is converted according to the time history variation of the rainfall process line to obtain the time history inflow process line at each boundary location.
[0022] In a preferred embodiment, S3 includes:
[0023] A two-dimensional shallow water wave equation based on the laws of mass conservation and momentum conservation is established. The source and sink terms include rainfall source terms, pavement infiltration terms, and storm drain discharge terms. The storm drain discharge rate is calculated using the weir flow formula or the orifice flow formula, and the friction source term is described using the Manning formula.
[0024] Unstructured triangular meshes are generated using the constrained Delaunay triangulation algorithm. Locally refined regions are set at the locations of rainwater inlets, and mesh quality is optimized through edge flipping and node smoothing operations.
[0025] The finite volume method is used to spatially discretize the governing equations, the interface numerical flux is calculated using the HLL approximate Riemann solver, the time step is driven by the second-order Runge-Kutta method, and the time step is adaptively determined according to the CFL condition.
[0026] Set a water accumulation threshold, count the number of grid cells with water depth exceeding the threshold to obtain the water accumulation area, count the maximum water depth value to obtain the maximum water accumulation depth, and count the duration of water accumulation.
[0027] In a preferred embodiment, S4 includes:
[0028] A neural network proxy model is established. The network input is a vector of design parameters, including the road cross slope, longitudinal slope, road crown height, and the spacing and size of the storm drains. The network output is a vector of performance indicators, including the maximum water accumulation depth, the maximum water accumulation area, and the drainage time.
[0029] A physical constraint loss function is constructed. The total loss function is a weighted combination of data fitting loss and physical constraint loss. The physical constraint loss includes three constraints: water depth is negatively correlated with cross slope, water area is positively correlated with distance between rainwater inlets, and drainage time is positively correlated with rainfall intensity. The partial derivatives of the network prediction results with respect to the parameters are calculated to determine whether the monotonicity conforms to physical laws. If the laws are violated, a penalty term is added.
[0030] An incremental sampling strategy is adopted. After training is completed, new sampling points are selected in the parameter space for CFD simulation calculation, and the new samples are added to the training set to retrain the network.
[0031] In a preferred embodiment, S4 further includes:
[0032] A multi-objective optimization mathematical model is established, with objective 1 being to minimize the average maximum water accumulation depth under multiple working conditions and objective 2 being to minimize the engineering modification cost.
[0033] Design an adaptive constraint penalty function to calculate the constraint violation amount for each constraint condition. The total constraint violation amount is equal to the sum of all constraint violations. After each generation of evolution, adjust the penalty coefficient based on the comparison results between the proportion of feasible solutions in the population and the preset first and second thresholds.
[0034] The non-dominated sorting genetic algorithm is executed. A tournament selection strategy is used to select parent individuals. A simulated binary crossover operator is used for crossover operations, and a polynomial mutation operator is used for mutation operations. The merged population is non-dominated sorted and the crowding distance is calculated. Individuals are selected to enter the next generation according to the principle of hierarchical priority and second-best crowding distance.
[0035] In a preferred embodiment, S5 includes:
[0036] The process of road water accumulation is divided into three stages: water accumulation rise stage, peak water depth stage and water accumulation recede stage. Three characteristic indicators are extracted: water accumulation rise rate, peak water depth and water accumulation recede rate.
[0037] Sensitivity analysis of parameters showed that the rate of water accumulation rise is sensitive to the roughness coefficient but not to the discharge coefficient, while the rate of water accumulation recedes is sensitive to the discharge coefficient but not to the roughness coefficient.
[0038] The roughness coefficient is identified by the water rise rate index. A lookup table of roughness coefficient and rise rate is established. The measured rise rate is substituted into the lookup table and linear interpolation is used to obtain the roughness coefficient.
[0039] Using a fixed roughness coefficient as the calibration value, the discharge coefficient is identified by the water receding rate index. A lookup table of discharge coefficient and receding rate is established. The discharge coefficient is obtained by substituting the experimentally measured receding rate into the lookup table and using linear interpolation.
[0040] In a preferred embodiment, S6 includes:
[0041] Each monitoring point is equipped with an ultrasonic water level gauge. The water level gauge probe is pointed downwards at the water line on the road surface. Recording will start automatically when the water depth exceeds the preset recording threshold and stop after the water depth remains below the preset recording threshold for a preset time.
[0042] Each monitoring point is equipped with video surveillance equipment, set to record event-triggered recording. When the water level gauge detects that the water depth exceeds the preset recording threshold, the camera will automatically start recording.
[0043] The monitoring equipment transmits data to the monitoring center server via a wireless network. The system sets an early warning threshold, and automatically sends an early warning message when the water depth exceeds the threshold.
[0044] In a preferred embodiment, S6 further includes:
[0045] Extract the maximum water depth and duration of water accumulation during rainfall events and compare them with the values under similar historical rainfall conditions before the renovation.
[0046] Social benefit assessment was conducted based on traffic flow monitoring data, which showed average vehicle speed and travel time during rainy weather.
[0047] Economic benefit assessment is conducted by calculating annual cost savings based on the average daily traffic volume of roads, the proportion of rainy days, the time saved during each rainfall event, and the time value of vehicles.
[0048] In a preferred embodiment, a pavement drainage structure optimization system based on fluid dynamics simulation is used to perform the steps of a pavement drainage structure optimization method based on fluid dynamics simulation, including:
[0049] The 3D modeling module is used to collect geographic information of the target road section and uses UAV oblique photography and point cloud data processing methods to obtain a refined 3D digital model of the road surface.
[0050] The boundary condition generation module is used to acquire historical rainfall data and use frequency analysis and design rainfall pattern generation methods to obtain a set of rainfall boundary conditions;
[0051] The simulation calculation module is used to obtain drainage performance indicators based on the refined three-dimensional digital model of the road surface and the set of rainfall boundary conditions, using the two-dimensional shallow water wave equation simulation method.
[0052] The optimization solution module is used to obtain drainage structure parameter schemes based on drainage performance indicators, using a neural network surrogate model and an improved genetic algorithm optimization method.
[0053] The parameter calibration module is used to obtain the calibrated simulation model parameters based on the drainage structure parameter scheme, using field tests and parameter calibration methods.
[0054] The effect verification module is used to obtain the engineering effect verification results based on the calibrated simulation model parameters and drainage structure parameter scheme, using engineering implementation and online monitoring methods.
[0055] The beneficial effects of this invention are as follows:
[0056] By constructing a refined three-dimensional digital model of the road surface through UAV oblique photography, and combining two-dimensional shallow water wave equation simulation and physical constraint neural network proxy model, a rapid and accurate evaluation of drainage performance was achieved. An adaptive constraint penalty genetic algorithm was used for multi-objective optimization, which can obtain a drainage structure parameter scheme that balances performance and cost under the condition of meeting engineering constraints.
[0057] By using a segmented parameter identification method based on water depth response characteristics, the sensitivity differences of roughness coefficient and discharge coefficient at different stages of the water accumulation process are utilized to achieve parameter decoupling identification. The model is calibrated by combining on-site artificial rainfall test data, which improves the reliability of simulation prediction. The implementation effect of the project is verified by an online monitoring system. Attached Figure Description
[0058] Figure 1This is a flowchart of a method for optimizing road drainage structure using fluid dynamics simulation according to the present invention;
[0059] Figure 2 This is a flowchart of a method for optimizing road drainage structure based on fluid dynamics simulation according to the present invention. Detailed Implementation
[0060] The subject matter described herein will now be discussed with reference to exemplary embodiments. It should be understood that these embodiments are discussed only to enable those skilled in the art to better understand and implement the subject matter described herein, and changes may be made to the function and arrangement of the elements discussed without departing from the scope of this specification. Various processes or components may be omitted, substituted, or added as needed in the examples. Furthermore, some features described in the examples may be combined in other examples.
[0061] At least one embodiment of the present invention discloses a method for optimizing road drainage structures based on fluid dynamics simulation, such as... Figures 1 to 2 As shown, it includes the following steps:
[0062] S1. Collect geographic information of the target road section and use UAV oblique photography and point cloud data processing methods to obtain a refined three-dimensional digital model of the road surface.
[0063] S11, Obtain geographic information of the target road segment by using UAV multi-view image acquisition method to obtain high-resolution oblique image dataset;
[0064] Aerial photogrammetry was conducted using a drone equipped with a high-resolution camera. Multiple parallel flight paths were planned, with appropriate spacing and flight altitude. The camera tilt angle was set to vertical and multiple tilt angles to achieve five-lens oblique photography. The image forward and lateral overlap were both set to be greater than 70%, ensuring that each ground point was covered by images from multiple different angles, achieving centimeter-level ground resolution. Several ground control points were established, and their 3D coordinates were determined using RTK surveying equipment, achieving centimeter-level accuracy. The drone then performed an automated flight mission along the planned flight paths, acquiring a large number of effective images to form a multi-view oblique image dataset covering the target road segment.
[0065] S12. Based on the high-resolution oblique image dataset, the motion recovery structure algorithm and the multi-view stereo dense matching method are used to obtain a high-density three-dimensional point cloud model of the road surface.
[0066] The acquired images were imported into photogrammetric software for processing. First, aerial triangulation was performed, and the software automatically extracted and matched image feature points. The incremental structure-of-motion (SfM) algorithm was then used to recover camera extrinsic parameters and sparse 3D point clouds. Combined with ground control points, a bundle adjustment algorithm was used for overall optimization, ensuring that the control point errors met accuracy requirements.
[0067] Multi-view stereo dense matching processing is performed. The software uses the semi-global matching (SGM) algorithm for stereo matching to generate a high-density point cloud with an average density of several thousand to tens of thousands of points per square meter. Statistical filtering is then applied to the point cloud to remove noise and outliers. This results in a high-quality, filtered 3D road surface point cloud model.
[0068] S13. Based on the high-density three-dimensional point cloud model of the road surface, a parametric road surface geometric model is obtained by using point cloud classification and geometric parameter extraction methods.
[0069] Point clouds are classified, and a threshold segmentation method based on elevation and intensity features is used to identify different road elements. The elevation values of the road surface point cloud are near the road design elevation, the elevation of the curb stone point cloud changes abruptly, and the point cloud density in the storm drain grate area is low. By setting elevation and density thresholds, the road surface point cloud is automatically extracted, and the locations of storm drain inlets and road boundaries are identified.
[0070] Geometric parameters are extracted from the road surface point cloud. The road is divided into multiple segments at certain intervals along its longitudinal direction. For each segment, a least-squares method is used to fit a plane, and the cross slope and longitudinal slope are calculated. The cross slope is the gradient of the fitted plane in the transverse direction, and the longitudinal slope is the gradient of the fitted plane in the longitudinal direction. The cross slope and longitudinal slope values of each segment are extracted, and the road crown height is calculated as the height difference between the road center and the curb. The initial value of the Manning roughness coefficient is determined based on the road surface material type. A parametric road surface geometric model containing slope parameters and roughness coefficient is obtained.
[0071] S14. Based on the parametric pavement geometry model and field survey data, a refined three-dimensional digital model of the pavement is obtained by using the drainage facility parameter input and model fusion method.
[0072] A field survey team was organized to conduct detailed measurements of the storm drains, recording parameters such as drain type, grate size, gap ratio, depth, and branch pipe diameter. The cross-sectional shape and dimensions of the drainage ditches were also recorded. GPS equipment was used to record the planar coordinates of each storm drain.
[0073] The geometric parameters of the storm drains and side ditches were entered into 3D modeling software to create 3D models of storm drain grates, collection wells, and collection ditches. The drainage facility model was spatially registered and fused with the road surface point cloud model, and topology checks were performed to correct gaps and overlaps. This resulted in a complete 3D model of the road drainage structure, including the road surface, storm drains, and collection ditches.
[0074] S2, acquire historical rainfall data, and use frequency analysis and design rainfall pattern generation methods to obtain a set of rainfall boundary conditions;
[0075] S21. Obtain long-term rainfall records and use the annual maximum value method and Pearson type III distribution fitting method to obtain the parameters of the rainstorm intensity formula;
[0076] Collect hourly rainfall data from long-term weather stations in the target road area, with a data period preferably greater than 30 years. Use a sliding time window method to extract the maximum rainfall over different durations, setting various rainfall durations such as 30 minutes, 60 minutes, and 120 minutes to obtain the annual maximum value sample.
[0077] Frequency analysis was performed on samples from each time period, and the sample mean and standard deviation were calculated. A Pearson Type III distribution was used as the theoretical distribution for fitting, and the coefficient of variation Cv and skewness coefficient Cs were calculated. The frequency factor Kp corresponding to different return periods, such as 5 years, 10 years, and 20 years, was obtained from tables, and the design rainfall for each return period was calculated.
[0078] Based on design rainfall data points with different durations and return periods, a least squares method is used to fit a rainfall intensity formula. The rainfall intensity formula expresses the relationship between rainfall intensity and return period and rainfall duration; rainfall intensity is positively correlated with return period and negatively correlated with rainfall duration. The formula contains several undetermined parameters, which are solved through nonlinear regression to obtain the rainfall intensity formula for the target road. The goodness of fit of the formula should meet the engineering accuracy requirements.
[0079] S22, based on the parameters of the rainstorm intensity formula, the Chicago rainfall pattern method was used to obtain the time-history rainfall intensity process line;
[0080] The Chicago rainfall pattern method is used to design the rainfall process curve. This method assumes that the rainstorm peak occurs at a certain point in the rainfall duration. A rainstorm peak location coefficient r is set, typically ranging from 0.3 to 0.5, and the time step is determined based on the simulation accuracy requirements. The rainfall intensity before and after the rainstorm peak decreases according to the rainstorm intensity formula. The rainfall intensity sequence is calculated step by step according to the set time step, and the sequence is normalized to make the total rainfall equal to the design value. The resulting time-history rainfall intensity process curve exhibits a single-peak distribution, with an initial increase, reaching its maximum at the rainstorm peak, and then decreasing in the later stages.
[0081] Using the same method, rainfall process curves for various design conditions with combinations of rainfall duration and return periods were generated. Multiple time-series rainfall intensity process curves were obtained, and each process curve was stored in time series format.
[0082] S23, obtain information on the surrounding catchment area and use the inference formula method to obtain the external inflow flow parameters;
[0083] Analyze the water catchment area around the target road section to identify the area, surface type, and runoff coefficient of the catchment area. The catchment area may include different types such as commercial building complexes, residential areas, and green spaces, each with a different runoff coefficient.
[0084] Runoff volume is calculated using a deductive formula, whereby runoff volume equals the comprehensive runoff coefficient multiplied by the rainfall intensity and then multiplied by the catchment area. Runoff volume is calculated separately for different catchment areas based on their parameters.
[0085] Determine the boundary locations of external inflows and identify the locations and quantities of rainwater discharged into roads within the catchment area. Calculate the inflow according to the time-history variation of the rainfall process curve, as the inflow is directly proportional to the rainfall intensity. Obtain the time-history inflow process curves for each boundary location.
[0086] S24. Based on the time-history rainfall intensity process line and external inflow parameters, a multi-condition combined design method is adopted to obtain the set of rainfall boundary conditions;
[0087] Multiple design conditions are formed by arranging various rainfall durations and return periods in a full range of combinations. These conditions should cover rainfall scenarios ranging from common to extreme, such as short-duration high-intensity rainfall and long-duration moderate-intensity rainfall.
[0088] For each working condition, the boundary conditions include a road surface rainfall source term and a lateral inflow term. The rainfall source term uses the time-history rainfall intensity process line generated in step S22, with rainfall acting uniformly across the entire road surface. The lateral inflow term uses the inflow rate calculated in step S23, applied to the boundary locations according to the time-history variation synchronized with the rainfall intensity. A complete boundary condition dataset for multiple working conditions is obtained, and the data is stored in a structured file format.
[0089] S3. Based on the refined three-dimensional digital model of the road surface and the set of rainfall boundary conditions, the drainage performance index is obtained by using the two-dimensional shallow water wave equation simulation method.
[0090] S31. Based on the physical process of road runoff, the two-dimensional shallow water wave control equations are obtained by using the Saint-Venant equations mathematical description method.
[0091] A two-dimensional shallow water wave equation for road surface runoff is established. The equation is based on the laws of conservation of mass and momentum, assuming a uniform vertical velocity distribution and that the pressure distribution satisfies the hydrostatic pressure assumption. The continuity equation describes the variation of water depth with time and space, while the momentum equation describes the variation of flow velocity with time and space. The equation includes convection, pressure, and source terms.
[0092] The source and sink terms include rainfall source terms, pavement infiltration terms, and storm drain discharge terms. Rainfall intensity is directly calculated using the rainfall process data from step S22. Infiltration rate is calculated using appropriate infiltration models for different pavement types, such as the Horton infiltration model or the Green Amter model, with model parameters determined based on pavement material properties. The storm drain discharge rate is non-zero only at the storm drain grid cells and is calculated using weir flow or orifice flow formulas. Parameters in the formulas, such as the discharge coefficient and water passage perimeter, are determined based on the actual size and shape of the storm drain.
[0093] The source terms also include the bottom slope source term and the friction source term. The bottom slope source term characterizes the driving effect of the ground slope on the water flow, and the ground elevation is provided by the pavement geometry model in step S1. The friction source term characterizes the resistance effect of pavement roughness on the water flow, described by the Manning formula, and the Manning roughness coefficient is determined based on the pavement material and texture characteristics. This yields a complete two-dimensional shallow water wave governing equation set, which is a hyperbolic partial differential equation set.
[0094] S32. Based on the refined three-dimensional digital model of the road surface, the Delaunay triangulation method is used to obtain the spatial discrete computational grid.
[0095] Import the 3D road surface model from step S1 into the mesh generation software. Extract the road surface boundary lines as the computational domain boundary, including the curb lines on both sides of the road and the cross-sectional lines at the start and end points of the road. Set mesh size control points on the boundary lines, and set different mesh sizes according to the importance of hydraulic characteristics in different areas. Set localized refinement areas at the locations of storm drains, with significantly reduced mesh sizes to precisely capture the local hydraulic characteristics at the storm drains. In other regular areas of the road surface, the mesh size can be appropriately increased to improve computational efficiency.
[0096] An unstructured triangular mesh is generated using the constrained Delaunay triangulation algorithm. The algorithm first generates uniformly distributed initial nodes within the computational domain, then performs local refinement or sparsification adjustments based on boundary constraints and size control. The algorithm automatically optimizes mesh quality by improving the regularity of triangular elements through edge flipping and node smoothing operations. A quality check is performed on the mesh, calculating quality indicators such as the minimum interior angle and area-to-length ratio of each triangle to ensure that the mesh quality meets computational requirements.
[0097] Interpolate the two-dimensional coordinates of the grid nodes with the road surface elevation model from step S1, assigning a ground elevation value to each grid node. Use inverse distance weighted interpolation or other appropriate interpolation methods. Obtain complete computational grid data containing node coordinates, node elevations, and cell topology relationships.
[0098] S33. Based on the two-dimensional shallow water wave control equations and the spatial discrete computational grid, the numerical solution scheme is obtained by using the finite volume method for discretization.
[0099] The finite volume method is used to spatially discretize the governing equations. For any triangular mesh element, the governing equations are integrated over the element, and the volume integral is transformed into a boundary integral according to the Gaussian divergence theorem. The discrete equations of the element are expressed as follows: the time derivative of the element's conserved variables equals the sum of the negative interface fluxes plus the source term.
[0100] The interface numerical flux was calculated using the HLL approximate Riemann solver. This solver divides the interface state into three regions based on two wave velocity estimates (left and right). The region to which the interface belongs is determined by the wave velocity sign, and the flux is calculated using the corresponding formula. The HLL solver is computationally simple and stable, and is suitable for solving shallow water wave equations.
[0101] In the discretization of the source terms, the bottom slope source term is calculated based on the elevation difference between the two sides of the element, the frictional source term is calculated according to Manning's formula, the rainfall source term is the rainfall intensity multiplied by the element area, and the inlet discharge source term is calculated according to the weir flow formula, only non-zero for elements containing inlets. The time stepping is performed using the second-order Runge-Kutta method in two steps. The time step size is adaptively determined based on the CFL condition, with the CFL number set to 0.5 to ensure numerical stability. A complete numerical solution scheme is obtained.
[0102] S34. Based on the numerical solution format and the set of rainfall boundary conditions, the GPU parallel computing method was used to obtain the road surface water depth field and flow velocity field data.
[0103] The numerical solver is written in CUDA and implements GPU parallel computing. The program runs on a workstation equipped with a high-performance GPU, which has a large number of CUDA cores and ample video memory. Mesh cells are allocated to thread blocks on the GPU, with each thread block containing several threads. Mesh data, elevation data, and conservation variables are stored in the GPU's global memory.
[0104] Write CUDA kernel functions to implement each step of the time-progression computation. The kernel functions are responsible for tasks such as calculating the interface numerical flux, accumulating the flux and calculating source terms to update conserved variables, and calculating the time step size based on CFL conditions. Each thread processes one or more grid cells, achieving large-scale parallel computing. The main program controls the time loop on the CPU, calling the corresponding kernel function sequentially at each time step.
[0105] Simulations were performed for various rainfall scenarios. Initial conditions were set as a dry road surface. Boundary conditions were set according to actual conditions, including free outflow boundaries, solid wall boundaries, and flow inflow boundaries. The simulation covered the period from the start of rainfall to the complete removal of accumulated water after rainfall. The time step was adaptively adjusted according to CFL conditions. Calculation results were output periodically, and water depth and flow velocity data were saved. GPU parallel computing achieved a speedup of tens of times compared to CPU serial computing.
[0106] This yields water depth and velocity field data at multiple output times. Each data set contains water depth and velocity component values for all grid cells. The data is stored in binary or other efficient formats.
[0107] S35, based on the road surface water depth field and flow velocity field data, uses spatiotemporal data analysis methods to obtain drainage performance indicators;
[0108] The output water depth field data is analyzed. A water accumulation threshold is set; areas with water depths exceeding this threshold are considered to have significant water accumulation. For each time-instance water depth field, the number of grid cells with water depths exceeding the threshold is counted, and multiplied by the average cell area to obtain the water accumulation area. The maximum water depth among all cells with water depths exceeding the threshold is then calculated to obtain the maximum water accumulation depth at that time.
[0109] This study analyzes the spatiotemporal evolution of water accumulation during rainfall. After rainfall begins, water gradually accumulates on the road surface, initially appearing at the lowest points of the longitudinal slope and along the cross slope's water catchment line. The area and depth of the accumulated water reach their maximum values during the peak of rainfall. After the rainfall ends, the water gradually recedes until it completely disappears.
[0110] Statistical analysis of drainage performance indicators includes maximum water depth, maximum water accumulation area, and water accumulation duration. Water accumulation duration is defined as the time from the appearance of water to its complete removal. The performance of the existing drainage structure is assessed to determine whether it meets the design objectives.
[0111] Simulation calculations and performance statistics were performed for various rainfall conditions to obtain drainage performance index data for each condition. The data were then compiled into a table to serve as a benchmark for subsequent optimization. Performance benchmark data for the existing drainage structure under multiple conditions were obtained.
[0112] S4. Based on the drainage performance indicators, a drainage structure parameter scheme is obtained by using a neural network surrogate model and an improved genetic algorithm optimization method.
[0113] S41, Obtain drainage structure parameters and engineering constraints, and use the design space definition method to obtain the mathematical model of the optimization problem;
[0114] This study analyzes key structural parameters affecting pavement drainage performance, including cross slope angle, longitudinal slope angle, road crown height, inlet spacing, and inlet size. It determines the current values and adjustable ranges for each parameter. Increasing the cross slope accelerates lateral drainage but affects driving comfort; the longitudinal slope is controlled by the road longitudinal profile design; the road crown height affects the lateral water collection path; the inlet spacing affects water collection capacity and project cost; and the inlet size corresponds to different specifications of standard grates.
[0115] Engineering constraints are set, including constraints on road surface elevation adjustment, storm drain location, traffic safety, and economic factors. Road surface elevation adjustment constraints require that elevation changes not exceed certain limits to meet the requirements for integration with existing ancillary facilities. Storm drain location constraints require that newly added or relocated storm drains avoid existing manholes and underground pipelines. Traffic safety constraints require that cross slopes and longitudinal slopes not exceed the limits specified in the standards. Economic constraints limit the number of new storm drains to control project costs.
[0116] Establish a parameterized optimization model. Define a decision variable vector containing key structural parameters. Determine the value range of each variable. Establish the mathematical expressions for the constraints. Define the objective function, which is a comprehensive evaluation of drainage performance indicators, including a weighted combination of maximum water accumulation depth, water accumulation area, and drainage time under multiple operating conditions. The objective is to minimize this comprehensive indicator. This yields the complete definition of the constrained optimization problem.
[0117] S42. Based on the mathematical model of the optimization problem, the parameter-performance mapping sample set is obtained by using Latin hypercube sampling and CFD simulation calculation methods.
[0118] The Latin hypercube sampling (LHS) method is used to generate initial sample points in the multidimensional parameter space. The LHS method divides the range of each parameter into several equal intervals, randomly samples a value within each interval, and ensures that the sample points are uniformly distributed in the parameter space through permutations and combinations. This generates a certain number of parameter samples, each containing the specific values of all parameters.
[0119] Each sample group's parameters are checked to determine if they meet the constraints. Samples that do not meet the constraints are removed, and new samples are collected to supplement the sample set. After constraint filtering, a valid sample set that meets the constraints is obtained.
[0120] CFD simulations were performed on the valid samples. For each sample group, the road geometry model from step S1 was adjusted based on the parameter values, modifying the cross slope, longitudinal slope, road crown, and storm drain layout parameters, and a new computational mesh was generated. Then, the method from step S3 was used to perform hydrodynamic simulations for various rainfall conditions, calculating the drainage performance indicators for each condition. Parallel computation using multiple GPU workstations significantly reduced the total computation time.
[0121] The parameters and performance mapping data of the samples are obtained. Each set of data includes input parameters and output performance indicators, which include the maximum water accumulation depth, maximum water accumulation area, and drainage time for various operating conditions. The data is stored in a structured format.
[0122] S43. Based on the parameter and performance mapping sample set, a physical constraint neural network training method is used to obtain the surrogate model;
[0123] To address the issue of long computation times in CFD simulations, a neural network surrogate model is constructed to replace simulation for rapid performance evaluation. To improve the generalization ability of the surrogate model under small sample conditions, hydraulic physical constraints are innovatively embedded in the neural network training.
[0124] A neural network proxy model is established. The network input is a design parameter vector x, which includes key drainage structure parameters such as road cross slope, longitudinal slope, road crown height, and the spacing and size of storm drain inlets. The network output is a performance index vector y, which includes performance indicators such as maximum water accumulation depth, maximum water accumulation area, and drainage time under various typical rainfall conditions. The network adopts a multi-layer fully connected structure, with the number of neurons in the hidden layers determined according to the input and output dimensions. The activation function used is the ReLU function or other nonlinear activation functions.
[0125] A physical constraint loss function is constructed. The total loss function consists of two parts: data fitting loss and physical constraint loss. These two parts are weighted and combined using weight coefficients, which are adjusted based on the training effect. The data fitting loss uses mean squared error to calculate the error between the network's predicted value and the true value of the CFD simulation sample. The physical constraint loss is constructed based on the hydraulic laws of drainage performance, including three constraints: First, the water depth is negatively correlated with the cross slope; the larger the cross slope, the faster the drainage and the shallower the water accumulation. Second, the water accumulation area is positively correlated with the distance between the rainwater inlets; the larger the distance, the weaker the water collection capacity and the more water accumulation. Third, the drainage time is positively correlated with the rainfall intensity; the stronger the rainfall, the longer the drainage time. These three physical laws are expressed as mathematical constraints. The partial derivatives of the network prediction results with respect to the parameters are calculated to determine whether the monotonicity conforms to the physical laws. If the laws are violated, a penalty term is added.
[0126] An incremental sampling strategy is employed to improve model accuracy. In the initial stage, the neural network is trained using the LHS samples generated in step S42, employing the Adam optimizer or other gradient optimization algorithms, with appropriate learning rates and training epochs. After training, several new sampling points are uniformly selected in the parameter space, and CFD simulations are performed on these points. The new samples are then added to the training set, and the network is retrained. This process is repeated several times, gradually expanding the size of the training set.
[0127] Verify the accuracy of the surrogate model. Several test samples were randomly selected, and performance metrics were calculated using both the surrogate model and CFD simulation. The predicted values were compared with the actual values. Statistical results show that the average relative error of all performance metrics predicted by the surrogate model is less than 10%, and the accuracy meets the optimization requirements. The surrogate model's prediction time per run is significantly shorter than that of CFD simulation, achieving a speedup of thousands to tens of thousands of times. A high-precision, fast surrogate model is thus obtained.
[0128] S44. Based on the surrogate model and the mathematical model of the optimization problem, the Pareto optimal solution set is obtained by using the adaptive constraint penalty genetic algorithm method.
[0129] To address the multi-objective and multi-constraint characteristics of the road drainage optimization problem, an improved non-dominated sorting genetic algorithm (NSGA-II) is employed for optimization. An innovative adaptive constraint penalty mechanism is designed, transforming hard constraints into penalty terms in the objective function to ensure that the optimization process converges within the feasible region.
[0130] A multi-objective optimization mathematical model is established. Two optimization objectives are defined: Objective 1 is to minimize the average maximum water accumulation depth under multiple working conditions, and Objective 2 is to minimize the engineering modification cost. The engineering cost includes the cost of adjusting the road surface elevation and the cost of modifying the storm drains. The cost of adjusting the elevation is proportional to the amount of adjustment, and the cost of modifying the storm drains is proportional to the number of new storm drains. The two objectives conflict, and a Pareto optimal balance needs to be sought.
[0131] Design an adaptive constraint penalty function. Define multiple types of constraints, including cross slope not exceeding the specification limit, longitudinal slope not exceeding the specification limit, road surface elevation adjustment not exceeding the set limit, and the number of newly added storm drains not exceeding the budget limit. For each constraint, calculate the constraint violation amount. If the constraint is satisfied, the violation amount is zero; otherwise, the violation amount is the excess. The total constraint violation rate equals the sum of all constraint violations. Design an adaptive penalty coefficient, setting an initial value. After each generation of evolution, if the proportion of feasible solutions in the population is lower than the first threshold, the penalty coefficient is decreased to relax the constraint; if the proportion of feasible solutions is higher than the second threshold, the penalty coefficient is increased to strengthen the constraint. The modified objective function equals the original objective function plus the penalty coefficient multiplied by the total constraint violation rate, transforming the constrained optimization problem into an unconstrained optimization problem.
[0132] Execute the NSGA-II optimization algorithm. Initialize the population; the population size is determined based on the parameter dimension. Randomly generate several sets of individual parameters, ensuring that the parameters are within their range. For each individual, call the surrogate model from step S43 to evaluate performance metrics, calculate the two objective function values and constraint violation degree, and obtain the corrected objective function value.
[0133] Perform genetic operations. A tournament selection strategy is used to select parent individuals. A simulated binary crossover (SBX) operator is used for crossover, with a high crossover probability. A polynomial mutation operator is used for mutation, with a low mutation probability. Crossover and mutation produce the offspring population.
[0134] The merged parent and offspring population is subjected to a non-dominated ordination process, dividing individuals into multiple non-dominated levels. The first level is the Pareto front, where individuals are not mutually dominant. For individuals within the same level, crowding distance is calculated, which characterizes the distribution density of individuals in the target space. Following the principle of level priority followed by crowding distance, a number of individuals equal to the population size are selected to enter the next generation.
[0135] The genetic operations and environmental selection are performed iteratively, co-evolving for several generations until the algorithm converges. The Pareto front of each generation is recorded. After the algorithm converges, the final Pareto front is extracted, which contains several non-dominated solutions.
[0136] The recommended solution is selected from the Pareto front. A comprehensive evaluation index is calculated for each solution; this index is a weighted sum of the two objectives, with weights determined based on decision preferences. The solution with the optimal comprehensive index is selected as the recommended solution. This solution achieves a balance between performance and cost while satisfying all constraints. The Pareto solution set for multi-objective optimization and the recommended balanced solution are obtained.
[0137] S45. Based on the Pareto optimal solution set, a multi-scheme comparison and evaluation method is used to obtain the recommended parameter configuration scheme;
[0138] A detailed analysis is performed on the Pareto front solution set obtained in step S44. The Pareto front contains several non-dominated solutions, each representing a set of parameter configurations and corresponding performance-cost indices. These solutions are plotted as a scatter plot in a two-dimensional target space, with the horizontal axis representing the average maximum water accumulation depth under multiple operating conditions and the vertical axis representing the engineering modification cost.
[0139] The distribution characteristics of the Pareto front are analyzed. Solutions at one end of the front have the best performance but the highest cost, while solutions at the other end have the lowest cost but poorest performance. Solutions in the middle achieve varying degrees of balance between performance and cost. Several representative solutions are selected for comparison, including performance-priority solutions, balanced solutions, and cost-priority solutions. The performance-priority solution achieves optimal drainage performance by increasing the cross slope and densely arranging storm drains, but this involves a large amount of engineering work and is costly. The cost-priority solution controls costs by reducing the amount of modification, but the performance improvement is limited. The balanced solution achieves a balance between the two.
[0140] Based on the actual engineering needs and budget constraints, the recommended scheme is selected from the Pareto front. This scheme significantly reduces the maximum water accumulation depth under common rainfall return period conditions, showing a marked improvement in performance compared to the current situation. The modification cost is within budget, making it economically feasible. The recommended optimal parameter configuration scheme is thus obtained.
[0141] S46. Based on the recommended parameter configuration scheme, the drainage structure parameter scheme is obtained by using CFD simulation verification method.
[0142] Accurate CFD simulations were performed to verify the candidate optimal solution. Based on the optimization parameters, the pavement geometry model from step S1 was modified, adjusting parameters such as cross slope, longitudinal slope, road crown height, and the spacing and size of storm drain inlets. The locations of newly added storm drain inlets were determined based on site surveys, avoiding prohibited areas.
[0143] The 3D model of the road surface and the computational mesh were regenerated based on the adjusted parameters. CFD simulations were performed for various rainfall conditions using the GPU parallel method described in step S3.
[0144] The performance metrics calculated by simulation were obtained. The reliability of the surrogate model was verified by comparing its predictions with those of the surrogate model. The simulation results were analyzed to evaluate the degree of performance improvement for each return period condition. For common return period conditions, the maximum water depth met the design target. For more extreme return period conditions, although the stringent target was not met, the water depth was significantly reduced compared to the current situation.
[0145] Considering the special characteristics of extreme operating conditions and the limitations imposed by underground pipeline capacity and engineering constraints, the current optimization scheme has achieved the best performance achievable under the given constraints. Taking into account performance, constraints, and economy, the current scheme is confirmed as the optimal balanced solution.
[0146] A drainage structure parameter scheme, verified by precise simulation, was obtained. This scheme enables the design target to be achieved under common return period conditions, and significantly improves performance under extreme return period conditions, thus meeting engineering requirements.
[0147] S5. Based on the drainage structure parameter scheme, the parameters of the calibrated simulation model are obtained by using field tests and parameter calibration methods.
[0148] S51, obtain representative road section information, and use on-site artificial rainfall test method to obtain measured data on water accumulation response;
[0149] A representative test section was selected from the target road segment. This section included the lowest point of the longitudinal slope and areas with severe water accumulation, making it representative of the target road. The current parameters of the test section were recorded.
[0150] A portable artificial rainmaking device was constructed. The device consisted of a water tank, a booster pump, pipes, and nozzles. A support frame was erected at a certain height above the test section, and horizontal branch pipes were installed. Several conical nozzles were mounted on the branch pipes, covering the entire road surface of the test section. Different rainfall intensities were achieved by adjusting the pump flow rate and the number of nozzles.
[0151] Multiple sets of rainfall test conditions were designed, including different types of rainfall such as moderate-intensity rainfall and heavy rainfall, to simulate different rainfall scenarios.
[0152] Water depth measuring points were set up on the road surface. Portable ultrasonic water level gauges were used, with measurement range and accuracy meeting the requirements. Several measuring points were set up longitudinally along the test section. The measuring points were located at the water catchment line on the cross slope of the road surface, i.e., the deepest point of the water accumulation. The water level gauges were fixed on temporary supports, with the probes pointing downwards at the road surface, and water depth data were recorded in real time.
[0153] Conduct various test conditions. Start the rainfall device and adjust the flow rate to stabilize the rainfall intensity at the set value. After rainfall begins, the water level gauges at the measuring points automatically record data. Rainfall continues for a certain period before stopping, but recording continues until the water recedes. The water flow pattern is recorded by a camera throughout the test.
[0154] Multiple sets of water depth time history curves were obtained for the experimental conditions. The data include information such as time, measuring point number, and water depth.
[0155] S52. Based on the measured data of water accumulation response, the parameters of the calibrated simulation model are obtained by using the segmented parameter identification method.
[0156] To address the uncertainty of key parameters such as pavement roughness coefficient and storm drain discharge coefficient in simulation models, an innovative piecewise parameter identification method based on water depth response characteristics is proposed. This method utilizes the differences in sensitivity to different parameters at different stages of the water accumulation process to achieve decoupled parameter identification, avoiding the parameter coupling problem of traditional global optimization methods.
[0157] Analyze the phased characteristics of the water accumulation process. Road surface water accumulation can be divided into three phases: the water rise phase, the peak maintenance phase, and the water receding phase. The water rise phase is mainly controlled by rainfall intensity and road surface roughness; the greater the roughness, the slower the water flow and the faster the water rises. The peak maintenance phase is mainly controlled by the balance between rainfall intensity and the drainage capacity of the storm drains. The water receding phase is mainly controlled by the drainage capacity of the storm drains; the greater the drainage coefficient, the faster the water recedes.
[0158] Extract water depth response indices at characteristic moments. For the experimental data in step S51, extract three characteristic indices: Indicator 1 is the water rise rate, defined as the average rise rate from the appearance of water to reaching its peak; Indicator 2 is the peak water depth, defined as the maximum water depth value during the entire process; Indicator 3 is the water receding rate, defined as the reciprocal of the time required for the water to recede from its peak to a set proportion. For multiple sets of experimental conditions at multiple measuring points, calculate these three indices respectively, obtaining several sets of characteristic data.
[0159] Sensitivity analysis of parameters was conducted. The single-factor perturbation method was used to analyze the influence of parameters on each index. With the discharge coefficient fixed, the roughness coefficient was varied within a reasonable range, and simulations were performed under typical working conditions to calculate the changes in the three indices. The results show that the water rise rate is highly sensitive to the roughness coefficient, with a sensitivity significantly higher than that of other indices; the peak water depth is moderately sensitive to the roughness coefficient; and the water receding rate is not sensitive to the roughness coefficient.
[0160] With the roughness coefficient fixed, the discharge coefficient was varied within a reasonable range, and the same analysis was performed. The results showed that the water receding rate was highly sensitive to the discharge coefficient, significantly more so than other indicators; the peak water depth was moderately sensitive to the discharge coefficient; and the water rise rate was not sensitive to the discharge coefficient.
[0161] Sensitivity analysis shows that the roughness coefficient mainly affects the water rise process, while the discharge coefficient mainly affects the water receding process. The effects of the two parameters at different stages are relatively independent, and decoupled identification can be achieved.
[0162] Perform segmented parameter identification. Identify the roughness coefficient. The water rise rate index is used for identification; this index is sensitive to the roughness coefficient but not to the discharge coefficient. Establish the relationship between the roughness coefficient and the rise rate. A lookup table is created using CFD simulation. The roughness coefficient is taken at certain intervals within a reasonable range, and the corresponding rise rate is calculated. Substitute the experimentally measured rise rate into the lookup table and use linear interpolation to obtain the roughness coefficient. Identify several sets of experimental data separately, obtaining several roughness coefficient values. Calculate the average value as the calibration value.
[0163] Identify the discharge coefficient. Using a fixed roughness coefficient as the calibration value in the first step, the discharge rate index is employed for identification. A relationship between the discharge coefficient and the discharge rate is established. A lookup table is created using CFD simulation. The discharge coefficient is taken at certain intervals within a reasonable range, and the corresponding discharge rate is calculated. The experimentally measured discharge rate is substituted into the lookup table, and linear interpolation is used to obtain the discharge coefficient. Several sets of experimental data are identified separately, yielding several discharge coefficient values. The average value is calculated as the calibration value.
[0164] Verification of calibration effectiveness. The experimental conditions were re-simulated using the calibrated parameters, and the root mean square error (RMSE) between the predicted and measured water depths was calculated. The error after calibration was significantly reduced compared to before calibration, with a reduction of over 50%. More importantly, the simulated curves after calibration showed good agreement with the measured curves during the rising, peak, and receding phases, indicating that the segmented identification method can accurately capture the hydraulic characteristics at different stages.
[0165] Compared with traditional global optimization methods, the piecewise parameter identification method has the following advantages: it achieves decoupled identification by utilizing differences in parameter sensitivity, avoiding the problem of multiple solutions caused by parameter coupling; based on the stage division of the physical process, the identification results have clear physical meaning; and it uses lookup tables and linear interpolation, resulting in high computational efficiency and eliminating the need for iterative optimization.
[0166] The simulation model parameters, calibrated by segment identification, are obtained. These parameters are applied to subsequent full-segment simulations to improve prediction reliability.
[0167] S53. Based on the calibrated simulation model parameters and drainage structure parameter scheme, the CFD simulation calculation method is used to obtain the performance prediction results of the optimized scheme.
[0168] The calibrated parameters were applied to the simulation model for the entire road section. The optimal solution obtained in step S4 was then recalculated using CFD simulations under various rainfall conditions. The pavement roughness coefficient and storm drain discharge coefficient in the model were set to the calibrated values, while all other settings remained unchanged.
[0169] Simulation results show that the water depth is slightly lower under calibrated parameters compared to the simulation using the initial parameters. Analysis reveals that the calibrated roughness coefficient is lower than the initial value, indicating that the actual road surface friction is lower than the model's preset value, resulting in smoother water flow and faster drainage. Conversely, the calibrated discharge coefficient is higher than the initial value, indicating that the actual water collection capacity of the storm drain is stronger than the model's preset value, leading to more efficient discharge.
[0170] High-reliability simulation results based on field test calibration were obtained. These results demonstrate that the optimized scheme outperforms the current scheme under all return periods. The reliability of the simulation results is supported by field test data, indicating high credibility.
[0171] S54. Based on the performance prediction results of the optimization scheme and the measured data of water accumulation response, a comparative verification method is used to obtain the reliability assessment conclusion of the optimization scheme.
[0172] A comparative analysis of the performance before and after optimization was conducted. Under typical operating conditions, the maximum water accumulation depth, maximum water accumulation area, and drainage time of the existing scheme were all higher than those of the optimized scheme. The optimization effect was significant, with all indicators decreasing by more than 60%.
[0173] Statistical analysis was conducted on the average performance under various operating conditions. The average maximum water accumulation depth, average maximum water accumulation area, and average drainage time of the current solution were all higher than those of the optimized solution. All three indicators showed significant improvement.
[0174] The robustness of the optimized solution was analyzed. The performance standard deviation under various operating conditions was calculated. The standard deviation of the optimized solution was lower than that of the current solution, indicating that the optimized solution's performance is more stable under different operating conditions. The optimized solution has better robustness and can adapt to various rainfall conditions.
[0175] The reliability of the simulation model was evaluated. Through parameter calibration via field tests, the simulation error was significantly reduced, and the accuracy of the calibrated model met the requirements for engineering applications. The performance indicators of the optimized scheme predicted by the simulation were verified and calculated multiple times, and the results were stable and repeatable. In conclusion, the simulation model is reliable, and the performance predictions of the optimized scheme are credible.
[0176] The reliability assessment of the optimized solution shows that it can control the maximum water accumulation depth within the design target under common return period rainfall conditions. Compared with the current situation, the performance is improved, the solution meets the engineering design objectives, is technically feasible, and is recommended for implementation.
[0177] S6. Based on the calibrated simulation model parameters and drainage structure parameter scheme, engineering implementation and online monitoring methods are used to obtain the engineering effect verification results;
[0178] S61. Based on the drainage structure parameter scheme, the engineering design document is obtained by using the construction drawing preparation method;
[0179] Based on the optimal parameters obtained in step S4 and the verification results in step S5, detailed construction drawings are prepared. Road surface elevation adjustments are designed using road design software. The cross slope, longitudinal slope, and camber height of the road surface are adjusted according to the optimized parameters by modifying the road cross-section and vertical curve parameters. The software automatically recalculates the design elevation of each point on the road surface. A road surface elevation adjustment table is exported, and the elevation changes are checked to ensure they meet the constraints.
[0180] Optimize the layout of storm drains. Adjust the spacing and size of storm drains based on optimization parameters to determine the locations of retained, removed, and newly added storm drains. The locations of new storm drains are determined based on the road surface model and site survey, avoiding existing manholes and underground pipelines. Mark the plane coordinates, elevation, type, and specifications of each storm drain in the design software. Mark the connection method between the storm drain and the storm drain branch pipe.
[0181] Construction drawings were prepared. These drawings included a road surface elevation adjustment plan, longitudinal section, cross section, storm drain inlet layout plan, and detailed drawings of the storm drains. The drawings were in CAD format and the scale met construction requirements. This resulted in the final construction drawings.
[0182] Develop a construction plan. The plan includes construction techniques, construction sequence, construction machinery, and quality control measures. Construction techniques include road milling, subgrade leveling and compaction, asphalt concrete paving, and storm drain inlet modification. The construction sequence is segmented and closed construction, determining the length and period of each segment. Construction machinery includes milling machines, pavers, road rollers, and excavators. Quality control measures include elevation verification measurement, road surface evenness testing, compaction testing, and storm drain inlet installation location inspection. Obtain the construction organization design document, which, after approval, serves as the basis for construction.
[0183] S62. Based on the engineering design documents, a segmented construction and quality inspection method was adopted to obtain the modified road drainage system.
[0184] Construction will be organized according to the construction plan. The construction will be divided into several stages, with each stage covering a certain length. Before construction begins, traffic will be diverted, construction barriers and signs will be set up, the construction lane will be closed, and the other lanes will remain open for traffic.
[0185] The construction process at each stage includes road milling, subgrade leveling and compaction, asphalt concrete paving, and storm drain modification. A milling machine is used to mill the road surface, with the milling depth controlled according to an elevation adjustment table. After milling, debris is removed, and the subgrade is leveled and compacted, with the compaction degree tested to ensure it meets requirements.
[0186] Total stations were used for paving elevation control. Elevation control stakes were set on both sides of the road section. During paving, the total station was used to measure the paving thickness in real time to ensure that the road surface elevation after paving met the design value. Asphalt concrete was paved, and then compacted with a road roller to ensure that the compaction degree met the standards. After paving was completed, the road surface smoothness was tested to ensure that it met the specifications.
[0187] Construction of the storm drain inlet renovation project was carried out. An excavator was used to excavate the storm drain inlet chamber, a concrete foundation was poured at the bottom of the chamber, and precast concrete manhole rings were installed on the side walls. A hole was made at the bottom of the chamber to connect to the storm drain branch pipe, and a storm drain grate was installed on the top of the chamber, with the top surface of the grate flush with the road surface. The elevation of the storm drain inlet was checked using a level instrument, and after it was deemed acceptable, backfilling was carried out around the chamber.
[0188] After each stage of construction is completed, a quality inspection will be conducted. The inspection will include checking the road surface elevation, smoothness, compaction, and the location and elevation of drainage inlets. Once the inspection results are satisfactory, the road will be opened to traffic, and the next stage of construction will begin.
[0189] After all construction is completed, a final acceptance inspection will be conducted. The inspection will involve measuring the overall elevation of the road surface, drawing an as-built cross-section diagram, comparing it with the design cross-section, and ensuring that all indicators meet the requirements. Upon successful final acceptance, the modified road drainage system will be obtained.
[0190] S63, based on the modified road drainage system, online monitoring technology is used to obtain drainage process monitoring data;
[0191] Online monitoring equipment was installed on the renovated road surface. Several monitoring points were selected, located in areas with historically severe waterlogging and upstream catchment areas of storm drains, and these monitoring points were numbered.
[0192] Each monitoring point is equipped with an ultrasonic water level gauge. The gauge's measurement range and accuracy meet the monitoring requirements. The gauge is mounted on a monitoring pole on the roadside and fixed to the curb using a fixing device. The gauge probe is pointed downwards towards the waterline on the road surface. The sampling interval is set according to monitoring needs; recording automatically begins when the water depth exceeds a set threshold and stops after the water depth remains below the threshold for a certain period.
[0193] Each monitoring point is equipped with a video surveillance system. The surveillance cameras are network high-definition cameras with resolution and frame rate meeting monitoring requirements. The cameras are mounted on roadside monitoring poles, with a top-down angle and field of view covering the road surface surrounding the monitoring point. The cameras are set to record event-triggered data; when the water level gauge detects that the water depth exceeds a set threshold, the camera automatically starts recording, and the recording is saved to local storage.
[0194] The monitoring equipment transmits data to the monitoring center server via a wireless network. The server deploys a data management system that receives and stores monitoring data in real time, providing data query and chart display functions. The system sets warning thresholds; when the water depth exceeds the warning value, it automatically sends a warning message to management personnel.
[0195] After the monitoring system was put into operation, it encountered several significant rainfall events during its operation. The monitoring data of these rainfall events were analyzed.
[0196] For rainfall events with different return periods, the monitoring points recorded the peak water depth and duration. Under common return period rainfall conditions, the water depth at the monitoring points was controlled within the design target range, and the drainage process was smooth. Video surveillance showed no obvious water accumulation areas on the road surface, and vehicles could pass normally.
[0197] Under rainfall conditions with a relatively large recurrence interval, the peak water depth at the monitoring points increased, but the water depth was significantly lower than that under similar historical rainfall conditions before the upgrade. Video surveillance showed slight water accumulation on the road surface, but the water depth was shallow, the distribution area was small, and vehicle traffic was basically unaffected.
[0198] Under extreme recurrence interval rainfall conditions, the peak water depth at the monitoring points was relatively high, but the actual water depth was still lower than before the upgrade, indicating a significant performance improvement. Video surveillance showed obvious water accumulation on the road surface, but the water receded quickly, and the water accumulation was basically eliminated within a short period after the rainfall ended, a substantial reduction compared to before the upgrade. Although there were brief periods of slow traffic, no traffic jams occurred.
[0199] Monitoring data for several rainfall events were obtained, including the time, amount, intensity, water depth time-history curves at each monitoring point, and video recordings. The data is stored in a server database and is available for querying and analysis.
[0200] S64. Based on the drainage process monitoring data, statistical analysis and comparative evaluation methods are used to obtain the engineering effect verification data;
[0201] Statistical analysis was performed on the monitoring data. Key performance indicators for several rainfall events were extracted, including maximum water depth and duration of water accumulation. Historical rainfall data before the renovation were compared, and the water depth and duration under similar rainfall conditions before the renovation were estimated using historical video footage and maintenance records.
[0202] For different rainfall events, the maximum water depth and duration of water accumulation after the modification were both lower than those under similar historical rainfall conditions before the modification. Average data from multiple rainfall events showed that both the maximum water depth and duration of water accumulation were significantly reduced after the modification compared to before. Monitoring data verified the effectiveness of the optimization scheme; the actual results were largely consistent with the simulation predictions, indicating high simulation reliability.
[0203] A social benefit assessment was conducted. The improved road section exhibited enhanced traffic capacity during rainfall, reduced traffic congestion time, and increased vehicle throughput. Traffic flow monitoring data showed that average vehicle speeds increased during rainfall, travel time was significantly shortened, and the number of traffic accidents caused by rainfall decreased substantially. Nearby businesses and residents reported a significant improvement in flooding issues and enhanced convenience of travel during heavy rain.
[0204] An economic benefit assessment was conducted. The total investment for this renovation project includes the costs of road milling and paving, drainage inlet renovation, and monitoring equipment. Annual cost savings were calculated based on parameters such as average daily traffic volume, the proportion of rainy days, travel time saved per rainfall event, and the time value of vehicles, leveraging the time cost reduction resulting from reduced traffic congestion. The project has a short payback period and good economic benefits.
[0205] A project effectiveness evaluation report was prepared. The report includes sections on project overview, monitoring plan, monitoring data, performance comparison analysis, social benefit assessment, economic benefit assessment, conclusions, and recommendations, and includes monitoring data charts and on-site photographs. The report concludes that the optimized plan achieved the expected results in the actual project, improving road surface waterlogging; drainage performance met design targets under common return-time rainfall conditions; and performance was significantly improved under extreme return-time conditions. The plan is technically feasible and economically reasonable, and its application in other similar road sections is recommended. Verification data of the project's effectiveness was obtained.
[0206] In one embodiment of this invention, the focus is on the application of road engineering drainage system optimization and renovation, particularly suitable for improving the drainage capacity of existing urban arterial roads, expressways, and urban center roads in the context of frequent rainstorms and flooding. Addressing practical constraints such as limited underground space in urban built-up areas, difficulties in upgrading stormwater pipe networks, and the inability to interrupt traffic operations for extended periods, this method optimizes surface drainage structural parameters such as road cross slope, longitudinal slope, road crown height, and stormwater inlet layout to achieve a refined improvement in road surface drainage performance without significantly altering the underground pipe network. This method can be applied to various engineering scenarios, including urban flood control, sponge city construction, renovation of old roads, and optimization of drainage design for new roads, providing technical support for the scientific design and precise renovation of road engineering drainage systems.
[0207] The following are examples of key data obtained in practical applications. Table 1 shows a comparison of pavement drainage structure parameters before and after optimization, including key parameters such as pavement cross slope, longitudinal slope, road crown height, and the spacing, size, and number of storm drains. Table 2 shows monitoring data of pavement drainage performance under various rainfall conditions, including rainfall characteristic parameters and drainage performance indicators under different return periods of rainfall.
[0208] Table 1 Comparison of pavement drainage structure parameters before and after optimization
[0209]
[0210] Table 2. Monitoring data of pavement drainage performance under various rainfall conditions.
[0211]
[0212] Note: Rainfall duration is in minutes, cumulative rainfall is in millimeters, average rainfall intensity is in millimeters per minute, maximum water depth is in centimeters, maximum water area is in square meters, and waterlogging duration is in minutes.
[0213] By applying the method of this invention, the road drainage structure parameters are optimized and adjusted. The cross slope of the road is appropriately increased to accelerate lateral drainage, the road crown height is increased to optimize the water collection path, and the spacing between rainwater inlets is reduced while their number is increased to improve water collection capacity. The optimized drainage system exhibits good drainage performance under various rainfall conditions, and key indicators such as maximum water accumulation depth, maximum water accumulation area, and water accumulation duration are effectively controlled, meeting the drainage design requirements of road engineering. This method provides an effective technical means for the scientific optimization of road engineering drainage systems.
[0214] The embodiments of the present invention have been described above. However, the embodiments are not limited to the specific implementation methods described above. The specific implementation methods described above are merely illustrative and not restrictive. Those skilled in the art can make more equivalent embodiments under the guidance of the present embodiments, and all of them are within the protection scope of the present embodiments.
Claims
1. A method for optimizing road drainage structure based on fluid dynamics simulation, characterized in that, Includes the following steps: S1. Collect geographic information of the target road section and use UAV oblique photography and point cloud data processing methods to obtain a refined three-dimensional digital model of the road surface. S2, acquire historical rainfall data, and use frequency analysis and design rainfall pattern generation methods to obtain a set of rainfall boundary conditions; S3. Based on the refined three-dimensional digital model of the road surface and the set of rainfall boundary conditions, the drainage performance index is obtained by using the two-dimensional shallow water wave equation simulation method. S4. Based on the drainage performance indicators, a drainage structure parameter scheme is obtained by using a neural network surrogate model and an improved genetic algorithm optimization method. S5. Based on the drainage structure parameter scheme, the calibrated simulation model parameters are obtained by using field tests and parameter calibration methods. S6. Based on the calibrated simulation model parameters and drainage structure parameter scheme, engineering implementation and online monitoring methods are used to obtain the engineering effect verification results.
2. The method for optimizing road drainage structure based on fluid dynamics simulation according to claim 1, characterized in that, S1 includes: Based on a high-resolution oblique image dataset, an incremental motion recovery structure algorithm is used to recover camera extrinsic parameters and sparse 3D point clouds. Combined with ground control points, a bundle adjustment algorithm is used for overall optimization. A semi-global matching algorithm is used to generate a high-density point cloud through stereo matching. Statistical filtering is then applied to the point cloud to remove noise points and outliers, resulting in a high-density three-dimensional point cloud model of the road surface. A threshold segmentation method based on elevation and intensity features is used to identify different road elements in the point cloud. By setting elevation and density thresholds, the point cloud of the road surface is automatically extracted, and the location of rainwater inlets and road boundary lines are identified. The road is divided into multiple segments along its longitudinal direction. For each segment, the least squares method is used to fit the plane to calculate the cross slope and longitudinal slope, resulting in a parametric pavement geometric model that includes slope parameters and roughness coefficients.
3. The method for optimizing road drainage structure based on fluid dynamics simulation according to claim 1, characterized in that, S2 includes: The maximum annual rainfall was obtained by extracting the maximum rainfall over different durations using the sliding time window method. The coefficient of variation and skewness coefficient were calculated by fitting the data using the Pearson type III distribution. The formula for rainstorm intensity was fitted using the least squares method. The Chicago rainfall pattern method was used to design the rainfall process line. The peak rainfall location coefficient was set, and the rainfall intensity sequence was calculated step by step according to the set time step. The sequence was normalized and adjusted so that the total rainfall was equal to the design value, and the time history rainfall intensity process line was obtained. Runoff volume is calculated using a reasoning formula. Runoff volume equals the comprehensive runoff coefficient multiplied by the rainfall intensity and then multiplied by the catchment area. The inflow volume is converted according to the time history variation of the rainfall process line to obtain the time history inflow process line at each boundary location.
4. The method for optimizing road drainage structure based on fluid dynamics simulation according to claim 1, characterized in that, S3 includes: A two-dimensional shallow water wave equation based on the laws of mass conservation and momentum conservation is established. The source and sink terms include rainfall source terms, pavement infiltration terms, and storm drain discharge terms. The storm drain discharge rate is calculated using the weir flow formula or the orifice flow formula, and the friction source term is described using the Manning formula. Unstructured triangular meshes are generated using the constrained Delaunay triangulation algorithm. Locally refined regions are set at the locations of rainwater inlets, and mesh quality is optimized through edge flipping and node smoothing operations. The finite volume method is used to spatially discretize the governing equations, the interface numerical flux is calculated using the HLL approximate Riemann solver, the time step is driven by the second-order Runge-Kutta method, and the time step is adaptively determined according to the CFL condition. Set a water accumulation threshold, count the number of grid cells with water depth exceeding the threshold to obtain the water accumulation area, count the maximum water depth value to obtain the maximum water accumulation depth, and count the duration of water accumulation.
5. The method for optimizing road drainage structure based on fluid dynamics simulation according to claim 1, characterized in that, S4 includes: A neural network proxy model is established. The network input is a vector of design parameters, including the road cross slope, longitudinal slope, road crown height, and the spacing and size of the storm drains. The network output is a vector of performance indicators, including the maximum water accumulation depth, the maximum water accumulation area, and the drainage time. A physical constraint loss function is constructed. The total loss function is a weighted combination of data fitting loss and physical constraint loss. The physical constraint loss includes three constraints: water depth is negatively correlated with cross slope, water area is positively correlated with distance between rainwater inlets, and drainage time is positively correlated with rainfall intensity. The partial derivatives of the network prediction results with respect to the parameters are calculated to determine whether the monotonicity conforms to physical laws. If the laws are violated, a penalty term is added. An incremental sampling strategy is adopted. After training is completed, new sampling points are selected in the parameter space for CFD simulation calculation, and the new samples are added to the training set to retrain the network.
6. The method for optimizing road drainage structure based on fluid dynamics simulation according to claim 1, characterized in that, S4 further includes: A multi-objective optimization mathematical model is established, with objective 1 being to minimize the average maximum water accumulation depth under multiple working conditions and objective 2 being to minimize the engineering modification cost. Design an adaptive constraint penalty function to calculate the constraint violation amount for each constraint condition. The total constraint violation amount is equal to the sum of all constraint violations. After each generation of evolution, adjust the penalty coefficient based on the comparison results between the proportion of feasible solutions in the population and the preset first and second thresholds. The non-dominated sorting genetic algorithm is executed. A tournament selection strategy is used to select parent individuals. A simulated binary crossover operator is used for crossover operations, and a polynomial mutation operator is used for mutation operations. The merged population is non-dominated sorted and the crowding distance is calculated. Individuals are selected to enter the next generation according to the principle of hierarchical priority and second-best crowding distance.
7. The method for optimizing road drainage structure based on fluid dynamics simulation according to claim 1, characterized in that, S5 includes: The process of road water accumulation is divided into three stages: water accumulation rise stage, peak water depth stage and water accumulation recede stage. Three characteristic indicators are extracted: water accumulation rise rate, peak water depth and water accumulation recede rate. Sensitivity analysis of parameters showed that the rate of water accumulation rise is sensitive to the roughness coefficient but not to the discharge coefficient, while the rate of water accumulation recedes is sensitive to the discharge coefficient but not to the roughness coefficient. The roughness coefficient is identified by the water rise rate index. A lookup table of roughness coefficient and rise rate is established. The measured rise rate is substituted into the lookup table and linear interpolation is used to obtain the roughness coefficient. Using a fixed roughness coefficient as the calibration value, the discharge coefficient is identified by the water receding rate index. A lookup table of discharge coefficient and receding rate is established. The discharge coefficient is obtained by substituting the experimentally measured receding rate into the lookup table and using linear interpolation.
8. The method for optimizing road drainage structure based on fluid dynamics simulation according to claim 1, characterized in that, S6 includes: Each monitoring point is equipped with an ultrasonic water level gauge. The water level gauge probe is pointed downwards at the water line on the road surface. Recording will start automatically when the water depth exceeds the preset recording threshold and stop after the water depth remains below the preset recording threshold for a preset time. Each monitoring point is equipped with video surveillance equipment, set to record event-triggered recording. When the water level gauge detects that the water depth exceeds the preset recording threshold, the camera will automatically start recording. The monitoring equipment transmits data to the monitoring center server via a wireless network. The system sets an early warning threshold, and automatically sends an early warning message when the water depth exceeds the threshold.
9. The method for optimizing road drainage structure based on fluid dynamics simulation according to claim 1, characterized in that, S6 further includes: Extract the maximum water depth and duration of water accumulation during rainfall events and compare them with the values under similar historical rainfall conditions before the renovation. Social benefit assessment was conducted based on traffic flow monitoring data, which showed average vehicle speed and travel time during rainy weather. Economic benefit assessment is conducted by calculating annual cost savings based on the average daily traffic volume of roads, the proportion of rainy days, the time saved during each rainfall event, and the time value of vehicles.
10. A road drainage structure optimization system based on fluid dynamics simulation, used to execute the steps in the road drainage structure optimization method based on fluid dynamics simulation as described in any one of claims 1-9, characterized in that, include: The 3D modeling module is used to collect geographic information of the target road section and uses UAV oblique photography and point cloud data processing methods to obtain a refined 3D digital model of the road surface. The boundary condition generation module is used to acquire historical rainfall data and use frequency analysis and design rainfall pattern generation methods to obtain a set of rainfall boundary conditions; The simulation calculation module is used to obtain drainage performance indicators based on the refined three-dimensional digital model of the road surface and the set of rainfall boundary conditions, using the two-dimensional shallow water wave equation simulation method. The optimization solution module is used to obtain drainage structure parameter schemes based on drainage performance indicators, using a neural network surrogate model and an improved genetic algorithm optimization method. The parameter calibration module is used to obtain the calibrated simulation model parameters based on the drainage structure parameter scheme, using field tests and parameter calibration methods. The effect verification module is used to obtain the engineering effect verification results based on the calibrated simulation model parameters and drainage structure parameter scheme, using engineering implementation and online monitoring methods.