CMSFlow:Nonequilibrium Sediment Transport
UNDER CONSTRUCTION
written by Alejandro Sanchez Last date modified: December 8, 2010
Totalload Transport Equation
The singlesized sediment transport model described in Sánchez and Wu (2011a) was extended to multiplesized sediments within CMS by Sánchez and Wu (2011b). In this model, the sediment transport is separated into current and waverelated transports. The transport due to currents includes the stirring effect of waves, and the waverelated transport includes the transport due to asymmetric oscillatory wave motion as well as steady contributions by Stokes drift, surface roller, and undertow. The currentrelated bed and suspended transports are combined into a single totalload transport equation, thus reducing the computational costs and simplifying the bed change computation. The 2DH transport equation for the currentrelated total load is

(1) 
for j= 1,2; k=1,2...,N, where N is the number of sediment size classes and
 = actual depthaveraged totalload sediment concentration [kg/m^{3}] for size class k defined as in which is the totalload mass transport
 = equilibrium depthaveraged totalload sediment concentration [kg/m^{3}] for size class k and described in the equilibrium concentration and transport rates section
 = totalload correction factor described in the TotalLoad Correction Factor section []
 = fraction of suspended load in total load for size class k described in fraction of suspended sediments section []
 = horizontal sediment mixing coefficient described in the horizontal sediment mixing coefficient section [m^{2}/s]
 = totalload adaptation coefficient described in the adaptation coefficient section []
 = sediment fall velocity [m/s].
In the above equation, the first term represents the temporal variation of ; the second term represents the horizontal advection; the third term represents the horizontal diffusion and dispersion of suspended sediments; and the last term represents the erosion and deposition. The equation may be applied to singlesized 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 units of sediment concentration used here are kg/m^{3} rather than dimensionless volume concentrations in order to avoid numerical precision errors at low concentrations.
In the above equations, it is assumed that the wave flux velocity is not included in the momentum equations. If the wave flux velocity is included, then the total flux velocity (V) should be used instead of the depthaveraged current velocity (U). The reason for this is because without a waveinduced sediment transport to counter the offshore directed transport due to the undertow, the model would predict excessive movement of sediment offshore.
Fraction of Suspended Sediment
In order to solve the system of equations for sediment transport implicitly, the fraction of suspended sediment must be determined explicitly. This is done by assuming

(2) 
where are the actual suspended and totalload transport rates and are the equilibrium suspended and totalload transport rates.
Adaptation Coefficient
The totalload adaptation coefficient is an important parameter in the sediment transport model. There are many variations of this parameter in literature (Lin 1984; Gallappatti and Vreugdenhil 1985; and Armanini and di Silvio 1986). CMS uses a totalload adaptation coefficient that is related to the totalload adaptation length (L_{t} ) and time (T_{t}) by

(3) 
where:
 = sediment fall velocity corresponding to the transport grain size for singlesized sediment transport or the median grain size for multiplesized sediment transport [m/s]
 U = depthaveraged current velocity [m/s]
 h = water depth [m].
The adaptation length or time is a characteristic distance or time for sediment to adjust from nonequilibrium to equilibrium transport. Because the total load is a combination of the bed and suspended loads, the associated adaptation length may be calculated as are the suspended and bedload adaptation lengths. The symbol is defined as

(4) 
in which are the adaptation coefficient lengths for suspended load. The adaptation coefficient can be calculated either empirically or based on analytical solutions to the pure vertical convectiondiffusion equation of suspended sediment. One example of an empirical formula is that proposed by Lin (1984):

(5) 
where u_{*} is the bed shear stress, and is the von Karman constant. Armanini and di Silvio (1986) proposed an analytical equation:
(6) 
where is the thickness of the bottom layer defined by is the zerovelocity distance from the bed. Gallappatti (1983) proposed the following equation to determine the suspendedload adaptation time:
(7) 
where u_{*} is the current related bottom shear velocity, .
The bedload 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 bedload adaptation length can be estimated as L_{b} = a_{b}h in which a_{b} is an empirical coefficient on the order of 510. 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.
Totalload correction factor
The totalload correction factor (β_{tk}) accounts for the vertical distribution of the suspended sediment concentration and velocity profiles as well as the fact that bed load travels at a slower velocity than the depthaveraged current velocity (Figure 1) (Wu 2007). By definition, β_{tk} is the ratio of the depthaveraged totalload and flow velocities.
Figure 1. Schematic of sediment and current vertical profiles.
In a combined bedload and suspendedload model, the correction factor is given by
(8) 
where u_{bk} is the bedload velocity, and β_{sk} is the suspendedload correction factor and is defined as the ratio of the depthaveraged suspended sediment and flow velocities. Since most sediment is transported near the bed, both the total and suspendedload correction factors (β_{tk} and β_{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 suspendedload correction factor (β_{sk}) may be obtained as (Sánchez and Wu 2011b)
(9) 
where:
A = a/h []
Z = z_{a}/h []
= vertical mixing coefficient [m^{2}/s]
a = reference height near the bed for the suspended load [m]
z_{a} = apparent roughness length [m]
(exponential integral).
The equation can be further simplified by assuming that the reference height is proportional to the apparent roughness length (e.g., a = 30z_{z}) so that . Figure 2 shows a comparison of the suspendedload correction factor based on the logarithmic velocity with exponential and Rouse suspended sediment concentration profiles. Both cases show similar behavior for the suspendedload correction factor (). For fine sediments (small fall velocity), is close to 1.0 and experiences smaller influences by the bottom roughness, while for course sediments, can be as low as 0.5 and is largely influenced by the bottom roughness. This is to be expected since course sediments are transported more closely to the bottom, compared to fine sediments.
Figure 2. Suspended load correction factors based on the logarithmic velocity profile and (a) exponential and (b) Rouse suspended sediment profile. The Rouse number is
The bedload velocity (u_{bk}) in Equation 8 is calculated using the van Rijn (1984a) formula with recalibrated coefficients from Wu et al. (2006):
(10) 
where:
s = specific gravity []
g = gravitational constant (~9.81 m/s^{2})
d_{k} = characteristic grain diameter for the k^{th} size class [m]
= grainrelated bed shear stress [Pa]
= grainrelated Manning’s roughness coefficient [s/m^{1/3}]
= critical bed shear stress for the k^{th} size class [Pa].
Bed Change Equation
The change in the water depth is calculated by
(11) 
where:
z_{b} = bed elevation with respect to the vertical datum [m]
p_{m}^{'} = bed porosity []
= sediment (material) density (~2650 kg/m^{3} for quartz sediment)
D_{s} = empirical bedslope coefficient (constant) []
is the bedload mass transport rate magnitude [kg/m/s].
The first term on the righthand side of the above equation represents the bed change due to sediment exchange near the bed. The last term accounts for the effect of the bed slope on bedload transport. The bedslope coefficient (D_{s}) is usually about 0.1 to 3.0. 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 changes for all size classes:
(12) 
Bed Sorting and Gradation (nonuniform sediments)
When simulating nonuniform sediments it is necessary to keep track of the vertical variation in bed composition. In order to do this, the bed is divided into discrete layers. The top layer is the mixing or active layer and is the layer of sediment which is actively being exchanged with the bed and suspended loads. The temporal variation of the bedmaterial in the mixing and second layers is calculated as (Wu 2004)

(11) 
where and are the thicknesses of the mixing and second layers respectively. for or otherwise. Figure 3 shows an example of the bed layer evolution over time for varying bed deposition and erosion.
Boundary Conditions
There are three types of boundary conditions in the sediment transport: Wetdry, Outflow and Inflow.
1. Wetdry interface.
 The interface between wet and dry cells has a zeroflux boundary condition. Both the advective and diffusive fluxes are set to zero at the wetdry interfaces. Note that avalanching may still occur between wetdry cells.
2. Outflow Boundary Condition
 Outflow boundaries are assigned a zerogradient boundary condition and sediments are allowed to be transported freely out of the domain.
3. Inflow Boundary Condition
 When flow is entering the domain, it is necessary to specify the sediment concentration. In CMSFlow, the inflow sediment concentration is set to the equilibrium sediment concentation. For some cases, it is desired to reduce the amount of sediment entering from the boundary such as in locations where the sediment source is limited (i.e. coral reefs). The inflow equilibrium sediment concentration may be adjusted by multiplying by a loading scaling factor and is specified by the Advanced Card NET_LOADING_FACTOR or SEDIMENT_INFLOW_LOADING_FACTOR.
Numerical Methods
The governing equations are discretized using the Finite Volume Method on a staggered, nonuniform Cartesian grid. Time integration is calculated with a simple explicit forward Euler scheme. Diffusion terms are discretized with the standard central difference scheme. Advection terms are discretized with either the first order upwind scheme or the second order Hybrid Linear/Parabolic Approximation (HLPA) scheme of Zhu (1991). The default advection scheme is HLPA but may be changed with the Advanced Card ADVECTION_SCHEME.
References
 Buttolph, A. M., C. W. Reed, N. C. Kraus, N. Ono, M. Larson, B. Camenen, H. Hanson, T. Wamsley, and A. K. Zundel. (2006). “Twodimensional depthaveraged circulation model CMSM2D: Version 3.0, Report 2: Sediment transport and morphology change.” Coastal and Hydraulics Laboratory Technical Report ERDC/CHL TR069. Vicksburg, MS: U.S. Army Engineer Research and Development Center, U.S.A.
 Camenen, B., and Larson, M. (2007). “A unified sediment transport formulation for coastal inlet application”. Technical Report ERDCCHL CR0701. Vicksburg, MS: U.S. Army Engineer Research and Development Center, U.S.A
 Parker, G., Kilingeman, P. C., and McLean, D. G. (1982). “Bed load and size distribution in paved gravelbed streams.” J. Hydr. Div., ASCE, 108(4), 544571.
 Soulsby, R. L. (1997). Dynamics of marine sands, a manual for practical applications. H. R. Wallingford, UK: Thomas Telford.
 Watanabe, A. (1987). “3dimensional numerical model of beach evolution”. Proc. Coastal Sediments ’87, ASCE, 802817.
 Wu, W. (2004).“Depthaveraged 2D numerical modeling of unsteady flow and nonuniform sediment transport in open channels”. J. Hydraulic Eng., ASCE, 135(10), 1013–1024.
 van Rijn, L. C. (1985). “Flume experiments of sedimentation in channels by currents and waves.” Report S 347II, Delft Hydraulics laboratory, Deflt, Netherlands.
 Zhu, J. (1991). “A low diffusive and oscillationfree convection scheme”. Com. App. Num. Meth., 7, 225232.
 Zundel, A. K. (2000). “Surfacewater modeling system reference manual”. Brigham Young University, Environmental Modeling Research Laboratory, Provo, UT.
External Links
 Aug 2006 TwoDimensional DepthAveraged Circulation Model CMSM2D: Version 3.0, Report 2, Sediment Transport and Morphology Change [1]
 Aug 2008 CMSWave: A Nearshore Spectral Wave Processes Model for Coastal Inlets and Navigation Projects [2]
Powerpoint presentation on NET