The bed sorting and layering model begins with the calculation of the layer thickness of the first and second layers. This is followed by the dry bulk density of the first layer and finally the grain fractions by weight of the first and second layers.

The active layer thickness is given by

1) \delta _{1i}^{n+1}=\min \left[\max \left(f_{1,90}d_{90,i},0.5\Delta ,\Delta z_{bi},\delta _{1,\min }\right),\delta _{1,\max }\right]

where

δ_{1,min} : user-specified minimum active layer thickness [L]

δ_{1,max}  : user-specified maximum active layer thickness [L]

f_{1,90} : user-specified active layer scaling factor corresponding to d90 (~ 1 to 10) [-]

d_{90,i} : 90th percentile diameter [L]

\Delta : bed form height [L]

\Delta z_{bi} : bed change [L]

If beforms are not being simulated, the bedform height is set to zero. The active layer thickness has a lower limit equal to the bed change so that in the case of deposition the boundary between he first and second layers does not go above the bed elevation of the previous time step. Basically, all deposition is occurs in the active layer and not deposition occurs in the active stratum.

The bed layering is computed at each subarea. However, here the subscripts for the subarea are dropped for simplicity. The minimum and maximum active layer thicknesses are used to avoid excessively small and large active layer thicknesses. It is noted that by setting the limits to the same value, then the above equation is equivalent to specifying a constant value. The thickness of the second layer is given by

2) \delta _{2i}^{n+1}=\delta _{2i}^{n}+\Delta \delta _{2i}

where

\Delta \delta _{2i}=\Delta z_{bi}-\Delta \delta _{1i} : change in second layer thickness [L]

\Delta z_{bi}=z_{bi}^{n+1}-z_{bi}^{n} : bed change [L]

\Delta \delta _{1i}=\delta _{1i}^{n+1}-\delta _{1i}^{n} : change in active layer thickness [L]

In order to avoid the second layer from becoming extremely thin or thick, a layer merging and splitting algorithm is implemented between layers 2 and 3. If the second layer is too thick, it is divided into two layers; thus, the previous third layer becomes the new fourth layer, and the last two bottom layers are merged into one. If the second layer is too thin, it is merged with the previous third layer to form a new second layer; thus, the previous fourth layer becomes the new third layer. To illustrate the bed layering process, the figure below shows an example of the temporal evolution of 7 bed layers during erosional and depositional regimes. Therefore, a minimum of 3 layers are required for the model. The bed layering model requires a minimum of 3 layers but does not have a maximum number of layers. Furthermore, the number of layers can vary from cell to cell. The number of layers is specified for each cell. All of the subregions within a cell have the same number of layers.

Figure 1. Schematic showing an example bed layer evolution. 


The bed sorting equations for the first and second bed layers are discretized explicitly as

3) m_{1ki}^{n+1}=\frac{\Delta M_{bki}+m_{1ki}^{n}\delta _{1i}^{n}-m_{*ki}^{n}\Delta \delta _{2i}}{\delta _{1i}^{n+1}}
4) m_{2ki}^{n+1}=\frac{m_{2ki}^{n}\delta _{2i}^{n}+m_{*ki}^{n}\Delta \delta _{2i}}{\delta _{2i}^{n+1}}

where

m_{1ki} = f_{1ki}ρ_{d1i} : fractional mass concentration in the first (active) layer [M/L3]

m_{2ki} = f_{2ki}ρ_{d2i} : fractional mass concentration in the in second layer [M/L3]

ρ_{d1i} = \sum_k m_{1ki} : dry bulk density of first layer [M/L3]

ρ_{d2i} = \sum_k m_{2ki} : dry bulk density of second layer [M/L3]

\Delta M_{bki} : fractional mass exchange with the bed in [M/L2]

\Delta \delta _{2i} : change in second layer thickness [L]

m_{*ki}^{n}=\left\{\begin{array} m_{1ki}^{n}\,\,\mathrm{for}\,\,\Delta \delta _{2i}\geq 0\\ m_{2ki}^{n}\,\,\mathrm{for}\,\,\Delta \delta _{2i}<0 \end{array}\right.

In the above equations, the superscripts indicate the time step level and the number subscripts indicate the bed layer. The above formulation takes into account the existing porosity of the first and second layers and the porosity of newly deposited material. The formulation allows for the bed to be eroded more than the active layer thickness in a single time step.

In order to close the system of equations, an approximation is required to estimate the porosity of the deposited and eroded material (referred to here as the exchange material) is computed based on the net volume erosion or deposition as

5) \phi _{b}=\left\{\begin{array} \phi _{0}\,\,\,\mathrm{for}\,\,\Delta z_{b}>0\\ \phi _{1}\,\,\,\,\mathrm{for}\,\Delta z_{b}\leq 0 \end{array}\right.

where

\Delta z_{b} : bed change [L]

ϕ_{0} : porosity of (newly) deposited material [-]

ϕ_{1} : porosity of layer 1 (active layer) [-]