Momentum Dispersion

Definition of Momentum Dispersion Stresses

The momentum dispersion stresses are defined as

1)	$\begin{array}{l}\displaystyle \boldsymbol{D} = D_{ij} = \frac{\rho}{h} \int_{z_b}^{z_s} \left( u_i(z) - U_i \right) \left( u_j(z) - U_j \right) dz\end{array}$

where

$\begin{array}{l}\boldsymbol{D} = D_{ij}\end{array}$ : dispersion stress tensor [M/L/T²]
$\begin{array}{l}U_i, U_j\end{array}$ : depth-averaged velocity components [L/T]
$\begin{array}{l}u_i(z), u_j(z)\end{array}$ : depth-varying velocity components [L/T]
$\begin{array}{l}\rho\end{array}$ : water density [M/L³]
$\begin{array}{l}h\end{array}$ : water depth [L]
$\begin{array}{l}z_b\end{array}$ : bed elevation [L]
$\begin{array}{l}z_s\end{array}$ : water surface elevation [L]

The first term in the above equation generates a momentum flux towards the outside of a bend.

Secondary Flow

In river bends, the balance between the outward centrifugal acceleration and the inward pressure gradient produces a secondary flow which flows towards the outer bank near the surface and towards the inner bank near the bed (Boussinesq 1868, Bowker 1988). The most accurate way of modeling the secondary currents at river bends is using a 3D model. However, the effects of the secondary flow can be approximated to some extent in depth-averaged 2D models by approximating the vertical current velocity profiles in the streamwise and transverse directions to estimate dispersive stresses (Duan 2004; Wu et al. 2004, 2005; Song et al. 2012). The streamwise velocity profile is typically assumed to be either power-law (e.g. . Odgaard 1986; Wu et al. 2004) or logarithmic (e.g. Duan 2004; De Vriend 1977). Several modeling approaches have been proposed to parameterize the vertical profile of the transverse current velocity which represents the secondary flow (e.g. Odgaard 1986; De Vriend 1977; 1981; Bernard and Schneider 1992). The most common vertical profile for the transverse current velocity is a linear one with zero mean (Odgaard 1986). The purpose here is not to accurately resolve the 3D velocity field but rather estimate the simple and efficient relations for the dispersion terms with reasonable accuracy. The linear transverse velocity profile has been found to produce good results compared to laboratory measurements of current velocity profiles and surface water surface elevations (Wu et al. 2004; 2005; Begnudelli et al. 2010; Lazzarin and Viero 2023). The magnitude of the linear transverse velocity is parameterized as a function of the local flow curvature.

The secondary/spiral flow intensity, $\begin{array}{l}I\end{array}$ , is defined such that it is related to the transversal velocity at the free-surface by the simple expression

2)	$\begin{array}{l}\displaystyle \boldsymbol{u}_{ns} = \boldsymbol{n}_n b_s I\end{array}$

where

$\begin{array}{l}\boldsymbol{u}_{ns}\end{array}$ : surface transverse secondary flow velocity [L/T]
$\begin{array}{l}\boldsymbol{n}_n\end{array}$ : unit vector perpendicular to the streamwise direction [-1,1] [-]
$\begin{array}{l}b_s\end{array}$ : an empirical non-dimensional coefficient [-]
$\begin{array}{l}I\end{array}$ : secondary flow intensity [L/T]

The spiral flow intensity is defined using the right-hand-rule such that positive values indicate counter-clockwise rotations and negative values indicate clockwise rotations.

Various formulations have been used in literature for the coefficient . Here the formulation of Kalwijk and Booij (1986) is adopted:

3)	$\begin{array}{l}\displaystyle b_s = \frac{3}{\kappa^2} \left(\frac{1}{2} - \frac{1}{m} \right)\end{array}$

The secondary flow intensity can be either assumed to be in equilibrium or non-equilibrium. The non-equilibrium approach requires solving a transport equation and is therefore much more computationally expensive.

The transverse or spanwise velocity is approximated by the linear profile

4)	$\begin{array}{l}\displaystyle \frac{u_n(\xi)}{u_{ns}} =2\xi - 1\end{array}$

where

$\begin{array}{l}\xi\end{array}$ : non-dimensional vertical coordinate between 0 and 1 [-]
$\begin{array}{l}u_{ns}=u_n(\xi=1)\end{array}$ : surface transversal (i.e. spiral or helical) velocity [L/T]
$\begin{array}{l}u_{nb}=u_n(\xi=0)\end{array}$ : near-bed transversal (i.e. spiral or helical) velocity [L/T]

The intensity of the secondary flow can therefore be described by the magnitude of velocity magnitude at the surface and bed. In general, there are two approaches for estimating the secondary flow intensity: (1) Non-equilibrium, and (2) equilibrium. In the equilibrium approach, the secondary flow intensity is estimated based solely on the local hydrodynamic conditions, and the production and dissipation of secondary flow are assumed to be equal. In the non-equilibrium approach, a transport equation is used to simulate the unsteady production, dissipation, advection, and diffusion of the secondary flow intensity. Here two models are implemented which parameterize the secondary flow intensity as a function of a transport variable.

Equilibrium Approach

In the equilibrium approach the secondary flow intensity is computed simply as

5)	$\begin{array}{l}\displaystyle I = I_e = \beta_I h \|\boldsymbol{U}\| K_c\end{array}$

where

$\begin{array}{l}I\end{array}$ : secondary flow intensity [L/T]
$\begin{array}{l}I_e\end{array}$ : equilibrium secondary flow intensity [L/T]
$\begin{array}{l}\beta_I\end{array}$ : user-specified non-dimensional calibration parameter approximately equal to 1 [-]
$\begin{array}{l}h\end{array}$ : water depth [L]
$\begin{array}{l}K_c\end{array}$ : flow curvature [1/L]
$\begin{array}{l}\boldsymbol{U}\end{array}$ : depth-averaged current velocity [L/T]

Non-equilibrium Approach

In the non-equilibrium approach, the secondary flow intensity is calculated by solving the following transport equation

6)	$\begin{array}{l}\displaystyle \frac{\partial(h I)}{\partial t} + \nabla \cdot \left( h \boldsymbol{U} I \right) = \nabla \cdot \left( \varepsilon_{Ih} h \nabla I\right) + \frac{h \| \boldsymbol{U}\| }{L_a} \left( I_e - I \right)\end{array}$

7)	$\begin{array}{l}\displaystyle I_e = \beta_I h \|\boldsymbol{U}\| K_c\end{array}$

in which

$\begin{array}{l}I\end{array}$ : secondary flow intensity [L/T]
$\begin{array}{l}\beta_I\end{array}$ : secondary flow intensity scale factor [-]
$\begin{array}{l}h\end{array}$ : water depth [L]
$\begin{array}{l}L_a\end{array}$ : adaptation length [L]
$\begin{array}{l}\boldsymbol{U}\end{array}$ : current velocity magnitude [L/T]
$\begin{array}{l}\varepsilon_{Ih}\end{array}$ : secondary flow intensity horizontal diffusion coefficient [L²/T]
$\begin{array}{l}K_c = 1/R_c\end{array}$ : local flow curvature [1/L]
$\begin{array}{l}R_c\end{array}$ : local radius of curvature [L]

Secondary Flow Adaptation Length

The adaptation length is components as Kalkwijk and Booij (1986)

8)	$\begin{array}{l}\displaystyle L_a = f_L \frac{(m-1)h}{2\kappa^2}\end{array}$

where

$\begin{array}{l}f_L\end{array}$ : user-specified scaling factor (default of 1)
$\begin{array}{l}m\end{array}$ : non-dimensional roughness parameter [-]
$\begin{array}{l}h\end{array}$ : water depth [L]
$\begin{array}{l}\kappa\end{array}$ : von Karman constant of 0.41

Secondary Flow Diffusion Coefficient

The secondary flow intensity diffusion coefficient is scaled by the turbulent eddy viscosity as

9)	$\begin{array}{l}\displaystyle L_a = f_L \frac{(m-1)h}{2\kappa^2}\end{array}$

where

$\begin{array}{l}\varepsilon_{Ih}\end{array}$ : secondary flow intensity diffusion coefficient [L²/T]
$\begin{array}{l}\nu_t\end{array}$ : turbulent eddy viscosity [L²/T]
$\begin{array}{l}f_{D}\end{array}$ : scaling factor O(1) [-]

Flow Curvature

The horizontal curvature of flow streamlines is represented by the flow curvature. The local flow curvature is estimated using the local depth-averaged velocities and gradients as (Theisel 1995)

10)	$\begin{array}{l}\displaystyle K_c = \frac{1}{R_c} = \frac{U^2\dfrac{\partial V}{\partial x} + UV\left( \dfrac{\partial V}{\partial y} - \dfrac{\partial U}{\partial x} \right) -V^2\dfrac{\partial U}{\partial y}}{\|\boldsymbol{U}\|^3}\end{array}$

where

$\begin{array}{l}K_c\end{array}$ : flow curvature [1/L]
$\begin{array}{l}R_c\end{array}$ : radius of curvature [L/T]
$\begin{array}{l}\boldsymbol{U} = (U,V)^T\end{array}$ : depth-averaged velocity [L/T]

The sign of the flow curvature indicates the direction in which the velocity is rotating along a streamline. HEC-RAS utilizes the right-hand convention in which positive flow curvatures indicate a counter-clockwise rotation, and negative flow curvatures indicate a clockwise rotation. With this convention, the flow curvature is positive for river bends which turn to the left and negative for river bends which turn to the right.

To improve stability the radius of curvature is limited as

11)	$\begin{array}{l}\displaystyle R_c > R_{c,\text{min}} = \min ( 10h,3 \sqrt{\Delta x} )\end{array}$

where

$\begin{array}{l}h\end{array}$ : water depth [L]
$\begin{array}{l}\Delta x\end{array}$ : local grid resolution estimated as the square root of the cell area [L]

Streamwise Velocity Profile

The streamwise component of the velocity can be assumed to have a power-law profile as

12)	$\begin{array}{l}\displaystyle \frac{u_s(\xi)}{\|\boldsymbol{U}\|} = \frac{m+1}{m} \xi^{1/m}\end{array}$

where

$\begin{array}{l}u_s(\xi)\end{array}$ : streamwise velocity [L/T]
$\begin{array}{l}\xi\end{array}$ : non-dimensional vertical coordinate [0, 1] [-]
$\begin{array}{l}m\end{array}$ : non-dimensional roughness parameter (approximately equal to 7) [-]
$\begin{array}{l}|\boldsymbol{U}|\end{array}$ : depth-averaged velocity magnitude [L/T]

The roughness parameter $\begin{array}{l}m\end{array}$ can be calculated from the Manning’s roughness coefficient by

13)	$\begin{array}{l}\displaystyle m = \frac{\kappa R^{1/6}}{g \sqrt{g}}\end{array}$

in which

$\begin{array}{l}\kappa\end{array}$ : von Karman constant [-]
$\begin{array}{l}n\end{array}$ : Manning’s roughness coefficient [T/L^1/3]
$\begin{array}{l}g\end{array}$ : gravity [L/T²]
$\begin{array}{l}R\end{array}$ : hydraulic radius [L]

Calculation of Momentum Dispersion Stresses

The momentum dispersion stresses can be calculated analytically assuming utilizing the streamwise and transverse velocity profiles described above (e.g. Wu and Wang 20024; Wu 2007). However, here the momentum dispersion stresses are calculated following the approach by Deltares (2024) as

14)	$\begin{array}{l}\displaystyle D_{xx} = - D_{yy} = 2 \rho \dfrac{UV}{\|\boldsymbol{U}\|} I \beta\end{array}$

15)	$\begin{array}{l}\displaystyle D_{xy} = D_{yx} = - \rho \dfrac{U^2 - V^2}{\|\boldsymbol{U}\|} I \beta\end{array}$

in which

$\begin{array}{l}I\end{array}$ : secondary flow intensity [L/T]
$\begin{array}{l}\beta = \beta_c \left( 5\alpha - 15\alpha^2 + 37.5\alpha^3 \right)\end{array}$ (Deltares 2024)
$\begin{array}{l}\beta_c \in [0,1]\end{array}$ : user-defined calibration parameter [-]
$\begin{array}{l}\alpha = 1/m = \dfrac{ n \sqrt g}{\kappa R^{1/6}}\end{array}$ [-]
$\begin{array}{l}\kappa\end{array}$ : von Karman constant [-]
$\begin{array}{l}n\end{array}$ : Manning’s roughness coefficient [T/L^1/3]
$\begin{array}{l}g\end{array}$ : gravity [L/T²]
$\begin{array}{l}R\end{array}$ : hydraulic radius [L]

The Deltares (2024) approach utilizes a reduced set of terms. The method is preferred here because it is generally more robust than other approaches and can be more readily calibrated, especially as implemented here with scaling factors for the spiral intensity, adaptation length, and diffusion coefficient. The dispersion terms are such that when the flow curvature is zero, the dispersion terms also go to zero, making model calibration simpler and easier to understand.

Although, the streamwise and transverse velocity profiles are not utilized for computing the dispersion stresses, the are utilizes for computing surface and near-bed velocities. The surface velocities are only utilized for output whereas the near-bed velocities are utilized for output and sediment transport calculations.