Influence of Dunes on Channel‐Scale Flow and Sediment Transport in a Sand Bed Braided River

Current understanding of the role that dunes play in controlling bar and channel‐scale processes and river morphodynamics is incomplete. We present results from a combined numerical modeling and field monitoring study that isolates the impact of dunes on depth‐averaged and near‐bed flow structure, with implications for morphodynamic modeling. Numerical simulations were conducted using the three‐dimensional computational fluid dynamics code OpenFOAM to quantify the time‐averaged flow structure within a 400 m × 100 m channel using digital elevation models (DEMs) for which (i) dunes and bars were present within the model and (ii) only bar‐scale topographic features were resolved (dunes were removed). Comparison of these two simulations shows that dunes enhance lateral flows and reduce velocities over bar tops by as much as 30%. Dunes influence the direction of modeled sediment transport at spatial scales larger than individual bedforms due to their effect on topographic steering of the near‐bed flow structure. We show that dunes can amplify, dampen, or even reverse the deflection of sediment down lateral bar slopes, and this is closely associated with 3‐D and obliquely orientated dunes. Sediment transport patterns calculated using theory implemented in depth‐averaged morphodynamic models suggest that gravitational deflection of sediment is still controlled by bar‐scale topography, even in the presence of dunes. However, improved parameterizations of flow and sediment transport in depth‐averaged morphodynamic models are needed that account for the effects of both dune‐ and bar‐scale morphology on near‐bed flow and sediment transport.

The relative simplicity and computational efficiency of 2-D models compared to 3-D approaches (those based on the 3-D Reynolds-averaged Navier-Stokes [RANS] equations or eddy-resolving models) means that simulation of decadal to centennial river evolution over large (10 2 -10 5 m) spatial scales is currently more practical using 2-D models (e.g., Baar et al., 2019;Nicholas, 2013) than 3-D. However, compared to 3-D models, such 2-D approaches inevitably involve a higher degree of abstraction of the governing physics (Lane et al., 1999;Sloff & Mosselman, 2012). For example, 2-D depth-averaged models must parameterize, rather than represent explicitly, phenomena such as secondary circulation, the effects of topographic features (e.g., dunes and ripples), and sediment transport processes (e.g., Furbish et al., 2012;Schmeeckle, 2014) that cannot be resolved by the model mesh. Consequently, depth-averaged models necessarily adopt a highly simplified treatment of many factors that control the linkages between the flow, sediment transport, and channel change within natural rivers.
These issues can be illustrated by considering the relationships between flow and sediment transport directions as represented within a typical 2-D morphodynamic model, for the simple case of a river with a transverse bed slope perpendicular to the depth-averaged fluid flow direction (Figure 1). When calculating sediment transport directions, 2-D depth-averaged models typically account for at least two processes that promote a deviation between the mean flow and sediment transport directions: first, the deviation between the depth-averaged and near-bed flow direction (labeled αFlow in Figure 1) resulting from secondary circulation of the flow (Kleinhans et al., 2008;Struiksma et al., 1985); and second, the deviation between the near-bed flow and sediment transport directions (labeled αSed in Figure 1), which results from the deflection of sediment particles down the transverse slope of the river bed due to the effects of gravity (Antoniazza et al., 2019;Apsley & Stansby, 2008;Baar et al., 2019;Engelund, 1974;Francalanci & Solari, 2008;Ikeda, 1982;Sekine & Parker, 1992;Struiksma et al., 1985).
It has been demonstrated that river morphology and dynamics simulated using 2-D models can be highly sensitive to the parameterization of the processes outlined above (Antoniazza et al., 2019;Baar et al., 2019;Mosselman & Binh Le, 2016;Nicholas et al., 2013;Schuurman et al., 2013;Sloff & Mosselman, 2012). In addition, while existing parameterizations of transverse slope effects on sediment transport have been developed and evaluated using laboratory data sets (Baar et al., 2018;Koch & Flokstra, 1981;Sieben & Talmon, 2011), their applicability at field scales remains uncertain. This is due, in part, to the greater complexity of field environments, which are characterized by multiple scales of topography (e.g., channels, bars, and bedforms) that influence flow and sediment transport. For example, dunes and ripples are ubiquitous in sand bed rivers and exert a significant control on hydraulic roughness and hence on the channel-scale spatial distribution of flow (Karim, 1995;Rodrigues et al., 2015;Sandbach et al., 2012Sandbach et al., , 2018. Flow structures and turbulence produced by dunes are also commonly 3-D (Baar et al., 2018;Best, 2005;Lefebvre, 2019;Omidyeganeh & Piomelli, 2013a, 2013b and have been proposed to break down, or suppress, secondary flow structures usually found in meander bends and braided rivers (Abad et al., 2013;Best, 2005;Blanckaert, 2010;Dietrich et al., 1979;Drake et al., 1988;Konsoer et al., 2016;Parsons et al., 2007). Moreover, dunes give rise to complex patterns of local bed slope that may have a significant influence on the gravitational deflection of sediment (Baar et al., 2018;Sieben & Talmon, 2011;Weij, 2012). How these processes interact with, and modify, patterns of flow and sediment transport at the scale of whole channels and bars has yet to be quantified.
The extent to which existing parameterizations of flow and sediment transport processes in depth-averaged morphodynamic models (e.g., as shown in Figure 1) require revision to account for the effects discussed above remains largely unknown. Furthermore, evaluating and improving these parameterizations in field settings is problematic due to the difficulty of obtaining contemporaneous flow, topography, and sediment transport data sets of sufficient resolution and quality. The work presented herein aims to address these Figure 1. Definition of the coordinate system and processes under investigation herein for the hypothetical case of depth-averaged flow aligned parallel to the x direction (downstream) and perpendicular to a transverse bed slope. Directions x, y, and z, correspond to velocities u, v, and w. αFlow is the deviation between the depth-averaged flow and the near-bed flow direction. αSed is the deviation between the near-bed flow and the bed sediment transport direction, resulting from gravitational deflection of sediment in the direction of the transverse bed slope (tanω∂z/∂y). Adapted from Sekine and Parker (1992). issues by (i) applying 3-D computational fluid dynamics (CFD) to quantify the effect of dunes on bar-scale flow characteristics within a sand bed river, and (ii) evaluating the implications for the parameterization of sand transport on transverse bed slopes within depth-averaged morphodynamic models.

Study Area
Field measurements and numerical modeling were conducted within a 400 m long reach of the South Saskatchewan River, near the town of Outlook (51°30′N, 107°03′W), Canada ( Figure 2). The South Saskatchewan River is a braided sand bed channel with a mean bed material grain size of approximately 0.3 mm and a mean gradient of approximately 0.23 m km −1 (Cant, 1976;Parker et al., 2013;Sambrook Smith et al., 2006). The discharge that inundates the largest nonvegetated sand braid bars is 280 m 3 s −1 (Cant & Walker, 1978;Lunt et al., 2013;Sambrook Smith et al., 2010), and the bankfull discharge (above which islands and the floodplain are inundated) is approximately 1,240 m 3 s −1 . The study reach is located approximately 15 km downstream of the Gardiner Dam (51°16′12″N, 106°52′15″W), which was completed in 1967 and traps >98% of suspended sediment in Lake Diefenbaker (Wiebe & Drennan, 1973). This sediment trapping promotes relatively clear water within the study reach, even at high flow stages, which is fundamental to the acquisition of the bathymetric data described below.

Field Data Collection
A survey was undertaken in September 2016 to measure the flow and morphology of the study reach (see Figure 2c). The reach was selected to include a straight channel with well-defined bar pool topography and simple inlet and outlet conditions and to be small enough to facilitate 3-D CFD modeling at a high (~0.08 m) spatial resolution. Average water depth in the study reach was 0.8 m, discharge at the time of measurement was 40 m 3 s −1 (River discharge was~80 m 3 s −1 ), and the reach average depth mean velocities were~0.5 m s −1 . A digital elevation model (DEM) for the reach was constructed from specially commissioned aerial imagery (pixel resolution 0.06 m) captured at a height of~1,500 m from a fixedwing airplane with an UltraCamXp sensor. Full details of the procedures used in DEM construction and associated error analysis are outlined in Strick et al. (2019). In brief, a 0.06 m (horizontal) resolution DEM was constructed using a combination of large format Structure-from-Motion (SfM) photogrammetry (using the commercial software Pix4D) and a regression model between water depth and image brightness (after Lane et al., 2010). SfM techniques were applied principally in emergent areas and within a 20-40 m wide section along the true left bank where the channel bed is composed predominantly of gravel. This generated a DEM in these areas with a horizontal spatial resolution of 0.24 m. A depth-brightness model was applied to the inundated areas of the reach, where the channel is sand bedded, to produce a DEM with a horizontal resolution of 0.06 m. This depth-brightness regression model was developed and validated for the aerial imagery using water depth measurements obtained with a Navisound© NS 215 system with a Reson© TC 2024 200 kHz high-resolution dual frequency single-beam echo sounder (SBES) operating at a sampling frequency of 10 Hz, with simultaneous geographical location provided by a Leica© 1,230 Real-Time Kinematic (RTK) dGPS system. The DEMs generated using SfM and the depth-brightness model were merged along the left bank of the study reach to produce a single DEM for use in CFD mesh generation. Dune geometry was measured on this DEM via zero-crossing analysis. Average dune lengths were approximately 2.5 m with a maximum of 6 m, and average dune heights were approximately 0.17 m with a maximum of 0.30 m.
Repeat flow surveys of the channel were undertaken along cross-channel transects to provide CFD model boundary conditions and data for model evaluation. A minimum of four repeat surveys (pairs in opposite directions) were conducted to measure the mean flow structure (Szupiany et al., 2007), with T1 having six repeats. Three-dimensional flow velocities were recorded using a Sontek© M9© acoustic Doppler current profiler (aDcp) mounted onto a small zodiac boat. The position of the aDcp was also recorded using the RTK dGPS. Velocities are calculated relative to the instrument head and corrected using the dGPS-calculated boat speed and heading. The error in aDcp measurements is recorded by the instrument as the difference in measured vertical velocity from pairs of beams (Mueller & Wagner, 2009;Muste et al., 2004). The mean instantaneous error velocity for these surveys was −0.074 m s −1 , and the aDcp manufacturer reports an accuracy of ±0.05 m s −1 for velocity measurements. The aDcp data were collected in pulse coherent mode, averaging seven pings per second, and were averaged to 1 Hz to improve the signalto-noise ratio. Three-dimensional flow velocities were measured on individual track lines and gridded to a resolution of 0.1 m (vertical bin size) and 0.25 m (horizontal bin size) using the VMT aDcp postprocessing software (Parsons et al., 2013).

CFD Modeling
The open source CFD package OpenFOAM (Weller et al., 1998) was applied to solve the RANS equations for incompressible flow (Equations 1 and 2). This procedure splits the flow variables into time-averaged (e.g., U) and fluctuating (e.g., U′) components: where p is pressure, ρ is water density, and the dynamic viscosity is μ.
The terms in Equation 2 are expressed as time-averaged quantities, with the Reynolds Stress terms originating from the product of the fluctuating components. These terms are modeled by a Boussinesq approximation (written in Tensor form): where U is the 3-D velocity vector U = [u v w], X is the Cartesian coordinate vector X = [x y z], δ ij is the Kronecker delta (δ ij = 1 if i = j, otherwise δ ij = 0), and k is the turbulent kinetic energy. The eddy viscosity (v t ) in Equation 3 was modeled using the Re-Normalization Group theory (RNG) k-ε turbulence closure (Yakhot & Orszag, 1986) that has been shown to perform well in riverine applications, where high strain and flow separation are common (e.g., Bradbrook et al., 2000;Hardy et al., 2011;Marjoribanks et al., 2017;Sandbach et al., 2012). Herein, the turbulent viscosity v t (=C μ k 2 /ϵ) is computed from the turbulent kinetic energy (k) and turbulent dissipation rate (ε). These are obtained by solving their own scalar transport equations: where is the production of turbulent kinetic energy through shear, C μ = 0.0845, σ k = σ ε = 0.7194, C 1ε = 1.42, and are constants which have been determined experimentally, and η À ¼ k ffiffiffiffiffiffiffiffiffiffiffi P k =v t p =ε Á is a dimensionless parameter (Yakhot & Orszag, 1986). A uniform near-bed cell thickness of 0.0035 m with an equivalent sand grain roughness (k s ) of 0.0017 m was implemented herein to ensure appropriate dimensionless wall distance (z + ) values in the near-bed grid cell. z + is defined as where z is the vertical coordinate, v is the kinematic viscosity, and u* is the friction velocity computed assuming an equilibrium boundary condition (law of the wall). For the near bed, k ¼ u*2= ffiffiffiffiffi ffi C μ p À Á and dis- where κ is the von Karmon constant (0.41) and Δz is the distance between the bed and the center of the near bed cell. The specific wall functions used in OpenFOAM were "epsilonWallFunction" for ε, "kqRWallFunction" for k, and "nutkRoughWallFunction" for v t , with a no-slip condition on velocity and zero-normal gradient for pressure. Mean z + values for both simulations were 55, with medians of 58 for the DEM without dunes and 57 for the DEM with dunes (described below), well above the minimum value of z + = 30 and below z + = 300, which are the bounds of the wall roughness model (i.e., inside the log law layer). In the present modeling, the wall roughness function represents the net effect of both grain roughness and small bedforms (e.g., ripples) that cannot be resolved by the CFD mesh (see below).
The OpenFOAM solver "simpleFOAM" (a steady-state solver for incompressible, turbulent flow) was used to compute the pressure using the Semi-Implicit Method for Pressure-Linked Equations "SIMPLE" (Patankar & Spalding, 1972) algorithm. The free surface was represented in the model with a rigid-lid approximation (e.g., Stoesser et al., 2008;Xie et al., 2013), adopted herein as the water surface was flat with a Froude number 0.15. Inlet conditions were defined using the aDcp surveys at the same channel discharge (40 m 3 s −1 ) at T1 (Figure 2). The depth-averaged velocity entering the domain at each point across the inlet section was determined from the aDcp data and used to run a 1-D vertical boundary layer profile at each point across the inlet channel to generate profiles of horizontal velocity, turbulent kinetic energy (k), and the dissipation rate (ε) for each boundary grid cell, while maintaining the unit discharge across the inlet equal to that of the aDcp data. Vertical velocity components at the inlet to the model domain were set to 0. A Neumann pressure condition was set at the outlet, with a uniform pressure of 0. The outlet was located 5 m downstream from the zone used to interpret model results.
To provide the minimum acceptable numerical dissipation for the mesh, second-order central differencing numerical schemes were used for gradients and second-order bounded central differencing for divergence, and an unbounded second-order deferred corrected scheme was employed for the Laplacian surface normal gradients (see Jasak, 1996;& Robertson et al., 2015 for details and sensitivity tests of the OpenFOAM numerics). Convergence criteria were tested iteratively and set to a tolerance of 1 × 10 −8 for pressure and 1 × 10 −10 for velocity, k and ε. Solver tolerances were set to 1 × 10 −10 for pressure and 1 × 10 −12 for velocity. Flow variables were recorded at 128 regularly spaced locations within the computational domain at each computational time step. Using this sampling, the convergence of the simulation was tested and confirmed by a consistent lack of variation in the values of flow variables. Convergence for both simulations carried out herein occurred before the preset criteria with pressure converging first, above the solver tolerance.

CFD Mesh Generation
To investigate the effect of dunes on the 3-D flow field, and the associated implications for the parameterization of sediment transport in 2-D morphodynamic models, CFD simulations were conducted using two meshes, one with the measured channel morphology and one with the dune-scale topography removed from the DEM. To remove dunes from the original 0.06 m resolution DEM without affecting the representation of bank or bar topography, the outlines of the banks and bar fronts were identified visually in ARCMAP©, and raster masks produced to separate each bar surface for individual filtering. River banks and the lee slopes of the unit bars were therefore not modified. To remove the dunes from the individual bar tops, a moving-average, weighted-mean filter was used, with a window size of 6 m × 1 m, with the longest axis in the downstream direction. Figure 3 shows the resulting bed topographies, referred to hereafter as the Dune DEM and Smoothed DEM. Following this process, the two surface meshes were compared quantitatively to confirm that the larger bar-scale morphology of the reach had not been altered by the removal of dune-scale topography (the mean and maximum differences between the two DEMs for topography at 10 m resolution were 0.08 and 1.4 mm). Therefore, the differences in mean flow structure reported below are only due to the removal of dunes, as opposed to changes in the larger-scale channel morphology. Both DEMs were down-sampled to a resolution of 0.24 m and then used to construct 0.08 m horizontal resolution CFD surface DEMs using 2-D bilinear interpolation. This procedure produced two 0.08 m resolution DEMs, which were each adjusted to limit CFD mesh gradients, cell skewness, and nonorthogonality. This was achieved by applying a gradient filter, which altered bed elevations over approximately 4% of each model domain by approximately 3 mm. Finally, two structured finite volume CFD meshes were constructed for CFD modeling, both comprising 4,948 × 1,172 × 20 cells in the downstream, cross-stream, and vertical directions, respectively. Horizontal mesh resolution was an average of 0.08 m with 20 cells in the vertical and a planar water surface. Cells in the vertical become finer toward the bed as a power law function of depth to ensure that the majority of the vertical grid cells were close to the bed in areas of greater velocity gradient. The vertical resolution of the boundary-fitted mesh therefore varied as a function of flow depth and height above the bed. The average vertical grid cell size is 0.04 m, with a standard deviation of 0.03 m. Grid independence tests were carried out by performing three simulations using grids with horizontal resolutions of 0.03, 0.08, and 0.12 m. Comparison of global and local flow characteristics for these simulations (see details in supporting information) demonstrated model results were relatively insensitive to grid resolution. All results generated for the study reach that are presented herein were generated using 0.08 m resolution meshes. Additional model validation was performed via comparison with the laboratory data of  and , and is reported in the supporting information.

Method of Analyzing Model Results
To assess the role that dunes might play in bed sediment transport within a depth-averaged morphodynamic models framework, quantification of the effect of dunes on flow at larger scales (i.e., than a single dune), are needed. Five metrics describing the relationships between the depth-averaged flow, near-bed flow, and sediment transport directions (as illustrated in Figure 1 and discussed in section 1) were calculated. Throughout the manuscript there are three symbols used to define particular variables: θ = a direction (in degrees); δ = the difference between a variable compared between simulations, and α = the difference between variables in a single simulation.
The effect of dunes on the mean flow was evaluated by calculating δMean, defined as the deviation between the depth-averaged flow directions in the simulations using the Dune DEM and the Smoothed DEM: where θ DAF dune and θ DAF smoothed are the depth-averaged angular planform flow directions for the simulations that used the Dune DEM and the Smoothed DEM, respectively, such that positive values indicate flow rotated toward the left bank in the Dune DEM simulation relative to the Smoothed DEM simulation. Similarly, the deviation between the near-bed flow directions in the two simulations (δBed) was calculated as where θ NBF dune and θ NBF smoothed are the near-bed angular planform flow directions for the simulations that use the Dune DEM and Smoothed DEM, respectively.
To investigate the influence of dunes on bed sediment transport, in the context of depth-averaged morphodynamic model parameterizations of these processes, near-bed flow velocities derived from the model simulations were used to calculate sediment transport directions. In order to relate the results of this analysis to past work on transverse slope effects on sediment transport (e.g., Baar et al., 2018;Sekine & Parker, 1992), sediment transport directions were quantified using an approach adopted in the Delft3D morphodynamic model: where θ Sed is the angular planform sediment transport direction, S 1 and S 2 are the sediment transport direction in the x and y directions (respectively) after accounting for gravitational deflection of sediment by transverse slope effects, and S 1 ′ and S 2 ′ are the sediment transport direction in the x and y directions (respectively) before correcting for transverse slope effects. τ and τ c are the bed shear stress and critical bed shear stress, respectively, ∂z /∂n is the bed slope normal to the near-bed flow direction, and β is an empirical constant. The unit sediment transport components prior to adjusting for transverse slope effects are specified by assuming that sediment transport moves in the direction of the near-bed flow, which obviates the need to determine the magnitude of the sediment transport rate. Although the parameter β is a function of many factors including grain size and bed roughness (see Baar et al., 2018), our aim in this first instance is to estimate an upper bound for the magnitude of transverse slope effects; therefore, calculations are based on the case where sediment transport is at the threshold of motion and assume β = 1.5 (the default value in Delft3D). The implications of this choice of value are discussed further below. It should also be noted that although sediment transport directions derived using Equations 11a-11d lead to deflection of sediment in the local downslope direction, these equations do not include an explicit representation of avalanching processes in the leeside of dunes, because the complexity of such processes precludes a simple treatment of their net effect on sediment transport (Drake et al., 1988;Nield et al., 2017;Pelletier et al., 2015;Sutton et al., 2013).
The deviation (δSed) between the sediment transport directions calculated using Equations 11a-11d, for the model simulations using the Dune DEM and the Smoothed DEM, was calculated as where θSed dune and θSed smoothed are the sediment transport directions derived for the Dune DEM simulation and the Smoothed DEM simulation, respectively. The difference (δ Flow) between the near-bed flow direction (θNBF) and the depth-averaged flow direction (θDAF) in the model results generated using the Dune DEM and Smoothed DEM was calculated as This metric quantifies the degree to which the deviation between the near-bed and depth-averaged flow directions changes in the presence of dunes. An equivalent metric (δTop) was derived to quantify the extent to which dunes influence the deviation between the near-bed flow direction and the sediment transport direction (i.e., including the effects of sediment deflection due to the transverse bed slope): Depth-averaged models are typically implemented at spatial resolutions that are coarser than individual dunes. Additionally, such models must parameterize the effects of the processes that are driven by dune-scale topography, rather than represent these processes explicitly. It is therefore relevant to consider whether differences in flow (and sediment transport) that are evident in the two model simulations at fine spatial resolution are also evident at coarser spatial scales that are more typical of those used in depth-averaged morphodynamic models. To evaluate this, the CFD model results, and associated metrics defined by Equations 9-14, are examined below at two spatial scales: (i) at the resolution of the 3-D CFD modeling (referred to hereafter as the fine-resolution results), and (ii) following spatial averaging of model results at a horizontal resolution of 5 m (referred to hereafter as the spatially averaged results). The spatial resolution of 5 m was selected to be large enough to ensure that most grid cells include at least one dune, while also being small enough to resolve the bar-scale morphology within the model domain. This is also a resolution at which it would be appropriate to apply a 2-D depth-averaged morphodynamic model within the study reach.

Model Comparison With Measured Velocities
Figures 4 and 5 compare modeled flow velocity magnitudes and horizontal directions with velocities surveyed using the aDcp at Transects T2 to T4 (i.e., downstream of the model inlet at T1, see Figure 2). The agreement between modeled and observed velocity directions ( Figure 5 and Table S3) is good in most locations. For example, both models replicate the transition from diverging flow in upstream areas (e.g., at Transect T2) to flow that is predominantly from right to left in downstream areas (e.g., at Transect T4).
Previous studies have demonstrated that perfect replication of observed velocities is unlikely at cm scales due to multiple factors, including flow measurement errors (such as inherent instrument error and geolocation error), unresolved turbulence structures, differences between the CFD mesh and natural topography, and the computational constraints on obtaining true grid independence (Marjoribanks et al., 2017;Sandbach et al., 2012). Consequently, the purpose of the comparison conducted herein is to (i) assess if the model replicates the measured average flow structure; (ii) confirm that differences between measured and modeled velocities do not invalidate the model as a quantitative description of the flow (Lane & Richards, 2001;Lane et al., 2005;French, 2010); and (iii) ascertain whether the inclusion of dunes within the model mesh produces flow fields that provide an improved fit to the field data. A statistical analysis is detailed in the supporting information summaries the overall qualitative (Figures 4 and 5) and quantitative (Table S3) comparisons between the model results and aDcp data sets. The analysis indicates that both models reproduce the macroscale structure of the flow, but the inclusion of dunes into the model DEM results in improved reproduction of the observed flow field. One important aspect of the difference in flow fields using the  (Table S3 and Figure S4). Flow is through the page (in the x direction).
Dune DEM is the increased cross-channel variability in the velocity direction and magnitudes compared to the Smoothed DEM simulations. This horizontal banding has been seen in many aDcp transects and has been attributed to field measurements capturing large-scale turbulence from form-roughness, and/or a lack of repeat measurements to produce satisfactory averages (Engel & Rhoads, 2016;Parsons et al., 2005;Szupiany et al., 2012). Yet herein, the Dune DEM simulation does replicate this effect, suggesting that it is not just an artifact of field sampling but is also an inherent aspect of the effect of dunes on the timeaveraged 3-D flow structure.

Effect of Dunes on Near-Bed Flow Structure
While Figures 4 and 5 illustrate that the model results generated using the two DEMs exhibit similar macroscale patterns of velocity magnitude and direction, differences between the results of the two simulations are greatest in the near-bed regions due to the influence of bedforms. Figure 6 shows an example of near-bed morphology and flow vector patterns over a 5 m × 6 m patch of bed located immediately upstream of a bar front. Dune leesides on the top of this bar tend to be orientated oblique to the dominant downstream flow direction, such that near-bed flow directions are rotated significantly relative to the Smoothed DEM case. Consequently, near-bed flow velocities exhibit significant variations in local magnitude and direction in the presence of dunes (Figure 6a) that are almost entirely absent in the model results generated using the Smoothed DEM (Figure 6b). Such flow structures in the leesides of oblique, or 3-D dunes, have been well documented for individual dunes (Allen, 1968;Omidyeganeh & Piomelli, 2013a, 2013bSieben & Talmon, 2011;Talmon et al., 1995;Walker, 1999;Walker & Shugar, 2013) but until now have proven difficult to measure or simulate across multiple bedforms.

Journal of Geophysical Research: Earth Surface
To assess the changes in flow direction over the dunes in the study reach, individual dune crests and their orientations were digitized by hand (Figure 7a), and the corresponding near-bed flow direction downstream of the crest was measured over a distance of 1 m (see supporting information for method). These measurements were made for dunes across the whole study reach (Figure 6c) and show that only~12% of dune crestlines are orientated orthogonal to the flow direction just upstream of the dune crest (±5°), whereas~67% of crests are oblique to the flow (between ±5°and ±50°) and~20% are highly oblique to the flow (± > 50°). In the dune troughs, the distribution of flow directions is heavily skewed, with a peak in the distribution at 50°, indicating that near-bed flow in the dune troughs has been deflected toward the orientation of the dune crests but that flow parallel to the crestline orientation is uncommon. Downstream of the troughs (the yellow histogram), the frequency distribution of flow directions peaks between the crest and trough distributions, indicating flow realignment toward the flow directions found at the dune crests. This pattern is

10.1029/2020JF005571
Journal of Geophysical Research: Earth Surface similar to those found in laboratory studies by Allen (1968) and Talmon et al. (1995). Clearly, the dune leesides have a large effect on the near-bed flow structure. How this translates to flow and sediment transport at the coarser scales typically employed by 2-D morphodynamic models is investigated in the following sections.  Figure 7a) tends to be highest over the bar tops but also where channel width is constricted (e.g., near the inlet and outlet). Both banks are  Figures 7b and 7d, which indicates that in the simulation that includes the effects of dunes, the deflection of flow toward the right bank in the upstream half of the reach and the movement of flow toward the left bank in downstream areas are both enhanced. The net influence of the dunes, which appears to be to enhance lateral flows and reduce velocities over bar tops, is likely to reflect the combined effect of enhanced roughness in shallow flow (bar tops) and topographic steering of the near-bed flow (see Figure 6).  Figure 8d) is deflected toward the left bank by approximately 4-5°more in the Dune DEM simulation than in the Smoothed DEM simulation, which represents an amplification of the deviation from the downstream direction by approximately 25%. These trends indicate a potential role for dune leeside topography in steering the near-bed flow at scales larger than individual bedforms. Figure 9 is the sediment transport directions (calculated using Equations 11a-11d) for the Smoothed DEM (Figure 9a) and Dune DEM (Figure 9b) and the deviation between the two simulations ( Figure 9c). While the overall directions of sediment transport across the DEM are comparable for the two simulations, there are subtle variations in transport direction when dunes are introduced. The close visual match between the patterns of near-bed flow direction (Figure 8a) and sediment transport direction (Figures 9a and 9b) suggest that gravitational deflection of sediment has a relatively small influence on sediment transport directions. Moreover, the similarity between patterns highlighted in the areas labeled iii, iv, and v in Figure 8d (for δBed) and the corresponding areas labeled i, ii, and iii in Figure 9c (for δSed) suggests that differences in the direction of sediment transport between the Dune DEM and Smoothed DEM are driven mainly by the influence of dunes on the near-bed flow (as opposed to the differences in gravitational deflection associated with dunes).

Displayed in
To investigate the effect of 3-D flow structure of sediment transport, the magnitude of the difference between the depth-averaged and near-bed flow directions is shown in Figure 10.

10.1029/2020JF005571
Journal of Geophysical Research: Earth Surface DEM simulation are strongly affected by the local planform orientation of bedform troughs ( Figure 6) it is likely that the degree to which dunes tend to amplify or reverse the deviation between the depth-averaged and near-bed flow directions should to be determined by the local planform orientation of bedform troughs. However, no systematic trend in crest orientations and the regions labeled in Figures 8-10 was found.
To observe the significance (or lack thereof) to the gravitational deflection of sediment by bed slopes, the magnitude of the difference between the near-bed flow and the sediment transport directions for the Smoothed DEM is shown in Figure 11a and Dune DEM (Figure 11b). The difference between Figures 11a and 11b (indexed by the metric δTop) is shown in Figure 11c and demonstrates the effect of including dune morphology. The patterns of topographic deflection of bed sediment in both simulations (Figures 11a and 11b) are controlled mainly by bar-scale morphology, with narrow zones of high topographic deflection corresponding directly to the location of bar fronts. The magnitude of the gravitational effect is similar for both simulations. However, despite the introduction of greater variability in the transverse bed slope in the presence of dunes, when spatially averaged at 5 m the magnitude of the difference in topographic deflections (Figure 11b) between the two simulations (indexed by δTop) is typically <1°, with no obvious spatial pattern or consistency (Figure 11c). Therefore, even in the presence of dunes, bar-scale topography still drives gravitational deflection of sediment.

Discussion
The results outlined in section 3 show that (i) the majority (>60%) of dune orientations are oblique to their upstream flow direction ( Figure 6); (ii) there is substantial steering of flow in dune troughs, with the majority of deviations (relative to the flow just upstream of the dune crest) peaking at~+50°( Figure 6) and oriented toward the left bank; and (iii) depth-averaged velocities are up to 30% lower over bar tops and 30% higher in thalwegs (Figure 7) for simulations that represent dunes compared with simulations where dunes are neglected. Results also illustrate that dunes have significant effects (identified by spatial averaging of model results) on flow and on inferred sediment transport at spatial scales larger than individual bedforms (Figure 8). With spatial averaging employed, spatially coherent trends can be seen when comparing between simulations (δ symbol variables in Figures 8-11) and when looking at individual variables of a simulation (θ symbol variables in Figures 8-11). For example, Figures 9 and 10 identify coherent regions where flow and sediment transport directions (relative to the dominant downstream flow direction) have been reversed or amplified by the presence of dunes. In contrast, Figure 11 demonstrates that dunes may have limited effects on the gravitational deflection of sediment transport. To relate these variables together within a single framework, the following section examines the relative importance of each variable to the overall difference in sediment transport directions between the two simulations (the measure δSed).

Interpretation of Results
The results presented above illustrate significant differences in flow and sediment transport directions in the two simulations conducted using the Dune DEM and Smoothed DEM. Equations 5-10 allow quantification of the relative contributions to δSed, resulting from differences in the depth-averaged and near-bed flow and the gravitational deflection of sediment in the presence and absence of dunes. For example, Figure 12 displays the cumulative frequency distribution of the difference between sediment transport directions calculated for the Dune DEM and Smoothed DEM simulations (δSed). Distributions are also shown for (i) differences in the depth-averaged flow direction between the simulations (indexed by δMean/δSed); (ii) differences in the deviation between the near-bed flow, and the depth-averaged flow directions between the simulations (indexed by δFlow/δSed); and (iii) differences in the strength of the transverse slope effect (indexed by δTop/δSed). Each point in these distributions is derived from the mean values of the model variables calculated within a single 5 m grid cell (i.e., following spatial averaging of the high resolution CFD results). It should also be noted that δSed = δMean + δFlow + δTop and that the metrics δFlow and δTop represent the difference between the results of the simulations carried out using the Dune and Smoothed DEMs for the two components of the deviation between the depth-averaged flow and sediment transport directions illustrated in Figure 1 (i.e., αFlow and αSed).
The distributions in Figure 12 show that δSed is positive over approximately 75% of the model domain, indicating that overall, the presence of dunes tends to deflect sediment toward the left bank relative to the sediment transport direction for the Smoothed DEM. The magnitude of δFlow/δSed is approximately 3 times the magnitude of δMean/δSed over much of the distribution, although δMean/δSed does make a positive contribution to δSed on average (δMean/δSed values are positive in 66% of locations). The distribution of δTop/δSed is approximately symmetrical with a large area of the domain characterized by low values. These distributions indicate that the differences in depth-averaged flow directions between the two simulations make a

10.1029/2020JF005571
Journal of Geophysical Research: Earth Surface relatively small contribution to the differences in sediment transport directions. Differences between near-bed and depth-average flow directions is the dominant control on sediment transport directions in the Dune DEM. The overall effects of dunes on topographic deflection of sediment are relatively small, except at the extreme tails of the distribution shown in Figure 12.
To compare the results of these simulations with values that are more indicative of conditions within 2-D depth-averaged morphodynamic models (that do not represent dunes explicitly), the net deviation between the sediment transport direction and depthaveraged flow direction results for the Smoothed DEM simulation was calculated. These values are plotted on Figure 13a as "αFlow Smoothed + αSed Smoothed ," against δSed. Several previous studies have highlighted that 2-D models are sensitive to the parameterization of secondary circulation effects and gravitational deflection of sediment Schuurman et al., 2013), both of which contribute to αFlow Smoothed + αSed Smoothed . Consequently, it is pertinent to consider the extent to which dunes tend to amplify or dampen these effects.  Figure 13b illustrates the spatial distribution of data points in these four quadrants. Approximately half (56%) of the study reach is associated with data in Q1, indicating that sediment transport in the presence of dunes deflects sediment toward the left bank and in the same direction as the steering of sediment transport (relative to the depth-averaged flow) in the Smoothed DEM simulation. In contrast, approximately 35% of the study reach is associated with situations where dunes deflect sediment in the opposite direction to the steering of sediment transport (relative to the depth-averaged flow) in the Smoothed DEM simulation. Figure 13b illustrates a strong spatial coherence to these tendencies. In the area labeled i in Figure 13b, sediment transport in the Smoothed DEM simulation is steered toward the left bank, yet dunes tend to steer sediment in the opposite direction and up the bar face (i.e., Q4). In contrast, the area labeled ii in Figure 13b is characterized by strong steering of the near-bed flow and/or topographic deflection of sediment toward the right bank in the Smoothed DEM, while the presence of dunes in this area tends to deflect sediment transport in the opposite direction. Moreover, in the case of all four quadrants, the magnitude of δSed is, on average, approximately twice the value of when dunes were not present (αFlow Smoothed + αSed Smoothed ). This difference suggests that failure of 2-D morphodynamic models to represent the effects of dunes on near-bed flow and sediment transport will lead to significant and systematic differences in sediment transport direction between model and reality.
Because 2-D morphodynamic models operate at relatively coarse spatial resolutions and do not represent variations in bed topography that are associated with dunes and other alluvial bedforms, they are likely to underestimate spatial variability in flow velocities. Moreover, this will have implications for the calculation of sediment transport rates over a range of scales (e.g., Ferguson, 2003;Furbish et al., 2012). To try and quantify these effects for our present example, a simple metric of sediment transport (velocity magnitude raised to the power five), based on the Engelund and Hansen (1967) sediment transport formulae, is used. The effects of neglecting dune-scale variations in flow can be quantified at a resolution relevant to 2-D morphodynamic modeling (e.g., 5 m in the case of the present study reach) by calculating the mean value of this sediment transport rate within each 5 m grid cell using the high resolution (8 cm) model output and dividing this by the equivalent metric derived using coarse-resolution (5 m) spatially averaged velocities. The spatial  Figure 14. The order of magnitude difference in the values shown in Figure 14a compared with Figure 14c demonstrates that neglecting variations in near-bed velocity associated with dunes has the potential to substantially affect sediment transport rates calculated (Figure 14a).
The effects of neglecting dunes can also be quantified by calculating the sediment transport metric using the mean velocity values within each 5 m grid cell derived from the Dune DEM simulation, divided by the metric derived from the Smoothed DEM simulation (Figures 14b and 14d). As 2-D morphodynamic models use sediment transport formulae (e.g., Engelund-Hansen) based on depth-averaged flow variables, on which the dune-scale topography has a more limited effect (e.g., in the 3-D simulations reported herein), the neglect of dune-scale topography has little effect on sediment transport rates (when measured using the simple metric employed here-see Figure 14c). In the context of 2-D morphodynamic modeling, the main effect of neglecting dunes on this sediment transport rate is associated with redistribution of the flow, which reduces velocity over bar tops and enhances flow and sediment transport in thalwegs and along the channel margins. This effect is equally significant in the context of both the depth-averaged and near-bed flow (Figure 14b vs. Figure 14d). . This scatter plot illustrates the magnitude of the net effect of dunes compared to magnitude of the deviation between the depth-averaged flow and sediment transport directions in the absence of dunes. Data are segmented into four quadrants. Points in Q1 and Q3 indicate areas where dunes deflect sediment transport in the same direction as the deviation between the mean flow and sediment transport in the absence of dunes (i.e., the addition of dunes can be inferred to amplify tendencies evident in the smoothed DEM simulation. (b) The planform distribution of data points for these four quadrants. Areas i and ii highlight regions of spatially consistent Q2, along a bar edge, and Q4, along a bar top near the left bank. Q1 dominates the latter third of the DEM, where the presence of dunes has amplified the deflection of flow from true right to left.

Implications
The results presented herein demonstrate that dunes can either amplify or suppress the deflection of sediment via steering of the near-bed flow, as a result of which it may be difficult to produce a simple method of accounting for dunes in 2-D morphodynamic models. Deviations between the depth-averaged and near-bed flow are incorporated in depth-averaged models as secondary circulation effects, often modeled as a function of streamline curvature (Struiksma et al., 1985). However, the simple (straight) morphology of the current study reach results in low planform curvature so that values of αFlow for the simulations reported herein are driven by topographic steering of the near-bed flow (e.g., by bar fronts or dunes aligned oblique to the mean flow). This may indicate that, particularly in the absence of strong channel planform Figure 14. Ratios of sediment transport rates. To assess how spatial variability of flow could affect sediment transport rates, near-bed velocity (a, b) and depth-averaged velocity (c, d) magnitudes were raised to the power 5 (Engelund & Hansen, 1967). For Panels (a) and (c) this sediment transport rate was calculated at the 8 cm CFD grid resolution, then averaged to 5 m, and divided by a sediment transport rate calculated from the corresponding 5 m averaged velocity magnitude. This comparison demonstrates an order of magnitude difference in sediment transport rate when using near-bed flow velocities (a) compared to using depth-averaged velocities (c). Panels (b) and (d) have the same calculation of sediment transport rate at the CFD resolution, but are divided by the equivalent metric from smooth DEM simulation. Panels (b) and (d) demonstrate similar magnitudes and spatial distributions, with the depth-averaged case (d) smoother than the near-bed flow case (b). Panels (b) and (d) therefore demonstrate that redistribution of the flow away from bar tops to thalwegs is equally significant for near-bed (b) and depth-averaged (d) flow.
curvature, there is a need to develop new parameterizations for use within depth-averaged models of the effects of near-bed topography on the vertical distribution of flow velocity.
There are approaches for predicting the spatial distribution of dune height within depth-averaged morphodynamic models (e.g., Giri et al., 2008). The results presented herein suggest that dune height (or more pertinently, roughness) will affect the spatial distribution of sediment transport rates through their effect on the mean flow (Figure 14c and 14d). When dunes are included, the majority of the differences in sediment transport directions is due to the effects of dunes on steering of near-bed flow (δFlow, Figure 12). To account for these effects, new approaches for predicting dune orientation are needed that will provide a basis for deriving improved parameterizations of the deviation between the depth-averaged and near-bed flow. However, this task is difficult to accomplish a priori (e.g., Weij, 2012), as dune orientation is often a product of the interplay between bed slope, spatially nonhomogeneous sediment transport rates (Rubin, 2012;Terwisscha van Scheltinga et al., 2020;Unsworth et al., 2016;Weij, 2012), and the interaction that dunes have on the mean flow field.
The relatively small contribution of gravitational deflection of sediment by local (dune-scale) bed slopes to the differences between sediment transport directions for the Dune DEM and Smoothed DEM simulations is noteworthy, particularly given that the analysis conducted using Equations 11a-11d was intended to maximize the potential magnitude of this effect. It is important to clarify that sediment transport directions estimated herein (using Equations 11a-11d) are only for bedload transport. While dunes are nearly always responsible for suspending some sand (Kuhnle & Wren, 2009;Naqshband et al., 2014;Wren et al., 2007), the effect of gravitational deflection (measured by δTop) on suspended material remains highly uncertain (Talmon & Wiesemann, 2006). Moreover, the conditions considered herein are more representative of the critical threshold of motion, since they evaluate the maximum likely contribution of transverse slope effects. At higher flow stages, and where a significant proportion of sediment moves in suspension, the effects of gravitational deflection are likely to be weaker than those reported herein (cf. Baar et al., 2018). This suggests that while parameterization of gravitational deflection due to transverse slope effects within depth-averaged models remains a key control on simulated channel morphology (Baar et al., 2019;Nicholas, 2013;Schuurman et al., 2013), improved parameterization of these effects may not necessitate explicit representation of the effects of dunes on δTop within such models. However, this conclusion may reflect, in part, the nature of the analysis conducted herein, which has been designed to represent the calculations performed by depth-averaged morphodynamics models, rather than the physical controls on sediment transport at fine (subdecimeter) spatial and temporal scales. Similar conclusions can be drawn in relation to the analysis presented in Figure 14 that indicates that the neglect of spatial variations in velocity associated with dunes will have a limited effect on sediment transport rates calculated in 2-D morphodynamic models, where these transport rates are based on formulae that utilize depth-averaged velocities. Full account within such models of the potential effects of dunes on sediment transport would necessitate the development of theory that accounts for the controls on near-bed turbulent sediment transport (Furbish et al., 2012) and its relationship to the depth-averaged flow and local bed morphology.
One of the primary limitations of the present analysis is that the sediment transport directions calculated herein do not include the effects of avalanching on the lee faces of bars or dunes. Flow reworking of sediment on the lee face is approximated by the present analysis but notably does not include the important contribution of turbulence-derived sediment transport events in the leeside (Nabi et al., 2015;Schmeeckle, 2014). Lee face reworking by the flow can form a substantial proportion of the sediment transport directions along the slip face, especially when the bedform crests and troughs are oblique to the mean flow direction (Sieben & Talmon, 2011;Terwisscha van Scheltinga et al., 2020;Unsworth et al., 2016). The contribution of lee face reworking has mostly only been quantified in the field for eolian environments (e.g., Nield et al., 2017;Pelletier et al., 2015) and has not yet received detailed field study in aqueous unidirectional flows where there are 2-D and 3-D bedforms migrating over lateral slopes. Cisneros et al. (2020) recently demonstrated that dune leesides in many of the world's larger rivers tend to have low-angle lee slopes. Such low-angle dunes produce smaller flow separation zones or no permanent flow separation at all (Best & Kostaschuk, 2002; resulting in reduced form drag (Kwoll et al., 2016). Crestline obliqueness effects on near-bed velocity directions is produced in the present work, but how leeside flow deflection is then changed with a lower leeside slope angles or composite leesides commonly found in larger rivers is an avenue of future research.

Conclusions
This paper reports on the application of a 3-D CFD model to simulate time-averaged flow patterns within a 400 m reach of the sand-bedded South Saskatchewan River, Canada. Two simulations were conducted using high-resolution DEMs with and without dunes. The differences between these simulations are discussed in relation to 2-D morphodynamic modeling. The main findings are the following: 1. The presence of dunes promotes redistribution of flow away from shallow bar tops toward the channel margins. This effect reduces bar top velocities by up to 30% and promotes a concomitant increase in bank edge velocities. 2. For depth-averaged flow, dunes amplify the dominant patterns of cross-stream flow induced by the bar-scale channel morphology. 3. Dunes promote substantial local variability in near-bed flow directions, even when averaged up to the scale at which 2-D morphodynamic models operate. In many areas, near-bed flow is steered by dune troughs which are oriented oblique to the downstream flow direction; this leads to differences in near-bed flow directions between simulations with and without dunes that exceed 45°in some locations. 4. Spatial averaging of high-resolution model results indicates that the effects of dunes on flow and sediment transport directions are evident at spatial scales coarser than the largest dunes and up to the scale of individual channel-scale bar forms. Parameterization of this subgrid-scale effect for morphodynamic models requires prediction of dune crestline orientation. 5. Significant differences between depth-averaged and near-bed flow directions are evident in simulations that use a Smoothed DEM (i.e., with dune bedforms removed). These differences reflect topographic steering of the flow (e.g., by bar fronts) in the absence of significant planform curvature, rather than the effects of secondary circulation driven by planform curvature. 6. Sediment transport patterns inferred using theory implemented in depth-averaged morphodynamic models suggests that gravitational deflection of sediment is controlled by bar-scale topography. Dune-scale topography promotes local variability in such gravitational deflection, but this effect is not associated with a coherent signal at spatial scales larger than individual dunes. The dominant effect of dunes on sediment transport directions is through their influence on near-bed flow directions, with a secondary effect linked to their impact on the depth-averaged flow.
These results have implications for the parameterization of flow and sediment transport processes in depth-averaged morphodynamic models. Specifically, they indicate that secondary circulation schemes based on streamline curvature may fail to represent controls on the deviation between depth-averaged and near-bed flow directions. Improvements in the realism of depth-averaged flow models may only be achieved through the development of new approaches for predicting the spatial distribution of bedform height and orientation.

Data Availability Statement
Project data are stored in, and available from, the UK Centre for Ecology and Hydrology (http://eidc. ceh.ac.uk).