- Email: [email protected]
A model for simulation of debris flow sedimentation in slit detention-dam reservoirs
Journal of Hydro-environment Research 27 (2019) 65–74
Contents lists available at ScienceDirect
Journal of Hydro-environment Research journal homepage: www.elsevier.com/locate/jher
A model for simulation of debris flow sedimentation in slit detention-dam reservoirs
T
⁎
Mohammad Ebrahim Banihabib , Farzaneh Nazarieh Department of Irrigation and Drainage Engineering, College of Aburaihan, University of Tehran, Pakdasht 20 km, Tehran, Iran
A R T I C LE I N FO
A B S T R A C T
Keywords: Slit dam Sedimentation Two-dimensional Debris flow
Debris flow is a natural phenomenon that requires appropriate structural measures such as detention dams to reduce hazards resulting from it. The objective of this paper is to present a two-dimensional numerical model to investigate debris flow deposition in the reservoir of the detention slit dams. The governing equations include shallow water equations considering the bed level variation, sediment continuity, and debris-flow transport equations. Shallow water equations were discretized using the MacCormack scheme and then Banihabib Scheme was applied for solving sediment continuity equation. The governing equations were solved in an iterative procedure and as a result, flow properties and bed level variations were obtained in silt dam reservoir. The model performance was evaluated against experimental and numerical models in two aspects including the outflow hydrograph and reservoir bed level variations or trapping efficiency. The model performed quite well in outflow hydrograph simulation, and the comparison of the bed level variations or trapping efficiency illustrated that the simulated values were acceptable. A sensitivity analysis of the parameters was carried out to analyze the influence of parameters variation on outflow and sedimentation volume. The results of this sensitivity analysis revealed sedimentation is more sensitive to the parameter variation than outflow hydrograph. Generally, the proposed numerical method can provide insights on the debris flow deposition in slit dam reservoir, which could be useful for slit dam design.
1. Introduction Debris flow leads to infrastructure destruction and loss of lives. It mostly occurs in mountainous areas where unstable sediments exist on steep slopes. Debris flows have received specific attention all over the world such as Japan (Han et al., 2015; Nakatani et al., 2008; Takahashi, 2007), China (Pan et al., 2013; Qiu et al., 2018), Taiwan (Tsai, 2007), the European Alps (Ghilardi et al., 2001; Michelini et al., 2017; Rickenmann et al., 2006) and other part of the world (Hunt, 1994; Paik, 2015; Pirulli and Pastor, 2012). Given the probable occurrence of debris flows and their potential destructiveness, they should be controlled by proper measures. There are different measures for flood control, including flood diversion channels, ring levees, river floodplain management, detention dams and etc (VanDine, 1996). Detention dam, which is a structural measure built on small or moderate rivers, can be quite satisfactory for flood control in short term. However, in the cases that debris and sediments fill the reservoir of the detention dam, the capability of the dam to alleviate floods may be decreased, and dam performance is challenged (Banihabib and Bahram, 2009; Hassan-Esfahani and
⁎
Banihabib, 2016). Thus, in rivers with a high concentration of sediment, instead of ordinary detention dams, the utilization of open detention dams is proposed (Armanini and Larcher, 2001). This type of dam, which is called slit dam is a sort of open dams having a slit opening in the central part of the dam from the base up to the crest of the dam and pass more sediments through their outlet (Campisano et al., 2014; Wehrmann et al., 2006). The function of slit dam is based on the backwater effect, which causes a portion of sediment deposit upstream of dam (Armanini and Larcher, 2001). The performance of debris flow mitigation measures is investigated in different studies, which can be categories in three groups. The first group of studies uses filed data in order to assess mitigation measures performance (Banihabib and Forghani, 2017; Catella et al., 2005; Hürlimann et al., 2003; Zeng et al., 2009). Difficulties in real-time field measurements of debris flow are the main problem in this method. Therefore, in the second group of studies, it is tried to simulate debris flow events in experimental setups (Banihabib and Bahram, 2009; Lin et al., 2017; Schwindt et al., 2016). Experimental setups could facilitate the possibility of assessment of mitigation measures in different conditions. Complementary to the field and experimental debris flow
Corresponding author. E-mail address: [email protected] (M.E. Banihabib).
https://doi.org/10.1016/j.jher.2019.09.001 Received 9 May 2018; Received in revised form 23 August 2019; Accepted 16 September 2019 Available online 19 September 2019 1570-6443/ © 2019 International Association for Hydro-environment Engineering and Research, Asia Pacific Division. Published by Elsevier B.V. All rights reserved.
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
level, ρt is specific density of debris flow and finally τbx and τby are shear stresses in x and y direction, respectively. The components of shear stresses can be calculated by the following equations:
studies, numerical models are developed (Hassan-Esfahani and Banihabib, 2016; Papa et al., 2018). Numerical models are useful for the analysis of the performance of mitigation measures before construction. Various research has focused on exploring slit dam efficiency. Armanini and Larcher (2001) proposed a theoretical model for designing opening of slit dams. This theoretical method was used to evaluate the efficiency of four slit dams located in Versilia River of Italy by Catella et al. (2005). A non-dimensional equation for the determination of the sedimentation profile in slit dam reservoirs was represented in the study of Banihabib and Bahram (2009). Moreover, Campisano et al. (2014) and Liu et al.,(2012) evaluated slit dam efficiency with one-dimensional numerical models. Although two-dimensional numerical debris flow simulation has been utilized for slot and sectional detention dams (Banihabib and Mokhtari, 2001; Chen et al., 2018), to the best of our knowledge, there is no research for simulating debris flow in slit dam two-dimensionally. The objective of this study is to propose a two-dimensional mathematical model for simulation of debris flow deposition in the slit dam reservoir. In accordance with this objective, the numerical method of Banihabib and Mokhtari (2001) for debris flow simulation in slot dams was upgraded for debris flow simulation in slit dam reservoir. The efficiency of the proposed numerical method is assessed using debris flow data obtained from the experimental setup and other numerical models. Finally, the sensitivity of the numerical model to changes in the selected parameters is investigated.
τbx =
ρt u u2 + v 2 ϕ2
(3)
τby =
ρt v u2 + v 2 ϕ2
(4)
where φ is the dimensionless velocity coefficient, which is the ratio of sectional mean velocity (U) to the shear velocity (u*) (Takahashi, 2007). This coefficient (φ = U/u*) is in the range of 1–10 for debris flows (Banihabib and Mokhtari, 2001). According to the velocity coefficient, flood volume can be calculated in each section based on Eq. (5): m
Vm =
2.1.2. Governing equations of solid phase The continuity equation of the solid phase can be derived by the mass conservation law. Considering all components of this equation, the entire sediment concentration variations in time and space are considered. This equation can be defined as follows:
∂ (ch) ∂Z ∂ (uhc ) ∂ (vhc ) + + =0 + Cs ∂t ∂t ∂x ∂y
Cs − c ∂ (Z ) ∂ (c ) ∂c ∂c + +u +v =0 h ∂t ∂t ∂x ∂y
where
⎡h⎤ U = ⎢uh ⎥, ⎢ ⎣ vh ⎥ ⎦ ⎡ uh ⎤ 2 E = ⎢ βx u h ⎥, ⎢ β uvh ⎥ ⎣ y ⎦ vh ⎡ ⎤ F = ⎢ βx uvh ⎥, ⎢ β v 2h ⎥ ⎣ y ⎦
c=
τ bx ρt τ by ρt
(7)
Different equations have been presented to calculate mean volumetric concentration including those by Banihabib (1997), Hashimoto and Hirano (1995) and Takahashi et al. (1992). In the equation of Takahashi et al. (1992) the static force is considered to estimate the debris flow transport rate. This equation employs the mean diameter of the sediment and the slope of the channel to estimate the sediment concentration. Hashimoto and Hirano (1995) considered shear stress, the physical characteristics of the sediment and also gradient for the estimation of the sediment transport rate. Investigations performed by Banihabib (1997) indicated that Takahashi equation may underestimate the sediment transport rate and Hashimoto and Hirano's formula overestimates the sediment transport rate (Banihabib, 1997). In this paper, the formula presented by Banihabib (1997) is used where the dimensionless sediment transport rate is dependent upon the relative depth and the dimensionless bed shear stress. Banihabib’s equation is defined as follows:
(1)
⎤ ⎥ − gh sin θx ⎥ ⎥ − gh sin θy ⎥ ⎦
(6)
where c is the mean volumetric concentration of sediment over the entire flow depth and Cs is the sediment concentration in the bed. The continuity equation of sediment phase (Eq. (6)) can be rewritten in the form of Eq. (7), considering the continuity equation of the liquid phase (Banihabib et al., 1996).
2.1.1. Governing equations of the liquid phase The shallow water equations, which are derived from the depthintegrated Navier-Stokes equations for an incompressible fluid, are widely used in computational fluid dynamics. Governing equations for hyper-concentrated sediment-laden flow can be acquired by considering the interaction of flow, sediment transport, and bed level variations in shallow water equations (Wu, 2007; Wu et al., 2012). Accordingly, the matrix form of the two-dimensional mass and momentum conservation equations for debris flow can be described as follows (Takahashi, 2007):
∂U ∂E ∂F + + ST = 0 + ∂t ∂x ∂y
(5)
where Ri is the hydraulic radius, S is channel bed slope, A is flow crosssection, i is time counter, and m is a number of time steps that flow depth is measured.
2.1. Governing equations
∂Z
gRi S Ai dt
i=1
2. Material and methods
⎡ ∂t ⎢ ∂ (Z + h) ST = ⎢ gh ∂x cos θx + ⎢ ⎢ gh ∂ (Z + h) cos θy + ∂y ⎣
∑ϕ
qs q
qs sgd3
(8) .65
h = α ⎛ ⎞ τ∗1.5 ⎝d⎠
(9)
where c is the definition of sediment transport concentration, q is the flow discharge per unit width, qs is the sediment discharge per unit width, α is a dimensionless coefficient, d is the sediment diameter, s is the relative submerged unit weight of the sediment and τ∗ is the dimensionless shear stress. The dimensionless shear stress is calculated from Eq. (10):
(2)
where h is the flow depth, u is the flow velocity in the x-direction, v is the flow velocity in the y-direction, g is the acceleration due to gravity, θx, and θy are inclination of the original bed surface, in the x-axis and yaxis directions, βx and βy are momentum correction factors, Z is bed 66
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
τ∗ =
(u2 + v 2) ϕ2gsd
predictor and a corrector. The flow variables at the k time step are known, and their values are calculated in time step k + 1. Thus, for grid points i, j, the backward (Eq. (13)) and forward (Eq. (14)) finite-difference approximation is used for the numerical integration of Eq. (1) in the predictor and corrector stages (Chaudhry, 2007).
(10)
The following equation can be acquired by substituting the shear stress of Eq. (10) into Eq. (9).
qs sgd3
1.5
h .65 u2 + v 2 ⎞ = α ⎛ ⎞ ⎛⎜ 2 ⎟ ⎝ d ⎠ ⎝ ϕ gsd ⎠
Ui∗, j = Uik, j −
(11)
2≤i≤N Δt Δt ∇x Eik, j − ∇y Fik, j − ΔtSik, j ⎧ ⎨ Δx Δy ⎩2 ≤ j ≤ M
(12)
1≥i≥N−1 Δt Δt Δx Ei∗, j − Δy Fi∗, j − ΔtSi∗, j ⎧ ⎨ Δx Δy ⎩1 ≥ j ≥ M − 1
(13)
The potential of sediment transport by debris flow per unit of width (qs) can be calculated by Eq. (11). Therefore, by substituting (qs) in Eq. (8), the mean volumetric concentration can be obtained.
k Ui∗∗ , j = Ui, j −
2.2. Numerical methods
where U* and U** are the midway values of U, Δx and Δy are the size of the spatial grid, in the x-and y-directions and ∇x and Δx are backward and forward differences, which are defined as follows:
One solution method that can be employed to solve the governing equations of this study is the decoupled model. This method assumes that change in one variable during a given time step is small enough that its effect on the other variables can be assumed negligible (Cao et al., 2002; Kassem and Chaudhry, 1998). Therefore, the flow equations and continuity equation of sediment phase can be solved separately. In the decoupled method, flow equations are solved at first, and subsequently, the flow variables obtained from solving flow equations are used to predict the bed level changes. According to above-mentioned method, Eq. (2) was solved by MacCormack scheme assuming bed level change rate is the same as its value obtained in the previous time step. MacCormack scheme was used in the various studies for debris flow simulation (Campisano et al., 2014; He et al., 2016; Ouyang et al., 2015). After calculating flow variables by Eq. (2), the bed level and sediment concentration were estimated using Banihabib et al. (1996) method. The previous steps were repeated iteratively during a computational time step until the changes between two continuous bed levels estimations became less than 0.01 maximum depth of debris flow. A flowchart, which explains the computational sequences of the numerical method, is shown in Fig. 1.
Δx Ui, j = U(i + 1, j) − U(i, j)
(14)
∇x Ui, j = U(i, j) − U(i − 1, j)
(15)
The new values of vector U at time k + 1 are attained by the following equation:
Uik, j+ 1 =
1 ∗ (Ui, j + Ui∗∗ ,j ) 2
(16)
By using a second-order finite-difference scheme like the MacCormack scheme, some numerical oscillations occur at the watersurface discontinuities. For smoothing these oscillations, some schemes are proposed such as Total Variation Diminishing Property (TVDP) schemes (Chakravarthy and Osher, 1985), Monotonic Upstream Scheme For Conservation Laws (MUSCL) (Van Leer, 1974) and compact schemes and Essentially Non-Oscillatory (ENO) schemes (Liu et al., 1994). In this research, to decrease these oscillations, the Jameson artificial viscosity (Jameson et al., 1981) was selected for its simplicity. The artificial viscosity dissipates fluctuations in the zones with high variations, but it does not affect quiescent zones. The amount of artificial viscosity can be regulated by adjusting dissipation constants (Kx and Ky) in x-and y-directions. The details of this procedure can be found in (Chaudhry, 2007). For the numerical stability of the MacCormack scheme, the CourantFriedrichs-Lewy criterion (CFL) should be equal or less than one at all of the simulation nodes (Chaudhry, 2007). Thus, the time step is restricted by the following condition:
2.2.1. Numerical integration of the flow equations The MacCormack scheme is composed of a two-stage chain of a
CFL =
Δt max[(|u| + min[dx , dy]
gh ), (|v| +
gh )] ≤ 1
(17)
In the simulation, CFL number was set to be less than the unit. 2.2.2. Determination of the bed level variation In the reservoir of dams in steep channels, both subcritical and supercritical flows can occur. Banihabib et al. (1996) proposed a scheme that considers both the subcritical and supercritical condition for integrating the continuity equation of the solid phase (Eq. (7)). In this scheme, the forward difference scheme is considered for the supercritical state, and the backward difference scheme is used for the subcritical state. According to this scheme, Eq. (7) is discretized as shown in Eq. (18).
Zik, j+ 1 = Zik, j +
Δt . hik, j+ 1
k+1 − c k i, j
⎡ ci, j ⎣
(cik, j+ 1 − Cs ) ⎢
Δt
+ AFx . uik, j+ 1
(c k + 1 − cik−+1,1j ) k + 1 i, j
+ (1 − AFx ). ui, j
(c k + 1 − c k + 1 ) k + 1 i, j + 1 i, j
AFy . vi, j Fig. 1. A flow chart showing the computational sequence of the two-dimensional numerical model for debris flow in slit dam reservoir.
Δy
Δx
Δx
+
+(1 − AFy ). vik, j+ 1
where AFx and AFy are calculated as: 67
(cik++1,1j − cik, j+ 1 )
(cik, j+ 1 − cik, j+−11 ) Δy
⎤ ⎥ ⎦
(18)
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
Fig. 2. Experimental setup; (a): plan and (b): lateral view of the flume.
⎧ 0 if AFx = ⎨ 1 if ⎩
u Fr 2 − 1 u Fr 2 − 1
≤0 >0
⎧ 0 if AFy = ⎨ 1 if ⎩
v Fr 2 − 1 v Fr 2 − 1
bed level, and p is a coefficient that considers the sediment distribution with depth (in percent). Considering that at the start point of each time step, cik, j+ 1 was not calculated for computation of the bed level variation (Eq. (18)), sediment concentrations of the last two time steps, k and k − 1, were substituted for those of time steps k and k + 1. Then the sediment concentration for time level k + 1 was computed using Eq. (20) using achieved bed level variation. The bed level variation which obtained in the previous step was used to adjust flow property in the current time step. Adjustment of flow property, bed level variation, and sediment concentrations for a given time step was repeated until bed level difference between two consecutive iterations became less than 0.01 maximum depth of debris flow (Fig. 1).
≤0 >0
(19)
Eq. (19) is an implicit equation of bed level Z and volumetric concentration of sediment c. A control volume equation was used for the computation the volumetric concentration of the sediment at different grid points (Banihabib and Mokhtari, 2001). According to the control volume, volumetric concentration of sediment at time step k + 1 was calculated by Eq. (20).
cik, j+ 1 =
k k (∑ qs+c.i, j dx . dy. hi, j ) − qz dx . dy. hik, j+ 1
(20)
Eq. (20) is taken from continuity equation, where ∑ qs is the total amount of sediment, entering or leaving the element at k + 1, c.ik, j dx . dy. hik, j is available sediment in the element at k, qz is the volume of sediment deposited or eroded in the element at k + 1, and dx . dy. hik, j+ 1 is the volume of the element at k + 1. To compute the net sediment input volume to the element due to the concentration gradient, Eq. (21) was used.
∑ qs = (qs (in) )x − (qs (out ) )x + (qs (in) ) y − (qs (out ) ) y
2.2.3. Initial and boundary conditions Debris flow simulation in slit dam reservoir is like a dam break simulated by Fennema and Chaudhry (1989) which employed the MacCormack Scheme for numerical modeling. Fennema and Chaudhry (1989) demonstrated that the MacCormack Scheme fails when the ratio of tail-water depth to reservoir depth being less than 0.25. Here, flow depth in the initial condition was assumed 0.25 of inflow, and the velocity was determined by Eq. (27)
(21)
where:
(qs (in) )x = (qs (out ) )x =
(qs (in) ) y = (qs (out ) ) y =
1 k (ci, j + cik− 1, j )(uik, j + uik− 1, j )(hik, j + hik− 1, j )Δy. Δt 8 1 k (ci, j + cik+ 1, j )(uik, j + uik+ 1, j )(hik, j + hik+ 1, j )Δy. Δt 8
1 k (ci, j + cik, j − 1 )(vik, j + vik, j − 1 )(hik, j + hik, j − 1 )Δx . Δt 8 1 k (ci, j + cik, j + 1 )(vik, j + vik, j + 1 )(hik, j + hik, j + 1 )Δx . Δt 8
u (i,1) = ϕ gRS
(22)
where φ is the dimensionless velocity coefficient and R is the hydraulic radius. Solid walls boundary was defined according to the reflection boundary condition. The details of reflection boundary for MacCormack scheme are given in Chaudhry (2007). The boundary at the slit of the dam depends on the ratio of the slit width to the channel width. Due to sudden narrowing of the channel by the slit dam, the flow energy is dispersed by creating vortexes upstream of the slit, on both sides. Armanini and Larcher (Armanini and Larcher, 2001) presented an equation to calculate the energy loss, as shown in Eq. (28).
(23) (24) (25)
The sediment fall velocities can be utilized to calculate the sedimentation volume in an element at k + 1, as shown in Eq. (26):
qz = (Zik. j+ 1 − Z ). dx . dy. Cs.
ws. Δt 1 . p hik, j
(27)
(26)
2 Fr 2 ΔE 2 = u ⎡1 − (Fru Rr )−2/3⎤ hu 2 ⎣ 3 ⎦
where ws is the sediment fall velocity, Cs is sediment concentration in 68
(28)
Journal of Hydro-environment Research 27 (2019) 65–74
0.6 2.6
0.9
Sediment concentration in bed (Cs) (-)
ghd
(29) (30)
The accuracy of the numerical model in this study was tested according to a laboratory experiment and alternative numerical models. The laboratory experiment was designed for simulating debris flow in slit dam reservoir, then the debris flow simulated by numerical method of this study. The other numerical models utilized for a comparison were Kanako model (Nakatani et al., 2008) and a model presented by Campisano et al. (2014). Kanako model applied by Liu et al. (2012) to analyze the effects of slit dam on discharge reduction. Campisano et al. (2014) proposed a numerical method for evaluating the slit dam trapping efficiency. The data obtained from numerical method of this study were compared with other models outputs using statistical indices, correlation coefficient (R2) and Root Mean Square Error (RMSE).
1
2.3.1. Experimental setup The experimental setup was composed of four major parts, including a tilting flume, sediment tank, a slit dam and ultrasonic water level gauges (Fig. 2). The slit dam was made of plexiglass and placed 140 cm from the flume end. During the experiment, three ultrasonic water level gauges were recording flow depth with a precision of 0.1 mm. Fig. 2 shows the plan and lateral views of the flume in this experiment. The flume slope was 43% and its height was 25 cm. The width of the slit was 5 cm, and the sediment was uniform sand with sediments mean size (D50) equals 1 mm and density 2.65 gr/cm3. The first, second and third ultrasonic water level gauges were at distances of 410, 210, and 60 cm from the downstream end of the flume, respectively. The first and second water level gauges were located upstream of the dam, and the third one was downstream of the dam reservoir. The main characteristics of the experiment are displayed in Table 1. In the experiments, first, sediment with a known weight was placed on the sediment tank (Fig. 2). Then, the sediment was slowly saturated with water until a little water covers sediment with the horizontal free surface. This volume of water and sediment formed the debris flood volume in the experiment. The upstream of the flume moved upward to a specific height until the required slope was achieved. Then, the gate of the sediment tank was opened suddenly, the flood moved toward the slit dam, and some of the sediment was deposited in the dam reservoir (Banihabib and Bahram, 2009). To have a precise measurement of the sediment thickness deposited in the dam reservoir, a digital point gauge was used. This point gauge was moved vertically along the flume and measured three points in each section (aa', bb' and cc'. sections Fig. 3). These three points were located in the middle, right and left of the section. The right and the left points were 9 cm away from the middle of the section. Fig. 3 shows the plan and lateral views of debris flow deposition in slit dam reservoir schematically. Debris flow occurred in experiment was fully developed according to the observations and criteria proposed by Takahashi (1991). Debris flow is fully developed when following ratio (Eq. (31)) is more than 0.15 (Takahashi 1991).
25
5
0.085 30 43
5
Sediments mean size D50 (mm) Sediment weight (kg) Flood volume (m3) Flume width (cm) Flume bed slope (%)
Slit width (cm)
2⎡ u 2 hu + u − ΔE⎤ ⎥ 3⎢ 2g ⎣ ⎦
2.3. Model test
Flood characteristic
Dam height (cm)
hd =
ud =
Dam characteristics
Table 1 Main characteristics of the experimental setup.
where ΔE is the energy loss, hu is the flow depth, Fru is Froude number at node u located upstream of the slit, and Rr is the ratio of the channel width to the slit width. Since critical velocity occurs at the slit, the flow depth in the slit is equal to two-thirds of the specific energy at the slit. Here, the energy at the slit was obtained from the difference of the energy upstream and the energy loss determined by Eq. (28). Thus, the boundary condition at the slit was determined as:
Sediment density (g/ cm3)
Sediment distribution in depth (p) (%)
M.E. Banihabib and F. Nazarieh
69
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
Fig. 3. Sketch of debris flow sedimentation in slit dam reservoir, (a) plan and (b) lateral views. 0.5
0.007
Discharge Volumetric concentration
Discharge (m^2/s)
0.005
0.4
0.3
0.004
0.003
0.2
0.002 0.1
Volumetric concentration (percent)
0.006
0.001
0.0
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 2.2 2.4 2.6 2.8 3.0 3.2 3.4 3.6 3.8 4.0 4.2 4.4 4.6 4.8 5.0 5.2 5.4 5.6 5.8 6.0
0.000
Time (s) Fig. 4. Flow and sediment concentration at second ultrasonic water level gauge used as the upstream boundary condition.
sin(S ) > 0.15 (ρs / ρ) − 1
2.3.2. Numerical test The Kanako model was used by Liu et al. (2012) for assessing the performance of slit dams in Tanjutani Gully, Kyoto City, Japan. The average width of this Gullys’ was 6 m, its slope at the investigated section was about 14.7%. and sediment mean size (D50) was 300 mm. Tested dams had three widths of slit including 0.4, 0.8 and 1.2 m. The upstream boundary conditions assume to be outflow hydrograph from the previous check dam. In the current study, slit dam examples in Tanjutani Gully, was simulated by the proposed numerical model considering 0.1 × 0.1 m gride size mesh and 0.004 s time step. Simulated outflow hydrograph downstream of slit dam were compared with the digitized data obtained from Kanako model in Liu et al. (2012) research. Campisano et al. (2014) proposed a numerical method for exploring slit dams trapping efficiency. The case study of Campisano et al. (2014) was characterized with a river with a rectangular section width of 10 m, length of 300 m, initial slope of 5% and sediment mean size (D50)
(31)
where S is channel bed slpe, ρs is the sediment material density, and ρ is water density. In the experiments of this study, the above mentioned ratio was more than 0.25 at downstream of the second ultrasonic water level gauge. The upstream boundary condition of the reservoir was defined according to the flow depths recorded by the second ultrasonic water level gauge. Inflow velocity and sediment concentration at upstream were calculated by Eqs. (27) and (8). The flow and sediment concentration at second ultrasonic water level gauge were determined, as shown in Fig. 4. The explained experiment was simulated by the proposed numerical method considering 0.01 m × 0.01 m grid size and time step 0.001 s.
70
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh 0.0030
Discharge (m^2/s)
0.0025
0.0020
0.0015
0.0010
Experimental data This study model
0.0005
0.0000 0
5
10
15
20
25
30
35
Time (s)
Fig. 5. Comparison of outflow hydrograph simulated by current study numerical model and the experimental data. 0.06
(a)
supposed to be 0.005 s.
Flume bed Experimental data This study model
0.05
3. Results and discussion
Elevation (m)
0.04
3.1. Validation
0.03
3.1.1. Validation based on experimental data Four parameters φ, α, Kx and Ky in the numerical model were calibrated based on the experimental data. Firstly, the dimensionless velocity coefficient (φ) was calibrated considering Eq. (5). The flood volume and flow depth in ultrasonic water level gauges locations are known, therefore, velocity coefficient (φ) is only unknown in Eq. (5) and could be calculated. This coefficient was calibrated as 5.6 for the location of second ultrasonic water level gauge. Secondly, the dimensionless coefficient α was calibrated to the amount of 1.5 by equating the transported sediment volume in both the experimental and mathematical model. Finally, for dissipation constants, Kx and Ky values of 0.49 and 0.3 were obtained by trial and error. To investigate the capability of the proposed numerical method to simulate debris flow passing the slit dam reservoir, simulated outflow hydrograph was compared with the experimental hydrograph. Fig. 5 shows comparison of hydrographs in third ultrasonic water level gauge location. This figure illustrates that the simulated flood peak discharge was estimated earlier and greater than the observed one. The difference between simulated flood peak and experimental flood peak may happen as a result of the larger time step used in the experimental data acquisition (0.2 s) in comparison to the fine time step applied in numerical method (0.001 s). This may cause losing real peak flow in the experimental data. It is inferred from Fig. 5 that the numerical model simulates the general trend of the outflow hydrograph correctly as values of statistical criteria, R2 was 0.93 and RMSE was 0.000182 m2/s. Numerical model capability in simulating sedimentation was assessed in aa', bb' sections which are represented in Fig. 2a. Comparisons of the simulated bed profile and measured data along these sections at the end of the experiment are shown in Fig. 6. Evaluation criteria at section aa' showed an acceptable bed profile simulation (R2 = 0.98 and RMSE = 0.0037 m). However, slight sedimentation was predicted at section aa' by the numerical model while sedimentation was not observed on experimental setup (Fig. 6a). The reason may be due to the solution method which could not fully address continuity of water and sediment simultaneously. Hence, the error of flow properties affects the simulation of sedimentation along section aa'. In Fig. 6(b), the bed profile of the numerical model and experimental data are compared along section bb'. The numerical model simulated the location of the sedimentation starting point and the bed profile accurately with R2 = 0.98 and RMSE = 0.0012 m. Hydraulic jump has a significant effect on deposition in slit dam
0.02
0.01
0.00 0.00
0.02
0.04
0.06
0.08
0.10
Distance from dam (m) 0.06
(b)
Flume bed Experimental data This study model
0.05
Elevation (m)
0.04
0.03
0.02
0.01
0.00 0.00
0.02
0.04
0.06
0.08
0.10
Distance from dam (m)
Fig. 6. Comparison of simulated bed profile and measured experimental data along (a) section aa' and (b) section bb'.
50 mm. It was assumed that slit dam was positioned at 250 m from the upstream of the river. Various contraction ratios of the slit width to the channel width were investigated. The slit widths used for comparison include 2, 3 and 4 m. Inflow was defined by triangular hydrograph with flood peak of 60 m3/s, time to peak of 3600 s and the total duration of 7200 s. Trapping efficiency E (t ) was calculated according to the following equation:
E (t ) =
Vdep (t ) Vsed (t )
(32)
where Vdep (t ) is the volume of sedimentation deposited in slit dam reservoir at time t and Vsed (t ) is volume of sediment entered to the reservoir until time t. Variation of trapping efficiency by time, depicted by Campisano et al. (2014), was digitized to compare with the result of the numerical method presented in the current study. The grid size used for the proposed numerical method was 0.5 × 0.5 m and time step was 71
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh 0.007
Discharge (m^2/s)
0.006
0.25
0.005
0.20
0.004
0.15
0.003
0.10
0.002
0.05
0.001
0.00
Distance of hydraulic jump from dam (m)
0.30
Inflow hydrograph Outflow hydrograph Hydraulic jump distance from dam
0.000 0
2
4
6
8
Time (s) Fig. 7. Variations of the simulated hydraulic jump location with time along section bb'. 50
50
(a)
Inflow hydrograph Outflow hydrograph (This study) Outflow hydrograph (Kanako)
30
20
10
0
Inflow hydrograph Outflow hydrograph (This study ) Outflow hydrograph (Kanako)
40
Discharge (m^3/s)
Discharge (m^3/s)
40
(b)
30
20
10
0
200
400
600
800
1000
1200
1400
1600
0
1800
Time (s)
200
400
600
800
1000
1200
1400
1600
1800
Time (s)
50
(c)
Inflow hydrograph Outflow hydrograph (This study) Outflow hydrograph (Kanako)
40
Discharge (m^3/s)
0
30
20
10
0
0
200
400
600
800
1000
1200
1400
1600
1800
Time (s) Fig. 8. Comparison of outflow hydrograph simulated by this study numerical model and Kanako model for the slit width of (a) 0.4 m, (b) 0.8 m and (c) 1.2 m.
reservoir, in fact by flow transition from supercritical to subcritical, flow transport capacity decrease and it causes a part of transported sediment deposit in the reservoir (Armanini and Larcher, 2001). The variation of hydraulic jump location along section bb' is shown versus time in Fig. 7. It should be noted that the hydraulic jump displacement along section bb' was the same as section aa'. The figure shows that, as the flood was filling the reservoir, the hydraulic jump moved upstream, and this process continued until most of the dam reservoir was emptied from the flood. Then, the hydraulic jump location moved downstream.
3.1.2. Validation based on other numerical models The outflow hydrographs generated by Kanako model in Liu et al. (2012) study and outflow hydrographs simulated with the numerical model of this study for the same example are compared in Fig. 8. These outflow hydrographs belong to the slit width of 0.4 m (Fig. 8a), 0.8 m (Fig. 8b) and 1.2 m (Fig. 8c). Calibrated amounts of φ, α, Kx and Ky parameters in this example were 3.7, 0.9, 0.48 and 0.35 respectively. The statistical criteria R2 and RMSE for the width of 0.4 m were 0.93 and 4.3, for the width of 0.8 were 0.98 and 1.8 and for the width of 72
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh 1.0
Trapping Efficiency (-)
0.8
(2014) numerical models. The statistical criteria, R2 and RMSE, for slit width of 2 m were 0.94 and 0.34, for slit width of 3 m were 0.98 and 0.045 and for slit width of 3 m were 0.99 and 0.038. It can be seen in Fig. 8 that Campisano et al. (2014) model underestimate the trapping efficiency during the decreasing branch of the hydrograph, which dam self-cleaning process happen in this period. It may be as a result of onedimensional simulation of flow behind slit dam in Campisano et al. (2014) model which can not simulate sediment trasport proccess in the corners of the reservoir properly. Considering statistical criteria in all evaluation tests, it can be concluded that there is a good agreement between the result of the proposed numerical model and other numerical models. Furthermore, calibrated parameters (φ, α, r Kx Ky) had restricted a range of variation in all tests. These calibrated values are also located in the range of acceptable values presented in other researches such as Chaudhry (2007) and Takahashi (2007), which can give some support to the credibility of the proposed numerical model.
d=2 m
0.6
d=3 m
0.4
d=4 m
0.2
Campisano Model This study Model
0.0 0
20
40
60
80
100
120
Time (min)
Fig. 9. Time-variation of trapping efficiency at slit dam simulated by this study and Campisano et al. (2014) numerical models at slit width of 2 m, 3 m and 4 m.
Changes in total outflow volume (%)
30
3.2. Sensitivity analysis
(a)
Four input parameters were used in the sensitivity analysis, including velocity coefficient (φ), dimensionless coefficient (α), and coefficients of dispersion Kx and Ky. Each input parameter was modified by ± 5, ± 10, ± 15 and ± 20 percent from its original value (calibrated values based on experimental data), while the other input parameters were kept constant at their calibration values. The numerical method was not stable in some cases, so some sensitivity curves are not complete for the entire range of parameter variations as shown in Fig. 10. Two different simulated outputs were assessed in the sensitivity analysis: the total outflow volume (Fig. 10a) and the total deposition volume trapped in slit dam reservoir (Fig. 10b). Sedimentation has a higher rate of variation (Fig. 10b) in comparison with outflow volume (Fig. 10a). It demonstrates that sedimentation volume was more sensitive to the parameter changes. Generally, the response of outflow and sedimentation volume to changes in parameters was reverse. In other words, when variation in parameters leads to increase in outflow volume, sedimentation volume decreased and vice versa. Sensitivity analysis in Fig. 10 also illustrates that raising in velocity coefficient (φ) leads to increase of outflow volume and decrease of sedimentation while raising in parameter (α) causes more sedimentation and less outflow volume. Changes in outflow and deposition volumes to variation in Kx and Ky seem to be similar. The curves related to these two parameters are so that any changes in them (regardless of increase or decrease of Kx and Ky), result in increase of outflow volume and decrease the volume of trapped sedimentation (Fig. 10).
20
10
0
-10 Velocity coefficient Dimensionless coefficient Coefficients of dispersion Kx
-20
Coefficients of dispersion Ky
-30
-20%
-15%
-10%
-5%
0%
5%
10%
15%
20%
Parameter changes (%) 40
Changes in total deposition volume (%)
(b) 30 Velocity coefficient Dimensionless coefficient Coefficients of dispersion Kx
20
10
Coefficients of dispersion Ky
0
-10
4. Conclusion -20
The MacCormack and Banihabib schemes were combined to simulate debris flow sedimentation in a slit detention dam reservoir. Performance of the proposed numerical methods was evaluated by comparing simulated data of this model with other experimental and numerical models, and the final conclusions can be summarized as follows:
-30 -20%
-15%
-10%
-5%
0%
5%
10%
15%
20%
Parameter changes (%)
Fig. 10. The sensitivity of (a) outflow and (b) sedimentation volume to the numerical model parameters variations.
• The MacCormack scheme and banihabib scheme could simulate the
1.2 m were 0.99 and 0.54. Furthermore, Fig. 8 illustrates that by increasing the width of slit, agreement between two models modifies. The time-variation of the trapping efficiency of the slit dam computed by the numerical model of this study, was compared with the result obtained from Campisano model (Campisano et al., 2014) in Fig. 9. The calibrated parameters in this case were 5.6 for φ, 1.3 for α, 0.42 for Kx and 0.3 for Ky. The results presented in Fig. 8 illustrates a generally good agreement between this study and Campisano et al.
• • 73
debris flow regime variation from supercritical to subcritical and the movement of the hydraulic jump in the reservoir of the slit dam. According to the comparison of outflow hydrographs in investigated models, the slit dam boundary condition was capable of simulating the general trend of the outflow hydrograph. The simulation of the bed level changes was fairly appropriate along the side longitudinal section, and it was acceptable along the central longitudinal section of the experimental setup. The existence of error in bed level simulation at the central section may be due to the
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
•
solution method, which could not fully address continuity equations. Sensitivity analysis revealed that simulated sedimentation is more sensitive to the parameter variation than outflow hydrograph. Generally, outflow and sedimentation showed an opposite response to the parameters changes.
channel, 50th annual conference of JSCE Japan, pp. 495–500. Hassan-Esfahani, L., Banihabib, M.E., 2016. The impact of slit and detention dams on debris flow control using GSTARS 3.0. Environ. Earth Sci. 75, 328. https://doi.org/ 10.1007/s12665-015-5183-z. He, S., Ouyang, C., Liu, W., Wang, D., 2016. Coupled model of two-phase debris flow, sediment transport and morphological evolution. Acta Geol. Sin. Engl. 90, 2206–2215. https://doi.org/10.1111/1755-6724.13031. Hunt, B., 1994. Newtonian fluid mechanics treatment of debris flows and avalanches. J. Hydraul. Eng. 120, 1350–1363. https://doi.org/10.1061/(ASCE)0733-9429(1994) 120:12(1350). Hürlimann, M., Rickenmann, D., Graf, C., 2003. Field and monitoring data of debris-flow events in the Swiss Alps. Can. Geotech. 40 (1), 161–175. https://doi.org/10.1139/ t02-087. Jameson, A., Schmidt, W., Turkel, E., 1981. Numerical solution of the Euler equations by finite volume methods using Runge Kutta time stepping schemes, 14th fluid and plasma dynamics conference, Palo Alto, California, pp. 1259. Kassem, A.A., Chaudhry, M.H., 1998. Comparison of coupled and semicoupled numerical models for alluvial channels. J. Hydraul. Eng. 128 (4), 794–802. https://doi.org/10. 1061/(ASCE)0733-9429(1998)124:8(794). Lin, X., Huo, M., Zhou wen, J., Cao, T., Yang Rong, F., Zhou wei, H., 2017. An experimental study on controlling post-earthquake debris flows using slit dams. Environ. Earth. Sci. 76 (22), 780. https://doi.org/10.1007/s12665-017-7140-5. Liu, J., Nakatani, K., Mizuyama, T., 2012. Hazard mitigation planning for debris flow based on numerical simulation using Kanako simulator. J. Mt. Sci. 9, 529–537. https://doi.org/10.1007/s11629-012-2225-9. Liu, X.-D., Osher, S., Chan, T., 1994. Weighted essentially non-oscillatory schemes. J. Comput. Phys. 115, 200–212. https://doi.org/10.1006/jcph.1994.1187. Michelini, T., Bettella, F., D’Agostino, V., 2017. Field investigations of the interaction between debris flows and forest vegetation in two Alpine fans. Geomorphology 279, 150–164. https://doi.org/10.1016/j.geomorph.2016.09.029. Nakatani, K., Wada, T., Satofuka, Y., Mizuyama, T., 2008. Development of “Kanako 2D (Ver. 2.00)”, a user-friendly one- and two-dimensional debris flow simulator equipped with a graphical user interface. Int. J. Eros. Control. Eng. 1, 10. Ouyang, C., He, S., Tang, C., 2015. Numerical analysis of dynamics of debris flow over erodible beds in Wenchuan earthquake-induced area. Eng. Geol. 194, 62–72. https:// doi.org/10.1016/j.enggeo.2014.07.012. Paik, J., 2015. A high resolution finite volume model for 1D debris flow. J. Hydro-environ. Res. 9, 145–155. https://doi.org/10.1016/j.jher.2014.03.001. Pan, H., Wang, R., Huang, J., Ou, G., 2013. Study on the ultimate depth of scour pit downstream of debris flow sabo dam based on the energy method. Eng. Geol. 160, 103–109. https://doi.org/10.1016/j.enggeo.2013.03.026. Papa, M., Sarno, L., Vitiello, F., Medina, V., 2018. Application of the 2D depth-averaged model, FLATModel, to pumiceous debris flows in the Amalfi Coast. Water 10 (9), 1159. https://doi.org/10.3390/w10091159. Pirulli, M., Pastor, M., 2012. Numerical study on the entrainment of bed material into rapid landslides. Geotechnique 62, 959–972. https://doi.org/10.1680/geot.10.P.074. Qiu, F., Huang, J., Li, Y., Han, Z., Wang, W., Chen, G., Qu, X., Su, B., 2018. Protecting highway bridges against debris flows using lateral berms: a case study of the 2008 and 2011 Cheyang debris flow events, China. Geomat. Nat. Haz. Risk 9, 196–210. https://doi.org/10.1080/19475705.2017.1414718. Rickenmann, D., Laigle, D., McArdell, B.W., Hübl, J., 2006. Comparison of 2D debris-flow simulation models with field events. Comput. Geosci. 10, 241–264. https://doi.org/ 10.1007/s10596-005-9021-3. English. Schwindt, S., De Cesare, G., Boillat, J.-L., Schleiss, A., 2016. Physical modelling optimization of a filter check dam in Switzerland. Conference Proceedings, INTERPRAEVENT 2016, Hazard and Risk Mitigation 828–836. Takahashi, T., 1991. Debris flow. In: IAHR Monograph Serie. Balkema, Rotterdam, pp. 1–165. Takahashi, T., 2007. Debris Flow: Mechanics, Prediction and Countermeasures, second ed. CRC Press, United States. Takahashi, T., Nakagawa, H., Harada, T., Yamashiki, Y., 1992. Routing debris flows with particle segregation. J. Hydraul. Eng. 118, 1490–1507. https://doi.org/10.1061/ (ASCE)0733-9429(1992) 118:11(1490). Tsai, Y., 2007. A debris-flow simulation model for the evaluation of protection structures. J. Mt. Sci. 4, 193–202. https://doi.org/10.1007/s11629-007-0193-2. English. Van Leer, B., 1974. Towards the ultimate conservative difference scheme. II. Monotonicity and conservation combined in a second-order scheme. J. Comput. Phys. 14, 361–370. https://doi.org/10.1016/0021-9991(74)90019-9. VanDine, D., 1996. Debris flow control structures for forest engineering, Research Branch, Ministry of Forests Research Program, Province of British Columbia. Wehrmann, H., Hübl, J., Holzinger, G., 2006. Classification of dams in torrential watersheds, Disaster Mitigation of Debris Flows, Slope Failures and Landslides, Universal Academy Press, Tokyo, Japan., pp. 829–838. Wu, W., 2007. Computational River Dynamics, first ed. CRC Press, United States. Wu, W., Marsooli, R., He, Z., 2012. Depth-averaged two-dimensional model of unsteady flow and sediment transport due to noncohesive embankment break/breaching. J. Hydraul. Eng. 138, 503–516. https://doi.org/10.1061/(ASCE)HY.1943-7900. 0000546. Zeng, Q.L., Yue, Z.Q., Yang, Z.F., Zhang, X.J., 2009. A case study of long-term field performance of check-dams in mitigation of soil erosion in Jiangjia stream, China. Eng. Geol. 58 (4), 897–911. https://doi.org/10.1007/s00254-008-1570-z.
Generally, the numerical method presented in this article can provide some insights on debris flow sedimentation upstream of the slit dams. However, because of limited availability of measured debris flow data, the model in this study was just compared with the result of an experimental data and one-dimensional numerical models. Further comparison with field data and other numerical models is required to enhance the reliability of the established model. Acknowledgments The authors appreciate Dr. Mary Kay Camarillo of the University of the Pacific and Dr. Seyyed-Mehdi Hashemi-Shahedani of University of Tehran for editing and reviewing this paper. Moreover, we thank Iran’s Ministry of Energy, Water, and Wastewater Standards and Projects Bureau for supporting this research through contract No. 720/1015/87. Appendix A. Supplementary data Supplementary data to this article can be found online at https:// doi.org/10.1016/j.jher.2019.09.001. References Armanini, A., Larcher, M., 2001. Rational criterion for designing opening of slit-check dam. J. Hydraul. Eng. 127, 94–104. https://doi.org/10.1061/(ASCE)07339429(2001) 127:2(94). Banihabib, M.E., 1997. Sediment transport of mud flow. In: 1st Iranian Hydraulic Conference. Iranian Hydraulic Association, Tehran, Iran, pp. 1–10. Banihabib, M.E., Bahram, E., 2009. Experimental Analyses of Sedimentation in the Slit Dam Reservoir. Great Rivers, Kansas City, Missouri, USA, World Environmental and Water Resources Congress, pp. 1–12. Banihabib, M.E., Forghani, A., 2017. An assessment framework for the mitigation effects of check dams on debris flow. CATENA 152, 277–284. https://doi.org/10.1016/j. catena.2017.01.018. Banihabib, M.E., Hirano, M., Harada, T., 1996. A mathematical model for hazard zone mapping of debris flow. Geo-Inf. 7, 87–90. https://doi.org/10.6010/ geoinformatics1990.7.1-2_87. Banihabib, M.E., Mokhtari, A.H., 2001. Numerical simulation in detention dams during a high concentrated flow. In: 3rd International Symposium on Environmental Hydraulics. Arizona, USA, pp. 1–6. Campisano, A., Cutore, P., Modica, C., 2014. Improving the evaluation of slit-check dam trapping efficiency by using a 1D unsteady flow numerical model. J. Hydraul. Eng. 140, 04014024. https://doi.org/10.1061/(ASCE)HY.1943-7900.0000868. Cao, Z., Day, R., Egashira, S., 2002. Coupled and decoupled numerical modeling of flow and morphological evolution in alluvial rivers. J. Hydraul. Eng. 128 (3), 306–321. https://doi.org/10.1061/(ASCE)0733-9429(2002)128:3(306). Catella, M., Paris, E., Solari, L., 2005. Case study: efficiency of slit-check dams in the mountain region of versilia basin. J. Hydraul. Eng. 131 (3), 145–152. https://doi.org/ 10.1061/(ASCE)0733-9429(2005)131:3(145). Chakravarthy, S., Osher, S., 1985. A new class of high accuracy TVD schemes for hyperbolic conservationlaws. In: 23rd Aerospace Sciences Meeting, Reno, Nevada, pp. 363. Chaudhry, M.H., 2007. Open-Channel Flow, Second ed. Springer, US. Chen, S.-C., Tfwala, S., 2018. Evaluating an optimum slit check dam design by using a 2D unsteady numerical model, E3S Web of Conferences. Fennema, R.J., Chaudhry, M.H., 1989. Implicit methods for two-dimensional unsteady free-surface flows. J. Hydraul. Res. 27, 321–332. https://doi.org/10.1080/ 00221688909499167. Ghilardi, P., Natale, L., Savi, F., 2001. Modeling debris flow propagation and deposition. Phys. Chem. Earth Part C 26, 651–656. https://doi.org/10.1016/S1464-1917(01) 00063-0. Han, Z., Chen, G., Li, Y., Tang, C., Xu, L., He, Y., Huang, X., Wang, W., 2015. Numerical simulation of debris-flow behavior incorporating a dynamic method for estimating the entrainment. Eng. Geol. 190, 52–64. https://doi.org/10.1016/j.enggeo.2015.02. 009. Hashimoto, H., Hirano, M., 1995. Sediment transport rate and resistance on a steep
74
Contents lists available at ScienceDirect
Journal of Hydro-environment Research journal homepage: www.elsevier.com/locate/jher
A model for simulation of debris flow sedimentation in slit detention-dam reservoirs
T
⁎
Mohammad Ebrahim Banihabib , Farzaneh Nazarieh Department of Irrigation and Drainage Engineering, College of Aburaihan, University of Tehran, Pakdasht 20 km, Tehran, Iran
A R T I C LE I N FO
A B S T R A C T
Keywords: Slit dam Sedimentation Two-dimensional Debris flow
Debris flow is a natural phenomenon that requires appropriate structural measures such as detention dams to reduce hazards resulting from it. The objective of this paper is to present a two-dimensional numerical model to investigate debris flow deposition in the reservoir of the detention slit dams. The governing equations include shallow water equations considering the bed level variation, sediment continuity, and debris-flow transport equations. Shallow water equations were discretized using the MacCormack scheme and then Banihabib Scheme was applied for solving sediment continuity equation. The governing equations were solved in an iterative procedure and as a result, flow properties and bed level variations were obtained in silt dam reservoir. The model performance was evaluated against experimental and numerical models in two aspects including the outflow hydrograph and reservoir bed level variations or trapping efficiency. The model performed quite well in outflow hydrograph simulation, and the comparison of the bed level variations or trapping efficiency illustrated that the simulated values were acceptable. A sensitivity analysis of the parameters was carried out to analyze the influence of parameters variation on outflow and sedimentation volume. The results of this sensitivity analysis revealed sedimentation is more sensitive to the parameter variation than outflow hydrograph. Generally, the proposed numerical method can provide insights on the debris flow deposition in slit dam reservoir, which could be useful for slit dam design.
1. Introduction Debris flow leads to infrastructure destruction and loss of lives. It mostly occurs in mountainous areas where unstable sediments exist on steep slopes. Debris flows have received specific attention all over the world such as Japan (Han et al., 2015; Nakatani et al., 2008; Takahashi, 2007), China (Pan et al., 2013; Qiu et al., 2018), Taiwan (Tsai, 2007), the European Alps (Ghilardi et al., 2001; Michelini et al., 2017; Rickenmann et al., 2006) and other part of the world (Hunt, 1994; Paik, 2015; Pirulli and Pastor, 2012). Given the probable occurrence of debris flows and their potential destructiveness, they should be controlled by proper measures. There are different measures for flood control, including flood diversion channels, ring levees, river floodplain management, detention dams and etc (VanDine, 1996). Detention dam, which is a structural measure built on small or moderate rivers, can be quite satisfactory for flood control in short term. However, in the cases that debris and sediments fill the reservoir of the detention dam, the capability of the dam to alleviate floods may be decreased, and dam performance is challenged (Banihabib and Bahram, 2009; Hassan-Esfahani and
⁎
Banihabib, 2016). Thus, in rivers with a high concentration of sediment, instead of ordinary detention dams, the utilization of open detention dams is proposed (Armanini and Larcher, 2001). This type of dam, which is called slit dam is a sort of open dams having a slit opening in the central part of the dam from the base up to the crest of the dam and pass more sediments through their outlet (Campisano et al., 2014; Wehrmann et al., 2006). The function of slit dam is based on the backwater effect, which causes a portion of sediment deposit upstream of dam (Armanini and Larcher, 2001). The performance of debris flow mitigation measures is investigated in different studies, which can be categories in three groups. The first group of studies uses filed data in order to assess mitigation measures performance (Banihabib and Forghani, 2017; Catella et al., 2005; Hürlimann et al., 2003; Zeng et al., 2009). Difficulties in real-time field measurements of debris flow are the main problem in this method. Therefore, in the second group of studies, it is tried to simulate debris flow events in experimental setups (Banihabib and Bahram, 2009; Lin et al., 2017; Schwindt et al., 2016). Experimental setups could facilitate the possibility of assessment of mitigation measures in different conditions. Complementary to the field and experimental debris flow
Corresponding author. E-mail address: [email protected] (M.E. Banihabib).
https://doi.org/10.1016/j.jher.2019.09.001 Received 9 May 2018; Received in revised form 23 August 2019; Accepted 16 September 2019 Available online 19 September 2019 1570-6443/ © 2019 International Association for Hydro-environment Engineering and Research, Asia Pacific Division. Published by Elsevier B.V. All rights reserved.
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
level, ρt is specific density of debris flow and finally τbx and τby are shear stresses in x and y direction, respectively. The components of shear stresses can be calculated by the following equations:
studies, numerical models are developed (Hassan-Esfahani and Banihabib, 2016; Papa et al., 2018). Numerical models are useful for the analysis of the performance of mitigation measures before construction. Various research has focused on exploring slit dam efficiency. Armanini and Larcher (2001) proposed a theoretical model for designing opening of slit dams. This theoretical method was used to evaluate the efficiency of four slit dams located in Versilia River of Italy by Catella et al. (2005). A non-dimensional equation for the determination of the sedimentation profile in slit dam reservoirs was represented in the study of Banihabib and Bahram (2009). Moreover, Campisano et al. (2014) and Liu et al.,(2012) evaluated slit dam efficiency with one-dimensional numerical models. Although two-dimensional numerical debris flow simulation has been utilized for slot and sectional detention dams (Banihabib and Mokhtari, 2001; Chen et al., 2018), to the best of our knowledge, there is no research for simulating debris flow in slit dam two-dimensionally. The objective of this study is to propose a two-dimensional mathematical model for simulation of debris flow deposition in the slit dam reservoir. In accordance with this objective, the numerical method of Banihabib and Mokhtari (2001) for debris flow simulation in slot dams was upgraded for debris flow simulation in slit dam reservoir. The efficiency of the proposed numerical method is assessed using debris flow data obtained from the experimental setup and other numerical models. Finally, the sensitivity of the numerical model to changes in the selected parameters is investigated.
τbx =
ρt u u2 + v 2 ϕ2
(3)
τby =
ρt v u2 + v 2 ϕ2
(4)
where φ is the dimensionless velocity coefficient, which is the ratio of sectional mean velocity (U) to the shear velocity (u*) (Takahashi, 2007). This coefficient (φ = U/u*) is in the range of 1–10 for debris flows (Banihabib and Mokhtari, 2001). According to the velocity coefficient, flood volume can be calculated in each section based on Eq. (5): m
Vm =
2.1.2. Governing equations of solid phase The continuity equation of the solid phase can be derived by the mass conservation law. Considering all components of this equation, the entire sediment concentration variations in time and space are considered. This equation can be defined as follows:
∂ (ch) ∂Z ∂ (uhc ) ∂ (vhc ) + + =0 + Cs ∂t ∂t ∂x ∂y
Cs − c ∂ (Z ) ∂ (c ) ∂c ∂c + +u +v =0 h ∂t ∂t ∂x ∂y
where
⎡h⎤ U = ⎢uh ⎥, ⎢ ⎣ vh ⎥ ⎦ ⎡ uh ⎤ 2 E = ⎢ βx u h ⎥, ⎢ β uvh ⎥ ⎣ y ⎦ vh ⎡ ⎤ F = ⎢ βx uvh ⎥, ⎢ β v 2h ⎥ ⎣ y ⎦
c=
τ bx ρt τ by ρt
(7)
Different equations have been presented to calculate mean volumetric concentration including those by Banihabib (1997), Hashimoto and Hirano (1995) and Takahashi et al. (1992). In the equation of Takahashi et al. (1992) the static force is considered to estimate the debris flow transport rate. This equation employs the mean diameter of the sediment and the slope of the channel to estimate the sediment concentration. Hashimoto and Hirano (1995) considered shear stress, the physical characteristics of the sediment and also gradient for the estimation of the sediment transport rate. Investigations performed by Banihabib (1997) indicated that Takahashi equation may underestimate the sediment transport rate and Hashimoto and Hirano's formula overestimates the sediment transport rate (Banihabib, 1997). In this paper, the formula presented by Banihabib (1997) is used where the dimensionless sediment transport rate is dependent upon the relative depth and the dimensionless bed shear stress. Banihabib’s equation is defined as follows:
(1)
⎤ ⎥ − gh sin θx ⎥ ⎥ − gh sin θy ⎥ ⎦
(6)
where c is the mean volumetric concentration of sediment over the entire flow depth and Cs is the sediment concentration in the bed. The continuity equation of sediment phase (Eq. (6)) can be rewritten in the form of Eq. (7), considering the continuity equation of the liquid phase (Banihabib et al., 1996).
2.1.1. Governing equations of the liquid phase The shallow water equations, which are derived from the depthintegrated Navier-Stokes equations for an incompressible fluid, are widely used in computational fluid dynamics. Governing equations for hyper-concentrated sediment-laden flow can be acquired by considering the interaction of flow, sediment transport, and bed level variations in shallow water equations (Wu, 2007; Wu et al., 2012). Accordingly, the matrix form of the two-dimensional mass and momentum conservation equations for debris flow can be described as follows (Takahashi, 2007):
∂U ∂E ∂F + + ST = 0 + ∂t ∂x ∂y
(5)
where Ri is the hydraulic radius, S is channel bed slope, A is flow crosssection, i is time counter, and m is a number of time steps that flow depth is measured.
2.1. Governing equations
∂Z
gRi S Ai dt
i=1
2. Material and methods
⎡ ∂t ⎢ ∂ (Z + h) ST = ⎢ gh ∂x cos θx + ⎢ ⎢ gh ∂ (Z + h) cos θy + ∂y ⎣
∑ϕ
qs q
qs sgd3
(8) .65
h = α ⎛ ⎞ τ∗1.5 ⎝d⎠
(9)
where c is the definition of sediment transport concentration, q is the flow discharge per unit width, qs is the sediment discharge per unit width, α is a dimensionless coefficient, d is the sediment diameter, s is the relative submerged unit weight of the sediment and τ∗ is the dimensionless shear stress. The dimensionless shear stress is calculated from Eq. (10):
(2)
where h is the flow depth, u is the flow velocity in the x-direction, v is the flow velocity in the y-direction, g is the acceleration due to gravity, θx, and θy are inclination of the original bed surface, in the x-axis and yaxis directions, βx and βy are momentum correction factors, Z is bed 66
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
τ∗ =
(u2 + v 2) ϕ2gsd
predictor and a corrector. The flow variables at the k time step are known, and their values are calculated in time step k + 1. Thus, for grid points i, j, the backward (Eq. (13)) and forward (Eq. (14)) finite-difference approximation is used for the numerical integration of Eq. (1) in the predictor and corrector stages (Chaudhry, 2007).
(10)
The following equation can be acquired by substituting the shear stress of Eq. (10) into Eq. (9).
qs sgd3
1.5
h .65 u2 + v 2 ⎞ = α ⎛ ⎞ ⎛⎜ 2 ⎟ ⎝ d ⎠ ⎝ ϕ gsd ⎠
Ui∗, j = Uik, j −
(11)
2≤i≤N Δt Δt ∇x Eik, j − ∇y Fik, j − ΔtSik, j ⎧ ⎨ Δx Δy ⎩2 ≤ j ≤ M
(12)
1≥i≥N−1 Δt Δt Δx Ei∗, j − Δy Fi∗, j − ΔtSi∗, j ⎧ ⎨ Δx Δy ⎩1 ≥ j ≥ M − 1
(13)
The potential of sediment transport by debris flow per unit of width (qs) can be calculated by Eq. (11). Therefore, by substituting (qs) in Eq. (8), the mean volumetric concentration can be obtained.
k Ui∗∗ , j = Ui, j −
2.2. Numerical methods
where U* and U** are the midway values of U, Δx and Δy are the size of the spatial grid, in the x-and y-directions and ∇x and Δx are backward and forward differences, which are defined as follows:
One solution method that can be employed to solve the governing equations of this study is the decoupled model. This method assumes that change in one variable during a given time step is small enough that its effect on the other variables can be assumed negligible (Cao et al., 2002; Kassem and Chaudhry, 1998). Therefore, the flow equations and continuity equation of sediment phase can be solved separately. In the decoupled method, flow equations are solved at first, and subsequently, the flow variables obtained from solving flow equations are used to predict the bed level changes. According to above-mentioned method, Eq. (2) was solved by MacCormack scheme assuming bed level change rate is the same as its value obtained in the previous time step. MacCormack scheme was used in the various studies for debris flow simulation (Campisano et al., 2014; He et al., 2016; Ouyang et al., 2015). After calculating flow variables by Eq. (2), the bed level and sediment concentration were estimated using Banihabib et al. (1996) method. The previous steps were repeated iteratively during a computational time step until the changes between two continuous bed levels estimations became less than 0.01 maximum depth of debris flow. A flowchart, which explains the computational sequences of the numerical method, is shown in Fig. 1.
Δx Ui, j = U(i + 1, j) − U(i, j)
(14)
∇x Ui, j = U(i, j) − U(i − 1, j)
(15)
The new values of vector U at time k + 1 are attained by the following equation:
Uik, j+ 1 =
1 ∗ (Ui, j + Ui∗∗ ,j ) 2
(16)
By using a second-order finite-difference scheme like the MacCormack scheme, some numerical oscillations occur at the watersurface discontinuities. For smoothing these oscillations, some schemes are proposed such as Total Variation Diminishing Property (TVDP) schemes (Chakravarthy and Osher, 1985), Monotonic Upstream Scheme For Conservation Laws (MUSCL) (Van Leer, 1974) and compact schemes and Essentially Non-Oscillatory (ENO) schemes (Liu et al., 1994). In this research, to decrease these oscillations, the Jameson artificial viscosity (Jameson et al., 1981) was selected for its simplicity. The artificial viscosity dissipates fluctuations in the zones with high variations, but it does not affect quiescent zones. The amount of artificial viscosity can be regulated by adjusting dissipation constants (Kx and Ky) in x-and y-directions. The details of this procedure can be found in (Chaudhry, 2007). For the numerical stability of the MacCormack scheme, the CourantFriedrichs-Lewy criterion (CFL) should be equal or less than one at all of the simulation nodes (Chaudhry, 2007). Thus, the time step is restricted by the following condition:
2.2.1. Numerical integration of the flow equations The MacCormack scheme is composed of a two-stage chain of a
CFL =
Δt max[(|u| + min[dx , dy]
gh ), (|v| +
gh )] ≤ 1
(17)
In the simulation, CFL number was set to be less than the unit. 2.2.2. Determination of the bed level variation In the reservoir of dams in steep channels, both subcritical and supercritical flows can occur. Banihabib et al. (1996) proposed a scheme that considers both the subcritical and supercritical condition for integrating the continuity equation of the solid phase (Eq. (7)). In this scheme, the forward difference scheme is considered for the supercritical state, and the backward difference scheme is used for the subcritical state. According to this scheme, Eq. (7) is discretized as shown in Eq. (18).
Zik, j+ 1 = Zik, j +
Δt . hik, j+ 1
k+1 − c k i, j
⎡ ci, j ⎣
(cik, j+ 1 − Cs ) ⎢
Δt
+ AFx . uik, j+ 1
(c k + 1 − cik−+1,1j ) k + 1 i, j
+ (1 − AFx ). ui, j
(c k + 1 − c k + 1 ) k + 1 i, j + 1 i, j
AFy . vi, j Fig. 1. A flow chart showing the computational sequence of the two-dimensional numerical model for debris flow in slit dam reservoir.
Δy
Δx
Δx
+
+(1 − AFy ). vik, j+ 1
where AFx and AFy are calculated as: 67
(cik++1,1j − cik, j+ 1 )
(cik, j+ 1 − cik, j+−11 ) Δy
⎤ ⎥ ⎦
(18)
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
Fig. 2. Experimental setup; (a): plan and (b): lateral view of the flume.
⎧ 0 if AFx = ⎨ 1 if ⎩
u Fr 2 − 1 u Fr 2 − 1
≤0 >0
⎧ 0 if AFy = ⎨ 1 if ⎩
v Fr 2 − 1 v Fr 2 − 1
bed level, and p is a coefficient that considers the sediment distribution with depth (in percent). Considering that at the start point of each time step, cik, j+ 1 was not calculated for computation of the bed level variation (Eq. (18)), sediment concentrations of the last two time steps, k and k − 1, were substituted for those of time steps k and k + 1. Then the sediment concentration for time level k + 1 was computed using Eq. (20) using achieved bed level variation. The bed level variation which obtained in the previous step was used to adjust flow property in the current time step. Adjustment of flow property, bed level variation, and sediment concentrations for a given time step was repeated until bed level difference between two consecutive iterations became less than 0.01 maximum depth of debris flow (Fig. 1).
≤0 >0
(19)
Eq. (19) is an implicit equation of bed level Z and volumetric concentration of sediment c. A control volume equation was used for the computation the volumetric concentration of the sediment at different grid points (Banihabib and Mokhtari, 2001). According to the control volume, volumetric concentration of sediment at time step k + 1 was calculated by Eq. (20).
cik, j+ 1 =
k k (∑ qs+c.i, j dx . dy. hi, j ) − qz dx . dy. hik, j+ 1
(20)
Eq. (20) is taken from continuity equation, where ∑ qs is the total amount of sediment, entering or leaving the element at k + 1, c.ik, j dx . dy. hik, j is available sediment in the element at k, qz is the volume of sediment deposited or eroded in the element at k + 1, and dx . dy. hik, j+ 1 is the volume of the element at k + 1. To compute the net sediment input volume to the element due to the concentration gradient, Eq. (21) was used.
∑ qs = (qs (in) )x − (qs (out ) )x + (qs (in) ) y − (qs (out ) ) y
2.2.3. Initial and boundary conditions Debris flow simulation in slit dam reservoir is like a dam break simulated by Fennema and Chaudhry (1989) which employed the MacCormack Scheme for numerical modeling. Fennema and Chaudhry (1989) demonstrated that the MacCormack Scheme fails when the ratio of tail-water depth to reservoir depth being less than 0.25. Here, flow depth in the initial condition was assumed 0.25 of inflow, and the velocity was determined by Eq. (27)
(21)
where:
(qs (in) )x = (qs (out ) )x =
(qs (in) ) y = (qs (out ) ) y =
1 k (ci, j + cik− 1, j )(uik, j + uik− 1, j )(hik, j + hik− 1, j )Δy. Δt 8 1 k (ci, j + cik+ 1, j )(uik, j + uik+ 1, j )(hik, j + hik+ 1, j )Δy. Δt 8
1 k (ci, j + cik, j − 1 )(vik, j + vik, j − 1 )(hik, j + hik, j − 1 )Δx . Δt 8 1 k (ci, j + cik, j + 1 )(vik, j + vik, j + 1 )(hik, j + hik, j + 1 )Δx . Δt 8
u (i,1) = ϕ gRS
(22)
where φ is the dimensionless velocity coefficient and R is the hydraulic radius. Solid walls boundary was defined according to the reflection boundary condition. The details of reflection boundary for MacCormack scheme are given in Chaudhry (2007). The boundary at the slit of the dam depends on the ratio of the slit width to the channel width. Due to sudden narrowing of the channel by the slit dam, the flow energy is dispersed by creating vortexes upstream of the slit, on both sides. Armanini and Larcher (Armanini and Larcher, 2001) presented an equation to calculate the energy loss, as shown in Eq. (28).
(23) (24) (25)
The sediment fall velocities can be utilized to calculate the sedimentation volume in an element at k + 1, as shown in Eq. (26):
qz = (Zik. j+ 1 − Z ). dx . dy. Cs.
ws. Δt 1 . p hik, j
(27)
(26)
2 Fr 2 ΔE 2 = u ⎡1 − (Fru Rr )−2/3⎤ hu 2 ⎣ 3 ⎦
where ws is the sediment fall velocity, Cs is sediment concentration in 68
(28)
Journal of Hydro-environment Research 27 (2019) 65–74
0.6 2.6
0.9
Sediment concentration in bed (Cs) (-)
ghd
(29) (30)
The accuracy of the numerical model in this study was tested according to a laboratory experiment and alternative numerical models. The laboratory experiment was designed for simulating debris flow in slit dam reservoir, then the debris flow simulated by numerical method of this study. The other numerical models utilized for a comparison were Kanako model (Nakatani et al., 2008) and a model presented by Campisano et al. (2014). Kanako model applied by Liu et al. (2012) to analyze the effects of slit dam on discharge reduction. Campisano et al. (2014) proposed a numerical method for evaluating the slit dam trapping efficiency. The data obtained from numerical method of this study were compared with other models outputs using statistical indices, correlation coefficient (R2) and Root Mean Square Error (RMSE).
1
2.3.1. Experimental setup The experimental setup was composed of four major parts, including a tilting flume, sediment tank, a slit dam and ultrasonic water level gauges (Fig. 2). The slit dam was made of plexiglass and placed 140 cm from the flume end. During the experiment, three ultrasonic water level gauges were recording flow depth with a precision of 0.1 mm. Fig. 2 shows the plan and lateral views of the flume in this experiment. The flume slope was 43% and its height was 25 cm. The width of the slit was 5 cm, and the sediment was uniform sand with sediments mean size (D50) equals 1 mm and density 2.65 gr/cm3. The first, second and third ultrasonic water level gauges were at distances of 410, 210, and 60 cm from the downstream end of the flume, respectively. The first and second water level gauges were located upstream of the dam, and the third one was downstream of the dam reservoir. The main characteristics of the experiment are displayed in Table 1. In the experiments, first, sediment with a known weight was placed on the sediment tank (Fig. 2). Then, the sediment was slowly saturated with water until a little water covers sediment with the horizontal free surface. This volume of water and sediment formed the debris flood volume in the experiment. The upstream of the flume moved upward to a specific height until the required slope was achieved. Then, the gate of the sediment tank was opened suddenly, the flood moved toward the slit dam, and some of the sediment was deposited in the dam reservoir (Banihabib and Bahram, 2009). To have a precise measurement of the sediment thickness deposited in the dam reservoir, a digital point gauge was used. This point gauge was moved vertically along the flume and measured three points in each section (aa', bb' and cc'. sections Fig. 3). These three points were located in the middle, right and left of the section. The right and the left points were 9 cm away from the middle of the section. Fig. 3 shows the plan and lateral views of debris flow deposition in slit dam reservoir schematically. Debris flow occurred in experiment was fully developed according to the observations and criteria proposed by Takahashi (1991). Debris flow is fully developed when following ratio (Eq. (31)) is more than 0.15 (Takahashi 1991).
25
5
0.085 30 43
5
Sediments mean size D50 (mm) Sediment weight (kg) Flood volume (m3) Flume width (cm) Flume bed slope (%)
Slit width (cm)
2⎡ u 2 hu + u − ΔE⎤ ⎥ 3⎢ 2g ⎣ ⎦
2.3. Model test
Flood characteristic
Dam height (cm)
hd =
ud =
Dam characteristics
Table 1 Main characteristics of the experimental setup.
where ΔE is the energy loss, hu is the flow depth, Fru is Froude number at node u located upstream of the slit, and Rr is the ratio of the channel width to the slit width. Since critical velocity occurs at the slit, the flow depth in the slit is equal to two-thirds of the specific energy at the slit. Here, the energy at the slit was obtained from the difference of the energy upstream and the energy loss determined by Eq. (28). Thus, the boundary condition at the slit was determined as:
Sediment density (g/ cm3)
Sediment distribution in depth (p) (%)
M.E. Banihabib and F. Nazarieh
69
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
Fig. 3. Sketch of debris flow sedimentation in slit dam reservoir, (a) plan and (b) lateral views. 0.5
0.007
Discharge Volumetric concentration
Discharge (m^2/s)
0.005
0.4
0.3
0.004
0.003
0.2
0.002 0.1
Volumetric concentration (percent)
0.006
0.001
0.0
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 2.2 2.4 2.6 2.8 3.0 3.2 3.4 3.6 3.8 4.0 4.2 4.4 4.6 4.8 5.0 5.2 5.4 5.6 5.8 6.0
0.000
Time (s) Fig. 4. Flow and sediment concentration at second ultrasonic water level gauge used as the upstream boundary condition.
sin(S ) > 0.15 (ρs / ρ) − 1
2.3.2. Numerical test The Kanako model was used by Liu et al. (2012) for assessing the performance of slit dams in Tanjutani Gully, Kyoto City, Japan. The average width of this Gullys’ was 6 m, its slope at the investigated section was about 14.7%. and sediment mean size (D50) was 300 mm. Tested dams had three widths of slit including 0.4, 0.8 and 1.2 m. The upstream boundary conditions assume to be outflow hydrograph from the previous check dam. In the current study, slit dam examples in Tanjutani Gully, was simulated by the proposed numerical model considering 0.1 × 0.1 m gride size mesh and 0.004 s time step. Simulated outflow hydrograph downstream of slit dam were compared with the digitized data obtained from Kanako model in Liu et al. (2012) research. Campisano et al. (2014) proposed a numerical method for exploring slit dams trapping efficiency. The case study of Campisano et al. (2014) was characterized with a river with a rectangular section width of 10 m, length of 300 m, initial slope of 5% and sediment mean size (D50)
(31)
where S is channel bed slpe, ρs is the sediment material density, and ρ is water density. In the experiments of this study, the above mentioned ratio was more than 0.25 at downstream of the second ultrasonic water level gauge. The upstream boundary condition of the reservoir was defined according to the flow depths recorded by the second ultrasonic water level gauge. Inflow velocity and sediment concentration at upstream were calculated by Eqs. (27) and (8). The flow and sediment concentration at second ultrasonic water level gauge were determined, as shown in Fig. 4. The explained experiment was simulated by the proposed numerical method considering 0.01 m × 0.01 m grid size and time step 0.001 s.
70
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh 0.0030
Discharge (m^2/s)
0.0025
0.0020
0.0015
0.0010
Experimental data This study model
0.0005
0.0000 0
5
10
15
20
25
30
35
Time (s)
Fig. 5. Comparison of outflow hydrograph simulated by current study numerical model and the experimental data. 0.06
(a)
supposed to be 0.005 s.
Flume bed Experimental data This study model
0.05
3. Results and discussion
Elevation (m)
0.04
3.1. Validation
0.03
3.1.1. Validation based on experimental data Four parameters φ, α, Kx and Ky in the numerical model were calibrated based on the experimental data. Firstly, the dimensionless velocity coefficient (φ) was calibrated considering Eq. (5). The flood volume and flow depth in ultrasonic water level gauges locations are known, therefore, velocity coefficient (φ) is only unknown in Eq. (5) and could be calculated. This coefficient was calibrated as 5.6 for the location of second ultrasonic water level gauge. Secondly, the dimensionless coefficient α was calibrated to the amount of 1.5 by equating the transported sediment volume in both the experimental and mathematical model. Finally, for dissipation constants, Kx and Ky values of 0.49 and 0.3 were obtained by trial and error. To investigate the capability of the proposed numerical method to simulate debris flow passing the slit dam reservoir, simulated outflow hydrograph was compared with the experimental hydrograph. Fig. 5 shows comparison of hydrographs in third ultrasonic water level gauge location. This figure illustrates that the simulated flood peak discharge was estimated earlier and greater than the observed one. The difference between simulated flood peak and experimental flood peak may happen as a result of the larger time step used in the experimental data acquisition (0.2 s) in comparison to the fine time step applied in numerical method (0.001 s). This may cause losing real peak flow in the experimental data. It is inferred from Fig. 5 that the numerical model simulates the general trend of the outflow hydrograph correctly as values of statistical criteria, R2 was 0.93 and RMSE was 0.000182 m2/s. Numerical model capability in simulating sedimentation was assessed in aa', bb' sections which are represented in Fig. 2a. Comparisons of the simulated bed profile and measured data along these sections at the end of the experiment are shown in Fig. 6. Evaluation criteria at section aa' showed an acceptable bed profile simulation (R2 = 0.98 and RMSE = 0.0037 m). However, slight sedimentation was predicted at section aa' by the numerical model while sedimentation was not observed on experimental setup (Fig. 6a). The reason may be due to the solution method which could not fully address continuity of water and sediment simultaneously. Hence, the error of flow properties affects the simulation of sedimentation along section aa'. In Fig. 6(b), the bed profile of the numerical model and experimental data are compared along section bb'. The numerical model simulated the location of the sedimentation starting point and the bed profile accurately with R2 = 0.98 and RMSE = 0.0012 m. Hydraulic jump has a significant effect on deposition in slit dam
0.02
0.01
0.00 0.00
0.02
0.04
0.06
0.08
0.10
Distance from dam (m) 0.06
(b)
Flume bed Experimental data This study model
0.05
Elevation (m)
0.04
0.03
0.02
0.01
0.00 0.00
0.02
0.04
0.06
0.08
0.10
Distance from dam (m)
Fig. 6. Comparison of simulated bed profile and measured experimental data along (a) section aa' and (b) section bb'.
50 mm. It was assumed that slit dam was positioned at 250 m from the upstream of the river. Various contraction ratios of the slit width to the channel width were investigated. The slit widths used for comparison include 2, 3 and 4 m. Inflow was defined by triangular hydrograph with flood peak of 60 m3/s, time to peak of 3600 s and the total duration of 7200 s. Trapping efficiency E (t ) was calculated according to the following equation:
E (t ) =
Vdep (t ) Vsed (t )
(32)
where Vdep (t ) is the volume of sedimentation deposited in slit dam reservoir at time t and Vsed (t ) is volume of sediment entered to the reservoir until time t. Variation of trapping efficiency by time, depicted by Campisano et al. (2014), was digitized to compare with the result of the numerical method presented in the current study. The grid size used for the proposed numerical method was 0.5 × 0.5 m and time step was 71
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh 0.007
Discharge (m^2/s)
0.006
0.25
0.005
0.20
0.004
0.15
0.003
0.10
0.002
0.05
0.001
0.00
Distance of hydraulic jump from dam (m)
0.30
Inflow hydrograph Outflow hydrograph Hydraulic jump distance from dam
0.000 0
2
4
6
8
Time (s) Fig. 7. Variations of the simulated hydraulic jump location with time along section bb'. 50
50
(a)
Inflow hydrograph Outflow hydrograph (This study) Outflow hydrograph (Kanako)
30
20
10
0
Inflow hydrograph Outflow hydrograph (This study ) Outflow hydrograph (Kanako)
40
Discharge (m^3/s)
Discharge (m^3/s)
40
(b)
30
20
10
0
200
400
600
800
1000
1200
1400
1600
0
1800
Time (s)
200
400
600
800
1000
1200
1400
1600
1800
Time (s)
50
(c)
Inflow hydrograph Outflow hydrograph (This study) Outflow hydrograph (Kanako)
40
Discharge (m^3/s)
0
30
20
10
0
0
200
400
600
800
1000
1200
1400
1600
1800
Time (s) Fig. 8. Comparison of outflow hydrograph simulated by this study numerical model and Kanako model for the slit width of (a) 0.4 m, (b) 0.8 m and (c) 1.2 m.
reservoir, in fact by flow transition from supercritical to subcritical, flow transport capacity decrease and it causes a part of transported sediment deposit in the reservoir (Armanini and Larcher, 2001). The variation of hydraulic jump location along section bb' is shown versus time in Fig. 7. It should be noted that the hydraulic jump displacement along section bb' was the same as section aa'. The figure shows that, as the flood was filling the reservoir, the hydraulic jump moved upstream, and this process continued until most of the dam reservoir was emptied from the flood. Then, the hydraulic jump location moved downstream.
3.1.2. Validation based on other numerical models The outflow hydrographs generated by Kanako model in Liu et al. (2012) study and outflow hydrographs simulated with the numerical model of this study for the same example are compared in Fig. 8. These outflow hydrographs belong to the slit width of 0.4 m (Fig. 8a), 0.8 m (Fig. 8b) and 1.2 m (Fig. 8c). Calibrated amounts of φ, α, Kx and Ky parameters in this example were 3.7, 0.9, 0.48 and 0.35 respectively. The statistical criteria R2 and RMSE for the width of 0.4 m were 0.93 and 4.3, for the width of 0.8 were 0.98 and 1.8 and for the width of 72
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh 1.0
Trapping Efficiency (-)
0.8
(2014) numerical models. The statistical criteria, R2 and RMSE, for slit width of 2 m were 0.94 and 0.34, for slit width of 3 m were 0.98 and 0.045 and for slit width of 3 m were 0.99 and 0.038. It can be seen in Fig. 8 that Campisano et al. (2014) model underestimate the trapping efficiency during the decreasing branch of the hydrograph, which dam self-cleaning process happen in this period. It may be as a result of onedimensional simulation of flow behind slit dam in Campisano et al. (2014) model which can not simulate sediment trasport proccess in the corners of the reservoir properly. Considering statistical criteria in all evaluation tests, it can be concluded that there is a good agreement between the result of the proposed numerical model and other numerical models. Furthermore, calibrated parameters (φ, α, r Kx Ky) had restricted a range of variation in all tests. These calibrated values are also located in the range of acceptable values presented in other researches such as Chaudhry (2007) and Takahashi (2007), which can give some support to the credibility of the proposed numerical model.
d=2 m
0.6
d=3 m
0.4
d=4 m
0.2
Campisano Model This study Model
0.0 0
20
40
60
80
100
120
Time (min)
Fig. 9. Time-variation of trapping efficiency at slit dam simulated by this study and Campisano et al. (2014) numerical models at slit width of 2 m, 3 m and 4 m.
Changes in total outflow volume (%)
30
3.2. Sensitivity analysis
(a)
Four input parameters were used in the sensitivity analysis, including velocity coefficient (φ), dimensionless coefficient (α), and coefficients of dispersion Kx and Ky. Each input parameter was modified by ± 5, ± 10, ± 15 and ± 20 percent from its original value (calibrated values based on experimental data), while the other input parameters were kept constant at their calibration values. The numerical method was not stable in some cases, so some sensitivity curves are not complete for the entire range of parameter variations as shown in Fig. 10. Two different simulated outputs were assessed in the sensitivity analysis: the total outflow volume (Fig. 10a) and the total deposition volume trapped in slit dam reservoir (Fig. 10b). Sedimentation has a higher rate of variation (Fig. 10b) in comparison with outflow volume (Fig. 10a). It demonstrates that sedimentation volume was more sensitive to the parameter changes. Generally, the response of outflow and sedimentation volume to changes in parameters was reverse. In other words, when variation in parameters leads to increase in outflow volume, sedimentation volume decreased and vice versa. Sensitivity analysis in Fig. 10 also illustrates that raising in velocity coefficient (φ) leads to increase of outflow volume and decrease of sedimentation while raising in parameter (α) causes more sedimentation and less outflow volume. Changes in outflow and deposition volumes to variation in Kx and Ky seem to be similar. The curves related to these two parameters are so that any changes in them (regardless of increase or decrease of Kx and Ky), result in increase of outflow volume and decrease the volume of trapped sedimentation (Fig. 10).
20
10
0
-10 Velocity coefficient Dimensionless coefficient Coefficients of dispersion Kx
-20
Coefficients of dispersion Ky
-30
-20%
-15%
-10%
-5%
0%
5%
10%
15%
20%
Parameter changes (%) 40
Changes in total deposition volume (%)
(b) 30 Velocity coefficient Dimensionless coefficient Coefficients of dispersion Kx
20
10
Coefficients of dispersion Ky
0
-10
4. Conclusion -20
The MacCormack and Banihabib schemes were combined to simulate debris flow sedimentation in a slit detention dam reservoir. Performance of the proposed numerical methods was evaluated by comparing simulated data of this model with other experimental and numerical models, and the final conclusions can be summarized as follows:
-30 -20%
-15%
-10%
-5%
0%
5%
10%
15%
20%
Parameter changes (%)
Fig. 10. The sensitivity of (a) outflow and (b) sedimentation volume to the numerical model parameters variations.
• The MacCormack scheme and banihabib scheme could simulate the
1.2 m were 0.99 and 0.54. Furthermore, Fig. 8 illustrates that by increasing the width of slit, agreement between two models modifies. The time-variation of the trapping efficiency of the slit dam computed by the numerical model of this study, was compared with the result obtained from Campisano model (Campisano et al., 2014) in Fig. 9. The calibrated parameters in this case were 5.6 for φ, 1.3 for α, 0.42 for Kx and 0.3 for Ky. The results presented in Fig. 8 illustrates a generally good agreement between this study and Campisano et al.
• • 73
debris flow regime variation from supercritical to subcritical and the movement of the hydraulic jump in the reservoir of the slit dam. According to the comparison of outflow hydrographs in investigated models, the slit dam boundary condition was capable of simulating the general trend of the outflow hydrograph. The simulation of the bed level changes was fairly appropriate along the side longitudinal section, and it was acceptable along the central longitudinal section of the experimental setup. The existence of error in bed level simulation at the central section may be due to the
Journal of Hydro-environment Research 27 (2019) 65–74
M.E. Banihabib and F. Nazarieh
•
solution method, which could not fully address continuity equations. Sensitivity analysis revealed that simulated sedimentation is more sensitive to the parameter variation than outflow hydrograph. Generally, outflow and sedimentation showed an opposite response to the parameters changes.
channel, 50th annual conference of JSCE Japan, pp. 495–500. Hassan-Esfahani, L., Banihabib, M.E., 2016. The impact of slit and detention dams on debris flow control using GSTARS 3.0. Environ. Earth Sci. 75, 328. https://doi.org/ 10.1007/s12665-015-5183-z. He, S., Ouyang, C., Liu, W., Wang, D., 2016. Coupled model of two-phase debris flow, sediment transport and morphological evolution. Acta Geol. Sin. Engl. 90, 2206–2215. https://doi.org/10.1111/1755-6724.13031. Hunt, B., 1994. Newtonian fluid mechanics treatment of debris flows and avalanches. J. Hydraul. Eng. 120, 1350–1363. https://doi.org/10.1061/(ASCE)0733-9429(1994) 120:12(1350). Hürlimann, M., Rickenmann, D., Graf, C., 2003. Field and monitoring data of debris-flow events in the Swiss Alps. Can. Geotech. 40 (1), 161–175. https://doi.org/10.1139/ t02-087. Jameson, A., Schmidt, W., Turkel, E., 1981. Numerical solution of the Euler equations by finite volume methods using Runge Kutta time stepping schemes, 14th fluid and plasma dynamics conference, Palo Alto, California, pp. 1259. Kassem, A.A., Chaudhry, M.H., 1998. Comparison of coupled and semicoupled numerical models for alluvial channels. J. Hydraul. Eng. 128 (4), 794–802. https://doi.org/10. 1061/(ASCE)0733-9429(1998)124:8(794). Lin, X., Huo, M., Zhou wen, J., Cao, T., Yang Rong, F., Zhou wei, H., 2017. An experimental study on controlling post-earthquake debris flows using slit dams. Environ. Earth. Sci. 76 (22), 780. https://doi.org/10.1007/s12665-017-7140-5. Liu, J., Nakatani, K., Mizuyama, T., 2012. Hazard mitigation planning for debris flow based on numerical simulation using Kanako simulator. J. Mt. Sci. 9, 529–537. https://doi.org/10.1007/s11629-012-2225-9. Liu, X.-D., Osher, S., Chan, T., 1994. Weighted essentially non-oscillatory schemes. J. Comput. Phys. 115, 200–212. https://doi.org/10.1006/jcph.1994.1187. Michelini, T., Bettella, F., D’Agostino, V., 2017. Field investigations of the interaction between debris flows and forest vegetation in two Alpine fans. Geomorphology 279, 150–164. https://doi.org/10.1016/j.geomorph.2016.09.029. Nakatani, K., Wada, T., Satofuka, Y., Mizuyama, T., 2008. Development of “Kanako 2D (Ver. 2.00)”, a user-friendly one- and two-dimensional debris flow simulator equipped with a graphical user interface. Int. J. Eros. Control. Eng. 1, 10. Ouyang, C., He, S., Tang, C., 2015. Numerical analysis of dynamics of debris flow over erodible beds in Wenchuan earthquake-induced area. Eng. Geol. 194, 62–72. https:// doi.org/10.1016/j.enggeo.2014.07.012. Paik, J., 2015. A high resolution finite volume model for 1D debris flow. J. Hydro-environ. Res. 9, 145–155. https://doi.org/10.1016/j.jher.2014.03.001. Pan, H., Wang, R., Huang, J., Ou, G., 2013. Study on the ultimate depth of scour pit downstream of debris flow sabo dam based on the energy method. Eng. Geol. 160, 103–109. https://doi.org/10.1016/j.enggeo.2013.03.026. Papa, M., Sarno, L., Vitiello, F., Medina, V., 2018. Application of the 2D depth-averaged model, FLATModel, to pumiceous debris flows in the Amalfi Coast. Water 10 (9), 1159. https://doi.org/10.3390/w10091159. Pirulli, M., Pastor, M., 2012. Numerical study on the entrainment of bed material into rapid landslides. Geotechnique 62, 959–972. https://doi.org/10.1680/geot.10.P.074. Qiu, F., Huang, J., Li, Y., Han, Z., Wang, W., Chen, G., Qu, X., Su, B., 2018. Protecting highway bridges against debris flows using lateral berms: a case study of the 2008 and 2011 Cheyang debris flow events, China. Geomat. Nat. Haz. Risk 9, 196–210. https://doi.org/10.1080/19475705.2017.1414718. Rickenmann, D., Laigle, D., McArdell, B.W., Hübl, J., 2006. Comparison of 2D debris-flow simulation models with field events. Comput. Geosci. 10, 241–264. https://doi.org/ 10.1007/s10596-005-9021-3. English. Schwindt, S., De Cesare, G., Boillat, J.-L., Schleiss, A., 2016. Physical modelling optimization of a filter check dam in Switzerland. Conference Proceedings, INTERPRAEVENT 2016, Hazard and Risk Mitigation 828–836. Takahashi, T., 1991. Debris flow. In: IAHR Monograph Serie. Balkema, Rotterdam, pp. 1–165. Takahashi, T., 2007. Debris Flow: Mechanics, Prediction and Countermeasures, second ed. CRC Press, United States. Takahashi, T., Nakagawa, H., Harada, T., Yamashiki, Y., 1992. Routing debris flows with particle segregation. J. Hydraul. Eng. 118, 1490–1507. https://doi.org/10.1061/ (ASCE)0733-9429(1992) 118:11(1490). Tsai, Y., 2007. A debris-flow simulation model for the evaluation of protection structures. J. Mt. Sci. 4, 193–202. https://doi.org/10.1007/s11629-007-0193-2. English. Van Leer, B., 1974. Towards the ultimate conservative difference scheme. II. Monotonicity and conservation combined in a second-order scheme. J. Comput. Phys. 14, 361–370. https://doi.org/10.1016/0021-9991(74)90019-9. VanDine, D., 1996. Debris flow control structures for forest engineering, Research Branch, Ministry of Forests Research Program, Province of British Columbia. Wehrmann, H., Hübl, J., Holzinger, G., 2006. Classification of dams in torrential watersheds, Disaster Mitigation of Debris Flows, Slope Failures and Landslides, Universal Academy Press, Tokyo, Japan., pp. 829–838. Wu, W., 2007. Computational River Dynamics, first ed. CRC Press, United States. Wu, W., Marsooli, R., He, Z., 2012. Depth-averaged two-dimensional model of unsteady flow and sediment transport due to noncohesive embankment break/breaching. J. Hydraul. Eng. 138, 503–516. https://doi.org/10.1061/(ASCE)HY.1943-7900. 0000546. Zeng, Q.L., Yue, Z.Q., Yang, Z.F., Zhang, X.J., 2009. A case study of long-term field performance of check-dams in mitigation of soil erosion in Jiangjia stream, China. Eng. Geol. 58 (4), 897–911. https://doi.org/10.1007/s00254-008-1570-z.
Generally, the numerical method presented in this article can provide some insights on debris flow sedimentation upstream of the slit dams. However, because of limited availability of measured debris flow data, the model in this study was just compared with the result of an experimental data and one-dimensional numerical models. Further comparison with field data and other numerical models is required to enhance the reliability of the established model. Acknowledgments The authors appreciate Dr. Mary Kay Camarillo of the University of the Pacific and Dr. Seyyed-Mehdi Hashemi-Shahedani of University of Tehran for editing and reviewing this paper. Moreover, we thank Iran’s Ministry of Energy, Water, and Wastewater Standards and Projects Bureau for supporting this research through contract No. 720/1015/87. Appendix A. Supplementary data Supplementary data to this article can be found online at https:// doi.org/10.1016/j.jher.2019.09.001. References Armanini, A., Larcher, M., 2001. Rational criterion for designing opening of slit-check dam. J. Hydraul. Eng. 127, 94–104. https://doi.org/10.1061/(ASCE)07339429(2001) 127:2(94). Banihabib, M.E., 1997. Sediment transport of mud flow. In: 1st Iranian Hydraulic Conference. Iranian Hydraulic Association, Tehran, Iran, pp. 1–10. Banihabib, M.E., Bahram, E., 2009. Experimental Analyses of Sedimentation in the Slit Dam Reservoir. Great Rivers, Kansas City, Missouri, USA, World Environmental and Water Resources Congress, pp. 1–12. Banihabib, M.E., Forghani, A., 2017. An assessment framework for the mitigation effects of check dams on debris flow. CATENA 152, 277–284. https://doi.org/10.1016/j. catena.2017.01.018. Banihabib, M.E., Hirano, M., Harada, T., 1996. A mathematical model for hazard zone mapping of debris flow. Geo-Inf. 7, 87–90. https://doi.org/10.6010/ geoinformatics1990.7.1-2_87. Banihabib, M.E., Mokhtari, A.H., 2001. Numerical simulation in detention dams during a high concentrated flow. In: 3rd International Symposium on Environmental Hydraulics. Arizona, USA, pp. 1–6. Campisano, A., Cutore, P., Modica, C., 2014. Improving the evaluation of slit-check dam trapping efficiency by using a 1D unsteady flow numerical model. J. Hydraul. Eng. 140, 04014024. https://doi.org/10.1061/(ASCE)HY.1943-7900.0000868. Cao, Z., Day, R., Egashira, S., 2002. Coupled and decoupled numerical modeling of flow and morphological evolution in alluvial rivers. J. Hydraul. Eng. 128 (3), 306–321. https://doi.org/10.1061/(ASCE)0733-9429(2002)128:3(306). Catella, M., Paris, E., Solari, L., 2005. Case study: efficiency of slit-check dams in the mountain region of versilia basin. J. Hydraul. Eng. 131 (3), 145–152. https://doi.org/ 10.1061/(ASCE)0733-9429(2005)131:3(145). Chakravarthy, S., Osher, S., 1985. A new class of high accuracy TVD schemes for hyperbolic conservationlaws. In: 23rd Aerospace Sciences Meeting, Reno, Nevada, pp. 363. Chaudhry, M.H., 2007. Open-Channel Flow, Second ed. Springer, US. Chen, S.-C., Tfwala, S., 2018. Evaluating an optimum slit check dam design by using a 2D unsteady numerical model, E3S Web of Conferences. Fennema, R.J., Chaudhry, M.H., 1989. Implicit methods for two-dimensional unsteady free-surface flows. J. Hydraul. Res. 27, 321–332. https://doi.org/10.1080/ 00221688909499167. Ghilardi, P., Natale, L., Savi, F., 2001. Modeling debris flow propagation and deposition. Phys. Chem. Earth Part C 26, 651–656. https://doi.org/10.1016/S1464-1917(01) 00063-0. Han, Z., Chen, G., Li, Y., Tang, C., Xu, L., He, Y., Huang, X., Wang, W., 2015. Numerical simulation of debris-flow behavior incorporating a dynamic method for estimating the entrainment. Eng. Geol. 190, 52–64. https://doi.org/10.1016/j.enggeo.2015.02. 009. Hashimoto, H., Hirano, M., 1995. Sediment transport rate and resistance on a steep
74