12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022

DIRECT NUMERICAL SIMULATIONS OF HYPORHEIC EXCHANGE
INDUCED BY ROUGH BEDFORMS

Guangchen Shen
Department of Mechanical Engineering
Michigan State University
428 S. Shaw Ln., East Lansing, MI, 48824
shenguan@msu.edu

Junlin Yuan
Department of Mechanical Engineering
Michigan State University
428 S. Shaw Ln., East Lansing, MI, 48824
junlin@egr.msu.edu

Mantha S. Phanikumar
Department of Civil and Environmental Engineering
Michigan State University
428 S. Shaw Ln., East Lansing, MI, 48824
phani@egr.msu.edu

ABSTRACT

same approach as for impermeable walls. From the literatureof impermeable-wall turbulence , wall roughness is known toenhance the momentum transfer toward the wall and dynamically modify the near-wall flow (Jiménez, 2004). The effects
become more important in non-equilibrium flows (Yuan & Piomelli, 2015; Mangavelli et al., 2021) subject to strongly varying longitudinal pressure gradients in space and time, which
can be induced by the presence of bedforms.
Some insights into how small-scale features of bed topography may play a role in the exchange were provided by a few
studies on multiscale or fractal bed topographies (Aubeneau
et al., 2015; Lee et al., 2020). They showed that the interfacial flux and bed storage (i.e., retention of surface flow in thebed) increase with small-scale bed features. Our earlier work
(Shen et al., 2020, 2022) based on pore-resolved DNS revealed
that, on a macroscopically flat bed (without bedform), grainscale roughness yields significant volumetric flux into the sediment. It also induces significant subsurface flow and storageas a result of multiscale pressure variations at the sedimentwater interface (SWI), a mechanism termed “advective pumping” by interfacial pressure, first observed for bedforms (Elliott& Brooks, 1997b). Given the multiscale nature of a realistic
stream bed, however, it is not clear whether the roughness effects described above are still significant in the presence of a
bedform.
The objective of this work is to identify to what extent
the grain-scale bed roughness affects the exchange of water
between the surface and subsurface flows induced by a duneshape bedform. To this end, pore-resolved DNS simulationsof turbulent half-channel flows bounded by a dune-shaped permeable bed with two different bed roughnesses are carried out
and compared.

To understand how different bed morphological scales affect stream hyporheic exchange, pore-resolved direct numerical simulations (DNS) of half channels at a Reτ of 1580 on a
porous symmetric dune with two different grain-scale roughnesses are conducted. The permeability Reynolds number is
ReK = 2.5, representative of streamflows over sandy beds. We
show that the dune-shaped bedform leads to large-scale variations of the interfacial pressure, while roughness superimposedon the bedform modifies the pressure variation, as well as thewall friction and shear penetration depth. Characteristics ofthe subsurface flow are more sensitive to the roughness texture when a bedform is present than when it is not. The results
highlight the importance of nonlinear interactions between theeffects of roughness and those due to bed morphological features at larger scales. DNS such as the ones reported hereincan be used to directly characterize the pore-scale dynamics
and to gain insights that could improve pore-unresolved models of hyporheic mixing and biogeochemical processes.

INTRODUCTION

The hyporheic zone, the region in and around a stream
channel where groundwater and stream water mix, plays animportant role in biogeochemical processes in stream ecosystems (Boano et al., 2014). Streamflows can be idealized as turbulent open channel flows bounded by permeable beds madeof packed sediment grains on the sides and the bottom. In
natural streamflows, the bed surface morphologies are multiscale and often self-affine fractals, consisting of large-scale
bedforms (e.g. ripples and dunes that are orders of magnitudelarger than individual grains) and smaller-scale features, down
to the bed roughness at the scale of a grain.
Bedforms can induce significant form resistance, production of turbulence and surface transport. Most analyses
and models on the role of bed morphology on the surfacesubsurface exchange focused on the scale of the bedforms.
Knowledge on the effect of small-scale bed features down
to the grain roughness is very limited. The bed roughness
has been either ignored or modeled heuristically based on the

PROBLEM FORMULATION
Governing Equations

We solve the equations of conservation of mass and momentum governing the incompressible flow of a Newtonian1

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022
Table 1. Summary of parameters. D is the grain diameter; Lxi is the simulation domain size in xi ; ReK = 2.5, Reτ = 1580, D+ = 79.
∆x+ , ∆y+ , and ∆z+ are DNS grid sizes in x, y and z, respectively, normalized using the viscous length scale ν/uτ .

Roughness

D/δ

Hs /δ

H/δ

λ /δ

Lc /λ

(Lx , Lz )/δ

+(∆x+ , ∆y+min , ∆z )Case 1

Regular

0.05

2

0.15

3

0.5

(3, 1)(3.7, 0.3, 4.1)

Case 2

Random

0.05

2

0.15

3

0.5

(3, 1)(3.7, 0.3, 4.1)

Figure 1. (a) Simulation domain and synthesized bedforms with either regular (Case 1) or random (Case 2) bed roughnesses shownin (c) colored by y. (b) Streamwise power spectral density of local roughness height fluctuations.
fluid:
∂ ui
= 0,
∂ xi
∂ u j ∂ ui u j∂P+=−+ ν$2 u j + Fj .
∂t
∂ xi
∂xj

(1)
(2)

Cases

xc /D

yc /D

θ

L/D

CD

Present

0.72

0.26

47.4

0.71

1.32

Mittal (1999)

−

−

47.2

0.72

1.26

48.6

−

1.26

Taneda (1956)

0.68

0.25

Table 2. Validation of spatial resolution of grains by comparing various characteristics of a uniform flow past a single
sphere at Re = 75 with results from DNS of Mittal (1999) andexperiments of Taneda (1956).

Here, x1 , x2 and x3 (or x, y and z) are, respectively, the streamwise, vertical and spanwise directions, and u j (or u, v andw) are the velocity components in those directions; P = p/ρ
is the modified pressure, where p is the static pressure andρ the density. The term Fj is a body force imposed by an
immersed boundary method (Yuan & Piomelli, 2014) to impose no-slip boundary conditions on the fluid-solid interface
in a Cartesian grid. Details of the Navier-Stokes simulations
were described by Shen et al. (2020). Symmetric boundary conditions are applied at both top and bottom boundariesof the simulation domain. Periodic conditions are applied
at x and z boundaries. A constant mean pressure gradientis used to drive the flow. A double-averaging (DA) decomposition of an instantaneous flow variable φ to differentiate
the turbulent fluctuation (φ ′ ) from the spatial fluctuation of
the time-mean value (φ!) (Raupach & Shaw, 1982) is applied:φ (!x,t) = 〈φ 〉(x, y) + φ!(!x) + φ ′ (!x,t), where φ is the temporal
average and 〈φ 〉 is the intrinsic line average along z.

ness (Case 1), in which the arrangement of grains in theuppermost layer is regular, and the other with random arrangements of uppermost-layer grains (Case 2). In bothcases, the same porosity of 0.4 is imposed in the bulk of
the sediment. The friction Reynolds number and permeability Reynolds number are √the same for both cases: Reτ =δ uτ /ν = 1580 and ReK = Kuτ /ν = 2.5, where K is the permeability. The grain diameter in wall units is D+ = Duτ /ν =79. The values of ReK and porosity match those in experimental measurements with a flat sediment bed by Voermans
et al. (2017). D+ is the same as the values used in Shenet al. (2020, 2022), which matched the values used by Voermans et al. (2017) in Case L12 therein. Here, the friction velocity uτ is obtained from the streamwise average ofthe peak magnitudes of the total shear stress profiles (sum of
the viscous, Reynolds and dispersive shear stresses): uτ ="$ %1/2! #(1/Lx ) oLx ν∂ 〈ū〉/∂ y − 〈u′ v′ 〉 − 〈!
u v!〉 dx
, where Lx isthe domain size in the x-direction. The number of grid points
are 1280, 1147 and 384 in x, y and z, respectively. The grain
geometry is resolved by at least 20 grid points in each direction. The y mesh is refined in the region bounded by the troughand the peak of the dune. The total simulation time used
for data collection is around 10 large-eddy turn-over times(LETOTs, defined as δ /uτ ).

Parameters

The simulation domain includes both the surface flow andsub-surface flow regions (Figure 1(a)). The elevation of y = 0
is set at the bedform trough elevation. δ is the half-channelheight measured from y = 0; Hs = 2δ is the sediment depth
measured downward from y = 0 to the bottom boundary of thesimulation domain. H = 0.15δ is the bedform height measured
from y = 0 to the dune crest. The bedform wavelength is λ =3δ . A symmetric dune is studied, with Lc = λ /2. The bedformgeometry matches that studied by Cardenas & Wilson (2007).
Details of the simulations and parameters are listed inTable 1. Two cases are considered, one with regular rough-2

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022

Figure 2. Profiles of time- and spanwise-averaged streamwise velocity:
Case 1,Case 2.
Local bedformheight.
To evaluate the adequacy of the grid size in resolving thepore flow, additional simulations of a single sphere, resolved
by 25 grid points in each direction, in a uniform unboundedflow was conducted at a Reynolds number UD/ν (where U is
the uniform freestream velocity and D is the sphere diameter)of 75, comparable to the local Reynolds number at the SWI
represented by D+ = 79. The simulation setup is similar tothat of Mittal (1999). Results are compared to both those of
Mittal (1999) and the experimental measurements of Taneda(1956) in Table 2. The results are 1-6% different for the flow
characteristics (separation angle θ , recirculation bubble size L
and recirculation center location (xc , yc )) and 5% for the dragcoefficient CD . The comparison gives confidence that the poreflow herein, especially the local mean shear layers and dragforce, are overall well captured by the spatial resolution.
The porous bedforms are immobile and modeled as
closely packed mono-disperse hard spheres. The locations of
the grains are determined based on molecular dynamics simulations (see details in Shen et al. (2020)). Two bedforms
with the same macroscopic dune geometry but different roughnesses at the uppermost layer are synthesized. One is the “regular” case (Case 1), formed by regular distribution in (x, z)of uppermost-layer grains, and the other is the “random” case
(Case 2) formed by random grain distribution of grains in both(x, z) and y. Statistically, both roughnesses are similar to those
imposed on a flat bed in Shen et al. (2020, 2022). Specifically,
the random roughness yields a larger Taylor microscale of theroughness height and a larger root-mean-square height (Shenet al., 2020).
Figures 1(c) compare the bed surfaces in Cases 1 and 2,
showing the differences in grain arrangement at the top of the
bedform. To quantify the differences between the two roughness geometries, the one-dimensional power spectra (Ek ) of
local roughness heights measured from the bedform surfaces
are compared in Figures 1(b) as functions of the streamwise
wavenumber (κ1 = 1/λ1 , where λ1 is the streamwise wavelength). As expected, similar spectral distributions as thoseseen on a macroscopically flat bed in Shen et al. (2022) are
observed here, as the same types of roughness are imposed.
Below the bed surface, a random distribution of grains is usedin both cases.

Figure 3. Profiles of (a) viscous and (b) Reynolds shear
stresses:
Case 1,Case 2.
Local bedformheight.

Figure 4. Streamwise variations of of the Clauser parameter:◦ Case 1, △ Case 2.

Velocity Statistics

The y profiles of the time- and spanwise-averaged streamwise velocity normalized by uτ are compared in Figure 2.
They show a thinner shear layer on the bed as the flow is accelerated toward the bedform crest, and a thicker shear layer
during flow deceleration after the crest. The random roughness leads to a lower bulk velocity in the half channel compared to the regular one and, consequently, results in a 17%
higher friction coefficient, C f = 2(uτ /Uδ )2 (where Uδ is themean velocity at the channel centerline), of 0.0109 compared
to 0.0093 in the regular-roughness case. The higher C f is dueto the larger effective length scale of the random roughness
enhancing the local wall friction on the bedform surface, consistent with the observations made on a macroscopically flat
bed by Shen et al. (2020).
Profiles of the viscous shear stress normalized by local
wall units (uτ and viscous length scale ν/uτ ) are shown in
Figure 3(a). Despite the flow deceleration on the lee side anda sharp edge of the bedform at its crest, no mean flow separation (or change of sign of the y-gradient of streamwise mean
velocity) is observed. This may be because that a permeable
wall promotes inner-outer interaction in a turbulent boundary
layer and, consequently, reduces the likelihood of pressuregradient-induced flow separation. The Reynolds shear stressprofiles are compared between the two cases in Figure 3(b).
They demonstrate deeper penetration of turbulence into thebed on the random roughness than the regular one, indicating
dependence of wall-normal mean momentum transfer on the

RESULTS

Results for Cases 1 and 2 are compared to understand how
different bed roughnesses may affect the SWI turbulence andthe subsurface flow in the presence of a bedform.

312th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022

Figure 5. Contours of time- and spanwise-averaged modified pressure (P) normalized by u2τ , in Case 1 (a) and Case 2 (b).
Bedform height. Streamwise distributions of interfacial pressure differences (with respect to the pressure value at dune trough, x = 0),
normalized based on friction velocity (c) or bulk velocity (d).
When normalized using u2τ (Figure 5 (c)), the pressure variation induced by the random roughness is significantly weaker(by almost a half in magnitude) than that generated by the regular one. However, the streamwise profiles almost collapse
when the pressure drop is normalized using u2b instead (Figure 5 (d)).
To explain the approximate scaling of interfacial pressure
variation on u2b , assume that the half-channel flow is steadyand inviscid, to which the Bernoulli’s equation applies along a
streamline (ignoring the hydrostatic pressure). Together withthe conservation of mass, one obtains: ∆P(x) = P(x) − Po =(1/2)Vo2 [1 − Ao /A(x)], where the subscript “o” denotes thereference point (i.e. x = 0), and V and A are the velocity magnitude and the cross-sectional fluid area of the half channel,
respectively. It follows that, for a given dune geometry (i.e. agiven Ao /A(x) ratio), ∆P varies as a quadratic function of velocity in the region of approximatley inviscid flow, which canbe represented using ub . Notice in Figure 5 (c) that the scaling based on ub does not collapse perfectly the two cases, as
the Bernoulli’s equation does not apply exactly. Nevertheless,these observations imply that the difference in pressure variation at the SWI between the two cases is mainly due to thehigher wall resistance (as shown by C f ) in the random caseleading to a reduced mass flux in the surface flow. According
to the pumping mechanism (Elliott & Brooks, 1997b), a difference in interfacial pressure is expected to change the subsurface flow, which is discussed next.

roughness details. On the lee side, different from the viscousshear stress which reaches the maximum values on the bed surface, the maximum magnitudes of the Reynolds shear stress
are reached consistently at the crest elevation of the dune,
probably due to strong TKE production near the dune crestand subsequent convection of new turbulence downstream.

Pressure Variations

To quantify the strength of half-channel flow accelerationand deceleration due to the bedform, the Clauser parameterβ (x) = (δ ∗ /|τ|max ) (d pδ /dx) is calculated and shown in Figure 4. Here, d pδ /dx is the streamwise gradient of the staticpressure at the channel centerline, δ ∗ (x) is the local displacement thickness, and |τ|max (x) is the local maximum total shear
stress magnitude at an x location. In both cases, β takes negative values (indicting favorable pressure gradients, FPG) on thestoss side and positive ones (adverse pressure gradient, APG)
on the lee side. The regular roughness (Case 1) yields larger
β magnitudes in both signs than does the random roughness(Case 2), due to both a stronger centerline pressure gradient(shown later) and a lower wall friction in Case 1. The β in the
APG region takes values up to around 2, indicating relatively
weak pressure gradient induced by the bedform. To investigate
the effect of wall permeability on APG flow and separation in
the context of hyporheic exchange, one needs to investigate
other bedform shapes, such as asymmetric bedforms, that generate stronger pressure gradients.
The distribution of the time- and spanwise-averaged pressure across the half channel and inside the permeable bed iscompared in Figure 5 (a,b). As expected, the bedform introduces negative streamwise pressure gradients on the stoss sideand positive gradients on the lee side. Similar to what was observed by Cardenas & Wilson (2007) using Reynolds-averaged
Navier-Stokes simulations of surface flows only, the minimumpressure is located at the crest of the bedform.
The pressure distribution along the SWI drives the flowin the sediment and controls the surface-subsurface transport.
The time-averaged interfacial pressure (measured using ∆P defined as the difference of local interfacial pressure from its
value at the dune trough) non-dimensionalized using either the
average bulk velocity ub or uτ is shown in Figures 5 (c) and (d).

Subsurface Particle Transit Times (τ )The role of hyporheic zones (i.e., subsurface flows) intransforming the chemical signature of stream water depends,
among other factors, on the exchange fluxes of water and solutes (nutrients, dissolved organic carbon and oxygen), as wellas the transport time scales which are often described using
distributions of transit times (Harman, 2015). Here, we focus
on the transit time, which is defined as the time spent by a parcel of water in the sediment, after its entry into the sedimentand before its exit from the sediment back to the surface flow.
The transit time is calculated based on a particle-tracking
method, which has been used with pore-scale simulation data
in a few studies (for example, Kim & Kang (2020)). The

4

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022

Figure 8. Probability density functions of transit times.
Fitted power law with slope indicated.
Table 3. The 50th and 95th percentiles of transit time distributions. All times are normalized by δ /uτ .
Figure 6.

Interface

τ50

τ95

Case 1

Regular3

385

Case 2

Random

0.09

7.47

Subsurface flow paths in Cases 1 (a) and 2 (b).

both cases, multiscale subsurface flow paths are induced, withthe longest ones characterized by roughly a half of the bedform length scale λ . The general pattern is consistent withobservations in existing numerical and experimental data ofsubsurface flows induced by bedforms (Cardenas & Wilson,
2007; Elliott & Brooks, 1997a). Here, the regular roughness isshown to induce a larger fraction of deep-reaching paths and,
consequently, overall longer transit times. This is probably due
to the higher-magnitude interfacial pressure variation observedin Figure 5 for this case. Figure 7 shows the probability density function (PDF) of the maximum depths (measured using
the minimum y value along a flow path) reached by the subsurface paths; the integral area equals one. It is evident that alarger fraction of the paths in the regular case reaches beneaththe trough of the bedform (i.e., y = 0) than in the random case.
−
The backward transit-time distributions, ←p Q (τ), definedas the PDF of parcel travel times in the bed (Harman, 2015)
are calculated based on the tracked paths for each case andcompared in Figure 8. The transit time values are normalizedusing δ /uτ . The distributions in both cases display powerlaw decays with a slope of around −1.5 at large transit times.
The existence of a heavy tail of this kind is consistent withexperimental observations by Aubeneau et al. (2015) and results of numerical simulations obtained by Lee et al. (2020)of exchanges induced by multiscale bed geometries, and also
indicates that the sediment domain in the DNS simulations
is sufficiently large (in both x and y) to capture deeper flow
paths. Compared to our earlier flat-bed results (Shen et al.,2020) which showed a power-law slope of around −1.3, theinclusion of a bedform does not appear to significantly modify
the probability of paths with very long transit times. However,
the addition of a bedform leads to a difference between thetwo cases that is opposite from the flat-bed results: the regular roughness produces higher probabilities of flow paths withvery long transit times compared to the random one, while results without the bedform Shen et al. (2020) showed that the

Figure 7. Probability density functions of the maximum
depths reached by individual flow paths (quantified using thelowest y values reached).
work of Shen et al. (2022) and the present work represent thefirst attempts in using the particle tracking approach based on
turbulent flow DNS. Specifically, a large number of tracked
fluid parcels are seeded uniformly on the bedform surface
and traced into the bed, yielding a large collection of threedimensional (3D) subsurface flow paths for each case. The
tracking is based on the convection of the 3D time-mean velocity only: dxi /dt = ui [xi (t)], where xi (t) is the parcel location and t is a fictitious time used for the tracking. Once a
subsurface flow path is determined, the transit time (denotedby τ) for the corresponding fluid parcel is calculated as thetotal time it spends along the path. In total, around o(106 )parcels are released for each case. The advection by instantaneous turbulent fluctuations u′i is not accounted for. This is because, at the present ReK of o(1) representative of sandy riverbeds, the magnitudes of normal Reynolds stress components
were shown to be significantly smaller than those of the forminduced stresses in the sediment below the Brinkman layer,
which is very thin (20% to 30% of grain diameter) (Shen et al.,2020). At a higher ReK range, however, the effect of u′i inside
the sediment on the transit times is likely to be more importantand the turbulent advection of parcels should not be neglected.
Figure 8 compares the tracked subsurface flow paths. For

5

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022
random roughness yields more long-transit-time paths instead.
The differences in transit time distributions between thetwo cases are quantified in Table 3, using the 50th and 95th
percentiles (denoted as τ50 and τ95 , respectively) of the transittime distributions. Here, τ50 represents the central tendency ofthe distribution, while τ95 is considered as a characteristic timeassociated with the rare, long-transit-time paths. The comparison shows that both indicators in the regular case are two orders of magnitude higher than in the random case. An important conclusion that follows is that the difference in grain-scaleroughness texture on top of a bedform affects the flow deep inthe sediment and significantly modifies the transit time distributions.

morphology produce fractal scaling of water storage times.
Geophys. Res. Lett. 42, 5309–5315.
Boano, F., Harvey, J. W., Marion, A., Packman, A. I., Revelli, R., Ridolfi, L. & Worman, A. 2014 Hyporheic flow
and transport processes: Mechanisms, models, and biogeochemical implications. Reviews of Geophysics 52, 603–679.
Cardenas, M. B. & Wilson, J. L. 2007 Dunes, turbulent eddies,and interfacial exchange with permeable sediments. Water
Resour. Res. 43, W08412–1–16.
Elliott, A. H. & Brooks, H. 1997a Transfer of nonsorbing solutes to a streambed with bed forms: Laboratory experiments. Water Resour. Res. 33, 137–151.
Elliott, A. L. & Brooks, H. 1997b Transfer of nonsorbing solutesto a streambedwith bed forms: Theory. Water Resour.
Res. 33, 123–136.
Harman, C. J. 2015 Time-variable transit time distributions
and transport: Theory and application to storage-dependenttransport of chloride in a watershed. Water Resour. Res. 51,
1–30.
Jiménez, J. 2004 Turbulent flows over rough walls. Annu. Rev.
Fluid Mech. 36, 173–196.
Kim, J. S. & Kang, P. K. 2020 Anomalous transport throughfree-flow-porous media interface: Pore-scale simulationand predictive modeling. Advances in Water Resources 135,
103467.
Lee, A., Aubeneau, A. F. & Cardenas, M. B. 2020 The
sensitivity of hyporheic exchange to fractal properties ofriverbeds. Water Resour. Res. 56, e2019WR026560–1–15.
Mangavelli, S. C., Yuan, J. & Brereton, G. J. 2021 Effects ofsurface roughness topography in transient channel flows. J.
Turbul. 22, 434–460.
Mittal, R. 1999 A Fourier–Chebyshev spectral collocation
method for simulating flow past spheres and spheroids. International Journal for Numerical Methods in Fluids 30 (7),921–937.
Raupach, M. R. & Shaw, R. H. 1982 Averaging procedures forflow within vegetation canopies. Bound.-Lay. Meteorol. 22,
79–90.
Shen, G., Yuan, J. & Phanikumar, M. S. 2020 Direct numerical simulations of turbulence and hyporheic mixing nearsediment–water interfaces. J. Fluid Mech. 892, A20.
Shen, G., Yuan, J. & Phanikumar, M. S. 2022 Quantifying theeffects of bed roughness on transit time distributions via direct numerical simulations of turbulent hyporheic exchange.
Water Resour. Res. 58, e2021WR030503–1—16.
Taneda, S. 1956 Experimental investigation of the wake behinda sphere at low Reynolds numbers. Journal of the PhysicalSociety of Japan 11 (10), 1104–1108.
Voermans, J. J., Ghisalberti, M. & Ivey, G. N. 2017 The variation of flow and turbulence across the sediment–water interface. J. Fluid Mech. 824, 413–437.
Yuan, J. & Piomelli, U. 2014 Numerical simulations of sinkflow boundary layers over rough surfaces. Phys. Fluids
26 (1), 015113.
Yuan, J. & Piomelli, U. 2015 Numerical simulation of a spatially developing accelerating boundary layer over roughness. J. Fluid Mech. 780, 192–214.

CONCLUSIONS

Pore-resolved DNS simulations of half-channel flows
bounded by immobile permeable dune-shaped walls are carried out at ReK = 2.5 and Reτ = 1580, to characterize the effects of two bed roughnesses, regular and random, on the interfacial and subsurface flows. The mean velocity profiles showthat the roughness modifies the flow velocity near the SWI, resulting in a difference in the friction coefficient depending on
the roughness texture. The random roughness yields a higher
friction coefficient and deeper shear penetration than the regular roughness, consistent with observations made in the absence of the bedform (Shen et al., 2022).
The presence of bedform is shown to induce static pressure variations along the bedform surface due to the acceleration and the subsequent deceleration of the half-channel
flow, consistent with observations from a large body of workon permeable and impermeable bedforms. Surprisingly, theroughness textures modify significantly the interfacial pressure
variation induced by the bedform: a larger pressure variation
is observed on the bedform with regular roughness than thatwith the random one. This is contrary to the observations
made without the bedform where the random roughness induces more intense pressure variations at the roughness lengthscales. This difference in hyporheic exchange between the flatbed and bedform cases is attributed to the different roles playedby the roughness: on a macroscopically flat bed, roughness
protuberances create roughness-scale pressure heterogeneity,
while on a bedform, the roughness leads to higher wall friction, reducing the mass flux through the half channel and consequently the pressure variations.
Correspondingly, contrary to the flat-bed flows where therandom roughness augments mass storage in the bed by inducing longer and deeper flow paths and longer transit times, withthe presence of bedform the opposite is observed. The results
show that grain-scale roughness superimposed on a bedform
does not simply act as a smaller-scale “bedform” with a lesser
effect on pumping. Instead, roughness modulates the bedfrominduced exchange, by changing the macroscopic pressure distribution along SWI. This work demonstrates that grain-scaleroughness on top of a bedform should not be ignored whenstudying hyporheic exchange, as shown by its influence on thepressure variation and transit times.

REFERENCES

Aubeneau, A. F., Martin, R. L., Bolster, D., Schumer, R., Jerolmack, D. & Packman, A. I. 2015 Fractal patterns in riverbed

6