Shallow water flow simulation method, apparatus, device, and medium

Abstract

The application discloses a shallow water flow simulation method, device, equipment and medium, relates to the technical field of shallow water dynamics, and comprises the following steps: dividing a plurality of river channel units according to to-be-simulated river channel parameters and determining an initial state of water flow; acquiring full-field and interface conservation variables of a current time step based on the state in the unit, determining a time step length according to a CFL condition; calculating a global maximum wave speed and completing flux splitting and wave speed estimation; reconstructing a WENO-JS3 first flux and a QUICK second flux with a three-point template, calculating corresponding GLF-HLL numerical flux according to the estimated wave speed, and weightedly fusing the GLF-HLL numerical flux into a mixed interface flux; discretely obtaining a time derivative of the conservation variable from the mixed flux and a flow source term by finite volume, and integrating to obtain full-field conservation variables at a next time; and iteratively updating until a stop condition is met, so that the shallow water flow simulation is completed. The adaptability, simulation precision and stability of a single numerical format in the shallow water flow simulation are improved.

Description

Technical Field

[0001] This invention relates to the field of shallow water dynamics technology, and particularly to shallow water flow simulation methods, apparatus, equipment and media. Background Technology

[0002] Existing numerical simulation methods for shallow water waves mainly employ a single high-resolution scheme for solution. While WENO-type schemes can stably capture strong discontinuities such as hydraulic jumps, dam-break shocks, and wet-dry interfaces, they often introduce significant numerical dissipation in smooth flow fields, leading to a decrease in the analytical capability for fine-scale flow structures. QUICK-type schemes, although exhibiting lower numerical dissipation and better preservation of flow details, are prone to non-physical oscillations in strong discontinuities, affecting computational stability.

[0003] Furthermore, existing hybrid schemes often employ reconstruction templates of varying lengths and rely on thresholds for scheme switching, which can increase the difficulty of flux matching and conservation processing, and may introduce additional numerical errors in the switching region. Meanwhile, adaptive mesh refinement (AMR) is gradually becoming an important technique for large-scale shallow water flow simulation, but wide template schemes are prone to generating cross-level data dependencies in the boundary region between coarse and fine meshes, increasing the complexity of algorithm implementation.

[0004] Therefore, there is an urgent need to provide a shallow water wave numerical simulation method that combines low dissipation characteristics with strong discontinuity capture capability, adopts a unified narrow template structure, and is applicable to AMR extension, so as to improve the accuracy, stability, and engineering applicability of complex flow simulation. Summary of the Invention

[0005] In view of this, the purpose of this invention is to provide a method, apparatus, device, and medium for simulating shallow water flow, thereby improving the adaptability, simulation accuracy, and stability of single numerical formats in shallow water flow simulation. The specific solution is as follows: In a first aspect, this application discloses a shallow water flow simulation method, comprising: Based on the state within the channel unit at the current time step of the channel to be simulated, obtain the global and interface conserved variables at the current time step, and obtain the step size of the current time step based on the CFL stability condition and the state within the channel unit. The global wave velocity parameters are obtained by calculating the global maximum wave velocity of the global conserved variables, and the estimated wave velocity is obtained by estimating the wave velocity based on the global wave velocity parameters after flux splitting. Based on the three-point template, candidate fluxes in WENO-JS3 format corresponding to each candidate template are constructed for the interface conservation variables. Smoothing indicators are generated using the candidate fluxes. Based on the smoothing indicators, the candidate fluxes are weighted and fused to obtain the first reconstructed flux in WENO-JS3 format. The second reconstructed flux in QUICK format is obtained based on the interface conservation variables after upwind interpolation. Based on the estimated wave velocity, the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux are obtained respectively, and the first numerical flux and the second numerical flux are weighted and fused to obtain the mixed interface flux. The global conserved variables are updated using the hybrid interface flux and the step size, and the simulation is carried out at the next time step based on the updated global conserved variables until the preset simulation stopping condition is met, thus completing the shallow water flow simulation.

[0006] Optionally, the state within the channel unit includes the water depth and water velocity of the channel to be simulated; the step of calculating the global maximum wave velocity of the global conserved variables to obtain global wave velocity parameters, and estimating the estimated wave velocity based on the global wave velocity parameters after flux splitting, includes: Traverse all channel units of the river to be simulated, and calculate the maximum propagation wave velocity in the entire field based on the water depth and water velocity of each channel unit to obtain the global wave velocity parameters. The flux is obtained based on the water depth and water velocity of the channel unit, and the flux is globally split into positive and negative fluxes based on the global wave velocity parameter to decompose the flux into split fluxes including positive and negative fluxes. The left and right traveling wave velocities corresponding to the splitting flux are estimated based on the shock wave intensity correction factor to obtain the estimated wave velocity.

[0007] Optionally, obtaining the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux based on the estimated wave velocity includes: When the estimated post-wave velocity satisfies that the leftward wave velocity is greater than or equal to zero, the first numerical flux and the second numerical flux are determined by the leftward flux contribution corresponding to the first reconstructed flux and the second reconstructed flux, respectively. When the estimated post-wave velocity satisfies that the right-moving wave velocity is less than or equal to zero, the first numerical flux and the second numerical flux are determined by the right-side flux contribution corresponding to the first reconstructed flux and the second reconstructed flux, respectively. When the estimated wave velocity satisfies that the left-moving wave velocity is less than zero and the right-moving wave velocity is greater than zero, the first numerical flux and the second numerical flux are determined by the contributions of the left and right sides of the first reconstructed flux and the second reconstructed flux, respectively.

[0008] Optionally, the step of constructing candidate fluxes in WENO-JS3 format for each candidate template based on the three-point template for the interface conservation variables, generating smoothing indicators using the candidate fluxes, and performing weighted fusion of the candidate fluxes based on the smoothing indicators to obtain the first reconstructed flux in WENO-JS3 format includes: Candidate fluxes in WENO-JS3 format corresponding to each candidate template are constructed from the interface conservation variables using a three-point template. A smoothing indicator characterizing the smoothness of the local flow field is constructed based on the difference between the conserved variables, candidate fluxes, or reconstructed fluxes of adjacent computational units within the candidate template; wherein the smoothing indicator in the smooth flow field region is smaller than the smoothing indicator in the region where the flow field has discontinuities or drastic changes; Nonlinear weights are assigned according to the smoothing indicator to obtain the weight coefficients corresponding to each candidate template; The candidate fluxes are weighted and fused using the weighting coefficients to obtain the first reconstructed flux in WENO-JS3 format.

[0009] Optionally, obtaining the second reconstructed flux in QUICK format based on the interface conservation variables after upwind interpolation includes: Based on the flow direction of water flow at the interface representing the state within the river unit, the interface conservation variables of the first preset number of river units adjacent to the windward side of the interface and the second preset number of river units on the leeward side of the interface are selected for secondary windward interpolation to obtain the interpolated interface conservation variables. The second reconstructed flux in QUICK format is obtained based on the interpolated interface conservation variables.

[0010] Optionally, the weighted fusion of the first numerical flux and the second numerical flux to obtain the mixed interface flux includes: A smoothness index is constructed based on the interface conservation variables; The continuous mixing weight is determined based on the smoothness index and the preset threshold; The second numerical flux and the first numerical flux are weighted and summed based on the continuous mixing weights to obtain the mixed interface flux.

[0011] Optionally, updating the global conserved variables using the hybrid interface flux and the step size includes: The riverbed elevation in the river channel parameters is balanced and discretized to obtain the discretized topographic source term; The friction source terms in the river channel parameters are semi-implicitly discretized to obtain the discretized friction source terms; Finite volume discretization is performed on the mixed interface flux, the discretized terrain source term, and the discretized friction source term to obtain the time derivative of the conserved variable. The time derivative of the conserved variable is integrated with the step size using the third-order TVD Runge-Kutta method to obtain the global conserved variable for the next time step.

[0012] Secondly, this application discloses a shallow water flow simulation device, comprising: The step size determination module is used to obtain the global conservation variables and interface conservation variables of the current time step based on the state of the channel unit of the channel to be simulated at the current time step, and to obtain the step size of the current time step based on the CFL stability condition and the state of the channel unit. The wave velocity estimation module is used to calculate the global maximum wave velocity of the global field conserved variables to obtain global wave velocity parameters, and to estimate the wave velocity based on the global wave velocity parameters after flux splitting to obtain the estimated wave velocity. The first flux reconstruction module is used to construct candidate fluxes in WENO-JS3 format corresponding to each candidate template based on the interface conservation variables using the three-point template, generate smoothing indicators using the candidate fluxes, and perform weighted fusion on the candidate fluxes based on the smoothing indicators to obtain the first reconstructed flux in WENO-JS3 format. The second flux reconstruction module is used to obtain the second reconstructed flux in QUICK format based on the interface conservation variables after upwind interpolation. The flux mixing module is used to obtain the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux respectively according to the estimated wave velocity, and to perform weighted fusion of the first numerical flux and the second numerical flux to obtain the mixed interface flux. The flow simulation module is used to update the global conserved variables using the mixed interface flux and the step size, so as to perform the simulation at the next time step based on the updated global conserved variables, until the preset simulation stopping condition is met, thus completing the shallow water flow simulation.

[0013] Thirdly, this application discloses an electronic device, including: Memory, used to store computer programs; A processor is used to execute the computer program to implement the steps of the aforementioned disclosed shallow water flow simulation method.

[0014] Fourthly, this application discloses a computer-readable storage medium for storing a computer program; wherein, when the computer program is executed by a processor, it implements the steps of the aforementioned disclosed shallow water flow simulation method.

[0015] The beneficial effects of this application are as follows: This application obtains the global conserved variables and interface conserved variables at the current time step based on the state within the channel unit of the channel to be simulated at the current time step; it obtains the step size of the current time step based on the CFL stability condition and the state within the channel unit; it calculates the global maximum wave velocity for the global conserved variables to obtain global wave velocity parameters, and estimates the wave velocity based on the global wave velocity parameters after flux splitting; it constructs candidate fluxes in WENO-JS3 format corresponding to each candidate template for the interface conserved variables based on the three-point template, generates smoothing indicators using the candidate fluxes, and applies the smoothing indicators to the candidate fluxes. Weighted fusion is performed to obtain the first reconstructed flux in WENO-JS3 format; the second reconstructed flux in QUICK format is obtained based on the interface conservation variables after upwind interpolation; the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux are obtained based on the estimated wave velocity, and the first numerical flux and the second numerical flux are weighted and fused to obtain the mixed interface flux; the global conservation variables are updated using the mixed interface flux and the step size, and the simulation is performed on the updated global conservation variables for the next time step until the preset simulation stopping condition is met, thus completing the shallow water flow simulation. Therefore, this application divides the channel into units and determines the initial state based on channel parameters, discretizes the complex channel space into computable units, and provides initial flow conditions. It also obtains the global and interface conserved variables at the current time step, providing fundamental data for subsequent flux calculations. Furthermore, by determining the current time step based on CFL stability conditions and the state within the channel units, the numerical calculations remain stable and do not diverge. Next, the global maximum wave velocity is calculated for the global conserved variables, and flux splitting and wave velocity estimation are performed, avoiding grid-by-grid eigenvalue decomposition, reducing computational load, and improving efficiency. Finally, WENO-JS3 candidate fluxes are constructed based on a three-point template, and smoothing indicators are generated for weighted fusion. This effectively suppresses numerical oscillations in discontinuous regions and ensures stable flow simulation. Simultaneously, QUICK reconstruction flux is obtained through upwind interpolation, which can reduce numerical dissipation and preserve flow details in smooth regions. Then, based on the estimated wave velocity, the GLF-HLL numerical fluxes corresponding to the two types of reconstruction fluxes are calculated and weighted and fused, which can achieve a smooth transition between the two formats and take into account both discontinuity capture and detail preservation. Subsequently, the mixed interface flux and flow source terms are discretized by finite volume and the time derivative of the discretized conserved variables is integrated with the time derivative of the discretized conserved variables at the stated step size to obtain the full-field conserved variables for the next time step. This can accurately reflect the influence of topography and friction on the flow and improve simulation accuracy. Finally, the full-field conserved variables are updated and the process is repeated until the stopping condition is met, which can continuously iterate the simulation of the flow evolution process and achieve high-precision, high-efficiency and stable shallow water flow simulation. Attached Figure Description

[0016] To more clearly illustrate the technical solutions in the embodiments of this application or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on the provided drawings without creative effort.

[0017] Figure 1 This is a flowchart of a shallow water flow simulation method disclosed in this application; Figure 2 This is a schematic diagram of a shallow water flow simulation device disclosed in this application; Figure 3 This is a structural diagram of an electronic device disclosed in this application. Detailed Implementation

[0018] The technical solutions of the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative effort are within the scope of protection of the present invention.

[0019] Shallow water wave equations (SWEs) are core mathematical models in hydraulic engineering, estuary and coastal engineering, and ecohydraulics, used to simulate engineering scenarios such as river evolution, tidal flow, and the evolution of wet-dry interfaces. Shallow water wave flow is characterized by multi-scale coupling. 1. Strong discontinuous structures exist, such as hydraulic jumps, dam-break shock waves, and wet and dry fronts, requiring numerical schemes to be oscillatory and highly robust; 2. The presence of fine turbulent structures, such as nearshore eddies and small-scale fluctuations related to sediment transport, requires numerical schemes with low dissipation and high detail preservation. 3. Large-scale engineering simulation requires adaptive mesh refinement (AMR): This requires narrow and localized numerical templates to avoid cross-mesh dependency conflicts.

[0020] Currently, the field of gas dynamics employs a hybrid compact-WENO scheme and a single WENO / QUICK scheme in shallow water wave domains, with the specific schemes as follows: 1. Gas dynamics hybrid scheme: It adopts a five-point wide template to fuse the compact low-dissipation scheme with the WENO discontinuous capture scheme. It is only applicable to compressible fluids, not optimized for shallow water wave equations, and cannot be compatible with dry-wet interfaces and bottom slope / friction source terms. 2. Shallow water wave single WENO scheme: The third-order / fifth-order WENO scheme can stably capture hydraulic jump and dam-break shock waves, but the numerical dissipation in the smooth region is too large, distorting the microstructure of turbulence, and is not suitable for sediment transport and ecological flow simulation. 3. Shallow water wave single QUICK format: Third-order low dissipation, which can preserve turbulent details in smooth regions, but Gibbs oscillations are generated at discontinuities, which can easily lead to calculation instability and misjudgment at the wet-dry interface; 4. Traditional HLL solvers still need to perform local wave velocity estimation and Riemann state construction at each interface, which limits computational efficiency.

[0021] The current solution has the following technical problems: 1. The switching function is a hard switch: when the indicator crosses the threshold, the format suddenly switches from low dissipation to shock capture. This discontinuous switch will introduce spurious fluctuations in the solution, and the choice of threshold is sensitive to the results.

[0022] 2. The indicator design is simple: it adopts local gradient or density difference, which has limited ability to distinguish between fine-scale fluctuations and strong discontinuities in complex shallow water flow fields. It is easy to falsely trigger WENO in turbulent regions, increasing dissipation.

[0023] 3. Inconsistent templates: The compact format and the WENO format use different template lengths, which leads to flux mismatch at the mixing boundary and affects global conservation.

[0024] 4. AMR compatibility not considered: Ren's solution is not designed for adaptive meshes. When AMR is used, the problem of incomplete templates for coarse and fine mesh interfaces will be further aggravated.

[0025] However, solutions that use WENO-JS3 or QUICK alone have inherent problems such as excessive heat dissipation or inability to handle interruptions.

[0026] Therefore, this application provides a shallow water flow simulation scheme to improve the adaptability, simulation accuracy, and stability of a single numerical format in shallow water flow simulation.

[0027] See Figure 1 As shown in the figure, this application discloses a shallow water flow simulation method, including: Step S11: Based on the state within the channel unit at the current time step of the channel to be simulated, obtain the global conservation variables and interface conservation variables at the current time step, and obtain the step size of the current time step based on the CFL stability condition and the state within the channel unit.

[0028] In this embodiment, the state within the channel unit includes the water depth and water velocity of the channel to be simulated. It is understood that before iteration, each channel unit needs to be assigned an initial flow state, i.e., an initial water depth h and an initial flow velocity (u, v), and the initial flow state of the channel is determined as the state within the channel unit at the current time step.

[0029] A minimum threshold ε is pre-set (e.g., 1e-6 meters). All channel units are traversed and the water depth of the current unit is determined. When the water depth of a unit is less than the minimum threshold ε, the channel unit is marked as a dry unit. At the same time, the velocity component of the water flow in the unit is forced to be set to zero. This avoids problems such as zero denominator, numerical overflow, or calculation divergence at small water depths in subsequent flux calculation, wave velocity solution, source term discretization, and conservation variable update. It ensures the computational robustness near the wet-dry interface and the continuity of the global simulation, and is suitable for numerical solution requirements of complex shallow water flow scenarios such as dam breaks, tides, and alternating wet and dry conditions in the channel.

[0030] Using the shallow water wave equation as the governing equation, the law of water flow conservation is described: ; ; ; In the formula, t is time; x and y are spatial coordinates; U is the conserved variable vector; F and G are the flux vectors in the x and y directions, respectively; For friction source term; This represents the topographic source term. h is the water depth; u and v represent the velocity components in the x and y directions, respectively; g is the acceleration due to gravity, taken as 9.81 m / s². 2 ; These are the riverbed shear stresses in the x and y directions, respectively; It is density; It is the riverbed elevation. The expression is as follows: ; Based on the geometric parameters (length, width, riverbed elevation) and hydraulic parameters (Manning roughness coefficient n, cross-sectional morphology) of the river channel to be simulated, the continuous river channel is discretized into a series of non-overlapping, full-coverage structured Cartesian grid elements (i.e., channel elements). The specific formula is as follows: ; In the formula, Indicates the grid index in the x / y direction. Indicates the grid step size. This represents the interface flux in the x-direction.

[0031] Each cell serves as the basic control volume for finite volume calculations, defining the cell center coordinates and grid step size. In other words, the continuous shallow water equations are transformed into discretely solvable finite volume equations, providing a computational grid foundation for subsequent flux calculations and source term discretization.

[0032] Based on the state within the channel unit, obtain the global and interface conserved variables for the current time step, and based on the CFL stability condition and the state within the channel unit, obtain the step size for the current time step.

[0033] Based on the state within each channel unit at the current time step, the conserved variable vector U of the shallow water wave equation is constructed, and the formula is: Where h is the water depth of the cell, and hu and hv are the unit width momentum in the x and y directions, respectively. By traversing all channel cells, the conserved variables of each cell are assembled according to the grid index to obtain the global conserved variables covering the entire computational domain. This unifies the expression of physical quantities in the flow field, providing a globally unified variable basis for subsequent flux calculation, wave velocity solution, and time iteration, and satisfying the conservation requirements of the finite volume method.

[0034] Based on the global conserved variables of adjacent channel units, state interpolation is performed on each grid interface (e.g., x-direction interface (i+1 / 2, j) and y-direction interface (i, j+1 / 2)) to obtain the conserved variables on the left side of the interface. and right-hand conserved variables That is, the interface conserved variable.

[0035] Based on the current state within the channel unit (water depth, flow velocity), calculate the maximum allowable step size for the current time step according to the CFL (Courant-Friedrichs-Lewy condition) stability criteria. The formula is as follows: ; In the formula, This represents the entire computational threshold, i.e., the entire region covered by the numerical simulation. This represents the acceleration due to gravity, taken as 9.81 m / s². 2 , Indicates water depth. Represents the velocity component in the x-direction. This indicates the grid step size.

[0036] Step S12: Calculate the global maximum wave velocity for the global conserved variables to obtain the global wave velocity parameters, and estimate the wave velocity based on the global wave velocity parameters after flux splitting to obtain the estimated wave velocity.

[0037] In this embodiment, the step of calculating the global maximum wave velocity of the global conserved variables to obtain global wave velocity parameters, and estimating the wave velocity based on the global wave velocity parameters after flux splitting to obtain the estimated wave velocity, includes: traversing all channel units of the river to be simulated, calculating the maximum propagation wave velocity across the entire field based on the water depth and water velocity of each channel unit to obtain global wave velocity parameters; obtaining the flux based on the water depth and water velocity of the channel units, and performing global positive and negative splitting on the flux based on the global wave velocity parameters to decompose the flux into split fluxes including positive and negative fluxes; estimating the left and right traveling wave velocities corresponding to the split fluxes based on the shock wave intensity correction factor to obtain the estimated wave velocity.

[0038] This step is the core process for constructing the HLL Riemann solver improved by global Lax-Friedrichs flux splitting. It aims to simplify interface flux calculations and ensure shock wave capture stability by unifying the wave velocity across the entire field. Specifically, the process involves: traversing all channel cells and calculating the maximum propagating wave velocity within the entire computational domain Ω based on the current water depth h and flow velocity u of each cell against the conserved variable U. ,Right now This is used to determine the global wave velocity parameters and avoid grid-by-grid eigenvalue decomposition.

[0039] Subsequently, physical fluxes were calculated based on the water depth and flow velocity of each channel unit, and global wave velocity parameters were utilized. Perform a global Lax-Friedrichs flux split, decomposing the original flux into positive and negative fluxes to obtain the split flux, the expression of which is: ; in, Represents a vector of conserved variables. Represents the physical flux vector. Represents the global wave speed parameter. and These represent positive flux and negative flux, respectively.

[0040] Subsequently, based on the flow field states on the left and right sides of the interface, a shock wave intensity correction factor q was introduced to estimate the traveling wave velocities on the left and right sides, resulting in the estimated wave velocities: ; in, , These represent the velocities of the traveling waves on the left and right sides of the interface, respectively. , These represent the velocity components on the left and right sides of the interface, respectively. Indicates the wave speed in shallow water. This represents the left and right shock wave intensity correction factors. Through the above processing, a unified wave velocity benchmark is provided for subsequent GLF-HLL interface flux calculations based on the first reconstructed flux in the WENO-JS3 format and the second reconstructed flux in the QUICK format, improving the stability and computational efficiency of complex shallow water flow simulations.

[0041] Step S13: Based on the three-point template, construct candidate fluxes in WENO-JS3 format corresponding to each candidate template for the interface conservation variables, generate smoothing indicators using the candidate fluxes, and perform weighted fusion on the candidate fluxes based on the smoothing indicators to obtain the first reconstructed flux in WENO-JS3 format.

[0042] In this embodiment, the process of constructing candidate fluxes in WENO-JS3 format for each candidate template based on the interface conservation variables using a three-point template, generating smoothing indicators using the candidate fluxes, and weighted fusing the candidate fluxes based on the smoothing indicators to obtain a first reconstructed flux in WENO-JS3 format includes: constructing candidate fluxes in WENO-JS3 format for each candidate template from the interface conservation variables using a three-point template; constructing smoothing indicators characterizing the smoothness of the local flow field based on the difference between the conservation variables, candidate fluxes, or reconstructed fluxes of adjacent computational units within the candidate templates; wherein the smoothing indicators in the smooth flow field region are smaller than the smoothing indicators in the region where the flow field has discontinuities or drastic changes; assigning nonlinear weights according to the smoothing indicators to obtain weight coefficients corresponding to each candidate template; and weighted fusing each candidate flux using the weight coefficients to obtain a first reconstructed flux in WENO-JS3 format.

[0043] This embodiment describes the process of high-precision, non-oscillating flux reconstruction based on a three-point compact template in the third-order WENO-JS3 format. The aim is to adapt to adaptive mesh refinement and accurately distinguish between smooth and discontinuous flow fields. The first reconstructed fluxes of the WENO-JS3 format in the x and y directions are obtained respectively. The following explanation uses the x-direction as an example to illustrate the process of obtaining the first reconstructed flux of the WENO-JS3 format: Using a three-point template, adjacent channel unit data are selected from the interface-conserved variables to construct candidate fluxes for each candidate template corresponding to the WENO-JS3 format. , The specific formula is as follows: ; This three-point template is a narrow template adapted to adaptive grid encryption design, which can avoid cross-grid data dependencies and ensure the numerical continuity of the interface between coarse and fine grids.

[0044] Subsequently, the sum of squared differences in candidate fluxes between adjacent channel cells is calculated to generate a smoothing indicator for quantifying the smoothness of the local flow field. , The specific formula is as follows: ; The smoothness indicator has a small value in the smooth flow field region, and a significant increase in the value in the flow field region with discontinuous or drastic changes such as hydraulic jump, dam break shock wave, dry-wet interface, etc., which can accurately identify the flow field characteristics.

[0045] Next, nonlinear weights are assigned to each candidate template based on the smoothing indicator, and the weight coefficients corresponding to each candidate template are calculated. The specific formula is shown below: ; In the formula, This represents the weight coefficient corresponding to each candidate template. Indicates a smoothing indicator. for , In this way, the weights of smooth regions are biased towards the higher-order precision template, and the weights of discontinuous regions are biased towards the stable, non-oscillating template.

[0046] Finally, the candidate fluxes are weighted and fused using the weight coefficients corresponding to each candidate template to obtain the first reconstructed flux in the WENO-JS3 format, which combines third-order accuracy in the smooth region with non-oscillatory characteristics in the discontinuous region. ,Right now This provides a high-precision and robust shock wave capture flux component for subsequent hybrid interface flux calculation.

[0047] Step S14: Obtain the second reconstructed flux in QUICK format based on the interface conservation variables after upwind interpolation.

[0048] In this embodiment, obtaining the second reconstructed flux in QUICK format based on the interface conservation variables after windward interpolation includes: selecting the interface conservation variables of a first preset number of adjacent channel units on the windward side of the interface and a second preset number of channel units on the leeward side of the interface according to the flow direction of the water flow at the interface representing the state within the channel unit, performing secondary windward interpolation on the interface conservation variables to obtain the interpolated interface conservation variables; and obtaining the second reconstructed flux in QUICK format based on the interpolated interface conservation variables.

[0049] Specifically, a three-point reconstruction template is constructed by selecting a first preset number of channel units on the windward side and a second preset number of channel units on the leeward side, based on the flow direction at the interface. The first preset number is 2, and the second preset number is 1. Subsequently, the interface conserved variables are interpolated twice on the windward side using the three-point template to construct the interface reconstruction state, and the second reconstruction flux in QUICK format is further calculated. Since a unified three-point narrow template is used throughout the process, it maintains third-order accuracy and low dissipation characteristics while adhering to the template structure consistent with the WENO-JS3 format, improving adaptability to adaptive grid densification environments.

[0050] This embodiment describes the core process of third-order low-dissipation QUICK flux reconstruction based on a three-point template. The aim is to effectively reduce numerical dissipation and finely preserve the microscopic flow structure of turbulence in smooth flow fields.

[0051] Based on the current state within the channel unit, the direction of water flow at the interface is determined. For the interface (i+1 / 2,j) in the x-direction, when the water flow at the interface is in the positive direction (from left to right), the interface conservation variables are selected from two adjacent channel units on the upstream side (i.e., the windward side) and one channel unit on the downstream side (i.e., the leeward side). A second windward interpolation is performed based on a three-point template to obtain the interface conservation variables after positive interpolation. The specific formula is as follows: ; When the water flow at the interface flows in the opposite direction (from right to left), the interface conservation variables are selected from two adjacent channel units upstream of the interface and one adjacent channel unit downstream of the interface. Similarly, based on the three-point template, a second upwind interpolation is performed to obtain the specific formula for the interface conservation variables after reverse interpolation: ; By using the above-mentioned secondary upwind interpolation of the flow direction, the QUICK scheme can be guaranteed to have third-order accuracy and low dissipation characteristics in a smooth flow field, avoiding the smoothing of the micro-flow structure by numerical dissipation.

[0052] Finally, based on the interface conservation variables obtained after forward or backward interpolation, the QUICK format interface reconstruction state is constructed, with the following formula: ; The physical flux corresponding to the interface reconstruction state is calculated to obtain the second reconstruction flux in QUICK format. The specific formula is as follows: ; The second reconstructed flux serves as the input for subsequent GLF-HLL interface numerical flux calculations. It characterizes the low-dissipation flux features of the smooth flow field region and complements the first reconstructed flux of the WENO-JS3 format, thus balancing the numerical simulation requirements of both shallow-water smooth flow and strong discontinuous flow.

[0053] It should be noted that replacing WENO-JS3 with WENO-Z3 and QUICK with the third-order upwind compact scheme can still achieve mixed solutions.

[0054] Step S15: Based on the estimated wave velocity, obtain the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux, respectively, and perform weighted fusion of the first numerical flux and the second numerical flux to obtain the mixed interface flux.

[0055] In this embodiment, obtaining the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux based on the estimated wave velocity includes: when the estimated wave velocity satisfies a leftward wave velocity greater than or equal to zero, determining the first numerical flux and the second numerical flux using the left-side flux contributions corresponding to the first reconstructed flux and the second reconstructed flux, respectively; when the estimated wave velocity satisfies a rightward wave velocity less than or equal to zero, determining the first numerical flux and the second numerical flux using the right-side flux contributions corresponding to the first reconstructed flux and the second reconstructed flux, respectively; and when the estimated wave velocity satisfies a leftward wave velocity less than zero and a rightward wave velocity greater than zero, jointly determining the first numerical flux and the second numerical flux using the left-side flux contributions corresponding to the first reconstructed flux and the second reconstructed flux, respectively.

[0056] This embodiment uses a global Lax-Friedrichs improved HLL Riemann solver and estimated wave velocity to construct corresponding GLF-HLL numerical fluxes based on the first reconstruction flux of the WENO-JS3 format and the second reconstruction flux of the QUICK format, respectively, to obtain the first and second numerical fluxes in the GLF-HLL format. The specific formulas are as follows: ; In the formula, , These represent the estimated wave velocities of the left and right traveling waves, respectively. This represents the flux contribution on the left side of the interface obtained by reconstructing the GLF forward splitting flux at cell i using the corresponding reconstruction format. This represents the flux contribution on the right side of the interface obtained by reconstructing the negative GLF split flux at cell i+1 using the corresponding reconstruction format.

[0057] Specifically, the direction of wave propagation is determined by estimating the wave velocity of the left and right traveling waves obtained earlier.

[0058] when When this occurs, it indicates that fluctuations on both sides of the interface propagate to the right, and the numerical flux at the interface is determined by the positive flux contribution from the left side. Therefore, we take... As a numerical flux of GLF-HLL; when When this occurs, it indicates that fluctuations on both sides of the interface propagate to the left, and the numerical flux at the interface is determined by the negative flux contribution from the right side. Therefore, we take... As a numerical flux of GLF-HLL; when This indicates the presence of a cross-wave propagation structure at the interface. The numerical flux at the interface is determined by the combined flux contributions from both the left and right sides. Therefore, the GLF-HLL numerical flux is constructed based on the combined flux contributions from both sides. In this embodiment, the following method is used. The calculation is performed in the form of .

[0059] The entire process employs global wave velocity simplification calculations, eliminating the need for interface-by-interface feature decomposition. It balances strong discontinuous robustness with computational efficiency while maintaining compatibility with three-point compact templates and adaptive mesh refinement.

[0060] In this embodiment, the weighted fusion of the first numerical flux and the second numerical flux to obtain the mixed interface flux includes: constructing a smoothness index based on the interface conservation variable; determining a continuous mixing weight based on the smoothness index and a preset threshold; and performing a weighted summation of the second numerical flux and the first numerical flux based on the continuous mixing weight to obtain the mixed interface flux.

[0061] This embodiment describes the process of smoothly switching between the QUICK low-loss format and the WENO-JS3 shock capture format in the x and y directions to obtain the hybrid interface flux. The specific execution process is described below using the x direction as an example (the y direction is similar).

[0062] First, the high-order smooth indicator ξ is calculated using the sixth-order central difference method for each component of the interface conserved variables, water depth and unit width momentum. The specific formula is as follows: ; In the formula, Represents the conserved variable components. This smoothness indicator can accurately distinguish between smooth and discontinuous regions of the flow field. The value is extremely small in the smooth turbulent region and significantly increases in discontinuous regions such as hydraulic jump, dam break shock wave, and dry and wet fronts.

[0063] Understandably, the sixth-order smoothness indicator can be replaced with a fourth-order / fifth-order finite difference indicator, and the smoothness of the flow field can still be determined.

[0064] The smoothing indicator is then logarithmically calculated with a preset threshold to generate a mixed weight that varies continuously between 0 and 1. The specific formula is as follows: ; In the formula, θ is a user-specified threshold parameter, with a typical value range of... ; Alternatively, the log-threshold function can be replaced with the sigmoid activation function, which can still achieve smooth subformat switching. That is, the sigmoid activation function and the smooth indicator are used to perform smooth format switching to generate continuous mixed weights.

[0065] Finally, based on this continuous hybrid weighting, the second numerical flux corresponding to the QUICK format and the first numerical flux corresponding to the WENO-JS3 format are weighted and summed to obtain a hybrid interface flux that combines the low dissipation characteristics of the smooth region with the non-oscillating robustness of the discontinuous region. In other words, when performing a weighted summation, the following situations exist: When ξ << θ Since the value is a large negative number, it is limited to 1 by min, β≈1 → QUICK is used; When ξ >> θ Since it is a positive number, it is restricted to 0 by max, β≈0 → use WENO-JS3; When ξ≈θ, β changes continuously between 0 and 1, achieving a smooth transition.

[0066] The specific formula is as follows: ; The entire process uses a three-point compact template, which is naturally compatible with adaptive mesh densification AMR and adaptable to complex flow scenarios in shallow water wave equations.

[0067] Step S16: Update the global conserved variables using the hybrid interface flux and the step size, and perform the simulation for the next time step based on the updated global conserved variables until the preset simulation stop condition is met, thus completing the shallow water flow simulation.

[0068] In this embodiment, updating the global conserved variables using the mixed interface flux and the step size includes: performing equilibrium discretization on the riverbed elevation in the channel parameters to obtain the discretized topographic source term; performing semi-implicit discretization on the friction source term in the channel parameters to obtain the discretized friction source term; performing finite volume discretization on the mixed interface flux, the discretized topographic source term, and the discretized friction source term to obtain the time derivative of the conserved variables; and integrating the time derivative of the conserved variables using the third-order TVD Runge-Kutta method with the step size to obtain the global conserved variables for the next time step.

[0069] The flux at the mixing interface and the flow source terms obtained based on the channel parameters are discretized by finite volume to obtain the time derivatives of the conserved variables. The time derivatives of the conserved variables are then integrated with the step size to obtain the global conserved variables for the next time step.

[0070] First, the riverbed elevation in the river channel parameters is processed using the Xing & Shu equilibrium discretization method to obtain the discretized topographic source term that satisfies the still water harmony, avoiding spurious flow. The splitting formula is shown below: ; In the formula, , Indicates the elevation of the riverbed after the split. This represents the maximum propagation wave speed, also known as the global wave speed parameter. This indicates the step size of the current time step.

[0071] The friction source term in the channel parameters is discretized using a semi-implicit scheme to obtain the discretized friction source term, which improves the calculation stability in shallow water depth regions. The specific formula is shown below: ; In the formula, Indicates water depth. , These represent the velocity components in the x and y directions, respectively. This represents the acceleration due to gravity, taken as 9.81 m / s². 2 , where n represents the Manning roughness coefficient. This represents the frictional attenuation coefficient per unit momentum.

[0072] Subsequently, the mixed interface flux, the discretized topographic source term, and the discretized friction source term are substituted into the semi-discrete scheme of the finite volume method for spatial discretization, yielding the time derivative of the conserved variable, as shown in the following formula: ; In the formula, Indicates the grid index in the x / y direction. Indicates the grid step size. This represents the interface flux in the x-direction.

[0073] Finally, using the current time step size, the time derivative of the conserved variables is integrated over time using the third-order TVD Runge-Kutta method to obtain stable and highly accurate full-field conserved variables for the next time step. The specific formula is as follows: ; This provides an updated flow field state for the next iteration.

[0074] Understandably, RK3-TVD can also be replaced with RK2-TVD / fourth-order Runge-Kutta scheme, and time integration can still be completed.

[0075] The global conserved variables for the next time step, obtained through third-order TVD Runge-Kutta time integration, are updated to the global conserved variables for the current time step, completing the single-step flow field iteration. The step size calculation, flux reconstruction, GLF-HLL interface flux calculation, mixed flux construction, and finite volume discretization processes are then re-executed until the preset simulation stopping conditions are met. These preset stopping conditions include the simulation duration reaching a set value, calculation results converging, the number of iteration steps reaching a preset upper limit, or the dry-wet interface evolution becoming stable. After the stopping conditions are met, the global water depth, velocity, and momentum distribution results are output, completing the shallow water flow simulation.

[0076] The beneficial effects of this application are as follows: This application obtains the global conserved variables and interface conserved variables at the current time step based on the state within the channel unit of the channel to be simulated at the current time step; it obtains the step size of the current time step based on the CFL stability condition and the state within the channel unit; it calculates the global maximum wave velocity for the global conserved variables to obtain global wave velocity parameters, and estimates the wave velocity based on the global wave velocity parameters after flux splitting; it constructs candidate fluxes in WENO-JS3 format corresponding to each candidate template for the interface conserved variables based on the three-point template, generates smoothing indicators using the candidate fluxes, and applies the smoothing indicators to the candidate fluxes. Weighted fusion is performed to obtain the first reconstructed flux in WENO-JS3 format; the second reconstructed flux in QUICK format is obtained based on the interface conservation variables after upwind interpolation; the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux are obtained based on the estimated wave velocity, and the first numerical flux and the second numerical flux are weighted and fused to obtain the mixed interface flux; the global conservation variables are updated using the mixed interface flux and the step size, and the simulation is performed on the updated global conservation variables for the next time step until the preset simulation stopping condition is met, thus completing the shallow water flow simulation. Therefore, this application divides the channel into units and determines the initial state based on channel parameters, discretizes the complex channel space into computable units, and provides initial flow conditions. It also obtains the global and interface conserved variables at the current time step, providing fundamental data for subsequent flux calculations. Furthermore, by determining the current time step based on CFL stability conditions and the state within the channel units, the numerical calculations remain stable and do not diverge. Next, the global maximum wave velocity is calculated for the global conserved variables, and flux splitting and wave velocity estimation are performed, avoiding grid-by-grid eigenvalue decomposition, reducing computational load, and improving efficiency. Finally, WENO-JS3 candidate fluxes are constructed based on a three-point template, and smoothing indicators are generated for weighted fusion. This effectively suppresses numerical oscillations in discontinuous regions and ensures stable flow simulation. Simultaneously, QUICK reconstruction flux is obtained through upwind interpolation, which can reduce numerical dissipation and preserve flow details in smooth regions. Then, based on the estimated wave velocity, the GLF-HLL numerical fluxes corresponding to the two types of reconstruction fluxes are calculated and weighted and fused, which can achieve a smooth transition between the two formats and take into account both discontinuity capture and detail preservation. Subsequently, the mixed interface flux and flow source terms are discretized by finite volume and the time derivative of the discretized conserved variables is integrated with the time derivative of the discretized conserved variables at the stated step size to obtain the full-field conserved variables for the next time step. This can accurately reflect the influence of topography and friction on the flow and improve simulation accuracy. Finally, the full-field conserved variables are updated and the process is repeated until the stopping condition is met, which can continuously iterate the simulation of the flow evolution process and achieve high-precision, high-efficiency and stable shallow water flow simulation.

[0077] The following section provides a detailed explanation using a simulation of shallow water flow in a strongly discontinuous flow scenario within a two-dimensional circular dam breach.

[0078] Construct a two-dimensional computational region as A circular dam is set at the center of the area, with the center coordinates being... ,radius The computational region is divided into Each river unit.

[0079] Initially, the water depth inside the dam was... The water depth outside the dam The initial flow velocity of each unit is .

[0080] After the simulation begins, the global conserved variables are first constructed based on the water depth and velocity in each channel unit, and the current time step is determined according to the CFL condition (CFL number is 0.1).

[0081] Then, calculate the maximum propagation wave speed across the entire field. Furthermore, a global Lax-Friedrichs flux splitting method was performed on the physical flux at the center of each unit to obtain positive and negative split flux components. Simultaneously, the traveling wave velocities on the left and right sides of the interface were obtained using a wave velocity estimation method. and .

[0082] For the WENO-JS3 sub-format: a three-point template is used to reconstruct the positive splitting flux and the negative splitting flux respectively, resulting in the flux contribution on the left side of the interface. and flux contribution on the right side of the interface .according to The symbol combination yields the numerical flux of the WENO-JS3 subformat.

[0083] ; For the QUICK sub-format: the same split flux components are used, but QUICK interpolation reconstruction is performed using a three-point upwind scheme to obtain the left-side contribution. and contributions on the right Similarly, the numerical flux of the QUICK sub-scheme is obtained by combining the wave velocity rules described above. .

[0084] Then, a higher-order smoothness index is constructed to generate continuous mixed weights. The numerical fluxes of the two sub-formats are weighted and fused to obtain the hybrid interface flux: ; The flux at the mixed interface and the source term are discretized together in a finite volume to obtain the time derivatives of the conserved variables. The third-order TVD Runge-Kutta method is then used to advance the time and update the global conserved variables to the next time step.

[0085] Repeat the above steps until the preset simulation termination time is reached, and finally obtain the spatiotemporal evolution results of the water depth field, velocity field and momentum field during the two-dimensional circular dam break.

[0086] The above embodiments demonstrate that this application can achieve stable numerical simulation of complex shallow water flow processes such as dam break shock wave propagation, while taking into account both the ability to capture discontinuities and the low dissipation characteristics of smooth regions.

[0087] The following section uses the simulation of shallow water flow in a parabolic terrain free surface oscillation scenario with friction as an example to elaborate on this application. To verify the ability of the shallow water flow simulation method of this application to handle terrain source terms, friction source terms, and wet-dry interface motion, a parabolic terrain free surface oscillation problem is used for simulation.

[0088] Construct a one-dimensional computational domain whose riverbed elevation satisfies: ; Where h0 is the height of the parabolic terrain feature; a is the horizontal scale of the terrain; and x is the river channel coordinate.

[0089] The initial free surface satisfies the conditions given by the analytical solution, and friction parameters are set. .

[0090] Based on the riverbed elevation, free surface elevation, and flow velocity information, the state within the river channel unit is constructed, and the overall conservation variables and interface conservation variables are further obtained.

[0091] Then, the time step is calculated based on the CFL stability condition, and the following steps are performed: 1. Calculate the maximum propagation velocity across the entire field and perform GLF flux splitting; 2. Estimate the velocities of the left and right traveling waves; 3. Construct the first reconstructed flux using the WENO-JS3 format; 4. Construct the second reconstruction flux using the QUICK format; 5. Calculate the corresponding numerical fluxes based on the GLF-HLL format; 6. Obtain the flux of the hybrid interface using continuous hybrid weights.

[0092] In the source term processing, the riverbed elevation term is balanced and discretized, the friction source term is semi-implicitly discretized, and the two terms are substituted together with the mixing interface flux into the finite volume discretization scheme to obtain the time derivative of the conserved variable.

[0093] Then, the third-order TVD Runge-Kutta method is used to advance the time, update the conserved variables for the next time step, and repeat the above steps until the preset simulation termination time is reached.

[0094] For friction parameters In this case, the analytical solution for the free surface is: ; The analytical solution for the flow velocity is: ; in ; Simulation results show that this application can accurately maintain the free surface oscillation characteristics under parabolic terrain conditions and stably handle the amplitude decay process caused by frictional dissipation. As time progresses, the free surface gradually approaches hydrostatic equilibrium, and the flow velocity gradually decreases to near zero, maintaining good consistency with the theoretical analytical solution.

[0095] Furthermore, the riverbed elevation balance discretization method used in this application can effectively maintain the numerical balance relationship between the source term and the flux term, and the semi-implicit discretization method of the friction source term can improve the computational stability, thereby ensuring the computational accuracy and robustness in the simulation process of complex shallow water flow.

[0096] See Figure 2 As shown in the figure, this application discloses a shallow water flow simulation device, including: The step size determination module 11 is used to obtain the global conservation variables and interface conservation variables of the current time step based on the state of the channel unit of the channel to be simulated at the current time step, and to obtain the step size of the current time step based on the CFL stability condition and the state of the channel unit. The wave velocity estimation module 12 is used to calculate the global maximum wave velocity of the global field conserved variables to obtain global wave velocity parameters, and to estimate the wave velocity based on the global wave velocity parameters after flux splitting to obtain the estimated wave velocity. The first flux reconstruction module 13 is used to construct candidate fluxes in WENO-JS3 format corresponding to each candidate template based on the interface conservation variables using the three-point template, generate smoothing indicators using the candidate fluxes, and perform weighted fusion on the candidate fluxes based on the smoothing indicators to obtain the first reconstructed flux in WENO-JS3 format. The second flux reconstruction module 14 is used to obtain the second reconstructed flux in QUICK format based on the interface conservation variables after upwind interpolation. The flux mixing module 15 is used to obtain the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux respectively according to the estimated wave velocity, and to perform weighted fusion of the first numerical flux and the second numerical flux to obtain the mixed interface flux. The flow simulation module 16 is used to update the global conserved variables using the mixed interface flux and the step size, so as to perform the simulation of the next time step based on the updated global conserved variables, until the preset simulation stop condition is met, and the shallow water flow simulation is completed.

[0097] The shallow water flow simulation device of this application is based on a three-point base template, dynamic weighted fusion, and an improved GLF-HLL flux solution. It is fully adapted to the shallow water wave equation and implemented in five steps: spatial discretization, flux calculation, sub-scheme reconstruction, dynamic mixing, and time integration. Among them, spatial discretization is based on finite volume discretization of the shallow water wave conservation equation using a uniform Cartesian grid, deriving a semi-discrete scheme of the element average conservation variables, and transforming the partial differential equation into an ordinary differential equation composed of element interface flux and source terms, laying the foundation for subsequent high-order reconstruction and time progression. Flux calculation adopts a global Lax-Friedrichs split improved HLL Riemann solver, first traversing the entire field to obtain the full... The maximum wave velocity is determined first, then global flux splitting and wave velocity estimation are performed, deriving a GLF-HLL flux expression that does not require element-by-element eigenvalue decomposition, balancing computational efficiency and strong discontinuity robustness. Sub-format reconstruction is based on a unified three-point template, deriving interface state reconstruction formulas for WENO-JS3 and QUICK formats respectively, ensuring consistent accuracy and compatibility with adaptive mesh refinement. Dynamic hybrid weights are constructed based on a sixth-order central difference to create a smooth indicator, and a continuously varying hybrid weight from 0 to 1 is derived through a logarithmic function, achieving a smooth transition between the smooth QUICK region and the discontinuous WENO-JS3 region, avoiding spurious fluctuations caused by hard switching. The time format uses a third-order TVD. The Runge-Kutta method derives a three-step strongly stable time-progression formula based on a semi-discrete scheme. The step size is determined by combining the CFL condition. At the same time, the topographic source term is balanced and discretized, and the friction source term is semi-implicitly processed to ensure the harmony of still water and the stability at shallow water depths. Finally, the semi-discrete scheme builds a bridge between spatial discretization and time integration, forming a complete derivation chain from the governing equations to the numerical solution, achieving high-precision, high-efficiency, and high-stability numerical solutions for shallow water waves.

[0098] In other words, this invention employs a three-point basis narrow template throughout, eliminating cross-grid data dependencies. This allows for natural adaptation to adaptive mesh refinement (AMR) technology, ensuring no numerical solution anomalies occur at AMR refinement boundaries and maintaining good conservation properties. It can stably support large-scale shallow water wave engineering simulations. Secondly, a dynamic weighted fusion strategy achieves synergistic format performance. In smooth flow regions, the low-dissipation QUICK format is automatically activated, significantly reducing numerical dissipation and substantially improving the accuracy of turbulent microstructure and sediment transport simulations. In discontinuous flow regions, it automatically switches to the WENO-JS3 format, effectively suppressing Gibbs oscillations and ensuring stable simulations of scenarios such as wet-dry interfaces and dam-break shock waves. Furthermore, no computational instability occurs. Next, an improved GLF-HLL Riemann solver is used, replacing traditional edge-by-edge eigenvalue decomposition with global flux splitting, effectively reducing computational load and significantly lowering the time consumption of large-scale simulations. In addition, specific optimizations are carried out for shallow water wave problems, which can adapt to bottom slope and friction source terms and satisfy accurate hydrostatic harmony C-property, ensuring distortion-free steady-state simulation and compatibility with all engineering scenarios such as tidal and wet-dry interface evolution. Furthermore, the invention has overall third-order computational accuracy, and the error converges in the third order during mesh refinement, obtaining high-precision numerical solutions at lower mesh resolutions, further reducing overall computational costs.

[0099] Furthermore, embodiments of this application also provide an electronic device. Figure 3 This is a structural diagram of an electronic device 20 according to an exemplary embodiment. The content of the diagram should not be construed as limiting the scope of this application.

[0100] Figure 3 This is a schematic diagram of the structure of an electronic device provided in an embodiment of this application. Specifically, it may include: at least one processor 21, at least one memory 22, a power supply 23, a communication interface 24, an input / output interface 25, and a communication bus 26. The memory 22 stores a computer program, which is loaded and executed by the processor 21 to implement the relevant steps in the shallow water flow simulation method performed by the electronic device disclosed in any of the foregoing embodiments.

[0101] In this embodiment, the power supply 23 is used to provide operating voltage for various hardware devices on the electronic device; the communication interface 24 can create a data transmission channel between the electronic device and external devices, and the communication protocol it follows can be any communication protocol applicable to the technical solution of this application, and is not specifically limited here; the input / output interface 25 is used to acquire external input data or output data to the outside world, and its specific interface type can be selected according to specific application needs, and is not specifically limited here.

[0102] The processor 21 may include one or more processing cores, such as a quad-core processor or an octa-core processor. The processor 21 may be implemented using at least one hardware form selected from DSP (Digital Signal Processing), FPGA (Field-Programmable Gate Array), and PLA (Programmable Logic Array). The processor 21 may also include a main processor and a coprocessor. The main processor, also known as a CPU (Central Processing Unit), is used to process data in the wake-up state; the coprocessor is a low-power processor used to process data in the standby state. In some embodiments, the processor 21 may integrate a GPU (Graphics Processing Unit), which is responsible for rendering and drawing the content to be displayed on the screen. In some embodiments, the processor 21 may also include an AI (Artificial Intelligence) processor, which is used to handle computational operations related to machine learning.

[0103] In addition, the memory 22, as a carrier for resource storage, can be a read-only memory, random access memory, disk or optical disk, etc. The resources stored on it include operating system 221, computer program 222 and data 223, etc., and the storage method can be temporary storage or permanent storage.

[0104] The operating system 221 manages and controls the various hardware devices and computer programs 222 on the electronic device to enable the processor 21 to perform calculations and processing on the massive amounts of data 223 in the memory 22. The operating system can be Windows, Unix, Linux, etc. The computer program 222, in addition to including a computer program capable of performing the shallow water flow simulation method executed by the electronic device as disclosed in any of the foregoing embodiments, may further include computer programs capable of performing other specific tasks. The data 223 may include data received by the electronic device from external devices, as well as data collected by its own input / output interface 25.

[0105] Furthermore, this application also discloses a computer-readable storage medium for storing a computer program; wherein, when the computer program is executed by a processor, it implements the aforementioned shallow water flow simulation method. Specific steps of this method can be found in the corresponding content disclosed in the foregoing embodiments, and will not be repeated here.

[0106] The various embodiments in this specification are described in a progressive manner, with each embodiment focusing on its differences from other embodiments. Similar or identical parts between embodiments can be referred to interchangeably. For the apparatus disclosed in the embodiments, since it corresponds to the method disclosed in the embodiments, the description is relatively simple; relevant parts can be referred to in the method section.

[0107] Those skilled in the art will further recognize that the units and algorithm steps of the various examples described in conjunction with the embodiments disclosed herein can be implemented in electronic hardware, computer software, or a combination of both. To clearly illustrate the interchangeability of hardware and software, the components and steps of the various examples have been generally described in terms of functionality in the foregoing description. Whether these functions are implemented in hardware or software depends on the specific application and design constraints of the technical solution. Those skilled in the art can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of this application. The steps of the methods or algorithms described in conjunction with the embodiments disclosed herein can be implemented directly in hardware, software modules executed by a processor, or a combination of both. The software module may be located in random access memory (RAM), memory, read-only memory (ROM), electrically programmable programmable read-only memory (EPROM), electrically erasable programmable read-only memory (EEPROM), register, hard disk, removable disk, CD-ROM (Compact Disc Read-Only Memory), or any other form of storage medium known in the art.

[0108] Finally, it should be noted that in this document, relational terms such as "first" and "second" are used only to distinguish one entity or operation from another, and do not necessarily require or imply any such actual relationship or order between these entities or operations. Furthermore, the terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Without further limitations, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes said element.

[0109] The present invention provides a detailed description of a shallow water flow simulation method, apparatus, device, and medium. Specific examples have been used to illustrate the principles and implementation methods of the present invention. The descriptions of the above embodiments are only intended to help understand the method and core ideas of the present invention. At the same time, those skilled in the art will recognize that, based on the ideas of the present invention, there will be changes in the specific implementation methods and application scope. Therefore, the content of this specification should not be construed as a limitation of the present invention.

Claims

1. A method for simulating shallow water flow, characterized in that, include: Based on the state within the channel unit at the current time step of the channel to be simulated, obtain the global and interface conserved variables at the current time step, and obtain the step size of the current time step based on the CFL stability condition and the state within the channel unit. The global wave velocity parameters are obtained by calculating the global maximum wave velocity of the global conserved variables, and the estimated wave velocity is obtained by estimating the wave velocity based on the global wave velocity parameters after flux splitting. Based on the three-point template, candidate fluxes in WENO-JS3 format corresponding to each candidate template are constructed for the interface conservation variables. Smoothing indicators are generated using the candidate fluxes. Based on the smoothing indicators, the candidate fluxes are weighted and fused to obtain the first reconstructed flux in WENO-JS3 format. The second reconstructed flux in QUICK format is obtained based on the interface conservation variables after upwind interpolation. Based on the estimated wave velocity, the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux are obtained respectively, and the first numerical flux and the second numerical flux are weighted and fused to obtain the mixed interface flux. The global conserved variables are updated using the hybrid interface flux and the step size, and the simulation is carried out at the next time step based on the updated global conserved variables until the preset simulation stopping condition is met, thus completing the shallow water flow simulation. The state within the channel unit includes the water depth and velocity of the channel to be simulated; the calculation of the global maximum wave velocity of the global conserved variables to obtain global wave velocity parameters, and the estimation of the estimated wave velocity based on the global wave velocity parameters after flux splitting, includes: Iterate through all channel units of the river to be simulated, and calculate the maximum propagating wave velocity across the entire field based on the water depth and velocity of each channel unit, thus obtaining global wave velocity parameters. Obtain the flux based on the water depth and velocity of each channel unit, and perform global positive and negative splitting on the flux based on the global wave velocity parameters to decompose the flux into split fluxes including positive and negative fluxes. The formula for global positive and negative splitting is: ； in, Represents a vector of conserved variables. Represents the physical flux vector. Represents the global wave speed parameter. and These represent positive flux and negative flux, respectively; The left and right traveling wave velocities corresponding to the splitting flux are estimated based on the shock wave intensity correction factor to obtain the estimated wave velocity; the formula for obtaining the estimated wave velocity is: ； in, , These represent the velocities of the traveling waves on the left and right sides of the interface, respectively. , c represents the velocity components on the left and right sides of the interface, respectively. L c R These represent the gravitational wave velocities of the left and right sides in shallow water, respectively. , These represent the left shock wave intensity correction factor and the right shock wave intensity correction factor, respectively. The weighted fusion of the first numerical flux and the second numerical flux to obtain the mixed interface flux includes: A smoothness index is constructed based on the interface conservation variables; a continuous mixing weight is determined according to the smoothness index and a preset threshold; the second numerical flux and the first numerical flux are weighted and summed based on the continuous mixing weight to obtain the mixed interface flux; The step of obtaining the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux based on the estimated wave velocity includes: When the estimated wave velocity satisfies a leftward wave velocity greater than or equal to zero, the first numerical flux and the second numerical flux are determined by the left-side flux contributions corresponding to the first reconstructed flux and the second reconstructed flux, respectively; when the estimated wave velocity satisfies a rightward wave velocity less than or equal to zero, the first numerical flux and the second numerical flux are determined by the right-side flux contributions corresponding to the first reconstructed flux and the second reconstructed flux, respectively; when the estimated wave velocity satisfies a leftward wave velocity less than zero and a rightward wave velocity greater than zero, the first numerical flux and the second numerical flux are determined by jointly using the left and right-side flux contributions corresponding to the first reconstructed flux and the second reconstructed flux, respectively. The formulas for obtaining the first and second numerical fluxes are as follows: ； In the formula, , These represent the estimated wave velocity of the left-traveling wave and the estimated wave velocity of the right-traveling wave, respectively. This represents the flux contribution on the left side of the interface obtained by reconstructing the GLF forward splitting flux at cell i using the corresponding reconstruction format. This represents the flux contribution on the right side of the interface obtained by reconstructing the negative GLF split flux at cell i+1 using the corresponding reconstruction format.

2. The shallow water flow simulation method according to claim 1, characterized in that, The process involves constructing candidate fluxes in WENO-JS3 format for each candidate template based on the interface conservation variables using a three-point template, generating smoothing indicators using these candidate fluxes, and then weighting and fusing the candidate fluxes based on the smoothing indicators to obtain the first reconstructed flux in WENO-JS3 format, including: Candidate fluxes in WENO-JS3 format corresponding to each candidate template are constructed from the interface conservation variables using a three-point template. A smoothing indicator characterizing the smoothness of the local flow field is constructed based on the difference between the conserved variables, candidate fluxes, or reconstructed fluxes of adjacent computational units within the candidate template; wherein the smoothing indicator in the smooth flow field region is smaller than the smoothing indicator in the region where the flow field has discontinuities or drastic changes; Nonlinear weights are assigned according to the smoothing indicator to obtain the weight coefficients corresponding to each candidate template; The candidate fluxes are weighted and fused using the weighting coefficients to obtain the first reconstructed flux in WENO-JS3 format.

3. The shallow water flow simulation method according to claim 1, characterized in that, The step of obtaining the second reconstruction flux in QUICK format based on the interface conservation variables after upwind interpolation includes: Based on the flow direction of water flow at the interface representing the state within the river unit, the interface conservation variables of the first preset number of river units adjacent to the windward side of the interface and the second preset number of river units on the leeward side of the interface are selected for secondary windward interpolation to obtain the interpolated interface conservation variables. The second reconstructed flux in QUICK format is obtained based on the interpolated interface conservation variables.

4. The shallow water flow simulation method according to any one of claims 1 to 3, characterized in that, The update of the global conserved variables using the hybrid interface flux and the step size includes: The riverbed elevation in the river channel parameters is balanced and discretized to obtain the discretized topographic source term; The friction source terms in the river channel parameters are semi-implicitly discretized to obtain the discretized friction source terms; Finite volume discretization is performed on the mixed interface flux, the discretized terrain source term, and the discretized friction source term to obtain the time derivative of the conserved variable. The time derivative of the conserved variable is integrated with the step size using the third-order TVD Runge-Kutta method to obtain the global conserved variable for the next time step.

5. A shallow water flow simulation device, characterized in that, The steps for implementing the shallow water flow simulation method as described in any one of claims 1 to 4 include: The step size determination module is used to obtain the global conservation variables and interface conservation variables of the current time step based on the state of the channel unit of the channel to be simulated at the current time step, and to obtain the step size of the current time step based on the CFL stability condition and the state of the channel unit. The wave velocity estimation module is used to calculate the global maximum wave velocity of the global field conserved variables to obtain global wave velocity parameters, and to estimate the wave velocity based on the global wave velocity parameters after flux splitting to obtain the estimated wave velocity. The first flux reconstruction module is used to construct candidate fluxes in WENO-JS3 format corresponding to each candidate template based on the interface conservation variables using the three-point template, generate smoothing indicators using the candidate fluxes, and perform weighted fusion on the candidate fluxes based on the smoothing indicators to obtain the first reconstructed flux in WENO-JS3 format. The second flux reconstruction module is used to obtain the second reconstructed flux in QUICK format based on the interface conservation variables after upwind interpolation. The flux mixing module is used to obtain the first numerical flux and the second numerical flux corresponding to the first reconstructed flux and the second reconstructed flux respectively according to the estimated wave velocity, and to perform weighted fusion of the first numerical flux and the second numerical flux to obtain the mixed interface flux. The flow simulation module is used to update the global conserved variables using the mixed interface flux and the step size, so as to perform the simulation at the next time step based on the updated global conserved variables, until the preset simulation stopping condition is met, thus completing the shallow water flow simulation.

6. An electronic device, characterized in that, include: Memory, used to store computer programs; A processor for executing the computer program to implement the steps of the shallow water flow simulation method as described in any one of claims 1 to 4.

7. A computer-readable storage medium, characterized in that, Used to store a computer program; wherein, when the computer program is executed by a processor, it implements the steps of the shallow water flow simulation method as described in any one of claims 1 to 4.

Images