FLOW ANALYSIS AND SCOUR HOLE COMPUTATION OF DYKE-BREACH 

 

Zhang Xiuzhong and Wang Guangqian

Key laboratory for Water and Sediment Sciences of Ministry of Education, Tsinghua University, Beijing, 100084, China

 

 

Abstract: Based on the 2D shallow water equations and the nonequilibrium suspended sediment transport equation, a mathematical model is presented to predict flow motion and scour hole evolution of an instantaneous dyke- or dam-breach. To simulate the complex flow phenomenon accurately and conveniently, a high-resolution finite element scheme is applied to capture the dam-break bore. Flow analysis and scour computation are conducted for a variety of breach widths, discharges and soil composition. The results obtained are used to provide hydrodynamic conditions for structural and stability analysis of the steel, timber and earth-rock combined dam. 

Keywords: dike breach, scour hole, high-resolution finite element, numerical simulation

1    INTRODUCTION

Levee or dam failures during floods often cause tremendous loss of lives and properties. For example, a 164-m long breach of the Putuo River in Hebei province of China resulted in flooding of 23 villages and 18 thousand hectares of farmland by the released 230 million cubic meter water in August 1996 [1]. Another failure of the anti-flood wall of the Jiujiang City by the Yangtze River in 1998 not only aggravated flood disaster but drove the flood control work into a severer situation. Whether the breach can be closed up quickly or not, and how high cost will be paid all depend on the hydraulic characteristics and the scour development of the breach. Therefore, further researches on problems of the dyke-breach play a vital role in preventing, reducing, and harnessing the flood disaster.

There are two major causes for a dike failure [2]. One is piping which causes collapse of levee and gives rise to overflow. The other is that flood stage exceeds the top of the levee to lead to overtopping. Both of the two kinds of flows are open channel flow, quite similar to the dam-break flow. Thus the dike-break flow can be described by the shallow water equations (also named the Saint-Venant equations for the 1D case) mathematically. After the dike fails with breaches, the bed at or near the breach will be continuously eroded to form scour hole. As the scour hole becomes bigger and deeper, it will further endanger the safety of the levee, and might cause a wider range of dyke breaching. So exactly predicting the flood routing and evolution of the scour hole downstream the breach will undoubtedly direct the implementation of flood prevention.

With the development of computer powers and the advance in numerical techniques, numerical simulation of dike-break problems is not only feasible but efficient. However, due to the convection-dominated and degenerative hyperbolic features of the shallow water equations, especially for the complex structure of the dyke-breach flow, with subcritical, supercritical flows, and even discontinuities either presenting simultaneously in different parts or occurring in sequence at different time. Thus how to overcome these numerical difficulties is the key to simulating the dyke-break flow.

So far, calculating results of the eroding process near the breach are rarely available in the literature. This paper is devoted to the study of flow and scour of the dyke breach. For that, a high-resolution FEM is adopted to track the dam-break bore [3]. Assuming that the bed material is noncohesive, and its motion satisfies the sediment continuity equations, so that the nonequilibrium sediment transport model can be applied to predict the motion of the sediment and the evolution of the scour hole.

2    MATHEMATICAL MODEL

2.1    Governing equations

It is commonly accepted to consider the dam-break flow as a shallow water phenomenon, and its solution with jump or shock discontinuities is said to constitute a weak solution of the conservation law. Under the hydrostatic approximation, the simplified N-S equations to describe incompressible free surface flow are as follows

                         (1)

              (2)

The nonequilibrium sediment transportation and bed deformation are described by

                   (3)

                         (4)

where ; and  = depth averaged velocity and unit discharge in  direction, respectively; C = the Chezy resistance coefficient; Eddy viscosity is estimated by , in which = Von Karman’s constant, = friction velocity; ; ; ; , and  = depth integrated concentration, sediment-carrying capability and setting velocity of the kth grain size fraction; = recovery coefficient; = dry bulk density of the bed material;  = sediment diffusive factor assumed to be the same as  in this paper; water surface elevation  is defined by water depth  and bed surface level .

2.2    Gradation of sediment-carrying capacity

In view of sediment either from upstream flow-carrying, or turbulence diffusion entering or the both, it is reasonable to calculate the gradation of sediment carrying capacity by considering both the flow condition and size distribution of the bed material simultaneously. So the treatment approach proposed by Li [4] is adopted.

2.3    Calculation of bed material gradation

Given the erosion or deposition depth of the kth fraction sediment , and the total erosion or deposition depth , the computation of the bed material gradation can be divided into two cases:

(1) , which means deposition occurring, the gradation of the mixed layer can be estimated by

where , = gradations of the mixed layer at the time level  and , respectively;  and  are the corresponding thickness of the mixed layer.

(2) , which means that erosion happened, the gradation of the mixed layer is then computed by

where = the average gradation in some sub-layers of the memory layer.

To simulate the armoring phenomenon of the bed material in degradation, two layers are incorporated in the model. The upper layer of the bed material is called mixed layer and the lower is named memory layer. The thickness of the mixed layer is , and its gradation is . The memory layer is divided again into many sub-layers, say n sub-layers, and the thickness and the gradation of each sub-layer are , , respectively. The number of the sub-layers of the memory layer will increase to n+1, and the gradation of the newly added sub-layer equals to that of the mixed layer in time level t, if the riverbed aggrades. Otherwise, if the riverbed degrades, the number of the sub-layers will subtract some layers, and the gradation of the memory layer should be adjusted correspondingly.

2.4    Finite element discretization and solution

The FEM is adopted due to its advantage to treat practical problem with complex geometries. Discretization for spatial variables and integration by part leads to

                           (5)

                       (6)

,

where M = lumped mass matrix; C = convective matrix; D = diffusive matrix; F = source term. They are expressed as follows.

 ;

 ;

;

Where N, W are interpolation and weight function, respectively.

When WN, it constitutes the Galerkin FEM which equals to central difference, so numerical oscillation will be created when solving convection-dominated problems. In this study, a high-resolution finite element scheme is incorporated to rearrange the convective term. Its main characteristic is that the scheme is obtained through the introduction of almost equal amount of diffusion and anti-diffusion, with limiter to preserve monotunity condition, so that both the requirement of high-order of accuracy and non-oscillatory condition are satisfied.

For the sake of convergence or fasting convergence, the source term F is transposed and linearized. After discretization for spatial variables, Eqs.(5), (6) are a system of ordinary differential equations, which can be solved by various explicit and implicit schemes. The Runge-Kutta time integration is used herein for unsteady flow. A sequential procedure is employed to solve the discretized equations, where Eq.(5) is a 2D convective equation for computing h; Eq.(6) is a 2D convection-diffusion transport equation to calculate u, v and s; Erosion depth can be obtained by Eq. (4).

2.5    Verification of the computational scheme

To verify accuracy, and allow comparisons with analytical solutions, the scheme is used to simulate a dam-break flow. The problem considered involves sudden and complete failure in a horizontal rectangular channel without friction. Initially the water is rest with an upstream depth of 10 m and a downstream of 1 m. Comparisons with stoker’s solution [5] at 120 second are displayed in Fig 1 and Fig 2.

3    RESULTS ANALYSIS FOR FLOW AND SCOUR ELUTION

Calculations concerning both dam-break and dike-breach are carried out. For the dam-break case, attention is primarily paid to considering flow and scour development near the breach under conditions that the inlet water surface elevation remain unchanged and all flow pass the breach. As for the case of dike-burst, focus is mainly put to observe flow and scour evolution under conditions that the channel diverts the flood and the upstream water stage is changeable.

3.1    Results analysis for dam-break

The physical region is shown in Fig 3. The horizontal channel is 2500m in length, 1000m in width. The dam is assumed to be 30m thick in the direction of flow. The Manning roughness is set to be 0.015, and the median size of the bed material is 0.3 mm. Initially the water is rest with an upstream depth of 7.5m and a downstream of 0.5m. As to the boundary conditions, water stage is specified to be 7.5m at the upstream, non-slip conditions are satisfied at and upstream the breach; at the downstream Neumann condition is applied to the variable u.

Three cases with the breach width of 50m, 100m, and 200m are simulated, respectively. Confined to the length of the paper, only some representative result figures are incorporated. More figures are available in the report [6]. With the breach width of 50m, Fig. 5 shows the varying process of the scour depth at the central line of the breach. At the same wide breach, Fig. 6 and 7 give the scour depth contours at different time. Fig. 8 and 9 present the velocity vector when the breach is 100m wide.

To demonstrate effects of the breach width on scour depth, and to investigate the relation between the sour depth and the breach width, the maximum scour depth and its location are listed in Table 1.

Table 1    Erosion results of different widths of the breach

 

maximum erosion depth (m)

location of maximum erosion (m)

width (m)

50

100

200

50

100

200

100s

0.091

0.16

0.14

96.0

120.0

165.0

600s

0.633

1.15

1.30

96.8

119.0

165.0

1200s

1.181

2.01

2.53

78.0

98.0

142.0

1800s

1.640

2.67

3.50

60.0

80.0

120.0

3600s

2.75

4.1

5.60

45.0

80.0

99.0

7200s

4.375

6.1

 

45.0

62.0

 

Notes: D50=0.3 mm, numbers in the first column is time.

It can be seen from Tab.1 that as the breach widens, the scour depth increases, and the distance between location of maximum erosion and the beach increases as well. The location of the maximum erosion depth gradually approaches the breach with time going on, till at a stable position at last.

From the scour depth contour we can see that the plane view of the scour hole is elliptical, and this trend is more obvious as time becomes longer. From Tab.1 and figure of scour hole elution we can also draw the following conclusions. The scour rate increased firstly, then dropped, and got zero at last. Compared to its downstream slope, the scour hole has a steeper upstream slope.

3.2    Results analysis for dike-break

As Fig.4 shows, the computational domain consists of the upper stream reach, the flood diversion area below, and the breach zone. The stream is 800m wide with a bed slope and a levee width being 0.0002 and 50m, respectively. Both the breach and the flood diversion area are horizontal. The Manning factor of the three portions is set to be 0.002. As boundary conditions of the stream, discharge process is given at the inlet, and the relation between water stage and discharge, which can be derived by the Chezy formula, is specified at the outlet. The initial conditions are determined by the computation of the channel flow.

At the same discharge and bed material, four cases with the breach width of 100m, 200m, 300m, and 500m are computed. With the same breach width, calculations of different discharges and bed material composition are also included. Similarly, only typical result figures are attached to conserve space. More figures can refer to [6]. With the breach width of 300m, Figs.10, 11 give the velocity vectors, Fig.12, 13 describe the evolution of the scour depth at the central line of the breach. The velocity process at the central line of a 300-meter wide breach is demonstrated in Fig. 14.

To compare the different cases and to explore some useful laws, the maximum scour depths and its locations at time of 7200s are listed in Tab. 2 and 3.

Table 2    Maximum erosion depth of different widths of the breach

Times

7200

Breach width (m)

100

200

300

500

Maximum erosion depth (m)

2.01

2.14

2.0

0.95

Location of maximum erosion depth (m)

50.

85.

86.

86.

Notes: D50=0.25 mm, Q=10000.0 m3/s.

Table 3    Effects of different discharge and sediment on scour depth

D50mm

0.25

0.1

dischargem3/s

10000.

15000.

15000.

Maximum erosion depth (m)

2.0

2.85

3.1

Location of maximum erosion depth (m)

86.

90.

120.

Notes: the breach width is 300 m, and the time is 7200 s.

From Tabs. 2, 3 as well as the scour depth evolution figure, It is not difficult to see that the maximum scour depth increases firstly, then decreases when the breach widens. At about 200m of the breach, the scour depth is deepest. At the same breach width, the larger the inlet discharge of the channel is, the greater will be the scour depth, and vice versa. Moreover for the same inlet discharge and breach width, the scour hole depth will be reduced if coarse bed material is used. In addition, the scour length is proportional to the scour depth. Fig.14 indicates that the velocity at the central line of the breach increases initially and then dropped.

4    CONCLUSIONS

The numerical scheme is validated by comparison with Stoker’s analytic solutions. As the present model is based on unstructured grid, it is very flexible to deal with geometrically complex problems that will occur in practical engineering.

The plane view of the scour hole is elliptic. Contrast to its downstream slope, the scour hole has a steeper upstream slope. There is a linear relation between the scour distance and the scour depth, which manifests that the geometrical similarity exists in the evolution process of the scour hole. The scour rate increases firstly, then drops, and gets zero at last. 

For the dam-break scour with upstream water stage remaining unchanged, the wider the breach is, the deeper the scour hole gets. For dyke-break scour with the stream inflow discharge being constant, the scour depth will have a maximum value when the breach has a certain width. In addition, the evolution characteristics of the scour hole is closely related to the stream discharge and the grain size of the bed material. The velocity is of magnitude near the breach, and it almost approached 10. m/s in the computation. It is a fact that the velocity increases in the beginning and then decreases.

Owing to the complexity of problems of the dike burst, some simplifications or assumptions are made in the model, so improvements and further researches are essentially required.

References

[1]    The 27th Army Group headquarters of Chinese People's Liberation Army. 1998, Report on application research of the breach-stopping technique for steel, timber and earth-rock combined dam. Rep. (in Chinese).

[2]    Gui, S. X., Zhang, R. D., and Xue, X. Z. 1998, Overtopping reliability models for river levee. J. Hydr. Engrg., ASCE, 124(12), 1227-1234.

[3]    Zhang, X. Z., Wang, G. Q., and Jin, S. 2001, Finite-element analysis for shallow water flows and high-resolution FEM scheme. J. Yangtze River Scientific Res. Inst., No.1. (in Chinese).

[4]    Li, Y. T. 1987, Preliminary study on the gradation of bed load in an equilibrium state. J. Sediment Res., No. 1. (in Chinese).

[5]    Stoker, J. J. 1957, Water waves. Interscience Publications, Wiley, New York.

[6]    Tsinghua University, the 27th Army Group. 2000, Hydraulic characteristics and scour hole computation of the dike-break. Rep., Beijing (in Chinese).

                

Fig. 1    Comparison of water surface elevation              Fig. 2    Comparison of velocity

         

Fig. 3   Definition sketch for partial dam break          Fig. 4   Definition sketch for Dyke breach

 

              

Fig. 5    Depth evolution of scour hole                    Fig. 6    Contour of scour hole at t=600 s

 

          

Fig. 7    Contour of scour hole at t=7200 s       Fig. 8    Velocity vectors at t=600 s

 

 

           

Fig. 9    Velocity vectors at t=3600 s             Fig. 10    Velocity vectors at t=100 s 

 

               

Fig. 11    Velocity vectors at t=7200 s         Fig. 12    Depth evolution of scour hole

 

 

 

 

 

        

Fig. 13    Depth evolution of scour hole         Fig. 14    Evolution of the velocity at the

                                               central line of the breach