International Journal of Hydraulic Engineering

p-ISSN: 2169-9771    e-ISSN: 2169-9801

2016;  5(1): 9-18

doi:10.5923/j.ijhe.20160501.02

 

Numerical Investigation on 1D and 2D Embankment Dams Failure Due to Overtopping Flow

Omid Saberi 1, Gerald Zenz 2

1PhD Student, Institute of Hydraulic Engineering and Water Resources Management, Stremayrgasse 10/2, A-8010 Graz, AUSTRIA

2Prof, Institute of Hydraulic Engineering and Water Resources Management, Stremayrgasse 10/2, A-8010 Graz, AUSTRIA

Correspondence to: Omid Saberi , PhD Student, Institute of Hydraulic Engineering and Water Resources Management, Stremayrgasse 10/2, A-8010 Graz, AUSTRIA.

Email:

Copyright © 2016 Scientific & Academic Publishing. All Rights Reserved.

This work is licensed under the Creative Commons Attribution International License (CC BY).
http://creativecommons.org/licenses/by/4.0/

Abstract

This article provides three numerical investigations on the overtopping failure of embankment dams which are modelled with non-cohesive fill material. The first investigation is based on a one-dimensional approach in order to calculate the outflow hydrograph during dam overtopping failure. The written program calculates the flow parameters by solving the one-dimensional Saint-Venant equation. Furthermore, the bed evolution is calculated by solving the Exner equation with the finite difference scheme. Here, the model accuracy is increased by dividing bed slopes into three categories: zero slope (0.0-1%), mild slope (1%-20%) and steep slope (>20%). At each time step based on the real bed slope the bed load transport is modelled. In the second investigation the two dimensional open-source TELEMAC software has been modified in order to model embankment dam failure. This modification is done by calibrating slope correction formula in the sediment part of this software. Finally, the influence of different geometrical parameters on the outflow hydrograph of dam failure is numerically modelled. Our modelling results show that the downstream slope has significant influence on the outflow hydrograph in comparison to the other geometrical parameters.

Keywords: Embankment dam, Overtopping failure, Numerical modelling, Steep slopes

Cite this paper: Omid Saberi , Gerald Zenz , Numerical Investigation on 1D and 2D Embankment Dams Failure Due to Overtopping Flow, International Journal of Hydraulic Engineering, Vol. 5 No. 1, 2016, pp. 9-18. doi: 10.5923/j.ijhe.20160501.02.

1. Introduction

Dam failure is a catastrophic event and can cause loss of life and property damage to the structures and properties below the dam. According to the International Commission on Large Dams (ICOLD [7]) a large number of dams worldwide have been failed as embankment dams in comparison to the other types of the dams. Within this type of dam, overtopping failure is one of the most common causes of failures in comparison to the other types of failures like piping and slope failure. Therefore, in hydraulic engineering, water resource management and risk management fields, overtopping failure is categorized as more dangerous failure in comparison to other types of failures. In this paper, dam failure due to overtopping is reassessed.
Generally, dam failure parameters can be calculated either by empirical or numerical methods. In empirical methods, dam failure parameters can be predicted by finding relationships between previous dam failures parameters by using regression analysis. The main disadvantage of this method is ignoring many site specific parameters like grain size, roughness, density, cohesion, velocity, stress, etc. Therefore, using these simplifications, can lead to errors and the necessity to include more appropriate failure parameters.
In the numerical methods dam failure parameters are calculated by solving the governing equations of fluid like Navier–Stokes, shallow water or Saint-Venant. The main advantage of this method is considering more reasonable parameters based on mechanical relationships like cohesion, friction, etc. This means if appropriate parameters are available, the numerical method is more accurate in comparison to the empirical method. Numerical modelling herein is performed in two steps: first to solve the hydrodynamic part, and second the erosion part. In the hydrodynamic part, flow parameters like water depth and velocity of water in downstream of the dam is determined by using one-dimensional Saint-Venant equation. In the erosion part by using the prior calculated flow parameters, a bed evolution during failure process is calculated with the sediment continuity equation. However, the final results of these two steps show that this is not accurate. The main problem arises from the formulation of the bed load equations. Available bed load equations are formulated for uniform flow with mild slopes. While in the embankment dam failure, flow regime is unsteady and the bed slope changes from steep to mild rapidly. This phenomena is investigated in this study and as a solution, a one-dimensional program is written in Fortran 90. Moreover, in the erosion part three different equations are used. The original form of Meyer-Peter & Muller [13] equation is used for modelling zero slopes, the SMART [18] equation is used for modelling slopes up to 20% and for modelling slopes more than 20% the modified form of SMART [18] equation is calibrated and used in this program.
The second investigation in this paper is done on the sediment part of the open-source TELEMAC software (SISYPHE software) to model embankment dam failure. This modification is done by calibrating slope correction formula in the SISYPHE software.
The accuracy of the 1-D program and the TELEMAC-2D software are validated and compared to each other by using L. Schmocker [17] small scale dam failure tests. In the last part of this paper the influences of geometrical parameters on the breach outflow hydrograph is numerically investigated and compared to each other. The final results show that the downstream slope of the embankment dams have more influence on the breach outflow hydrograph in comparison to the other geometrical parameters of the embankment dams.

2. One-Dimensional Dam Failure Modelling

Dam failure parameters as failure time and peak outflow discharge are determined by using either empirical or numerical methods. Empirical methods use the data from the failure dam histories to predict dam failure parameters. However, in numerical methods by solving governing equations of flow and considering a wider range of scenarios and parameters results based on physical processes are gained. Furthermore, the numerical approach is normally performed in two steps:
- hHydrodynamic part
- Erosion part
In the following of this paper these two steps are explained in more detail.

2.1. Hydrodynamic Part

The hydrodynamic part in the numerical approach for dam failure modelling concerns on one-dimensional flow parameters like water depth and its velocity. These parameters can be determined by solving the Saint-Venant equation. This one-dimensional equation is valid under four main assumptions:
1. The flow velocity and water depth change only in the direction of flow.
2. Vertical acceleration is neglected and hydrostatic pressure is considered along the channel.
3. The density of fluid is constant along the channel and fluid is incompressible.
4. Slope of the channel is small.
The conservation form of this equation can be described as the following equations (Cunge [3]):
(1)
Here, is the hydrostatic pressure. However, this complex form can be simplified if a rectangular channel with constant width along the channel are assumed. This means that can be written as equation (2):
(2)
There are different schemes to solve this equation like MACCORMACK Scheme, FTSC Scheme, Lax Scheme and Leap-Frog scheme. In this study the revised MACCORMACK scheme (Garcia Navarro [4]) is selected for solving the Saint-Venant equation. The main advantages of this method in comparison to the other methods are:
1. The MACCORMACK scheme has two steps, predictor step and corrector step which is capable of capturing the discontinuities in the flow.
2. The MACCORMACK scheme has higher accuracy because of using two differential equations both in space and time.
3. In the MACCORMACK scheme the primary results which are determined during the predictor part, are used during the corrector part as shown in the Figure 1.
Figure 1. The MACCORMACK Scheme
In general, in the MACCORMACK scheme, the Saint-Venant equation is solved in two steps:
- Predictor step
- Corrector step
In the predictor step the forward finite difference method normally is used to calculate the primary flow parameters as shown in the equations (3) and (4):
(3)
(4)
Here, are bed slope and friction slope and are defined as the following equations:
(5)
In the corrector step the backward finite difference method is used to calculate the secondary flow parameters as shown in the equations (6) and (7):
(6)
(7)
Here, are bed slope and friction slope and are defined as the following equations:
(8)
Finally, the value for each hydraulic parameter is determined by averaging between the results from these two steps (predictor and corrector):
(9)
The above methodology, describes a procedure to solve the Saint-Venant equation, the stability of this calculation should be controlled by using the Courant-Friedrich-Lewy (CFL) formula. In general, calculation get unstable when the CFL number become larger than 1.0. Therefore the value of the CFL number must be less than 1.0 or 1.0 to have stable calculation. The CFL number is expressed as equation (10):
(10)
(11)
Here, is the wave velocity, is the flow velocity, A (m2) is the cross section area and Z (m) is the bed elevation.
Now, the only unknown parameters remaining in the Equation (9) are the manning number and the boundary condition values. In the following sections these parameters are determined and explained in more detail.
2.1.1. Boundary Condition
In the MACCORMACK method like other explicit methods, all computational nodes can be calculated except two nodes, the first node (n=1) and the last node (n=N). These nodes represent in total four unknown parameters which are related to the height of water and water velocity at points 1 and N. In order to calculate these unknown parameters a new method has been used in this paper. This method was suggested by Garcia Navarro [4]. In this method unknown parameters are estimated by using characteristics and linear interpolation approaches.
2.1.2. Roughness Computation
In open channel, the total flow resistance results from the interaction between different elements which are located on the bed of the channel. Among these elements, some of them have more influence on flow resistance like grain size, vegetation, bed slope, bed aliment and obstruction in the channel. In fact, roughness coefficient shows the effect of these parameters in stream flow. The importance of calculating roughness coefficient is related to define exact value for velocity in the manning equation. In this investigation the following equations are used for calculating roughness coefficient:
The Ghani [5] equation, for mild slope
(12)
The Jarrett [8] equation for steep slope
(13)
Here, R (m) is the hydraulic radius, H (m) is the water depth and d50 (m) is the mean diameter.
2.1.3. Validate Hydraulic Part
In this part, in order to verify the accuracy and performance of the MACCORMACK scheme, the Goutal and Maurel [6] test is employed. This test was carried out in a frictionless rectangular channel with 25m length and 1.0m width. In the middle of this channel one bump has been considered. The bed topography for this bump is illustrated in the Equations (14) and (15).
(14)
(15)
In this part, the model accuracy is demonstrated by using the transcritical Flow over bump (supercritical with hydraulic jump) and the final results are compared against the analytical results.
Figure 2. Geometry of bump test
Transcritical Flow over bump
In this test the downstream water level is kept constant at 0.33 m and inflow discharge per unit width is kept constant at 0.18 (m3/s)(1/m). The final results are shown in the figure 3.
Figure 3. Supercritical with hydraulic jump
Figure 3, shows the MACCORMACK scheme has acceptable accuracy for modelling hydraulic jump in the open channels.

2.2. Erosion Part

The second step in the numerical modelling of the embankment dam failure is calculating the surface erosion. In this step the bed elevation is defined by solving the one-dimensional sediment continuity equation. This equation is expressed as the following equation.
(16)
Equation (16) is known as the Exner equation. This equation can be discretized by using the Modified –Lax scheme as shown in the following equations:
(17)
(18)
By substituting equations 17 and 18 into equation 16, the bed elevation for time step (j+1) can be determined as:
(19)
Here, is the porosity and is the stability parameter and can be determined, based on Vreugdenhil [20] equation:
(20)
(21)
Now, the only unknown parameter remaining in the equation (19) is the bed load parameter In order to find this parameter many different experimental equations with different accuracy are available. Generally, all available bed load equations are formulated for uniform flow with mild slopes. While in the embankment dam failure, flow regime is unsteady and the bed slope changes from steep to mild rapidly. Consequently, using these equations for modelling embankment dam failures may show excessive errors in the final results. In this study to overcome this deficiency the bed slope parameter in the SMART [18] equation is calibrated. Furthermore, in order to reduce errors during calculation different bed load equations are selected and applied for slopes less than 20%. The applied equations for different slopes are discussed briefly in the following:
For zero slopes, original form of the Meyer-Peter & Muller [13] equation is used:
(22)
For slopes up to 20% slopes original form of the SMART [18] equation is used:
(23)
To model beds with slopes more than 20%, the calibrated form of the SMART [18] equation is selected:
(24)
Here, S0 is the bed slope, is the relative density, θ is the shield parameter, is the critical shield parameter and ω is the calibrated parameter.
In this paper to calibrate ω parameter, the Chinnarasri [1] dam failure tests are selected. The final value for ω parameter is shown in the Table 1.
Table 1. Calibrated
      parameter (for slope more than 20%)
     

3. TELEMAC MASCARETE

The second investigation in this paper is related to updating the sediment part of the TELEMAC MASCARETE software. This system has been developed by the Department Laboratorie National d´Hydrauliquie (LNH) at the Electronic de France Direction des Etudes et Recherché (EDF-DER).The software system is owned by EDF-R&D and this software available as the open-source software which is designed for modelling free surface fluid, flood wave propagation, ground water flow and sediment transport in the open channels and rivers. The TELEMAC MASCARETE system is able to solve the Shallow Water Equation for two dimensional modelling (TELEMAC-2D) and the Navier–Stokes Equation for three dimensional modelling (TELEMAC-3D). For modelling erosion and sediment transport in the rivers the SISYPHE software is used in TELEMAC system.
The SISYPHE software can simulate sediment transport and bed elevation in the complex morphology same as coastal, rivers, lakes and estuaries with the different discharge rate, different sediment grain size and different sediment transport equations. The SISYPHE software can be easily coupled with the TELEMAC-2D and TELEMAC-3D software. In this coupling, at each time step TELEMAC-2D or 3D send the calculated hydrodynamic parameters like height of the water (H) and water velocity (U, V) to the SISYPHE software. The SISYPHE software model bed elevation by solving the two-dimensional sediment continuity equation which is called the Exner equation as shown in the Equation 25 (TASSI [21]).
(25)
Here, is the bed load transport, z (m) is the bed elevation and p is the bed porosity.
Furthermore, SISYPHE software use the Meyer-Peter-Müller, Einstein-Brown and Hunziker equations to determine bed load transport To improve the accuracy of the final results, the SISYPHE software uses two correction formula which are called Koch and Flokstra [9] and Soulsby [10] equations. In this paper, in order to model sediment transport over the steep slopes, the βTel parameter in the Koch and Flokstra [9] formula is calibrated by using the Chinnarasri [1] dam failure tests as shown in the equations (26) and (27):
Koch and Flokstra formula (1981):
(26)
(27)
Here, is the correction formula, is an empirical parameter, (m3/s) is the bed load transport and (m3/s) is the corrected bed load transport.
The calibrated value for βTel parameter is shows in the Table 2.
Table 2. Calibrated βTel parameter
     

4. Validate 1-D Program and TELEMAC MASCARETE Software

In order to verify the accuracy and performance of the 1-D program and TELEMAC-2D software, two experimental tests are selected from the Schmocker [17] dam failure tests. The Schmocker tests are carried out in a glass-side flume with 8.0m length, 0.2m width and 0.7m height. The small scale dike is installed at 4.0m distance from the channel intake with 0.2m height and 0.1m crest width. The upstream and downstream slopes of this dam are fixed at 1V:2H. Furthermore, these dam failure tests have been done based on two different dam materials:
- Homogeneous sand with mean sediment diameter 2.0mm.
- Homogenous sand with mean sediment diameter 0.31mm.
In both materials the sediment is non-cohesive with density of 2650 kg/m3. The inflow discharge also was kept constant at 6.0 l/s during failure process.The results for 1-D program and TELEMAC-2D software are shown in the Figures 4 and 5.
Figure 4. Comparison between measurement, 1-D and TELEMAC-2D results for d=2.0 mm
Figure 5. Comparison between measurement, 1-D and TELEMAC-2D results for d=0.31 mm
These figures confirm that the applied calibration approach in both software have a good accuracy for modelling slopes more than 20%.

5. Influence of Initial Breach

One of the major differences between one and two dimensional numerical modelling in embankment dam failure is consideration of the initial breach. In order to find the influence of the mentioned issue on the final results, two different scenarios have been made on the Morris [14] dam failure test.
- First scenario, initial breach is considered in the center of the crest.
- Second scenario, no initial breach is considered on the dam crest (Figure 6).
Figure 6. Influence of initial breach on the outflow hydrograph
These two scenarios are modelled numerically with the TELEMAC-2D software and 1-D program and the final results are compared against each other.

5.1. IMPACT Project (Test No.2)

The IMPACT project (Morris [14]) is one of the well-known large scale dam-break projects which has been done in Norway in 2002. This embankment dam was built mainly from non-cohesive soil with D50=4.75mm. The main purpose of this test was to have a better understanding of the failure mechanism in the homogeneous non-cohesive embankment dams which are failed by the overtopping flow. Table 3 shows the more details about this test.
Table 3. List of geometrical parameters
     

5.2. Compare Results of 1-D Program and TELEMAC-2D Software

In this part, the modelling results of the one-dimensional program and the TELEMAC-2D software, are compared against to the measurement results, these comparisons are shown in the Figure 7.
Figure 7. Influence of initial breach on outflow hydrograph in 1-D and 2-D software
Figure 7 shows the outflow hydrograph without initial breach gives higher value in peak outflow discharge and lower value in failure time compared with the outflow hydrograph with initial breach.

6. Influence of Upstream Slope, Downstream Slope and Crest Width

In this part, the focus mainly is on the influence of dam’s geometrical parameters like upstream slope, downstream slope and crest width on the outflow hydrograph calculation during dam failure process. In order to define which geometrical parameters have a significant influence on the breach outflow hydrograph the Chinnarasri [1] dam failure test is chosen. In the following, this dam failure test is explained in more detail.

6.1. Chinnarasri Dam Failure Tests (Test No.1)

This experimental test was performed in a flume with dimensions 35m×1m×1m. A small scale dam is considered in the middle of this flume with 0.8m in height and 0.3m crest width. The upstream slope of this dam is fixed at 1V:3H and downstream slope of this dam is fixed at 1V:2H. For all tests, the initial water level in the reservoir is kept constant at 0.83m and downstream water level is kept constant at 0.03m. The soil porosity is 0.35 and soil density is 2.65 103 kg/m3. This dam failure test is numerically modelled by modifying shape parameters as expressed in Table 4 by using TELEMAC-2D software. Herein, the final results are compared against each other within the given dataset and are shown in the Figures 8, 9 and 10.
Table 4. List of geometrical parameters
     
Figure 8. Influence of different downstream slope on outflow hydrograph
Figure 9. Influence of different upstream slope on outflow hydrograph
Figure 10. Influence of different crest width on outflow hydrograph
Figures 8, 9 and 10 show, the downstream slope of the embankment dam has more influence on the outflow hydrograph in comparison to the other geometrical parameters.

7. Discussion

Our sensitivity analysis in the previous part (section 6) reveals that the downstream slope in the embankment dams has more influence on the breach outflow hydrograph in comparison to the other geometrical parameters like upstream slope and crest width. The main reason is related to the high influence of the downstream slope on water velocity and bottom shear stress during failure process as being illustrated in the Figure 11.
Figure 11. Water velocity and energy line in downstream slope
From the above figure, it can be concluded by decreasing downstream slope, water velocity and thus the erosion rate of the dam's material are decreased. This means that the risk of dam failure is reduced significantly.
The results of this investigation can be validated with other work studies. The most important one is the Minimum energy loss (MEL) method which is based on the work of Mackay [23]. In the following MEL method is explained in more detail.

7.1. Minimum Energy Loss (Mackay [23])

The Minimum energy loss (MEL) for the first time was designed and introduced by Mackay [23] and is one of the unusual overtopping protection designs. It was first built on the Chinchilla weir in Australia and it is still in use without any damage. The Minimum energy loss (MEL) is working based on the reduce energy of large flood during passing on the embankment dams to prevent any erosion and damage by decreasing downstream slope as shown in the Figures 12 and 13.
Figure 12. Side view of Minimum energy loss (MEL)
Figure 13. Plan view of Minimum energy loss (MEL)
Figure 12 reveals that this method is in agreement with our conclusion which was mentioned in the previous section (section 6).
(28)
(29)
Here:
(H1-z) (m) is the upstream head above spillway crest
is the crest width (Figure 1-7)
is the design upstream head
Z1 (m) is the any elevation above the toe
B (m) is the width of the channel at Z1 elevation

8. Conclusions

In this paper a one-dimensional program is introduced and developed in Fortran 90 to model non-cohesive embankment dam which fail by overtopping flow. Additionally the two-dimensional open-source TELEMAC software has been applied and modified in order to model embankment dam failure. This modification is done by calibrating the slope correction formulas in the sediment part of the TELEMAC software. During validation of these two software, some investigations have been done on geometrical parameters to find out the influence of these parameters on the breach flow hydrograph of embankment dam failure.
In the following the most important results are listed:
- Computational stability in 1-D program is depending on the Courant number. Best result can be reach if the Courant number is less than 1.0.
- In 1-D program two basic assumptions are considered. The first assumption is that whole length of the crest is eroded at same time. In the second assumption, rectangular shape with constant width is considered to model the reservoir of a dam.
- Change of the downstream slope in the embankment dam failure has more influence on the outflow hydrograph in comparison to the other geometrical parameters.

Notation

The following symbols are used in this paper

References

[1]  Chinnarasri, C. and Tingsanchali, T. and Weesakuli S. and Wongwises, S. (2003). Flow patterns and damage of dike overtopping. International Journal of Sediment Research, 301-309.
[2]  Costa, J. E. (1985). “Floods from dam failures.” U. S. Geological Survey Open-File Report 85 560, Denver, CO, 54 p.
[3]  Cunge, J. A., Holly, F. M., Verwey, A. (1980). Practical Aspects of Computational River. London: Pitman Publishing Ltd.
[4]  Garcia-Navarro, P. and Saviron, J., M. (1992). McCormack's method for the numerical simulation of onedimensional discontinuous unsteady open channel flow. Journal of Hydraulic Research, 95-105.
[5]  Ghani, A. A. (2007). Revised equations for Manning's coefficient for sand-bed rivers. International Journal of River Basin Management, Vol. 5, No. 4, 329-346.
[6]  Goutal, N. and Maurel, F. (1997). Proceeding of the second workshop on dam-break wave simulation. Laboratoire National d' Hydraulique Chatou, Electricité de France.
[7]  J.F.A. SILVEIRA. (2011). SMALL DAMS Design, Surveillance and Rehabilitation. International Committee on Large Dams (ICOLD 2011), 30-34.
[8]  Jarrett, R. (1984). Hydraulics of high-gradient streams. American Society of Civil Engineers, Journal of Hydraulic Engineering, HY11, v. 110, p. 1519-1539.
[9]  Joe D. Hoffman, Steven Frankel. (2001). Numerical Methods for Engineers and Scientists. New York: Marcel Dekker.
[10]  Koch, F.G. and Flokstra, C. (1981). Bed level computations for curved alluvial channels. XIXth Congress of the International Association for Hydraulic Research, New Delhi India.
[11]  Lax, P, D, and Wendroff, B. (1964). Difference schemes for hyperbolic equations with high order of accuracy. Comm. Pure appl. Math.
[12]  Mac Cormack, R. (1969). The Effect of Viscosity in Hypervelocity Impact Cratering. AIAA Paper, 354.
[13]  Meyer‐Peter, E. and R. Müller. (1948). Formulas for bed load transport. paper presented at 2nd Meeting, International Association for Hydraulics research, 39-64.
[14]  Morris, M. (2006). “IMPACT Project Field Tests Data Analysis.” HR Wallingford, UK.
[15]  P. Garcia-Navarro and J. M. Saviron. (2010). McCormack's method for the numerical simulation of one dimensional discontinuous unsteady open channel flow. Journal of Hydraulic Research, 95-105.
[16]  Parker, G., and Paola, C., and Leclair, S. (2000). Probabilistic Exner Sediment Continuty Equation for Mixtures with No Active Layer. Jounral of Hydraulic engineering, 818-826.
[17]  Schmocker, L.r; Hager, W. (2012). Plane dike-breach due to overtopping: effects of sediment, dike height and. IAHR, 576-587.
[18]  Smart, G. (1984). Sediment Transport Formula For Steep Channels. ASCE, 267-276.
[19]  Soulsby, R. (1997). Dynamics of Marine Sands. London: Thomas Thelford Edition.
[20]  Vreugdenhil, C. B. and Wijbenga, J. H. A.. (1982). COMPUTATION OF FLOW PATTERNS IN RIVERS. ASCE, 87-90.
[21]  TASSI, P. (2014). Sisyphe v6.3 User's Manual. EDF R&D
[22]  Chanson, H. (2013). Embankment overtopping protection systems. Springer-Verlag Berlin Heidelberg, 305-318.
[23]  MacKAY, G. (1971). Design of Minimum Energy Culverts. Brisbane, Australia: Dept of Civil Eng., Univ. of Queensland.