Linear stability analysis of a boundary layer with rotating wall-normal cylindrical roughness elements

Abstract Cylinder stubs of finite length, mounted in a wall-normal manner on a flat plate, exert considerable control on their wake when rotating around their central axis. The present paper investigates these effects with linear stability theory and a direct numerical simulation. Three configurations are considered: evenly spaced corotating roughness elements, as well as positive and negative counter-rotating roughness pairs. Here, ‘positive’ and ‘negative’ are defined in accordance with the induced high- and low-speed streaks, respectively. A primary feature of the rotating-cylinder-induced wake is a ‘dominating inner vortex’ (DIV), which intensifies the lift-up effect and creates high-amplitude streaks. Linear stability analysis shows that the modified streaky flow is capable of effectively stabilizing Tollmien–Schlichting (TS) modes. The mechanism of TS mode stabilization, as found by a perturbation kinetic energy analysis, is attributed to the reduction of the wall-normal perturbation production. On the other hand, an inviscid inflectional instability mode related to the presence of the roughness appears which destabilizes the boundary-layer flow, primarily due to an increase in wall-normal perturbation energy production, but also due to increasing spanwise energy production, depending on the case. The inflectional instability roughness mode is more amplified with thicker cylinders since the induced DIV tends to support longer-living inflection points. Regarding an imaginable laminar–turbulent transition delay, positive rotating thinner roughness pairs would be preferable.


Introduction
Surface roughness elements of boundary layer size and smaller are ubiquitous and have a crucial influence on flow instability and laminar-turbulent transition. The velocity streaks generated by the so-called lift-up effect (Landahl 1980) had long been believed to trigger earlier transition, until the work of  and Fransson et al. (2006) demonstrated that the streaky flow produced by roughness arrays is also capable of delaying laminar-turbulent transition. This motivated additional investigations of transition delay methods with roughness elements in both two-dimensional (2-D) (Fransson et al. 2006) and three-dimensional (3-D) (Saric, Carpenter & Reed 2011) boundary layers, due to the fact that a successful transition delay technology has a great benefit in the reduction of operational costs for moving vehicles and engineering parts or processes.
The laminar-turbulent transition process in a Blasius boundary layer can take several paths depending on the level of external disturbances (Morkovin 1994), like free stream turbulence, surface roughness or controlled excitation. At low levels of excitation, linear stability theory (LST) reveals the contribution of Tollmien-Schlichting (TS) waves, which are 2-D, viscous and grow exponentially. As the amplitude of the TS-wave grows to 1 % of the free stream velocity, secondary instability sets in and it quickly breaks down into turbulence (Kachanov 1994). At higher external disturbance levels, very weak streamwise vortices are able to push high-speed fluid towards the wall and the low-speed fluid to the opposite direction, i.e. the lift-up mechanism (Landahl 1980), thus creating elongated streamwise streaks. In such circumstances, a transient algebraic increase of the perturbation, which is due to the non-normality of the controlling stability operator (Trefethen et al. 1993), followed by viscous decay can be expected (Reshotko 2001). If the external disturbance level is relatively low, the non-modal growth falls back to the modal growth of the primary instability. With moderate levels of external disturbances, the generated velocity streaks are able to support inviscid inflectional instability. A sinuous type of instability sets in when the streak amplitude exceeds 26 % of the free stream velocity u ∞ , while a varicose type is expected when the streak amplitude reaches 37 % of u ∞ (Andersson et al. 2001). The amplification of inviscid inflectional modes activates generation of higher harmonics, leading to the breakdown of the streaks and the rise of turbulent spots (Bakchinov et al. 1995;Brandt, Schlatter & Henningson 2004). If the boundary layer flow is exposed to even higher levels of external disturbances, the non-modal growth of disturbances can evolve directly into the turbulent state, bypassing the primary instability mechanism (Morkovin 1985).
Regarding surface-roughness-induced laminar-turbulent transition, the investigations date back to the 1950s, for example Tani & Sato (1956) for 2-D roughness elements and Gregory & Walker (1956) for 3-D isolated roughness elements. Two-dimensional roughness elements, such as gaps and steps, can be understood as an amplifier of 2-D TS waves (Ergin & White 2006), where recirculation zones behind the 2-D roughness elements provide a destabilization mechanism for the primary instability and thus enhance the growth rate of TS waves (Klebanoff & Tidstrom 1972;Goldstein 1985). Understanding of boundary layer instability with 3-D roughness elements is, however, much more complicated. Gregory & Walker (1956) were among the first to investigate a boundary layer with a 3-D isolated roughness element, where a horseshoe vortex that wraps around the roughness element and two following elongated counter-rotating streamwise vortex legs were visualized with smoke flow. Another distinct feature which is well known today is a system of downstream velocity streaks of alternating high-and low-speed fluid generated by the counter-rotating vortex legs. Because of the lift-up effect (Landahl 1980), this streaky flow can give rise to transient growth which can be strong enough to trigger laminar-turbulent transition (Joslin & Grosch 1995). In this regard, it is believed that such spanwise periodic steady streamwise vortices are the most dangerous disturbances according to optimal transient growth theory (Butler & Farrell 1992;Luchini 2000). Furthermore, an optimal streak as predicted by Andersson, Berggren & Henningson (1999) is one with a non-dimensioned spanwise wavenumber β = 0.45. However, good agreement has only been found in experiments with free stream, turbulence-generated streaks (Matsubara & Alfredsson 2001). The transient growth with surface-roughness-elements-generated streaks is only suboptimal (Ergin & White 2006;Rizzetta & Visbal 2007), i.e. it reaches its maximum far downstream of the optimal prediction. At high streak amplitude, secondary instability arises either from wall-normal or spanwise inflectional velocity profiles which finally leads to breakdown of the streaks and ends up in rapid laminar-turbulent transition (Bakchinov et al. 1995;Andersson et al. 2001;Brandt et al. 2004).
In contrast to the above-mentioned streak-induced non-modal instability, which is found to trigger earlier transition, a stabilization effect with steady and unsteady streaks were experimentally observed by Kachanov & Tararykin (1987) and Boiko et al. (1994), respectively. The stabilization mechanism as revealed by numerical simulation and linear stability analysis (Cossu & Brandt 2002 is that the extra perturbation production energy from the spanwise shear brought in by the streaks turns to become negative, which, together with the viscous dissipation, outweigh the wall-normal production term, thus resulting in an overall stabilization effect. Inspired by this theoretical prediction, cylindrical roughness-element-induced streaks are observed to delay transition experimentally (Fransson et al. 2006). This approach needs a careful set-up to avoid the fast growing inflectional secondary instability. It is also found that the stabilization effect prefers higher amplitude streaks which produce stronger spanwise shear and accordingly greater stabilization effects. However, the global instability arising from the wake behind big cylinders prevents the use of large cylindrical roughness elements (Loiseau et al. 2014). Therefore, winglet-type miniature-vortex generators (MVG) were used to generate high-amplitude streaks and smaller recirculation zones behind the MVG (Shahinfar et al. 2012;Siconolfi, Camarri & Fransson 2015), so as to avoid the global instability. Accordingly, the streak-amplitude threshold is pushed from 12 % of u ∞ with cylindrical roughness elements (Fransson et al. 2004) to 32 % with MVG .
The present paper considers a novel active method to control velocity streaks by simply rotating the roughness elements along their axis, which is oriented normal to the flat-plate surface. Existing investigations which can be found in literature deal with infinite static and rotating circular cylinders in a uniform free stream flow, see Williamson (1996) and Rao et al. (2015) for a comprehensive review. By rotating an infinite cylinder, the 2-D Bénard-von Kármán vortex shedding can be efficiently suppressed (Mittal & Kumar 2003). Despite the vast number of investigations of flow instability with rotary infinite cylinders, there is no study of boundary layer stability with rotating finite-length cylindrical elements apart from ours.
Before the relevant instability mechanisms were revealed, a correlation between the observed transitional roughness Reynolds number Re kk = u(k)k/ν and the roughness element aspect ratio η = D/k was developed (Klebanoff, Schubauer & Tidstrom 1955;Tani 1969). Here, u(k) is the unperturbed velocity at roughness height k and D the roughness diameter. The transition diagram published by Von Doenhoff & Braslow (1961) indicates that the critical transitional Reynolds number Re kk scales with the aspect ratio proportional to η 2/5 . At supercritical Re kk , the laminar-turbulent transition is expected shortly behind it. From this correlation diagram, the critical Reynolds number for a roughness element with aspect ratio η = 1 is 484 < Re kk < 882, which increases to 610 < Re kk < 1017 for a roughness element with η = 0.5. This means that thin roughness elements are, therefore, less likely to trigger an early transition. With the present method, higher amplitude streaks can be obtained by rotating thinner roughness elements, rather than employing thicker roughness elements. The benefits of this method are two-fold: TS-wave attenuation with higher amplitude streaks and prevention of early transition with thinner roughness elements.
In this work, the instability of the streaky flow induced by rotating roughness elements is studied with LST. In § 2, the set-up of the roughness elements studied and the relevant numerical methods are introduced. Next, the resulting baseflow and its linear stability properties are discussed in § 3. The underlying mechanisms are identified with a perturbation kinetic energy (PKE) analysis in § 3.3 which is followed by a direct numerical simulation (DNS) study in § 3.4 which excludes the threat of a possible bypass transition by the presence of the roughness elements. The results are summarized and concluded in § 4.

Set-up
The numerical set-ups of the boundary layer study with embedded corotating and counter-rotating cylindrical roughness elements are illustrated in figures 1(a) and 1(b), respectively. The roughness elements with height k and diameter D are placed at the location x k from the flat plate leading edge, under zero-pressure gradient in the streamwise direction. At this station, the non-dimensional displacement thickness is δ * /k = 0.6883. A roughness pair consists of two rotating roughness elements with spacing λ, as shown in figure 1. The spanwise spacing between roughness pairs is Λ. The origin of the local x, y, z coordinate system is placed at the location of the roughness element x k for the corotating case and is placed at the centre of the roughness pair for the counter-rotating case. Hereafter, all length scales are non-dimensionalized by the dimensional roughness heightk. The incoming free stream velocity is u ∞ . Depending on the rotation sense of the roughness pair, different types of downstream streaks are created. For the corotating roughness case, the generated downstream vortices rotate in the same direction, which will counteract the momentum exchange effect of neighbouring vortices. In contrast, this momentum exchange effect is intensified in the case of counter-rotating roughness elements. The strength of the generated downstream vortices depends on the rotation velocity of the roughness elements. This can be expressed in non-dimensional form by the ratio of the tangential velocity induced by the cylindrical roughness element to the incoming local velocity at its top, i.e.
where Ω is the angular velocity of the cylinder and u(k) is the unperturbed Blasius velocity at the upper edge of the cylinder. The configurations considered in the present study are summarized in table 1. The counter-rotating cases C1, C2 and C4 differ by their aspect ratio η, the spacings (λ, Λ) and the positive and negative rotation sense as indicated in the table. In the positive   Table 1. Parameters of simulation, non-dimensionalized with respect to roughness heightk = 0.01 m and free stream velocityū ∞ = 0.937 m s −1 . a +: rotation direction with high-momentum creation in the centre. b −: rotation direction with low-momentum creation in the centre. counter-rotating case, a high-momentum fluid area is induced in the middle of the roughness pair. In the negative rotation case, a low momentum fluid is induced in the middle of the roughness pair. The roughness spacings (λ, Λ) in case C1 follow the case C01 from Siconolfi et al. (2015), which demonstrated an effective delay of laminar-turbulent transition induced by TS waves. The roughness elements in the corotating case C3 are equally spaced, to simulate a regular roughness array. A small and a moderate rotation rate are considered for the counter-rotating case, i.e. Ω u = 0.2 and 0.46, respectively. The static roughness array case, i.e. Ω u = 0 in case C3, is taken as reference.

Base flow computation
In LST, the flow quantities q(x, t) = {u, v, w, p}(x, t) are split into a steady baseflow q 0 (x) and an unsteady perturbation q (x, t) similar to the Reynolds decomposition. In this study, the steady baseflows have been computed with the method of selective frequency damping (SFD) (Åkervik et al. 2006). In the SFD approach, a filtered stateũ is introduced and its evolutionary equation is solved together with the Navier-Stokes equations. The following non-dimensional governing equations are obtained: where Re k =ū ∞k /ν is the Reynolds number based on free stream velocity, ω c is the filter cutoff circular frequency and χ is the feedback control coefficient. An overbar denotes dimensional values, hence t =tū ∞ /k and p =p/(ρū 2 ). The additional forcing term on the right-hand side of the momentum equation works as a temporal low-pass filter. From theoretical analysis it is known that the cutoff frequency ω c should be lower than the frequency of the lowest instability, while the feedback control coefficient χ should be higher than the growth rate of that instability (Åkervik et al. 2006). The above governing equations are implemented and solved with the open source OpenFOAM solver icoFoam, which solves the incompressible Navier-Stokes equations using the pressure implicit with splitting of operators (known as PISO) algorithm. The pressure equation is solved by the geometric algebraic multigrid (known as GAMG) solver and the velocity equation by the preconditioned biconjugate gradient (known as PBiCG) solver. The integration domain has a streamwise extent of L x = 320 (−20 ∼ x ∼ 300) and a wall-normal extent of L y = 40. The spanwise extent L z is determined by the roughness pair spacing Λ, as given in table 1. For computational efficiency (Shrestha & Candler 2019), the Blasius boundary layer velocity profile is prescribed to the inlet according to the distance from the leading edge of the flat plate. No-slip wall boundary conditions are applied at the bottom wall. The cylinder can be rotated with constant Ω about its vertical axis leading to a constant tangential velocity u w = ΩD/2 at the cylinder wall. The following direction mixed boundary condition is imposed at the top of the integration domain: Here, n denotes surface norm. At the outlet, any possible reflection is eliminated by using the following advective boundary condition together with a sponge zone with gradually increasing damping strength towards the outlet: A structured grid is used to discretize the above governing partial differential equations. Following the grid convergence study in , 80 grid points are used to resolve the boundary layer, which is sufficient for a laminar flow computation. An expansion (with ratio r = 1.02-1.2) of the grid spacing toward the far field is applied to reduce the total grid points, and 300 equidistant grid points in the circumferential direction are used around the cylinder.

Biglobal linear stability analysis
In this paper, the instability of the boundary layer is studied with the so-called biglobal LST (Theofilis 2011). With this method, the asymptotic behaviour of infinitesimally small 3-D perturbations superimposed on a 2-D base flow can be analysed. Transient growth is not considered here. The analysed base flow is assumed to be parallel in x, i.e. ∂q 0 /∂x ≈ 0.
The perturbations q (x, t) are modelled by the normal mode ansatz, whereq is a 2-D complex amplitude function, α the real streamwise wavenumber, ω = ω r + iω i the complex frequency, and c.c. the complex conjugate. The real part of the eigenvalue ω r corresponds to the angular frequency of the eigenmode/instability, while the imaginary part ω i determines its exponential growth. The instability wave amplifies in case of positive ω i and it decays with negative ω i . Inserting (2.8) into the linearized Navier-Stokes equations and recasting the coefficients results in the following generalized eigenvalue problem for temporal stability analysis: where L and M are the coefficient matrices, given in Appendix A.
The generalized eigenvalue problem is discretized on a y-z plane at consecutive streamwise stations x = constant, which are interpolated from the numerical simulation with a cubic spline. Equation (2.9) is discretized by the Fourier spectral method in the wall-normal direction and the summation-by-parts method (Mattsson & Nordström 2004) in the spanwise direction. The grid is clustered toward the wall following the mapping successfully used by Staudenmeyer, Schnoebel & Rist (2019) for streamwise corner-flow instability with the same eigenvalue solver, with the coefficients b 1 = 0.5, b 2 = −1.2 and the cutoff wall-normal extent y max ≈ 6-10. A grid convergence study has been performed to identify the necessary discretization: N y = 65, and N z = 110. At the wall surface, the boundary condition for (2.9) isq = 0. The von Neumann boundary condition ∂q/∂y = 0 is specified at the top boundary, and for the spanwise boundary a periodic boundary condition is used. The ARPACK (Lehoucq, Sorensen & Yang 1998) routine which uses an implicitly restarted Arnoldi method is employed to solve the generalized eigenvalue problem, (2.9). The solution procedure is implemented in Python code, the validation of which has been performed by comparing the TS-wave amplitude and its growth rate with numerical simulation . Another successful application of the code can be found in Puckert, Wu & Rist (2020).

PKE analysis
For an analysis of the instability mechanisms, the PKE analysis following  is used. The basic idea is to derive the evolution equation for the PKE e = (u 2 + v 2 + w 2 )/2 under the assumption of periodic waves in the streamwise direction. By integrating the PKE in the y-z plane and over a wavelength in the x direction, the following reduced Reynolds-Orr equation is obtained: where E is the total PKE and D the viscous dissipation energy. Here T y and T z are the production energies from the interaction of Reynolds stresses u v and u w with wall-normal shear ∂u 0 /∂y and spanwise shear ∂u 0 /∂z, respectively. The above identity means that the temporal growth rate of the instability mode is composed of three parts: the wall-normal, the spanwise production and the viscous dissipation. Based on the normal-mode ansatz of the perturbation q , the perturbation energy balance terms can be expressed in the form (E, D, T y , where ω i is the temporal growth rate from the normal mode ansatz in (2.8). Substituting this expression back into (2.11), leads to where the terms on the right-hand side are as follows: (2.16) Hereτ xy ,τ xz denote wall-normal and spanwise Reynolds stress components calculated from the instability wave amplitude function, respectively, andd is the viscous dissipation calculated from the perturbation vorticity vector.

Results and discussion
3.1. Base flow In figure 2, vortex structures induced by an isolated static and an isolated rotating cylindrical roughness element are illustrated by means of λ 2 isosurfaces (Jeong & Hussain 1995). This serves the purpose of illustrating the influence of rotation on the vortex structure. While two pairs of equally strong counter-rotating vortices are excited by the static roughness element, the rotating case shows differences with respect to the symmetry and relative strength of its vortices. The naming convention for the induced vortices follows that of Groskopf & Kloker (2016), due to the similarity of the generated vortex structure between the rotating roughness element and their oblique roughness element. For the typical static roughness element (cf. figure 2a), spanwise vorticity is created in front of the roughness element as the boundary layer is blocked by it. The vorticity rolls up and wraps around the cylindrical roughness element, and then develops into a pair of HV legs in the streamwise direction. Another weaker pair of IV is generated directly behind the roughness element, beginning from the downstream reverse flow region. For the rotating case, the vortices are no longer symmetric, see figure 2(b). Instead, one IV gets strengthened and becomes the DIV in the downstream area. On the other hand, the other IV is weakened and rotates around the DIV, hence named SIV. The HV pair is also weakened and vanishes in the near-wake region. The formation of the DIV is well illustrated in figure 3, and can be summarized as follows.  The single leg of the streamwise DIV is stronger compared with those of the HV in the static case, thus it is more effective in pulling high-speed fluid towards the wall and pushing low-speed fluid to the outer region of the boundary layer. In other words, the lift-up effect (Landahl 1990) is stronger in the rotating roughness case compared with its static counterpart, therefore it is able to create stronger streamwise streaks. The locally accelerated flow at the bottom also creates a slight cross-flow in the downstream wake of the cylinder.
It turns out that the characteristics of the DIV created by the isolated roughness element can be used to control velocity streaks. It is the purpose of this paper to investigate the linear stability property of a boundary layer with different combinations of rotating roughness elements, i.e. corotating roughness arrays and counter-rotating roughness pairs. Figure 4 presents slices of the steady base flow obtained from the SFD solver with the identified vortices, high shear regions (by means of I 2 = constant, see Meyer (2003)) and the velocity gradients (∂u/∂y and ∂u/∂z). Two slices are shown for each case, one at x = 5 (in the two top rows) the other at x = 40 (in the two lower rows). In the near-wake region (x = 5), the location of the HVs on the sides of the DIV is observable for all cases. Since the DIV rotates in a clockwise sense (see figure 4c for instance), a stronger high velocity streak is created on its right-hand side and a stronger low velocity streak is created on its left-hand side compared with its static counterpart. This feature is more obvious at the downstream slice, i.e. at x = 40. This alternating high-low-high velocity streak pattern would attenuate the strength of the DIV, if the neighbouring corotating roughness elements were too close to each other. Contrary to this attenuation effect, the counter-rotating roughness pair intensifies either the high-speed or the low-speed velocity streak, depending on the rotation direction. In figure 4(a), the counter-rotating roughness pair rotates in the direction which intensifies the high-speed velocity streak and the low-speed velocity streaks on the side are dissipated by viscosity. It is typical for streaky flow to have a high-shear region above the low-speed velocity streak. At slice x = 5, it is obvious that a high-shear region (I 2 = −450) composed of wall-normal (∂u/∂y) and spanwise (∂u/∂z) terms exists above the DIV. This structure belongs to an inflectional velocity profile, where inviscid inflectional instability could be expected. At slice x = 40, this high-shear region with inflection point still exists for corotating case C3 due to the long-living DIV, while for the counter-rotating case C1 only the intensified high-speed streak persists. The latter happens to be a preferable velocity profile with respect to the stabilization of TS-waves (Siconolfi et al. 2015). Figure 5 shows the streamwise evolution of the velocity gradient maxima in y-z planes for all three types of streaky flows by corotation (case C3), positive counter-rotation (case C2) and negative counter-rotation (case C4) of the roughness elements. It is clear that the primary components of shear are the ∂u/∂y and ∂u/∂z terms in both static and rotating cases. The other terms are several magnitudes lower. Comparing figure 5(b) with figure 5(a), it is clear that the rotation effect not only promotes these two terms, but also slightly increases the ∂w/∂y and ∂v/∂z terms. However, the streamwise gradient terms ∂/∂x remain negligibly low except for the very short near-wake region (x < 5), i.e. the reverse-flow region. This observation validates use of the parallel flow assumption as in Groskopf & Kloker (2016), for instance. Since the production of the PKE is proportional to the shear magnitude, it can be inferred that the associated instability mode is fundamentally determined by the wall-normal and spanwise productions which justifies using the reduced Reynolds-Orr equation, i.e. (2.11) for perturbation energy analysis.
To quantify the amplitude of velocity streaks, the following definition of Groskopf & Kloker (2016) is used: where u is the spanwise mean value, which represents the 3-D baseflow deformation caused by the low-and high-speed streaks. In figure 6, the streamwise evolution of the streak amplitudes for the cases from table 1 are compared. A rapid growth of the streak at the roughness element (x = 0) can be observed. A subsequent transient growth of the streak amplitude begins from x ≈ 30 for case C1 with Ω u = 0.1. With increasing rotation rate, this starting position moves forward until x ≈ 5. Such tendency is also observed for other cases, especially case C3. For cases C1 and C3 with Ω u = 0.46, the streaks reach a peak (u st = 0.32) at around x ≈ 55 and thereafter decay gradually. With increasing rotation rate, the streak amplitudes of these two cases evolve similarly with an almost constant offset. The streak evolution for cases C2 and C4 exhibit some differences. Here, the velocity streak maxima occur closer to the roughness for lower rotation rates (Ω u = 0.2). For higher rotation rate (Ω u = 0.46), an additional peak can be observed at x ≈ 100 for case C4 and x ≈ 200 for case C2. Both postpone the decay of the streaks leading to a higher streak amplitude.

Linear stability analysis
Linear stability analysis is performed at the same streamwise location x = 100, where the base flow is quasi-parallel, for all cases presented in table 1. Two types of modes are identified, the viscous TS-like mode and the inviscid inflectional mode, both clearly shown in figure 7 for four cases. The TS-like mode contains a 3-D distortion relative to the 2-D TS mode of the Blasius boundary layer, due to the influence of the velocity streaks. The TS-like mode's amplitude maximum is located close to the wall and the mode travels downstream with a phase speed c ph = ω r /α ≈ 0.3-0.4. Due to the influence of the streaks, the TS-like mode exhibits alternating signs in spanwise direction, as shown in figures 7(a) and 7(b). This amplitude modulation is mainly induced by the low-speed streak, which is clearest in figure 7(b) at z = ±2. The inviscid inflectional modes in figures 7(c) and 7(d) are found to reside on the high-shear region which is identified by the I 2 criterion. Since the high-shear region is induced by the roughness element, the mode is hereafter called the roughness mode. The roughness modes display a primary structure in the high-shear region and secondary structures which are in antiphase to the primary at the sides of the high-speed streaks. The shown roughness modal structures for cases C2 and C4 are similar to each other. If the two low-speed regions in C2 move closer to each other, the shown modal structures would coalesce to the shape of case C4. The amplitude of thev,ŵ components of the shown TS-like mode is quite small compared with the correspondinĝ u components, whereas thev component of the shown roughness mode is similar to the amplitude of theû component. What is more, theû,v components always appear in antiphase. It is also to be noted that the typical symmetric sense of the roughness modes as found in the static roughness case is lost, the mode shown in this paper is always the most amplified one.
The effect of different streaks on the boundary layer stability is studied by computing the eigenvalue problem (2.9) at streamwise location x = 100 for a set of wavenumbers α. The temporal amplification rate ω i and corresponding phase velocity c ph are shown in figure 8 for both TS-like mode and roughness mode. The TS mode of the undisturbed Blasius boundary layer, which is slightly amplified around α = 0.35, is evaluated at the same streamwise location for reference. The two types of modes differ by their phase velocity c ph , such that the TS-like mode travels at a lower speed (c ph ≈ 0.3-0.4) than the roughness mode (c ph ≈ 0.6-0.8). The identified phase speed for the roughness mode is similar to the results of Di Giovanni & Stemmer (2018), who found that a roughness-induced mode travels at c ph = 0.68u ∞ in a hypersonic boundary layer. For the TS-like mode, all streaks are capable of attenuating the TS-like mode except for the case C4 with high rotation rate (Ω u = 0.46), where the roughness mode gets mostly amplified with a relatively high amplification rate. Only the roughness modes for cases C1 and C3 with low rotation rate (Ω u = 0.2) remain damped. While at lower rotation rate (Ω u = 0.2), the amplified mode peaks at a wavenumber around α = 0.5, the roughness mode with higher rotation rate (Ω u = 0.46) peaks at a larger wavenumber around α = 1.0. It is also found that the amplification rate for roughness modes increases with rotation rate. At the same time, the most amplified wavenumber moves to higher values. However, among all cases, the TS-wave attenuation effect is only observed to increase with increasing rotation rate in configuration C1, i.e. positive counter-rotating roughness elements. Such tendency is not that obvious for other cases. This indicates that only thin roughness elements with low rotation rate are suitable for boundary layer instability attenuation or possible transition delay, otherwise the additional roughness mode is always a threat leading to a fast-growing inflectional instability. In figure 9, temporal stability diagrams obtained by tracking eigenmodes in the (α, x)-space are presented. The TS-like modes of case C1 with Ω u = 0.1 and 0.2 are compared with the TS mode flat-plate boundary layer in the upper row. The x position corresponding to the critical Reynolds number for the unperturbed TS-mode is at x = 85 with wavenumber α = 0.37. With a slight rotation (Ω u = 0.1), the critical position moves forward to x = 38 and the wavenumber reduces to α = 0.34. However, the overall amplification rate ω i is greatly reduced, as the maximum is decreased by 40 % from 2 × 10 −3 to 1.2 × 10 −3 . With further increase of the rotation rate to Ω u = 0.2 the TS-like mode is fully damped in the studied streamwise range. The shown contour levels are hence all negative in figure 9(c). The TS-like mode for even higher rotation rate Ω u = 0.46, which is not shown here for convenience, is also fully damped in the studied streamwise range.
The stability diagram of the most amplified roughness mode for case C2 is compared with the static roughness array case from C3 in the lower row of figure 9. The roughness mode seems to be mainly amplified directly behind the roughness element. However, as indicated in the velocity gradients (cf. figure 5), the parallel flow assumption most probably fails in the direct vicinity of the roughness wake, this must be treated with caution. Hence, the roughness mode is only tracked from x = 5 downstream where the initial transients from figure 5 have decayed. Compared with the TS-like mode, the roughness mode is amplified more strongly and within a broader range of wavenumbers. Since the roughness elements impose an abrupt geometric discontinuity to the boundary layer flow, extreme high shear can be found in the near-wake region, where multiple unstable modes can be found in the LST analysis. While the amplified region of the static reference case only extends to x = 25, the rotation effect extends the amplified region to the full streamwise region for both cases. Besides, there are three notable features in the stability diagram of the roughness mode. The first is that the most amplified instability moves to a higher wavenumber. The second is that the highly amplified near-wake region moves closer to the roughness elements. The third feature is that a second peak appears at x = 60 for higher rotation rate in figure 9( f ). In fact, this peak already emerges at x = 100 for the case with Ω u = 0.2 in figure 9(e). The rotation effect moves the second peak to higher wavenumber and forward to the roughness elements as well. This observation is consistent with the observation of a local velocity streak amplitude maximum in figure 6(b).
3.3. Instability mechanisms Perturbation kinetic energy analysis is performed to get a deeper understanding of the linear instability mechanisms. Results are shown in figure 10. The integral energy terms (T y ,T z ) are normalized with respect to their corresponding dissipation termsD. In case C1, the amplified reference TS mode is also plotted together with the TS-like mode. As the TS mode is 2-D, the spanwise production termT z is negligibly small. Since the wall-normal production termT y exceeds the dissipation termD, the TS mode is amplified as a result of extracting energy from the baseflow. At very low rotation rate (Ω u = 0.1), theT y term is slightly lower than that of the TS mode and theT z term is found to be negative. This negative spanwise production term is considered to be the responsible mechanism for TS-wave attenuation and transition delay with static roughness elements Fransson et al. 2006). The reference case with static roughness (Ω u = 0) is plotted with C3 as well. With this set-up, the spanwise production energy is also found to be negative, but it is not strong enough to surpass the wall-normal production energy. With increasing Ω u , the negative spanwise production termT z becomes positive in case C1. On the other hand, the wall-normal production termT y reduces dramatically. Therefore, the sum of production terms becomes smaller than the dissipation term, resulting in an overall attenuation of the TS-like mode. Such mechanisms, i.e. reduction of wall-normal production and amplification of spanwise production perturbation energy, are also observed for the TS-like mode in cases C2, C3 and C4. In cases C2 and C4, the wall-normal production term is even reduced to a negative value at high rotation rate (Ω u = 0.46), although the negative side is the overshoot of the spanwise production term over the dissipation, see case C4 with Ω u = 0.46 in figures 10 and 8. Perturbation kinetic energy analysis results for the roughness mode are shown in the right column of figure 10. In case C1, both production terms are relatively low when the rotation rate is small (Ω u = 0.1, 0.2). However, the wall-normal termT y rises significantly at high rotation rate (Ω u = 0.46). Such notable increase of the wall-normal production term is also evident for all other cases, especially in cases C2 and C4 where the wall-normal productions rise up to four and six times the dissipation term, respectively. Although the increase ofT y in case C3 is relatively small, it still reaches 1.5 times the dissipation term. In addition, the spanwise production termT z is also increased to 1.8 times the dissipation term. This unique feature of the roughness mode is due to the fact that case C3 is the only corotation set-up in which both wall-normal shear ∂u/∂y and spanwise shear ∂u/∂z are significantly increased. To conclude, the roughness mode can be destabilized by the streaks induced by rotating roughness elements and the according destabilization effect is much more significant than the stabilization effect of the TS-like mode, in agreement with the amplification rates in figure 8. In § 3.4 we shall investigate whether this destabilization leads to immediate laminar-turbulent transition at the roughness elements or not. The above mechanisms can be further explained by examining the spatial distribution of the PKE terms for each type of mode. Contour plots of the PKE terms for the TS-like mode from case C1 are shown in figure 11. BothT y andT z are the product of Reynolds shear stresses (τ xy ,τ xz ) and baseflow velocity gradients (∂u/∂y, ∂u/∂z), both of which can be positive or negative. The corresponding signs of the involved Reynolds shear stresses and the baseflow velocity gradients determine the sign of the PKE production. As mentioned above, the spanwise production term for the TS mode is zero, since there is no spanwise velocity gradient. However, neither the spanwise velocity gradient nor the spanwise Reynolds shear stress is zero in streaky flows, the multiplication of the two terms results in unavoidable non-negative spanwise production. Figure 11(a) shows the spatial distribution of the PKE terms for case C1 with low rotation rate (Ω u = 0.1), which is a case closest to the result of the TS mode. As already explained, the spanwise production term for the TS mode is zero and the wall-normal production is positive throughout the boundary layer. Due to the effect of the velocity streaks, the otherwise constant sign of the wall-normal production is broken, see the top subpanels of figure 11. The positive τˆx y · (∂u/∂y) τˆx y · (∂u/∂y) τˆx y · (∂u/∂y) τˆx y · (∂u/∂y) and negative production terms cancel each other in the spanwise mean, leading to the reduction of the wall-normal productionT y for the TS-like mode as shown in figure 10. The spanwise production is also composed of parts with opposite signs, which inhibit the growth of the integral value as rotation rate increases. This explains why the spanwise productionT z of the TS-like mode for case C1 remains almost unchanged for Ω u = 0.1 and Ω u = 0.2 in figure 10. With higher rotation rate (Ω u = 0.46), the negative parts of the spanwise production become negligibly low, making the positive parts dominant, which explains the rapid rise of the integral value. The rise of the spanwise production can be so strong that it again destabilizes the TS-mode, see the TS-like mode of case C4 in figures 10 and 8. The spatial distributions of the PKE terms for the roughness mode of cases C2 and C3 are shown in figure 12. For the inflectional instability, the productions are mostly concentrated above the low-speed streak where inflection points can be found. The negative part of the wall-normal production is lower compared with the spanwise production for the case shown. As rotation rate increases, the wall-normal production intensifies and results in a higher integral value. However, the robust negative part in the spanwise production restricts the rise of the integral value. This explains the behaviour of spanwise productionT z for the roughness mode in figure 10. Nevertheless, the high wall-normal production has already destabilized the roughness mode. The rise ofT z in case C3 is due to the unique feature of the velocity streak induced by the corotating roughness elements, that both wall-normal and spanwise production are intensified, see figures 12(c) and 12(d).

Exclusion of a possible bypass transition
The results presented so far have identified two opposing effects: stabilization of TS-modes and destabilization of roughness modes. Especially, the latter appears to be strongest Amplitude Figure 13. Snapshot (y = 1) of disturbance field u at t/T 0 = 4.1 for case C2 with Ω u = 0.46 and disturbance strip amplitudeÂ = 0.5u 0 . The disturbance strip is marked by vertical lines ahead of the roughness pairs. Its spatial distribution and temporal modulation are shown above. Two spanwise periods Λ are shown. Note that the visualization is presented with a logarithmic scale so as to exaggerate the otherwise imperceptible instability in the near-wake.
towards the roughness elements where the local stability analysis is no longer valid. Therefore, in order to exclude both, a possible bypass transition and a possible global instability, a DNS is conducted for case C2 with pairs of counter-rotating roughness elements and Ω u = 0.46, for which the previous LST predicts that the amplification rate of the roughness mode is 100 times larger than that of the T-S mode. The steady-state base flow is used as the initial field, and the nonlinear Navier-Stokes equations are marched in time. A blowing-suction disturbance strip is placed ahead of the roughness elements from x = −11 to −10 to trigger disturbances with a broad spectrum. The wall-normal velocity component mapped to the wall within the disturbance strip is specified as (cf. Kurz & Kloker 2016) v(x * , y * = 0, z * , t) =Â sin(πx * ) exp(−(πx * ) 2 / √ 2) 1 50 50 n=1 cos(2π f 0 nt), x * ∈ [−1; 1], (3.2) whereÂ is the amplitude and f 0 = 0.25 Hz is the fundamental frequency. Frequencies with n = 1-50 are superimposed to produce a periodic pulse disturbance which will excite a broad spectrum of instability waves. The frequencies of this spectrum are chosen such that its lower part covers the dominant TS-modes and its higher part covers the dominant roughness modes. Their different frequencies can be inferred from their different wavenumbers in figure 9. In addition, the spatial distribution and temporal modulation of the disturbance strip is switched on and switched off after four fundamental periods T 0 = 1/ f 0 , as shown by the dashed envelope curve together with a snapshot of the disturbance field u = u − u 0 at wall distance (y = 1) in figure 13.
In the latter, two spanwise periods are presented such that the impact of both sides of the roughness-element pairs can be seen, the acceleration in-between two cylinders and the deceleration between two pairs. Note the different distances between the four cylinders shown. Initially, a wave package observable at x = 2-20 is triggered when the pulse passes the roughness elements. Its wavefront is not only interrupted by the individual wakes of each cylinder but it oscillates in antiphase, depending on its position relative to the cylinder pairs. Oscillating structures start to occur at x ≈ 25 as a result of the dispersive and convective nature of the wave package. Farther down they lead to braid-like structures due to amplification of the roughness mode. A disturbance maximum occurs between x ≈ 100 and x ≈ 170. Here, we notice that the near-wake of the roughness elements is dominated by instability waves instead of a bypass transition, despite the large disturbance amplitudeÂ = 0.5u 0 and the broad band of forcing.
On the other hand, the near-wake region is reported to be the 'wavemaker' for self-sustained global instability (Loiseau et al. 2014). This means that the near-wake can be the starting region of a global instability, although immediate laminar-turbulent transition may not occur. However, this hypothesis can be excluded as well for the present case, as is shown in figure 14 which presents x-t diagrams of disturbances sampled at ( y, z) = (1, 0) over time. Two disturbance strip amplitudesÂ = 0.05u ∞ andÂ = 0.5u ∞ are shown to present their influence on the transition front. In both cases, the disturbances are observed to grow gradually while they are convected downstream. Bypass transition does not occur in the near-wake. Increasing the disturbance amplitude 10 times does not change the convective nature of the instabilities either, but the transition front moves to approximately x ≈ 200. The decisive test is to switch the forcing off at t/T 0 = 4. Figure 14(b) shows that all perturbations are then convected away by the flow. This is the final proof that no global instability exists.

Conclusions
The present investigations have shown that rotation has a decisive influence on the wake behind wall-mounted, finite-length cylindrical roughness elements in a laminar boundary layer. Local acceleration by the circumferential velocity shifts the separation zone at the wall to one side of the cylinder. This process prevents the otherwise typical formation of a HV at the bottom of the cylinder. Instead, a single vortex (a DIV) is formed, similar to the case of roughness elements in a boundary layer with cross-flow or oblique roughness elements in a 2-D boundary layer, see Groskopf & Kloker (2016). Compared with the HV in the non-rotating case, the DIV is stronger and extends farther downstream. In addition to arrays with corotating roughness elements, cases with roughness element pairs have also been investigated. For these, the distances and the sense of rotation with respect to each other have been varied as well. All these investigations show that the deformation of the boundary layer by flow-induced streaks can be controlled considerably by the rotational speed of the roughness elements.
More than 10 cases are presented which have been studied in numerical simulations of the steady laminar base flow and a subsequent biglobal linear stability analysis. The latter has been performed at constant downstream distances from the roughness elements but in a wide range that begins at five roughness heights behind their centre. Use of the parallel-flow approximation for LST can be justified by comparing the much smaller gradients in streamwise direction with respect to other directions (see figure 5). Two unstable or least stable modes are observed: the classical TS-mode of boundary layer instability and a roughness mode due to the presence of inflection points and dominating wall-normal as well as spanwise shear in the streamwise velocity component of the laminar base flow. The first is deformed three-dimensionally by the roughness or its wake. The latter occurs initially only in the immediate vicinity of the roughness. The amplification maxima for both modes occur at different wavenumbers and hence different frequencies, with the roughness mode being 5 to 10 times higher.
As the cylinders are rotated, the TS mode is attenuated and finally damped, while the amplification rate of the roughness mode continues to increase. An analysis of the PKE has revealed the specific contributions of the wall-normal and spanwise shear components to the overall amplification or attenuation of the different modes under the influence of rotation. The investigations also show that the control of the roughness mode is less dangerous for slim roughness elements compared with thicker ones. The risk of a possible bypass transition due to a strong amplification of the roughness mode in the immediate vicinity of the roughness elements, where the local stability analysis could not be performed, was excluded by a DNS study.
Overall, the present method complements earlier investigations by Fransson et al. on TS-mode attenuation with a series of cylindrical roughness elements and shows the possibility of increasing their effect by rotation. Based on the observation that some modes can be further destabilized by the control, one may propose that the current control method could also be used to trigger laminar-turbulent transition instead of delaying it. Construction of a device with rotating roughness elements appears feasible, at least to the present authors. In addition to what has been shown here, we have already observed in another DNS study that rotating cylindrical roughness elements could also be considered in a turbulent flow for the purpose of delaying boundary layer separation. where L T = iαu 0 + v 0 (∂/∂y) + w 0 (∂/∂z) − 1/Re(∂ 2 /∂y 2 + ∂ 2 /∂z 2 − α 2 ).