1D bed evolution model

Flow model

The shallow water flow model (St. Venant Eq.) is used to calculate the 1D flow field. The governing equation is discretized on a Staggered grid system.

Governing equation

$$
\frac{\partial h}{\partial t}+\frac{\partial uh}{\partial x} = 0
$$
$$
\frac{\partial u}{\partial t}+u\frac{\partial u}{\partial x} = -g\frac{\partial H}{\partial x} -\frac{gn_m^2u|u|}{h^{4/3}}
$$
where, $t$: time, $x$: downstream cordinate, $u$: flow velocity, $h$: water depth, $H$: water surface elevation, $g$: gravitational acceleration, and $n$: Manning's roughness coefficient.

Morphodynamics

The model deals with bed- and suspended-load transport. The bed evolutaiton is then calculated using Exner equation.

Governing equation

$$
(1-\lambda)\frac{\partial\eta}{\partial t}+\frac{\partial q_b}{\partial x} = D-E
$$
where, $\eta$: bed elevation $q_b$: bedload flux, $\lambda$: porosity of bed, $E$: entrainment rate, $D$: deposition rate.
The bedload transport rate is calculated by a modified version of Meyer, Peter and Muller Eq. ( Wong and Parker, 2006 ).
$$
q_b = 4(\theta-\theta_c)^{3/2}\sqrt{Rgd^3}
$$
where, $\theta$: Shields number, $\theta_c$: critical Shields number, $R$: specific weight of sediment in fluid, $d$: grain size.
The critical Shields number, $\theta_c$ is obtained by Iwagaki's formula ( Iwagaki, 1956 ).
$$
\begin{eqnarray*}
u_{*c}^2 & = & 0.14Rgd &\quad (R_{ep} < 2.14)\\
& = & (0.1235Rg)^{25/32}\nu^{7/16}d^{11/32} &\quad (2.14\leq R_{ep} < 54.2)\\
& = & 0.034Rgd &\quad (54.2\leq R_{ep} < 162.7)\\
& = & (0.01505Rg)^{25/22}\nu^{-3/11}d^{31/22} &\quad (162.7\leq R_{ep} < 671)\\
& = & 0.05Rgd &\quad (671 < R_{ep} )
\end{eqnarray*}
$$
$$
R_{ep} =\frac{\sqrt{Rgd^3}}{\nu}
$$
$$
\theta_c =\frac{u_{*c}^2}{Rgd}
$$
where, $\nu$: kinematic viscosity coefficient.
The depth-averaged suspended sediment concentration, $c$, is calculated by the following advection equation.
$$
\frac{\partial ch}{\partial t}+\frac{\partial cuh}{\partial x} = E-D
$$
The entrainment rate is calculated by Garcia and Parker relation ( Garcia and Parker, 1991 ),
$$
E = w_f E_{*}\\
E_* =\frac{AZ}{1+AZ/0.3}\\
Z =\frac{u_{*}}{w_f}R_{ep}^{0.6}\\
A = 1.3\times 10^{-7}
$$
The fall velocity of sediment in fluid is calculated by the relation of Dietrich ( Dietrich, 1982 )
$$
\begin{equation}
w_f = exp[-b_1+b_2 lnR_{ep} - b_3 [lnR_{ep}]^2 - b_4 [lnR_{ep}]^3 + b_5 [lnR_{ep}]^4]\sqrt{Rgd}
\end{equation}
$$
where, $b_1=2.891394, b_2=0.95296, b_3 = 0.056835, b_4=0.002892, b_5=0.000245$.

Discretization

The governing equations are discretized on the staggered grid system below,

meaning that the scalar variables like water depth, bed elevation, suspended sediment concentration are defined at the center of cells (i.e., at the open circle in the above figure), and the vector variables like flow velocity and bedload flux are defined at the side of the cells (i.e., open triangle in the figure).
For discretizing the momentum equation, the 1st order upwind scheme is used for the advection term and the roughness term is discretized by using semi-implicit method to avoid Vasiliev instability. In the following equations, the $n$ denotes the discretized time step.
$$
\begin{equation}
u_i^{n+1} = u_i^{n}-\left[u_i^n\frac{u_i^n-u_{i-1}^n}{\Delta x}
+g\left(\frac{\eta_{i+1}-\eta_i}{\Delta x}+\frac{h_{i+1}-h_i}{\Delta x}
+\frac{gn_m^2u_i^nu_i^{n+1}}{\bar{h}_i^{4/3}}\right)\right]\Delta t
\end{equation}
$$
$$
\begin{equation}
u_i^{n+1} =\frac{u_i^{n}-\left[u_i^n\frac{u_i^n-u_{i-1}^n}{\Delta x}
+g\left(\frac{\eta_{i+1}-\eta_i}{\Delta x}+\frac{h_{i+1}-h_i}{\Delta x}
\right)\right]\Delta t}{1+gn_m^2u_i^nu_i^{n+1}\Delta t/\bar{h}_i^{4/3}}
\end{equation}
$$
where, $\bar{h}_i=(h_i+h_{i+1})/2$.
The discretized continuity equation is;
$$
h_i^{n+1} = h_i^n-\frac{q_i-q_{i-1}}{\Delta x}\Delta t\\
q_i = u_i\bar{h}_i
$$
The Exner equation is discretized as;
$$
\eta_i^{n+1} =\eta_i^n -\frac{1}{1-\lambda}\frac{q_{bi}-q_{bi-1}}{\Delta x}\Delta t +\frac{1}{1-\lambda}(D_i - E_i)
$$
The transport equation of suspended sediment concentration is discretized by using the donner cell scheme for advection term as;
$$
(ch)_i^{n+1} = (ch)_i^{n} -\frac{q_i c_i-q_{i-1} c_{i-1}}{\Delta x}\Delta t+(E_i-D_i)\Delta t
$$

Example1: 1D bed evolution after a cutoff

Discharge = 200 m3/s
Channel width = 200 m
Manning's n = 0.03
Bed slope = 0.002 (mild reach), 0.01 (steep reach)
Grain size = 2 mm

Example2: Sedimentation in a dam

Discharge = 200 m3/s
Channel width = 200 m
Manning's n = 0.03
Bed slope = 0.002
Grain size = 0.2 mm
Downstream WL = 3 m (fixed)