Abstract
The silted-up sediment in the reservoir may have a significant influence on the propagation of dam-break flows. In this paper, a three-dimensional numerical simulation of the silted-up dam-break flow is carried out. In this paper, simulations of three-dimensional silted-up dam-break flow are carried out. A kind of Eulerian–Eulerian two-fluid model (TFM), coupled level set and volume of fluid (CLSVOF) methods, is presented. In order to calculate the motions of the air–water interface and the sediment simultaneously, kinetic particle theory (KPT) and computational fluid dynamics (CFD) are combined. The rheology-based constitutive equations of sediment are also considered to simulate scouring and deposition. In addition, a partial-slip boundary condition (BC) for the velocity of the sediment phase at stationary walls is implemented. The simulation results of the benchmark case demonstrate that the proposed model can effectively simulate the silted-up dam-break flow while taking into account multi-interface capturing problems. Subsequently, the simulations of the silted-up dam-break flow over dry are investigated numerically in a three-dimensional long channel. The simulated results reveal that, in the dam-break flows, the silted-up sediment height has a significant influence on wave propagation, dynamic pressure loads, sediment transport, and sediment deposition.
HIGHLIGHTS
Three-dimensional simulation of silted-up dam-break flow is carried out.
The air–water interface movement is captured by CLSVOF method.
A Eulerian–Eulerian multiphase model coupling kinetic particle theory and computational fluid dynamics is used.
The rheology-based constitutive equations of sediment are considered.
The effect of the silted-up sediment height during a dam-break flow is investigated.
INTRODUCTION
Dam-break flood is triggered by the sudden structural damage of dams, which leads to a surge in downstream water level and generates unstable flood waves. In general, the investigation of dam-break flows could be roughly divided into air–water interface, dynamic pressure loads, and sediment transport. For air–water interface, Kocaman & Ozmen-cagatay (2012) explored the effect of lateral channel contraction on dam-break waves. Kocaman & Ozmen-cagatay (2015) conducted a similar dam-break flow experiment to investigate water level fluctuations and free surface evolution. Ye & Zhao (2017) used the two-liquid volume of fluid (VOF) approach to analyze the evolution of water–water interface between reservoir water and downstream water in wet-bed dam-break flows. By using records of water surface profiles, stage hydrographs, and cross-sectional mean velocities, Wang et al. (2020) investigated the similarity of dam-break flows on dry beds. In addition, three-dimensional simulations of dam-break flows striking rigid structures have also been presented by Arote et al. (2021) and Alloush et al. (2021). These studies enable a thorough understanding of flood front propagation and water stage changes of clear-water dam-break flows.
For dynamic pressure loads in dam-break flows, Lobovský et al. (2014) conducted a series of experiments to study the pressure loads on the downstream vertical wall during dam breaks. Furthermore, Mohamed et al. (2018) performed an experimental investigation of the dam-break flow impacting on a vertical cylinder positioned over a horizontal dry bed. Impacting on downstream walls or inflexible obstacles would usually lead to violent free surface deformations and surging dynamic pressure loads in dam-break flows (Meng et al. 2021a; He et al. 2022; Maghsoodi et al. 2022). The impacting pressure of dam-break floods on a solitary structure or various types of obstacles in the downstream region has also been numerically studied (Meng et al. 2021a; He et al. 2022; Maghsoodi et al. 2022), which provides accurate fluctuation profiles of air–water interface and time evolution of dynamic pressure on the rigid-body surface.
However, the propagation of air–water–sediment multiphase dam-break flows may be significantly different from that of air–water two-phase flow. The influence of sediment on dam-break wave propagation may be of great importance. In recent years, studies on the scouring process and the meandering effect have also shown that sediment plays an important role in dam-break flows (Bihs & Olsen 2011; Pu & Lim 2014; Zheng et al. 2016). In addition, turbulence has been shown to be another key effect to be investigated (Pu et al. 2012a, 2012b, 2014; Pu 2015), especially in dam-break flows. Moreover, sediment transport in dam-break flows may result in turbulence and affect the dam-break flow process. Previous studies mainly focused on clear-water dam-break flows. There is still a lack of quantitative investigations of the effect of silted-up sediment height on dam-break-induced flood. However, there has been no numerical investigation on silted-up dam-break flows because it remains a challenge for existing solvers to simulate rheology-based sediment movement with complicated fluctuations in air–water interface. Therefore, to fill the gap and improve numerical simulation methods, this study aims to investigate the effects of silted-up sediment height on a dam-break flow by proposing an Eulerian–Eulerian multiphase model. Different from Lagrangian methods such as SPH (Chang & Chang 2013; Pu et al. 2013; Shi et al. 2019), this model adopts the Eulerian approach. Compared with the one-fluid Eulerian model such as SWEs (Pu et al. 2012a, 2012b), this study employs the two-fluid Eulerian multiphase method to simulate both air–water interface and sediment movement.
The innovations of this study can be summarized as follows:
- (1)
The proposed Eulerian–Eulerian two-fluid model (TFM), coupled with the level set and volume of fluid (CLSVOF) method, is a combination of the kinetic particle theory (KPT) and computational fluid dynamics (CFD) to calculate air–water interface evolution and sediment movement.
- (2)
In the proposed multiphase model, the rheology-based constitutive equations and a partial-slip boundary condition (BC) of sediment are implemented to simulate scouring and deposition.
- (3)
Unlike previous clear-water dam-break flows that only considered the water interface, this study also considers the effect of silted-up sediment, providing a new insight into the emerging silted-up dam-break flow problem.
EULERIAN–EULERIAN CLSVOF METHOD
Interface capturing


















Momentum equations

























Rheology-based constitutive equations of sediment
















Sediment–fluid interaction model
In silted-up dam-break flows, the carrier phases (air and water) are responsible for promoting the motion of the solid phase (sediment), and the sediment will in turn affect the fluid movement as well. Since drag force plays a dominant role in dam-break flows with high sediment concentration (Shi et al. 2019; Meng et al. 2022), only the drag force () is considered in this study.












Partial-slip boundary condition for the sediment phase











Numerical scheme and procedure
For the simulations performed in this study, the Euler implicit scheme is selected for time discretization. For spatial terms: Gauss vanLeer for ; Gauss linear for other divergence terms; Gauss linear for gradient terms and Laplacian terms. The interpolation scheme is set to linear interpolation. The PCG solver and the smooth solver are chosen to solve the equations of the pressure field and the velocity field, respectively. The tolerance of each iteration of every variable is set to be less than 1.0 × 10−6. Except for the sediment phase velocity, other velocity boundary conditions are set as no-slip. The pressure BC of walls is set as a fixed flux pressure boundary. The boundary conditions of walls for other scalars (
,
and
) are all set to zero gradient.
The governing equations are embedded into the PIMPLE scheme to form a new multiphase flow solver. The new overall solution algorithm for the proposed multiphase flow model is summarized as follows:
(Step 1) Determine and input the constant parameters in the simulation, such as gravity vector, sediment particle size, sediment density, tension coefficient, etc.
(Step 2) Initialize the fields and set the boundary conditions.
(Step 3) Use the CLSVOF solver for the new
and
fields by Equations (1) and (2).
(Step 4) Update the physical properties of the fluid (including air and water) phase in the whole computational domain by Equation (3) and calculate the surface tension term.
(Step 5) Compute and update sediment intergranular pressure
, sediment kinematic viscosity
, drag force
, viscous stress
, and
(k=f, s).
(Step 6) Solve Equation (5) for the intermediate predicted fluid (including air and water) velocity.
(Step 7) Re-solve Equation (5) for the fluid pressure and then correct the final fluid velocity based on the new pressure field solved in this step.
(Step 8) Solve Equations (6) and (22) for the sediment phase velocity and correct the velocity at wall boundaries.
(Step 9) Solve Equation (4) to obtain the new sediment concentration distribution.
(Step 10) If
< 0, set
= 0; if
>
, set
= 0.995
.
(Step 11) Repeat Step 3 to Step 10 until the end of the simulation.
MODEL VERIFICATION
The simulation of a three-dimensional silted-up dam-break flow is employed to demonstrate the solver's capability. In the experiment carried out by Duarte et al. (2011), the initial sediment column (1.5 m in length × 0.42 m in width × hs0 in height) was set in a reservoir (1.5 m in length × 0.42 m in width × hw0 in height). In their study, eight sets of experiments were conducted in a rectangular tank with the size of 5.5 m in length × 0.42 m in width × 0.42 m in height. Two sets of their experiments are taken in this study for simulation: (1) hw0 = 0.41 m, hs0 = 0.43 m, = 1,565 kg/m3, and (2) hw0 = 0.40 m, hs0 = 0.22 m,
= 1,671 kg/m3. Since the initial volume fraction of sediment was not specified in their experiments (Duarte et al. 2011), we take the initial volume fraction as
= 0.533 according to a similar experiment by Wang et al. (2017). The density and viscosity of water are 1,000 kg/m3 and 1.0 × 10−6 m2/s, respectively. The density and viscosity of the air phase are 1 kg/m3 and 1.0 × 10−5 m2/s, respectively. The model parameters and material constant values are summarized in Table 1. Among them, K = 1.0 × 105 Pa is used because an extremely high Young's modulus causes numerical instability and waste of computational cost (Lee 2021; Meng et al. 2022). The maximum Courant–Friedrichs–Lewy (CFL) number in high-concentration regions should be less than 0.005, because a too high CFL number will lead to non-convergence results, especially in sediment siltation areas. Therefore, the time step ranges from 1.0 × 10−6 to 1.0 × 10−4 s in this study.
Model parameters adopted in this study
Parameter . | ![]() | ![]() | ![]() | b . | a . | ![]() | ![]() | ![]() | ![]() |
---|---|---|---|---|---|---|---|---|---|
Value | 0.1 | 5 | 0.6 | 2 | 0.1 | 0.585 | 1.5 | 0.3 | 0.8 |
Parameter | ![]() | ![]() | ![]() | ![]() | ![]() | ![]() | K (Pa) | ![]() | ![]() |
Value | 0.48 | 0.62 | 0.85 | 0.01 | 1.000 | 3.6 | 1.0 × 105 | 0.33 | 288 |
Parameter . | ![]() | ![]() | ![]() | b . | a . | ![]() | ![]() | ![]() | ![]() |
---|---|---|---|---|---|---|---|---|---|
Value | 0.1 | 5 | 0.6 | 2 | 0.1 | 0.585 | 1.5 | 0.3 | 0.8 |
Parameter | ![]() | ![]() | ![]() | ![]() | ![]() | ![]() | K (Pa) | ![]() | ![]() |
Value | 0.48 | 0.62 | 0.85 | 0.01 | 1.000 | 3.6 | 1.0 × 105 | 0.33 | 288 |
Comparisons between simulated and experimental results (Duarte et al. 2011): (a) air–water–sediment interface and (b) displacements of the water wave front (XW) and sediment front (XS).
Comparisons between simulated and experimental results (Duarte et al. 2011): (a) air–water–sediment interface and (b) displacements of the water wave front (XW) and sediment front (XS).
NUMERICAL RESULTS AND DISCUSSION
Geometry and conditions
Simulation conditions in this study
No. . | hs (m) . | hs/H (%) . |
---|---|---|
B1 | 0 | 0 |
B2 | 0.02 | 13.33 |
B3 | 0.03 | 20.0 |
B4 | 0.045 | 30.0 |
B5 | 0.06 | 40.0 |
No. . | hs (m) . | hs/H (%) . |
---|---|---|
B1 | 0 | 0 |
B2 | 0.02 | 13.33 |
B3 | 0.03 | 20.0 |
B4 | 0.045 | 30.0 |
B5 | 0.06 | 40.0 |
To better investigate water level fluctuations, four virtual depth probes are positioned along the center line of the channel (the x direction). Positions of these virtual probes are shown in Figure 2 and Table 3. In addition, four pressure probes are set to analyze the dynamics of the silted-up dam-break problem. Figure 2 (right) and Table 3 show the probe layout. To monitor the dynamic pressure during sediment collapse, the first centerline probe (P1) is positioned directly at the bottom of the channel. The other three pressure probes of the centerline are arranged on the downstream wall and the corresponding locations are 10 mm (P2), 20 mm (P3), and 40 mm (P4) above the channel bed.
Positions of water level and dynamic pressure probes
Probe . | Position (m) . | Probe . | Position (m) . |
---|---|---|---|
H1 | (0.1,0.05,0) | P1 | (0.1,0.05,0) |
H2 | (0.2,0.05,0) | P2 | (1,0.05,0.01) |
H3 | (0.4,0.05,0) | P3 | (1,0.05,0.02) |
H4 | (0.6,0.05,0) | P4 | (1,0.05,00.4) |
Probe . | Position (m) . | Probe . | Position (m) . |
---|---|---|---|
H1 | (0.1,0.05,0) | P1 | (0.1,0.05,0) |
H2 | (0.2,0.05,0) | P2 | (1,0.05,0.01) |
H3 | (0.4,0.05,0) | P3 | (1,0.05,0.02) |
H4 | (0.6,0.05,0) | P4 | (1,0.05,00.4) |
Simulation results
Temporal evolution of silted-up dam-break flows in B1–B5 cases. t = 0.2 s (a), 0.9 s (b), 1.6 s (c), and 2.4 s (d).
Temporal evolution of silted-up dam-break flows in B1–B5 cases. t = 0.2 s (a), 0.9 s (b), 1.6 s (c), and 2.4 s (d).
The present model has successfully simulated typical dam-break structures, such as sediment being deposited upstream, suspended sediment moving downstream, and sediment being exposed to air. The sediment front interacts with the surrounding water as the sediment column slides, making the water head more complicated than the initial stage. The silted-up sediment column and water flow are sheared by gravity in the early stages, which explains why the collapse occurs at the same time. The turbulence in the water flow is strong enough to suspend the sediment particles in the middle and later stages as the flow develops, causing a small amount of separated sediment particles to go downstream. However, the accumulated sediment does not move significantly, which is quite different from the suspended sediment.
Effect of silted-up sediment height
Initial collapse
Hydrodynamics
Similarly, the measured peaks at P2–P4 suggest that the dam-break wave has reached and impacted the downstream wall. The main difference among different cases appears at the pressure peak, including the peak response time and the peak value. In clear-water dam-break flow (B1), the peak pressure value is always the highest. In silted-up cases, as the silted-up sediment height increases, the peak value of the pulsating pressure also increases.
Wave propagation and sediment transport
Propagation of the positive water wave (a) and sediment interface (b) in all cases.
Propagation of the positive water wave (a) and sediment interface (b) in all cases.
In addition, as shown in Figure 7(b), in all the cases, the propagation of sediment front over time is also traced while the reverse wave is neglected. The trajectories are relatively steep before it reaches the downstream wall at about 0.6 s, and then gradually become gentler beyond that time. As the silted-up sediment height rises, the slope of the sediment propulsion curve becomes milder. Furthermore, the time for the sediment to reach the downstream wall also differs. On the whole, the sediment propagation curves of B3 and B4 are relatively similar, but they are observably different from those of B2 and B5. It is indicated that the silted-up sediment height has a significant influence on its propagation. Theoretically, different silted-up sediment heights determine whether the dam-break process is dominated by water–sediment drag force or by stress among sediment particles. This is the decisive factor that affects the speed of dam-break flow and wave propagation. Only the propagation of the positive wave is considered in this study. The reflection and attenuation of positive and reverse waves is an interesting but quite challenging topic, which deserves our attention in the future.
Figure 8 also shows the final deposition state (the results between about 1.4 and 2.4 s) of upstream silted-up sediment for cases B2–B5 if the downstream wall boundary is ignored. This is because, after about 1.4 s, the longitudinal features of the sediment deposited upstream no longer change until the reverse wave reaches the upstream at about 2.5 s. When the reverse wave reaches the upstream, the reverse wave washes the deposited sediment, thereby changing its deposition feature. It can be concluded that in silted-up dam-break flows, the silted-up sediment height will have a significant impact on the sediment collapse process and final deposition feature.
CONCLUSION
Based on the Eulerian–Eulerian multiphase model, a three-dimensional numerical study has been conducted to investigate silted-up dam-break flows with various silted-up sediment heights and ambient water depths. The improved CLSVOF approach and the proposed two-liquid method are used to calculate both the air–water interface evolution and sediment transport. In addition, by coupling the dynamic mechanism of sediment particles with CFD, rheology-based constitutive equations are used to simulate the deposition and erosion of sediment. The results of benchmark cases of submarine sediment collapse demonstrate that the proposed model can effectively simulate evolutions of both air–water and sediment interfaces, suggesting that the solver can reproduce gravity-dominated motion in multi-interface capturing problems. Then, the effects of silted-up sediment height on wave propagation, dynamic pressure loads, sediment transport, and sediment deposition are numerically investigated by the proposed multiphase flow model. The focus of the analysis is on the effects of silted-up sediment height on wave propagation, dynamic pressure loads, sediment transport, and sediment deposition. As the height of silted-up sediment height grows, more gravitational potential energy is consumed by intergranular stresses of the sediment phase, thereby taking longer for the dam-break flow to reach the downstream wall. However, the reverse wave reaches the upstream wall earlier. In terms of dynamic pressure loads on the downstream wall, the peak pressure values in clear-water dam-break flows are always larger than those in silted-up cases. However, in the silted-up cases, the peak value of the pulsating pressure rises with increasing silted-up sediment height. For the silted-up sediment upstream, as the silted-up sediment height declines, the time when the leftmost sediment begins to drop that indicates the beginning of a full collapse will be postponed. It demonstrates that raising the silted-up sediment height accelerates the collapse. The silted-up sediment height also has an influence on the sediment collapse process and final deposition appearance. This study provides a new insight into the emerging silted-up dam-break flow problem by proposing a new numerical model. However, the model still has some shortcomings, such as high computational time cost and less stability than a two-phase flow solver. The parameters of the model are more complex, which need to be optimized. In the future, we can study the differences in using various models (such as SPH models) for simulating silted-up dam-break flows. Turbulence might also be considered in three-dimensional silted-up dam-break flows in the future, which is a challenging but interesting problem.
ACKNOWLEDGEMENTS
This work was supported by the Natural Science Foundation of Zhejiang Province, China (No. LY20D010009, No. LHY22E080004), Open Fund of State Key Lab of Hydraulics and Mountain River Engineering, China (No. SKHL2015), Open Fund of Key Laboratory of Flood & Drought Disaster Defense, the Ministry of Water Resources, China (No. 2021-01), National Natural Science Foundation of China (No. 51779162), National Natural Science Foundation of China (No. 51979178), National Natural Science Foundation of China (No. 52122904).
DATA AVAILABILITY STATEMENT
All relevant data are included in the paper or its Supplementary Information.
CONFLICT OF INTEREST
The authors declare there is no conflict.