3.1. Identifying the critical parameter governing secondary flow polarity

In this subsection, we aim to study parameters in the multicolumn roughness configuration that govern the polarity of secondary flows. The roughness element arrangement considered involves multiple length scales of interest, as shown in figure 2. Here, the length scale $s_l$ refers to the local spanwise gap between adjacent elements of the same wider column, while $s_a$ refers to the spanwise gap between adjacent elements of the neighbouring wider columns. Since the elements are staggered in the streamwise direction as shown in figure 1(a), $s_a$ remains constant across all rows. The length scale $s_w$ refers to the total width of a column, defined as the maximum spanwise gap between elements within the same column. Although the left-hand column appears wider than the right-hand column in figure 2, $s_w$ is identical for both columns due to staggered arrangement, which can be appreciated from the mean flow field visualized in figure 1(b). The length scale $s_y$ refers to the spanwise distance between the centres of the wider columns, while $H$ is the height of the half-channel. In this paper, the lateral length scales confined to a single wider column ( $s_l$ and $s_w$ ) are normalized with $D$ , while the length scales corresponding to different wider columns ( $s_a$ and $s_y$ ) are normalized with $H$ . The streamwise spacing between element rows, denoted by $s_x$ , is normalized by $h$ .

Based on variations in these length scales, four different arrangements of roughness elements are obtained which are outlined in table 1. The case names reflect changes in only the second symbol, which corresponds to the length scale $s_l/D$ . As shown in the table, altering $s_l/D$ also affects parameters $s_w/D$ and $s_a/H$ , due to their inherent interdependence. Since $s_y/H$ and $s_x/h$ have previously been identified as critical parameters influencing the polarity of secondary flows (Yang & Anderson Reference Yang and Anderson2018; Joshi & Anderson Reference Joshi and Anderson2022), they are fixed at a constant value across all four cases to eliminate their influence from the current analysis.

Table 1. Suite of simulations for secondary flow reversal in § 3.1. The naming scheme for cases is defined in (2.3) and length scales ( $s_l$ , $s_w$ , $s_a$ and $s_y$ ) are shown in figure 2. Here $s_x$ is the streamwise gap between the element rows.

Figure 2. Schematic of length scales considered in the roughness element arrangement. Here $s_l$ , spanwise gap between adjacent elements of the same wider column; $s_w$ , width of the wider column; $s_a$ , spanwise gap between adjacent elements of different wider columns; $s_y$ , spanwise gap between centres of different wider columns.

Figure 3 presents the horizontally averaged first- and second-order statistics for the four configurations listed in table 1. In rough-wall flows, the mean streamwise velocity profile in the logarithmic region is often characterized using the aerodynamic roughness length $z_0$ and the displacement height $d$ , as

(3.1) \begin{equation} \frac {\langle \overline {u} \rangle }{u_\tau } = \frac {1}{\kappa } \ln\left (\frac {z-d}{z_0}\right)\! , \end{equation}

where $\kappa = 0.384$ is the von Kármán constant used in this study. Owing to the pronounced spanwise heterogeneities introduced by secondary flows, as will be demonstrated later through flow visualizations, prior studies have shown that horizontally averaged profiles in such flows generally deviate from the classical log-law scaling (Castro et al. Reference Castro, Kim, Stroh and Lim2021; Sathe & Giometto Reference Sathe and Giometto2024). Surprisingly, however, all four cases examined here exhibit a discernible logarithmic region, suggesting that the porous and elevated nature of the roughness elements helps preserve near-logarithmic scaling in the mean profile. The corresponding aerodynamic roughness lengths $z_0$ and displacement heights $d$ are listed in table 1. For reference, the equivalent sand-grain roughness height, computed as $k_s = z_0 e^{8.5\kappa }$ , is also reported in the same table (Womack et al. Reference Womack, Volino, Meneveau and Schultz2022). Furthermore, while the normalized velocity defect profiles collapse across all four cases, the total streamwise fluctuation profiles, comprising both temporal and dispersive contributions, exhibit a noticeable departure for case S2-5-2.

Figure 3. (a) Spatially averaged mean streamwise velocity, (b) velocity defect and (c) root mean squared velocity profiles for the cases considered in table 1. The black dashed line in (a) denotes the log-law slope. Here $\hat {z} = (z-d)/(H-d)$ .

Figure 4. Pseudocolour plot of vertical velocity for different cases mentioned in table 1, taken at a streamwise location coinciding with the elements. Here (a) S2-5-2, (b) S2.75-5-2, (c) S3.5-5-2, (d) S4-5-2. The black arrows indicate the vectors of spanwise and vertical velocity. The green line indicates contour of 95 % of horizontally averaged maximum mean streamwise velocity.

Figure 4 shows the pseudocolour plot of mean vertical velocity for these cases. The mean flow field is averaged over repeating patterns of roughness elements in the streamwise direction to improve symmetry. Since the structure of the mean secondary flows remains largely unchanged as the flow evolves within each repeating unit, we present only a representative $Y$ – $Z$ slice taken at the streamwise location coinciding with the element rows. Here, the obstacles appear as distinct elements with their own local updraughts and downdraughts divided at the middle of the rotor disk. This is due to the non-rotational nature of the actuator disk model used in this study.

It can be clearly seen that as the local spanwise gap between the elements ( $s_l/D$ ) is varied systematically, the response of the flow due to the roughness arrangement changes significantly. For the case S2-5-2, where elements are placed close to each other, half-channel height-scale secondary flows (here onwards referred to as delta-scale secondary flows) are distinctly seen. Such delta-scale secondary flows are also observed in numerical studies of flow over modelled roughness (Willingham et al. Reference Willingham, Anderson, Christensen and Barros2014; Anderson et al. Reference Anderson, Barros, Christensen and Awasthi2015a ; Chung, Monty & Hutchins Reference Chung, Monty and Hutchins2018), resolved roughness of different types (Hwang & Lee Reference Hwang and Lee2018; Yang & Anderson Reference Yang and Anderson2018; Stroh et al. Reference Stroh, Schäfer, Frohnapfel and Forooghi2020) and in experimental studies (Vanderwel & Ganapathisubramani Reference Vanderwel and Ganapathisubramani2015; Kevin et al. Reference Kevin, Monty, Bai, Pathikonda, Nugroho, Barros, Christensen and Hutchins2017), when $s_y/H \gtrsim 1.5$ . However, different polarity of these delta-scale secondary flows is also reported in these studies. In figure 4(a), clear downdraughts are seen over the elements for case S2-5-2. Consequently, strong updraughts are observed in the valley between the elements due to the secondary flows. As $s_l/D$ is systematically increased, the flow response for cases S2.75-5-2 and S3.5-5-2 differs significantly. The delta-scale vortices disappear in the mean flow for these configurations, and as a consequence, no prominent updraughts or downdraughts are observed coinciding with both the element columns symmetrically. As $s_l/D$ is further increased, the delta-scale secondary flows reappear, but this time, strong updraughts are seen over the elements for the case S4-5-2, accompanied by strong downdraughts in the valley. Observations of the vertical velocity at the domain centre show that the extent of the updraught in the valley, originating at the wall, decreases monotonically with increasing $s_l/D$ . Figure 4 also includes the 95 % contour of the horizontally averaged maximum mean streamwise velocity, highlighting the spanwise variations in the mean streamwise flow. For the case S2-5-2, the contour shows significant modulation of mean streamwise velocity in a way such that the higher momentum fluid flows over the elements and the flow slows down significantly in the valley between the elements. Significant reduction in the modulation of mean streamwise velocity is observed for S2.75-5-2, whereas a reversed modulation starts appearing for S3.5-5-2, where the flow slows down over the elements and speeds up in the valley. This reversed modulation is strengthened for the case S4-5-2, where the LMPs align with the updraughts over the elements. This clearly showcases that as the parameter $s_l/D$ is changed, the mean secondary rolls rearrange themselves to change the locations of LMPs and HMPs with respect to the roughness arrangement.

Figure 5 shows the TKE at Y–Z plane coinciding with the element row for the cases discussed above. For all the cases, maximum TKE lies in the vicinity of the roughness compared with the valley, regardless of the polarity of secondary flows. However, the way this TKE, generated near the roughness elements, is subsequently advected throughout the flow depends strongly on the polarity of the secondary circulation. For the case S2-5-2, due to strong downdraughts occurring over the element rows, the TKE generated by the elements is advected down and towards the valley. This can be appreciated by looking at the contour of 20 % of maximum TKE. Despite having two large sources in the spanwise direction that generate TKE, the 20 % contour does not show any significant hills or valleys for this case, indicating that the TKE is well-mixed in the spanwise direction. As $s_l/D$ is increased, hills in the contour of TKE start appearing coinciding with the element rows. This is expected, as the secondary rolls rearrange themselves to favour updraughts over the elements, which results in the generated TKE by elements being advected upwards. This creates significant spanwise heterogeneity in the magnitude of TKE which is sustained deep in the outer layer.

Figure 5. Pseudocolour plot of TKE for different cases mentioned in table 1, taken at a streamwise location coinciding with the elements. Here (a) S2-5-2, (b) S2.75-5-2, (c) S3.5-5-2, (d) S4-5-2. The green line indicates contour of 20 % of maximum TKE on the visualized plane.

While a reversal in secondary flow polarity is observed as $s_l/D$ is changed, it does not necessarily imply that $s_l/D$ is the critical parameter governing the polarity of secondary flows. From table 1, it is seen that changing $s_l/D$ also changes parameters $s_w/D$ and $s_a/H$ due to their interdependence, and each of these parameters could be critical for the flow response observed in figure 4. Furthermore, these observations are made for a locally staggered arrangement, which was chosen so that the parameter $s_l/D$ can be increased without the formation of local valleys within the same column. However, in a staggered arrangement, as in contrast with an aligned arrangement, shear emanating from the upstream element may freely impinge on the downstream element, which can have a critical impact on the value of production term in the TKE budget equation. As the imbalance between production, dissipation and transport terms in the TKE budget is shown to drive the secondary flow polarity (Hinze Reference Hinze1967; Anderson et al. Reference Anderson, Barros, Christensen and Awasthi2015a ; Hwang & Lee Reference Hwang and Lee2018; Joshi & Anderson Reference Joshi and Anderson2022), it is unclear, thus far, whether the flow response observed is specific to the local staggered arrangement used in the study. Hence, further investigation is needed to accurately isolate these parameters and determine their effect on the secondary flow polarity. The subsequent analysis in this subsection only focuses on the cases S2-5-2 and S4-5-2 as they show clear reversal in secondary flow polarity. While the long-time-averaged flow fields for cases S2.75-5-2 and S3.5-5-2 may give the impression that delta-scale secondary flows are absent, this is not an accurate reflection of the instantaneous flow dynamics. In these cases, no single pathway remains dominant over time; instead, frequent switching between pathways leads to a weakening of the secondary flow signatures in the mean visualization. A detailed discussion of this behaviour is provided in § 4.

To investigate the parameter critical to the reversal in secondary flow polarity, we use the energy transport tubes introduced in Meyers & Meneveau (Reference Meyers and Meneveau2013) to identify pathways taken by the kinetic energy of the mean flow (MKE) which is extracted by a given roughness element. Similar to the concept of classical mass transport tubes, energy transport tubes are defined such that no fluxes of energy exist over the tube mantle. The flux vector of MKE can be defined based on its transport equation,

(3.2) \begin{equation} \frac {\partial }{\partial x_j}\left (\overline {E}_{j}\right)=-\frac {1}{\rho }\frac {\partial \overline {u}_i \overline {p}_\infty }{\partial x_i}+ \overline {u_i^{\prime } u_j^{\prime }} \frac {\partial \overline {u}_i}{\partial x_j}-2 \nu \overline {S}_{\textit{ij}} \overline {S}_{\textit{ij}}+\frac {1}{\rho }\overline {u}_i \overline {F}_i \ , \end{equation}

where

(3.3) \begin{equation} \overline {E}_{j}= \overline {K} \ \overline {u}_j+\left (\frac {\overline {p}}{\rho } \delta _{\textit{ij}} + \overline {u_i^{\prime } u_j^{\prime }}-2 v \overline {S}_{\textit{ij}}\right) \overline {u}_i, \end{equation}

is the mean kinetic energy flux vector per unit mass, $\overline {K} = ( \overline {u}_i \overline {u}_i / 2)$ is the MKE, $\overline {S}_{\textit{ij}}$ is the mean strain rate tensor and $\overline {p}$ represents the part of pressure responsible for pressure gradient other than the externally applied pressure gradient ( $\overline {p}_\infty)$ . Equation (3.3) shows that the energy transport tube differs from classical mass transport tube due to the pressure, turbulence and viscous transport terms. As the Reynolds number for all the cases under investigation is very high, the contribution of viscosity is neglected. From (3.2), it is evident that the right-hand side contains both sources and sinks of energy. Therefore, unlike mass transport tubes, the total energy within an energy transport tube may not be conserved (Meyers & Meneveau Reference Meyers and Meneveau2013).

Before analysing the mean energy transport tubes for cases S2-5-2 and S4-5-2, we first consider simpler aligned cases to understand the preferential pathways for MKE entrainment for the range of $s_a/H$ and $s_y/H$ values under investigation, which are detailed in table 1. Meyers & Meneveau (Reference Meyers and Meneveau2013, see their figure 10) showcased variations in the shapes of energy tubes with varying $s_y/D$ . In their investigation for aligned configuration, they considered different values of $s_y/D$ , ranging from 4.54 to 10.5. They showcased that as $s_y/D$ is progressively increased, the elements favour spanwise entrainment of MKE compared with the vertical entrainment. When normalized by the half-channel height, i.e. considering $s_y/H$ instead of $s_y/D$ , this parameter range corresponds to 0.454–1.05 for their study. However, it is well known that as $s_y/H$ is increased further, delta-scale secondary flows appear and strong inhomogeneities in Reynolds stresses are observed associated with these secondary motions (Willingham et al. Reference Willingham, Anderson, Christensen and Barros2014; Vanderwel & Ganapathisubramani Reference Vanderwel and Ganapathisubramani2015; Yang & Anderson Reference Yang and Anderson2018). Thus, considering the contribution of advection and Reynolds stresses in (3.3), it is unclear whether the preferential spanwise entrainment of MKE with increasing $s_y/H$ continues for $s_y/H \gt 1$ . To investigate this, we run two cases with $s_y/H = 1 \ \text{and} \ 1.5$ , as detailed in table 2. These cases have aligned configuration with only one element per column in the spanwise direction to mimic the cases considered in Meyers & Meneveau (Reference Meyers and Meneveau2013). The streamwise gap between the elements is kept equal to the cases shown in table 1 to facilitate comparison with S2-5-2 and S4-5-2.

Table 2. Additional suite of simulations for energy tube evolution. The naming scheme is defined in (2.3). The length scales are the same as in table 1.

Figure 6 presents the pseudocolour plot of vertical velocity and modulation of mean streamwise velocity for these cases. For case A0-3-3, which corresponds to $s_y/H=1$ , weak downdraughts are observed over the elements and updraughts are observed in the valley. However, the effect of these circulations on the mean streamwise velocity is minimal. On the other hand, for the case A0-2-2, which corresponds to $s_y/H=1.5$ , strong delta-scale secondary flows are observed with HMPs and corresponding downdraughts coinciding with the elements. The effect of these circulations on the mean energy transport tubes is depicted in figure 7. For the case A0-3-3, the energy tube expands in both the vertical and the spanwise direction equally. We note that this shape differs from the case with $s_y/H=1.05$ considered in Meyers & Meneveau (Reference Meyers and Meneveau2013) due to the comparatively smaller streamwise gap considered in our study. For A0-3-3, initially, the growth rate is higher in the spanwise direction. However, as the upstream distance increases, the lower-half of the energy tube starts converging to a single line indicating diminishing spanwise growth, whereas the vertical growth is still noticeable until $x=-34s_xD$ . However, in contrast, the energy tube for A0-2-2 is noticeably narrower and taller. By analysing the maximum vertical and spanwise distance of the location of the MKE which is entrained by the elements, which is shown in figure 7(b,d), it is clear that the case A0-2-2 favours vertical entrainment of MKE over the spanwise entrainment compared with A0-3-3. This showcases that the observed trend of increasing spanwise entrainment with increasing $s_y/H$ for $s_y/H \lesssim 1$ is reversed due to secondary circulations as $s_y/H$ increases beyond unity. This is expected, as in this case, the secondary circulations align the higher-momentum fluid with the region above the elements, and strong downdraughts facilitate its downward transport. From this analysis, we can see that the shape and evolution of the mean energy transport tube encodes the information of polarity of secondary circulations. This is because as we move away from the elements, the contribution of pressure transport and turbulent transport in (3.3) decreases. Thus, the location of extremities of the tube far away from the elements indicates the direction of advection. Specifically, a narrower and taller tube indicates downdraught over the elements, which forms a circulation such that the flow is advected from the elements towards the valley (see, e.g. figure 4 a). In contrast, a shorter and wider tube indicates spanwise advection of MKE, which forms a circulation such that the flow is advected from the valley towards the elements, which results in updraughts over the elements (see, e.g. figure 4 d). For the cases such as A0-3-3 where the energy tube is equally wider and taller, formation of HMPs and LMPs is not observed which is reflected in weak modulation of mean streamwise velocity shown in figure 6(a). It should also be noted that, for both A0-3-3 and A0-2-2, the two parameters $s_a/H$ and $s_y/H$ are identical. Thus, this change in entrainment pathway of MKE with changing $s_y/H$ can also be interpreted as a change due to changing $s_a/H$ , which denotes the gap between adjacent elements in the spanwise direction.

Figure 6. Pseudocolour plot of vertical velocity for different cases mentioned in table 2, taken at a streamwise location coinciding with the elements. Here (a) A0-3-3, (b) A0-2-2. The black arrows indicate the vectors of spanwise and vertical velocity. The green line indicates contour of 95 % of maximum mean streamwise velocity.

Figure 7. Mean energy transport tubes for (a) A0-3-3 and (c) A0-2-2. Element location is shown by the black dashed line. Tube outlines are plotted at different upstream locations $x=-ns_xD$ , with $n=2, 4, \ldots , 34$ . The outline colour gradient indicates the upstream distance from the element, with lighter colours representing locations farther upstream. (b) Maximum height and (d) maximum spanwise distance at which the MKE is entrained by the element relative to the element centre. Here purple, A0-3-3; red, A0-2-2.

Figure 8 shows the mean energy transport tubes for cases S2-5-2 and S4-5-2 along with the maximum vertical and spanwise distance at which the MKE is entrained by the elements. To enable direct comparison, the two layouts are overlaid, and we focus on the column segment containing three elements, where the central elements coincide and the edge elements are offset according to their respective $s_l/D$ values. To maintain visual clarity, the transport tubes associated with the central element are omitted, as the focus here is on the entrainment behaviour of the edge elements. The figure shows a clear contrast in the entrainment mechanisms between the two cases. In case S2-5-2, the elements predominantly entrain MKE from above, indicating a vertically dominant entrainment pathway. In contrast, case S4-5-2 exhibits stronger spanwise entrainment. Interestingly, this shift in entrainment direction occurs despite both configurations sharing the same $s_y/H$ value, as their central elements are aligned at $y/D = 0$ . This indicates that a parameter other than $s_y/H$ is responsible for the different entrainment pathways observed for these cases. As discussed in the previous paragraph, a change in $s_a/H$ – which differs from $s_y/H$ in these configurations – may also influence the MKE entrainment pattern. Specifically, we noted that the cases may fall into different $s_a/H$ regimes, wherein the elements tend to favour either vertical or lateral entrainment. Table 1 highlights this distinction, showing that $s_a/H$ equals 1.25 for S2-5-2 and 1.00 for S4-5-2, which may place them in separate $s_a/H$ regimes. This shift in $s_a/H$ appears to be the driving factor behind the distinct entrainment patterns observed in S2-5-2 and S4-5-2, which ultimately result in the reversal of mean secondary flow polarity seen in figure 4. Thus, we hypothesize that $s_a/H$ is the critical parameter in a multicolumn roughness configuration that governs the polarity of secondary flows. In subsequent subsections, we will disprove $s_l/D$ , staggered configuration and $s_w/D$ as being crucial to the polarity reversal observed in figure 4, thereby reinforcing our hypothesis.

Table 3. Suite of simulations for § 3.2. The naming scheme for cases is defined in (2.3). The length scales are the same as in table 1.

Figure 8. (a) Mean energy transport tubes for S2-5-2 (purple) and S4-5-2(blue). Elements are shown by the black lines with locations for S2-5-2 and S4-5-2 being $y/D=\{-2, 0, 2\}$ and $y/D=\{-4, 0, 4\}$ , respectively. Tube outlines are plotted at different upstream locations $x=-ns_xD$ , with $n=3, 5, \ldots , 37$ . The outline colour gradient indicates the upstream distance from the element, with lighter colours representing locations farther upstream. (b) Maximum height and (d) maximum spanwise distance at which the MKE is entrained by the element relative to the element centre. Here purple, S2-5-2; blue, S4-5-2.

3.2. Disproving the local spanwise gap ( $s_l/D$ ) as the critical parameter

In this subsection, we provide further evidence supporting our conclusion that $s_a/H$ is the critical parameter by showing that $s_l/D$ does not play a decisive role. At a first glance, based on the plots shown in figure 4(a,d), one might speculate that increasing $s_l/D$ , which results in a significantly larger local spanwise gap between the elements in case S4-5-2 compared with S2-5-2, could influence the preferential fluid pathways locally. A narrower spanwise gap (S2-5-2) may require vertical entrainment of the flow, resulting in a global downdraught over the elements, which may not be necessary when a wider spanwise gap is available for the flow (S4-5-2). To disprove this intuition, we present four new cases along with cases S2-5-2 and S3.5-5-2 considered previously, as outlined in table 3.

Figure 9. Pseudocolour plot of vertical velocity for different cases mentioned in table 3, taken at a streamwise location coinciding with the elements. Here (a) S2-3-2, (b) S3.5-3-2, (c) S2-5-2, (d) S3.5-5-2, (e) S2-7-2, (f) S3.5-7-2. The black arrows indicate the vectors of spanwise and vertical velocity. The green line indicates contour of 95 % of maximum mean streamwise velocity.

Figure 9 presents the pseudocolour plot of vertical velocity for the cases listed in table 3. We first focus on the cases where $s_l/D = 2$ . Although these cases share the same $s_l/D$ , they have different $s_a/H$ values due to the increasing number of local columns. For cases S2-3-2 and S2-5-2, distinct downdraughts are observed above the elements. Additionally, the modulation of mean streamwise velocity shows formation of HMPs over the elements and LMPs in the valley. This behaviour contrasts with case S2-7-2, where intermediate secondary flows occur, where the symmetry of the secondary flows breaks down. The modulation of mean streamwise velocity is also markedly different for S2-7-2, as no distinct HMPs or LMPs coincide symmetrically with the elements. Thus, a clear transition is observed: from HMPs over the elements in S2-3-2 and S2-5-2 to intermediate secondary flows in S2-7-2, as the value of $s_a/H$ decreases, while keeping $s_l/D$ constant. Similarly, for cases with $s_l/D = 3.5$ , we observe a reversal in flow polarity. The flow transitions from HMPs over the elements in S3.5–3-2, to intermediate secondary flows in S3.5-5-2, and to LMPs over the elements in S3.5-7-2, as $s_a/H$ decreases. These six cases can together be ranked in decreasing order of $s_a/H$ as $\text{S2-3-2} \gt \text{S3.5-3-2} \gt \text{S2-5-2} \gt \text{S3.5-5-2} \approx \text{S2-7-2} \gt \text{S3.5-7-2}$ . In this sequence, we observe a progressive transition from HMPs over the elements, to intermediate secondary flows, to LMPs over the elements, independent of changes in the local spanwise gap. Thus, from the observations of different flow response for the same $s_l/D$ (cases S3.5-3-2, S3.5-5-2 and S3.5-7-2) and similar flow response for different $s_l/D$ (cases S2-3-2 and S3.5-3-2), we conclude that $s_l/D$ is not the critical parameter driving the reversal in secondary flow polarity observed in § 3.1. The variations in flow response observed with decreasing $s_a/H$ provide further reinforcement for our conclusion that $s_a/H$ is the key parameter governing this behaviour.

3.3. Disproving the local element alignment as the critical parameter

In this subsection, we showcase that the secondary flow polarity reversal observed in § 3.1 is not limited to the staggered configuration within a wider column. To prove this, three new cases are simulated with aligned configuration as outlined in table 4. In these cases, the variation in the parameter $s_l/D$ is kept modest compared with the cases considered in § 3.1 to prevent the formations of local valleys within a wider column. However, by including four elements per wider column, different regimes of $s_a/H$ are obtained despite the limited variation in $s_l/D$ . If the polarity reversal observed in § 3.1 depends solely on $s_a/H$ and is independent of the staggered configuration used, a corresponding shift in the secondary flow response is anticipated for these aligned cases as well.

Table 4. Suite of simulations for § 3.3. The naming scheme for cases is defined in (2.3). The length scales are the same as in table 1.

Figure 10 shows the flow response for these three cases. For the case A1.5-8-2, clear downdraughts are observed over the elements and updraughts are observed within the valley. However, these circulations have only a modest impact on the modulation of mean streamwise velocity, with the formation of weak HMPs over the elements and weak LMPs in the valley. As the ratio $s_a/H$ decreases, intermediate secondary flows are observed for case A2-8-2, where no significant updraughts or downdraughts coincide with the element locations. With a further reduction in $s_a/H$ , pronounced updraughts are observed over the elements, accompanied by the formation of LMPs over the elements and HMPs within the valley. Although fully symmetric secondary flows are not observed in this case, the flow response between cases A1.5-8-2 and A2.5-8-2 is markedly different. These findings indicate that as $s_a/H$ decreases, there is a transition from downdraughts to updraughts over the elements for the aligned configuration, consistent with the observations reported in §§ 3.1 and 3.2 for the staggered configuration. This demonstrates that the observed flow response in § 3.1 is not limited to the local alignment of the elements within a wider column. However, it is important to ensure that the alignment does not create significant local valleys within the wider column. If such valleys form, the very concept of a ‘wider column’ may break down, which may lead to a disruption in the expected flow behaviour.

Figure 10. Pseudocolour plot of vertical velocity for different cases mentioned in table 4, taken at a streamwise location coinciding with the elements. Here (a) A1.5-8-2, (b) A2-8-2, (c) A2.5-8-2. The black arrows indicate the vectors of spanwise and vertical velocity. The green line indicates contour of 95 % of maximum mean streamwise velocity.

3.4. Disproving the width of the wider column ( $s_w/D$ ) as the critical parameter

In this subsection, we discuss why the secondary flow polarity reversal observed in § 3.1 is not due to the width of wider columns. For this, we first refer the reader to the cases considered in table 2 and the corresponding energy tube visualization in figure 7. For these cases, there is no change in the width of the columns, as each column contains only a single element. Nevertheless, distinct patterns of MKE entrainment emerge when the parameter $s_a/H$ is varied. This demonstrates that even in the absence of changes in $s_w/D$ , different flow responses can be observed by varying $s_a/H$ .

Next, we discuss the nature of coupling between parameters $s_l/D$ and $s_w/D$ . For cases with an equal number of elements per row, as outlined in table 1, there exists a direct coupling between $s_l/D$ and $s_w/D$ . Specifically, any modification in $s_l/D$ necessitates a corresponding adjustment in $s_w/D$ , precluding independent variation of these parameters. Therefore, $s_w/D$ can be regarded as a local parameter, concerned only with a given wider column. In § 3.2, it is demonstrated that $s_l/D$ does not play a crucial role in determining secondary flow polarity. Given that $s_l/D$ influences only the local interactions within a wider column, this implies that, within the studied spanwise sparsity range (i.e. $s_y/H=1.5$ ), the elements’ local spanwise columnar interactions are not crucial to the global secondary flow response. This suggests that the different MKE entrainment patterns observed in figure 8 are not critically affected by these local interactions. Since $s_w/D$ also governs local interactions, it can be inferred that the distinct MKE entrainment pathways observed in this figure are not attributable to changes in $s_w/D$ . Therefore, $s_w/D$ is not a critical parameter governing the secondary flow polarity reversal observed in § 3.1. We further add that the coupling between $s_l/D$ and $s_a/H$ is of different nature and is dependent through the parameter $s_y/H$ . Due to this dependence on $s_y/H$ , $s_a/H$ governs interactions between neighbouring wider columns, rather than local interactions within a single wider column. Independent variations in $s_l/D$ and $s_a/H$ can be achieved by changing $s_y/H$ while maintaining the same number of elements per row. Consequently, while it is reasonable to conclude that $s_w/D$ is not a critical parameter – given that $s_l/D$ is non-critical – the same reasoning cannot be applied to $s_a/H$ .

Overall, the analysis presented in this section clearly demonstrates that variations in the parameter $s_a/H$ are the primary driver behind the reversal in secondary flow polarity observed in figure 4. This finding highlights that, in multicolumn configurations, roughness elements located at the edges of the columns can play a decisive role in setting the polarity of secondary flows. While the polarity reversal is evident in the long-time-averaged flow field, the extent to which these structures persist in time remains an open question. In particular, it is unclear whether the presence (figure 4 a,d) or absence (figure 4 b,c) of counter-rotating delta-scale secondary flows in the mean reflects fundamentally different instantaneous dynamics. We address this aspect in the following section.

4.1. Visualization of instantaneous flow field

In this subsection, we examine the instantaneous flow field behaviour of case S3.5-5-2. Figure 4(c) presents the long-time-averaged flow field for this case, and the prominent observations are described briefly here. Unlike the cases in figure 4(a,d), which exhibit symmetrical secondary flows, this case displays a pronounced asymmetry, characterized by a strong updraught in the left-hand column and a significantly weaker updraught in the right-hand column. Consequently, symmetrical delta-scale counter-rotating vortices are absent. The time-averaged field also reveals that an updraught remains confined near the wall at the domain centre, with a weak downdraught positioned above it. Additionally, the mean streamwise velocity contour indicates that an HMP aligns with the domain centre, while an LMP is positioned over the roughness elements.

Figure 11. (a,b) Pseudocolour plots of normalized streamwise velocity $(\tilde {u}/u_\tau)$ at different time steps, taken at a streamwise plane intersecting the element locations for the case S3.5-5-2. White circles indicate element positions. (c,d) Visualization of streamwise vortical structures using isosurfaces of signed swirling strength (4 % of maximum $\lambda _{ci}^2$ ). Panels (c) and (d) correspond to the same time steps as panels (a) and (b), respectively.

Figure 11(a,b) presents two realizations of the instantaneous streamwise velocity for this case, taken approximately 12 large eddy turnovers ( $H/u_\tau$ ) apart. The colour bar is capped so that the velocity values exceeding 80 % of the maximum appear in a uniform light shade, enhancing the visibility of HMRs. In figure 11(a), a dominant HMR is clearly visible at the centre of the domain, while low-momentum fluid remains near the roughness elements and extends upward into the outer layer. However, in figure 11(b), an LMR occupies the centre of the domain, penetrating deep into the outer layer. While the low-momentum fluid is always present near the roughness elements due to the drag they impose, these two realizations illustrate that instantaneous secondary circulations can organize themselves to either advect this low-momentum fluid upward or, at a different time step, switch circulation direction and transport it laterally at the centre – despite the roughness arrangement remaining unchanged. Thus, while an HMP appears at the domain centre in the long-time-averaged visualization, this does not imply a continuous presence of HMR at the centre throughout the simulation. Instead, the regions of low- and high-momentum alternate, and what we observe in the long-time-averaged field is the superposition of these alternating structures. Depending on their relative strength and persistence over time, one pathway may appear dominant in the long-time-averaged visualization.

Figure 11(c,d) shows the isosurfaces of two-dimensional signed swirl strength, which indicates the strength of local swirling in the streamwise direction. The swirl strength ( $\lambda _{ci}$ ) is computed as the imaginary component of the complex eigenvalue of the two-dimensional velocity gradient tensor, computed in the spanwise-wall-normal plane, and is signed based on the direction of vorticity (Zhou et al. Reference Zhou, Adrian, Balachandar and Kendall1999; Wu & Christensen Reference Wu and Christensen2006; Anderson et al. Reference Anderson, Barros, Christensen and Awasthi2015a ). The figure displays isosurfaces corresponding to 4 % of the maximum $\lambda _{ci}^2$ , a threshold selected to enhance the visualization of stronger circulations while minimizing background noise. In figure 11(c), most of the vorticity is concentrated near the elements, while the valley between the arrays exhibits a prominent region of low circulation. In contrast, figure 11(d) shows a noticeable increase in vorticity within the valley. Taken together with the streamwise velocity visualizations, these observations indicate that LMRs are associated with strong streamwise circulations, whereas HMRs tend to suppress them. This observation aligns with the findings in the literature, where it is reported that LMPs exhibit stronger random fluctuations, while HMPs remain comparatively calmer (Barros & Christensen Reference Barros and Christensen2014; Anderson et al. Reference Anderson, Barros, Christensen and Awasthi2015a ; Kevin et al. Reference Kevin, Monty, Bai, Pathikonda, Nugroho, Barros, Christensen and Hutchins2017; Vanderwel et al. Reference Vanderwel, Stroh, Kriegseis, Frohnapfel and Ganapathisubramani2019). Additional analysis of instantaneous vortices (not shown) reveals that, across all cases listed in table 1, the largest instantaneous vortex is approximately an order of magnitude smaller than the delta-scale vortex. This analysis also indicates that larger instantaneous circulations generally tend to be weaker, while stronger circulations are typically confined to smaller spatial regions. These findings confirm that delta-scale vortices do not exist instantaneously but instead arise as a consequence of time-averaging.

The observations presented in this subsection give rise to three key questions.

1. (i) For how long does a given region of low- or high-momentum persist at the centre?

2. (ii) How frequently do the regions switch between one another?

3. (iii) Is this behaviour universal, or does it depend on the roughness arrangement?

Exploring these questions is essential to understanding why some cases exhibit symmetrical secondary flows, while others appear disrupted. The next subsection investigates these aspects in greater detail, shedding light on the underlying flow dynamics.

4.2. Conditional averaging of instantaneous flow field

In this subsection, we investigate the conditionally averaged flow fields for the cases listed in tables 1 and 3. As discussed in the previous subsection, the HMR and the LMR can both occur at the domain centre at different instances, indicating reversals in the swirling direction of the secondary flows. Since downdraughts and updraughts serve as key indicators for the presence of HMRs and LMRs, respectively, we categorize the instantaneous flow into two distinct events, based on the relative occurrence of updraughts and downdraughts at the domain centre. For this conditional averaging, velocity data from eight spanwise locations are considered: four located closest to the centre of the domain, and two each near the edges in the spanwise direction. These edge locations also represent valleys between elements due to the periodic boundary condition. At each of these locations, vertical velocity values from the entire streamwise and vertical extent are included, resulting in approximately 885 000 grid points at which updraught and downdraught occurrences are evaluated for each time step. Consequently, two distinct averaged velocity fields are generated: one representing the event when updraughts outnumber downdraughts, and the other corresponding to the event when downdraughts predominate.

To formalize this classification, let $w(\boldsymbol{x}_i,t)$ be the instantaneous vertical velocity at grid point $\boldsymbol{x}_i$ and time $t$ , with $\mathcal{D}$ denoting the set of diagnostic grid points. At each time instant $t$ , the number of updraught and downdraught points is computed as

(4.1) \begin{equation} \begin{aligned} N_{\textit{up}}(t) &= \# \{ \boldsymbol{x}_i \in \mathcal{D} \, | \, w(\boldsymbol{x}_i,t) \gt 0 \}, \\N_{\textit{do}w\textit{n}}(t) &= \# \{ \boldsymbol{x}_i \in \mathcal{D} \, | \, w(\boldsymbol{x}_i,t) \leqslant 0 \}. \end{aligned} \end{equation}

The flow field at time $t$ is classified as updraught-dominated if $N_{\textit{up}}(t) \gt N_{\textit{do}w\textit{n}}(t)$ and as downdraught-dominated otherwise. The corresponding conditional averages are then defined as

(4.2) \begin{equation} \begin{aligned} \overline {u}_{i_{c+}}(\boldsymbol{x}) &= \frac {1}{T_{\textit{up}}} \sum _{\substack {t \, : \, N_{\textit{up}}(t) \gt N_{\textit{do}w\textit{n}}(t)}} \tilde {u}_i(\boldsymbol{x},t)\ , \\ \overline {u}_{i_{c-}}(\boldsymbol{x}) &= \frac {1}{T_{\textit{do}w\textit{n}}} \sum _{\substack {t \, : \, N_{\textit{do}w\textit{n}}(t) \ge N_{\textit{up}}(t)}} \tilde {u}_i(\boldsymbol{x},t)\ , \end{aligned} \end{equation}

where $T_{\textit{up}}$ and $T_{\textit{do}w\textit{n}}$ are the number of time instants satisfying each condition.

Figure 12. Conditionally averaged pseudocolour plots of vertical velocity at a streamwise location aligned with the elements for the cases listed in table 1: (a,b) S2-5-2, (c,d) S2.75-5-2, (e,f) S3.5-5-2, (g,h) S4-5-2. Panels (a), (c), (e) and (g) correspond to events where updraughts outnumber downdraughts, while panels (b), (d), (f) and (h) represent events dominated by downdraughts. The vertical bar on the left-hand side of each row indicates the percentage of time updraughts were dominant (red, with numerical value at the top) and the percentage of time downdraughts were dominant (blue, with numerical value at the bottom). Black arrows indicate the vectors of spanwise and vertical velocity. The green line indicates contour of 95 % of maximum mean streamwise velocity.

Figure 12 presents the conditionally averaged vertical velocity fields for the cases listed in table 1. Figure 12(a,c,e,g), corresponds to events where updraughts outnumber downdraughts at the domain centre, while figure 12(b,d,f,h), depicts events dominated by downdraughts. The vertical bar on the left-hand side of each row indicates the percentage of simulation time for which each specific event occurs, with the extent of red colour denoting dominant updraughts and blue representing dominant downdraughts. It is evident that in all cases – including those where symmetrical secondary flows with a dominant polarity were observed – the updraughts and downdraughts do not persist for the entire duration of the simulation. For instance, in case S2-5-2, an updraught exists at the domain centre for 93 % of the simulation time, while downdraughts dominate 7 % of the time. As a result, the conditionally averaged field shown in figure 12(a) closely resembles the long-time averaged field in figure 4(a), as the events with opposite polarity lack sufficient temporal persistence to leave a significant imprint on the long-time averaged flow field. As the parameter $s_a/H$ is decreased, moving to the case S2.75-5-2, events with dominant downdraughts become substantially more frequent, constituting 46.6 % of the simulation time steps, while updraughts persist for the remaining 53.4 % of the time. In periods where downdraughts dominate, the streamwise velocity clearly exhibits strong spanwise modulation, with HMPs appearing at the domain centre and LMPs aligned over the roughness elements. In contrast, when updraughts dominate, the modulation is notably weaker, with weak LMPs near the domain centre and weak HMPs near the elements. This observation is significant because, in the long-time averaged field shown in figure 4(b), an HMP appears at the domain centre, despite the updraught (and its associated LMP) persisting longer. This highlights that the long-time average may not always reflect the flow structures that are most frequently present throughout the simulation. Another interesting observation for this case is the reversal of flow pathways at the roughness elements. Specifically, figure 12(c) shows downdraughts clearly aligned with the elements, while figure 12(d) shows updraughts aligned over these arrays. This behaviour differs from the findings of Vanderwel et al. (Reference Vanderwel, Stroh, Kriegseis, Frohnapfel and Ganapathisubramani2019), who reported that for flow over streamwise-aligned ridges, LMRs remained anchored to elevated surface regions, extending deep into the boundary layer as narrow and meandering streaks. We attribute this discrepancy to differences in the physical nature of roughness elements and the presence of significant streamwise gap in our roughness configuration. Overall, despite both the conditionally averaged flow fields showcasing the presence of delta scale vortices, as no event has dominant presence throughout the simulation, the long-time averaged visualization shows seemingly disrupted secondary flows. Further decreasing $s_a/H$ results in case S3.5-5-2, where downdraughts become even more frequent, dominating approximately 65.5 % of the simulation duration compared with 34.5 % for updraughts. Similar to the previous case, downdraught-dominated events strongly modulate streamwise velocity, while updraught-dominated events exhibit only weak modulation. Interestingly, all delta-scale vortices observed in the conditional averages for this case rotate in a anticlockwise direction. This explains why a prominent anticlockwise delta-scale vortex appears clearly in the long-time averaged visualization (figure 4 c), even though the rest of the flow appears somewhat disrupted due to the lack of persistent dominance of one particular pathway. This is likely due to the dominance of one-sided large-scale rotational motions in this configuration, a feature previously identified by Kevin et al. (Reference Kevin, Monty, Bai, Pathikonda, Nugroho, Barros, Christensen and Hutchins2017). Finally, in case S4-5-2, where clear delta-scale secondary flows and dominant pathways appear in the long-time averaged visualization, downdraughts dominate at the domain centre for approximately 73 % of the simulation. However, the reversed secondary flow pattern remains less symmetric than in case S3.5-7-2 (figure 9 f), where downdraughts dominate approximately 95 % of the simulation time, resulting in a symmetrical long-time averaged flow pattern. Overall, as $s_a/H$ is progressively decreased, we observe a consistent shift towards downdraught dominance at the domain centre. The configurations transition from exhibiting a strong central updraught to a pronounced central downdraught, eventually reaching a level where this reversal is clearly evident in the long-time-averaged flow field.

Figure 13. Temporal persistence of updraught and downdraught events for cases: (a) S2-5-2, (b) S2.75-5-2, (c) S3.5-5-2, (d) S4-5-2, (e) S3.5-7-2. The red colour denotes a net updraught at the domain centre, while the blue colour denotes a net downdraught. The vertical bar on the right-hand side of each row indicates the percentage of time updraughts were dominant (red, with numerical value at the top) and the percentage of time downdraughts were dominant (blue, with numerical value at the bottom).

From the analysis so far, it is abundantly clear that for most of the secondary flow simulations, the updraughts and downdraughts, and the associated LMRs and HMRs, do not remain stable and may interchange with one another. However, it remains unclear whether these reversals occur in a periodic manner or follow a more chaotic, non-periodic pattern. To investigate this, figure 13 marks the starting time step where a particular updraught or downdraught event occurs and the marker’s extent on y axis indicates the duration for which that event persists. The analysis is conducted for the cases mentioned in table 1 along with the case S3.5-7-2 in table 3 to emphasize dynamics of a case where downdraughts dominate at the centre of the domain for the majority of the simulation.

For case S2-5-2, we observe that the updraughts at the domain centre are relatively stable and can persist continuously for extended durations. In contrast, downdraughts occur infrequently and in a non-periodic manner, and they are typically short-lived. For the intermediate cases shown in figure 13(b,c), the average temporal extent of both updraughts and downdraughts is significantly reduced. The vertical axis scale is adjusted individually for each subplot in the figure to enhance visual clarity. In these cases, the relative occurrence of updraughts and downdraughts at the domain centre appears roughly balanced, with the exception of two long-duration updraughts events in case S2.75-5-2. Due to the reduced persistence of individual events, these intermediate configurations undergo frequent, non-periodic reversals. As a result, the secondary flow patterns appear disrupted in the long-time-averaged visualizations. As $s_a/H$ decreases further, in cases such as S4-5-2 and S3.5-7-2, we find that downdraughts at the domain centre begin to stabilize and persist for longer durations compared with updraughts. With downdraughts sustaining for extended periods, reversals occur less frequently, and the long-time-averaged flow fields exhibit greater symmetry, with clearly organized LMPs and HMPs.

Figure 14. (a) Percentage occurrence of updraught and downdraught events and (b) average persistence time for the cases outlined in tables 1 and 3. The red colour denotes a net updraught at the domain centre, while the blue colour denotes a net downdraught. Circles indicate cases exhibiting a dominant polarity in long-time-average visualization, whereas triangles indicate cases with disrupted secondary flows in long-time-average visualization. The cases are arranged in decreasing order of $s_a/H$ .

Figure 14 provides a quantitative assessment of the visual observations discussed above. This figure includes all cases listed in tables 1 and 3, arranged in decreasing order of $s_a/H$ . The percentage existence of updraught and downdraught events refers to the fraction of the simulation during which a particular event is dominant. The average time of existence denotes the mean duration for which a given event persists continuously. To avoid the influence of very short-lived events on the average time of existence, only the events persisting for $t^*\gt 1$ are considered. As shown in figure 14(a), cases in which one event exists for more than 2.5 times the duration of the other tend to exhibit a dominant polarity in the long-time averaged flow field, accompanied by clear delta-scale vortices. In contrast, intermediate cases show roughly equal presence of both polarities, with no single event significantly dominating the other. A similar trend is observed in the average time of existence. For the intermediate cases, both events persist for nearly equal durations, whereas in cases with a dominant polarity, the dominant event lasts at least twice as long as the non-dominant one.

In conclusion, the analysis presented in this subsection reveals a clear transitional behaviour as the parameter $s_a/H$ decreases from case S2-5-2 to case S3.5-7-2. For the cases outlined in tables 1 and 3, it is observed that when $s_a/H \gtrsim 1.2$ , delta-scale secondary flows emerge with LMPs and associated updraughts aligned at the domain centre. These updraughts persist for the majority of the simulation time, and although reversals can occur, they are brief and contribute negligibly to the long-time averaged flow field. When the cases lie in the regime $1.2 \gtrsim s_a/H \gtrsim 1$ , the updraughts at the domain centre become increasingly unstable and do not sustain for longer periods during the simulations. As a result, the reversal occurs much more frequently, and both events appear to be short-lived. As no particular event significantly dominates the other in this regime, the secondary flows appear to be disrupted in the long-time average visualization. However, for these cases, delta-scale vortices can still be observed in the conditionally averaged flow field, as shown in figure 12(c,d). When $s_a/H$ is reduced further such that the cases lie in the regime $s_a/H \lesssim 1$ , the HMPs and associated downdraughts at the domain centre begin to dominate and sustain for longer durations. As a result, the propensity of the system for reversals decreases, and we again observe delta scale secondary flows with a prominent alignment of pathways in the long-time average visualization.

4.3. Characterization of instantaneous vertical velocity distribution

In the previous subsection, we classified events based on the sign of the vertical velocity at grid points located at the domain centre, thereby distinguishing between updraughts and downdraughts. While this approach provided valuable insight into the switching behaviour between these events at the domain centre, it did not capture variations in the magnitude of vertical velocity within individual events. As observed in figure 13, different roughness configurations influence the stability of the events sustained at the domain centre. This raises the question of whether the time series of vertical velocity exhibits distinct signatures depending on whether a dominant event – either updraught or downdraught – is sustained at the domain centre, or whether the flow frequently switches between the two. Therefore, in this subsection, we analyse the temporal evolution of vertical velocity to gain a deeper understanding of its dynamics during updraught and downdraught events.

Figure 15. (a) Probability density functions of the time series of vertical velocity fluctuations for the cases listed in table 1. Data are sampled at two streamwise positions: the location of the first row of elements (solid lines), and the midpoint between the first and second element rows (circles). In both cases, the spanwise and vertical location of the sampling point corresponds to the centre of the Y–Z plane, i.e. $y = L_y/2$ and $z = H/2$ . (b) Time series of vertical velocity fluctuations for case S2-5-2, sampled at a grid point aligned with the first element row. (c) Low-pass filtered version of the signal shown in (b). (d) Signal sampled and averaged over every $x$ location at $y = L_y/2$ and $z = H/2$ .

Figure 15(a) presents the PDFs of the time series of vertical velocity fluctuations for the cases listed in table 1. To assess the influence of streamwise heterogeneity introduced by the roughness arrangement, vertical velocity time signals are sampled at two streamwise locations: one coinciding with the first row of elements, and the other situated at the midpoint between the first and second element rows. The spanwise position of the sampling point is fixed at the domain centre (i.e. $y = L_y/2$ ), and the vertical position is located at half the half-channel height (i.e. $z = H/2$ ). As seen from the figure, the PDFs at the two streamwise locations exhibit a high degree of similarity, suggesting that the streamwise heterogeneity in the roughness layout does not significantly affect the statistical properties of vertical velocity at this particular spanwise location and height.

In this figure, the case S2-5-2 exhibits a distinctly broader PDF – indicating a higher standard deviation – while the PDFs for the other cases appear remarkably similar. The similarity between the intermediate cases (S2.75-5-2 and S3.5-5-2) and S4-5-2 is intriguing, as it suggests that the frequent pathway reversals characterizing the intermediate cases do not leave a distinct imprint on the temporal statistics of vertical velocity, compared with cases with a persistent central pathway and infrequent reversals. Instead, the broader distribution observed for S2-5-2, along with the monotonic (albeit modest) decrease in standard deviation across the other cases, implies that the intensity of fluctuations is primarily governed by how long the flow locally maintains an updraught—typically associated with stronger turbulent activity – or a downdraught, which tends to exhibit weaker fluctuations.

Figure 15(b) shows the time series of vertical velocity fluctuations for case S2-5-2. While the signal exhibits strong turbulent fluctuations, it also reveals an underlying low-frequency modulation. To isolate and examine this low-frequency behaviour, we apply the maximal overlap discrete wavelet transform using a Coiflet2 wavelet. The maximal overlap discrete wavelet transform smooth component acts as a zero-phase, shift-invariant low-pass filter, preserving the temporal alignment of the low-frequency content with the original signal (Percival & Mofjeld Reference Percival and Mofjeld1997). The support of the scaling function used for filtering spans approximately 45 000 time steps ( $t^* = tu_\tau /h \approx 9.1$ ), applied to a signal with a total duration of 7 000 000 time steps ( $t^* = 1417.5$ ). The resulting low-pass filtered signal is shown in figure 15(c). Figure 15(d) presents the streamwise-averaged vertical velocity fluctuations, obtained by averaging the signal across all the streamwise grid points at $y = L_y/2$ and $z = H/2$ at each time step. The resulting time series closely resembles the low-pass filtered signal derived from a single point, indicating that the streamwise averaging operation at this spanwise and vertical location effectively preserves the long-time variations in the vertical velocity fluctuations observed during the simulation. Therefore, to mitigate sensitivity to sampling location, we hereafter focus on the streamwise-averaged vertical velocity signal.

Figure 16. Time series of vertical velocity fluctuations that are sampled and averaged over every $x$ location at $y = L_y/2$ and $z = H/2$ , for cases: (a) S2-5-2, (b) S2.75-5-2, (c) S3.5-5-2, (d) S4-5-2. The red colour denotes a positive vertical velocity $(\langle \overline {w} \rangle _x + \langle w^\prime \rangle _x \gt 0)$ while the blue colour denotes a negative vertical velocity.

Figure 16 shows the streamwise-averaged vertical velocity fluctuations at $y = L_y/2$ and $z = H/2$ for the cases listed in table 1. To enable direct comparison with figure 13, which illustrates the time-dependent switching of updraught and downdraught events, the signal is colour-coded such that red indicates regions where the vertical velocity is positive (i.e. $\langle \overline {w} \rangle _x + \langle w' \rangle _x \gt 0$ ), while blue indicates regions where it is negative. From the figure, it is evident that all four cases exhibit similar low-frequency, periodic modulations in the vertical velocity, with the standard deviation remaining within 10 % across all cases. Although conditional averaging previously indicated that the updraught at the domain centre is relatively stable in case S2-5-2, persisting for long durations with infrequent reversals, this stability of the event does not necessarily correspond to a vertical velocity signal with lower temporal variability. As seen in figure 16(a), the vertical velocity remains predominantly positive for extended periods, but its magnitude continues to undergo quasiperiodic oscillations, similar to those observed in the other cases (figure 16 b–d). However, because these oscillations mostly occur within the positive regime, they do not manifest as reversals in the conditionally averaged statistics. In contrast, for the intermediate cases (S2.75-5-2 and S3.5-5-2), the vertical velocity oscillations frequently cross zero, leading to repeated sign changes. This is primarily due to the mean vertical velocity at this location being close to zero. Thus, the more frequent switching observed in figure 13(c), relative to figure 13(a), reflects the near-neutral mean state rather than an inherently more chaotic fluctuation pattern. In fact, as noted earlier in the discussion of figure 15, updraught-dominated regions are associated with stronger turbulent fluctuations. Therefore, case S2-5-2, which sustains more persistent updraughts, exhibits higher overall variance in vertical velocity compared with the intermediate cases, despite appearing less dynamic in terms of switching between updraughts and downdraughts.

Figure 17. Probablity density function of time series of vertical velocity fluctuations that are sampled and averaged over every $x$ location at $y = L_y/2$ and $z = H/2$ , for cases: (a) S2-5-2, (b) S2.75-5-2, (c) S3.5-5-2, (d) S4-5-2. The red region denotes a positive vertical velocity $(\langle \overline {w} \rangle _x + \langle w^\prime \rangle _x \gt 0)$ while the blue region denotes a negative vertical velocity.

Figure 17 shows the PDFs of the streamwise-averaged vertical velocity fluctuations discussed above. The PDFs are colour-coded such that red corresponds to regions where the vertical velocity is positive (i.e. $\langle \overline {w} \rangle _x + \langle w^\prime \rangle _x \gt 0$ ), and blue corresponds to negative vertical velocity. In contrast to the unimodal distributions observed in figure 15, the streamwise-averaged PDFs in this figure reveal a clear multimodal character across all cases. For cases S2-5-2 and S4-5-2, which exhibit a dominant polarity in the long-time-averaged flow field, the prominent modes lie predominantly on one side of the mean: positive for S2-5-2 and negative for S4-5-2. However, for the intermediate cases (S2.75-5-2 and S3.5-5-2), the modes are distributed on both sides of the mean, reflecting frequent transitions between updraught- and downdraught-dominated states. These observations imply that secondary circulations are inherently unsteady, regardless of whether the time-averaged flow exhibits a persistent polarity or not. In all cases, the vertical velocity signal fluctuates between multiple preferred states, as evidenced by the modal structure of the PDFs. This behaviour persists even in flows that do not visibly exhibit frequent switching, indicating that temporal variability in the vertical momentum transport is a fundamental feature of the secondary flow dynamics.