Sediment Transport 1: Difference between revisions

From CIRPwiki
Jump to navigation Jump to search
No edit summary
(Undo revision 13527 by Rdchllbm (talk))
Tag: Undo
 
(51 intermediate revisions by 2 users not shown)
Line 1: Line 1:
'''Sediment Transport'''
= Sediment Transport Overview =
 
'''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: <math> q_{tk*} = q_{sk*} + q_{bk*}</math>.
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: <math> q_{tk*} = q_{sk*} + q_{bk*}</math>.
Line 15: Line 13:
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.  
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'''
=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.
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.
Line 44: Line 42:
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).
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'''
= 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.
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.
Line 85: Line 83:
::<math>c_{a*}</math> = equilibrium (capacity) sediment concentration at an elevation a above the bed [kg/m<sup>3</sup>]
::<math>c_{a*}</math> = equilibrium (capacity) sediment concentration at an elevation a above the bed [kg/m<sup>3</sup>]


Bed Change Equation
= Bed Change Equation =


The bed change is calculated as  
The bed change is calculated as  
Line 103: Line 101:
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).
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'''
= Nonequilibrium Total-Load Transport Model =


Total-load Transport Equation
Total-load Transport Equation
Line 137: Line 135:
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.  
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
=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  
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  
Line 146: Line 144:
where <math>q_{sk}</math>  and <math>q_{tk}</math>  are the actual fraction of suspended- and total-load transport rates and <math>q_{sk*}</math>  and <math>q_{tk*} </math>  are the equilibrium fraction of suspended- and total-load transport rates.
where <math>q_{sk}</math>  and <math>q_{tk}</math>  are the actual fraction of suspended- and total-load transport rates and <math>q_{sk*}</math>  and <math>q_{tk*} </math>  are the equilibrium fraction of suspended- and total-load transport rates.


Adaptation Coefficient
= Adaptation Coefficient =


The total-load adaptation coefficient, <math>\alpha_t  </math>  , 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 <math>\alpha_t  </math>  that is related to the total-load adaptation length <math>L_t</math>  and time <math>T_t</math> by
The total-load adaptation coefficient, <math>\alpha_t  </math>  , 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 <math>\alpha_t  </math>  that is related to the total-load adaptation length <math>L_t</math>  and time <math>T_t</math> by
Line 195: Line 193:
The bed-load adaptation length, <math>L_b</math>  , 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 <math>L_b = a_b h</math>  in which <math>a_b</math>  is an empirical coefficient on the order of 5-10. Fang (2003) found that <math>L_b</math>  of approximately two or three times the grid resolution works well for field applications.  Although limited guidance exists on methods to estimate <math>L_b</math>, the determination of <math>L_b</math> 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.
The bed-load adaptation length, <math>L_b</math>  , 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 <math>L_b = a_b h</math>  in which <math>a_b</math>  is an empirical coefficient on the order of 5-10. Fang (2003) found that <math>L_b</math>  of approximately two or three times the grid resolution works well for field applications.  Although limited guidance exists on methods to estimate <math>L_b</math>, the determination of <math>L_b</math> 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
= Total-Load Correction Factor =


The correction factor, <math>\beta_{tk}</math> , 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, <math>\beta_{tk}</math> is the ratio of the depth-averaged total-load and flow velocities.
The correction factor, <math>\beta_{tk}</math> , 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, <math>\beta_{tk}</math> is the ratio of the depth-averaged total-load and flow velocities.
Line 225: Line 223:


where s is the specific gravity, g  is the gravitational constant, <math>d_{50}</math>  is the median grain size diameter, <math>\tau_b^'</math>  is the bed shear stress related to the grain roughness and is determined by <math>\tau_b^' = (n^' / n)^{3/2}\tau_b</math>  where <math>n^' = d_{50}^{1/6} / 20</math>  is the Manning’s coefficient corresponding to the grain roughness and <math>\tau_{crk}</math>  is the critical bed shear stress.
where s is the specific gravity, g  is the gravitational constant, <math>d_{50}</math>  is the median grain size diameter, <math>\tau_b^'</math>  is the bed shear stress related to the grain roughness and is determined by <math>\tau_b^' = (n^' / n)^{3/2}\tau_b</math>  where <math>n^' = d_{50}^{1/6} / 20</math>  is the Manning’s coefficient corresponding to the grain roughness and <math>\tau_{crk}</math>  is the critical bed shear stress.
= Bed Change Equation =
The fractional bed change is calculated as
::<math>
\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)
</math> (2-57)
where
::<math>z_b</math>  = bed elevation with respect to the vertical datum [m]
::<math>p_m^'</math> = bed porosity [-]
::<math>f_{morph}</math>  = morphologic acceleration factor [-]
::<math>\rho_s</math>  = sediment density [~2650 kg/m<sup>3</sup> for quartz sediment]
::<math>C_{tk}</math> = depth-averaged total-load sediment mass concentration for size class k defined as <math>C_{tk} = q_{tk} /(Uh)</math> in which <math>q_{tk}</math>  is the total-load mass transport [kg/m<sup>3</sup>]
::<math>C_{tk*}</math> = depth-averaged total-load equilibrium sediment mass concentration for size class k and described in the Equilibrium Concentration and Transport Rates section [kg/m<sup>3</sup>]
::<math>\alpha_t</math>  = total-load adaptation coefficient described in the Adaptation Coefficient section [-]
::<math>\omega_{sk}</math>  = sediment fall velocity [m/s]
::<math>D_s</math> = empirical bed-slope coefficient (constant) [-]
::<math>q_{bk} = hUC_{tk}(1 - r_{sk})</math> 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
::<math>
\frac {\partial z_b}{\partial t} = \sum_k \left(\frac {\partial z_b}{\partial t}  \right)_k
</math> (2-58)
The purpose of the morphologic acceleration factor <math>f_{morph}</math>  is to speed-up the bed change so that the simulation time <math>t_{sim}</math>  represents approximately the change that would occur in <math>t{morph} = f_{morph}t_{sim} </math> . 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)
::<math>
\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)
</math>  (2-59)
where <math>\delta_1</math>  is the thickness of the first layer. <math>p_{1k}^*</math>  is equal to <math>p_{1k}</math>  for <math>\partial z_b / \partial t - \partial \delta_1 / \partial t \geq 0</math> , and equal to the bed material gradation in the second sediment layer for <math>\partial z_b / \partial t - \partial \delta_1 / \partial t < 0 </math> . The bed-material sorting in the se-cond layer is calculated as
::<math>
\frac {\partial (\delta_2 p_{2k})}  {\partial t} = - p_{1k}^*
\left(\frac {\partial \delta_1}        {\partial t} - \frac {\partial z_b} {\partial t} \right)
</math> (2-60)
where <math>\delta_2</math>  is the thickness of the second layer, and <\math>p_{2k}</math>  is the fraction of the <math>k^{th}</math>  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
::<math>
\delta_1 = min [max(\delta_{min},2d_{50},0.5\Delta),\delta_{max}]
</math>
where <math>\Delta</math>  is the beform height, and <math>\delta_{min}</math>  and <math>\delta_{max}</math>  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
::<math>
\delta_{1,hb} = min(\delta_1 , z_b - z_{hb})
</math>(2-62)
where <math>z_{hb} </math> 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 <math>C_{tk*}</math>  must be estimated from an empirical formula. The depth-averaged equilibrium concentration is defined as
::<math>
C_{tk*} = \frac {q_{}tk*}{Uh}
</math>(2-63)
where <math>q_{tk*}</math>  is the total-load transport for the <math>k^{th} </math> sediment size class esti-mated from an empirical formula.
For convenience, <math> C_{tk*} </math> is written in general form as
::<math>
C_{tk*} = p_{1k} C_{tk}^*
</math>(2-64)
where <math>p_{1k}</math>  is the fraction of the sediment size k in the first (top) bed layer and <math>C_{tk}^*</math>  is the potential equilibrium total-load concentration. The potential concentration  <math> C_{tk}^*</math> can be interpreted as the equilibrium concentration for uniform sediment of size <math>d_k</math> . 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
::<math>
\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)
</math>(2-65)
::<math>
\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]
</math> (2-66)
where
::<math>q_b^*</math>  = potential equilibrium bed load transport [kg/m/s]
::<math>q_s^*</math>  = potential suspended load transport [kg/m/s]
::<math>d_{50}</math>  = median grain size [m]
::<math> s </math> = sediment specific gravity or relative density [-]
::<math>g</math>  = gravitational constant (9.81 m/s<sup>2</sup>)
::<math>\rho_s</math>  = sediment density (~2650 kg/m<sup>3</sup>)
::<math>\Theta_C</math>  = Shields parameters due to currents [-]
::<math>\Theta_{cw,m}</math>  = mean Shields parameters due to waves and currents [-]
::<math>\Theta_cw</math>  = maximum Shields parameters due to waves and currents [-]
::<math>\Theta_cr</math>  = critical Shields parameter [-]
::<math>\epsilon</math>  = vertical sediment diffusivity [m<sup>2</sup>/s]
::<math>c_R</math>  = reference bed concentration [kg/m<sup>3</sup>]
::<math>f_b</math>  = bed-load scaling factor (default 1.0) [-]
::<math>f_s</math>  = suspended-load scaling factor (default 1.0) [-]
The reference concentration is given by
::<math>c_R = \rho_s A_{cR}\Theta_{CW,m}exp\left(-4.5\frac{\Theta_{cr}}{\Theta_{cw}}  \right)</math> (2-67)
where the coefficient <math>A_{cR}</math>  is determined by the following relationship
::<math>
A_{cR} = 0.0035 exp (-0.3d_*)
</math> (2-68)
The current-related shear stress is calculated as
::<math>
\tau_c = \rho c_b U^2
</math> (2-69)
where <math>\rho</math>  is the water density,<math>c_b</math>  is the bed friction coefficient, and <math>U</math>  is the current velocity magnitude. The drag coefficient is calculated as
::<math>
c_b = \left [\frac{\kappa}{ln(h/z_0 )-1}    \right]^2
</math>(2-70)
where <math>\kappa</math>  is the von Karman constant (0.4), <math>h</math>  is the total water depth, and <math>z_0</math>  is the roughness length calculated as <math>z_0 = k_s / 30</math>  , where <math>k_s</math>  is the total bed Nikuradse roughness. The total bed roughness is assumed to be a linear sum of the grain-related roughness <math>k_{s,d}</math>, form-drag (ripple) roughness <math>k_{s,r}</math>, and sediment-related roughness <math>k_{s,s}</math> . Bed forms are also separated into current and wave-related bed forms. The current- and wave-related total roughness is then
::<math>
k_{s,c\mid w} = k_{sg} + k_{sr,c \mid w} + k_{ss,c \mid w}
</math>(2-71)
where the subscript c|w indicates either the current (c) or wave (w) related component. The grain-related roughness is estimated as <math>k_{sg} = 2d_{50}</math>.
The ripple roughness is calculated as (Soulsby 1997)
::<math>
k_{sr,c \mid w} = 7.5 \frac{H_{r,c \mid w}^2 }{L_{r,c \mid w}}
</math> (2-72)
where <math>H_{r,c \mid w}</math>  and <math>L_{r,c \mid w}</math>  are either the current- or wave-related ripple height and length respectively. The current-related ripple height and length are calculated as
::<math>
L_{r,c} = 1000d_{50}
</math> (2-73)
::<math>
H_{r,c} = L_{r,c}/7
</math>(2-74)
The wave-related ripple height and length are calculated using the expressions proposed by van Rijn (1984b, 1989)
::<math>
H_{r,c} = \left\{
\begin{align}
&0.22 A_w                                  &\text{for} \psi_w < 10 \\
&2.8 \times 10^{-13}(250 - \psi_w)^5 A_w    &\text{for} 10\leq \psi_w < 250 \\
&0                                        &\text{for} 250 \leq \psi_w
\end{align}
\right.
</math>(2-75)
::<math>
L_{r,w} = \left\{
\begin{align}
&1.25 A_w                                  &\text{for} \psi_w < 10 \\
&1.4 \times 10^{-6}(250 - \psi_w)^{2.5} A_w &\text{for} 10 \leq \psi_w < 250 \\
&0 \quad                                    &\text{for} 250 \leq \psi_w
\end{align}
\right.
</math>(2-76)
where <math>A_W</math>  is the semi-orbital excursion and <math>\psi_W</math>  is the wave mobility parameter. The semi-orbital excursion is defined as
::<math>
A_W = \frac {u_w T}{2\pi}
</math><2-77)
in which <math>u_w</math>  is the peak bottom orbital velocity and <math>T</math>  is the wave period. For random waves, <math> u_w = \sqrt u_{rms}</math>  and <math> T = T_p</math>  are used as representative values. The wave mobility parameter,<math> \psi_w</math>  , is defined as
::<math>
\psi_W = \frac {u_W^2}{(s-1)gd_{50}}
</math>(2-78)
For the Lund-CIRP sediment transport equations, a wave bottom shear stress is calculated as
::<math>
\tau_W = \frac {1}{2} \rho f_w u_w^2
</math>(2-79)
where <math>f_w</math>  is the wave friction factor calculated using the expression of Swart (1976)
::<math>
f_w = \left\{
\begin{align}
&exp(5.21 r^{-0.19} - 6.0) \text{for} r > 1.57 \\
&0.3 \quad \quad\quad \quad\quad \quad\quad \quad \text{for} r\leq 1.57
\end{align}
\right.
</math>(2-80)
where <math>r</math>  is the relative roughness defined as <math>r = A_W/k_{sg|}</math>.
The current- and wave-related sediment roughnesses are estimated as
::<math>
k_{ss,c\mid w} = 5d_{50}\Theta_{c\mid w}
</math>(2-81)
Equation 2-73 above must be solved simultaneously with the expressions for the bottom shear stress because the roughness depends on the stress. The exact solution is approximated using explicit polynomial fits in order to avoid time-consuming iterations in calculating the bed shear stress.
The critical Shields parameter is estimated using the formula proposed by Soulsby (1997)
::<math>
\Theta_{cr} = \frac {0.3}{1+1.2 d_*} + 0.055[1-exp(-0.02d_*)]
</math> (2-82)
where <math>d_*</math>  is the dimensionless grain size
::<math>
d_* = d_{50}\left[\frac {(s-1)g}{v_2}\right]^{1/2}
</math>(2-83)
where <math> v</math>  is the kinematic viscosity.
The sediment fall velocity is calculated using the formula by Soulsby (1997)
::<math>
\omega_s = \frac {v}{d}\left[(10.36^2 + 1.049d_*^3)^{1/2} - 10.36    \right]
</math>(2-84)
where <math>d</math>  is the grain size. The vertical sediment diffusivity is calculated as
::<math>
\epsilon = h \left(\frac {D_e}{\rho} \right)^{1/3}
</math>(2-85)
where <math>D_e</math>  is the total effective dissipation given by
::<math>
D_e = k_b^3 D_b + k_c^3 D_c + k_w^3 D_w
</math>(2-86)
in which <math>k_b, k_c,  \text{and} k_w </math>are coefficients, <math>D_b</math> is the wave breaking dissipation, and <math>D_c \text{and} D_w</math> are the bottom friction dissipation due to currents and waves, respectively. The dissipation from bottom friction due to current, <math>D_c</math>, and the dissipation from bottom friction due to waves, <math>D_w</math> , are expressed as
::<math>
D_{c\mid w} = \tau_{c\mid w}u*_{c\mid w}
</math>(2-87)
where the subscript c|w indicates either the current (c) or wave (w) relat-ed component, and <math>u*_{c\mid w}</math>  and <math>\tau_{c\mid w} </math>  are the current- or wave-related bed shear velocity and stress, respectively. The coefficient <math>k_b</math>  =0.017, and <math>k_c \text{and} k_w</math> are function of the Schmidt number:
::<math>
k_{c\mid w} = \frac{\kappa}{6} \sigma_{c\mid w}
</math>(2-88)
where <math>\sigma_{c\mid w}</math>  is either the current or wave-related Schmidt number calculated from the following relationships (Camenen and Larson 2008):
::<math>
\sigma_c = \left\{
\begin{align}
&0.7 + 3.6 sin^{2.5} \left(\frac {\pi \omega_s}{2u*_c}    \right) & \text{for} \frac {\omega_s}{u*_c} \leq 1 \\
&1 + 3.3 sin^{2.5} \left(\frac{\pi u*_c}{2\omega_s}      \right) & \text{for} \frac {\omega_s}{u*_c} > 1
\end{align}
\right.
</math>
(2-89)
::<math>
\sigma_w = \left\{
\begin{align}
&0.09 + 1.4 sin^{2.5}\left(\frac{\pi\omega_s}{2u*_w}\right) & \text{for} \frac{\omega_s}{u*_w}\leq 1 \\
&1 + 0.49 sin^{2.5}\left(\frac{\pi u*_w}{2 \omega_s}\right) 
& \text{for} \frac{\omega_s}{u*_w}>1
\end{align}
\right.
</math>
(2-90)
In the case of coexisting or combined waves and currents, the wave-current Schmidt number is estimated as
::<math>
\sigma_{cw} = X^5_v \sigma_c + (1 - X^5_v)\sigma_w
</math>
(2-91)
where <math>X_v</math>  is a weighting factor equal to <math>X_v = U / (U + u_w)</math> .
For multiple-sized (nonuniform) sediments, the fractional equilibrium sediment transport rates are calculated as
::<math>
\frac{q_{bk*}}{\sqrt{(s-1)gd_k^3}} = f_b \zeta_k^{\zeta-1}p_{bk}\rho_s 12 \sqrt{\Theta_c}\Theta_{cw,m}exp \left(-4.5 \frac{\Theta_{crk}}{\Theta_{cw}}  \right)
</math>
(2-92)
::<math>
\frac{q_{bk*}}{\sqrt{(s-1)gd_k^3}} = f_s \zeta_k^{\zeta-1}p{bk}\rho_s c_{Rk}U\frac{\epsilon_k}{\omega_{sk}}\left[1 - exp \left( -\frac{\omega_{sk}h}{\epsilon}  \right)  \right]
</math>
(2-93)
where the subscript k indicates variables which are calculated based only on the sediment size class k. <math>\zeta_k</math>  is the hiding and exposure coefficient de-scribed in Hiding and Exposure.
==Van Rijn==
The van Rijn (1984a,b) sediment transport equations for bed load and suspended load are used with the recalibrated coefficients of van Rijn (2007a,b), as given by
::<math>
q_b^* = f_b \rho_s 0.015Uh
\left(\frac{U_e - U_{cr}}{\sqrt{(s-1)gd_{50}}}    \right)^{1.5} \left( \frac{d_{50}}{h} \right)^{1.2}
</math>(2-94)
::<math>
q_s^* = f_s \rho_s 0.012Ud_{50} \left(\frac{U_e - U_{cr}}{\sqrt{(s-1)gd_{50}}}  \right)^{2.4}d_*^{-0.6}
</math> (2-95)
where
::<math>q_b^*</math> = equilibrium bed-load transport [kg/m/s]
::<math>q_s^*</math>  = equilibrium suspended-load transport [kg/m/s]
::<math>d_{50}</math>  = median grain size [m]
:: s = sediment specific gravity or relative density [-]
:: g  = gravitational constant (9.81 m/<sup>2</sup>)
::<math>\rho_s</math>  = sediment density (~2650 kg/m<sup>3</sup>)
:: U = depth-averaged current velocity [m/s]
::<math>U_{cr}</math>  = critical depth-averaged velocity for incipiant motion [m/s]
::<math>U_e</math>  = effective depth averaged velocity [m/s]
::<math>f_s</math>  = suspended-load scaling factor (default 1.0) [-]
::<math>f_b</math>  = bed-load scaling factor (default 1.0) [-]
The effective depth-averaged velocity is calculated as <math>U_e = U + \gamma u_w</math>  with <math>\gamma</math> =0.4 for random waves and <math>\gamma</math> =0.8 for regular waves. Here <math>u_w</math>  is the bottom wave orbital velocity based on linear wave theory and the signifi-cant wave height The critical depth-averaged velocity is estimated as <math>U_{cr} = \beta U_{crc} + (1 - \beta)u_{crw}</math> where <math>\beta</math>  is a blending factor and <math>U_{crc}</math>  and <math>u_{crw}</math> are the critical velocities for currents and waves, respectively. As described in van Rijn (2007), the critical velocity for currents is
::<math>
U_{cr} =
\left\{
\begin{align}
&0.19 d_{50}^{0.1} log_{10} \left(\frac{4h}{d_{90}}\right), &for 0.1 \leq d_{50} \leq 0.5 mm\\
&8.5 d_{50}^{0.6} log_{10} \left(\frac{4h}{d_{90}}\right),
&for 0.5 \leq d_{50} \leq 2.0 mm
\end{align}
\right.
</math>(2-96)
where <math>d_{90}</math> is the sediment grain size in meters of 90<sup>th</sup> percentile. The critical velocity for waves is based on Komar and Miller (1975),
::<math>
U_{crw} =
\left\{
\begin{align}
&0.24[(s-1)g]^{0.66} d_{50}^{0.33}T_p^{0.33}, &for 0.1 \leq d_{50} \leq 0.5 mm \\
&0.95[(s-1)g]^{0.57} d_{50}^{0.43}T_p^{0.14}, &for 0.5 \leq d_{50} \leq 2.0 mm
\end{align}
\right.
</math>(2-97)
where <math>T_p</math> is the peak wave period. According to van Rijn (2007), the bed load transport formula predicts transport rates by a factor of 2 for velocities higher than 0.6 m/s, but underpredicts transport rates by a factor of 2-3 for velocities close to the initiation of motion.
The near-bed equilibrium sediment concentration for the van Rijn (formula is given by (van Rijn 1985)
::<math>
c_{a*} = f_s \rho_s 0.015 \frac {d_{50}}{a} \left(\frac {\tau^'_{max} - \tau_{cr}} {\tau_{cr}}  \right)^{1.5} d_*^{-0.3}
</math>(2-98)
where <math>r^1_{max}</math>  is the maximum skin shear stress, <math>r_{cr}</math>  is the critical shear stress, and <math>\alpha</math>  is reference height above the bed given by
::<math>
\alpha = max(0.5H_r, 0.01h)
</math>(2-99)
where <math>H_r</math>  is the bed ripple height.
The van Rijn formula (1984a,b, 2007a,b) were originally proposed for well-sorted sediments. The sediment availability is included by multiplication of transport rateswith the fraction of the sediment size class in the upper bed layer. The hiding and exposure is considered by in a correction factor which multiples to the critical velocity. When applied to multiple-sized sediments, the fractional equilibrium transport rate is calculated as
::<math>
q_{bk*} = f_b \rho_s p _{1k}0.015Uh \left(\frac{U_e - \zeta_k^{1/2}U_{crk}}{\sqrt{(s-1)gd_k}} \right)^{1.5} \left(\frac {d_k}{h} \right)^{1.2}
</math>(2-100)
::<math>
q_{sk*} = f_s \rho_s p_{1k}0.012Uh \left(\frac{U_e - \zeta_k^{1/2}U_{crk}}{\sqrt{(s-1)gd_k}}  \right)^{2.4} \left(\frac{d_k}{h}  \right) d_{*k}^{-0.6}
</math> (2-101)
where <math>p_{bk}</math>  is the fractional bed composition and <math>\zeta_k</math>  is the hiding and ex-posure coefficient. The subscript k indicates values which are calculated based on the k<sup>th</sup> sediment size class.
==Soulsby-van Rijn==
Soulsby (1997) proposed the following equation for the total load sediment transport rate under action of combined current and waves,
::<math>
q_{b*} = f_b \rho_s 0.005Uh \left(\frac{U_e - U_{cr}}{\sqrt{(s-1)gd_{50}}}\right)^{2.4} \left(\frac{d_{50}}{h} \right)^{1.2}</math>(2-102)
::<math>
q_{s*} = f_b \rho_s 0.012Uh \left(\frac{U_e - U_{cr}}{\sqrt{(s-1)gd_{50}}}\right)^{2.4} \left(\frac{d_{50}}{h} \right)d_*^{-0.6}</math>(2-103)
where
::<math>q_{t*}</math> = equilibrium total-load transport [kg/m/s]
::<math>\rho_s</math>  = sediment density (~2650 kg/m<sup>3</sup>)
:: U = depth-averaged current velocity [m/s]
::<math>U_e</math>  = effective velocity defined as <math>= \sqrt{U^2 + \frac{0.018}{C_d}u^2_{rms}}</math>  [m/s]
::<math>u{rms}</math>  = peak bottom wave orbital velocity based on the root-mean-squared wave height [m/s]
::<math>C_d</math> = drag coefficient due to currents alone [-]
::<math>U_{cr}</math>  = critical depth-averaged velocity for initiation of motion [m/s]
::<math>\alpha_{sv}</math>  = empirical bed-slope coefficient, approximately equal to 1.6
::<math>m_s</math>  = bed-slope [-]
The coefficients A<sub>sb</sub> and A<sub>ss</sub> are related to the bed and suspended loads, respectively, and given by
::<math>
A_{sb} = f_b \frac{0.005h(d_{50}/h)^{1.2}} {[(s-1)gd_{50}]^{1.2}}, A_{ss} = f_s \frac{0.012d_{50}d_*^{-0.6}} {[(s-1)gd_{50}]^{1.2}}
</math>(2-104)
Note that the bed and suspended load scaling factors (<math>f_b</math>  and <math>f_s</math>  ) are in-cluded in the calculations for the coefficients <math>A_{sb}</math> and <math>A_{ss}</math>.
The current drag coefficient is calculated as
::<math>
C_d = \left[\frac{\kappa}{ln(h/z_0) - 1}    \right]^2
</math>(2-105)
where <math>z_0</math>  is the bed roughness length (equal to 0.006 m).
The Soulsby-van Rijn formula is modified for multiple-sized sediments similarly to the van Rijn formula in the previous section with the equation
::<math>
q_{bk*} = f_b \rho_s p_{1k}0.005Uh \left(\frac{U_e - \zeta_k^{1/2}U_{crk} } {\sqrt{(s-1)gd_k}}  \right)^{2.4} \left(\frac{d_k}{h}  \right)^{1.2} (1 - \alpha_{sv}tan \beta)
</math>(2-106)
::<math>
q_{sk*} = f_b \rho_s p_{1k}0.012Uh \left(\frac{U_e - \zeta_k^{1/2}U_{crk} } {\sqrt{(s-1)gd_k}}  \right)^{2.4} \left(\frac{d_k}{h}\right) d_{*k}^{-0.6} \left(1 - \frac{e_s}{\omega_s}tan \beta  \right)
</math>(2-107)
where
::<math>
U_e = \sqrt {U^2 + \frac{0.018}{C_d} u^2_{rms}}
</math>
The subscript k indicates that the value is calculated based only on the size class k and not the median grain size. The availability of sediment fractions is included through <math>p_{1k}</math>, while hiding and exposure of grain sizes is accounted for by modifying the critical velocity. it is noted that the Soulsby-van Rijn (1997) formula are very similar to the van Rijn (1984ab, 207ab) except for the definition of the effective velocity and the recalibration of the bed-load formula coefficients in 2007a. The proposed changes should be verified with measurements or numerical simulations for multiple-sized sediment transport.
==Watanabe==
The equilibrium total load sediment transport rate  is determined by Watanabe (1987) as
::<math>
q_{t*} = A_{Wat}U \left(\frac{\tau_{bmax} - \tau_{cr}}  {\rho g}  \right)</math> (2-108)
where
::<math>q_{t*}</math>  = equilibrium total-load transport [kg/m/s]
::<math>\tau_{bmax}</math>  = combined wave-current maximum shear stress [Pa]
::<math>\rho</math>  = water density (~1025 kg/m<sup>3</sup>)
:: g = gravitational constant (9.81 m/s<sup>2</sup>)
::<math>\tau_{cr}</math>  = critical shear stress of incipient motion in [Pa]
::<math>A_{Wat}</math>  = empirical coefficient typically ranging from 0.1 to 2
The critical shear stress is determined from the Shields diagram. The maximum bed shear stress <math>\tau_{bmax}</math> is calculated as (Soulsby 1997)
::<math>
\tau_{bmax} = \sqrt{(\tau_b + \tau_w cos \varphi )^2 +(r_w sin \varphi)^2                    }
</math>(2-109)
where <math>\tau_b</math>  is the mean shear stress by waves and current over a wave cycle, <math>\tau_w</math>  is the mean wave bed shear stress, and <math>\varphi</math>  is the angle between the waves and current. The wave bed shear stress is calculated as
::<math>
\tau_w = \frac {1}{2}\rho f_w u_w^2
</math>(2-110)
where <math>f_w</math>  is a wave friction factor and <math>u_w</math>  is the bottom wave orbital velocity amplitude. The wave friction factor is calculated with the expression by Nielson (1992)
::<math>
f_w = exp(5.5 r^{-0.2} - 6.3)
</math>(2-111)
where r is the relative roughness defined as <math>r = A_w / k_sg</math>. <math>A_w</math>  is the semi-orbital excursion defined as <math>A_w = u_w T /(2\pi)</math>.
'''Hiding and Exposure'''
When the bed material is composed of multiple grain sizes, larger grains have a greater probability of being exposed to the flow while smaller particles have a greater probability of being hidden to the flow. Figure 2.6 shows an example of a sediment grain <math>d_k</math>  being exposed to the flow by an exposure height <math>\Delta_e</math> , and sediment grain <math> d_j</math>  being hidden by <math>d_k</math>.
Figure 2.6. Schematic of the exposure height of bed sediment grains
For the van Rijn (Reference(s)), Soulsby-van Rijn (Reference) and Watanabe (Reference) transport formulas, the hiding and exposure mechanism is considered by correcting the critical shear stress or velocity using a hiding and exposure correction function,<math>\zeta_k</math>  . For the Lund-CIRP transport formula, an alternate approach is required due to the way in which the Shields number and grain size are included in the formulation; thus, the hiding and exposure correction function is directly used to multiply the transport rate. Two methods are used to calculate <math>\zeta_k</math>  , depending on whether the sediment transport model is run with a single sediment size or with multiple sediment sizes and the methods are described in the following sections.
Single-sized sediment transport
In some applications, the coastal bed material is dominated by a single sediment size with patches of other sediment sizes or materials (e.g. shell hash) that may not contribute significantly to morphology change in the areas of interest; however, they may modify the sediment transport through hiding and exposure. For example, it is possible for the bed material to consist of mostly uniform sand with patches of shell fragments (bimodal distribution) in some regions. The shell material is difficult to model numerically because it is usually poorly sorted and its hydraulic properties are unknown. For such regions, sediment transport models commonly estimate excessive erosion due to the absence of the hiding effect of the coarser shell material (e.g. Cayocca 2001). A better and more physical plausible approach is to use the local bed composition along with a correction function for hiding and exposure to account for the effects of the uniform sand with the patches of coarser shell material. For single-sized sediment transport, the correction function for hiding and exposure is calculated following Parker et al. (1982) and others as,
::<math>
\zeta_k = \left(\frac {d_{50}}{d_k} \right)
</math>(2-112)
where <i>i</i>  is an empirical coefficient between 0.5-1.0. The aforementioned sediment transport equations are implemented by using the transport grain size <math>d_k</math>  rather than the bed material <math>d_{50}</math> . A single and constant transport size <math>d_k</math>  is used, while the bed material <math>d_{50}</math>  varies spatially. The spatial distribution of  can be obtained from field measurement data and for simplicity is assumed constant during the model simulation time. This is a significant assumption and may not be reasonable for some applications. However, this method provides a simple conceptual mechanism for considering an important process in the proposed single-sized sediment transport model. The approach has been successfully applied to Shinnecock Inlet, NY to simulate morphology change at a coastal inlet (Sánchez and Wu 2011a). A more accurate and complex approach is to simulate the transport and sorting of multiple sediments.
Multiple-sized sediment transport
The hiding and exposure of the each sediment size class is considered by modifying the critical shields parameter <math>\Theta_{ck}</math>  for each sediment size class based on Wu et al. 2000
::<math>
\zeta_k = \left (\frac{P_{ek}}{P_{hk}}    \right)^{-m}
</math>(2-113)
where ''m''  is an empirical coefficient that varies for each transport formula, approximately equal to 0.6-1.0. <math>P_{ek}</math>  and <math>P_{hk}</math>  are the total hiding and exposure probabilities and are calculated as
::<math>
P_{hk} = \sum^N_{j=1}p_{1j}\frac {d_j}{l_k + d_j} \quad P_{ek} = \sum^N_{j=1}p_{1j}\frac{d_k}{d_k + d_j}
</math> (2-114)
where N is the number of grain size classes.
'''Horizontal Sediment Mixing Coefficient'''
The horizontal sediment mixing coefficient, v<sub>s</sub>, represents the combined effects of turbulent diffusion and dispersion due to nonuniform vertical profiles. In CMS, the horizontal sediment mixing coefficient is assumed to be proportional to the turbulent eddy viscosity as
::<math>
v_s = v_t / \sigma_s
</math>(2-115)
where <math>\sigma_s</math>  is the Schmidt number and <math>v_t</math> is…. There are many formulas to estimate the Schmidt number. The Schmidt number is set by default equal to 1.0 in CMS, but may be modified by the user with cards.
'''Boundary Conditions'''
At the interface between wet and dry cells, the sediment transport rate is set to zero. The inflow boundary condition requires a given sediment concentration at the boundary. However, for most coastal applications, the actual sediment concentration is not available and the model implements the equilibrium concentration capacity (Dirichlet boundary condition). CMS also requires the size distribution of the inflow transport. For stability reasons, the inflow size distribution is assumed to be equal to the initial size distribution of the bed at the boundary. Inflow sediment transport rates are specified either as a total sediment transport rate, <math>Q_{sed,tot}</math>  in kg/s, or as a fractional sediment transport rate,<math>Q_{sed,k}</math>  in kg/s, for each sediment size class. If only the total sediment transport rate is specified, then the fractional sediment transport rate is calculated as <math>Q{sed,k}  = p_{1k}Q_{sed.tot}</math> where <math>p_{1k}</math>  is the bed material fraction of the first layer at the boundary. The fractional sediment transport <math>q_{tk,B} </math> at boundary cell ''B'' in kg/m/s is then calculated along the cell string according to
::<math>
q_{tk,B} = \frac {q^*_{tk,B}} {\sum \Delta l_f Q^*_{tk,N}} Q_{sed,k}
</math>(2-116)
where <math>q^*_{tk,B}</math>  is the potential sediment transport rate at boundary cell ''B'' and <math>\Delta l_f</math>  is the inflow cell face width. If the flow is directed outward of the do-main, a zero-gradient boundary condition is used for sediment concentration (Neumann boundary condition).
'''Salinity Transport
:Transport Equation'''
CMS calculates the salinity transport based on the following 2DH salinity conservation equation
::<math>
\frac {\partial(hC_{sal})} {\partial t} + \frac {\partial(hV_j C_{sal})} {\partial x_j} = \frac{\partial}{\partial x_f} \left(v_{sal}h \frac{\partial C_{sal}} {\partial x_j}  \right) + S_{sal}
</math>(2-117)
where
::<math>C_{sal}</math>  = depth-averaged salinity [ppt or parts per thousand]
::''h''  = water depth [m]
::<math>V_j</math>  = total flux velocity [m/s]
::<math>v_{sal}</math>  = horizontal mixing coefficient <math>v_{sal} = v_t / \sigma_s \quad  [m^2/s]</math>
::<math>S_{sal} </math> = source/sink term due to precipitation, evaporation, and structures (e.g. culverts) [ppt m/s]
Equation (2-106) represents the horizontal fluxes of salt in water bodies and is balanced by exchanges of salt via diffusive fluxes. Major processes that contribute to the salinity are: freshwater inflows from rivers, vertical fluxes of freshwater by precipitation and evaporation at the water surface, and groundwater fluxes.  These processes can be specified as the surface and bottom boundary conditions in the equation.
'''Boundary Conditions'''
The boundary condition for the salinity transport equation is dependent on the flow direction. At cell faces between wet and dry cells, a zero-flux boundary condition is applied. The cell face velocity <math>U_f</math>  is zero between wet and dry cells and this eliminates advection transport. In addition, the salinity horizontal mixing coefficient at the cell face, <math>\bar v_{sal,f}</math>  , is set to zero and this eliminates the diffusive flux. If the flow is directed inward of the modeling domain, then a user-specified salinity concentration is pre-scribed (Direchlet boundary condition). If the flow is directed outward of the modeling domain, then a zero-gradient boundary condition is applied (Neumann boundary condition).  Overall comments for this section: 1) A figure would be very HELPFUL showing inward and outward directed flow with the modeling domain and 2) Direchlet and Neumann BCs should be explained with reference citations.
'''Surface Roller'''
As a wave transitions from non-breaking to fully-breaking, some of the energy is converted into momentum that goes into the aerated region of the water column known.  This phenomenon is known as the surface roller. Under the assumption that the surface roller moves in the mean wave direction, the evolution and dissipation of the surface roller energy is calculated by a steady-steate energy balance equation (Stive and De Vriend 1994, Ruessink et al. 2001)
::<math>
\frac {\partial (2E_{sr}cw_j)} {\partial x_j} = -D_{sr} + f_e D_{br}
</math>(2-118)
where
::<math>E_{sr}</math> = surface roller energy density [N/m]
::<math>c</math>  = roller propagation speed [m/s]
::<math>w_j = (cos \theta_m, sin \theta_m)</math>  is the wave unit vector [-]
::<math>\theta_m</math> = mean wave direction [deg]
::<math>D_{sr}</math>  = roller energy dissipation [N/m/s]
::<math>D_{br}</math>  = wave breaking dissipation (from the wave model) [N/m/s]
The surface roller speed is calculated using the long-wave approximation <math>c = \sqrt{gh}</math> , where ''h'' is the local water depth. The surface roller dissipation is approximated as
::<math>
D_{sr} = \frac {g2E_{sr}\beta_D} {c}
</math>(2-119)
where <math>\beta_D</math>  is the surface roller dissipation coefficient approximately equal to 0.05-0.1. The surface roller contribution to the wave stresses, <math>R_{ij} </math> , is given by
::<math>
R_{ij} = 2E_{sr}w_i w_j
</math>(2-120)
One effect of the surface roller is to push the peak alongshore current ve-locity closer to shore.  Another side effect of the surface roller is to improve model stability (Sánchez et al. 2011a). The influence of the surface roller on the mean water surface elevation is relatively minor (Sánchez et al. 2011a).

Latest revision as of 15:32, 6 December 2022

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: .

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.

(2-42)

for , where N is the number of sediment size classes and

t = time [s]


h = total water depth [m]
=Cartesian coordinate in the jth direction [m]
= equilibrium total-load transport rate [kg/m/s]
= bed elevation with respect to the vertical datum [m]
= bed porosity [-]
= morphologic acceleration factor [-]
= sediment density [~2650 kg/m3 for quartz sediment]
= 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

(2-43)

where

t = time[s]
h = water depth [m]
= Cartesian coordinate in the
entrainment or pick-up function [kg/m2/s]
deposition or settling function [kg/m2/s]

The entrainment and deposition functions are calculated as

(2-44)

where

z = vertical coordinate from the bed [m]
= vertical sediment mixing coefficient [m2/s]
c = local sediment concentration [kg/m3]
= sediment fall velocity [m/s]
= calculated sediment concentration at an elevation a above the bed [kg/m3]
= equilibrium (capacity) sediment concentration at an elevation a above the bed [kg/m3]

Bed Change Equation

The bed change is calculated as

(2-46)

where

= bed elevation with respect to the vertical datum [m]
= bed porosity [-]
= morphologic acceleration factor [-]
= sediment density [~2650 kg/m3 for quartz sediment]
= empirical bed-slope coefficient (constant) [-]
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

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]
= Cartesian coordinate in the jth direction [m]
= depth-averaged total-load sediment mass concentration for size class k defined as in which is the total-load mass transport [kg/m3]
= depth-averaged total-load sediment mass concentration for size class k and described in the Equilibrium Concentration and Transport Rates section [kg/m3]
= total-load correction factor described in the Total-Load Cor-rection Factor section [-]
= fraction of suspended load in total load for size class k and is described in Fraction of Suspended Sediments section [-]
= horizontal sediment mixing coefficient described in the Hori-zontal Sediment Mixing Coefficient section [m2/s]
= total-load adaptation coefficient described in the Adaptation Coefficient section [-]
= 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

(2-48)

where and are the actual fraction of suspended- and total-load transport rates and and are the equilibrium fraction of suspended- and total-load transport rates.

Adaptation Coefficient

The total-load adaptation coefficient, , 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 that is related to the total-load adaptation length and time by


(2-49)

where

= 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 or , where Ls and Lb are the suspended- and bed-load adaptation lengths. Ls is defined as

(2-50)

in which and are the adaptation coefficient lengths for suspended load for the adaptation coefficient 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)


(2-51)

where is the bed shear stress, and is the von Karman constant. Armanini and di Silvio (1986) proposed an analytical equation

(2-52)

where is the thickness of the bottom layer defined by and is the zero-velocity distance from the bed. Gallapatti (1983) proposed the following equation to determine the suspended load adaptation time

(2-53)

where is the current related bottom shear velocity, , and .

The bed-load adaptation length, , 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 in which is an empirical coefficient on the order of 5-10. Fang (2003) found that of approximately two or three times the grid resolution works well for field applications. Although limited guidance exists on methods to estimate , the determination of 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, , 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, 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

(2-54)

where is the bed load velocity and 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 ( and ) 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 may be obtained as (Sánchez and Wu 2011)

(2-55)

where in which is the sediment fall velocity for size class k, is the vertical mixing coefficient, a is a reference height for the suspended load, h is the total water depth, is the apparent roughness length, and is the exponential integral. The equation can be further simplified by assuming that the reference height is proportional to the roughness height (e.g. ), so that . 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 .

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

(2-56)

where s is the specific gravity, g is the gravitational constant, is the median grain size diameter, is the bed shear stress related to the grain roughness and is determined by where is the Manning’s coefficient corresponding to the grain roughness and is the critical bed shear stress.

Bed Change Equation

The fractional bed change is calculated as

(2-57)

where

= bed elevation with respect to the vertical datum [m]
= bed porosity [-]
= morphologic acceleration factor [-]
= sediment density [~2650 kg/m3 for quartz sediment]
= depth-averaged total-load sediment mass concentration for size class k defined as in which is the total-load mass transport [kg/m3]
= depth-averaged total-load equilibrium sediment mass concentration for size class k and described in the Equilibrium Concentration and Transport Rates section [kg/m3]
= total-load adaptation coefficient described in the Adaptation Coefficient section [-]
= sediment fall velocity [m/s]
= empirical bed-slope coefficient (constant) [-]
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

(2-58)

The purpose of the morphologic acceleration factor is to speed-up the bed change so that the simulation time represents approximately the change that would occur in . 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)

(2-59)

where is the thickness of the first layer. is equal to for , and equal to the bed material gradation in the second sediment layer for . The bed-material sorting in the se-cond layer is calculated as

(2-60)

where is the thickness of the second layer, and <\math>p_{2k}</math> is the fraction of the 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

where is the beform height, and and 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

(2-62)

where 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 must be estimated from an empirical formula. The depth-averaged equilibrium concentration is defined as

(2-63)

where is the total-load transport for the sediment size class esti-mated from an empirical formula.

For convenience, is written in general form as

(2-64)

where is the fraction of the sediment size k in the first (top) bed layer and is the potential equilibrium total-load concentration. The potential concentration can be interpreted as the equilibrium concentration for uniform sediment of size . 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

(2-65)
(2-66)

where

= potential equilibrium bed load transport [kg/m/s]
= potential suspended load transport [kg/m/s]
= median grain size [m]
= sediment specific gravity or relative density [-]
= gravitational constant (9.81 m/s2)
= sediment density (~2650 kg/m3)
= Shields parameters due to currents [-]
= mean Shields parameters due to waves and currents [-]
= maximum Shields parameters due to waves and currents [-]
= critical Shields parameter [-]
= vertical sediment diffusivity [m2/s]
= reference bed concentration [kg/m3]
= bed-load scaling factor (default 1.0) [-]
= suspended-load scaling factor (default 1.0) [-]

The reference concentration is given by

(2-67)

where the coefficient is determined by the following relationship

(2-68)

The current-related shear stress is calculated as

(2-69)

where is the water density, is the bed friction coefficient, and is the current velocity magnitude. The drag coefficient is calculated as

(2-70)

where is the von Karman constant (0.4), is the total water depth, and is the roughness length calculated as , where is the total bed Nikuradse roughness. The total bed roughness is assumed to be a linear sum of the grain-related roughness , form-drag (ripple) roughness , and sediment-related roughness . Bed forms are also separated into current and wave-related bed forms. The current- and wave-related total roughness is then

(2-71)

where the subscript c|w indicates either the current (c) or wave (w) related component. The grain-related roughness is estimated as .

The ripple roughness is calculated as (Soulsby 1997)

(2-72)

where and are either the current- or wave-related ripple height and length respectively. The current-related ripple height and length are calculated as


(2-73)
(2-74)

The wave-related ripple height and length are calculated using the expressions proposed by van Rijn (1984b, 1989)

(2-75)
(2-76)


where is the semi-orbital excursion and is the wave mobility parameter. The semi-orbital excursion is defined as

<2-77)

in which is the peak bottom orbital velocity and is the wave period. For random waves, and are used as representative values. The wave mobility parameter, , is defined as


(2-78)

For the Lund-CIRP sediment transport equations, a wave bottom shear stress is calculated as

(2-79)

where is the wave friction factor calculated using the expression of Swart (1976)

(2-80)

where is the relative roughness defined as .

The current- and wave-related sediment roughnesses are estimated as

(2-81)

Equation 2-73 above must be solved simultaneously with the expressions for the bottom shear stress because the roughness depends on the stress. The exact solution is approximated using explicit polynomial fits in order to avoid time-consuming iterations in calculating the bed shear stress.

The critical Shields parameter is estimated using the formula proposed by Soulsby (1997)

(2-82)

where is the dimensionless grain size

(2-83)


where is the kinematic viscosity.

The sediment fall velocity is calculated using the formula by Soulsby (1997)

(2-84)

where is the grain size. The vertical sediment diffusivity is calculated as

(2-85)

where is the total effective dissipation given by

(2-86)

in which are coefficients, is the wave breaking dissipation, and are the bottom friction dissipation due to currents and waves, respectively. The dissipation from bottom friction due to current, , and the dissipation from bottom friction due to waves, , are expressed as

(2-87)

where the subscript c|w indicates either the current (c) or wave (w) relat-ed component, and and are the current- or wave-related bed shear velocity and stress, respectively. The coefficient =0.017, and are function of the Schmidt number:

(2-88)

where is either the current or wave-related Schmidt number calculated from the following relationships (Camenen and Larson 2008):

(2-89)


(2-90)

In the case of coexisting or combined waves and currents, the wave-current Schmidt number is estimated as

(2-91)

where is a weighting factor equal to .

For multiple-sized (nonuniform) sediments, the fractional equilibrium sediment transport rates are calculated as


(2-92)

(2-93)

where the subscript k indicates variables which are calculated based only on the sediment size class k. is the hiding and exposure coefficient de-scribed in Hiding and Exposure.

Van Rijn

The van Rijn (1984a,b) sediment transport equations for bed load and suspended load are used with the recalibrated coefficients of van Rijn (2007a,b), as given by


(2-94)
(2-95)

where

= equilibrium bed-load transport [kg/m/s]
= equilibrium suspended-load transport [kg/m/s]
= median grain size [m]
s = sediment specific gravity or relative density [-]
g = gravitational constant (9.81 m/2)
= sediment density (~2650 kg/m3)
U = depth-averaged current velocity [m/s]
= critical depth-averaged velocity for incipiant motion [m/s]
= effective depth averaged velocity [m/s]
= suspended-load scaling factor (default 1.0) [-]
= bed-load scaling factor (default 1.0) [-]

The effective depth-averaged velocity is calculated as with =0.4 for random waves and =0.8 for regular waves. Here is the bottom wave orbital velocity based on linear wave theory and the signifi-cant wave height The critical depth-averaged velocity is estimated as where is a blending factor and and are the critical velocities for currents and waves, respectively. As described in van Rijn (2007), the critical velocity for currents is

(2-96)

where is the sediment grain size in meters of 90th percentile. The critical velocity for waves is based on Komar and Miller (1975),

(2-97)

where is the peak wave period. According to van Rijn (2007), the bed load transport formula predicts transport rates by a factor of 2 for velocities higher than 0.6 m/s, but underpredicts transport rates by a factor of 2-3 for velocities close to the initiation of motion.

The near-bed equilibrium sediment concentration for the van Rijn (formula is given by (van Rijn 1985)

(2-98)

where is the maximum skin shear stress, is the critical shear stress, and is reference height above the bed given by

(2-99)

where is the bed ripple height.

The van Rijn formula (1984a,b, 2007a,b) were originally proposed for well-sorted sediments. The sediment availability is included by multiplication of transport rateswith the fraction of the sediment size class in the upper bed layer. The hiding and exposure is considered by in a correction factor which multiples to the critical velocity. When applied to multiple-sized sediments, the fractional equilibrium transport rate is calculated as

(2-100)
(2-101)

where is the fractional bed composition and is the hiding and ex-posure coefficient. The subscript k indicates values which are calculated based on the kth sediment size class.

Soulsby-van Rijn

Soulsby (1997) proposed the following equation for the total load sediment transport rate under action of combined current and waves,

(2-102)
(2-103)

where

= equilibrium total-load transport [kg/m/s]
= sediment density (~2650 kg/m3)
U = depth-averaged current velocity [m/s]
= effective velocity defined as [m/s]
= peak bottom wave orbital velocity based on the root-mean-squared wave height [m/s]
= drag coefficient due to currents alone [-]
= critical depth-averaged velocity for initiation of motion [m/s]
= empirical bed-slope coefficient, approximately equal to 1.6
= bed-slope [-]

The coefficients Asb and Ass are related to the bed and suspended loads, respectively, and given by

(2-104)

Note that the bed and suspended load scaling factors ( and ) are in-cluded in the calculations for the coefficients and .

The current drag coefficient is calculated as

(2-105)

where is the bed roughness length (equal to 0.006 m).

The Soulsby-van Rijn formula is modified for multiple-sized sediments similarly to the van Rijn formula in the previous section with the equation

(2-106)
(2-107)

where

The subscript k indicates that the value is calculated based only on the size class k and not the median grain size. The availability of sediment fractions is included through , while hiding and exposure of grain sizes is accounted for by modifying the critical velocity. it is noted that the Soulsby-van Rijn (1997) formula are very similar to the van Rijn (1984ab, 207ab) except for the definition of the effective velocity and the recalibration of the bed-load formula coefficients in 2007a. The proposed changes should be verified with measurements or numerical simulations for multiple-sized sediment transport.

Watanabe

The equilibrium total load sediment transport rate is determined by Watanabe (1987) as

(2-108)

where

= equilibrium total-load transport [kg/m/s]
= combined wave-current maximum shear stress [Pa]
= water density (~1025 kg/m3)
g = gravitational constant (9.81 m/s2)
= critical shear stress of incipient motion in [Pa]
= empirical coefficient typically ranging from 0.1 to 2

The critical shear stress is determined from the Shields diagram. The maximum bed shear stress is calculated as (Soulsby 1997)

(2-109)

where is the mean shear stress by waves and current over a wave cycle, is the mean wave bed shear stress, and is the angle between the waves and current. The wave bed shear stress is calculated as

(2-110)

where is a wave friction factor and is the bottom wave orbital velocity amplitude. The wave friction factor is calculated with the expression by Nielson (1992)

(2-111)

where r is the relative roughness defined as . is the semi-orbital excursion defined as .

Hiding and Exposure

When the bed material is composed of multiple grain sizes, larger grains have a greater probability of being exposed to the flow while smaller particles have a greater probability of being hidden to the flow. Figure 2.6 shows an example of a sediment grain being exposed to the flow by an exposure height , and sediment grain being hidden by .

Figure 2.6. Schematic of the exposure height of bed sediment grains

For the van Rijn (Reference(s)), Soulsby-van Rijn (Reference) and Watanabe (Reference) transport formulas, the hiding and exposure mechanism is considered by correcting the critical shear stress or velocity using a hiding and exposure correction function, . For the Lund-CIRP transport formula, an alternate approach is required due to the way in which the Shields number and grain size are included in the formulation; thus, the hiding and exposure correction function is directly used to multiply the transport rate. Two methods are used to calculate , depending on whether the sediment transport model is run with a single sediment size or with multiple sediment sizes and the methods are described in the following sections.

Single-sized sediment transport

In some applications, the coastal bed material is dominated by a single sediment size with patches of other sediment sizes or materials (e.g. shell hash) that may not contribute significantly to morphology change in the areas of interest; however, they may modify the sediment transport through hiding and exposure. For example, it is possible for the bed material to consist of mostly uniform sand with patches of shell fragments (bimodal distribution) in some regions. The shell material is difficult to model numerically because it is usually poorly sorted and its hydraulic properties are unknown. For such regions, sediment transport models commonly estimate excessive erosion due to the absence of the hiding effect of the coarser shell material (e.g. Cayocca 2001). A better and more physical plausible approach is to use the local bed composition along with a correction function for hiding and exposure to account for the effects of the uniform sand with the patches of coarser shell material. For single-sized sediment transport, the correction function for hiding and exposure is calculated following Parker et al. (1982) and others as,

(2-112)

where i is an empirical coefficient between 0.5-1.0. The aforementioned sediment transport equations are implemented by using the transport grain size rather than the bed material . A single and constant transport size is used, while the bed material varies spatially. The spatial distribution of can be obtained from field measurement data and for simplicity is assumed constant during the model simulation time. This is a significant assumption and may not be reasonable for some applications. However, this method provides a simple conceptual mechanism for considering an important process in the proposed single-sized sediment transport model. The approach has been successfully applied to Shinnecock Inlet, NY to simulate morphology change at a coastal inlet (Sánchez and Wu 2011a). A more accurate and complex approach is to simulate the transport and sorting of multiple sediments.

Multiple-sized sediment transport

The hiding and exposure of the each sediment size class is considered by modifying the critical shields parameter for each sediment size class based on Wu et al. 2000

(2-113)

where m is an empirical coefficient that varies for each transport formula, approximately equal to 0.6-1.0. and are the total hiding and exposure probabilities and are calculated as

(2-114)


where N is the number of grain size classes.

Horizontal Sediment Mixing Coefficient

The horizontal sediment mixing coefficient, vs, represents the combined effects of turbulent diffusion and dispersion due to nonuniform vertical profiles. In CMS, the horizontal sediment mixing coefficient is assumed to be proportional to the turbulent eddy viscosity as

(2-115)

where is the Schmidt number and is…. There are many formulas to estimate the Schmidt number. The Schmidt number is set by default equal to 1.0 in CMS, but may be modified by the user with cards.

Boundary Conditions

At the interface between wet and dry cells, the sediment transport rate is set to zero. The inflow boundary condition requires a given sediment concentration at the boundary. However, for most coastal applications, the actual sediment concentration is not available and the model implements the equilibrium concentration capacity (Dirichlet boundary condition). CMS also requires the size distribution of the inflow transport. For stability reasons, the inflow size distribution is assumed to be equal to the initial size distribution of the bed at the boundary. Inflow sediment transport rates are specified either as a total sediment transport rate, in kg/s, or as a fractional sediment transport rate, in kg/s, for each sediment size class. If only the total sediment transport rate is specified, then the fractional sediment transport rate is calculated as where is the bed material fraction of the first layer at the boundary. The fractional sediment transport at boundary cell B in kg/m/s is then calculated along the cell string according to

(2-116)

where is the potential sediment transport rate at boundary cell B and is the inflow cell face width. If the flow is directed outward of the do-main, a zero-gradient boundary condition is used for sediment concentration (Neumann boundary condition).

Salinity Transport

Transport Equation

CMS calculates the salinity transport based on the following 2DH salinity conservation equation

(2-117)

where

= depth-averaged salinity [ppt or parts per thousand]
h = water depth [m]
= total flux velocity [m/s]
= horizontal mixing coefficient
= source/sink term due to precipitation, evaporation, and structures (e.g. culverts) [ppt m/s]

Equation (2-106) represents the horizontal fluxes of salt in water bodies and is balanced by exchanges of salt via diffusive fluxes. Major processes that contribute to the salinity are: freshwater inflows from rivers, vertical fluxes of freshwater by precipitation and evaporation at the water surface, and groundwater fluxes. These processes can be specified as the surface and bottom boundary conditions in the equation.

Boundary Conditions

The boundary condition for the salinity transport equation is dependent on the flow direction. At cell faces between wet and dry cells, a zero-flux boundary condition is applied. The cell face velocity is zero between wet and dry cells and this eliminates advection transport. In addition, the salinity horizontal mixing coefficient at the cell face, , is set to zero and this eliminates the diffusive flux. If the flow is directed inward of the modeling domain, then a user-specified salinity concentration is pre-scribed (Direchlet boundary condition). If the flow is directed outward of the modeling domain, then a zero-gradient boundary condition is applied (Neumann boundary condition). Overall comments for this section: 1) A figure would be very HELPFUL showing inward and outward directed flow with the modeling domain and 2) Direchlet and Neumann BCs should be explained with reference citations.

Surface Roller

As a wave transitions from non-breaking to fully-breaking, some of the energy is converted into momentum that goes into the aerated region of the water column known. This phenomenon is known as the surface roller. Under the assumption that the surface roller moves in the mean wave direction, the evolution and dissipation of the surface roller energy is calculated by a steady-steate energy balance equation (Stive and De Vriend 1994, Ruessink et al. 2001)

(2-118)

where

= surface roller energy density [N/m]
= roller propagation speed [m/s]
is the wave unit vector [-]
= mean wave direction [deg]
= roller energy dissipation [N/m/s]
= wave breaking dissipation (from the wave model) [N/m/s]

The surface roller speed is calculated using the long-wave approximation , where h is the local water depth. The surface roller dissipation is approximated as

(2-119)

where is the surface roller dissipation coefficient approximately equal to 0.05-0.1. The surface roller contribution to the wave stresses, , is given by

(2-120)

One effect of the surface roller is to push the peak alongshore current ve-locity closer to shore. Another side effect of the surface roller is to improve model stability (Sánchez et al. 2011a). The influence of the surface roller on the mean water surface elevation is relatively minor (Sánchez et al. 2011a).