3.2.1 DA velocity and friction

Figure 8.

(a) Diagnostic function for fitting the logarithmic profile and (b) DA velocity profiles for regular (– – –) and random (——) interfaces; - - - (red) fitted logarithmic profiles.

The location of the zero-plane displacement, $d$, and the von Kármán constant, $\unicode[STIX]{x1D705}$, of the wall-bounded flow are determined by fitting the streamwise DA velocity profile to the logarithmic law,

(3.1)$$\begin{eqnarray}U^{+}=\frac{1}{\unicode[STIX]{x1D705}}\log (y+d)^{+}+B,\end{eqnarray}$$

where $B$ is the intercept of the logarithmic profile. The fitting procedure is similar to the method used by Breugem et al. (Reference Breugem, Boersma and Uittenbogaard2006) and Suga et al. (Reference Suga, Matsumura, Ashitaka, Tominaga and Kaneda2010) and is demonstrated in figure 8(a). The fitted velocity profiles are shown in figure 8(b). The fitted values of $\unicode[STIX]{x1D705}$ and $d$ are listed in table 2; these values are close to the results ($\unicode[STIX]{x1D705}=0.31{-}0.32$ and $d/\unicode[STIX]{x1D6FF}=0.06{-}0.1$) reported by Breugem et al. (Reference Breugem, Boersma and Uittenbogaard2006), Suga et al. (Reference Suga, Matsumura, Ashitaka, Tominaga and Kaneda2010) and Manes et al. (Reference Manes, Poggi and Ridolfi2011) at similar $Re_{K}$. It is well established that $\unicode[STIX]{x1D705}$ for turbulence bounded by a porous wall is smaller than the corresponding value for an impermeable wall, which is usually within $5\,\%$ of 0.41 for both smooth and rough walls. Despite sharing similar $\unicode[STIX]{x1D705}$ values, Cases 2 and 3 differ significantly in $d$. Recalling that $d$ is measured from the sediment crest height, the higher $d$ in Case 3 is probably due to the larger variance in local sediment heights associated with the random interface. In the following plots we offset $y$ by the amount of $d$ to effectively collapse the logarithmic regions between the two cases.

Figure 9.

DA velocity along $y$ in linear scaling for regular (– – –) and random (——) interfaces. Inset magnifies the distribution at $y\approx -d$.

The DA velocity profiles in linear scaling (figure 9) show the differences in the velocities near the interface. The velocity in the vicinity of the regular interface is negative with a small magnitude. This can be related to the organized recirculation regions induced by each sphere at the top layer (discussed in § 3.2.3). This is an important observation with implications for interfacial exchange of solutes. Also, simpler models, especially those that ignore inertial effects, may not be able to capture these features. Such a flow pattern is also described in the conceptual model of Pokrajac & Manes (Reference Pokrajac and Manes2009). For lower $y$, the velocity decreases and eventually reaches a constant, positive value inside the sediment. For the random interface, the velocity variation is monotonic. Both velocity profiles exhibit an inflection point near the sediment crest, consistent with observations of Manes et al. (Reference Manes, Poggi and Ridolfi2011), Voermans et al. (Reference Voermans, Ghisalberti and Ivey2017).

The depth of shear-layer penetration, $\unicode[STIX]{x1D6FF}_{b}$, is defined as the distance from the $y=-d$ to the $y$ location separating the constant-velocity region in the bed and the shear layer above; $\unicode[STIX]{x1D6FF}_{b}$ is also referred to as the Brinkman-layer thickness. Here, $\unicode[STIX]{x1D6FF}_{b}$ is calculated in the same way as in Voermans et al. (Reference Voermans, Ghisalberti and Ivey2017), $\langle \overline{u}\rangle _{y+d=-\unicode[STIX]{x1D6FF}_{b}}=0.01(U_{i}-U_{p})+U_{p}$, where $U_{i}$ is the DA velocity at $y=-d$ and $U_{p}$ is the constant velocity deep in the sediment; $\unicode[STIX]{x1D6FF}_{b}$ is the penetration depth measured from $y=-d$. The value measured from the crest, $\unicode[STIX]{x1D6FF}_{b}^{\ast }=\unicode[STIX]{x1D6FF}_{b}+d$, for the two cases is shown in table 2, equalling $0.138\unicode[STIX]{x1D6FF}$ (regular case) and $0.165\unicode[STIX]{x1D6FF}$ (random case). They are of the same order as the grain diameter, as also observed by Goharzadeh, Khalili & Jørgensen (Reference Goharzadeh, Khalili and Jørgensen2005).

The friction coefficient is defined as

(3.2)$$\begin{eqnarray}C_{f}=2\left(\frac{u_{\unicode[STIX]{x1D70F}}}{U_{\unicode[STIX]{x1D6FF}}}\right)^{2}.\end{eqnarray}$$

Table 2 shows that the friction coefficient on the random interface is 70 % higher than the regular one. A higher $C_{f}$ for the random interface is expected as the larger interface length scales ($\unicode[STIX]{x1D70E}_{h}^{+}$ and $\unicode[STIX]{x1D706}^{+}$) lead to higher total drag, similar to the effects of an impermeable wall with a larger roughness length scale.

3.2.2 Turbulent fluctuations

Figure 10(a,b) compares the Reynolds-stress tensor and its anisotropy between the two cases. The Reynolds-stress anisotropy tensor is defined as

(3.3)$$\begin{eqnarray}b_{ij}=\frac{\langle \overline{u_{i}^{\prime }u_{j}^{\prime }}\rangle }{\langle \overline{u_{k}^{\prime }u_{k}^{\prime }}\rangle }-\frac{\unicode[STIX]{x1D6FF}_{ij}}{3}.\end{eqnarray}$$

It is shown that Townsend’s wall-similarity hypothesis for a wall-bounded flow applies to the outer layer ($y/\unicode[STIX]{x1D6FF}>0.2$). As $y$ decreases from the sediment crest, for both cases all components of the Reynolds-stress tensor are damped rapidly, together with a decrease of tensor anisotropy. In the interface region, the random interface results in a lower Reynolds-stress anisotropy with a higher fraction of TKE residing in $v^{\prime }$ and $w^{\prime }$ motions.

Figure 10.

Components of (a) the Reynolds-stress tensor and (b) the anisotropy tensor for regular (– – –) and random (——) interfaces.

Figure 11.

Contours of instantaneous ${u^{\prime }}^{+}$ in the $(x,z)$ plane at the respective $\langle \overline{u^{\prime 2}}\rangle$-peak elevations for the (a) regular and (b) random interfaces.

The more isotropic turbulence near the random interface is linked to the augmented disturbance from the sediment obstruction to the near-wall turbulent coherent structures. Figure 11 compares the contours of ${u^{\prime }}^{+}$ at the respective locations of the $\langle \overline{u^{\prime 2}}\rangle$ peaks in the two cases. The low-speed streaks near the random interface are significantly shorter in wall units than those near the regular interface due to the physical blockage from the spheres serving as roughness elements.

The penetration depth of turbulent shear stress into the bed, $\unicode[STIX]{x1D6FF}_{p}$, is obtained from $\langle \overline{u^{\prime }v^{\prime }}\rangle _{y=-\unicode[STIX]{x1D6FF}_{p}}=0.01\langle \overline{u^{\prime }v^{\prime }}\rangle _{i}$, following Voermans et al. (Reference Voermans, Ghisalberti and Ivey2017), and shown in table 2; $\langle \overline{u^{\prime }v^{\prime }}\rangle _{i}$ is the value of $\langle \overline{u^{\prime }v^{\prime }}\rangle$ at $y=-d$; $\unicode[STIX]{x1D6FF}_{p}^{\ast }$ is the same depth measured from the crest. The flow over the random interface gives significantly deeper turbulence penetration, as shown by a 47 % higher $\unicode[STIX]{x1D6FF}_{p}^{\ast }$ and 20 % higher $\unicode[STIX]{x1D6FF}_{p}$ compared to the regular interface.

Figure 12.

Integral length scales of $u^{\prime }$ motions along the streamwise direction for regular (– – –) and random (——) interfaces.

To quantitatively compare turbulent scales, the integral length scales are shown in figure 12. It is defined as

(3.4)$$\begin{eqnarray}L_{ij,x_{k}}(y)=\int _{0}^{\infty }R_{ij}(y,\unicode[STIX]{x0394}r_{k})d(\unicode[STIX]{x0394}r_{k}),\end{eqnarray}$$

where

(3.5)$$\begin{eqnarray}R_{ij}(y,\unicode[STIX]{x0394}r_{k})=\frac{\langle \overline{u_{i}^{\prime }(y,x_{k})u_{j}^{\prime }(y,x_{k}+\unicode[STIX]{x0394}r_{k})}\rangle }{\unicode[STIX]{x1D70E}_{u_{i}}(y)\unicode[STIX]{x1D70E}_{u_{j}}(y)}\end{eqnarray}$$

(no summation over repeated indices). Here, $\unicode[STIX]{x0394}r_{k}$ is the separation along the $x_{k}$ direction and $\unicode[STIX]{x1D70E}_{u_{i}}(y)$, $\unicode[STIX]{x1D70E}_{u_{j}}(y)$ are the root-mean-square of turbulent fluctuations. The integration is carried out to the value of $\unicode[STIX]{x0394}r_{k}$ at which the correlation coefficient first crosses 0.3. A threshold value between 0.3 and 0.5 is usually used to calculate the integral lengths from two-point correlation coefficients (see, for example, Krogstad & Antonia (Reference Krogstad and Antonia1994), Christensen & Wu (Reference Christensen and Wu2005) and Volino, Schultz & Flack (Reference Volino, Schultz and Flack2011)). It has been checked that varying this threshold from 0.2 to 0.5 would not noticeably affect the comparison. It is shown that, for the random interface, at $y=-d$ the coherent motions are more extensive in $x$ due to a deeper flow penetration. At the sediment crest, however, the coherent motions are noticeably shorter in $x$ for the random interface, due to the disturbances of the larger roughness height.

Farther away from the interface ($(y+d)/\unicode[STIX]{x1D6FF}>0.15$), the difference in the $x$ integral length between the two cases systematically increases with $y$, indicating a lack of Townsend’s wall-similarity hypothesis for this flow quantity. One reason is the relatively high roughness. Jiménez (Reference Jiménez2004) has proposed that Townsend’s hypothesis applies when the roughness height is less than $0.02\unicode[STIX]{x1D6FF}$. In the present study, the roughness height (measured from $y=-d$) is at least $0.06\unicode[STIX]{x1D6FF}$. Another possible reason is the limited Reynolds number used in the present DNS. The experimental study of Krogstad & Antonia (Reference Krogstad and Antonia1994) showed a lack of similarity in the outer layer $L_{11}$ in $x$ between a smooth-wall and a rough-wall turbulent boundary-layer flow. However, a repeated experiment of boundary layer flows by Krogstad & Efros (Reference Krogstad and Efros2012) at a higher Reynolds number indicated reduced discrepancy in the outer layer $L_{11}$.

3.2.3 Form-induced fluctuations

Figure 13.

Components of the form-induced stress tensor: (a) streamwise, (b) wall-normal and (c) spanwise components, for regular (– – –) and random (——) interfaces.

Figure 14.

(a) Time-averaged streamlines for the regular interface. (b) Streamlines on Slice A ($(x,y)$ plane across the crests of the grains). (c) Streamline distribution on Slice B ($(x,y)$ plane across the valley between two rows of grains). Dashed lines indicate $y=-d$.

Figure 15.

(a) Time-averaged streamlines for the random interface. Panels (b) and (c) show streamlines on two $(x,y)$ planes, Slice C and Slice D. Dashed lines indicate $y=-d$.

The form-induced stresses are shown in figure 13. The peak values of all components are much smaller than the corresponding Reynolds stresses, consistent with observations of Manes et al. (Reference Manes, Pokrajac, McEwan and Nikora2009b), Voermans et al. (Reference Voermans, Ghisalberti and Ivey2017) and Fang et al. (Reference Fang, Han, He and Dey2018). The streamwise form-induced stress near the regular interface region is concentrated in a narrower layer near the crest elevation with a much higher peak value. This is due to the vertical alignment of the wake regions of the grains, as shown in figure 14(b,c). Large values of $\langle \widetilde{u}^{2}\rangle$ are found in the troughs of the top sediment surface and in the wake regions of the surface protuberances. Since these regions are located in a narrow layer near the crest, the distribution of $\langle \widetilde{u}^{2}\rangle$ is narrower with a larger peak. In contrast, figure 15 shows that, for the random interface, the mean recirculation regions can reach much larger sizes and that the penetration of the time-mean flow is much deeper. In addition, $\langle \widetilde{v}^{2}\rangle$, $\langle \widetilde{w}^{2}\rangle$ and $\langle \widetilde{u}\,\widetilde{v}\rangle$ are higher near the random interface. This may be explained by the higher spatial variation of $\widetilde{u}$ associated with the disordered blockage, which leads to more intense $\widetilde{v}$ and $\widetilde{w}$ through continuity.

Figure 16.

Shear stresses for cases with (a) regular and (b) random interfaces: - - - (blue) Reynolds shear stress, – – – viscous shear stress and —— (red) form-induced shear stress.

To analyse the importance of form-induced stresses to the mean momentum balance, we compare the magnitudes of the various shear stresses in figure 16. For the regular interface, the form-induced shear stress is negligible; the viscous shear stress contributes significantly to the total shear stress at the sediment crest. However, for the random interface the viscous shear is weaker due to the milder DA velocity gradient associated with the thicker shear layer; a much stronger form-induced stress is present.

It is important to note that, in the vicinity of the random interface, the form-induced shear stress reaches a similar magnitude of the Reynolds shear stress (figure 16b). The fact that the form-induced stress potentially affects the DA velocity near the interface and that it strongly depends on the interface geometry calls for a systematic study of the generation of the form-induced stresses; however, this is beyond the scope of the present work.

3.2.4 Reynolds-stress budgets

The budget equations of the normal Reynolds stresses, $\langle \overline{u_{\unicode[STIX]{x1D6FC}}^{\prime 2}}\rangle _{s}$, can be written as (Raupach & Shaw Reference Raupach and Shaw1982; Mignot et al. Reference Mignot, Bartheleemy and Hurther2009; Yuan & Aghaei Jouybari Reference Yuan and Aghaei Jouybari2018)

(3.6)$$\begin{eqnarray}\displaystyle 0 & = & \displaystyle \underbrace{-2\langle \overline{u_{\unicode[STIX]{x1D6FC}}^{\prime }v^{\prime }}\rangle _{s}\frac{\unicode[STIX]{x2202}\langle \overline{u}_{\unicode[STIX]{x1D6FC}}\rangle }{\unicode[STIX]{x2202}y}}_{P_{s}}\underbrace{-2\langle \overline{u_{\unicode[STIX]{x1D6FC}}^{\prime }u_{j}^{\prime }}\rangle _{s}\left\langle \frac{\unicode[STIX]{x2202}\widetilde{u}_{\unicode[STIX]{x1D6FC}}}{\unicode[STIX]{x2202}x_{j}}\right\rangle }_{P_{m}}\underbrace{-2\left\langle \widetilde{u_{\unicode[STIX]{x1D6FC}}^{\prime }u_{j}^{\prime }}\frac{\unicode[STIX]{x2202}\widetilde{u}_{\unicode[STIX]{x1D6FC}}}{\unicode[STIX]{x2202}x_{j}}\right\rangle _{s}}_{P_{w}}\underbrace{-\left\langle \frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x_{j}}\widetilde{u_{\unicode[STIX]{x1D6FC}}^{\prime }u_{\unicode[STIX]{x1D6FC}}^{\prime }}\widetilde{u}_{j}\right\rangle _{s}}_{T_{w}}\nonumber\\ \displaystyle & & \displaystyle \underbrace{-\left\langle \frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x_{j}}\overline{u_{\unicode[STIX]{x1D6FC}}^{\prime }u_{\unicode[STIX]{x1D6FC}}^{\prime }u_{j}^{\prime }}\right\rangle _{s}}_{T_{t}}\underbrace{+\unicode[STIX]{x1D708}\left\langle \frac{\unicode[STIX]{x2202}^{2}\overline{u_{\unicode[STIX]{x1D6FC}}^{\prime 2}}}{\unicode[STIX]{x2202}x_{j}\unicode[STIX]{x2202}x_{j}}\right\rangle _{s}}_{T_{\unicode[STIX]{x1D708}}}\underbrace{-2\left\langle \overline{u_{\unicode[STIX]{x1D6FC}}^{\prime }\frac{\unicode[STIX]{x2202}P^{\prime }}{\unicode[STIX]{x2202}x_{\unicode[STIX]{x1D6FC}}}}\right\rangle _{s}}_{\unicode[STIX]{x1D6F1}}\underbrace{-2\unicode[STIX]{x1D708}\left\langle \overline{\frac{\unicode[STIX]{x2202}u_{\unicode[STIX]{x1D6FC}}^{\prime }}{\unicode[STIX]{x2202}x_{j}}\frac{\unicode[STIX]{x2202}u_{\unicode[STIX]{x1D6FC}}^{\prime }}{\unicode[STIX]{x2202}x_{j}}}\right\rangle _{s}}_{\unicode[STIX]{x1D716}}.\end{eqnarray}$$

On the right-hand side, the first three terms are, respectively, shear production $(P_{s})$ and additional productions due to the form-induced shear ($P_{w}$ and $P_{m}$). The following three terms are the transport terms due to the form-induced (wake) fluctuations $(T_{w})$, turbulent fluctuations $(T_{t})$ and viscous diffusion $(T_{\unicode[STIX]{x1D708}})$; $\unicode[STIX]{x1D6F1}$ is the pressure work and $\unicode[STIX]{x1D716}$ is the viscous dissipation. Summing over the equations of the three Reynolds-stress components yields the budget equation of the TKE; the sum of the $\unicode[STIX]{x1D6F1}$ terms yields the pressure transport, $T_{p}$. The production terms and viscous dissipation of the TKE are shown in figure 17(a), whereas the transport terms are shown in figure 17(b). The residuals are around $3\,\%$ of $P_{s}$ peak values.

For $(y+d)/\unicode[STIX]{x1D6FF}\geqslant 0.2$, the shear production is balanced by the viscous dissipation (not shown). The $P_{s}$ terms reach their maximum values close to the sediment crest. The total amounts of form-induced production ($P_{m}+P_{w}$) are significant near the interface, becoming dominant a short distance below the crest. Fang et al. (Reference Fang, Han, He and Dey2018) also observed non-negligible form-induced TKE production for regular interface particle distribution with $Re_{K}=0{-}27$.

In the vicinity of the regular interface, $T_{\unicode[STIX]{x1D708}}$, $T_{t}$ and $T_{p}$ remove TKE from the high-TKE region slightly above the crest and transfer it into the low-TKE region in the bed; $T_{\unicode[STIX]{x1D708}}$ dominates other TKE transports due to high $\unicode[STIX]{x2202}U/\unicode[STIX]{x2202}y$, while $T_{w}$ is negligible. In contrast, for the random interface $T_{\unicode[STIX]{x1D708}}$ is negligible and the augmented $\widetilde{v}$ fluctuations (figure 13) lead to a non-negligible form-induced transport, which works against the other transport processes by moving TKE upward from low-TKE region in the bed to the crest region corresponding to the TKE peak. Such process contributes to a more equilibrium turbulence (with more energy dissipated at the location where it is generated).

Figure 17.

Terms in the TKE budget equation. (a) Production and viscous dissipation for regular (symbols and thin lines) and random interface (thick lines) cases; $P_{s}$ (——, ▵), $P_{m}+P_{w}$ (—— (red), ○ (red)) and $\unicode[STIX]{x1D716}$ (– – –, ▫). Transport terms for (b) regular and (c) random interfaces; $T_{w}$ (—— (red)), $T_{t}$ (- - -), $T_{\unicode[STIX]{x1D708}}$ (——) and $T_{p}$ (– – –). All terms are normalized by $u_{\unicode[STIX]{x1D70F}}^{3}/\unicode[STIX]{x1D6FF}$.

Figure 18.

(a) Form-induced production and (b) pressure–strain-rate term in regular (– – –) and random interface (——) cases for the budgets of $\langle \overline{u^{\prime 2}}\rangle _{s}$ (○), $\langle \overline{v^{\prime 2}}\rangle _{s}$ (▫) and $\langle \overline{w^{\prime 2}}\rangle _{s}$ (▵); normalization is done with $u_{\unicode[STIX]{x1D70F}}$ and $\unicode[STIX]{x1D6FF}$.

Figure 19.

Time-averaged spanwise vorticity $\overline{\unicode[STIX]{x1D714}}_{z}$ normalized by $u_{\unicode[STIX]{x1D70F}}/\unicode[STIX]{x1D6FF}$ at (a) Slice B of the regular interface and (b) Slice D of the random interface.

The total form-induced production ($P_{w}+P_{m}$) for individual normal Reynolds stresses are compared in figure 18(a). For both cases below the crest, it predominantly contributes to $\langle \overline{u^{\prime 2}}\rangle$ other than the wall-normal and spanwise components, more so for the regular interface. $P_{m}$ can be rewritten as $-2\langle \overline{u_{\unicode[STIX]{x1D6FC}}^{\prime }v^{\prime }}\rangle \langle \overline{u}_{\unicode[STIX]{x1D6FC}}\rangle \text{d}\unicode[STIX]{x1D703}/\text{d}y$. Thus, $P_{m}$ remains positive for $\langle \overline{u^{\prime 2}}\rangle$ and is zero for $\langle \overline{v^{\prime 2}}\rangle$ and $\langle \overline{w^{\prime 2}}\rangle$. The negative $P_{w}+P_{m}$ in both cases in the vicinity of the crest is then due to negative $P_{w}$, indicating conversion of TKE to the kinetic energy of form-induced fluctuations, since this term also exists in the form-induced stress budgets with a positive sign (Raupach & Shaw Reference Raupach and Shaw1982). Such conversion is probably due to the work of particle drag against the turbulent motions of scales larger than the particles, generating form-induced fluctuations at the scale of the particles. Deeper into the porous bed, the scales of turbulent motions are rapidly reduced and are smaller than that of the particles (as shown by the integral scales in figure 12); here, the particle-scale form-induced shear layers generates turbulence, represented by positive $P_{w}$.

The pressure work for a normal Reynolds stress is decomposed into the pressure transport and a pressure–strain-rate term, ${\mathcal{R}}=2\langle \overline{P^{\prime }\unicode[STIX]{x2202}u_{\unicode[STIX]{x1D6FC}}^{\prime }/\unicode[STIX]{x2202}x_{\unicode[STIX]{x1D6FC}}}\rangle _{s}$. The latter is compared in figure 18(b); ${\mathcal{R}}$ distributes energy among different normal components of the Reynolds-stress tensor and forms the primary source for the $v^{\prime }$- and $w^{\prime }$-energy budgets that balances the viscous dissipation. The main difference is that, at the crest of the regular interface, the $u^{\prime }$ energy converts mostly to that of $w^{\prime }$, with ${\mathcal{R}}_{22}\approx 0$, while ${\mathcal{R}}_{22}$ and ${\mathcal{R}}_{33}$ are of comparable magnitudes at the crest for the random case. Liu & Katz (Reference Liu and Katz2018) observed experimentally that, in a shear layer formed over an open cavity, the magnitude of ${\mathcal{R}}_{11}$ is an order of magnitude higher than ${\mathcal{R}}_{22}$. The similar observation herein may also be attributed to the local mean shear layers formed above the topmost-layer particles, shown in figure 19 using time-mean spanwise vorticity in an $(x,y)$ plane. Such shear layers are stronger and aligned in $y$ above the regular interface, while the random distribution of particle for the random interface results in a spread in such effect on energy redistribution.

3.3.1 Variation of apparent permeability

As a convention of the oil and gas industry (Edwards et al. Reference Edwards, Shapiro, Bar-Yoseph and Shapira1990; Huang & Ayoub Reference Huang and Ayoub2006; Javadpour Reference Javadpour2009), the apparent permeability, $K^{\ast }$, is used to establish a relationship between mean velocity and macroscopic pressure gradient in a finite Reynolds number flow inside a porous medium. Here, the local value of $K^{\ast }$ is obtained based on the local volume-averaged velocity and the gradient of the local volume-averaged pressure,

(3.7)$$\begin{eqnarray}\left\langle \overline{u}\right\rangle _{s}(x,y,z)=-\frac{K^{\ast }(x,y,z)}{\unicode[STIX]{x1D707}}\frac{\unicode[STIX]{x2202}\langle \overline{p}\rangle (x,y,z)}{\unicode[STIX]{x2202}x},\end{eqnarray}$$

where $\langle \cdot \rangle _{s}$ and $\langle \cdot \rangle$ indicate, respectively, the superficial and intrinsic volume averaging performed locally (within a volume much larger than a grain size). Normalized by $u_{\unicode[STIX]{x1D70F}}$ and $\unicode[STIX]{x1D6FF}$, equation (3.7) becomes

(3.8)$$\begin{eqnarray}\langle \overline{u}\rangle _{s}^{+}=-Re_{\unicode[STIX]{x1D70F}}\left(\frac{K^{\ast }}{\unicode[STIX]{x1D6FF}^{2}}\right)\left(\frac{\unicode[STIX]{x2202}\langle \overline{P}\rangle ^{+}}{\unicode[STIX]{x2202}x/\unicode[STIX]{x1D6FF}}\right),\end{eqnarray}$$

keeping in mind that $P=p/\unicode[STIX]{x1D70C}$. For convenience, from here on we write the dimensionless $K^{\ast }/\unicode[STIX]{x1D6FF}^{2}$ as $K^{\ast }$. Note that $K^{\ast }$ is not the Darcy permeability, but a macroscopic descriptor of the complex flow at the interface, including the influences of both diffusion of the mean momentum (traditionally modelled using the Brinkman correction (Brinkman Reference Brinkman1949)) and microscopic inertial and drag forces (traditionally modelled using Forchheimer’s correction (Hassanizadeh & Gray Reference Hassanizadeh and Gray1987)).

Figure 20.

Horizontal variation of the apparent permeability $K^{\ast }$ in the Brinkman layer in the random interface case.

Since the Brinkman layer (measured from crest) has a thickness ($\unicode[STIX]{x1D6FF}_{b}^{\ast }$) of the order of a sphere diameter, $K^{\ast }$ in this layer is considered constant along $y$, but varying in $x$ and $z$. To calculate such a $K^{\ast }(x,z)$, we divide the Brinkman layer horizontally into small segments with a size of $2D\times 2D$ in $x$ and $z$. Enlarging the segment size to $4D\times 4D$ or $8D\times 8D$ does not fundamentally change the observation. The velocity and pressure values in (3.7) are spatially averaged in each small segment; the pressure gradient for each segment is obtained using a second-order central difference. The $(x,z)$ distribution of $K^{\ast }$ for the random interface is shown in figure 20; the variation for the regular interface is between $1\times 10^{-3}$ and $2\times 10^{-3}$, barely visible using the same colour bar and is thus not shown. The random interface leads to a much more heterogeneous permeability distribution with a higher mean value. Based on its natural logarithm, the mean, $\unicode[STIX]{x1D707}_{ln(K^{\ast })}$, and the variance, $\unicode[STIX]{x1D70E}_{ln(K^{\ast })}^{2}$, are calculated and summarized in table 2. For the regular and random interfaces, respectively, the $\unicode[STIX]{x1D707}_{ln(K^{\ast })}$ values are $-6.5$ and $-5.4$ (or twice higher in mean apparent permeability for the random case) and the $\unicode[STIX]{x1D70E}_{ln(K^{\ast })}^{2}$ values are $0.005$ and $0.055$. Since the square root of permeability has the physical meaning of the effective pore size (Breugem et al. Reference Breugem, Boersma and Uittenbogaard2006), the results here show that the random interface increases the effective pore size at the interface associated with a larger apparent permeability. This is consistent with the larger streamwise integral length scales of $u^{\prime }$ motions at $y\approx -d$, as shown in figure 12.

The $K^{\ast }$ value is also calculated below the Brinkman layer using (3.8). Here, $K^{\ast }$ becomes the bulk permeability, $K$, associated with laminar flow inside the porous medium, which is typically predicted using a semi-empirical model such as the Kozeny–Carman equation. The value of $\text{ln}(K^{\ast })$ obtained using the DNS data is almost a constant of $-10.9$, regardless of location. Such value indicates a $K^{\ast }$ that is almost two orders of magnitude smaller than the value in the Brinkman layer. In comparison, the Kozeny–Carman equation predicts a $\text{ln}(K)$ value of $-10.0$ inside the bulk of the sediment, overestimating the bulk permeability by almost 1.5 times.

In the next section, it will be shown that the spatial variation of apparent permeability in the Brinkman layer of the random interface is correlated with the roughness geometry. As $K^{\ast }$ in the Brinkman layer incorporates the effect of both macroscopic diffusion and microscopic inertia on the macroscopic (volume-averaged) velocity, we focus on the spatial variation of these two mechanisms as affected by the roughness geometry.

3.3.2 Momentum transport mechanisms and dependence on roughness geometry

Figure 21.

Terms in the DA $u$-momentum equation (3.9), normalized by $u_{\unicode[STIX]{x1D70F}}^{2}/\unicode[STIX]{x1D6FF}$ for (a) regular and (b) random interfaces: Reynolds-stress term(—— (red)), viscous stress term (– – – (red)), total drag (——), form-induced stress term (– – –) and pressure-gradient term (- - -).

The double-averaged $u$-momentum equation can be written as (Raupach & Shaw Reference Raupach and Shaw1982)

(3.9)$$\begin{eqnarray}-\frac{\unicode[STIX]{x2202}\langle \bar{P}\rangle _{s}}{\unicode[STIX]{x2202}x}+\underbrace{\left(\unicode[STIX]{x1D708}\frac{\unicode[STIX]{x2202}^{2}\langle \bar{u}\rangle _{s}}{\unicode[STIX]{x2202}y^{2}}-\frac{\unicode[STIX]{x2202}\langle \overline{u^{\prime }v^{\prime }}\rangle _{s}}{\unicode[STIX]{x2202}y}\right)}_{{\mathcal{D}}}-\underbrace{\left(\frac{\unicode[STIX]{x2202}\langle \tilde{u} \tilde{v}\rangle _{s}}{\unicode[STIX]{x2202}y}+F_{D}\right)}_{{\mathcal{I}}}=0,\end{eqnarray}$$

where $F_{D}$ is the total drag per unit mass (sum of both pressure drag and viscous drag around the grains). According to Whitaker (Reference Whitaker1996), the viscous and Reynolds shear stress terms represent the effects due to diffusion of the mean momentum, usually modelled using the Brinkman correction; the sum of these two terms is denoted herein as ${\mathcal{D}}$, representing ‘diffusion’. The terms related to the form-induced fluctuations (i.e. the form-induced shear stress and the total drag) represent the effects of fluid inertia at the grain scale on the mean momentum, usually modelled by the Forchheimer correction. Their sum is denoted as ${\mathcal{I}}$, representing (microscopic) ‘inertia’.

Figure 21 shows the balance of equation (3.9) across the interface, with the spatial averaging performed using area averaging in the $(x,z)$ plane of the whole domain. The residual is less than $2.3\,\%$ of the peak magnitude of the total force term for the two cases. For both cases, the spatial-averaged values of ${\mathcal{D}}$ and ${\mathcal{I}}$ are both significant inside the Brinkman layer, approximately balancing each other as the pressure gradient is comparatively weak. The total drag is the main contributor to ${\mathcal{I}}$ for both cases. For the regular interface, ${\mathcal{D}}$ is contributed to almost equally by both Reynolds and viscous stresses, while it is predominantly due to the Reynolds stress for the random one; this is consistent with the observations made for shear stress distribution (figure 19). The total drag peaks at the crest for the regular interface. This is because the main contributor of $F_{D}$ – the pressure drag – is caused by recirculation regions, which are distributed uniformly near the crest for this case. For the random interface, recirculation regions of various sizes are formed due to various scales of roughness protuberances, leading to a peak of $F_{D}$ below the crest.

Figure 22.

Horizontal variation of local volume-averaged ${\mathcal{D}}$ (a) and $|{\mathcal{I}}|$ (b), both normalized by $u_{\unicode[STIX]{x1D70F}}^{2}/\unicode[STIX]{x1D6FF}$. (c) Sediment surface height fluctuations $h^{\prime }$ normalized by $\unicode[STIX]{x1D6FF}$ for the random interface case. (d) Correlations of ${\mathcal{D}}^{\prime }$ (red) and $|{\mathcal{I}}|^{\prime }$ (black) with $h^{\prime }$.

We now explore potential correlation between the geometry of the random interface and the $(x,z)$ variation of locally averaged ${\mathcal{D}}$ and ${\mathcal{I}}$ magnitudes inside the Brinkman layer. The same local volume-averaging method described in § 3.3.1 is used. The results for the regular interface are not shown, as the variations are barely visible using the same range of contour levels. Figures 22(a) and 22(b) compare the $(x,z)$ distribution of locally averaged ${\mathcal{D}}$ and $|{\mathcal{I}}|$, respectively. The spatial average is similar for both contours, a result of a relatively weak pressure gradient discussed earlier. The spatial variation is much more intense for $|{\mathcal{I}}|$ than for ${\mathcal{D}}$. This indicates that, although roughness leads to local augmentation of $\overline{u^{\prime }v^{\prime }}$ (predominantly through wake-induced production), such spatial variation is rather weak and does not significantly affect the $(x,z)$ pattern of wall-normal $\langle \overline{u}\rangle _{s}$ transport. In contrast, the pressure drag is generated by the small-scale wake of each grain. Therefore, $|{\mathcal{I}}|$ is more sensitive to the local roughness height. Visual comparison with the roughness height distribution $h^{\prime }(x,z)$ in figure 22(c) shows that the distributions of $|{\mathcal{I}}|$ and $h^{\prime }$ are somewhat correlated. Indeed, the scatter plot (figure 22d) of $h^{\prime }$ and the spatial fluctuation of $|{\mathcal{I}}|$ (departure from the $x$–$z$ mean) displays a positive correlation (with preferred distribution in quadrants I and III). This is consistent with the expectation that a high roughness protuberance tends to generate a large (in $y$) recirculation region and, consequently, a higher-than-average value of pressure drag at the corresponding $(x,z)$ location.

As $K^{\ast }$ describes the global effect of both ${\mathcal{I}}$ and ${\mathcal{D}}$, the $K^{\ast }$ contour in figure 20 also displays correlation with $h^{\prime }(x,z)$ distribution; the correlation is negative as a higher local drag serves as a stronger momentum sink.

3.3.3 Mechanisms producing vertical mass flux

The hyporheic flux $(Q_{H})$ quantifies the magnitude of mass exchange across the interface during a period of time; it is defined as the wall-normal volumetric flow rate at $y=-d$ directed out of the permeable bed (Elliott & Brooks Reference Elliott and Brooks1997),

(3.10)$$\begin{eqnarray}Q_{H}(T)=\frac{1}{A}\int _{(x,z)}\overline{\left.v(t,\boldsymbol{x})\right|_{v>0}}^{T}\,\text{d}A,\end{eqnarray}$$

where $\overline{()}^{T}$ denotes running time averaging within a window size $T$ and $A$ is the fluid area of the $(x,z)$ plane at $y=-d$. Figure 23(a) shows the variation of the vertical flux, $Q_{H}$, as a function of $T$. The error bars quantify the variance of the values obtained from the DNS data for each $T$; the variance is zero for $T$ equal to the total simulation time and increases with decreasing $T$. For $T$ much larger than the large-eddy turn-over time, $\unicode[STIX]{x1D6FF}/u_{\unicode[STIX]{x1D70F}}$, the form-induced wall-normal velocity ($\widetilde{v}$) is the only contributor of $Q_{H}$, while for $T<\unicode[STIX]{x1D6FF}/u_{\unicode[STIX]{x1D70F}}$ the turbulence fluctuation ($v^{\prime }$) also contributes to the exchange. Evidently, the contribution of $\widetilde{v}$ is more significant than that of $v^{\prime }$ for both interfaces. Also, $Q_{H}$ is around five times larger for the random interface compared to the other, consistent with the more intense Reynolds and form-induced stresses near $y=-d$.

Figure 23.

(a) Wall-normal flux at $y=-d$ normalized by $u_{\unicode[STIX]{x1D70F}}\unicode[STIX]{x1D6FF}^{2}$ with different averaging window sizes $Tu_{\unicode[STIX]{x1D70F}}/\unicode[STIX]{x1D6FF}$ for regular (– – –) and random (——) interfaces. (b) Pressure work term of $\langle \widetilde{v}^{2}\rangle _{s}$ budget normalized by $u_{\unicode[STIX]{x1D70F}}^{3}/\unicode[STIX]{x1D6FF}$ for regular (– – –, ○) and random (——) interfaces.

The mechanisms governing $v^{\prime }$ energy transport have been discussed in § 3.2.4 using the stress budget. For the transport of $\widetilde{v}$ energy, studies on impermeable flows over rough surfaces (Raupach & Shaw Reference Raupach and Shaw1982; Yuan & Piomelli Reference Yuan and Piomelli2014b) showed that the process is governed by the work done by the form-induced pressure, $\widetilde{\unicode[STIX]{x1D6F1}}=-2\langle \widetilde{v}\unicode[STIX]{x2202}\widetilde{P}/\unicode[STIX]{x2202}y\rangle _{s}$, as the main source, the wake production predominantly as a sink. The wall-normal variation of $\widetilde{\unicode[STIX]{x1D6F1}}$ is compared in figure 23(b); the much higher overall value for the random interface is due to the more intense $\widetilde{P}$ fluctuations (shown in figure 24 in an $(x,y)$ slice for each case) leading to significant $\unicode[STIX]{x2202}\widetilde{P}/\unicode[STIX]{x2202}y$ near the interface.

Thus, the main difference between the $Q_{H}$ generation in the two cases is that the random particle arrangement leads to larger-scale and more intense time-mean pressure variations which, by working against the time-mean velocity, generate more intense $\widetilde{v}$ near the interface and promote vertical momentum exchange. Meanwhile, $v^{\prime }$ is also intensified by the pressure work near the interface, also contributing to $Q_{H}$, although with a smaller magnitude.

Figure 24.

Time-averaged pressure variation and three-dimensional mean streamlines for (a) regular and (b) random interfaces. Blue arrows indicate examples of hyporheic flow paths.

The hyporheic flow paths are also of interests due to their impact on important macroscopic exchange parameters such as the residence-time distribution and the penetration depth. Below the region of turbulence penetration ($1D$ below the crest), the flow is steady and, consequently, the mean streamlines shown in figure 24 can be interpreted as the hyporheic flow paths of particles released from the upstream.

For the regular interface, the spatial variation of $\widetilde{P}$ shows distributions with alternating negative and positive values with a length scale of one sphere diameter, representing organized mean recirculation regions that are shown in figure 14(b,c). The main mechanism of vertical time-mean exchange is due to these organized recirculation regions at the interface. The relatively low $K^{\ast }$ at the uppermost layer of the regular interface limits communication between the flows inside the sediment and in the region above; thus, in the sediment, there is little vertical exchange and the lateral flow paths are long.

In contrast, for the random interface, the hyporheic flow paths display a multiscale pattern with short flow paths near the interface and longer flow paths in the shape of circular arcs reaching deeper into the bed. These paths are similar to those visualized through dye transport for a flat bed in the experiments of Packman, Salehin & Zaramella (Reference Packman, Salehin and Zaramella2004). Near the bed, the strong $\widetilde{P}$ gradients of larger coherence length induce flow infiltration near the bed, reminiscent of bedform-induced advective pumping (Packman et al. Reference Packman, Salehin and Zaramella2004). This also results in more intense $\widetilde{v}$, increasing form-induced mixing. The streamwise lengths of the flow paths appear to be reduced as the wall-normal communication is strengthened.