Sediment Transport 1

Sediment Transport

Overview

For sand transport, the wash-load (i.e. sediment transport which does not contribute to the bed-material) can be assumed to be zero, and therefore, the total-load transport is equal to the sum of the bed- and the suspended-load transports: $q_{tk*}=q_{sk*}+q_{bk*}$.

There are currently three sediment transport models available in CMS:

(1) Equilibrium total load

(2) Equilibrium bed load plus non-equilibrium suspended load, and

(3) Non-equilibrium total-load.

The first two models are single-size sediment transport models and are only available with the explicit time-stepping schemes. The third is multiple-sized sediment transport model and is available with both the explicit and implicit time-stepping schemes.

Equilibrium Total-load Transport Model

In this model, both the bed load and suspended load are assumed to be in equilibrium. The bed change is solved using a simple mass balance equation known as the Exner equation.

${\frac {\rho _{s}(1-\rho _{m}^{'})}{f_{morph}}}{\frac {\partial z_{b}}{\partial t}}={\frac {\partial q_{t*j}}{\partial x_{j}}}+{\frac {\partial }{\partial x_{j}}}\left(D_{s}\mid q_{t*}\mid {\frac {\partial z_{b}}{\partial x_{j}}}\right)$ (2-42)

for $j=1,2;k=1,2,.....,N$ , where N is the number of sediment size classes and

t = time [s]

h = total water depth [m]

$x_{j}$ =Cartesian coordinate in the j^th direction [m]

$q_{t*}$ = equilibrium total-load transport rate [kg/m/s]

$z_{b}$ = bed elevation with respect to the vertical datum [m]

$p_{m}^{'}$ = bed porosity [-]

$f_{morph}$ = morphologic acceleration factor [-]

$\rho _{s}$ = sediment density [~2650 kg/m³ for quartz sediment]

$D_{s}$ = empirical bed-slope coefficient (constant) [-]

Because the model assumes that both the sediment transport is equilibri-um, it only recommended for coarse grids with resolutions larger than 50-100 m where the assumption of equilibrium sediment transport is more appropriate. As mentioned above the equilibrium total-load sediment transport model is a single-size sediment transport model and is only available with the explicit time-stepping scheme. For more information on the equilibrium sediment transport model, the reader is referred to Buttolph et al. (2007).

Equilibrium Bed-load plus Nonequilibrium Suspended Load Transport Model

Calculations of suspended load and bed load are conducted separately. The bed load is assumed to be in equilibrium and is included in the bed change equation while the suspended load is solved through the solution of an advection-diffusion equation. Actually the advection diffusion equation is a non-equilibrium formulation, but because the bed load is assumed to be in equilibrium, this model is referred to the "Equilibrium A-D" model.

Suspended-load Transport Equation

The transport equation for the suspended load is given by

${\frac {\partial (hC)}{\partial t}}+{\frac {\partial (hV_{j}C)}{\partial x_{j}}}={\frac {\partial }{\partial x_{j}}}\left(v_{s}h{\frac {\partial C}{\partial x_{j}}}\right)+E_{b}-D_{b}$ (2-43)

where

t = time[s]

h = water depth [m]

$x_{j}$ = Cartesian coordinate in the

$E_{b}$ entrainment or pick-up function [kg/m²/s]

$D_{b}$ deposition or settling function [kg/m²/s]

The entrainment and deposition functions are calculated as

$E_{b}=-\varepsilon _{s}{\frac {\partial c}{\partial z}}\left|_{z=a}\right.=\omega _{s}c_{a*}$ (2-44)

$E_{b}=\omega _{s}c_{a}(2-45)$

where

z = vertical coordinate from the bed [m]

$\varepsilon _{s}$ = vertical sediment mixing coefficient [m²/s]

c = local sediment concentration [kg/m³]

$\omega _{s}$ = sediment fall velocity [m/s]

$c_{a}$ = calculated sediment concentration at an elevation a above the bed [kg/m³]

$c_{a*}$ = equilibrium (capacity) sediment concentration at an elevation a above the bed [kg/m³]

Bed Change Equation

The bed change is calculated as

${\frac {\rho _{s}(1-p_{m}^{'})}{f_{morph}}}{\frac {\partial z_{b}}{\partial t}}={\frac {\partial q_{b*j}}{\partial x_{j}}}+{\frac {\partial }{\partial x_{j}}}\left(D_{s}\mid q_{t}*\mid {\frac {\partial z_{b}}{\partial x_{j}}}\right)+E_{b}-D_{b}$ (2-46)

where

$z_{b}$ = bed elevation with respect to the vertical datum [m]

$p_{m}^{'}$ = bed porosity [-]

$f_{morph}$ = morphologic acceleration factor [-]

$\rho _{s}$ = sediment density [~2650 kg/m³ for quartz sediment]

$D_{s}$ = empirical bed-slope coefficient (constant) [-]

$q_{bk}=hUC_{tk}(1-r_{sk})$ is the bed load mass transport rate [kg/m/s]

As mentioned above the equilibrium total-load sediment transport model is a single-size sediment transport model and is only available with the explicit time-stepping scheme. For more information on the equilibrium sediment transport model, the reader is referred to Buttolph et al. (2007).

Nonequilibrium Total-Load Transport Model

Total-load Transport Equation

The single-sized sediment transport model described in Sánchez and Wu (2011a) was extended to multiple-sized sediments within CMS by Sánchez and Wu (2011b). In this model, the sediment transport is separated into current- and wave-related transports. The transport due to currents includes the stirring effect of waves; and the wave-related transport includes the transport due to asymmetric oscillatory wave motion as well as steady contributions by Stokes drift, surface roller, and undertow. The current-related bed and suspended transports are combined into a single total-load transport equation, thus reducing the computational costs and simplifying the bed change computation. The two-dimensional horizontal (2DH) transport equation for the current-related total load is

${\frac {\partial }{\partial t}}\left({\frac {hC_{tk}}{\beta _{tk}}}\right)+{\frac {\partial (hV_{j}C_{tk})}{\partial x_{j}}}={\frac {\partial }{\partial x_{j}}}\left[v_{s}h{\frac {\partial (r_{sk}C_{tk})}{\partial x_{j}}}\right]+\alpha _{t}\omega _{sk}(C_{t*k}-C_{tk})(2-47)$

for j=1,2; k=1,2,......N , where N is the number of sediment size classes and

t = time [s]

h = water depth [m]

$x_{j}$ = Cartesian coordinate in the j^th direction [m]

$c{tk}$ = depth-averaged total-load sediment mass concentration for size class k defined as $c_{tk}=q_{tk}/(Uh)$ in which $q_{tk}$ is the total-load mass transport [kg/m³]

$c_{tk*}$ = depth-averaged total-load sediment mass concentration for size class k and described in the Equilibrium Concentration and Transport Rates section [kg/m³]

$\beta _{tk}$ = total-load correction factor described in the Total-Load Cor-rection Factor section [-]

$r_{sk}$ = fraction of suspended load in total load for size class k and is described in Fraction of Suspended Sediments section [-]

$v_{s}$ = horizontal sediment mixing coefficient described in the Hori-zontal Sediment Mixing Coefficient section [m²/s]

$\alpha _{t}$ = total-load adaptation coefficient described in the Adaptation Coefficient section [-]

$\omega _{sk}$ = sediment fall velocity [m/s]

The above equation may be applied to single-sized sediment transport by using a single sediment size class (i.e. N=1). The bed composition, however, does not vary when using a single sediment size class. The sediment mass concentrations are used rather than volume concentrations in order to avoid precision errors at low concentrations.

Fraction of Suspended Sediments

In order to solve the system of equations for sediment transport implicitly, the fraction of suspended sediments must be determined explicitly. This is done by assuming

$r_{sk}={\frac {q_{sk}}{q_{tk}}}\Box {\frac {q_{sk*}}{q_{tk*}}}$ (2-48)

where $q_{sk}$ and $q_{tk}$ are the actual fraction of suspended- and total-load transport rates and $q_{sk*}$ and $q_{tk*}$ are the equilibrium fraction of suspended- and total-load transport rates.

Adaptation Coefficient

The total-load adaptation coefficient, $\alpha _{t}$ , is an important parameter in the sediment transport model. There are many variations of this parameter in literature (e.g. Lin 1984, Gallappatti and Vreugdenhil 1985, and Armanini and di Silvio 1986). CMS uses a total-load adaptation coefficient $\alpha _{t}$ that is related to the total-load adaptation length $L_{t}$ and time $T_{t}$ by

$L_{t}={\frac {Uh}{\alpha _{t}\omega _{s}}}=UT_{t}$ (2-49)

where

$\omega _{s}$ = sediment fall velocity corresponding to the transport grain size for single-sized sediment transport or the median grain size for multiple-sized sediment transport [m/s]

U = depth-averaged current velocity [m/s]

h = water depth [m]

The adaptation length (time) is a characteristic distance (time) for sedi-ment to adjust from non-equilibrium to equilibrium transport. Because the total load is a combination of the bed and suspended loads, the associated adaptation length may be calculated as $L_{t}=r_{s}L_{s}+(1-r_{s})L_{b}$ or $L_{t}=max(L_{s},L_{b})$ , where L_s and L_b are the suspended- and bed-load adaptation lengths. L_s is defined as

$L_{s}={\frac {Uh}{\alpha \omega _{s}}}=UT_{s}$ (2-50)

in which $\alpha$ and $T_{s}$ are the adaptation coefficient lengths for suspended load for the adaptation coefficient $\alpha$ can be calculated either empirically or based on analytical solutions to the pure vertical convection-diffusion equation of suspended sediment. One example of an empirical formula is that proposed by Lin (1984)

$\alpha =3.25+0.55ln\left({\frac {\omega _{s}}{\kappa u_{*}}}\right)$ (2-51)

where $u_{*}$ is the bed shear stress, and $\kappa$ is the von Karman constant. Armanini and di Silvio (1986) proposed an analytical equation

${\frac {1}{\alpha }}={\frac {\delta }{h}}+\left(1-{\frac {\partial }{h}}\right)exp\left[-1.5\left({\frac {\delta }{h}}\right)^{-1/6}{\frac {\omega _{s}}{u_{*}}}\right]$ (2-52)

where $\delta$ is the thickness of the bottom layer defined by $\delta =33z_{0}$ and $z_{0}$ is the zero-velocity distance from the bed. Gallapatti (1983) proposed the following equation to determine the suspended load adaptation time

$T_{s}={\frac {h}{u_{*}}}exp\left[{\begin{aligned}&(1.57-20.12u_{r})\omega _{*}^{3}+(326.832u_{r}^{2.2047}-0.2)\omega _{*}^{2}\\&+(0.1385lnu_{r}-6.4061)\omega _{*}+(0.5467u_{r}-2.1963)\end{aligned}}\right]$ (2-53)

where $u*_{c}$ is the current related bottom shear velocity, $u_{r}=u*/U$ , and $\omega _{*}=\omega _{s}/u_{*}$ .

The bed-load adaptation length, $L_{b}$ , is generally related to the dimension of bed forms such as sand dunes. Large bed forms are generally proportional to the water depth and therefore the bed load adaptation length can be estimated as $L_{b}=a_{b}h$ in which $a_{b}$ is an empirical coefficient on the order of 5-10. Fang (2003) found that $L_{b}$ of approximately two or three times the grid resolution works well for field applications. Although limited guidance exists on methods to estimate $L_{b}$, the determination of $L_{b}$ is still empirical and in the developmental stage. For a detailed discussion of the adaptation length, the reader is referred to Wu (2007). In general, it is recommended that the adaptation length be calibrated with field data in order to achieve the best and most reliable results.

Total-Load Correction Factor

The correction factor, $\beta _{tk}$ , accounts for the vertical distribution of the suspended sediment concentration and velocity profiles, as well as the fact that bed load travels a slower velocity than the depth-averaged current velocity (see Figure 2.3). By definition, $\beta _{tk}$ is the ratio of the depth-averaged total-load and flow velocities.

fig_2_3.png

Figure 2.3. Schematic of sediment and current vertical profiles.

In a combined bed load and suspended load model, the correction factor is given by

$\beta _{tk}={\frac {1}{r_{sk}/\beta _{sk}+(1-r_{sk})U/u_{bk}}}$(2-54)

where $u_{bk}$ is the bed load velocity and $\beta _{sk}$ is the suspended load correction factor and is defined as the ratio of the depth-averaged sediment and flow velocities. Since most sediment is transported near the bed, both the total and suspended load correction factors ($\beta _{tk}$ and $\beta _{sk}$ ) are usually less than 1 and typically in the range of 0.3 and 0.7, respectively. By assuming logarithmic current velocity and exponential suspended sediment concentration profiles, an explicit expression for the suspended load correction factor $\beta _{sk}$ may be obtained as (Sánchez and Wu 2011)

$\beta _{sk}={\frac {\int _{a}^{h}uc_{k}dz}{U\int _{a}^{h}c_{k}dz}}={\frac {E_{1}(\phi _{k}A)-E_{1}(\phi _{k})+ln(A/Z)e^{-\phi _{k}A}-ln(1/Z)e^{-\phi _{k}}}{e^{-\phi _{k}A}[ln(1/Z)-1][1-e^{-\phi _{k}(1-A)}]}}$ (2-55)

where $\phi _{k}=\omega _{sk}h/\epsilon ,A=a/h,z=z_{a}/h,$ in which $\omega _{sk}$ is the sediment fall velocity for size class k, $\epsilon$ is the vertical mixing coefficient, a is a reference height for the suspended load, h is the total water depth, $z_{a}$ is the apparent roughness length, and $E_{1}$ is the exponential integral. The equation can be further simplified by assuming that the reference height is proportional to the roughness height (e.g.$a=30z_{a}$ ), so that $\beta _{sk}=\beta _{sk}(z,\phi _{k})$. Figure 2.4 shows a comparison of the suspended load correction factor based on the logarithmic velocity with exponential and Rouse suspended sediment concentration profiles.

Figure 2.4. Suspend load correction factors based on the logarithmic velocity profile and (a) exponential and (b) Rouse suspended sediment profile. The Rouse number is $r=\omega _{s}/(\kappa u_{*})$.

The bed load velocity, $u_{bk}$ , is calculated using the van Rijn (1984a) formula with re-calibrated coefficients from Wu et al. (2006)

$u_{bk}=1.64\left({\frac {\tau _{b}^{'}}{\tau _{crk}}}-1\right)^{0.5}{\sqrt {(s-1)gd_{k}}}$ (2-56)

where s is the specific gravity, g is the gravitational constant, $d_{50}$ is the median grain size diameter, $\tau _{b}^{'}$ is the bed shear stress related to the grain roughness and is determined by $\tau _{b}^{'}=(n^{'}/n)^{3/2}\tau _{b}$ where $n^{'}=d_{50}^{1/6}/20$ is the Manning’s coefficient corresponding to the grain roughness and $\tau _{crk}$ is the critical bed shear stress.

Bed Change Equation

The fractional bed change is calculated as

${\frac {\rho _{s}(1-\rho _{m}^{'})}{f_{morph}}}\left({\frac {\partial z_{b}}{\partial t}}\right)_{k}=\alpha _{t}\omega _{sk}(C_{tk}-C_{tk*})+{\frac {\partial }{\partial x_{j}}}\left(D_{s}\mid q_{bk}\mid {\frac {\partial z_{b}}{\partial x_{j}}}\right)$ (2-57)

where

$z_{b}$ = bed elevation with respect to the vertical datum [m]

$p_{m}^{'}$ = bed porosity [-]

$f_{morph}$ = morphologic acceleration factor [-]

$\rho _{s}$ = sediment density [~2650 kg/m³ for quartz sediment]

$C_{tk}$ = depth-averaged total-load sediment mass concentration for size class k defined as $C_{tk}=q_{tk}/(Uh)$ in which $q_{tk}$ is the total-load mass transport [kg/m³]

$C_{tk*}$ = depth-averaged total-load equilibrium sediment mass concentration for size class k and described in the Equilibrium Concentration and Transport Rates section [kg/m³]

$\alpha _{t}$ = total-load adaptation coefficient described in the Adaptation Coefficient section [-]

$\omega _{sk}$ = sediment fall velocity [m/s]

$D_{s}$ = empirical bed-slope coefficient (constant) [-]

$q_{bk}=hUC_{tk}(1-r_{sk})$ is the bed load mass transport rate [kg/m/s]

The sediment density is required in the previous equation (2-51) since mass concentrations are used. For a detailed derivation of the above equation, the reader is referred to Sanchez and Wu (2011a). The total bed change is calculated as the sum of the bed change for all size classes

${\frac {\partial z_{b}}{\partial t}}=\sum _{k}\left({\frac {\partial z_{b}}{\partial t}}\right)_{k}$ (2-58)

The purpose of the morphologic acceleration factor $f_{morph}$ is to speed-up the bed change so that the simulation time $t_{sim}$ represents approximately the change that would occur in $t{morph}=f_{morph}t_{sim}$ . This factor should be used with caution and only for idealized cases or time periods which are periodic (mainly tidal). If time-varying winds or waves are important pro-cesses for driving sediment transport, then it is recommended to use re-duced or representative wind and wave conditions. Since the CMS runs relatively fast, it is generally recommended to not use the morphologic acceleration factor when validating the sediment transport model using hindcast measurements. If good initial and boundary conditions are available should be available and therefore it. The morphologic acceleration factor is useful however when simulating idealized cases or analyzing project alternatives.

Bed material sorting and layering

Bed sorting is the process in which the bed material changes composition (fraction of each grain size class). The bed is descritized into multiple lay-ers each with a uniform bed composition. The fraction of each size class is then calculated and stored in each layer. The sorting of sediments is then calculating using the mixing or active layer concept (Hirano 1971; Karim and Kennedy 1982; and Wu 1991). The active layer is the top layer of the bed which exchanges material directly with the sediment transport.

The temporal variation of the bed-material gradation in the first (mixing or active) layer is calculated as (Wu 2007)

${\frac {\partial (\delta _{1}p_{1k})}{\partial t}}=\left({\frac {\partial z_{b}}{\partial t}}\right)_{k}+p_{1k}^{*}\left({\frac {\partial \delta _{1}}{\partial t}}-{\frac {\partial z_{b}}{\partial t}}\right)$ (2-59)

where $\delta _{1}$ is the thickness of the first layer. $p_{1k}^{*}$ is equal to $p_{1k}$ for $\partial z_{b}/\partial t-\partial \delta _{1}/\partial t\geq 0$ , and equal to the bed material gradation in the second sediment layer for $\partial z_{b}/\partial t-\partial \delta _{1}/\partial t<0$ . The bed-material sorting in the se-cond layer is calculated as

${\frac {\partial (\delta _{2}p_{2k})}{\partial t}}=-p_{1k}^{*}\left({\frac {\partial \delta _{1}}{\partial t}}-{\frac {\partial z_{b}}{\partial t}}\right)$ (2-60)

where $\delta _{2}$ is the thickness of the second layer, and <\math>p_{2k}</math> is the fraction of the $k^{th}$ sediment size in the second layer. It is noted that there is no material exchanged between the sediment layers below the second layer.

The sediment transport, bed change, and bed gradation equations are solved simultaneously (coupled), but are decoupled from the flow calculation at a given time step. To illustrate the bed layering process, Figure 2.5 shows an example of the temporal evolution of 7 bed layers during erosional and depositional regimes.

Figure 2.5. Schematic showing an example bed layer evolution. Colors indicate layer number and not bed composition.

Mixing Layer Thickness

The mixing layer thickness is calculated as

$\delta _{1}=min[max(\delta _{min},2d_{50},0.5\Delta ),\delta _{max}]$

where $\Delta$ is the beform height, and $\delta _{min}$ and $\delta _{max}$ are user-specified mini-mum and maximum layer thicknesses, respectively. At the beginning of each time step, the mixing layer thickness is calculated. For cell with a hard (non-erodable) bottom, the mixing layer is calculated as

$\delta _{1,hb}=min(\delta _{1},z_{b}-z_{hb})$(2-62)

where $z_{hb}$ is the elevation of the hard bottom. A hard bottom is a nonerodable bed surface such as bed rock or a coastal structure.

Equilibrium Concentrations and Transport Rates

In order to close the system of equations describing the sediment transport, bed change, and bed sorting equations, the fractional equilibri-um depth-averaged total-load concentration $C_{tk*}$ must be estimated from an empirical formula. The depth-averaged equilibrium concentration is defined as

$C_{tk*}={\frac {q_{}tk*}{Uh}}$(2-63)

where $q_{tk*}$ is the total-load transport for the $k^{th}$ sediment size class esti-mated from an empirical formula.

For convenience, $C_{tk*}$ is written in general form as

$C_{tk*}=p_{1k}C_{tk}^{*}$(2-64)

where $p_{1k}$ is the fraction of the sediment size k in the first (top) bed layer and $C_{tk}^{*}$ is the potential equilibrium total-load concentration. The potential concentration $C_{tk}^{*}$ can be interpreted as the equilibrium concentration for uniform sediment of size $d_{k}$ . Equation 2-57 above is essential for the coupling of sediment transport, bed change, and bed sorting equations.

Equilibrium Transport and Concentration Formulas

Lund-CIRP

Camenen and Larson (2005, 2007, 2008) developed general sediment transport formulas for bed and suspended loads under combined waves and currents. These are referred to as the Lund-CIRP transport formulas. The general transport formulas can be used for both symmetric and asymmetric waves, but for simplicity the waves are assumed to be symmetric. The current-related bed- and suspended-load transport with wave stirring is given by

${\frac {q_{b}^{*}}{\sqrt {(s-1)gd_{50}^{3}}}}=f_{b}\rho _{s}12{\sqrt {\Theta _{c}}}\Theta _{cw,m}exp\left(-4.5{\frac {\Theta _{cr}}{\Theta _{cw}}}\right)$(2-65)

${\frac {q_{b}^{*}}{\sqrt {(s-1)gd_{50}^{3}}}}=f_{s}\rho _{s}C_{R}U{\frac {\epsilon }{\omega _{s}}}\left[1-exp\left(-{\frac {\omega _{s}h}{\epsilon }}\right)\right]$ (2-66)

where

<math> = potential equilibrium bed load transport [kg/m/s]

<math> = potential suspended load transport [kg/m/s]

<math> = median grain size [m]

<math> = sediment specific gravity or relative density [-]

<math> = gravitational constant (9.81 m/s2)

<math> = sediment density (~2650 kg/m3)

<math> = Shields parameters due to currents [-]

<math> = mean Shields parameters due to waves and currents [-]

<math> = maximum Shields parameters due to waves and currents [-]

<math> = critical Shields parameter [-]

<math> = vertical sediment diffusivity [m2/s]

<math> = reference bed concentration [kg/m3]

<math> = bed-load scaling factor (default 1.0) [-]

<math> = suspended-load scaling factor (default 1.0) [-]