Boundary layer transition due to distributed roughness: effect of roughness spacing

Abstract The influence of roughness spacing on boundary layer transition over distributed roughness elements is studied using direct numerical simulation and global stability analysis, and compared with isolated roughness elements at the same Reynolds number $Re_h=U_eh/\nu$ ($U_e$ is the boundary layer edge velocity, h is roughness height and $\nu$ is the kinetic viscosity of the fluid). Small spanwise spacing ($\lambda _z=2.5h$) inhibits the formation of counter-rotating vortex pairs (CVP) and, as a result, hairpin vortices are not generated and the downstream shear layer is steady. For $\lambda _z=5h$, the CVP and hairpin vortices are induced by the first row of roughness, perturbing the downstream shear layer and causing transition. The temporal periodicity of the primary hairpin vortices appears to be independent of the streamwise spacing. Distributed roughness leads to a lower critical roughness Reynolds number for instability to occur and a more significant breakdown of the boundary layer compared with isolated roughness. When the streamwise spacing is $\lambda _x=5h$, the high-momentum fluid barely moves downward into the cavities and the wake flow has little impact on the following roughness elements. The leading unstable varicose mode is associated with the central low-speed streaks along the aligned roughness elements, and its frequency is close to the shedding frequency of hairpin vortices in the isolated case. For larger streamwise spacing ($\lambda _x=10h$), two distinct modes are obtained from global stability analysis. The first mode shows varicose symmetry, corresponding to the primary hairpin vortex shedding induced by the first-row roughness. The high-speed streaks formed in the longitudinal grooves are also found to be unstable and interact with the varicose mode. The second mode is a sinuous mode with lower frequency, induced as the wake flow of the first-row roughness runs into the second row. It extracts most energy from the spanwise shear between the high- and low-speed streaks.


Introduction
When laminar boundary layers encounter rough surfaces, the flow field can be greatly modified by surface roughness, and transition to turbulence can occur.Understanding roughness-induced transition is important since it leads to increase in friction drag and affects the performance of aeronautical and naval applications.Three-dimensional (3-D) surface roughness can be generally categorized into isolated and distributed roughness elements, both of which are commonly involved in engineering applications.While an isolated roughness element represents a single protuberance or a trip on the surface, it can also be considered as a foundational model to be extended to distributed surface roughness.
The effects of isolated roughness elements on boundary layers have been studied in past literature (Baker 1979).The streamwise vortices induced by an isolated roughness element create longitudinal streaks downstream whose unstable nature plays a crucial role in the transition process (Reshotko 2001;Fransson et al. 2004Fransson et al. , 2005)).The linear stability properties of isolated-roughness-induced transition have been investigated both computationally and experimentally (De Tullio et al. 2013;Loiseau et al. 2014;Citro et al. 2015;Bucci et al. 2021;Ma & Mahesh 2022).They are found to depend on the combined effects of parameters such as the ratio of roughness height h to the displacement boundary layer thickness δ * , aspect ratio η = d/h (d is roughness width) and Re h = U e h/ν (U e is the boundary layer edge velocity and ν is the kinetic viscosity of the fluid).Compared to the isolated roughness element, distributed roughness elements display phenomena not present in the isolated case.Fewer studies have been performed on how 3-D distributed roughness elements affect the stability properties and flow structures in transitional boundary layers.
Global linear stability analysis (Theofilis 2011) provides insight into the temporal disturbance growth during the early stages of transition, and is useful for non-parallel flows such as roughness wakes.It has been used to capture and understand the leading unstable modes in boundary layer transition due to isolated roughness (Loiseau et al. 2014;Citro et al. 2015;Bucci et al. 2021;Ma & Mahesh 2022).The potential of global stability analysis to detect the leading unstable modes induced by distributed surface roughness remains relatively unexplored.In this work, we therefore combine DNS and global stability analysis to investigate transition due to distributed surface roughness with varying streamwise and spanwise spacing.
Past studies of transition over distributed surface roughness have mainly focused on the effects of roughness height.Corke et al. (1986) studied the effects of distributed roughness on transition and suggested that transition is most likely to be triggered by the few highest peaks.They also found that the low-inertia fluid in the valleys between roughness elements could be influenced more easily by free-stream disturbances.For roughness with small amplitudes, transition is induced through a linear amplification of temporal disturbance growth followed by secondary instabilities and breakdown to turbulence (Reshotko 2001).In contrast, large-amplitude roughness creates local separations, leading to strong 3-D disturbances that can develop into turbulence directly.This is the socalled "bypass" transition which means that the linear instability processes such as Tollmien-Schlichting waves are bypassed.Vadlamani et al. (2018) conducted numerical investigations on transition of a subsonic boundary layer over sinusoidal roughness elements with different height.They suggested that secondary instabilities on the streaks promote transition to turbulence for roughness elements that are inside the boundary layer, and that the instability wavelengths appear to be governed by the fixed streamwise and spanwise spacings between roughness elements.For roughness elements that are higher than the boundary layer, they found that the scale of instability is shorter and the shedding from the obstacles leads to transition.von Deyn et al. (2020) investigated the influence of distributed roughness and freestream turbulence on bypass transition.They found that the bypass transition occurs as a result of secondary instability of low-speed streaks and suggested that the streak spacing does not vary with different roughness density at constant roughness height.
While roughness height is an important parameter to affect transition, roughness spacing can also potentially modify the transitional flow behavior and associated in-stability mechanisms.Investigation of spanwise proximity of roughness elements has commonly involved spanwise arrays of roughness elements.Ergin & White (2006) studied transitional flow behind a spanwise array of cylindrical roughness elements and suggested that transition results from a competition between unsteady disturbance growth and steady disturbance decay.Choudhari & Fischer (2006) performed numerical simulations of a flat plate boundary layer past a spanwise array of cylindrical roughness elements (spaced three diameters away from each other).They observed self-sustained and unsteady vortical structures that resemble the flow structures during natural transition for roughness with large height.Muppidi & Mahesh (2012) investigated supersonic boundary layer flow over distributed surface roughness and showed that the counter-rotating vortex pairs induced by roughness perturb the shear layer and triggers transition.With closely packed roughness elements, there is little upstream spacing to generate a strong horseshoe vortex and little spanwise spacing to produce counter-rotating vortices, which can prevent transition from occurring.The influence of spanwise spacing on instability and transition has been less-explored in past work.Our present study aims to quantitatively evaluate the effects of spanwise spacing of roughness elements on the vortical structures and associated instability characteristics in transition.
The role of streamwise proximity of roughness elements is also a crucial factor of transition due to distributed roughness.Choudhari et al. (2010) performed simulations of a high-speed boundary layer past an isolated diamond trip as well as two trips aligned with the flow direction (spaced three times the trip width).They found that the introduction of an additional trip could amplify the streak amplitude and cause transition even at a smaller trip height, while the case with larger trip height presents a weaker augmentation of the streak amplitude since the flow has not recovered as much as for smaller trip height.In a turbulent boundary layer, the classification of d-type and k-type behaviors (Jiménez 2004) is related to the streamwise roughness spacing and might be made based on the ratio of roughness pitch to height being equal to 3 (Tani 1987).Perry et al. (1969) and Ikeda & Durbin (2007) further claimed that the difference between k-type and d-type roughness is related to the state of vortex shedding: for d-type roughness, stable vortices form within the grooves and no eddy sheds into the flow above the crests; for k-type roughness, eddies with length scale of order h shed into the flow above the crests.Vanderwel & Ganapathisubramani (2015) suggested that distributed roughness with streamwise gaps less than 4-5h would act like continuous strips in turbulent boundary layers, while more than 5h would act like 3-D distributed roughness.Although the dependence of flow behaviors on different streamwise spacing has been discussed extensively in turbulent boundary layers, less is studied on how streamwise spacing affects transition.This work therefore investigates transition due to distributed surface roughness with streamwise spacing 5h and 10h, and compares the results to the isolated roughness case.
The paper is organized as follows.The numerical methodology is introduced in §2.The flow configuration and simulation details are discussed in §3.In §4, the influence of distributed roughness on the transitional boundary layer is studied, the formation mechanisms of hairpin vortices in different roughness distributions are examined, and the base flow and global stability results are analyzed and compared to the isolated case.Finally, the paper is summarized in §5.

Numerical methodology
The governing equations and numerical method are briefly summarized.An overview of linear stability and details regarding the iterative eigenvalue solver are provided.

Direct numerical simulation
The incompressible Navier-Stokes (N-S) equations are solved using the finite volume algorithm developed by Mahesh et al. (2004): where u i and x i are the i-th component of the velocity and position vectors respectively, p denotes pressure divided by density, ν is the kinematic viscosity of the fluid and K i is a constant pressure gradient (divided by density).Note that the density is absorbed in the pressure and K i .The algorithm is robust and emphasizes discrete kinetic energy conservation in the inviscid limit which enables it to simulate high-Re flows without adding numerical dissipation.A predictor-corrector methodology is used where the velocities are first predicted using the momentum equation, and then corrected using the pressure gradient obtained from the Poisson equation yielded by the continuity equation.The Poisson equation is solved using a multigrid pre-conditioned conjugate gradient method (CGM) using the Trilinos libraries (Sandia National Labs).
The DNS solver has been validated for a variety of problems on related wall-bounded flows, including: realistically rough superhydrophobic surfaces (Alamé & Mahesh 2019), random rough surfaces in turbulent channel flow (Ma et al. 2021) and boundary layer transition with an isolated roughness element (Ma & Mahesh 2022).

Linear stability analysis
Modal linear stability analysis is the study of the dynamic response of a base state subject to infinitesimal perturbations (Theofilis 2011).In the present work, the incompressible Navier-Stokes equations (2.1) are linearized about a base state, u i and p.The flow can be decomposed into a base state subject to a small O( ) perturbation ũi and p.The linearized Navier-Stokes (LNS) equations are obtained by subtracting the base state from equations (2.1) and can be written as follows: The same numerical schemes as the N-S equations are used to solve the LNS equations.
The LNS equations are subject to the following boundary conditions: where S is the boundary of the spatial domain.
The LNS equations can be rewritten as a system of linear equations, where A is the LNS operator and ũi is the velocity perturbation.The solutions to the linear system of equations (2.4) are: where ûi is the velocity coefficient, and both ûi and ω can be complex.This defines σ = Re(ω) as the growth/damping rate and ω a = Im(ω) as the temporal frequency of ûi .The linear system of equations can then be transformed into a linear eigenvalue problem: Table 1.Simulation parameters for the isolated and distributed roughness cases at Re h = 600.
where ω j = diag(Ω) j is the j-th eigenvalue and ûj i = U i [j, :] is the j-th eigenvector.For global stability analysis, solving the eigenvalue problem using direct methods is very computationally expensive.Instead, a matrix-free method, the implicitly restarted Arnoldi method (IRAM) implemented in the PARPACK library is used to efficiently solve for the leading eigenvalues and eigenmodes in the present work.The ratio of the first-row roughness height to the displacement boundary layer thickness h/δ * is 2.86, consistent with the isolated cube in Ma & Mahesh (2022).The roughness height is h = 1, the reference length in the simulations.The Blasius laminar boundary layer solution is specified at the inflow boundary, and convective boundary conditions are used at the outflow boundary.Periodic boundary conditions are used in the spanwise direction.No-slip boundary conditions are imposed on the flat plate and the roughness surfaces.The boundary conditions U e = 1, ∂v/∂y = 0 and ∂w/∂y = 0 are used at the

Downstream flow separation
A laminar boundary layer undergoes flow separation downstream of roughness elements.The wake flow downstream of roughness elements is visualized by instantaneous streamwise velocity in the symmetry plane in figure 2. Figure 2(a) shows that for Case (5h, 2.5h), the reverse flow is strong in the first cavity and becomes weak from the second cavity onwards.Wall-normal momentum transfer hardly occurs and the boundary layer above the roughness layer remains laminar.For Case (5h, 5h), the unsteadiness of the reverse flow is observed in the first cavity in figure 2(b).The length scale of flow separation is comparable to the streamwise spacing λ x = 5h.The high-momentum fluid above roughness tips hardly penetrates into the cavities.For a larger streamwise spacing λ x = 10h, as shown in figure 2(c), the high-momentum fluid transfers downwards into the first cavity, and impinges on the second-row roughness, which could possibly induce different instability modes.The penetration of high-momentum fluid into the cavities becomes weaker with increasing downstream distance.

Boundary layer integral parameters
The effect of streamwise roughness spacing on transitional boundary layers is examined through boundary layer integral parameters for Cases (5h, 5h) and (10h, 5h), with a comparison to the isolated cube.The streamwise variation of the displacement boundary layer thickness (δ * ), momentum boundary layer thickness (θ) and shape factor (H) is presented in figure 3(a)-(c).The parameters are computed from a time-averaged flow field and spanwise averaging is performed across the span.Comprehensive spatial averaging where the averaging volume includes both the fluid and solid area is used to ensure continuity in the profiles (Xie & Fuka 2018).To understand the spatial inhomogeneity of the flow field caused by distributed roughness, the local integral boundary layer parameters are also shown in figure 4.
In figure 3(a), an increase in δ * is seen due to the presence of roughness elements compared to the laminar Blasius solution.A more significant increase is observed for distributed roughness than the isolated roughness.The increase is more pronounced for distributed roughness with smaller streamwise spacing.Figure 4(a) shows that the region with lower values of local δ * corresponds to the location of lateral high-speed streaks induced by roughness.Figure 3(b) shows negative deviation from the Blasius solution in the profiles of θ for the distributed roughness cases.This is related to the strong reverse flow that occurs closely behind the first two rows of roughness elements, as demonstrated by figure 4(b).A steeper increase in the profiles of θ indicates that the breakdown of boundary layers is more significant for distributed surface roughness compared to the isolated roughness.
The shape factor in figure 3(c) shows significant increase compared to the isolated roughness.The high values of H are associated with inflection points and reveals locally the instability of the streaks induced by roughness elements, as shown in figure 4(c).The steep drop in H begins at approximately x = 5h for the three cases, suggesting that the streaks with high amplitudes break down and transition happens at similar downstream positions for three different roughness distributions.The shape factor converges farther downstream after the breakdown.The values are higher for Case (5h, 5h) when compared to Case (10h, 5h) and the isolated case, resulting from a stronger blockage effect of a denser roughness distribution.

Formation of hairpin-shaped vortices
Packets of hairpin-shaped vortices are key structures in roughness-induced transition.They are associated with global instability as known for isolated roughness.Cohen et al. (2014) proposed a model consisting of three key elements required for the formation of hairpin vortices: simple shear, counter-rotating vortices and two-dimensional vortex sheet, and highlights that this combination is unstable.The influence of roughness spacing on the generation of hairpin vortices is therefore important to understand and is investigated in this section for distributed surface roughness.

Counter-rotating vortex pairs (CVP)
As known for isolated roughness, the spanwise vortices formed upstream wrap around the roughness element and give birth to the streamwise vortices downstream of the roughness, and the streamwise counter-rotating vortices are known to play a critical role in the generation of hairpin vortices (Iyer & Mahesh 2013;Bucci et al. 2021).As the baseline, figure 5 depicts the characteristics of streamwise vortices induced by an isolated cube.The symmetry plane vortices (SP) located very close to the roughness have an upwash between them.They lift up and move towards each other with increasing downstream distance, and are dissipated at x = 25h.The off-symmetry plane vortices (OSP) located farther from the symmetry plane are the continuation of the vortex tubes from upstream.They have a central downwash which keeps them away from each other with increasing downstream distance.
The effects of spanwise spacing on CVP are investigated for Case (5h, 2.5h) in figure 6.At x = 0, the OSP vortices are observed in the groove between two adjacent cubes, similar as the isolated case.However, the SP vortices do not appear on either side of obstacles, in contrast to the isolated case.The generation mechanism of the SP vortices is examined using the mean streamlines in figure 6(a) and illustrated in figure 6(b).The secondary flow close to the cube sides moves downward due to the motion of the OSP vortices (from a to b), then moves toward the cube due to a positive spanwise pressure gradient (from b to c), and moves upward for mass conservation (from c to d).With small spanwise spacing, the OSP vortices in the groove are located closer toward each other, which strengthens the upward fluid motion in the groove and weakens the centrifugal forces.The last step under the effects of centrifugal forces for the generation of SP vortices (from d to a) is therefore missing.As a result, a weak CVP is observed at the roughness tip location at x = 4h, and is dissipated within a short downstream distance.The OSP vortices in the groove are also diminished with increasing downstream distance and mostly disappears beyond the second row of roughness elements, as shown in figure 6(c).The effects of streamwise spacing on CVP are investigated in figure 7 for Cases (5h, 5h) and (10h, 5h).For both cases, the SP and OSP vortices are observed and behave similarly as the isolated roughness at x = 0. Figure 7(a) shows that for Case (5h, 5h), the CVP grows with increasing downstream distance, and is distorted and pushed away from each other by the following roughness as observed at x = 4h.At x = 10h, the SP vortices developed from the front obstacles are weakened due to the stagnation effects of the following obstacles, but a new pair of SP vortices is induced again on either side of the cube, and the OSP vortices on the lateral sides are strengthened.Figure 7(c) shows that the CVP is dissipated beyond x = 30h, similarly as the isolated roughness case.
For Case (10h, 5h), the streamwise spacing is much larger than the streamwise length scale of the separation region.The behavior of CVP is similar as that for an isolated roughness at x = 0 and 4h, as shown in figure 7(b).Instead of being distorted by the roughness as Case (5h, 5h), both the SP and OSP vortices move closer towards the symmetry plane, enhancing the momentum exchange within the cavities.At x = 10h, they impinge onto the second row of roughness elements, and break down into small vortical structures.A new pair of SP vortices is generated and the OSP vortices are strengthened on the lateral sides.Figure 7(d) shows that after the second-row cubes, the SP vortices maintain in the cavities while the OSP vortices develop in the longitudinal grooves.The CVP persists farther than Case (5h, 5h), and is dissipated beyond x = 40h.

Perturbation to the shear layer
The CVP examined in §4.2.1 perturbs the shear layer in their vicinity, and their size and strength determine the nature of perturbation to the shear layer.The perturbation to the shear layer is visualized by instantaneous spanwise vorticity at the symmetry plane in figure 8.The shear layer downstream of the isolated cube is shown as the baseline in figure 8(a).The breakdown of the shear layer is observed downstream of the cube due to the perturbation of vortex shedding, and the vortex heads are dissipated at x = 20h.
For distributed surface roughness, figure 8(b) shows that for small spanwise spacing λ z = 2.5h, the downstream shear layer formed right above the roughness tips remains steady and no hairpin vortex shedding is produced in the wake flow.The inflection point in the mean velocity profiles is examined in figure 9.It is a necessary condition for instability in shear flows, and is corresponding to the peak location of spanwise vorticity.Although the inflection points formed by wall-normal shear can be identified for Case (5h, 2.5h) in figure 9, the absence of CVP due to the small spanwise spacing results in a failure of hairpin vortex generation.
With larger spanwise spacing λ z = 5h, the breakdown of the shear layer is seen in figures 8(c) and 8(d).Although the localized shear layers are induced by each roughness tip in Cases (5h, 5h) and (10h, 5h), the primary vortex shedding behaves similarly as the isolated case, and the vortex heads are dissipated at around x = 20h.The length scale of localized shear layers is equivalent to 5h, and few vortices penetrate into cavities for Case (5h, 5h).In contrast, for Case (10h, 5h), more complex vortical structures are evolved into the cavities from the second-row cubes and the strongest localized mean shear is demonstrated at the roughness location in figure 9.The periodicity of the primary vortex shedding seems to be similar as that for the isolated roughness and independent with different streamwise roughness spacing.The associated instability mechanisms will be further examined and discussed in §4.3.
The evolution of vortical structures is examined using the Q criterion (Hunt et al. 1988) in figure 10 for Case (5h, 5h).Both the SP and OSP vortices are observed in the vicinity of the first-row cubes.They interact with the shear layer, leading to the hairpin vortex shedding downstream of the first-row roughness elements.The packets of hairpintype structure with small legs labelled 1 can be identified in the top-left inset of figure 10.The shorter streamwise extent of the vortical motions is due to the blockage effects of the closely distributed cubes.As the vortex legs are inclined upward and undergo stretching by a positive wall-normal velocity gradient, they are cut off by the following cubes and break down into smaller vortical structures with low momentum within the roughness layer.The bottom-right inset of figure 10 shows that the interactions between vortical structures in the longitudinal grooves occur at approximately x = 15h.As the hairpin-type vortices break down into smaller structures at x = 20h, the small scales are still organized in spanwise coherent structures, labelled 2. The arch-shaped spanwise coherent structures are observed farther downstream.
For Case (10h, 5h), both SP and OSP vortices are observed, and the hairpin vortices are generated from the first-row cubes, as shown in figure 11.The top-left inset shows that the hairpin vortices labelled 1 behave similarly as the isolated roughness case.As the high-momentum fluid impinges onto the second row of roughness, another spanwise vortex wraps around the second-row cubes, and the vortical structures break down into small scales downstream of the second-row cubes.The CVP sheds from the second-row roughness tips, and a separation of the vortex heads labelled 2 occurs.This indicates that a different unstable mode might be induced by the second row of cubes.The bottomright inset shows that the vortical structures in the longitudinal grooves interact with each other from x = 15h.The spanwise coherent structures are less prominent than Case (5h, 5h).The longitudinal vortical structures labelled 3 are possibly related to the streak interactions in the grooves.

Global stability analysis
The hairpin vortices induced by roughness elements are inherently unstable as the CVP perturbs the shear layer, forming the inflection point in the base-flow velocity profiles.The perturbation of the shear layer is inhomogeneous in the spanwise direction and therefore the global stability analysis is useful to provide insights into the instability mechanisms associated with the distributed roughness wakes.

Base flow computation
The base flow for global instability analysis is computed for the distributed roughness cases.For Cases (5h, 2.5h) and (5h, 5h), the selective frequency damping (SFD) method ( Åkervik et al. 2006) is used to artificially settle the flow towards a steady equilibrium.An encapsulated formation of the SFD method developed by Jordi et al. (2014) is applied in the present work.
For Case (10h, 5h), it is found that the SFD method is unable to damp the oscillations due to the unsteady part of the solutions, even though careful choices are made for the control coefficient χ and the filter width ∆.The possible reason is that the SFD method is unable to get the steady state when there are multiple unstable modes, as indicated by Casacuberta et al. (2018).They discussed the effectivity of SFD for systems with more than one unstable eigenmode where the most unstable eigenvalue is µ c and other unstable eigenvalues are denoted by µ k .They concluded that SFD is unable to drive the system towards the base flow when µ k with large values of Im(µ k )/Re(µ k ) is present close to the origin.As discussed in Ma & Mahesh (2022), using the time-averaged mean flow as the base state for global stability analysis is able to capture the temporal frequency and associated mode shape of the primary vortical structures for roughnessinduced transition.The time-averaged mean flow is therefore considered as an alternate base flow in the present work to investigate global instability for Case (10h, 5h).

High-and low-speed longitudinal streaks
The high-and low-speed streaks are visualized in figure 12 by the streamwise deviation of the base flow from the Blasius solution for the three distributed roughness cases.Figure 12(a) shows that for Case (5h, 2.5h), the central and lateral low-speed streaks are merged with each other and form a shear layer above the roughness layer.The high-speed streaks are only induced at the first row of roughness elements and disappear within a short downstream distance due to the dense spanwise roughness distribution.In this case, the instability mechanism might be Kelvin-Helmholtz instability rather than streak instability.For Case (5h, 5h), both the central and lateral streaks are observed in figure 12 streaks develop away from the symmetry plane with increasing downstream distance, due to the presence of the following roughness elements.
As the streamwise spacing is increased to λ x = 10h, different behavior of high-and lowspeed streaks is seen in figure 12(c).In contrast to Case (5h, 5h), the central low-speed streak only forms within the vicinity downstream of the cubes.The high-speed streaks induced by the first-row cubes move close to each other and collide onto the following obstacles, which induces the high-speed streaks from the second-row cubes.The highspeed streaks grow and interact in the longitudinal grooves farther downstream.

Eigenspectra and eigenmodes
The leading eigenvalues for Cases (5h, 5h) and (10h, 5h) at different Re h are plotted in figure 13.For Case (5h, 5h), compared to the isolated case at the same Re h , the growth rate is larger and the temporal frequency is lower.This indicates that the distributed roughness elements lead to lower critical Re h for linear instability to occur compared to the isolated roughness.The associated eigenmode for Case (5h, 5h) is examined in figure 14.The varicose mode is observed along the central low-speed streak, similar to the varicose mode observed in the isolated case (Ma & Mahesh 2022).The dominant production terms of disturbance kinetic energy P y and P z are examined in figure 15.It shows that the unstable mode extracts energy from both the wall-normal and spanwise shear of the central low-speed streak as well as the localized shear layer induced by the cubes.
For Case (10h, 5h), two leading eigenvalues are captured in figure 13.One is the leading eigenvalue whose temporal frequency is close to that of the isolated case and Case (5h, 5h).This leading eigenvalue corresponds to the primary hairpin vortex shedding and is marginally stable, consistent with the state of marginal stability of mean flow for the isolated roughness suggested by Ma & Mahesh (2022).The corresponding eigenmode in figure 16 is associated with the high-speed streaks in the longitudinal grooves as well as the entire shear layer formed above the roughness tips.The production results in figure 17 indicate that the mode extracts most of energy from the localized shear caused by obstacles and the high-speed streaks in the grooves farther downstream.In contrast to the isolated case and Case (5h, 5h), the mode hardly extracts the energy from the central low-speed streak since the central streak is diminished for larger streamwise spacing.The other unstable leading eigenvalue with a lower frequency is also obtained for Case (10h, 5h) in figure 13. Figure 18 shows that the associated eigenmode displays sinuous symmetry, it is induced by the second row of roughness elements and fades away within short downstream distance.This indicates that when the streamwise spacing is larger than the region of flow separation, an additional unstable mode is generated as the wake flow from the first-row roughness impinges on the following roughness.The production results in figure 19 demonstrate that this sinuous mode mainly extracts its energy from the spanwise shear induced between the high-and low-speed streaks in the grooves.

Summary
The effects of roughness spacing on boundary layer transition due to distributed surface roughness are investigated.Both streamwise and spanwise proximities of roughness elements are considered.The transitional flow behavior and the primary vortical structures due to distributed surface roughness are examined by performing direct numerical simulations.Global stability analysis is performed to study the stability properties of the flow field modified by the distributed roughness, with a comparison to an isolated roughness element with the same geometry (Ma & Mahesh 2022).
When the spanwise spacing is small (λ z = 2.5h), the OSP vortices in the groove enhance the upward fluid motion and weaken the centrifugal forces, inhibiting the formation of SP vortices.With the absence of CVP, the hairpin vortices are not generated and the downstream shear layer remains steady at Re h = 600.The flow might be unstable at higher Re h due to Kelvin-Helmholtz instability.When the spanwise spacing increases to λ z = 5h, the wake flow becomes unsteady and the effects of streamwise spacing on transition are investigated for the fixed λ z = 5h.The steeper and larger increase in the streamwise variations of boundary layer thickness indicates that the breakdown of boundary layer is more significant for distributed surface roughness.The shape factor profiles suggest that transition occurs at similar downstream locations as the isolated roughness.
When the streamwise spacing is comparable to the length scale of flow separation (λ x = 5h), the high-momentum fluid above the roughness layer barely penetrates into the cavities, and the primary hairpin vortices with shorter legs shed at the similar frequency as the isolated roughness case.The global stability analysis indicates that the leading unstable eigenvalue is close to that of the isolated case, while the distributed roughness results in a lower critical Re h for instability to occur.The unstable eigenmode presents varicose symmetry, and extracts the energy from the central low-speed streaks and the localized shear induced by the cubes.
When the streamwise spacing increases to λ x = 10h, the high-momentum fluid transfers downward into the cavities.The CVP and hairpin vortices induced by the first row of roughness break down into small vortical structures when they run onto the second-row roughness.Two distinct eigenvalues are obtained from global stability analysis.One corresponds to the primary hairpin vortex shedding induced by the firstrow cubes, whose shedding periodicity is independent on different streamwise spacing.It is also associated with the high-speed streaks developed in the longitudinal grooves farther downstream.The other leading unstable eigenmode with lower frequency presents sinuous symmetry.It is induced as the wake flow of the first row of roughness impinges onto the second row of roughness, and extracts its energy from the spanwise shear between the high-and low-speed streaks.

Figure 1 .
Figure 1.Schematic of the flow configuration and roughness distribution.x, y and z coordinates are the streamwise, wall-normal and spanwise directions.The distance from the inlet of the computational domain to the first row of roughness elements remains constant as l = 15h.The streamwise and spanwise spacing are denoted by λx and λz respectively.

Figure 1
Figure 1 depicts the flow configuration and roughness distribution.At the inflow, a laminar Blasius boundary layer profile is prescribed.Cuboids with aspect ratio of width to height η = d/h = 1 are aligned in both the streamwise and spanwise directions.The ratio of the first-row roughness height to the displacement boundary layer thickness h/δ * is 2.86, consistent with the isolated cube inMa & Mahesh (2022).The roughness height is h = 1, the reference length in the simulations.The Blasius laminar boundary layer solution is specified at the inflow boundary, and convective boundary conditions are used at the outflow boundary.Periodic boundary conditions are used in the spanwise direction.No-slip boundary conditions are imposed on the flat plate and the roughness surfaces.The boundary conditions U e = 1, ∂v/∂y = 0 and ∂w/∂y = 0 are used at the

Figure 3 .
Figure 3. Streamwise variations of (a) displacement boundary layer thickness, (b) momentum boundary layer thickness and (c) shape factor of Cases (5h, 5h) and (10h, 5h) at Re h = 600, with comparison to the isolated roughness and laminar Blasius solution. (

Figure 5 .Figure 6 .
Figure 5. Mean streamwise vorticity ωx = ±0.5 with mean streamlines in the cross-flow planes (top) and the plane of y = 0.5h (bottom) from the DNS at Re h = 600 for the isolated roughness.

Figure 7 .
Figure 7. Mean streamwise vorticity ωx = ±0.5 with mean streamlines in the cross-flow planes and the plane of y = 0.5h from the DNS at Re h = 600 for (a, c) Case (5h, 5h) and (b, d) Case (10h, 5h).The contour levels are the same as figure 5.

Figure 10 .
Figure 10.Instantaneous vortical structures for Case (5h, 5h) in perspective and top-down views, visualized by isocontours of Q = 0.1U 2 e /h 2 and colored with instantaneous streamwise velocity.Plotted in the x-y and z-y planes are the contours of instantaneous spanwise vorticity, in the range of −1 (black) and 0 (white).

Figure 11 .
Figure 11.Instantaneous vortical structures for Case (10h, 5h) in perspective and top-down views, visualized by isocontours of Q = 0.1U 2 e /h 2 and colored with instantaneous streamwise velocity.The contour levels are the same as figure 10.

Figure 12 .
Figure 12.Top-down views of high-and low-speed streaks, visualized by isosurfaces (top) and contour plots at the plane y = 0.5h (bottom) of the streamwise velocity deviation of the base flow from the Blasius boundary layer solution, u d = U b − u bl , for (a) Case (5h, 2.5h), (b) Case (5h, 5h) and (c) Case (10h, 5h).The contour levels in (c) are the same as (a) and (b).

Figure 13 .Figure 14 .
Figure 13.Leading eigenvalues for Cases (5h, 5h) and (10h, 5h) at different Re h , with a comparison to the isolated roughness case.The vertical dash-dot line denotes the Strouhal number St = ωah/(2πu h ) of the primary hairpin vortices, where u h is the Blasius velocity at roughness tips.

Figure 17 .Figure 18 .
Figure 17.Contour plots of Py (left) and Pz (right) in the cross-flow planes at (a) x = 10h and (b) x = 20h for the leading varicose mode of Case (10h, 5h) at Re h = 600.The contour levels are the same as figure 15.

Figure 19 .
Figure 19.Contour plots of Py (left) and Pz (right) in the cross-flow planes at (a) x = 15h and (b) x = 20h for the leading sinuous mode of Case (10h, 5h) at Re h = 600.The contour levels are the same as figure 15.