5.1. Sediment ridge evolution

In the following, we provide more rigorous quantitative analysis of the relation between bedforms and the turbulent flow.

As mentioned in § 3, all simulations are initiated with an initially macroscopically flat sediment bed and ridges form solely under the action of the turbulent structures after a few bulk time units. Figure 4 shows the time evolution of the fluctuations of the streamwise-averaged sediment bed ${h_f^{\prime \prime }}(z,t)$. Note that the sediment ridges in the time period of interest are essentially parallel to and without relevant perturbations in the streamwise direction such that they can be considered as essentially statistically invariant in the mean flow direction.

Figure 4. Space–time plot of the streamwise-averaged sediment bed height fluctuations ${h_b^{\prime \prime }}(z,t)/D$. Blue and red regions refer to troughs and crests of the streamwise-averaged fluid–bed interface profiles, respectively. Cases (a) S250, (b) M250, (c) M850, (d) L250.

First, it is concluded that ridges form at different spanwise locations more or less simultaneously and hence independently from each other right after the release of the particles. It is interesting to see that in the early phase, say the first $20$ bulk time units of each simulation, the lateral spacing between adjacent crests and troughs is clearly lower than the conventionally found values of $(1-2){H_f}$ (Wolman & Brush Reference Wolman and Brush1961; Ikeda Reference Ikeda1981; Colombini Reference Colombini1993). Advancing in time, however, a coalescence of the bedforms is observed that causes a reduction of the number of individual ridges. After approximately $40$ bulk time units, merging or splitting events between ridges are less frequent but still observable and the remaining bedforms now feature a lateral spacing approximately in the expected range and a more pronounced amplitude throughout all simulations. Note that an exact match between the observed/estimated spacings and our results is not to be expected, as the time scales at which the former are observed differ by several orders of magnitude and represent quasi-asymptotic states whereas the currently observed state is highly transient. Further recall that the majority of the available experimental datasets have been conducted under the influence of lateral sidewalls in low to moderate aspect ratio flumes and underlie the sidewall-induced secondary currents. The larger and further developed ridges appear to be rather immobile in the lateral direction over the observed time period, i.e. the mean spanwise position of these ridges is quite stable. This and the regular spacing of sediment ridges seem to be no effect of laterally narrow domain sizes, as ridge crests follow the same straight vertical space–time lines in the sufficiently wide channel of case L250.

We further conclude that the narrow boxes with a lateral domain period ${L_z}/{H_f} \approx 3$ are capable of accommodating one or two ridges, whereas the large domain of L250 exhibits nine to ten ridge units. We can therefore consider the small to medium domains as close to minimal in the context of the number of available ridges, which shall be favourable for the subsequent analysis in a sense that individual ridges and their relation to turbulent coherent structures can be investigated excluding possible merging or splitting effects between individual bedforms. The large domain simulation L250, on the other hand, contains a sufficient number of individual ridges to allow statements on the collective behaviour of the bedforms.

Figure 5 provides the time evolution of the mean bedform geometry and the particle transport rate. Figure 5(a,b) shows the development of the mean pattern height expressed by the root mean square of the sediment bed height fluctuations ${\sigma _{h,z}}$ and that of the lateral spacing in terms of the mean bedform wavelength ${\lambda _{h,z}}$ (cf. Appendix B.1 for definitions). In accordance with the discussed space–time plots, it is seen that the bedform amplitude globally increases with time during the initial approximately $40$ bulk time units. In the first approximately $10$ bulk time units, low-amplitude disturbances first form, before they merge in the subsequent phase between $t=10{T_b}$ and $t=20{T_b}$, leading to a net reduction of the number of ridges. The bedforms continue growing in amplitude predominantly during the time interval between $t=20{T_b}$ and $t=40{T_b}$. The lower growth rate in the first $5$ bulk time units of the three low Reynolds number cases S250, M250 and L250 is a consequence of the lower Shields numbers when compared with the high Reynolds number case M850, in which the high Shields number causes a stronger erosion and thus a faster pattern growth. For the sake of comparison, figure 5(a) additionally provides the evolution of the spanwise-averaged root mean square of the sediment bed height fluctuations ${\sigma _{h,x}}$, that is a measure for the amplitude of possibly evolving transverse ripple-like bed features. It is thereby verified that transverse bed features are indeed not relevant in this initial time frame, consistent with the instantaneous bed snapshots in figure 3. The variation of the number of individual patterns is revealed by the oscillations of the mean wavelength ${\lambda _{h,z}}$ during the first $40$ bulk time units. For $t>40{T_b}$, then, ${\lambda _{h,z}}$ settles without further strong oscillations in all cases, attaining final values of $1.47{H_f}$ (S250), $1.14{H_f}$ (M850) and $1.44{H_f}$ (L250) that are of comparable size to the dominant wavelength ${\lambda _{h,z}}=1.3{H_f}$ obtained by means of linear stability analysis (Colombini Reference Colombini1993). The clearly lower attained final value of $0.76{H_f}$ in case M250 is in agreement with the corresponding space–time plot in figure 4 which shows that the bed in the final phase indeed features three ridge patterns.

Figure 5. Time evolution of fluid–bed interface dynamics. (a) Root mean square of the fluctuation of the streamwise-averaged $\sigma _{h,z}/D$ (solid lines) and spanwise-averaged fluid–bed interface $\sigma _{h,x}/D$ (dashed lines). (b) Mean lateral ridge wavelength $\lambda _{h,z}/{H_f}$. The black solid line marks the most amplified wavelength determined in the linear stability analysis of Colombini (Reference Colombini1993), while horizontal dashed lines of matching colour mark the lateral domain size ${L_z}/{H_f}$. (c) Streamwise particle flux ${{\langle {{q_{p,x}}} \rangle _{xz}}}/q_{ref}$ as a function of time. The reference particle flux $q_{ref}$ is computed based on the Wong & Parker (Reference Wong and Parker2006) version of the classical formula of Meyer-Peter & Müller (Reference Meyer-Peter and Müller1948), that is, $q_{ref}/({u_{g}}\, D)=4.93({\theta }(t)-{\theta _c})^{1.6}$, with the critical Shields number ${\theta _c}=0.034$ (Soulsby & Whitehouse Reference Soulsby and Whitehouse1997). (d) Mean particle flux density ${\langle {\phi {u_p}} \rangle _{xzt}}/(q_{ref}/D)$ as a function of the bed/wall distance. The inset shows the same quantity in the near-bed region, with the wall distance scaled with the particle diameter $D$. S250 (blue), M250 (green), M850 (yellow), L250 (red).

Figure 5(c) contains the temporal evolution of the streamwise particle flux ${{\langle {{q_{p,x}}} \rangle _{xz}}}(t)$ (cf. (B4a) in Appendix B) showing that the particle flow rate asymptotically approaches a quasi-steady state that is well estimated based on the empirical formula of Meyer-Peter & Müller (Reference Meyer-Peter and Müller1948) in the modified form of Wong & Parker (Reference Wong and Parker2006), viz.

(5.1)\begin{equation} q_{ref}/({u_{g}} D)=4.93({\theta}(t)-{\theta_c})^{1.6}, \end{equation}

where we have chosen a critical Shields number ${\theta _c}=0.034$ (Soulsby & Whitehouse Reference Soulsby and Whitehouse1997). This further supports earlier findings in the context of transverse patterns by Kidanemariam & Uhlmann (Reference Kidanemariam and Uhlmann2017) that relation (5.1) reasonably well describes the mean particle flow in the presence of sediment bedforms, even though it was actually developed for flow over a flat sediment bed. It should be stressed that the apparent ’overshoot’ of ${{\langle {{q_{p,x}}} \rangle _{xz}}}(t)/q_{ref}$ during the initial 20 bulk time units is not a feature of the particle flux evolution ${{\langle {{q_{p,x}}} \rangle _{xz}}}(t)$ which increases in an essentially monotonic way in this time window, before it settles to a quasi-stationary plateau (not shown). Instead, it is observed that relation (5.1) underpredicts the actual particle flux in this initial transient regime owing to the different growth rates of the particle flux and the non-dimensional wall-shear stress, while it provides a good approximation in the quasi-steady interval in which both quantities vary only weakly.

All discussed simulations belong to the bedload-dominated regime, that is, particle transport focusses on a layer of thickness ${O}(D)$ above the bed in which particles are transported by rolling, sliding and small jumps (saltation) without losing the contact with the bed for longer time intervals (van Rijn Reference van Rijn1984). This is verified in figure 5(d) which shows the wall-normal profile of the mean particle flux density ${\langle {\phi {u_p}} \rangle _{xzt}}(y)$ (cf. (B6)), highlighting that most of the transported particles indeed remain in the near bed region. It becomes evident from the wall-normal expansion of the flux density that the higher Shields number ${\theta }$ in case M850 leads to a thickening of the bedload layer by around $1D$ compared with the remaining cases. We conclude that in this latter case, a layer of $0.25{H_f}$ thickness above the interface is characterized by a non-negligible particle flux.

5.2. Turbulent mean flow

The presence of the mobile sediment severely modifies the mean fluid flow characteristics. Figure 6(a) shows the modification of the mean fluid velocity profile ${{\langle {\boldsymbol {u}_f} \rangle _{xzt}}}^{+}(y)$ compared with that of the single-phase cases. It is seen that the mean shear $\mathrm {d}{{\langle {\boldsymbol {u}_f} \rangle _{xzt}}}/\mathrm {d} y$ reduces to negligible value when approaching the mean fluid–bed interface at ${\tilde {y}}=0$, whereas in the smooth-wall case the mean shear peaks at the wall and is responsible for the entire wall-shear stress. In the multi-phase cases, on the other hand, the shear at the location of the virtual wall is clearly dominated by the contributions stemming from particle–fluid interactions (Kidanemariam & Uhlmann Reference Kidanemariam and Uhlmann2017). The mean velocity ${{\langle {\boldsymbol {u}_f} \rangle _{xzt}}}$ increases significantly slower with increasing distance to the bed in all particulate cases, leading to a velocity deficit compared with the smooth-wall cases, an essential feature arising in turbulent flows over rough walls due to an enhanced friction (Flores & Jiménez Reference Flores and Jiménez2006; Chan-Braun, García-Villalba & Uhlmann Reference Chan-Braun, García-Villalba and Uhlmann2011). The mean shear consequently increases further away from the bed compared with the smooth-wall cases in order to maintain the constant bulk flow rate ${q_f}$, resulting in a different slope of the velocity profile in the log and outer layers. For the low Reynolds number cases, the deviation is rather weak and most of the effect of the sediment bed can be compensated by the conventional shift of the log-layer velocity profile ${{\langle {u_f} \rangle _{xzt}}}^{+}=1/0.41 \ln ({\tilde {y}^{+}}) + 5.2 - {\Delta U^{+}}$ (Jiménez Reference Jiménez2004) by an offset ${\Delta U^{+}}=3.59$ (M250) and ${\Delta U^{+}}=3.73$ (L250), respectively. For the high Reynolds number case M850, on the other hand, the deviation of the slope is more pronounced such that a simple vertical shift will not lead to a satisfactory match of the profiles.

Figure 6. (a) Wall-normal profiles of the mean velocity ${{\langle {u_f} \rangle _{xzt}}}^{+}$ as a function of ${\tilde {y}^{+}}$ (semi-logarithmic scaling). The dashed-dotted line in (a) shows the standard predictions for the linear viscous sublayer and logarithmic layer, respectively, ${{\langle {u_f} \rangle _{xzt}}}^{+}={\tilde {y}^{+}}$ and ${{\langle {u_f} \rangle _{xzt}}}^{+}=1/\kappa \ln ({\tilde {y}^{+}})+B$ ($\kappa =0.41$, $B=5.2$). (b) Wall-normal profiles of the mean turbulence intensities $u_{rms}$ (solid lines), $v_{rms}$ (dashed lines), $w_{rms}$ (dotted lines). Here, M250 (green), M850 (yellow), L250 (red), M650_smooth (grey), L250_smooth (black).

The modifications of the mean velocity profile have direct implications for the intensity and distribution of the Reynolds stresses, which are depicted in figure 6(b) exemplary for cases M250, M850 and M650_smooth, respectively. Recall that the presented statistics in the multi-phase cases cannot be taken as fully converged as the observation interval is quite short and in addition highly transient. While we can assume to have captured enough statistics for the short-living small scales in the near-wall region, deviations between single- and multi-phase simulations close to the free surface are expected to be a consequence of the restricted averaging interval ${T_{obs}}$ with respect to the lifetime of the largest scales. Figure 6(b) indicates that in the low Reynolds number particle-laden cases, the typical buffer-layer peak of the streamwise normal stress ${{\langle {{\overline {u}}_f{\overline {u}}_f} \rangle _{xzt}}}$ is visibly reduced but still detectable at a comparable wall-normal position (data for L250 not shown), whereas it is completely absent in the high Reynolds number case M850 as a consequence of the intense particle motion in the latter case which destroys the near wall regeneration cycle in a similar way as in the case of flows past fully rough walls (Jiménez Reference Jiménez2004; Flores & Jiménez Reference Flores and Jiménez2006; Mazzuoli & Uhlmann Reference Mazzuoli and Uhlmann2017). In the low Reynolds number cases, on the other hand, particles are smaller by a factor of three in terms of wall units ($D\approx 10{\delta _{\nu }}$), which has been observed to be sufficiently small for particles to be transported by the buffer-layer structures rather than to destroy the near-wall cycle (Kidanemariam et al. Reference Kidanemariam, Chan-Braun, Doychev and Uhlmann2013). Deviations for the remaining normal stresses ${{\langle {{\overline {v}}_f{\overline {v}}_f} \rangle _{xzt}}}$ and ${{\langle {{\overline {w}}_f{\overline {w}}_f} \rangle _{xzt}}}$ are also observable, but are clearly weaker and visually restricted to a layer of thickness ${\tilde {y}}\approx 0.11{H_f}$ (${\tilde {y}}\approx 3D$) above the bed in the low Reynolds number cases and ${\tilde {y}}\approx 0.17{H_f}$ (${\tilde {y}}\approx 4.5D$) in the high Reynolds number case, respectively.

In the following, we scrutinize the impact of the mobile sediment bed on the distribution of the kinetic energy among the different scales away from the wall. To this end, the premultiplied streamwise energy spectra of cases M250 and M850 are analysed in figure 7 in the centre of the clear-fluid region (${\tilde {y}/{H_f}}=0.5$) for smooth-wall and particle-laden simulations. The instantaneous streamwise energy spectra at a given wall-normal location $y$ is introduced as

(5.2)\begin{equation} {\phi_{uu}}({k_x},y,{k_z},t)={\hat{u}_f}{\hat{u}_f^{*}}, \end{equation}

where ${\hat {u}_f}({k_x},y,{k_z},t)=\mathcal {F}({u_f^{\prime }}({\boldsymbol {x}},t))$ is the Fourier transform of the fluctuating field in the two wall-parallel directions, with the asterisk indicating complex conjugation. ${k_x}=2{\rm \pi} /{\lambda _x}$ and ${k_z}=2{\rm \pi} /{\lambda _z}$ are the streamwise and spanwise wavenumber–wavelength pairs. Velocity spectra for the remaining velocity components are defined accordingly.

Figure 7. Time-averaged premultiplied streamwise energy spectra ${k_x} {k_z} {{\langle {\phi _{uu}} \rangle _{t}}}({\lambda _x},y,{\lambda _z})\,{H_f}^{2}/{u_\tau}^{\!2}$ at the reference wall-normal location ${\tilde {y}/{H_f}} = 0.5$ for cases (a) M250, (b) M850. Coloured isolines are 0.2(0.2)0.6 times the maximum value of the respective energy spectrum, while grey-shaded areas represent the same quantity evaluated for the smooth-wall reference simulations (a) L250_smooth and (b) M650_smooth, respectively. The streamwise and spanwise domain periods ${L_x}$ and ${L_z}$ of the particle-laden simulations are highlighted by dashed lines in the respective colours.

Apart from some expected fluctuations due to the limited statistics, one may conclude a good agreement between the spectral patterns in both flow configurations. As expected, the medium domains cover only part of the full spectra, but this part is again observed to reasonably well collapse with that of the smooth-wall reference simulations. The good match between spectra in single- and multi-phase simulations is an indication that the evolving sediment bed does not significantly alter the streamwise energy spectra away from the bed and that the large-scale streaks over ridges are in fact essentially the same structures as those observed over smooth walls. While it is established for flow over fixed roughness elements that the roughness effect remains restricted to a layer of several multiples of their height (Jiménez Reference Jiménez2004,  Reference Jiménez2018; Mazzuoli & Uhlmann Reference Mazzuoli and Uhlmann2017) and does not strongly affect the organization of the large scales outside the roughness, a similar statement was not evident for the current case.

The wall-normal variation of the energy spectra is clearly identifiable in figures 8 and 9, where the premultiplied streamwise-integrated energy spectra are presented as functions of the lateral wavelength ${\lambda _z}$ and the wall-normal distance ${\tilde {y}}$ in inner and outer scales, respectively. Note that we have scaled each spectrum by its maximum value over all wall-normal locations to highlight how the energy distribution over the wall-normal direction is modified through the mobile sediment. The presented energy spectra are consistent with our claim that the relative particle size ${D^{+}}$ together with the high particle transport rate in the bedload layer of case M850 cause the immediate breakdown of the buffer-layer regeneration cycle. We observe that the near-wall peak in the streamwise spectra at ${\lambda _z}^{+}\approx 100$ associated with the buffer-layer streaks is still present for the low Reynolds number case M250 and kinetic energy is concentrated around this region whereas it is absent in case M850, in which the kinetic energy peaks further away from the wall outside the buffer layer. As indicated by the dashed line, most of the energy agglomerates in this latter case in regions above the bedload layer with less intense particle transport. A direct comparison with the smooth-wall reference case shows that this region is connected to the already existing second energetic peak at ${\lambda _z}/{H_f}\approx 1-1.5$ (${\lambda _z}^{+}\approx 1000$) that reflects the large-scale streaks at higher Reynolds numbers.

Figure 8. Time-averaged and streamwise-integrated premultiplied energy spectra ${k_z} \int {{\langle {\phi _{uu}} \rangle _{t}}}({\lambda _x},y,{\lambda _z}) \,\mathrm {d}{k_x} {H_f}^{2}/{u_\tau }^{2}$ as a function of the inner-scaled spanwise wavelength ${\lambda _z}^{+}$ and wall distance ${\tilde {y}^{+}}$, respectively. Cases (a) M250, (b) M850. Coloured isolines are 0.2(0.2)0.6 times the maximum value of the respective energy spectra, while grey-shaded areas indicate the same quantities determined for the smooth-wall reference simulations (a) L250_smooth and (b) M650_smooth, respectively. The mean fluid height ${H_f}^{+}={Re_\tau }$ and the spanwise domain period ${L_z}^{+}$ of the particle simulations are marked by coloured dashed lines. The dashed black line refers to the wall-normal distance at which the mean particle flux density ${\langle {\phi {u_p}} \rangle _{xzt}}$ attains its maximum (cf. figure 5d).

Figure 9. Time-averaged and streamwise-integrated premultiplied energy spectra for cases M850 and M650_smooth as in figure 8, but with the lateral wavelength ${\lambda _z}/{H_f}$ and the wall-normal distance ${\tilde {y}/{H_f}}$ scaled in outer units. An additional reference height is added in form of the crest height of the mean fluid–bed interface averaged over $t/{T_b}\,\in \,[40,85]$ (black solid line); (a) streamwise ${{\langle {\phi _{uu}} \rangle _{t}}}$ and (b) wall-normal spectra ${{\langle {\phi _{vv}} \rangle _{t}}}$.

Turning our attention to the wall-normal energy spectra depicted in figure 9(b), we first conclude that the medium domains considered here are just wide enough to accommodate the major fraction of the wall-normal energy spectra and hence are at the limit of being minimal near the free surface (Jiménez Reference Jiménez2013a), which further supports our observations in § 3. In the low Reynolds number case M250, the spectra of single-phase and particle laden simulations again almost collapse, highlighting the fact that the distribution of energy over the scales and the wall-normal layers is essentially the same as over the smooth wall (plot not shown). In the high Reynolds number case M850, on the other hand, the wall-normal kinetic energy is significantly reduced throughout the bedload layer and the region of intense vertical motion is found above the latter. However, the general ellipsoidal shape and the inclination angle are maintained. The contours of spectral energy density almost collapse in the vicinity of the free surface and for large wavelengths, whereas the intense energy region in the bulk is compressed, restricted by the bedload layer that damps the wall-normal energy even before the bed is reached, quite similar to the streamwise spectra.

5.3. Large-scale flow organization

Let us now turn to the dynamics of the large-scale flow structures in physical space. Figure 10 shows the space–time evolution of the streamwise velocity fluctuations with respect to the streamwise average, ${u_f^{\prime \prime }}(z,t)$. Note that streamwise averaging filters out lateral meandering of structures while retaining the footprint of the large streamwise-elongated structures.

Figure 10. Space–time plot of the streamwise-averaged fluctuations of the streamwise velocity component ${u_f^{\prime \prime }}/{u_\tau }$ extracted at ${\tilde {y}/{H_f}}=0.5$. Blue and red regions refer to streamwise-averaged zones of low and high streamwise velocity, respectively. Cases (a) M250, (b) M650_smooth, (c) M850, (d) L250.

In fact, the streamwise-averaged high and low-speed regions are found to develop comparably straight in the space–time plane, indicating that there is only weak lateral migration of these zones. It is noteworthy that there is no qualitative difference in this point between the ridge-featuring and smooth-wall cases during the current limited time interval, showing that even laterally homogeneous smooth channels are able to maintain a substantial spanwise modulation of the mean flow profile over an intermediate time scale as a consequence of the long-living large-scale streaks (Jiménez Reference Jiménez2013a), which is expected to disappear for sufficiently long time averaging intervals leading to statistical homogeneity. We shall discuss the role of this intermediate time scale dynamics in the creation of ridges below. Note that a comparison between the medium domain simulations and case L250 in figure 10(d) verifies that the observations are not a consequence of the limited box width ${L_z}/{H_f}\approx 3$, but that comparable regular organization of the large-scale streaks over intermediate time intervals similarly occur in domains as wide as ${L_z}/{H_f}\approx 16$. Also, it shows that the large-scale streaks are of significant length, as they are clearly detectable even in streamwise averages over long streamwise extents ${L_x}/{H_f}\approx 12$.

Intermittently, the amplitude of the high- and low-speed regions is seen to reduce as, for instance, around $t=40{T_b}$ in case M650_smooth (figure 10b) or in a period $10\lesssim t/{T_b} \lesssim 20$ in case M850. It appears likely that at least some of these ‘events’ are related to regularly occurring bursting events in which the large-scale streaks bend and eventually break, just as their viscous counterparts in the buffer layer do (Flores & Jiménez Reference Flores and Jiménez2010). Interestingly, there are events after which large-scale streaks recover at practically the same lateral position, while in other situations as after $t=45 {T_b}$ in case M250 (figure 10a), high- and low-speed regions reorganize and partly even change the number of large-scale streaks, as it is here the case. We will investigate the implications of such ‘events’ on sediment ridges in figure 12 below.

Throughout the previous sections, it has become clear that the large-scale streaks over ridges are no other structures than those found in smooth-wall turbulence, and we have in fact observed these structures throughout the entire ridge evolution interval. For homogeneous roughness, Flores, Jiménez & Del Álamo (Reference Flores, Jiménez and Del Álamo2007) similarly report that the self-similar vortex clusters and the associated velocity streaks are not significantly modified outside the roughness layer. We have, however, not explicitly discussed the question of cause and effect: based on the results seen so far, it seems reasonable that the spanwise heterogeneity of the flow in form of the well-organized large-scale low and high-speed streaks generates a laterally varying bed-shear stress that leads accordingly to spanwise heterogeneous erosion from the very beginning, hence implying that the large-scale streaks trigger the first ridges and not vice versa. Note that this does not rule out the possibility of a feedback of developed ridges on the mobility or meandering tendency of the large-scale streaks in later stages of their lifetime.

Large-scale streaks are entirely turbulent structures and as such not as smooth as their viscous counterparts in the buffer layer (Jiménez Reference Jiménez2013a). Also, large-scale instantaneous quasi-streamwise rotations are almost impossible to detect by the classical gradient-based vortex detection methods such as the $\lambda _2$-criterion of Jeong & Hussain (Reference Jeong and Hussain1995) as the velocity gradients are rather weak, as opposed to the near-wall quasi-streamwise vortices. Hence, for the understanding of the large-scale processes it is often customary to study fields in which small-scale motions are filtered out, either indirectly as in the context of ensemble averaging (see, for instance, Del Álamo et al. Reference Del Álamo, Jiménez, Zandonade and Moser2006; Lozano-Durán et al. Reference Lozano-Durán, Flores and Jiménez2012) and in some sense in our own streamwise-averaging procedure, or directly by means of instantaneous flow field filtering, for instance, with a Gaussian filter as in Motoori & Goto (Reference Motoori and Goto2019,  Reference Motoori and Goto2021). The advantage of direct filtering is that there is no information loss for all scales larger than the filter widths. In the current study, we therefore extend our analysis by applying a Gaussian cutoff filter to the flow field. The low-pass filtered field of the $i$th velocity component is obtained by the following convolution of the field with an anisotropic Gaussian kernel (Lozano-Durán, Holzner & Jiménez Reference Lozano-Durán, Holzner and Jiménez2016):

(5.3)\begin{equation} {{\widetilde{u}}_{f,i}}({\boldsymbol{x}}) = \iiint_{V} G \cdot {u_{f,i}}({\boldsymbol{x}}-{\boldsymbol{x}}^{\prime}) \exp \left( - \left( \dfrac{{\rm \pi} x^{\prime}}{{\varDelta_x}} \right)^{2} - \left( \dfrac{{\rm \pi} y^{\prime}}{{\varDelta_y}} \right)^{2} - \left( \dfrac{{\rm \pi} z^{\prime}}{{\varDelta_z}} \right)^{2} \right)\mathrm{d} x^{\prime} \,\mathrm{d} y^{\prime} \,\mathrm{d}z^{\prime}. \end{equation}

Therein, $\varDelta _i$ ($i=x,y,z$) is the filter cutoff width in the three Cartesian directions and ${\boldsymbol {x}}^{\prime }$ is the inner-convolution coordinate. The boundary conditions are treated similarly as in Lozano-Durán et al. (Reference Lozano-Durán, Holzner and Jiménez2016), i.e. the flow field is periodically extended in the wall-parallel directions, while it is mirrored vertically at the bottom wall, reversing the sign of the wall-normal velocity component to take care of the fundamental symmetries. The free surface can be treated in a completely analogous way as the symmetries along this plane resulting from the free-slip boundary condition do not differ from that at a smooth wall. Here, $V$ is the filter volume and the constant coefficient $G$ is normalized to ensure that the integral of the kernel over $V$ is unity. Note that the filtered flow field in a layer of thickness ${\varDelta _y}/{H_f}$ will be affected by the velocity field inside the bed. In the following, however, we will study only filtering results in regions sufficiently away from this wall-normal region or for single-phase flow simulations, in which the above effect is naturally absent. If not otherwise stated, we conventionally use filter widths $[{\varDelta _x},{\varDelta _y},{\varDelta _z}]/{H_f}=[0.6,0.2,0.3]$ that obey the typical aspect ratio of the vortex clusters as reported by Del Álamo et al. (Reference Del Álamo, Jiménez, Zandonade and Moser2006). In the following, we will also use a two-dimensional analogue of (5.3) in wall-parallel $(x,z)$-planes whenever we are mainly interested in the structure's wall-parallel extensions, using filter widths $[{\varDelta _x},{\varDelta _z}]/{H_f}= [0.6,0.3]$.

In the previous paragraphs, it has been discussed that the lateral organization of large-scale high- and low-speed velocity streaks and the associated ejection and sweep events may cause a laterally varying instantaneous shear stress along the bed or the lower wall, respectively, as recently discussed in the hydraulic context by Bagherimiyab & Lemmin (Reference Bagherimiyab and Lemmin2018). Underlying this hypothesis is today's understanding that the near-wall region is not exclusively populated by short scales of order ${O}(100{\delta _{\nu }})$, but that it also features large scales of order ${O}({H_f})$ that scale in outer units and that are correlated with the structures away from the wall (Del Álamo & Jiménez Reference Del Álamo and Jiménez2003). Even though the near-wall regeneration cycle itself is capable of functioning autonomously (Jiménez & Pinelli Reference Jiménez and Pinelli1999), large-scale influences have nonetheless been observed to modulate regions even close to the wall (Jiménez, Del Álamo & Flores Reference Jiménez, Del Álamo and Flores2004). Toh & Itano (Reference Toh and Itano2005) stated that there is a frequent interplay between small near-wall structures and large structures away from it, the former necessarily agglomerating below the large-scale ejections due to continuity (Jiménez Reference Jiménez2018).

In figure 11, the organization of buffer-layer structures with dimensions ${O}(100{\delta _{\nu }})$ into larger essentially ${H_f}$-scaled streaks and the corresponding outer large-scale streaks at ${\tilde {y}/{H_f}}=0.3$ can be observed for the smooth-wall case M650_smooth in visualizations of unfiltered and corresponding spatially filtered fields. Note that qualitatively similar results are found for the flow organization over the sediment beds just before the particle release. To avoid confusion with the mean wall shear ${\tau _w}$, let us introduce the instantaneous bottom shear stress ${\tau _b}(x,z,t)={\nu _f} \,\mathrm {d}{u_f}/\mathrm {d}{\tilde {y}}|_{{\tilde {y}}=0}$ for the following analysis. In figure 11, the high-speed buffer-layer streaks whose position can be inferred by their induced zones of locally increased shear are seen to cluster in two streamwise-elongated streak-like structures, where one spans the entire box length in a slightly meandering way, while the other is roughly half as long and laterally shifted by an offset somewhat larger than ${H_f}$. These two large-scale structures are found to lay essentially below the corresponding large-scale high-speed streaks at the upper end of the logarithmic layer, ${\tilde {y}/{H_f}}=0.3$, and are of comparable lateral extension, even though the streamwise extension differs in particular for the patchy structure with centre at $z/{H_f}=0.5$. In order to quantify the correlation of the boundary shear with the flow field at each height of the channel in M650_smooth in the streamwise-averaged framework, let us define the instantaneous two-point correlation coefficient as follows:

(5.4)\begin{equation} {\rho_{u,\tau}}(y,t)={\langle {{\tau^{\prime\prime}_b}(z,t) {u_f^{\prime\prime}}(y,z,t)} \rangle_{z}}/[ {\langle {({\tau^{\prime\prime}_b}(z,t))^{2}} \rangle_{z}} \, {\langle {({u_f^{\prime\prime}}(y,z,t))^{2}} \rangle_{z}}]^{1/2}. \end{equation}

It turns out that, in accordance with the previous findings, regions of high and low momentum over the entire logarithmic layer are highly correlated with the corresponding wall shear regions along the bottom wall, with a mean value of the correlation in the channel centre of ${\langle {{\rho _{u,\tau }}} \rangle _{t}}({\tilde {y}/{H_f}}=0.5) \approx 0.5$. Note that the restricted periodic domain length makes it hard to visually detect a streamwise phase shift between the wall-shear object and the log-layer streaks which is, however, expected to exist as a consequence of the varying mean shear and thus propagation speed between the two wall-normal locations.

Figure 11. (a,b) Instantaneous wall-shear stress fluctuations ${\tau ^{\prime }_b}(x,z,t)/({\rho _f} {u_\tau }^{2})$ for the smooth-wall simulation M650_smooth (a) without applied filtering and (b) after a two-dimensional Gaussian cutoff filter has been applied in the two homogeneous directions with filter widths ${\varDelta _x}=3y_{ref}$ and ${\varDelta _z}=1.5y_{ref}$, respectively. The reference value has been set to the height $y_{ref}=0.3{H_f}$. (c,d) Corresponding streamwise velocity fluctuation field $({u_f^{\prime }})^{+}$ extracted at ${\tilde {y}/{H_f}}=0.3$ (c) non-filtered and (d) filtered with the same filter size as in (b). The time at which all snapshots have been extracted corresponds to the second marker point in figure 14(b), $t/{T_b}=83$.

The observed correlation between log-layer streaks and the organization of the bottom shear stress is of direct relevance for the current simulations, as it provides a mechanism that couples the erosion rate along the bed with the large-scale structures. While it appears inevitable that large scales drive the particle motion in case M850 due to the absence of the near-wall cycle, the observed link between large-scale structures and the local wall shear might explain why the same ridge spacing is found for the low Reynolds number cases in which small scales are, in general, still able to transport particles (Kidanemariam et al. Reference Kidanemariam, Chan-Braun, Doychev and Uhlmann2013).

5.4. Large-scale streak–ridge interaction

The above arguments imply that particle erosion and hence bed topography are linked with the large-scale structures via their near-bed effects and the local wall shear. Indeed, it appears that the number of adjacent alternating high- and low-speed regions agrees reasonably well with the spanwise wavelength ${\lambda _z}/{H_f} \approx 1.1-1.3$ of the outer energetic peak determined from the streamwise energy spectra, which, in turn, is of striking similarity to the observed mean pattern wavelength ${\lambda _{h,z}}/{H_f}$. That, indeed, regions of high and low streamwise velocity are spatially correlated with the bed contour visualized in figure 12, in an exemplary manner for case M850. Here, the space–time evolution of the sediment bed fluctuation ${h_b^{\prime \prime }}(z,t)$ and the respective fluctuation of the streamwise velocity component at ${\tilde {y}/{H_f}}=0.5$ are repeated for convenience, supplemented by similarly presented data for the streamwise and spanwise particle flux components ${q_{p,x}^{\prime \prime }}$ and ${q_{p,z}^{\prime \prime }}$, respectively. The observed dependences can be summarized as follows: regions of high momentum which are the streamwise-averaged analogue of the large-scale streaks lead to spanwise alternating zones of locally increased erosion due to turbulent sweep structures that are located in the large-scale high-speed streaks. The wall-shear stress is locally increased and with it the local Shields number. Higher shear stress and erosion directly imply an increased streamwise particle flux ${q_{p,x}^{\prime \prime }}$ in these regions. The effect is further supported by the action of the lateral particle transport ${q_{p,z}^{\prime \prime }}$, which is predominantly directed in such sense that sediment is transported from the regions of high streamwise particle transport (${q_{p,x}^{\prime \prime }}>0$) to areas with lower streamwise particle flux, i.e. ${q_{p,x}^{\prime \prime }}<0$, as indicated by the red and blue lines in figure 12(d). The consequence is a local decrease of the sediment bed height in the form of a local trough. Corresponding opposite relations are found for low-speed regions, ridges and regions of low streamwise particle transport ${q_{p,x}^{\prime \prime }}<0$.

Figure 12. Space–time evolution in case M850 of (a) the sediment bed height fluctuation of the streamwise-averaged bed (identical to figure 4c); (b) the streamwise-averaged streamwise velocity fluctuation ${u_f^{\prime \prime }}/{u_\tau }$ at ${\tilde {y}/{H_f}}=0.5$ (identical to figure 10c); (c) the streamwise-averaged particle flux fluctuation ${q_{p,x}^{\prime \prime }}/q_{x,rms}$, with $q_{x,rms}(t)=({\langle {{q_{p,x}^{\prime \prime }}{q_{p,x}^{\prime \prime }}} \rangle _{z}})^{1/2}$. (d) Same data as in (c) are shown in the background as grey map, supplemented with red (blue) dots marking regions of lateral bed growth (decrease) that are regions with vanishing lateral particle flux ${q_{p,z}^{\prime \prime }}=0$ and negative (positive) lateral gradient $\partial _z {q_{p,x}^{\prime \prime }}(z,t)$.

The lateral particle flux is found to be two orders of magnitude smaller than the corresponding streamwise particle flux ${q_{p,x}^{\prime \prime }}$ in all simulations. It is therefore conceived that the initial sediment ridges readily appearing after the onset of particle motions are mostly a consequence of the laterally varying particle erosion in the streamwise direction, resulting in the comparably stronger streamwise particle flux. Lateral particle transport from the troughs to the crests which is due to the weaker lateral fluid motion then further supports the growth of these initial sediment ridges once the particle fluxes have reached a quasi-stationary level (cf. figure 5c).

The visually identified correlation between large-scale structures and the bed contour in the multi-phase simulations can be quantified with the aid of a two-point cross-correlation as defined in (5.4), replacing the wall-shear stress ${\tau ^{\prime \prime }_b}$ by the sediment bed height fluctuation ${h_b^{\prime \prime }}$. In agreement with the previous observations, it turns out that, over most of the simulation time, the bed contour shows a high instantaneous correlation with the flow field over most of the channel height reaching values of more than 0.6 at the bulk flow centreline ${\tilde {y}/{H_f}}=0.5$ (figure omitted). This is, however, not the case in the phase between $t=10{T_b}$ and $t=30{T_b}$, in which the bed contour is found to be essentially uncorrelated with the flow structures away from the bed. The missing correlation between the large-scale structures and the location of sediment ridges is also observed in figure 12(a,b); between $t=10{T_b}$ and $t=20{T_b}$, the amplitude of the streamwise-averaged velocity fluctuations reduces and the previously and subsequently observed intense regions temporarily ‘diffuse’ into several narrower regions. We have argued earlier that the large-scale streaks break up in this phase. Two large-scale streaks recover between $t=20{T_b}$ and $t=25{T_b}$, but at somewhat different lateral positions. Interestingly, the bed contour is seen to adapt to the new spanwise organization of the large-scale streaks, laterally migrating in time towards these spanwise positions, however, with a noticeable time delay.

Figure 13. Two-point cross-time correlations between ${u_f^{\prime \prime }}(y,z,t+{\delta t})$ and: (a) the sediment bed height fluctuations, $-{\rho _{u,h_b}^{t}}(y,{\delta t})$; (b) the velocity fluctuations at the reference height ${u_f^{\prime \prime }}({\tilde {y}/{H_f}}=0.5,z,t)$, ${\rho _{uu}^{t}}(y,{\delta t})$. The reference time and reference wall-normal position are indicated by a white cross. Red (blue) regions represent strong (weak) correlation of the compared quantities. Black dashed lines connect the maximum correlation values at each wall-normal distance ${\tilde {y}/{H_f}}$. Contours separating the coloured areas indicate values of 0(0.1)1 in all panels. Data are for case M850.

This time lag between the bed evolution and the large-scale streak organization is instructive as it indicates how information is propagating across the channel. In order to determine the time delay between the dynamics of the large-scale streaks and the deformation of the sediment bed contour, we introduce the two-time cross-correlation between the streamwise-averaged flow field and an arbitrary physical quantity $a$ as

(5.5)\begin{equation} {\rho_{ua}^{t}}(y,{\delta t})= {\langle {{u_f^{\prime\prime}}(y,z,t+{\delta t}) a^{\prime\prime}(z,t)} \rangle_{z}}/ [ {\langle {({u_f^{\prime\prime}}(y,z,t+{\delta t}))^{2}} \rangle_{z}} {\langle {(a^{\prime\prime}(z,t))^{2}} \rangle_{z}}]^{1/2}. \end{equation}

Figure 13(a) shows the cross-time correlation between the flow field and the bed contour, i.e. $a^{\prime \prime }={h_b^{\prime \prime }}(z,t)$, while figure 13(b) provides the cross-time correlation between the flow field at varying wall-normal positions and that at a reference height chosen as ${\tilde {y}/{H_f}}=0.5$, i.e. $a^{\prime \prime }={u_f^{\prime \prime }}({\tilde {y}/{H_f}}=0.5,z,t)$. Note that the correlation between bed and flow field is presented premultiplied with a global factor $-1$ to take account of the fact that a locally higher velocity implies a local decrease of the sediment bed height. In figures 13(a) and 13(b), a time lag is clearly identifiable. Figure 13(a) verifies our earlier observations that ridges align themselves with the location of the large-scale streaks with a time lag of around $10{T_b}$ (equivalent to approximately $0.9{H_f}/{u_\tau }$). A very similar time delay is observed in figure 13(b) between the dynamics of the flow field at the channel centre and that close to the bed. The observed time lag implies that the generation of the initial sediment ridges is indeed predominantly dictated by the organization of the large-scale streaks in the channel centre which interact in a ‘top–down mechanism’ with the fluid layers closer to the bed.

It also indicates that the existence of dominant large-scale velocity streaks in the core of the channel is of direct relevance to the stability and position of the underlying sediment ridges and troughs, respectively. In the following, we shall analyse such events during which a clear signature of high- and low-speed regions is missing, avoiding streamwise averaging for the moment to show that the modulation of ${u_f^{\prime \prime }}$ in such phases is indeed due to a break-up of the large-scale streaks. Figure 14 shows the time evolution of the most energetic Fourier modes of the streamwise velocity fluctuations, ${\phi _{uu}}({k_x},y,{k_z},t)$, in wall-parallel planes at ${\tilde {y}/{H_f}}=0.5$. For the sake of readability, let us refer to modes with wavelength pairs ${\lambda _x}={L_x}/i$ and ${\lambda _z}={L_z}/k$ ($i,k \in \mathbb {N}_0$) as the $(i,k)$-mode in the remainder. Note further that the shown curves are normalized by the total streamwise kinetic energy contained in the fluctuating field at that time and wall-normal position, respectively, such that all contributions sum up to unity. This is particularly noticeable regarding the transient nature of the currently studied systems: the insets show the absolute value of the streamwise ${\langle {{u_f^{\prime }}{u_f^{\prime }}} \rangle _{xz}}$ (black) and wall-normal energy ${\langle {{v_f^{\prime }}{v_f^{\prime }}} \rangle _{xz}}$ (red) contained in the same plane, indicating that in all sediment-laden cases the turbulent kinetic energy globally increases with time (at least initially) as a consequence of the increasing friction along the bed due to particle transport and bedform development.

Figure 14. Time evolution of the energy distribution between individual modes of the streamwise velocity spectrum ${\phi _{uu}}({k_x},{k_z},t)$ evaluated at a wall-normal position ${\tilde {y}/{H_f}}\approx 0.5$. Note that the spectrum at each time has been normalized such that the contributions of all individual modes sum up to unity. Note that the most dominant modes (those modes that carry more than $20\,\%$ of the total energy at least once during the observation time) are highlighted using the following colour scheme: blue line $(0,1)$-mode, red line $(0,2)$-mode, purple line $(1,2)$-mode. The smaller insets show the time evolution of the total streamwise and wall-normal fluctuation energy at the same wall-normal position, (black solid line) ${\langle {{u_f^{\prime }}{u_f^{\prime }}} \rangle _{xz}}$ and (red solid line) ${\langle {{v_f^{\prime }}{v_f^{\prime }}} \rangle _{xz}}$, normalized by their time-averaged mean. Cases (a) M250, (b) M650_smooth, (c) M850, (d) L250.

In the three medium-sized boxes, the range of wavelengths is limited by the box size and only a few large-scale modes can be accommodated, such that most of the time a large fraction of the total kinetic energy is contained in these few dominant modes that are highlighted in different colours. Here, we have classified those modes as dominant that contain more than $20\,\%$ of the total perturbation energy at least once during the observation interval. From the four dominant highly energetic mode pairs that are of relevance in the medium boxes, three feature the streamwise zero mode, i.e. the corresponding structures are of infinite streamwise length and do not appear in the premultiplied energy spectra.

These modes are in fact the alternating large-scale high- and low-speed streaks identified earlier, spanning the entire streamwise domain and are thus clearly detectable in the streamwise average. In the large box L250, on the other hand, there are no specially exposed modes, but energy is rather uniformly distributed among all modes. Comparing the dynamics of the individual modes in the medium boxes with the space–time plots in figure 10, the times during which the clear signature of high- and low-speed regions disappears correspond to instances in which there is temporarily no clear dominant mode in the spectra, in particular the energy in the infinitely long modes is markedly reduced. In case M250 (figure 14a), the earlier discussed change from two to only one pair of large-scale high- and low-speed streaks is clearly visible around $t\approx 50{T_b}$, accompanied by a reduction in the total energies ${\langle {{u_f^{\prime }}{u_f^{\prime }}} \rangle _{xz}}$ and ${\langle {{v_f^{\prime }}{v_f^{\prime }}} \rangle _{xz}}$, suggesting that at this time, the large-scale streaks indeed break up. Figure 14(b) shows a comparable behaviour for the smooth-wall case M650_smooth with the difference that after the breakdown of the infinitely long two-streak mode, the same mode recovers instead of being replaced by a single pair of large-scale streaks as in the previous case. The two states without and with a clear dominant mode indicated by the two symbols in figure 14(b) are visualized in physical space in terms of fluctuating streamwise velocity fields in figure 15(a,b), which also demonstrates that there is no strong lateral meandering. In the high Reynolds number sediment-laden case M850 provided in figure 14(c), however, the streamwise zero modes carry only slightly more energy than the remaining modes at the beginning of the simulation, but shortly after the particle release (between $t=10{T_b}$ and $t=20{T_b}$) even this slight plus in energy is lost and they cannot be differentiated from the remaining modes anymore. Again, this is in line with the space–time plots and the flow field visualization in figure 15(c) underlines the absence of a clear dominating large structure. Later, however, the mode that represents two pairs of elongated large-scale streaks increases first mildly, then more strongly until it settles at $t\approx 40{T_b}$ and dominates the spectra for the rest of the simulation. The dominance of this mode is also seen in figure 15(d).

Figure 15. Wall-parallel plane of the streamwise velocity fluctuation field $({u_f^{\prime }})^{+}$ extracted at ${\tilde {y}/{H_f}}\approx 0.5$ for cases (a,b) M650_smooth (c,d) M850 at two different times. The flow fields refer to the instances marked by filled circles in figures 14(b) and 14(c), respectively.

5.5. Large-scale streaks and mean secondary flow

Up to this point, we have seen that quasi-streamwise large-scale streaks characterize the appearance of the flow field in particular in the medium boxes where a considerable fraction of the energy resides in the infinitely long streamwise modes. Some authors have hypothesized that sediment ridges may ‘lock’ the spanwise position of these large-scale streaks such that they leave a footprint in the streamwise time average as up- and downwelling of the mean velocity profile together with mean secondary currents (Nezu Reference Nezu2005).

Having shown that, at least over intermediate time intervals, the large-scale streaks indeed propagate little in the lateral direction and observing that instantaneous spatial meandering appears to be rather weak, the assumption that mean secondary currents are the statistical footprints of large-scale streaks appears conclusive. Indeed, depth-spanning vortices are seen to populate the cross-section of all four periodic channels in figure 16, which provides mean flow field data averaged over short time intervals (between $20$ and $80$ bulk time units) and over the streamwise direction. The secondary mean flow is expressed in terms of the streamfunction ${\langle {\psi } \rangle _{xt}}(y,z)$, which is defined by the requirement that $\nabla _{2D} {\langle {\psi } \rangle _{xt}}(y,z)= (\partial _y {\langle {\psi } \rangle _{xt}},\partial _z {\langle {\psi } \rangle _{xt}})^{\textrm {T}} = (-{{\langle {w_f} \rangle _{xt}}},{{\langle {v_f} \rangle _{xt}}})^{\textrm {T}}$. For the particle-laden simulations, the conventionally discussed configuration of primary flow upwelling over ridge crests and downwelling over the troughs is observed with accordingly oriented secondary currents (Nezu & Nakagawa Reference Nezu and Nakagawa1993), representing the short-time averaged equivalent of the LMPs/HMPs studied by, for instance, Mejia-Alvarez & Christensen (Reference Mejia-Alvarez and Christensen2013). A comparison between the cases over mobile sediment beds and the smooth-wall simulation M650_smooth in figure 16(b) shows that, over adequately chosen intermediate time intervals of the order of the large-scale streaks’ lifetime, even statistically spanwise homogeneous simulations do generate a mean secondary flow pattern that is similar to the secondary currents next to developed ridges in figure 16(a,c,d).

Figure 16. Mean primary and secondary flow patterns for cases (a) M250, (b) M650_smooth, (c) M850, (d) L250. Note that for the sake of comparison, in panel (d) we have arbitrarily chosen a sub-domain of the entire cross-section of case L250 with a lateral width similar to ${L_z}/{H_f}$ in the medium domain cases. Isolines of the primary flow ${{\langle {u_f} \rangle _{xt}}}/{u_b}$ are shown in intervals 0(0.1)1.2, while the secondary flow pattern $({{\langle {v_f} \rangle _{xt}}},{{\langle {w_f} \rangle _{xt}}})^{\textrm {T}}/{u_b}$ is shown in terms of the secondary mean flow streamfunction ${\langle {\psi } \rangle _{xt}}$. Clockwise (counter-clockwise) secondary flow rotation is indicated by red (blue) contours. The time-averaged fluid–bed interface profile is indicated by a red curve. The time window over which data have been accumulated can be seen in figure 17. For the shown contours, the range $[\min _{y,z}{\langle {\psi } \rangle _{xt}},\max _{y,z}{\langle {\psi } \rangle _{xt}}]$ is divided into 20 equally spaced subintervals.

Figure 17. Instantaneous secondary flow intensity ${u_{\perp }}/{u_b}$ as a function of time (solid lines) and the corresponding value determined for the time-averaged fields in figure 16 (dashed lines) following the definition in (5.6). The length of the dashed lines indicates the time window over which the fields shown in the respective panels of figure 16 have been averaged: M250 (green), M650_smooth (grey), M850 (yellow), L250 (red).

In particular, the secondary flow cells in all cases share roughly the same spanwise spacing, which is equivalent to the preferential spanwise wavelength of the instantaneous large-scale streaks, ${\lambda _z}/{H_f}\approx 1-1.5$. An exception is case M250 for which in the observed time interval intermittently only one pair of large-scale low- and high-speed streaks exists, as has been seen earlier. Local mean up- and downflows are the statistical footprints of ejections and sweeps, occurring in the respective large-scale low- and high-speed streaks.

Also, we do not observe any significant quantitative difference in the amplitude of the secondary currents from case to case. Figure 17 shows the time-dependent secondary flow intensity ${u_{\perp }}(t)$ defined from the cross-plane-averaged kinetic energy contained in the cross-flow field $({\langle {{v_f}} \rangle _{x}},{\langle {{w_f}} \rangle _{x}})^{\textrm {T}}$, viz.

(5.6)\begin{equation} {u_{{\perp}}}(t) = \left[ \dfrac{1}{{L_z} {\langle {h_f} \rangle_{xz}}(t)} \displaystyle \int_{0}^{{L_z}} \displaystyle \int_{{\langle {h_b} \rangle_{xz}}(t)}^{{L_y}} ( {\langle {{v_f}} \rangle_{x}}^{2} + {\langle {{w_f}} \rangle_{x}}^{2} )\, \mathrm{d} y \,\mathrm{d}z\right]^{1/2}. \end{equation}

Note that a comparable formulation has been originally presented by Sakai (Reference Sakai2016) for quantification of sidewall-induced mean secondary flow in open ducts and has been adapted here to the multi-phase configuration. It is seen that ${u_{\perp }}(t)$ oscillates between $1.5\,\%$ and $2.5\,\%$ of the bulk velocity for all cases, with slightly lower values in case L250 as a consequence of the larger averaging ensemble in a sense that the wider domain allows accommodation of a larger number of secondary flows cells (by a factor of three to four). Note that the secondary flow intensity of the time-averaged cross-flow $({\langle {{v_f}} \rangle _{xt}},{\langle {{w_f}} \rangle _{xt}})^{\textrm {T}}$ in figure 16 is, as expected, systematically lower than their instantaneous counterparts. Interestingly, the latter values are of comparable size to the usually reported strength of sidewall-induced mean secondary flows in open ducts. For instance, Sakai (Reference Sakai2016) determined the intensity of secondary flow within a distance $z=1{H_f}$ from the sidewalls of a smooth-wall open duct in a range $1.1\,\%$ to $1.3\,\%$ of the bulk velocity.

Note that the observed secondary flow pattern disappears for conventional time and streamwise averaging over sufficiently long time intervals in the smooth-wall case, and it is only recovered in the context of conditional averaging over individual log- or outer-layer streaks (Lozano-Durán et al. Reference Lozano-Durán, Flores and Jiménez2012; Kevin et al. Reference Kevin, Monty and Hutchins2019a). Long-time prognosis of the situation over ridge-covered beds is only possible for the shortest simulation S250 (cf. figure 18 in Appendix A) where the observation time is not restricted through the development of transverse bedforms. Here, observations over longer time intervals of roughly $600{T_b}$ suggest that the observed ridges undergo continuous variations in form of merging, splitting or short lateral propagations, probably reacting to changes in the structure of the large-scale streak organization as seen before in the initial phase. The global trend, however, shows that the lateral positions at which ridge crests are found do not vary strongly in this narrow computational domain. In figure 18 it can be seen that the dominant ridge recovers at essentially the same lateral position even after an intermediate disappearance during the period between $t=400{T_b}$ and $t=500{T_b}$.