## Abstract

This paper investigates the cascading dam-break flood propagation on the downstream sloping channel and reservoir using the shallow water equations (SWEs) and the Reynolds-average Navier-Stokes equations (RANS). The calculated surface profiles, stage hydrographs and maximum run-up heights for 24 sets of initial conditions are elaborately compared with the experimental measurements, which show the SWEs reproduce the wave oscillation evolution and the maximum run-up height inaccurately. The maximum run-up heights calculated by the SWEs are all smaller than those by the RANS and the measured results, the maximum errors are within −10% to −25%, which may predict delay of the downstream dam-break. However, the maximum errors calculated by the RANS are within ±10%. So the RANS predict more accurate results than the SWEs. Additionally, the generation of short waves must be below a certain upstream-to-downstream ‘depth ratio’, roughly the ‘depth ratio’ <2/3 in this study. If the ratio is too high, it is difficult to form a wavy push due to air entrainment and turbulence. The SWEs predict more accurate results for shallow initial depths than deep initial depths. Therefore, the advantage of the RANS can be more obvious for deep initial depths.

## INTRODUCTION

In the past decades, the dam-break problem has been a very important topic for hydraulic community due to its catastrophic destruction to the downstream, so numerous studies (Whitham 1955; Katopodes & Strelkoff 1978; Fraccarollo & Toro 1995) have focused on the dam-failure problem and have proved that the shallow water equations (SWEs) are a suitable model for large-scale flood that can well balance calculation accuracy and efficiency. However, most of the previous researches have paid attention to dam-break floods in a single reservoir (Fraccarollo & Toro 1995; Shigematsu *et al.* 2004; Ozmen-Cagatay & Kocaman 2010; Cantero-Chinchilla *et al.* 2016; Castro-Orgaz & Chanson 2017), and very little attention has been focused on dam-break floods in cascade reservoirs. At present, due to the needs of social development and technological advancement, cascade reservoirs have been constructed in many river basins in China. Once the upstream dam fails and the flood falls into the downstream reservoir, the flood wave can have a strong impact on the downstream dam and a very high wave run-up on the downstream dam, which can cause the downstream reservoir to crash. This process is repeated downstream to form a cascading dam-break, which can lead to much more serious consequences than a single dam-break. For example, in 1975, ‘75.8’ torrential rain caused more than 62 cascading dam-breaks due to overtopping, including Banqiao, Shimantan, Tiangang and Zhugou *et al.* in Henan province, China, leading to more than 20,000 deaths, a considerable amount of property damage and ecological damage (Chen 2012). Therefore, it is of great significance to study the cascading dam-break. At present, only a small number of experimental results (Xue *et al.* 2011; Xu *et al.* 2013; Chen *et al.* 2014) and numerical simulation results for cascading dam-break (Sun *et al.* 2014; Zhou *et al.* 2014) are available. These numerical results are mainly concerned with macroscopic flood-risk analysis and have paid little attention to the detailed evolution process of the upstream dam-break flood. However, the most basic and core problem of cascading dam-break remains the single upstream dam-failure, but because connected intact reservoirs are downstream of the failure, the upstream dam-break wave overtops the downstream dams. Only after the mechanism of the single dam-break and the flood evolution process in the downstream reservoir are clearly understood, can the subsequent and much more serious cascading dam-break process be truly understood. Hence, this study focuses on the evolution process of the flood wave in the downstream reservoir in the case of the instantaneous failure of the upstream dam.

Many researchers have studied the single dam-break problem based on the SWEs (Fraccarollo & Toro 1995; Cao *et al.* 2011, 2014; Tabrizi & Hessaroeyeh 2015). At the same time, it is found that the SWEs cannot accurately reproduce the negative rarefaction wave for dry bed and the wave-like oscillation evolution process in the initial stage of the dam-break wave for wet bed due to the hydrostatic pressure assumption (Shigematsu *et al.* 2004; Ozmen-Cagatay & Kocaman 2010; Wang *et al.* 2016; Castro-Orgaz & Chanson 2017). Therefore, some scholars have studied the dam-break wave with a Boussinesq-type model (Cantero-Chinchilla *et al.* 2016; Wang *et al.* 2016; Castro-Orgaz & Chanson 2017) or full three-dimensional Reynolds-averaged Navier-Stokes equations (RANS) coupled with various turbulent closures (Shigematsu *et al.* 2004; Ozmen-Cagatay & Kocaman 2010), which are not limited by hydrostatic pressure assumption used in the SWEs. Additionally, some scholars have found solving the dam-break problem on the steeply sloping channel using the SWEs (Denlinger & O'Connell 2008; Castro-Orgaz *et al.* 2015) deviates from the measured results, due to the influence of the convection acceleration term leading to nonhydrostatic pressure on the steep slope. Based on the above reasons, in view of the fact that the cascade reservoirs are built on relatively steep rivers in the mountainous areas, even in steep valleys prone to mudslides, in order to better understand the evolution process of the cascading dam-break flood in the sloping channel, both the RANS and the SWEs models are used for comparative study.

The objective of this study is to clearly and accurately understand the evolution process of the dam-break wave in the downstream reservoir and lay the foundation for further understanding of the cascading dam-failure. The numerical simulation results of water surface profiles, stage hydrographs and run-up heights on the downstream dam for different initial depths were detailed and compared with the published experimental results by Xu *et al.* (2013).

## MATHEMATICAL MODELS

### SWEs

*t*is time,

*x*and

*y*are the spatial directions,

*h*is the water depth,

*U*and

*V*are the depth-averaged velocity in the

*x*and

*y*directions,

*g*is the gravitational acceleration, is the bottom elevation and

*n*is the Manning coefficient.

### RANS

*x*,

*y*and

*z*directions, can be written as where

*t*is time,

*x*,

*y*and

*z*represent directions, , and are fractional areas open to flow in the

*x*,

*y*and

*z*directions, is the fractional volume of the flow in the calculation cell,

*u*,

*v*and

*w*are the velocity in the

*x*,

*y*and

*z*directions, is the fluid density,

*p*is the pressure, , and are the body accelerations in the

*x*,

*y*and

*z*directions, , and are the face acceleration in the

*x*,

*y*and

*z*directions due to the viscosity stress and Reynolds stress.

*F*= 0 stands for no water in the calculated cell,

*F*= 1 stands for water filled in the calculated cell, and

*F*= 0–1 stands for water partially filled in the calculated cell. The free surface can be calculated by the following equation

## VALIDATION TEST AND MODEL SETUP

Recently, Xu *et al.* (2013) implemented a series of cascading dam-break experiments, and those published experimental results are used as benchmarks in this study. The experiment includes a 4-degree positive slope channel connected by the upstream and downstream reservoirs. The upstream reservoir is 2 m long, 0.4 m wide and 0.95 m high. The downstream channel is 10 m long, 0.4 m wide and 0.7 m high; the detailed geometry is presented in Figure 1. A total of four probe points are arranged along the sloping channel spaced 0.6 m from the downstream dam, four wave-height meters are used to record the stage hydrographs at the four probe points after the upstream dam-break. The upstream initial water depths are = 0.10 m, 0.15 m, 0.20 m and 0.25 m, and the downstream initial water depths are = 0 m, 0.10 m, 0.15 m, 0.20 m, 0.25 m and 0.30 m.

In this study, four different sizes of uniform cells, 0.1 m, 0.06 m, 0.02 m and 0.01 m, were examined respectively to conduct sensitivity analysis, with the global refinement ratio of 10 (=0.1 m/0.01 m) above of the recommended minimum value of 1.3 by Celik *et al.* (2008). Comparisons show the 0.02 m and 0.01 m cell sizes almost calculate the same free surface level and wave arrival time, but the calculation results of the above two cell sizes are very different from those of the other two cell sizes. Therefore, the calculation accuracy and efficiency of 0.02 m cell size are relatively high and satisfactory. Finally, an orthogonal grid in which cell size = 0.02 m is selected. The RANS coupled with RNG k-*ɛ* turbulence model solver in Flow-3D are used to solve the 3-D dam-break field. The detailed configurations of geometry, mesh, physical model, boundary conditions, initial conditions and numerical methods are summarized in Table 1. The SWEs are solved by the module of MIKE 21 Flow Model HD FM; its geometry, mesh, physical model, initial conditions, numerical methods, wetting and drying and bottom friction are also detailed in Table 1.

Models . | Items . | Configurations . | Notes . |
---|---|---|---|

RANS | Geometry | 11.9756 × 0.4 × 1.5 | Length × Width × Height |

Cell size | 0.02 × 0.02 × 0.02 | Length × Width × Height | |

Mesh | 599 × 20 × 75 | 898,500 cells hexahedral | |

Physics model | RNG k-ɛ + air entrainment | Default coefficients | |

Top boundary | Symmetry | – | |

Bottom, upstream, downstream, left and right side boundary | wall | = 0.0001 m | |

Initial conditions | = 0.10 m, 0.15 m, 0.20 m, 0.25 m = 0 m, 0.10 m, 0.15 m, 0.20 m, 0.25 m, 0.30 m | 24 runs | |

Numerics | Implicit GMRES pressure solver Explicit viscous stress/advection/free surface pressure solver Momentum advection first order | – | |

SWEs | Geometry | 11.9756 × 0.4 | Length × Width |

Cell size | 0.02 × 0.02 | Length × Width | |

Mesh | 599 × 20 | 11,980 cells square hexahedral | |

Physics model | SWEs | – | |

Initial conditions | = 0.10 m, 0.15 m, 0.20 m, 0.25 m = 0 m, 0.10 m, 0.15 m, 0.20 m, 0.25 m, 0.30 m | 24 sets | |

Solution technique | Time integration/high order Space discretization/high order Minimum time step = 0.001 s Maximum time step = 0.02 s CFL = 0.8 | – | |

Flood and dry | Drying depth = 5e - 005 m Flooding depth = 0.0005 m Wetting depth = 0.001 m | – | |

Bed resistance | Manning number = | 117.9 |

Models . | Items . | Configurations . | Notes . |
---|---|---|---|

RANS | Geometry | 11.9756 × 0.4 × 1.5 | Length × Width × Height |

Cell size | 0.02 × 0.02 × 0.02 | Length × Width × Height | |

Mesh | 599 × 20 × 75 | 898,500 cells hexahedral | |

Physics model | RNG k-ɛ + air entrainment | Default coefficients | |

Top boundary | Symmetry | – | |

Bottom, upstream, downstream, left and right side boundary | wall | = 0.0001 m | |

Initial conditions | = 0.10 m, 0.15 m, 0.20 m, 0.25 m = 0 m, 0.10 m, 0.15 m, 0.20 m, 0.25 m, 0.30 m | 24 runs | |

Numerics | Implicit GMRES pressure solver Explicit viscous stress/advection/free surface pressure solver Momentum advection first order | – | |

SWEs | Geometry | 11.9756 × 0.4 | Length × Width |

Cell size | 0.02 × 0.02 | Length × Width | |

Mesh | 599 × 20 | 11,980 cells square hexahedral | |

Physics model | SWEs | – | |

Initial conditions | = 0.10 m, 0.15 m, 0.20 m, 0.25 m = 0 m, 0.10 m, 0.15 m, 0.20 m, 0.25 m, 0.30 m | 24 sets | |

Solution technique | Time integration/high order Space discretization/high order Minimum time step = 0.001 s Maximum time step = 0.02 s CFL = 0.8 | – | |

Flood and dry | Drying depth = 5e - 005 m Flooding depth = 0.0005 m Wetting depth = 0.001 m | – | |

Bed resistance | Manning number = | 117.9 |

## RESULTS AND DISCUSSION

### Stage hydrographs at Gauge ‘A’, ‘B’, ‘C’ and ‘D’ for h _{1} = 0.2 m and h _{2} = 0.1 m

Figure 2 shows that, in general, the results calculated by the RANS are more accurate than those calculated by the SWEs in terms of both wave height and depth time process. As far as the maximum wave heights (= the maximum run-up height – downstream initial depth) are concerned, only the error of point ‘A’ calculated by the SWEs is large, and the errors of the remaining calculated results using the SWEs and RANS are small. The calculated water depths of both methods finally tending to steady state agree well with the experimental results. The arrival time of the flood has different degrees of deviation at the ‘C’ and ‘D’ points where the water depths are shallow. The main reason may be that the upstream flood wave interacts strongly with the downstream reservoir here; both air entrainment similar to hydraulic jump and strong turbulent mixing occurs. It can be seen that the calculation result of point ‘C’ by the RANS is more accurate than that of the SWEs, because the RANS model is coupled with the air entrainment model and RNG k-*ε* turbulent closure. There are two stage-like water-depth time processes at the three points ‘B’, ‘C’ and ‘D’. At point ‘C’, a local small peak appears at 4 s, and then a second stage increases steeply at 6 s. The main reason is that the flood pushes the downstream stationary water body forward, leading to the first small peak/stage. The second large stage is mainly caused by the downstream water body that hits the downstream dam and then is reflected back to the upstream. The water depth of point ‘A’ rises at 5 s, while the second water-depth stages of the three points ‘B’, ‘C’ and ‘D’ increase at 5.2 s, 6.2 s and 7.7 s, respectively. It can be seen that the second water-depth stage rise time of points ‘B’, ‘C’ and ‘D’ are delayed later in turn; that is, the reflected wave reaches the three points ‘B’, ‘C’ and ‘D’ in turn. This process can be clearly seen from the evolution of the surface profile in Figure 3. Since the stage hydrographs of the four points ‘A’, ‘B’, ‘C’ and ‘D’ for other upstream and downstream initial water depths are similar to present results, they are omitted to avoid redundancy. In this case, the 3D RANS in Flow-3D takes 638 s and the 2D SWEs in MIKE takes 14 s. So the SWEs greatly save computation time, but the accuracy is not as good as the RANS, such as the maximum run-up height of Gauge ‘A’ shown in Figure 2(a). On the contrary, the calculation accuracy of the RANS is higher than the SWEs, but it is extremely time consuming. Therefore, it is best to use the RANS when local details are considered.

### Water surface profile for h _{1} = 0.15 m and h _{2} = 0.30 m

The dam-break wave reaches the downstream reservoir at 3.6 s and a small hump forms, shown in Figure 3(a). The upstream flood continues to fall into the downstream reservoir, the small hump becomes higher and longer at 4.8 s, shown in Figure 3(b). And then the flood wave reaches the downstream dam at 6 s, shown in Figure 3(c). The calculated flood wave by the RANS runs to the maximum height on the downstream dam at 6.4 s, shown in Figure 3(d). It can be seen that the calculated results by the SWEs changes severely and steeply with discontinuity, while the results calculated by the RANS changes moderately and progressively. The calculated maximum run-up/wave height by the RANS is obviously higher than that calculated by the SWEs. Figure 3(e) shows that the dam-break wave is reflected back to the upstream with water surface oscillation at 8.4 s. However, the wavy reflection evolution processes have not been correctly calculated by the SWEs because of the limitation of the shallow hydrostatic pressure assumption. Due to the large reflection flow rate from the downstream dam, the discontinuous reflection shock wave front forms at the tail of the reservoir. The reflected wave continues to advance upstream shown in Figure 3(f), in which the RANS still show weak water surface oscillation. The reflected wave runs to the upstream sloping channel, and the maximum run-up surface elevation has exceeded the downstream surface elevation. It can be foreseen that the water body can continue to return to the downstream reservoir after it runs to the maximum height on the sloping channel. After repeating the above flood-evolution process several times, the downstream reservoir can tend to be stationary. In general, it can be seen that the surface profile between the two methods is slightly different at the negative rarefaction wave to the upstream reservoir and significantly different in the downstream reservoir, especially the run-up height on the downstream dam at 6.4 s. The maximum run-up height difference shown in Figure 3(d) is about 0.05 m, which can be also seen in Figure 5(f). However, the difference on the slope is very weak, indicating that it is sufficient to simulate the flow on a 4-degree mild slope using the SWEs but insufficient in the downstream reservoir.

Stage hydrographs at Gauge ‘A’ for h _{1}= 0.10 m and h _{1} = 0.15 m

h _{1} = 0.10 m

If the flood wave overtops the downstream dam, it can most likely lead to the downstream dam-break caused by serious erosion of the downstream earth-rock dam. Therefore, the run-up height of the flood wave on the downstream dam is the main indicator of the cascading dam-break problem. This section mainly analyzes the stage hydrographs at Gauge ‘A’, because it is closely related to the run-up height on the downstream dam under different initial conditions. The stage hydrographs at Gauge ‘A’ of the downstream dry reservoir calculated by the SWEs is the closest to the measured value and the RANS calculation result among the six runs for = 0.10 m, shown in Figure 4. The stage hydrographs calculated by the SWEs differ significantly from the measured values for the other five runs in terms of both the maximum wave heights and the water-depth time processes. The maximum absolute error of maximum wave heights between the results of SWEs and measured ones for all runs is nearly 0.05 m; the water-depth time process of RANS and measurement are wave-like oscillation processes and that of SWEs is a straight and discontinued change process. The maximum wave heights calculated by the SWEs are smaller than those calculated by the RANS and the measured values, but both calculated water levels tending to be stable are basically the same. The RANS predict the wave oscillation evolution process and the maximum wave height better than the SWEs for all the runs. The time that the wave reaches the maximum height calculated by the SWEs is always ahead of that by the RANS; however, the predicted arrival time by the RANS coincides reasonably well with experimental results. The higher the downstream initial depth is, the later the wave reaches point ‘A’ due to the stronger resistance to wave propagation from the downstream water in the reservoir.

h _{1} = 0.15 m

For the downstream initial dry bed, the deviation between the maximum wave height predicted by the SWEs and the measured value for upstream initial depth = 0.15 m (Figure 5(a)) is larger than that for upstream initial depth = 0.10 m (Figure 4(a)), which may be due to the greater impact of the upstream flood wave on the downstream. The higher downstream water depth leads to the more serious the water surface oscillation of the downstream reservoir, shown in Figure 5(e) and 5(f) and Figure 4(d)–4(f). The two methods basically predict the same arrival time of flood waves, shown in Figure 5. The RANS results are once again expected to agree with the measurements better than the SWE results. The calculation results of = 0.20 m and 0.25 m are basically consistent with the results of = 0.10 m and 0.15 m, and are not described here for brevity.

### Comparison of stage hydrographs at Gauge ‘A’ for different initial depths

#### The same upstream depth with different downstream depths

For the same upstream initial depth, the greater downstream initial depth results in the slightly smaller maximum wave height and the later wave-arrival time, shown in Figure 6(a)–6(d). Therefore, the downstream water body has a weak water-blocking effect in both the streamwise direction and the vertical direction. For the upstream initial depth = 0.10 m, the corresponding wave peak number is 3, 4 and 3, respectively, for different ratios = 0.1/0.2, 0.1/0.25 and 0.1/0.3 in Figure 6(h)–6(j), indicating that the wave oscillation frequency increases firstly and then decreases as the increases. The possible reason is, for the same upstream initial depth, the downstream initial depth must be some ‘critical depth’, the energy dissipation of the dam-break wave is minimal. For the case of upstream initial depth = 0.10 m, the downstream initial ‘critical depth’ may be = 0.25 m. When the downstream initial depth is greater than the ‘critical depth’, the energy dissipation increases due to the water-blocking effect, or the energy of the approaching flow is smaller relative to the downstream still water, leading to smaller wave oscillation frequency and duration. When the downstream initial depth is smaller than the ‘critical depth’, the dam-break wave directly forms a flow similar to a hydraulic jump or a breaking wave, and the energy dissipation is still larger than that for the ‘critical depth’ condition, so the wave oscillation frequency is still small. In fact, it can be seen from Figure 6(c) and 6(g) that as increases to more than 2/3, the wave's oscillations are more and more difficult to form, which indirectly indicates that the greater upstream initial depth causes the greater impact with the side wall and channel bed. The greater the energy dissipation caused by the breaking flow is, the more difficult it is to form a wave oscillation. Therefore, wave oscillation must be formed in a certain upstream-to-downstream ratio range. From the current comparison, it can be locked at < 2/3.

#### The same downstream initial depth with different upstream initial depths

From Figure 6(e)–6(j), for the same downstream initial depth, the greater upstream initial depth leads to the higher maximum wave height, the earlier wave arrival time and the higher final stable water level. The maximum wave height is always greater than the upstream initial depth. For the downstream initial depth = 0.30 m, shown in Figure 6(j), the smaller the upstream-to-downstream depth ratio is, the larger the wave oscillation frequency is, and the smaller the maximum wave height is; that is, the smaller upstream initial depth results in the longer wave oscillation process but smaller wave height. In contrast, the larger the upstream initial depth is, the larger the maximum wave height is, and the shorter the oscillation process lasts. Basically, it can be seen that a distinct wave oscillation process has not been formed in Figure 6(j) for = 0.25 m due to > 2/3. The possible reason is that the stronger upstream inflow impacts the downstream reservoir, the turbulent mixing is stronger, and the energy dissipation rate is higher, which in turn causes faster energy dissipation and faster wave attenuation. In Figure 6(g)–6(i), the same results as above can be clearly observed.

#### Run-up heights

In Figure 7, for the same upstream initial depth, the downstream initial depth increases, the maximum run-up height (i.e. the highest water level that the wave can reach from the bottom datum) increases, because the water cushion provided by the downstream initial water allows the flood wave to run up the dam higher. On the other hand, for the same downstream initial depth, the higher upstream initial depth leads to the higher run-up height because of the bigger impact force on the downstream. The maximum run-up heights are positively correlated with the upstream and downstream initial depths. Obviously, in Figure 7, the run-up heights calculated by the RANS agree better with the measured results than those calculated by the SWEs; all the results calculated by the SWEs are smaller than those calculated by the RANS and measured results.

*R*were nondimensionalized with the upstream initial depths, and the relationship between the dimensionless run-up heights and the downstream-to-upstream depth ratios is shown in Figure 8(a), the relationship can be expressed as Obviously, the results calculated by the SWEs are all smaller than those calculated by the RANS and measured results in Figure 8(a), which are consistent with the results shown in Figure 7. In addition, the errors between the dimensionless run-up heights calculated by the RANS and the measured values are within ±10%; however, the errors between the dimensionless run-up heights calculated by the SWEs and the measured values are within −10% to −25%, which are presented in Figure 8(b). So all the results shown in Figures 7 and 8 indicate that the SWEs underestimate the run-up heights with maximum error of −25% and the RANS estimate the maximum run-up heights relatively reasonably.

## CONCLUSION

This paper investigates the cascading dam-break flood evolution process on the sloping channel connected by upstream and downstream dams when the upstream dam collapses instantaneously. The water surface profiles, stage hydrographs and maximum run-up heights calculated by the RANS and SWEs for different upstream and downstream initial depths are compared with the experimental measurements detailed. The comparison of water surface profiles shows that the short wave is difficult to be reproduced by the SWEs due to the hydrostatic pressure assumption. Only the average change trend can be predicted. From the results of stage hydrographs, the RANS predict more accurate results than the SWEs in terms of both the maximum wave heights and wave-like oscillation evolution processes. Only the final stable water levels calculated by the SWEs are in good agreement with the RANS calculation results. The SWEs predicted results correspond more favorably with the RANS and the measured results for the smaller downstream and upstream initial depths. A greater upstream initial depth leads to a higher maximum wave height, earlier wave-arrival time and higher final stable water level. A greater downstream initial depth results in a slightly smaller maximum wave height and a later wave-arrival time. It can be seen that the formation of short waves must be below a suitable upstream-to-downstream initial depth ‘critical ratio’. If the ratio is higher than the ‘critical ratio’, it is difficult to form a wavy push due to the strong breaking, air entrainment and turbulence. The results of this study can be roughly seen as < 2/3. The calculation results show that the maximum run-up height can increase with the increase of the initial depth in both the upstream and downstream. The maximum run-up/wave heights calculated by the SWEs are lower than both the measured results and those calculated by the RANS. The errors of the dimensionless run-up heights calculated by the RANS are within ±10%, while the errors by the SWEs are between −10% and −25%. Therefore, the SWEs must lead to prediction delay of the downstream overtopping dam-failures. This paper especially focuses on the difference between the calculation results of RANS and SWEs, and fills in the widely ignored knowledge gap that a large deviation can occur when the traditional SWEs are used to calculate the run-up height of cascading dam-failure. This study also contributes to the understanding of the evolution mechanism of cascading dam-break and the analysis of flood risks, and it has direct guiding significance for the downstream overtopping dam-break. However, the more complicated and catastrophic downstream cascading dam-failures are not considered, which can be further tried in the future.

## ACKNOWLEDGEMENTS

This work was supported by the National Natural Science Foundation of China (Grant No. 51739011, 51909185), the National Key Research and Development Plan (2016YFC 0402707-03) and China Postdoctoral Science Foundation (No. 2019M652550).