Mixed mode transition in boundary layers: Helical instability

Recent (Bose \&Durbin, \textit{Phys. Rev. Fluids}, 1, 073602, 2016) direct numerical simulations (DNS) of adverse- and zero-pressure-gradient boundary layers beneath moderate levels of free stream turbulence ($Tu$ $\le 2\%$) revealed a \textit{mixed mode} transition regime, intermediate between orderly and bypass routes. In this regime, the amplitudes of the Klebanoff streaks and instability waves are similar, and potentially can interact. Three-dimensional visualizations of transitional eddies revealed a helical pattern, quite distinct from the sinuous and varicose forms seen in pure bypass transition. This raises the fundamental question of whether the helical pattern could be attributed to a previously unknown instability mode. Two-dimensional stability analyses are performed herein for base flows extracted from DNS flow fields. The three-dimensional structure of the eigenfunction of the most unstable mode indeed reveals a helical pattern. These instability modes are obtained significantly earlier than helical structures are seen in flow visualizations, even for innocuous base flow structures. It is inferred that the instability mode is the underlying cause of helical patterns that emerge in the transition region. The streak configuration leading to the formation of the helical mode instability is different from those leading to sinuous and varicose modes in previous studies of pure bypass transition. Thus, the mixed mode precursor is the distinctive cause.


I. INTRODUCTION
Direct numerical simulations (DNS) have provided a detailed description of the laminar to turbulent transition process in attached subsonic boundary layers.Both orderly and bypass transition have been investigated in detail.The theory of orderly transition is well established.It starts with exponential growth of instability waves, followed by a much faster secondary instability process forming Λ-vortices.Finally, lifted Λ-vortices break down into patches of turbulence -that final stage is impervious to theory.When transition does not emerge from growing instability waves it is classified as bypass transition.In most of the literature, bypass transition is synonymous with transition due to incident free-stream turbulence (FST).Over time various theoretical concepts have emerged [8].They explain the occurence of streamwise elongated perturbations, called Klebanoff streaks, as precursors to transition.These disturbances amplify via the lift up mechanism, or algebraic growth, until their amplitudes saturate.Klebanoff streaks undergo secondary instability.The nature of that secondary instability is less clear cut than the Λ-vortices of orderly transition.Analysis have characterized it in terms of sinuous and varicose modes, and outer and inner instabilities [1,9].
Experimental evidence suggests the existence of flow regimes in which both instability waves and streaks are involved in the transition process [13,14].Especially, this can occur in adverse pressure gradient (APG) under low to moderate levels of free-stream turbulence.
A recent DNS study [4][5][6] revealed an intriguing transition process: beneath moderate levels of FST (intensity T u ≤ 2%), and in the presence of instability waves at super-critical Reynolds numbers, a truly mixed mode transition process was reported.Three-dimensional visualizations of the perturbation fields revealed a features quite distinct from the sinuous and varicose modes of pure bypass transition.The flow structures resembled a helical pattern.A deep lying critical layer was reported for these helical eddies, and therefore, these were classified as inner instabilities-in the terminology of Vaughan and Zaki [12]-but, again, distinct from what had been seen previously.
Those empirical observations pose the fundamental questions, does a helical, secondary streak instability exist?and, is it the cause of the patterns seen in flow visualization?Those questions prompted the present study.We analyze the secondary stability of boundary layers in the mixed mode transition regime.Both zero and APG boundary layers are considered.
The approach parallels that of Hack and Zaki [9]: observed helical patterns are tracked back in time to identify their precursors.These provide base flow fields to which stability analysis is applied.An iterative technique is used to extract those eigenmodes having the highest growth rates.
As will be shown, the eigenfunction solution that is obtained is indeed of helical form.These eigenmodes are obtained even for apparently innocuous base streaks, significantly earlier than the helical instability is seen in DNS.This substantiates the hypothesis that the helices that were seen in computer simulations can be attributed to a new type of instability.
The streak configuration leading to the inner helical mode is distinct from the precursors to previous solutions of secondary instabilities.Outer sinuous and inner varicose streak instabilities have been attributed to the inflection points in the spanwise and wall-normal profiles of the streamwise velocity, respectively [1,9].The inner helical mode appears at the intersection of low-and high-speed streaks, where spanwise and wall-normal shear are both important.Such unique orientation of the streaks is due to the high growth rate of the mixed mode precursors, formed as a consequence of the instability waves interacting with the Klebanoff streaks.As is inevitable, the new instability, discovered herein, is a consequence of the precursor differing from previous studies.

A. Direct Numerical Simulations
The direct simulations were performed in the mixed mode transition regime for both ZPG and APG boundary layers in order to create the precursor fields.A detailed description of the numerical methods, boundary conditions, and the method for the generation of inlet turbulence can be found in [3,6] and are not repeated here.For all simulations, inlet FST with a prescribed intensity (T u ≡ 1 3 u ′ i u ′ i /U ∞ = 1% and 2%) is superimposed on a Blasius boundary layer [7,10] at a super-critical Reynolds number, Re b = U ∞ x 0 /ν ≈ 398.In ZPG simulations, the inlet additionally includes a TS-wave with amplitude A = 0.5%U ∞ .
The domain is rectangular for the ZPG simulations with size 320δ 0 × 24δ 0 × 24δ 0 (δ 0 is the boundary-layer thickness at inlet) in the streamwise, wallnormal, and spanwise directions.
The grid has 1, 537×161×241 points.The computational domain and grid specifications for the APG simulations are the same as the ZPG simulations except that the APG is induced by a curved upper boundary characterized by the Hartree parameter, β H = −0.14[6].The largest grid spacings in viscous units are reported in table I. Despite the use of different grids, the mean skin-friction curves from the present simulations are in good agreement with those from previous studies [3,5,6].Once a statistically stationary state was reached in the simulations, instantaneous snapshots were stored at constant time intervals of 1.8δ 0 /U ∞ , over a period spanning at least three flow-through times, for stability analysis.

B. Stability analysis
Base flow profiles are extracted in planes normal to the streamwise direction, from instantaneous DNS fields.Only the streamwise component of the velocity field is included, where, k x is the streamwise wavenumber, also extracted from the DNS.These idealizations can be formally justified, as the streaky base flow is predominantly in the streamwise direction, and changes far more slowly than the growth rate of the perturbations in x and t [1,9].
Consequently, the following set of coupled linear equations is obtained for the perturbations, . 0 indicates a zero in a matrix element of the discretized equations.
The perturbations satisfy the following boundary conditions, û2 , v2 , ŵ2 = 0 at y = 0 (3) As in the DNS, a periodic boundary condition is applied in the spanwise direction.
Following Barkley and Henderson [2], the set of autonomous equations may be written Denoting e Ct = A(t), formal integration of (5) yields The eigenspectrum of the matrix exponential, A(t), dictates the time evolution of the eigenfunctions of A(t).The real and imaginary parts of ω provide the frequency and growth rate of disturbances, respectively.Instead of an explicit computation of the eigenvalues of the large, two-dimensional stability problem, an iterative method is used following Sorensen [11].In this method, eigenvalues with the highest growth rates, and also those that yield phase speeds similar to that extracted from DNS, are extracted with the implicitly restarted Arnoldi method (IRAM).The Francis QR factorization technique is used to extract eigenvalues.A fuller description of the iterative scheme is provided in Hack and Zaki [9], and is not repeated here.
Formally, the iterative eigenvalue solver applies the linear vector differential operator A(t) repeatedly to an initial disturbance field as in (6).In practice, this is equivalent to repeatedly integrating ( 1) and ( 2) in time, with the boundary conditions ( 3) and ( 4).Time integration is with a fractional step method, similar to the solution scheme used for the DNS: advection uses the second order Adams Bashforth scheme, and diffusion uses the implicit Crank Nicolson scheme [5].The spatial discretizations in both wall-normal and spanwise directions is by second-order central differencing.For both the momentum equation and the pressure correction steps, periodicity in the spanwise direction is exploited by solving for coefficients of spanwise Fourier modes, followed by an inverse transform.161 and 482 grid points are used to discretize the equations in the wall-normal and spanwise directions, respectively.
An important consideration is the calculation of the frequency of the instability mode (equation 6) where, θ is the phase of a complex eigenvalue of A(t), and ∆t is the interval of time integration.The iterative scheme is based on most unstable mode dominating at large times.
Care is needed in choosing ∆t, as an inverse trigonometric function is required to compute the modal frequency.∆t was chosen in conjunction with applying A(t) multiple times, to ensure |ω r ∆t| ∼ π/4 which is within the range of the inverse trigonometric function used.states is extracted from DNS by streamwise displacement of the sampling plane at the phase speed of the disturbance observed in DNS.In this way, evolution of a helical disturbance from quasi-parallel stability analysis can be compared, in an approximate sense, to DNS.

III. RESULTS & DISCUSSION
In the mixed mode transition regime, the primary flow field is a juxtaposition of Kelbanoff streaks and instability waves on the boundary layers.In the APG simulations, the instability waves arise spontaneously [6], while in ZPG, an unstable TS wave is injected at the inlet of the simulation domain.on a high-amplitude, mushroom-like, low-speed streaky structure.Unlike the inner varicose instability described by Hack and Zaki [9], which forms when a high-speed streak climbs up a low-speed streak from behind, the inner helical instability forms when a low-speed structure overlaps a high-speed structure from the side (also see figure 2 in [4]).The low phase speed, c ≈ 0.546 of this instability mode confirms that the helical mode is an inner instability.The critical layer for the instability is shown in figure 2(b).
x, y-planes through a helical mode [4,5] reveal a Kelvin-Helmholtz type instability which is similar in description to the inner varicose mode.However, unlike the inner varicose modes [9], figure 2(b) shows that the low-speed structure positions on top of the high-speed structure for a helical mode.
Figure 3 is a three-dimensional view of velocity components of the most unstable eigenmode, over two streamwise wavelengths.This clearly reveals a helical pattern.It is quite distinct from the previously reported sinuous and varicose type instabilities: contrast to figure 20 of Hack and Zaki [9].
The disturbance field from the stability calculations is in qualitative agreement with the helical mode observed in the DNS of figure 1(c).The chirality of the helix depends on the direction of the side-wise overlap of the low-and the high-speed streaks.For this particular instance of the helical mode, the overlap is from the right of the low-speed streak, and therefore, looking downstream, the helix wraps clockwise around the low-speed streak.In cases where it is on the left, the helical instability rotates counter-clockwise.
Although the most unstable helical eigenmode (hereafter called mode 1) is shown in figures 2 and 3, the stability analysis predicts multiple unstable modes for base states like that shown in figure 2. The phase speed, c, and the temporal growth rate, ω i , for the three most unstable modes of the two base states extracted at time instants t 0 + 16.2 (figure 2) and t 0 + 27 (at this instant, the helical instability is about to become prominent in DNS) are plotted against the streamwise perturbation wavelength, λ x /δ 0 , in figure 4. The second most unstable mode (mode 2) has phase speed similar to mode 1.However, its growth rate is lower for all λ x /δ 0 .For mode 1, highest value of ω i is obtained for λ x /δ 0 ∼ 1.25 at t 0 + 16.2, and for λ x /δ 0 ∼ 1.55 at t 0 + 27.Therefore, the stability analysis predicts an increase in the streamwise wavelength of the helical instability with downstream distance.This is qualitatively consistent with the DNS.The third most unstable mode (mode 3) has a significantly lower phase speed than modes 1 and 2, and is stable for all λ x /δ 0 at inception of the helical instability, t 0 + 27.At the earlier instant, t 0 + 16.2, this mode is unstable for λ x /δ 0 > 1.5.The mode shapes for modes 2 and 3 are presented in appendix A.
The stability analysis was performed on a series of base states extracted from DNS following the evolution of the helical mode.The streamwise wavelength for the stability analysis is λ x /δ 0 = 1.55, for which ω i was highest at x/δ 0 ∼ 185.8 and t = t 0 + 16.2 (see figure 4), and is approximately the streamwise average wavelength of the helical mode in DNS.
The results (phase speed and growth rate) are compared with those extracted from DNS in figure 5. Phase speed is extracted from the DNS data by following the helical mode precursor backward in time and applying a finite difference formula.The spatio-temporal growth rate is extracted by following the peak of v ′ corresponding to the instability in a wall-parallel plane located at y/δ 0 = 1 (the plane is close to the critical layer obtained for the inner helical mode shown in figure 2) and using the formula In figure 5(a), phase speeds extracted from DNS show four discrete levels, indicating a superposition of instabilities propagating downstream.The phase speed of the unstable modes from stability analysis remains approximately constant with downstream distance.
The superposition of the unstable modes obtained from stability analysis has some elements of consistency with the DNS.The growth rate, ω i for the helical mode is significantly larger than those extracted from DNS.The temporal growth rate of mode 3 is in good agreement with the DNS; but, this mode has a smaller phase speed than those extracted from the DNS.Mode 2 is initially stable, and only becomes unstable at x/δ 0 ∼ 184.2.Note that, as shown in figure 2, the helical mode is the most unstable well before it becomes clearly seen, at t 0 + 27 and x/δ 0 ≈ 192.Its high growth rate means that the time interval between its inception and breakdown is very short.In DNS's, the breakdown of the helical mode after its genesis is also very short; it can only be clearly identified in iso-surface plots at t 0 + 27, and it breaks down at ∼ t 0 + 36.
In figure 5, the growth rate of the most unstable helical mode does not match with those extracted from DNS.There are two reasons: firstly, the stability analysis is linear, disturbances are monochromatic, and the growth rate is only temporal.The DNS is nonlinear and the growth rate is spatio-temporal.Secondly, and perhaps more importantly, the length and time scales of the helical instability are an order of magnitude smaller than the base-state instability waves (the scales of the Klebanoff streaks are even larger).Therefore, at the beginning of the transition region, the scales corresponding to the helical instability have negligible energy.Once they become unstable, their energy grows at exponential rate.
Hence, although the helical mode is the most unstable, at first, the primary instability has most of the disturbance energy and dominates the growth rate observed in the DNS.As the disturbance energy at the helical-mode scales become of similar magnitude as the primary instability, they show up in the flow fields, but quickly break down to turbulence.This has been previously verified by plotting spectra,, e.g., figure 17 in [6].
A. Helical mode in ZPG Bose and Durbin [4] also observed the helical mode in a ZPG boundary layer, exposed to free-stream turbulence and a Tolmein-Schlichting wave.The stability results for this case is presented in figure 6.It can be hypothesized that an inflectional base flow is not necessary for the inception of helical instability.
Figure 6(a) shows the base flow configuration at time t 0 in a case for which a helical instability was noticed at t 0 + 10.8 in the DNS-see figure 6(b).The cross-stream plane in figure 6(c) shows that the base streak configuration leading to the genesis of a helical mode is same as that previously described for APG (figure 2).The most unstable helical eigenfunction, 6(c), is also similar to the APG case.The eigenmode is obtained for an input λ x /δ 0 = 1.6 which is similar to the average wavelength of the helical instability seen in figure 6(b).The eigenmode from the stability analysis corresponds to c ≈ 0.624, and ω i ≈ 0.39.The phase speed and growth rates extracted from DNS at this streamwise station are c ∼ 0.63 and ω i ∼ 0.01.These results are also consistent with the analysis presented for APG. Figure 6(d) shows the top view of surfaces of iso-u ′ 2 values depicting the helical mode in ZPG.

IV. CONCLUSIONS
The fundamental question addressed in this paper is whether the helical patterns that have been seen in flow visualizations of mixed mode transition can be traced to an underlying instability.That has been shown to be the case.A new instability was found, that is distinct from the sinuous and varicose secondary instability modes that have previously been cited in studies of bypass transition.The helical secondary instability is peculiar to mixed mode transition-that is, when precursors consist of a combination of Klebanoff The stability analyses predict either two or three unstable modes, of which the helical, inner mode has the highest growth rate.Helical instability occurred significantly earlier than where helices were detected in DNS, indicating that it is the underlying cause.
The growth rate predicted by the stability analysis for the helical mode is larger compared to that extracted from DNS.This too is evidence that the instability gives rise to the observed vortical structures.The initial disturbance grows rapidly, but has small amplitude.Thus, appears to be similar to that of a helix with a higher azimuthal wavenumber than mode 1.
v ′ 2 and w ′ 2 appears to be symmetric w.r.t. the base streaky structure.In the DNS, the flow structure is likely to be a superposition of the three unstable modes, although the helical mode dominates over a long time (figure 1) due to its much higher growth rate.

Figure 1 (FIG. 2 :
Figure1(a) depicts the primary disturbance field at three time instants, from the APG simulation with T u = 1%.Instantaneous iso-surfaces of streamwise perturbation velocity, u ′

Figure 2 (
a) shows the full spanwise extent of the domain; a zoomed-in view of the precursor on which the helical instability develops is presented in figure 2(b).

Figure 2 (
Figure 2(b) includes contours of the streamwise perturbation eigenfunction corresponding to the eigenmode with the highest growth rate computed via (2) for this base flow configuration.k x is an input parameter of the stability calculation.It is estimated from the helical instability in the snapshot shown in figure 1(c).It is evident that the instability develops

FIG. 6 :
FIG. 6: Evolution of a helical instability in ZPG (T u = 1% and A = 0.5%U ∞ ): (a) Isosurfaces of u ′ = ±0.1 at t 0 and (b) at t 0 + 10.8.(c) Contours of −0.2 ≤ u ′ ≤ 0.2 (dark for negative values) at time t 0 and x/δ 0 ∼ 160.The yellow contours show the real part of u ′ 2 for the most unstable eigenmode (dashed lines indicate negative values).The grey solid line and the horizontal yellow dashed line mark the local critical layer and δ 99 .(d) Top view of the iso-surfaces of u ′ 2 = ±0.025(red and green represent positive and negative perturbation, respectively) for this mode plotted over two wavelengths.
streaks and instability waves.Stability analyses were performed for base states extracted from DNS of mixed mode transition in both APG and ZPG boundary layers.Those base states are unstable to modes of a helical form.Why these base states occur in the mixed mode regime remains an open question.

FIG. 7 :FIG. 9 :
FIG. 7: Contours of real part of the streamwise perturbation eigenfunction (dashed lines indicate negative values) for the (a) second and (b) third most unstable modes at time t 0 + 16.2 as in figure 2(b) at the streamwise location x/δ 0 ∼ 185.8.The grey solid line and the horizontal yellow dashed line mark the local critical layer and δ 99 , respectively.