A multi-scale wave-current-sediment coupled shore evolution prediction model method

Abstract

The present invention discloses a multi-scale wave-current-sediment coupled shore evolution prediction model method, which relates to the technical field of tidal flat landform evolution simulation. The method includes integrating tidal current dynamics, wave motion, and sediment dynamics modules and solving them to predict the landform evolution of tidal flat areas; establishing a coordinate system based on measured tidal flat topography data, obtaining tidal current, wave, and sediment data based on hydrological stations and importing them into the model as boundary conditions; constructing a shore evolution prediction model based on the integration of tidal current dynamics, wave motion, and sediment dynamics modules, and setting the required parameters for each module; using the finite difference method to solve the fluid continuity equation, the wave energy conservation equation, the one-dimensional vertical average suspended sediment convection-diffusion equation, and the bottom bed renewal equation to obtain tidal flat profile topography at different stages, thereby accurately simulating the tidal flat profile evolution process under the action of waves and currents. Compared with traditional tidal current and sediment models, the simulation accuracy is higher and the computational efficiency can be significantly improved.

Description

Technical Field

[0001] The present invention relates to the technical field of tidal flat landform evolution simulation technology, in particular to a multi-scale wave-current-sediment coupled shore evolution prediction model method. Background Art

[0002] Tidal flats serve as a buffer zone between humans and the ocean. Tidal flat profile development is controlled by factors such as tides, waves, sediment, and organisms. The complex feedback mechanisms between these multiple dynamic factors complicate the prediction of tidal flat profile morphology. Under the positive feedback influence of hydrodynamics and morphology, tidal flats are in a long-term dynamic equilibrium between erosion and siltation. Numerical models play a crucial role in analyzing and predicting the evolution of tidal flat morphology.

[0003] However, such models require complex boundary conditions and input parameters, generally in two-dimensional or three-dimensional form. The simulation and prediction of long-term landform evolution usually requires a long calculation time. In addition, the near-bottom sediment flux terms in the existing sediment module control equations are calculated using the shear stress method. The suspension and sedimentation states of sediment are usually controlled by the critical scour shear stress and critical erosion shear stress, while the critical erosion shear stress and erosion rate generally need to be determined in the laboratory. The confirmation of parameters is difficult, especially for cohesive sediment. Many studies have used arbitrary and unreasonable values. Summary of the Invention

[0004] In order to solve the above technical problems, the present invention provides a multi-scale wave-current-sediment coupled shore evolution prediction model method, comprising the following steps:

[0005] S1. Based on the given offshore tidal flat topographic data, a coordinate system is established, and then the tidal, wave, and sediment conditions of the area are obtained based on the observation station;

[0006] S2. Construct a wave-current-sediment coupled shore evolution prediction model. The wave-current-sediment coupled shore evolution prediction model includes a tidal dynamics module, a wave dynamics module, and a sediment dynamics module. The tidal dynamics module is constructed based on the fluid continuity equation with zero free surface gradient. The wave propagation module is generated based on the wave energy conservation equation. Starting from the sediment entrainment force formula, the near-bottom sediment flux F required for the one-dimensional average suspended sediment convection-diffusion equation is provided. s , combined with the bottom bed renewal equation dominated by the near-bottom sediment flux, a sediment dynamics module is constructed;

[0007] S3, arranging the data obtained in step S1 and using them as boundary conditions to import into the wave-current-sediment coupled shore evolution prediction model;

[0008] S4. Use the finite difference method to iteratively solve the wave-current-sediment coupled beach evolution prediction model to obtain the dynamic evolution process of tidal flat landform morphology.

[0009] The technical solution further defined in the present invention is:

[0010] Furthermore, in step S2, the tidal power module is based on the fluid continuity equation. It is assumed that at any time, the gradient of the free surface is 0, that is, there is no oscillation of the liquid surface. The momentum equation is ignored and the control equation is as follows:

[0011]

[0012] Where u represents the average velocity in the depth of the tidal flat section, t represents time, and x represents the distance across the shore; water depth h = η - z b , η represents the water surface height, z b Indicates bed elevation;

[0013] The water surface height η function satisfies the following formula:

[0014]

[0015] Where A represents the tidal amplitude, T c represents the tidal cycle, T c =12h.

[0016] As described above, in a multi-scale wave-current-sediment coupled shore evolution prediction model method, in step S2, a wave power module is formed based on the wave energy conservation equation, which is as follows:

[0017]

[0018] Where E represents wave energy, ρ represents the density of seawater, g represents the acceleration due to gravity, H rms Represents the root mean square wave height; S tot represents the total wave energy term, c g represents the group velocity, S fric represents the wave bottom friction term, S break represents the wave breaking source term;

[0019] Group velocity c g The propagation speed of wave energy is expressed as follows based on linear wave theory:

[0020]

[0021] Where c = σ / k represents the phase velocity of the wave, σ represents the wave frequency, k represents the wave number, and the direction of the wave phase velocity is fixed to the shoreward direction;

[0022] Wave bottom friction term S fric As shown in the following formula:

[0023]

[0024] Among them, α f represents the wave energy dissipation coefficient, H represents the wave height, T w represents the wave period;

[0025] When the water depth gradient is large and the wave height exceeds a certain proportion of the water depth, the waves break, as shown in the following formula:

[0026]

[0027] Where α represents a constant of 1; Qb represents the probability of wave breaking determined by the discriminant index b. H represents wave height, the maximum wave height is H max =0.78h.

[0028] As described above, a multi-scale wave-current-sediment coupled beach evolution prediction model method is proposed. In step S2, starting from the sediment entrainment force formula, the near-bottom sediment flux F required by the one-dimensional average suspended sediment convection-diffusion equation is provided. s , and combined with the bottom bed renewal equation dominated by the near-bottom sediment flux to form a sediment dynamics module. The sediment model control equation adopts the one-dimensional vertical average suspended sediment convection-diffusion equation, which is as follows:

[0029]

[0030] F s =αω(C * -C)

[0031] Where h represents water depth, C represents vertical average sediment concentration, t represents time, u represents depth average velocity in the tidal flat section, x represents cross-shore distance, ε represents horizontal suspended sediment diffusion coefficient, and F represents the vertical average sediment concentration. s represents the suspended sediment flux near the bottom, α represents a constant of 1, ω represents the sediment settling velocity, C * Indicates the sand-carrying force;

[0032] The near-bottom sediment flux will cause scouring or silting of the bottom bed. The bed surface renewal equation is expressed by the following formula:

[0033]

[0034] Where p represents the porosity of the bed sediment, ρ s represents the sediment density, z b Indicates the bed surface elevation.

[0035] As mentioned above, a multi-scale wave-current-sediment coupled beach evolution prediction model method, the sediment carrying force C * Same dimension as vertical average sediment concentration C:

[0036]

[0037] Among them, C D represents the drag coefficient, ρ and ρ s represent the density of water and sediment respectively, U represents the dynamic factor, and g represents the acceleration due to gravity;

[0038] Drag coefficient C D It is related to the Reynolds number and the bed surface roughness, and is expressed by the Manning formula:

[0039]

[0040] The calculation formula of power factor U is:

[0041]

[0042] Among them, n represents the roughness, u c represents the tidal flow velocity, u w represents the wave velocity, f c and f w represent the tidal and wave resistance coefficients f respectively c =gn 2 / h 1 / 3 ,, ξ represents the displacement amplitude of the bottom particles, ξ=u w T w / 2π,T w represents the wave period, k s represents the Nicolaz function value;

[0043] The formula of sand holding force is summarized as follows:

[0044] As described above, in a multi-scale wave-current-sediment coupled beach evolution prediction model method, in step S3, the boundary conditions required by the wave-current-sediment coupled beach evolution prediction model are the land boundary and the seaward boundary. The land boundary adopts a closed boundary condition, and the seaward open boundary adopts different boundary conditions according to different equations.

[0045] The tidal equation adopts the Dirichlet boundary condition of flow, and the boundary physical quantity is a constant; the wave equation adopts the Dirichlet boundary condition of root mean square wave height; the sediment equation adopts the Dirichlet boundary condition at high tide and the Neumann boundary condition at low tide, and the gradient of the boundary physical quantity is a constant.

[0046] As described above, in the multi-scale wave-current-sediment coupled shore evolution prediction model method, step S4 specifically includes the following sub-steps:

[0047] S4.1. Based on the provided boundary conditions, solve the tidal dynamics module and the wave propagation module using the finite difference method to calculate the flow velocity and wave velocity distribution on the tidal flat terrain;

[0048] S4.2. Substitute the velocity and wave velocity distribution obtained from the wave dynamics module and the tidal dynamics module into the sediment entrainment force formula to calculate the sediment entrainment force. Discretely solve the average suspended sediment convection-diffusion equation and the sediment entrainment force formula to obtain the suspended sediment concentration distribution on the tidal flat terrain.

[0049] S4.3. Based on the solved suspended sediment concentration distribution, substitute it into the bottom bed update equation and use the forward Euler method to perform iterative solution to obtain the tidal flat profile topography at different simulation stages.

[0050] As described above, a multi-scale wave-current-sediment coupled shore evolution prediction model method, in step S4.1, the tidal dynamic equation includes the time derivative term Where h represents water depth, t represents time, and the time derivative term By directly taking the derivative of the water surface function with respect to time, the equation is transformed into an ordinary differential equation:

[0051]

[0052] Where x represents the cross-shore distance and q represents the flow rate per unit width, in m 2 / s, q = hu, q is discretized using the first-order forward difference, u represents the depth-averaged velocity in the tidal flat section direction; i represents the spatial scale, Δx represents the spatial step size, q i+1 Indicates the single-width flow q value at spatial position i+1, q i Represents the single-width flow q value at spatial position i; when the boundary flow q is determined bd =∫Δh(x, t)dx / Δt, the equation is directly solved by iteration, where Δh represents the change in water level and Δt represents the time step.

[0053] As described above, in a multi-scale wave-current-sediment coupled shore evolution prediction model method, in step S4.2, the wave energy conservation equation in the wave dynamics module includes the convection term The first-order forward difference method is used to solve the spatial equation, and the backward Euler method is used to implicitly solve the time term:

[0054]

[0055] Among them, c g represents the wave group velocity, E represents the wave energy, x represents the cross-shore distance, Δx represents the spatial step size, and j represents the time scale; It represents the wave energy at spatial position i and time scale j+1; It represents the wave energy at spatial position i+1 and time scale j+1; represents the wave energy at spatial position i and time scale j; represents the group velocity at spatial position i+1 and time scale j+1; represents the group velocity at spatial position i and time scale j+1; It represents the total wave energy at spatial position i and time scale j.

[0056] As described above, a multi-scale wave-current-sediment coupled beach evolution prediction model method, in step S4.3, the suspended sediment convection and diffusion equation in the sediment dynamics module includes the convection term Where h represents the water depth, u represents the depth-averaged velocity in the tidal flat section, C represents the vertical average sediment concentration, and x represents the cross-shore distance. The discretization of the convection term adopts the upwind format, that is, when u>0 at low tide, first-order backward difference is used, and when u<0 at high tide, first-order forward difference is used. Time discretization also adopts the implicit backward Euler method. The source term less than 0 is implicitly linearized by -αωC, where α represents a constant of 1 and ω represents the sediment settling velocity. The specific discretization format is as follows:

[0057] Time item:

[0058]

[0059] Diffusion term:

[0060]

[0061] Convection term:

[0062]

[0063] in, represents the water depth at spatial position i and time scale j; represents the water depth at spatial position i-1 and time scale j; represents the water depth at spatial position i+1 and time scale j; It represents the suspended sediment concentration at spatial position i and time scale j+1; represents the suspended sediment concentration at spatial position i and time scale j; It represents the suspended sediment concentration at the spatial position i-1 and the time scale j+1; It represents the suspended sediment concentration at the spatial position i+1 and the time scale j+1; represents the flow velocity at spatial position i+1 and time scale j; represents the flow velocity at spatial position i and time scale j; represents the flow velocity at spatial position i-1 and time scale j;

[0064] Substitute the suspended sediment concentration distribution into the bed update equation and use the forward Euler method to iteratively solve the problem. The discrete format is:

[0065]

[0066] in, It represents the sand-carrying force at spatial position i and time scale j; It represents the sand-carrying force at spatial position i and time scale j+1;

[0067] Obtain the simulated tidal flat profile topography and analyze the tidal flat evolution law under the action of waves and currents.

[0068] The beneficial effects of the present invention are:

[0069] (1) In terms of computational efficiency, the present invention appropriately simplifies the fluid continuity equation and wave energy conservation equation under long-period simulation. Compared with the traditional tidal sediment model Delft3D, under the same grid parameters and simulation conditions, the computation time is only 1 / 5 of that of the traditional tidal sediment model, which significantly improves the computational efficiency of the evolution of coastal tidal flat landforms.

[0070] (2) In terms of calculation accuracy, the present invention is compared with the Delft3D model. The annual terrain curves of the present invention are almost identical. Moreover, as the dynamic conditions gradually adapt to the terrain over time, the error has a significant decreasing trend. The calculation accuracy is almost consistent with that of Delft3D, and the simulation accuracy is relatively high. BRIEF DESCRIPTION OF THE DRAWINGS

[0071] Figure 1 It is a schematic diagram of the overall process of the present invention;

[0072] Figure 2 Schematic diagram of the measured terrain and the model initial terrain in an embodiment of the present invention;

[0073] Figure 3 Schematic diagram of the landform evolution of a tidal flat within 5 years when wave conditions are not considered in an embodiment of the present invention;

[0074] Figure 4 This is a schematic diagram of the landform evolution of a tidal flat over five years when wave conditions are considered and the significant wave height is 0.3 m in an embodiment of the present invention;

[0075] Figure 5 This is a schematic diagram of the change in the magnitude of the sand-carrying force of the tidal flat within 5 years when the wave conditions are not considered in the embodiment of the present invention;

[0076] Figure 6Schematic diagram comparing the results of the model method according to the embodiment of the present invention and the traditional tidal sediment model. DETAILED DESCRIPTION

[0077] This embodiment provides a multi-scale wave-current-sediment coupled shore evolution prediction model method, such as Figure 1 As shown, the following steps are included:

[0078] S1. Based on the given offshore tidal flat sea area topography data, a coordinate system is established, and then the hydrological conditions (tidal, wave) and sediment conditions of the area are obtained based on the observation station.

[0079] S2. Construct a wave-current-sediment coupled shore evolution prediction model, which is implemented by a tidal dynamics module, a wave dynamics module, and a sediment dynamics module. The tidal dynamics module is constructed based on the fluid continuity equation with zero free surface gradient. The wave propagation module is generated based on the simplified wave energy conservation equation. Starting from the general sediment entrainment force formula, the near-bottom sediment flux F required for the one-dimensional average suspended sediment convection-diffusion equation is provided. s , combined with the bottom bed renewal equation dominated by the near-bottom sediment flux, a sediment dynamics module is constructed.

[0080] The tidal dynamics module is based on the fluid continuity equation. It assumes that at any time, the gradient of the free surface is 0 (i.e., there is no oscillation of the liquid surface) and ignores the momentum equation. The governing equation is as follows:

[0081]

[0082] Where u represents the average velocity in the depth of the tidal flat section, t represents time, and x represents the distance across the shore; water depth h = η - z b , η represents the water surface height, z b Indicates the bed surface elevation.

[0083] The water surface height η function satisfies the following formula:

[0084]

[0085] Where A represents the tidal amplitude, T c represents the tidal cycle, T c =12h.

[0086] The wave power module adopts a simplified phase average model and is composed of a wave power module based on a simplified wave energy conservation equation. The wave energy conservation equation is as follows:

[0087]

[0088] Where E represents wave energy, ρ represents the density of seawater, g represents the acceleration due to gravity, H rmsRepresents the root mean square wave height; S tot represents the total wave energy term, c g represents the group velocity, S fric represents the wave bottom friction term, which is the main source of wave energy dissipation near the shore, S break represents the wave breaking source term.

[0089] Group velocity c g The propagation speed of wave energy is expressed as follows based on linear wave theory:

[0090]

[0091] Wherein, c=σ / k represents the phase velocity of the wave, σ represents the wave frequency, and k represents the wave number. The direction of the wave phase velocity in this embodiment is fixed to the shoreward direction.

[0092] Wave bottom friction term S fric As shown in the following formula:

[0093]

[0094] Among them, α f represents the wave energy dissipation coefficient, H represents the wave height, T w Represents the wave period.

[0095] When the water depth gradient is large and the wave height exceeds a certain proportion of the water depth, the waves break, as shown in the following formula:

[0096]

[0097] Where α represents a constant of 1; Qb represents the probability of wave breaking determined by the discriminant index b. H represents wave height, the maximum wave height is H max =0.78h.

[0098] Based on the general sediment entrainment formula, the near-bottom sediment flux F required for the one-dimensional average suspended sediment convection-diffusion equation is provided. s , and combined with the bottom bed renewal equation dominated by the near-bottom sediment flux to form a sediment dynamics module. The sediment model control equation adopts the one-dimensional vertical average suspended sediment convection-diffusion equation, which is as follows:

[0099]

[0100] F s =αω(C * -C)

[0101] Where h represents water depth, C represents vertical average sediment concentration, t represents time, u represents depth average velocity in the tidal flat section, x represents cross-shore distance, ε represents horizontal suspended sediment diffusion coefficient, and F represents the vertical average sediment concentration. srepresents the suspended sediment flux near the bottom, α represents a constant of 1, ω represents the sediment settling velocity, C * Indicates sand-holding force.

[0102] Sand holding capacity C * Same dimension as vertical average sediment concentration C:

[0103]

[0104] Among them, C D represents the drag coefficient, ρ and ρ s They represent the density of water and sediment respectively, U represents the dynamic factor, and g represents the acceleration due to gravity.

[0105] Drag coefficient C D It is related to the Reynolds number and the bed surface roughness, and is expressed by the Manning formula:

[0106]

[0107] The calculation formula of power factor U is:

[0108]

[0109] Among them, n represents the roughness, u c represents the tidal flow velocity, u w represents the wave velocity, f c and f w represent the tidal and wave resistance coefficients f respectively c =gn 2 / h 1 / 3 ,, ξ represents the displacement amplitude of the bottom particles, ξ=u w T w / 2π,T w represents the wave period, k s represents the Nicolaz function value;

[0110] The formula of sand holding force is summarized as follows:

[0111] The near-bottom sediment flux will cause scouring or silting of the bottom bed. The bed surface renewal equation is expressed by the following formula:

[0112]

[0113] Where p represents the porosity of the bed sediment, ρ s represents the sediment density, z b Indicates the bed surface elevation.

[0114] S3. Arrange the data obtained in step S1 and use them as boundary conditions to import into the wave-current-sediment coupled beach evolution prediction model. The boundary conditions required by the model are the land boundary and the seaward boundary. The land boundary adopts a closed boundary condition, and the seaward open boundary adopts different boundary conditions according to different equations.

[0115] The tidal equation adopts the Dirichlet boundary condition of flow (the boundary physical quantity is a constant); the wave equation adopts the Dirichlet boundary condition of root mean square wave height; to prevent excessive sedimentation at the boundary, the sediment equation adopts the Dirichlet boundary condition at high tide and the Neumann boundary condition at low tide (the gradient of the boundary physical quantity is a constant).

[0116] S4. Use the finite difference method to iteratively solve the wave-current-sediment coupled shore evolution prediction model to obtain the dynamic evolution process of the tidal flat landform, which specifically includes the following steps:

[0117] S4.1. During the solution process, based on the provided boundary conditions, the tidal dynamics module and the wave propagation module are solved using the finite difference method to calculate the flow velocity and wave velocity distribution on the tidal flat terrain.

[0118] The tidal dynamic equation includes the time derivative term Where h represents water depth, t represents time, and the time derivative term By directly taking the derivative of the water surface function with respect to time, the equation is transformed into an ordinary differential equation:

[0119]

[0120] Where x represents the cross-shore distance and q represents the flow rate per unit width, in m 2 / s, q = hu, q is discretized using the first-order forward difference, u represents the depth-averaged velocity in the tidal flat section direction; i represents the spatial scale, Δx represents the spatial step size, q i+1 Indicates the single-width flow q value at spatial position i+1, q i Represents the single-width flow q value at spatial position i; when the boundary flow q is determined bd =∫Δh(x, t)dx / Δt, the equation is directly solved by iteration, where Δh represents the change in water level and Δt represents the time step.

[0121] S4.2. Substitute the velocity and wave velocity distribution obtained from the wave dynamics module and the tidal dynamics module into the sediment entrainment force formula to calculate the sediment entrainment force. Discretely solve the average suspended sediment convection-diffusion equation and the sediment entrainment force formula to obtain the suspended sediment concentration distribution on the tidal flat terrain.

[0122] The wave energy conservation equation in the wave power module includes a convection term The first-order forward difference method is used to solve the spatial equation. In terms of time terms, the backward Euler method is used to implicitly solve the equation. It is unconditionally stable at different time steps and improves the model operation efficiency, as shown in the following formula:

[0123]

[0124] Among them, c g represents the wave group velocity, E represents the wave energy, x represents the cross-shore distance, Δx represents the spatial step size, and j represents the time scale. This format requires solving a system of linear equations at each time step. The coefficient matrix is ​​a tridiagonal matrix with a dominant principal diagonal, which can be solved by the Thomas algorithm. It represents the wave energy at spatial position i and time scale j+1; It represents the wave energy at spatial position i+1 and time scale j+1; represents the wave energy at spatial position i and time scale j; represents the group velocity at spatial position i+1 and time scale j+1; represents the group velocity at spatial position i and time scale j+1; It represents the total wave energy at spatial position i and time scale j.

[0125] S4.3. Based on the solved suspended sediment concentration distribution, substitute it into the bottom bed update equation and use the forward Euler method to perform iterative solution to obtain the tidal flat profile topography at different simulation stages.

[0126] Substitute the flow velocity and wave velocity distribution obtained by the wave and tidal dynamics module into the sediment entrainment formula Calculate the sediment entrainment force and solve the suspended sediment concentration distribution based on the embedded average suspended sediment convection-diffusion equation and the sediment entrainment force.

[0127] Among them, the convection term is required to solve the suspended sediment convection diffusion equation. Special processing is performed. Since the tidal current has opposite directions during high and low tides, in order to avoid numerical oscillation, the discretization of the convection term adopts the upwind format, that is, the first-order backward difference is used when the flow velocity u>0 (low tide), and the first-order forward difference is used when the flow velocity u<0 (high tide); the time discretization also adopts the implicit backward Euler method; to ensure strict dominance of the main diagonal, the source term part less than 0 -αωC is linearized implicitly, where α represents the constant 1 and ω represents the sediment settling velocity; the specific discretization format is as follows:

[0128] Time item:

[0129]

[0130] Diffusion term:

[0131]

[0132] Convection term:

[0133]

[0134] in, represents the water depth at spatial position i and time scale j; represents the water depth at spatial position i-1 and time scale j; represents the water depth at spatial position i+1 and time scale j; It represents the suspended sediment concentration at spatial position i and time scale j+1; represents the suspended sediment concentration at spatial position i and time scale j; It represents the suspended sediment concentration at the spatial position i-1 and the time scale j+1; It represents the suspended sediment concentration at the spatial position i+1 and the time scale j+1; represents the flow velocity at spatial position i+1 and time scale j; represents the flow velocity at spatial position i and time scale j; represents the flow velocity at spatial position i-1 and time scale j;

[0135] Substitute the suspended sediment concentration distribution into the bed update equation and use the forward Euler method to iteratively solve the problem. The discrete format is:

[0136]

[0137] in, It represents the sand-carrying force at spatial position i and time scale j; It represents the sand-carrying force at spatial position i and time scale j+1;

[0138] Obtain the simulated tidal flat profile topography and analyze the tidal flat evolution law under the action of waves and currents.

[0139] In this embodiment, the measured section of the tidal flat near the Wanggang River estuary in Jiangsu Province is taken as an example, the cross-section terrain is simplified, and the initial slope of the model is taken as 0.2%. Figure 2 The figure shows a schematic diagram of the measured topography and the initial topography of the model. The elevation range is -6 to 4 m. The tidal ranges of M2 and S2 in this area are 3.4 m and 1.2 m respectively. The simplified tidal range is 4 m. Two working conditions, wave and no wave, are considered. The wave period is set to 3 s, and the effective wave height is set to 0 m and 0.3 m respectively. The sediment part provides parameters such as boundary suspended sediment concentration C, diffusion coefficient ε, porosity p, sedimentation velocity ω, and density.

[0140] The initial slope of the model is 0.2%, and the calculation grid spacing is 1m, so the bed elevation zb Stored in a 5000*1 grid, the water surface function is obtained based on the tidal range of 4m. Water depth h = η - z b , which constitutes the basic condition for solving the tidal power module.

[0141] wave Given the initial wave heights of 0 and 0.5 m, the results are as follows Figures 3 and 4 As shown, wave bottom friction Wave energy dissipation coefficient α f At a value of 0.015, the wave breaking term The wave breaking probability Qb is determined by the discriminant index b, Constitute the wave power module.

[0142] The sediment transport capacity module is a one-dimensional vertical average suspended sediment convection diffusion equation and a bottom bed renewal equation. According to the sediment conditions, the boundary suspended sediment concentration C is set to 0.3 kg / m 3 , the diffusion coefficient ε is 10m 2 / s, porosity p = 0.5, α is 1, sediment settling velocity ω is 0.0005m / s, sediment density and water density are 1650kg / m 3 and 1000kg / m 3 , thus forming a sediment power module.

[0143] In the tidal dynamic equation It can be obtained by directly derivatizing the water surface function with respect to time, where q is the flow rate per unit width (m 2 / s), q=hu, so the equation is transformed into:

[0144]

[0145] Discrete first-order forward difference iteration to solve variable u c .

[0146] The wave energy conservation equation includes the convection term The first-order forward difference is used to solve the space, and the backward Euler method is used to implicitly solve the equation in the time term to obtain the variable u w .

[0147]

[0148] Calculate the flow velocity u using the wave and tidal dynamic model c , wave speed u w The distribution of the sediment parameters provided is substituted into the sediment carrying force formula. Calculate the sand holding force, the result is as follows Figure 5As shown in the figure, the suspended sediment concentration distribution is calculated based on the average suspended sediment convection-diffusion equation and the obtained sediment entrainment force. Here, the space is discretized, and the first-order backward difference format is used at low tide, and the first-order forward difference format is used at high tide. The implicit backward Euler method is used for time discretization.

[0149] Finally, the suspended sediment concentration distribution obtained above is substituted into the bottom bed update equation, and the forward Euler method is used to iteratively solve the tidal flat profile topography at different simulation stages, such as Figure 6 The figure shows a comparison diagram of the results of the model method of the embodiment of the present invention and the traditional tidal sediment model. Compared with the Delft3D model, the annual terrain curves are almost identical, and the error is significantly reduced over time, which has a higher simulation accuracy.

[0150] This embodiment uses the sediment carrying capacity formula as the near-bottom sediment flux. The sediment carrying capacity refers to the sediment content of the water body under certain water flow and wave conditions. When the sediment content in the water body is different from the sediment carrying capacity of the water flow, the sediment returns to consistency (suspended or silted) at a certain rate, so that the sediment content in the water body is adapted to the hydrodynamic conditions.

[0151] This embodiment provides a multi-scale wave-current-sediment coupled shore evolution prediction model method. Based on given offshore tidal flat water area topography data, a coordinate system is established. Then, based on observation stations, tidal, wave, and sediment data in the area are obtained. The data are organized and used as boundary conditions to import into the wave-current-sediment coupled shore evolution prediction model. A tidal flat prediction model is constructed based on the integration of tidal dynamics, wave motion, and sediment dynamics modules, and the required model parameters of each module are set. The finite difference method is used to solve the fluid continuity equation, wave energy conservation equation, one-dimensional vertical average suspended sediment convection and diffusion equation, and bottom bed update equation in the model to obtain tidal flat profile topography at different simulation stages and analyze the shore evolution law of the tidal flat profile under the action of waves and currents.

[0152] This embodiment has a sufficient theoretical basis and can more accurately simulate the evolution process of tidal flat profile morphology under the combined action of waves and currents. Compared with the traditional tidal sediment model, it significantly improves the calculation efficiency and simulation accuracy of the evolution of coastal tidal flat morphology.

[0153] In addition to the above embodiments, the present invention may also have other implementations. Any technical solution formed by equivalent replacement or equivalent transformation falls within the scope of protection required by the present invention.

Claims

1. A multi-scale wave-current-sediment coupled shore evolution prediction model method, characterized by: The following steps are involved: S1. Based on the given offshore tidal flat topographic data, a coordinate system is established, and then the tidal, wave, and sediment conditions of the area are obtained based on the observation station; S2. Construct a wave-current-sediment coupled shore evolution prediction model, which includes a tidal current dynamics module, a wave dynamics module, and a sediment dynamics module. The tidal current dynamics module is constructed based on the fluid continuity equation with zero free surface gradients. The wave dynamics module is generated based on the wave energy conservation equation. Based on the sediment entrainment formula, the near-bottom sediment flux F required for the one-dimensional average suspended sediment convection-diffusion equation is provided. s , combined with the bottom bed renewal equation dominated by the near-bottom sediment flux, a sediment dynamics module is constructed; The tidal dynamics module is based on the fluid continuity equation. It assumes that at any time, the gradient of the free surface is 0, that is, there is no oscillation of the liquid surface. The momentum equation is ignored, and the control equation is as follows: Where u represents the average velocity in the depth of the tidal flat section, t represents time, and x represents the distance across the shore; water depth h = η - z b , η represents the water surface height, z b Indicates bed elevation; The water surface height η function satisfies the following formula: Where A represents the tidal amplitude, T c represents the tidal cycle, T c =12 hours; The wave power module is composed based on the wave energy conservation equation, which is as follows: Where E represents wave energy, ρ represents the density of seawater, g represents the acceleration due to gravity, H rms Represents the root mean square wave height; S tot represents the total wave energy term, c g represents the group velocity, S fric represents the wave bottom friction term, S break represents the wave breaking source term; Group velocity c g The propagation speed of wave energy is expressed as follows based on linear wave theory: Where c = σ / k represents the phase velocity of the wave, σ represents the wave frequency, k represents the wave number, and the direction of the wave phase velocity is fixed to the shoreward direction; Wave bottom friction term S fric As shown in the following formula: Among them, α f represents the wave energy dissipation coefficient, H represents the wave height, T w represents the wave period; When the water depth gradient is large and the wave height exceeds a certain proportion of the water depth, the waves break, as shown in the following formula: Where α represents a constant of 1; Qb represents the probability of wave breaking determined by the discriminant index b. H represents wave height, the maximum wave height is H max =0.78h; Based on the sediment entrainment formula, the near-bottom sediment flux F required for the one-dimensional average suspended sediment convection-diffusion equation is provided. s , and combined with the bottom bed renewal equation dominated by the near-bottom sediment flux to form a sediment dynamics module. The sediment model control equation adopts the one-dimensional vertical average suspended sediment convection-diffusion equation, which is as follows: F s =αω(C * -C) Where h represents water depth, C represents vertical average sediment concentration, t represents time, u represents depth average velocity in the tidal flat section, x represents cross-shore distance, ε represents horizontal suspended sediment diffusion coefficient, and F represents the vertical average sediment concentration. s represents the suspended sediment flux near the bottom, α represents the constant 1, ω represents the sediment settling velocity, and C* represents the sediment entrainment force; The near-bottom sediment flux will cause scouring or silting of the bottom bed. The bed surface renewal equation is expressed by the following formula: Where p represents the porosity of the bed sediment, ρ s represents the sediment density, z b Indicates bed elevation; S3, arranging the data obtained in step S1 and using them as boundary conditions to import into the wave-current-sediment coupled shore evolution prediction model; S4. Use the finite difference method to iteratively solve the wave-current-sediment coupled beach evolution prediction model to obtain the dynamic evolution process of tidal flat landform morphology.

2. The multi-scale wave-current-sediment coupled shore evolution prediction model method according to claim 1, characterized in that: The sediment entrainment capacity C* has the same dimension as the vertical average sediment concentration C: Among them, C D represents the drag coefficient, ρ and ρ s represent the density of water and sediment respectively, U represents the dynamic factor, and g represents the acceleration due to gravity; Drag coefficient C D It is related to the Reynolds number and the bed surface roughness, and is expressed by the Manning formula: The calculation formula of power factor U is: Among them, n represents the roughness, u c represents the tidal flow velocity, u w represents the wave velocity, f c and f w represent the tidal and wave resistance coefficients f respectively c =gn 2 / h 1 / 3 ,, ξ represents the displacement amplitude of the bottom particles, ξ=u w T w / 2π,T w represents the wave period, k s represents the Nicolaz function value; The formula of sand holding force is summarized as follows:

3. The multi-scale wave-current-sediment coupled shore evolution prediction model method according to claim 1, characterized in that: In step S3, the boundary conditions required by the wave-current-sediment coupled beach evolution prediction model are the land boundary and the seaward boundary. The land boundary adopts a closed boundary condition, and the seaward open boundary adopts different boundary conditions according to different equations. The tidal equation adopts the Dirichlet boundary condition of flow, and the boundary physical quantity is a constant; the wave equation adopts the Dirichlet boundary condition of root mean square wave height; the sediment equation adopts the Dirichlet boundary condition at high tide and the Neumann boundary condition at low tide, and the gradient of the boundary physical quantity is a constant.

4. The multi-scale wave-current-sediment coupled shore evolution prediction model method according to claim 1, characterized in that: The step S4 specifically includes the following sub-steps: S4.

1. Based on the provided boundary conditions, solve the tidal dynamics module and the wave dynamics module using the finite difference method to calculate the flow velocity and wave velocity distribution on the tidal flat terrain; S4.

2. Substitute the velocity and wave velocity distribution obtained from the wave dynamics module and the tidal dynamics module into the sediment entrainment force formula to calculate the sediment entrainment force. Discretely solve the average suspended sediment convection-diffusion equation and the sediment entrainment force formula to obtain the suspended sediment concentration distribution on the tidal flat terrain. S4.

3. Based on the solved suspended sediment concentration distribution, substitute it into the bottom bed update equation and use the forward Euler method to perform iterative solution to obtain the tidal flat profile topography at different simulation stages.

5. The multi-scale wave-current-sediment coupled shore evolution prediction model method according to claim 4, characterized in that: In step S4.1, the tidal current dynamic equation includes the time derivative term Where h represents water depth, t represents time, and the time derivative term By directly taking the derivative of the water surface function with respect to time, the equation is transformed into an ordinary differential equation: Where x represents the cross-shore distance and q represents the flow rate per unit width, in m 2 / s, q = hu, q is discretized using the first-order forward difference, u represents the depth-averaged velocity in the tidal flat section direction; i represents the spatial scale, Δx represents the spatial step size, q i+1 Indicates the single-width flow q value at spatial position i+1, q i Represents the single-width flow q value at spatial position i; when the boundary flow q is determined bd =∫Δh(x, t)dx / Δt, the equation is directly solved by iteration, where Δh represents the change in water level and Δt represents the time step.

6. The multi-scale wave-current-sediment coupled shore evolution prediction model method according to claim 4, characterized in that: In step S4.2, the wave energy conservation equation in the wave power module includes the convection term The first-order forward difference method is used to solve the spatial equation, and the backward Euler method is used to implicitly solve the time term: Among them, c g represents the wave group velocity, E represents the wave energy, x represents the cross-shore distance, Δx represents the spatial step size, and j represents the time scale; It represents the wave energy at spatial position i and time scale j+1; It represents the wave energy at spatial position i+1 and time scale j+1; represents the wave energy at spatial position i and time scale j; represents the group velocity at spatial position i+1 and time scale j+1; represents the group velocity at spatial position i and time scale j+1; It represents the total wave energy at spatial position i and time scale j.

7. The multi-scale wave-current-sediment coupled shore evolution prediction model method according to claim 4, characterized in that: In step S4.3, the suspended sediment convection diffusion equation in the sediment dynamics module includes the convection term Where h represents the water depth, u represents the depth-averaged velocity in the tidal flat section, C represents the vertical average sediment concentration, and x represents the cross-shore distance. The discretization of the convection term adopts the upwind format, that is, when u>0 at low tide, first-order backward difference is used, and when u<0 at high tide, first-order forward difference is used. Time discretization also adopts the implicit backward Euler method. The source term less than 0 is implicitly linearized by -αωC, where α represents a constant of 1 and ω represents the sediment settling velocity. The specific discretization format is as follows: Time item: Diffusion term: Convection term: in, represents the water depth at spatial position i and time scale j; represents the water depth at spatial position i-1 and time scale j; represents the water depth at spatial position i+1 and time scale j; It represents the suspended sediment concentration at spatial position i and time scale j+1; represents the suspended sediment concentration at spatial position i and time scale j; It represents the suspended sediment concentration at the spatial position i-1 and the time scale j+1; It represents the suspended sediment concentration at the spatial position i+1 and the time scale j+1; represents the flow velocity at spatial position i+1 and time scale j; represents the flow velocity at spatial position i and time scale j; represents the flow velocity at spatial position i-1 and time scale j; Substitute the suspended sediment concentration distribution into the bed update equation and use the forward Euler method to iteratively solve the problem. The discrete format is: in, It represents the sand-carrying force at spatial position i and time scale j; It represents the sand-carrying force at spatial position i and time scale j+1; Obtain the simulated tidal flat profile topography and analyze the tidal flat evolution law under the action of waves and currents.

Images