2.1. Experimental set-up and instrumentation

The experiments were conducted in the Arizona Low Speed Wind Tunnel (ALSWT) situated in the Department of Aerospace and Mechanical Engineering at the University of Arizona. The closed-loop wind tunnel has a test section of 0.9 m $\times$ 1.2 m $\times$ 3.65 m (height x width x length, 3 ft $\times$ 4 ft $\times$ 12 ft). A Pitot tube mounted 0.4 m downstream of the test section entry at the tunnel sidewall extends into the free stream to acquire total and static pressures to measure the flow speed. The maximum flow speed in the test section is $U_{\infty } \approx 80\,\mathrm {m\,s^-{^1}}$ (262.5 ft s^–1). Uniformity of the mean flow over the test section is at or better than $\pm 0.5\,\%$ . The turbulence intensity is less than Tu = 0.05 % in the range of 1 Hz to 10 kHz and considerably lower (e.g. Tu = 0.035 %) for the flow speed of the experiments discussed in this paper (7 m s^–1, see § 2.2). A more extensive investigation of the free-stream properties can be found in Borgmann et al. (Reference Borgmann, Hosseinverdi, Little and Fasel2020). The temperature inside the tunnel is regulated by a heat exchanger with a chilled water supply. Throughout the experiments, temperature was held within the range of $\pm$ 0.44 $^{\circ }$ C (1 $^{\circ }$ F) of 22.2 $^{\circ }$ C (72 $^{\circ }$ F).

A LSB is generated on a flat plate by an inverted wing: a NACA $64_3-618$ airfoil with a chord length of $c = {8}\,\mathrm {in.}$ (203.2 mm) shown in figure 2. The favourable-to-adverse pressure gradient along the suction side of the inverted wing leads to the formation of a LSB for a wide range of Reynolds numbers. Projection of the pressure gradient onto a flat plate neglects the influence of surface curvature to the transition process, analogous to early work by Gaster (Reference Gaster1967). In addition to the flexibility in the mounting location and incidence angle of the inverted wing, the relatively small streamwise extent and high aspect ratio provide a quasi-2-D LSB across a large spanwise extent (more than four times the time-averaged separation length). In addition, the low total blockage in the wind tunnel provides a swift pressure recovery downstream of the LSB. This contrasts with the large extent of APG generally observed downstream of reattachment when using a contoured ceiling (Watmuff Reference Watmuff1999; Diwan & Ramesh Reference Diwan and Ramesh2009). Varying the distance and the incidence angle of the wing relative to the flat plate allows for a wide range of possible free-stream pressure distributions acting on the flat-plate boundary layer. For the investigation discussed here the inverted wing is placed at a distance of 3.2 in. (81.3 mm) from the plate surface; 17.5 in. (0.44 m) downstream of the leading edge; at an incidence angle of ${2}^{\circ }$ (relative to the inverted airfoil, see drawing in figure 2). The wing is equipped with boundary-layer suction along $50\,\%{-}60\,\%$ -chord and $85\,\%{-}90\,\%$ -chord, on the suction side and a zigzag turbulator tape to trip the flow close to the leading edge at the pressure side to avoid separation on either surface of the inverted wing. Careful distribution of suction holes in a PVC pipe inside the plenum of the hollow downstream half of the inverted wing ensure even suction across the span. Suction is provided by LAMB AMETEK 117500-12 3-STAGE vacuum motor, connected to both ends of the inverted wing.

Figure 2.

Flat-plate model and inverted airfoil used for LSB transition experiments: (a) isometric view showing coordinate system (red), pressure tap locations (orange) and suction direction (blue); (b) cross-section of airfoil showing PVC pipe outline, plenum in the hollow wing and holes along the suction side; planform view shows suction distribution along the airfoil suction side where red lines highlight the active sections; (c) side view, highlighting the incidence angle between the inverted wing and flat plate (blue text); Coloured planes in (a) and (c) mark the PIV planes in the x-y (green) and x-z planes (orange); The origin of the coordinate system shown in (a) and (c) places $z = 0$ at the centre of the span and $y = 0$ at the flat-plate surface.

The flat-plate model shown in figure 2 is made from an aluminium honeycomb structure, which is covered on both sides by a 1 mm aluminium skin. The width of the plate is 47.75 in. (1.21 m), spanning nearly the entire width of the test section. Gaps between the sides of the plate and the walls of the test section were sealed with a soft expanding foam. The plate is 78 in. (1.98 m) long with a total thickness of 0.5 in. (13 mm) and a flatness of $\pm {0.008}\,\rm {in\, ft^-{^1}}$ ( $\pm {0.665}\,\rm {mm\, m^-{^1}}$ ). Streamlined support structures elevate the plate to 8.25 in. (210 mm) above the wind-tunnel floor to avoid wall effects. The interchangeable leading edge is configured as a 1 : 20 super-ellipse, see Lin, Reed & Saric (Reference Lin, Reed and Saric1992). Its total length of 9.5 in. (241 mm) and smooth transition to the plate top surface combined with the adjustable flap at the trailing edge of the plate ensure consistent and repeatable inflow conditions, as the flap compensates for blockage and allows for positioning of the stagnation point along the leading edge. The flow along the bottom side of the plate is tripped just downstream of the leading edge to reduce unsteadiness in the flow underneath the plate. The origin of the coordinate system used throughout the investigation is located at the centre of the leading edge, with the x axis in the streamwise, y axis in the vertical and z axis in the spanwise direction. The plate surface is considered zero in the vertical (y) direction. All results in this work are presented in dimensionless form where the length scales, $x$ , $y$ and $z$ , are scaled with a reference length of $L_\infty ^\ast = {0.0254}\,\mathrm {m}$ ( ${1}\,\mathrm {inch}$ ), in accordance with the earlier experiments by Gaster (Reference Gaster1967) and velocities are scaled with the free-stream velocity $U_\infty ^\ast$ , respectively

Pressure taps are located along two parallel lines in the streamwise direction, 1.5 in. (38 mm) offset to either sides of the centreline. The streamwise spacing is 12.5 mm (0.5 in.) in the vicinity of the LSB and 25 mm (1 in.) elsewhere. The spacing is gradually refined towards the leading edge with the first tap at 15 mm (0.6 in.) and an initial spacing of 5 mm 0.2 in.). Scanivalve Corp. ZOC33 pressure scanners (full scale range [FD] and accuracy of 2.49 k Pa (10 in. $\text {H}_2\text {O}$ ) $\pm 0.15\, \%$ FS) in combination with an ERAD Remote A/D module are used to record surface pressure. Pressure data are sampled at 625 Hz and averaged over 96 s.

Time-resolved velocity measurements were recorded using constant temperature anemometry (CTA). A Dantec Dynamics StreamlinePro frame was operated with a 55P11 miniature wire probe and a 55H21 support. The hot-wire probe was attached to a two axis traverse along the wind-tunnel ceiling providing variable positioning of the sensor along the x-y plane. The CTA probe was calibrated using a StreamLine Pro Automatic Calibrator and referenced to a temperature probe positioned above the CTA in the free stream on the traverse.

The CTA data in the free stream and along the LSB were recorded at 40 kHz for 20 s and analogue filtered between 1 Hz and 10 kHz. Spectral content was calculated using Welch’s method, with a resolution of 0.61 Hz averaged over 21 windows. Free-stream turbulence intensity was calculated as the root-mean-square of the velocity fluctuations scaled by the local free-stream velocity (2.1) over the same frequency range (1 Hz–10 kHz).

(2.1) \begin{equation} Tu = \frac {\sqrt {|u'|^2}}{U_\infty } .\end{equation}

In addition to hot-wire data, a LaVision GmbH particle image velocimetry (PIV) system was used to obtain spatially resolved velocity data in x-y planes of the wake of the inverted wing and in the LSB (PIV plane shown in figure 2). A Quantel Evergreen HP dual-head Nd:YAG laser illuminates submicron di-ethyl-hexyl-sebacat (DEHS) seed particles. In the wake of the inverted wing, 300 image pairs were acquired at 10 Hz with a 5-megapixel 16-bit Imager sCMOS camera, incorporating a 2x teleconverter lens, 50 mm lenses and narrow band-pass optical filters. For each preprocessed and time-filtered image pair, subregions were cross-correlated using decreasing window sizes ( $64^2$ – $24^2$ pixel) and multipass processed with $50\, \%$ overlap followed by a median filter and smoothing function. To increase the resolution in the LSB, two cameras with the same optics set-up were placed in parallel. The data in the x-y plane along the LSB consist of 1500 image pairs, recorded at 14 Hz. Images were preprocessed and time filtered, and subregions were cross-correlated with window sizes decreasing from $64^2$ – $16^2$ pixel with $50\, \%$ overlap. The final images were stitched and smoothed by a single pass with a median filter (size 3 × 3). The error based on a 95 % confidence interval for the velocity magnitude relative to the free-stream velocity (7 m s^–1) is ${\lt}3\,\%$ in the wake of the inverted wing and ${\lt}2\,\%$ in the separation bubble. The larger number of images recorded for the time-averaged data in the separation bubble is motivated by convergence of second-order statistical quantities (e.g. Reynolds stress, turbulent kinetic energy, etc.). Assuming a normal distribution of all statistical quantities the relative error in the turbulent kinetic energy based on a 95 % confidence interval is ${\lt}5.5\,\%$ in the LSB.

Stereoscopic PIV was used in the x-z plane of the LSB. Two 5-megapixel 16-bit Imager sCMOS cameras incorporating scheimpflug mounts, 85 mm lenses and narrow band-pass optical filters are used to record 800 image pairs at 14 Hz. Using the same processing routine described for the x-y plane data, image pairs are processed separately, using a final window size of $32^2$ pixel. The resulting velocity fields are further postprocessed to remove vectors with correlation peak ratios of less than 1.5 and correlation coefficients below 0.1. After applying the stereo calibration, a 3 $\times$ 3 Gaussian smoothing algorithm and polynomial filters were employed for the final data. The time-averaged results have less than 3.5 % statistical uncertainty (based on free stream) using a 95 % confidence interval. The spatial resolution is 1.6 mm (0.064 in.).

2.2. Free-stream conditions

The free-stream flow in the test section with and without the model (airfoil and flat plate) installed was measured upstream of the location of the inverted wing (when installed), at the mid point (vertically) between the flat plate and tunnel ceiling in both cases ( $x = 9$ , $y = 13$ , $z = 0$ ). The tunnel operating conditions were identical with respect to fan frequency and fan blade angle, resulting in a slight increase in free-stream velocity ( ${\approx}4\, \%$ ) when the model was installed due to the blockage effect of the airfoil and flat plate. This generates a slightly elevated noise floor in the PSD with the model installed. Prior to mounting the model in the wind-tunnel test section, assessment of the free-stream of the empty test section showed low levels of turbulence intensity (Tu $\leqslant 0.034\,\%$ , in the range of 1 Hz to 10 kHz). Data were collected at a Reynolds number of $Re_C = 90\,000$ , based on the chord of the inverted wing ( $c = {203.2}\,\mathrm {mm}$ ). Major contributions to the turbulence intensity are observed in the low frequency range ( $\lt {5}\,\mathrm {Hz}$ ), resulting in Tu $\leqslant 0.015\,\%$ if the signal is band pass filtered between 5 Hz to 10 kHz). Therefore, the PSD of the streamwise velocity is shown for reference in figure 3. With the flat-plate model and inverted wing installed, a minimal increase to Tu $\leqslant 0.035\, \%$ is observed in the range of 1 Hz to 10 kHz (or Tu $\leqslant 0.017\, \%$ in the range of 5 Hz to 10 kHz), while spectral content remains nearly unchanged. Peaks in the spectrum at 36 Hz and 72 Hz can be attributed to the fan blade passing frequency and its harmonic.

Figure 3.

Power spectral density (PSD) of the free-stream velocity in the empty test section and at test conditions for the LSB with the model installed, $Re_C = 90\,000$ .

2.3. Flow around the inverted wing

At low chord based Reynolds number ( $Re_C = 90\,000$ ) of the inverted wing and low incidence angles ( ${2}^{\circ }$ ) between the airfoil and the flat plate (figure 2), laminar separation was observed at the suction (bottom) and pressure (top) side of the NACA $64_3{-}618$ airfoil. Velocity magnitude contours in figure 4 show the necessity of flow control to avoid separation leading to unsteadiness in the free stream. Without control the flow separates on the wing just downstream of the maximum thickness at $x = 21.45$ ( $35\, \% c$ ) on the suction side facing the flat plate, figure 4(a). The boundary layer does not reattach to the wing and causes a significant wake. In addition, a separation bubble forms on the pressure side between $x = 22.05$ and $x = 24.8$ (figures 4 a and 4 c). The same was observed in precursor computational fluid dynamics (CFD) calculations shown in figures 4 b and 4 d. Steady 2-D boundary-layer suction through the mechanism described in § 2.1 was successful in preventing separation on the suction side. The method is different from the jet entrainment scheme at the trailing edge of the inverted wing used in the seminal work by Gaster (Reference Gaster1967, Reference Gaster2006). Successful application of quasi-continuous boundary-layer suction (plenum suction behind a perforated section of the wing) in a similar set-up in a water tunnel experiment was previously shown by Radi & Fasel (Reference Radi and Fasel2010) and Jagadeesh & Fasel (Reference Jagadeesh and Fasel2013). Figures 4(c) and 4(d) demonstrate the effectiveness of suction to suppress boundary-layer separation as the flow remained attached along the entire chord.

Figure 4.

Velocity magnitude $|U|$ along the inverted wing in the experiment (a, c, e) and from the base flow calculations in ANSYS Fluent withthe $\gamma -Re_\theta -SST$ transition model (b) and (d); (a) experiment with no boundary-layer suction and no trip on pressure side; (b) ANSYS fluent without suction; (c) experiment with suction, but no trip on pressure side; (d) ANSYS Fluent calculation with suction; (e) experiment with suction and pressure side tripped; (f) PSDs from measurement in the wake of the pressure side with and without trip ( $x = 25.5$ ; $y = 4.75$ ).

The small separation bubble at the pressure side had negligible influence on the time-averaged flow field between the wing and flat plate. Spectral analysis of CTA measurements taken above the trailing edge downstream of the separated region reveals two frequency peaks between 200 and 300 Hz (figure 4 f). To avoid interference between these periodic disturbances emanating from the pressure side of the inverted wing and the LSB on the flat plate in the experiment, the boundary layer on the pressure side was tripped just downstream of the leading edge. Then the turbulent boundary layer no longer separated, the frequency peaks in the power spectra were suppressed (figure 4 f) and the separation on the pressure side was no longer present in the time-averaged flow field along the inverted wing (figure 4 e). With boundary-layer tripping (pressure side) and suction (suction side), the flow remained attached along the entire chord on both sides of the wing (figure 4 e).

2.4. Laminar separation bubble on the flat plate

The pressure-gradient-induced LSB at the flat plate was recorded for the flow field shown in figure 4(e), for which the application of flow control to fully attach the flow along the entire inverted wing. Pressure measurements on the flat plate in figure 5 confirm the effects of suction at the wing to the surface pressure at the flat plate. In figure 5 the downstream $C_p$ development corresponds to figures 4(a) and 4(e), without and with flow control (suction and trip). As expected, the separation without flow control at the wing (figure 4 a) result in a weak favourable-to-adverse pressure gradient in the time-averaged sense. With suction, a significant plateau of constant $C_p$ between $x = 22.5$ and $x = 26$ indicates the presence of a LSB. The effect of boundary-layer tripping at the pressure side of the inverted wing was found to be negligible in the time-averaged pressure measurements (not shown here). The current experimental data compare very well with the pressure gradients in the experiments by Gaster (Reference Gaster1967) (Series I, Case 6 in figure 5) for a similar free-stream velocity 21.8 ft s^–1 (6.64 m s^–1)) compare well with the current experiment. To compare the pressure gradients of the current work with results from Gaster (Reference Gaster1967), the boundary layer along the flat plate was tripped, thus resulting in a quasi-inviscid pressure development. The resulting turbulent boundary layer stays attached and surface pressure measurements resemble the inviscid case (red and grey line in figure 5). In both cases (laminar and turbulent boundary layer along the flat plate), the NACA $64_3-618$ causes a more gradual favourable pressure gradient compared with the inverted wing used by Gaster (Reference Gaster1967). However, the APG and separation location match remarkably well, with a slightly larger LSB in the current study. Results from Gaster (Reference Gaster1967) (Series I, Case 6 in figure 5) were scaled in the streamwise direction to compensate for the smaller chord length of $c_{Gaster} = {139.7}\,\mathrm {mm}$ (5.5 in., chord-based Reynolds number of $63\,000$ ). The boundary-layer momentum thickness at separation was lower in Gaster’s experiments, $Re_{s,\theta }=218$ (Series I, Case 6 – 21.8 ft s^–1 (6.64 m s^–1)).

In figure 5(b) spanwise pressure measurements at the flat plate are provided for select downstream locations, for the most unsteady regions close to separation and reattachment of the LSB. The results in figure 5(b) indicate small variations support the assumption of essentially 2-D conditions in a time-averaged sense.

Figure 5.

Flat-plate surface $C_p$ development in the experiment. (a) Chordwise $C_p$ with and without boundary-layer suction at the inverted wing, pressure side tripped in both cases. Comparison with Series I, Case 6 (Gaster Reference Gaster1967), for quasi-inviscid (tripped flat-plate boundary layer). (b) Spanwise $C_p$ at select chordwise locations for suction and tripping at the inverted wing.

An overview of time-averaged quantities for the LSB is provided in figure 6. Normalized velocities ( $u/U_{\infty }$ and $v/U_{\infty }$ ) are shown next to the root-mean-square (r.m.s.) of the respective fluctuating quantities. In the region of APG, the incoming laminar boundary layer separates at $x=21.8$ . Mean reattachment was found at $x=27.6$ resulting in a bubble length of $L=5.8$ and a maximum bubble height of $H=0.35$ at $x = 26$ . Separation and reattachment are determined based on the dividing streamline extracted from the time-averaged flow field (figure 1). Velocity vectors along the wall ( $y\lt 0.04$ ) were omitted to avoid the impact of laser reflections from PIV. Vector fields were interpolated to ensure no-slip conditions at the wall. The dividing streamline is identified as the streamline between the time-resolved separation and reattachment location, closing the separated region (Kurelek et al. Reference Kurelek, Tuna, Yarusevych and Kotsonis2020). The fluctuating quantities in figure 6 were obtained from the r.m.s. of the $u'$ and $v'$ -velocity fluctuations. Low levels of streamwise velocity fluctuations near the separation location and along the dividing streamline upstream of maximum thickness agree with results in Zaman, Mckinzie & Rumsey (Reference Zaman, Mckinzie and Rumsey1989); Michelis & Kotsonis (Reference Michelis and Kotsonis2016) and appear to be related to low frequency ’flapping’, which was observed in separation bubbles at similar Reynolds numbers. Significant levels of streamwise and wall-normal fluctuations near the maximum bubble height and downstream of it are related to formation of strong periodic vortex shedding, as a result of the inviscid instability in the shear layer developing between the free stream and the region of reverse flow. Based on the bubble topology, the two parameter bursting criterion proposed by Gaster (Reference Gaster1967) was calculated as $\overline {P} = -0.206$ , based on the momentum thickness at separation $\theta = 0.024$ (0.613 mm), the dynamic viscosity $\nu = 1.658 \times 10^{-5}$ , the velocity decrease across the separation ( $\triangle U = {1.335}\,\mathrm {m\,s^-{^1}}$ ) and the separation length ( $L_s$ ). The resulting value is indicative of a short LSB and is almost identical to the LSB in Gaster (Reference Gaster1967) Series I, Case 6 as shown in figure 5 ( $\overline {P}_{Gaster,I,6} = -0.199$ ) at similar free-stream conditions.

Figure 6.

Time-averaged velocity components (a) $U/U_{\infty }$ and (b) $V/U_{\infty }$ , r.m.s. of velocity fluctuations (c)  $u'_{rms}/U_{\infty }$ and (d) $v'_{rms}/U_{\infty }$ across the separation bubble derived from PIV. Streamline overlay for reference. Dividing streamline (dashed line), separation ( $\blacktriangleleft$ ) and reattachment location ( $\blacktriangle$ ).

3.1. Set-up for stability investigations and DNS

The numerical investigations, consisting of stability analysis and DNS, were carried out to complement the experimental work in order to extract and better explain the relevant physics. The set-up was guided by the wind-tunnel experiments, as discussed in § 2. The computational approach involves two steps (see figure 7 a): in step $I$ , a ‘precursor calculation’ of the entire flow field based on the set-up of the wind-tunnel experiments including the wind-tunnel walls, inverted wing and the top of the flat plate was performed using ANSYS Fluent. The ‘precursor calculation’ uses a transitional $\gamma -Re_{\theta }-SST$ (shear stress transport) model and a slip wall along the suction side of the inverted wing. In step $II$ , high resolution DNS of the flow along the flat plate downstream of the leading edge was performed. Inflow and free-stream boundary conditions for the DNS are extracted from the flow field computed in the precursor calculations.

Figure 7.

Schematic of the computational set-up for precursor calculations, step $I$ , according to ALSWT experiments (a). Also shown are contours of streamwise velocity. The area indicated by dashed lines corresponds to the computational domain for DNS, step $II$ (b).

The computational domain for the presented stability analysis and DNS was defined as $10\leqslant x\leqslant 36$ (figure 7 b), the domain height is $l_y=2.8$ and the domain width in the spanwise direction is $l_z=2$ . At the outflow boundary, $x=x_o$ , all second derivatives in the streamwise direction are set to zero. In addition, a buffer domain in region $x_b\lt x\lt x_o$ smoothly dampens out the fluctuations generated inside the domain and prevents upstream contamination (figure 7 b). No-slip and no-penetration conditions are enforced on the surface of the flat plate.

For all DNS results presented in this paper, the same computational grid with $n_x\times n_y\times n_z=2601\times 200\times 220\approx 115$ million grid points was used. The grid spacing is uniform in the streamwise direction, and the grid points are exponentially displaced from the wall in the wall-normal direction to improve the resolution near the wall. In the spanwise direction, 139 Fourier modes (220 collocation points) were employed for the domain width of $l_z=2$ . The spanwise domain width was selected based on the maximum bubble height (Jones, Sandberg & Sandham Reference Jones, Sandberg and Sandham2008; Balzer & Fasel Reference Balzer and Fasel2016; Hosseinverdi & Fasel Reference Hosseinverdi and Fasel2019). It was found that, to resolve the largest spanwise structures, a domain width of at least four times the maximum bubble height is necessary. The ratio of domain width to the maximum bubble height was found to be approximately six. The grid resolution in wall units within the separation bubble and in the redeveloping turbulent boundary layer downstream of the reattachment point was then $\triangle x^+\leqslant 5.8$ , $\triangle z^+\leqslant 5.2$ and $\triangle y^+_w\leqslant 0.5$ .

Direct numerical simulation of transition requires numerical methods with low numerical dispersion and dissipation errors. An extensively validated high-order accurate Navier–Stokes solver developed in the University of Arizona CFD Laboratory was employed for the present numerical simulations (Meitz & Fasel Reference Meitz and Fasel2000; Hosseinverdi, Balzer & Fasel Reference Hosseinverdi, Balzer and Fasel2012; Hosseinverdi & Fasel Reference Hosseinverdi and Fasel2017b ). In this code, the 3-D unsteady incompressible Navier–Stokes equations are solved in vorticityvelocity formulation

(3.1) \begin{align} \frac {\partial \boldsymbol {\omega }}{\partial t} &=\boldsymbol {\nabla }\times (\boldsymbol {u}\times \boldsymbol {\omega })+\frac {1}{Re}\triangle \boldsymbol {\omega }, \end{align}

(3.2) \begin{align} \triangle \boldsymbol {u}&=\boldsymbol {\nabla }\times \boldsymbol {\omega }, \end{align}

where the vorticity vector, $\boldsymbol {\omega }$ , is defined as the negative curl of the velocity vector $\boldsymbol {\omega }=-\boldsymbol {\nabla }\times \boldsymbol {u}$ . In (3.1), the global Reynolds number is defined as $Re=U_\infty ^\ast L_\infty ^\ast /\nu ^\ast$ , where $\nu ^\ast$ is the kinematic viscosity. The asterisk denotes dimensional quantities. Here, $U_\infty ^\ast$ and $L_\infty ^\ast$ are reference velocity and length scales, respectively.

The governing equations are solved in a 3-D Cartesian coordinate system, where the streamwise, wall-normal and spanwise directions are denoted by $x$ , $y$ and $z$ , respectively (figure 7), and the corresponding velocity and vorticity components, ( $u$ , $v$ , $w$ ) and ( $\omega _x$ , $\omega _y$ , $\omega _z$ ). The governing equations are integrated in time using an explicit fourth-order RungeKutta scheme. All derivatives in streamwise and wall-normal directions are approximated with fourth-order accurate compact finite differences. Note that the finite-difference approximations for the derivatives with respect to $y$ are constructed for a non-equidistant grid instead of using a coordinate transformation. The flow field is assumed to be periodic in the spanwise direction, therefore, it can be expanded in Fourier cosine and sine series with a pseudo-spectral treatment of the nonlinear terms. The velocity Poisson equations are solved by a combination of a fourth-order standard compact difference in the wall-normal direction and Fourier sine transform in the streamwise direction. For details of the numerical method see Meitz & Fasel (Reference Meitz and Fasel2000), Hosseinverdi et al. (Reference Hosseinverdi, Balzer and Fasel2012)and Hosseinverdi & Fasel (Reference Hosseinverdi and Fasel2017b ).

3.2. Direct numerical simulations without free-stream turbulence

Contrary to the wind-tunnel environment, without the addition of FST and other artificial disturbances in the computational domain, the low numerical noise of the high-fidelity DNS leads to quasi-zero FST conditions in the calculations, providing ideal conditions to investigate the natural transition in the LSB. Time-averaged velocity profiles from the DNS without FST in figure 8(a) show the streamwise component in the boundary layer along the flat plate. The results are in good agreement with those extracted from PIV data in the experiment (figure 8 a), suggesting near identical boundary-layer development upstream of laminar separation between experiment and numerical calculations in a time-averaged sense, independent of the difference in FST. As expected, the APG decelerates the flow in the boundary layer downstream of the apex of the inverted wing ( $x = 20.42$ ), resulting in an inflection point in the boundary layer in the vicinity of laminar separation. Experimental data outside the symbols in figure 8(a) were interpolated, enforcing no-slip conditions at the wall.

Figure 8.

(a) Normalized boundary-layer profiles between the onset of APG and laminar separation, extracted from PIV and DNS. (b) Instantaneous flow visualizations for 3-D DNS for $Tu=0\, \%$ . Shown are contours of spanwise vorticity (averaged in spanwise direction) together with the mean dividing streamline identified by black dashed lines.

Instantaneous contours of spanwise vorticity (averaged in the $z$ -direction) obtained from 3-D DNS without FST are presented in figure 8(b). The location and size of the bubbles were obtained from the time- and spanwise-averaged data to extract the dividing streamlines in figure 8(b). The inflection of the velocity profile in the shear layer between free-stream and recirculating flow causes growth of disturbances upstream of the maximum bubble height due to the K–H instability. Subsequent ‘roll-up’ and shedding of vortical structures were observed near the location of the maximum bubble height. The shear-layer roll-up produces large wall-normal momentum which facilitates entrainment and aids in the reattachment process. Small vortical structures that are generated near the wall due to the breakdown of these large-scale structures further contribute to the reattachment process. It should be noted that transition to turbulence appears to be self-sustained for the DNS even without FST, i.e. no external forcing or disturbances were required to initiate or sustain the transition process. This observation agrees with the bounded range of globally unstable spanwise wavenumbers identified in § 3.3. This is in contrast to zero-pressure-gradient boundary layers and implies that the development of 3-D disturbances is a direct result of an absolute/ global instability mechanism.

Comparison of the mean flow via streamlines and time-averaged surface pressure coefficient shows the importance of FST for the topology of the LSB (9). The development of $C_p$ calculated using the transitional $\gamma -Re_{\theta }-SST$ model in the Reynolds-averaged Navier–Stokes (RANS) calculations which is known to struggle with LSBs in low Reynolds number flows, predicts – as expected – a shorter bubble than the experiments. The high-fidelity DNS is in good agreement with the experiment for $x\leqslant 25.5$ . However, streamlines and surface pressure indicate that reattachment occurs farther downstream compared with the experiments. The presence of minute ‘numerical noise’ in the DNS (with no added FST) reduces the amplitudes of the disturbances in the upstream boundary layer. The dominant 2-D mode in the shear layer requires extended growth downstream before reaching critical amplitudes (compared with the case with FST), necessary for the interaction with the global mode in the LSB to cause transition to turbulence. This leads to a later transition along the shear layer and results in significant deviation from the wind-tunnel experiments for $x \geqslant 25.5$ , despite the low levels of FST observed in the wind-tunnel experiments (§ 2.2). The influence of FST in the transition process, especially within the LSB, is well documented (Jacobs & Durbin Reference Jacobs and Durbin2001; Balzer & Fasel Reference Balzer and Fasel2016; Hosseinverdi & Fasel Reference Hosseinverdi and Fasel2019; Simoni et al. Reference Simoni, Lengani, Ubaldi, Zunino and Dellacasagrande2017; Istvan & Yarusevych Reference Istvan and Yarusevych2018) and results in this section corroborate that meaningful comparison of numerical results with experimental efforts will have to account for the disturbances present in the respective wind-tunnel or flight environment (§§ 3.4 and 4).

Figure 9.

Comparison of mean results between experiments and DNS. Plotted in (a) are streamlines obtained from experiments and DNS without FST. (b) Comparison of flat-plate surface pressure coefficient, $C_p$ , between experiment, DNS without FST as well as RANS calculations.

3.3. Linear local and global stability analyses

In addition to DNS, local and global linear stability analyses were carried out to investigate the initial stages of the transition process, and to help interpret observations from the experiments and DNS (§ 4) such as for example the self-sustained transition in the DNS without FST as discussed in the previous section. First, linear stability theory (LST) analysis was performed based on the mean flow field (time- and spanwise averaged) obtained from DNS without FST (§ 3.2) in order to map out the stability behaviour with respect to 2-D disturbance waves. The stability diagram from the solution to the OrrSommerfeld equation for the mean streamwise velocity profiles (see figure 8 a) is presented in figure 10, which shows contours of the spatial growth rates, $\alpha _i$ , in the $St$ -x plane. Since the spatial stability problem was considered ( $\alpha$ -complex), the local velocity profiles are either convectively stable ( $\alpha _i\gt 0$ ) or unstable ( $\alpha _i\lt 0$ ). The solid black contour line corresponds to the neutral curve ( $\alpha _i=0$ ). LST results reveal that the LSB is unstable with respect to disturbances in a broad range of frequencies.

Figure 10.

Linear stability theory result for the time- and spanwise-averaged $u$ -velocity of the DNS without FST. Shown are the contours of the amplification rate, $\alpha _i$ , together with contour lines of $N$ -factor for 2-D distributed waves in the x- $St$ plane. The solid contour line represents the neutral curve. The dominant shedding frequency observed in the experiment ( $St_{EXP}$ ) and DNS ( $St_{DNS}$ ) are indicated by horizontal dotted lines. Here, $x_{apg}$ and $x_s$ correspond to the $x$ -location of the onset of APG and separation point, respectively.

For the linear global stability analysis (LGSA), the linearized Navier–Stokes equations (LNSE) with a mean flow obtained from DNS were directly solved as an initial value problem

(3.3a) \begin{align} \frac {\partial \boldsymbol {\omega }^{\prime }}{\partial t}&={\nabla }\times \left (\boldsymbol {u}^{\prime }\times \boldsymbol {\Omega } + \boldsymbol {U}\times \boldsymbol {\omega }^{\prime }\right )+\frac {1}{Re}\triangle \boldsymbol {\omega }^{\prime }, \end{align}

(3.3b) \begin{align} \triangle \boldsymbol {u}^{\prime }&={\nabla }\times \boldsymbol {\omega }^{\prime }, \end{align}

where $\boldsymbol {U}(x,y)=(U,\,V,\,0)^T$ and $\boldsymbol {\Omega }(x,y)=(0,\,0,\,\Omega _z)^T$ correspond to the the steady baseflow, and $\boldsymbol {u}^\prime (\boldsymbol {x},t)$ and $\boldsymbol {\omega }^\prime (\boldsymbol {x},t)$ denote the disturbance flow. This framework allows us to capture temporally decaying/growing global modes in the context of convective/absolute instability mechanisms. Note that this approach fully takes into account all the non-parallel effects with respect to the baseflow and the form of the disturbances.

Since the base flow is only a function of $x$ and $y$ and homogeneous in the third direction, 3-D (oblique) disturbance waves are further decomposed in the $z$ -direction as

(3.4a) \begin{align} \left [u^\prime ,v^\prime , \omega _z^\prime \right ]^{T}(\boldsymbol {x},t)&=\left [ \hat {u}, \hat {v}, \hat {\omega }_{z}\right ]^T (x,y,t)\cos (\gamma z)\,, \end{align}

(3.4b) \begin{align} \left [w^\prime ,\omega _x^\prime , \omega _y^\prime \right ]^{T}(\boldsymbol {x},t)&=\left [ \hat {w}, \hat {\omega }_{x}, \hat {\omega }_{y}\right ]^T (x,y,t)\sin (\gamma z)\,, \end{align}

where $\gamma$ is the spanwise wavenumber. In the LGSA, disturbances with a specified $\gamma$ were perturbed through a narrow blowing and suction slot on the flat-plate surface by pulsing the wall-normal velocity of the following form:

(3.5) \begin{equation} \hat {v}_k(|x-x_f|\leqslant b/2,y=0,t)= S(x) \exp \left ( \left [\frac {t-t_p}{\delta t_p}\right ]^2 \right ), \end{equation}

where $x_f$ and $b$ correspond to the centre and streamwise width of the slot, respectively, and are set to $x_f=19$ and $b=0.2$ . The shape function, $S(x)$ , is a polynomial which is zero outside the forcing slot such that smooth derivatives at boundaries of the suction/blowing slot are obtained (Hosseinverdi & Fasel Reference Hosseinverdi and Fasel2019). In (3.5), the temporal function is realized by a Gaussian distribution, where $t_p=10 \triangle t$ and $\delta t_p=0.25\triangle t$ , for which the corresponding frequency spectrum is broadband with approximately the same amplitude levels for a frequency range up to $St = f L_\infty ^\ast / U_\infty = 100$ .

For a specified $\gamma$ , the temporal growth rate of the least stable mode (or the most unstable one in case of an absolute instability) will determine the long-time behaviour of the flow with respect to 2-D/3-D disturbances. Figure 11(a) shows the temporal evolution of the maximum disturbance kinetic energy, $\hat {q}_{k}$ , for a few selected values of $\gamma$ for the mean flow obtained from the DNS with zero FST, where $\hat {q}_{k}$ is computed as $\hat {q}^2_k= (\hat {u}^2 + \hat {v}^2 + \hat {w}^2 )$ . The linear response of the mean flow to the wave packet disturbances can be classified into two regions: (i) a short-time behaviour for $0\lt t\lt 30$ , which is associated with the convective instability characteristics of the base flow (cf. LST) and (ii) a long-time behaviour which indicates the global stability of the mean bubble, i.e. being globally stable or unstable.

Figure 11.

(a) Time evolution of maximum disturbance $\hat {q}_k$ in the logarithmic scale for DNS without FST. (b) Temporal growth rate as a function of spanwise wavenumber.

In figure 11(a), the decay/growth of the maximum $\hat {q}_k$ for 2-D/3-D modes indicate that the mean LSB obtained from the DNS without FST is globally stable/unstable with respect to the 2-D/3-D disturbances. Analysis for 2-D and 3-D disturbances in the range of $0\leqslant \gamma \leqslant 10\pi$ for the base flow from the DNS and the temporal growth rate associated with each $\gamma$ is then computed based on $\tilde {\sigma }=\log (\hat {q}_{k,max})/\triangle t$ . The distribution of the temporal growth rate as a function of $\gamma$ is displayed in figure 11(b). It can be observed that the LSB is asymptotically stable with regard to $\gamma =0$ (2-D mode) and becomes globally unstable for a range of $\gamma$ , i.e. $\pi \leqslant \gamma \leqslant 9\pi$ , where the maximum $\tilde {\sigma }$ is attained for $\gamma =3\pi$ . The presence of this global instability explains the self-sustaining transition to turbulence in the LSB, consistent with findings by Rodríguez et al. (Reference Rodríguez, Gennaro and Souza2021).

Subsequently, using the time-dependent data available from the LNSE, dynamic mode decomposition (DMD) was employed to extract the relevant unstable modes. When a discrete dataset of the flow field is taken at $N$ snapshots separated by a constant sampling time interval $\triangle t$ , DMD allows for an expansion of the velocity data in the form

(3.6) \begin{equation} \boldsymbol {u}(\boldsymbol {x},t)=\sum _{j=1}^{N} \boldsymbol {\psi }_j(\boldsymbol {x}) \exp {\left ([\sigma _j +i\omega _j] t\right )} , \end{equation}

where $\boldsymbol {\psi }_j$ are the spatial DMD modes and $\sigma _j$ and $\omega _j$ represent the growth rates and oscillatory radial frequencies, respectively.

For example in figure 12(a) the DMD eigenvalue spectrum in theSt- $\sigma$ plane is displayed for $\gamma =3\pi$ (i.e. most unstable $\gamma$ found, see figure 11(b) based on the mean flow from the DNS without FST. The DMD eigenvalues were found to be in complex conjugate pairs, therefore, the spectrum is symmetric about $St=0$ . Note that the data in the asymptotic state ( $45\lt t\lt 120$ ) were used for the DMD. The first three most unstable modes are highlighted as solid (red) symbols. It is worth noting that the largest growth rate obtained from the DMD ( $\sigma \approx 0.2$ ) is very close to $\tilde {\sigma }$ computed based on $\hat {q}_{k,max}$ for $\gamma =3\pi$ (see figure 11 b). The spatial eigenfunctions of the most unstable modes are plotted in figure 12(b), which indicate the trapped modes inside the separated region of the mean bubble.

Figure 12.

Dynamic mode decomposition analysis of the LNSE results for $\gamma =3\pi$ based on the mean flow obtained from the DNS with zero FST: (a) eigenvalue spectrum and (b) eigenfunctions of the most unstable modes. Shown are contour lines of $\Re (v)$ -velocity.

A similar analysis was performed for 2-D disturbances ( $\gamma =0$ ) using the mean flow obtained from the DNS with zero FST. The resulting DMD eigenvalue spectrum is displayed in figure 13. Consistent with the observation in figure 11(b), for this case the mean separation bubble is asymptotically stable ( $\sigma \lt 0$ ). Of special interest is the frequency associated with the least stable mode ( $St=0.67-0.68$ ), which is highlighted in figure 13. The relevance of this will be re-adressed later in § 4.

Figure 13.

Dynamic mode decomposition eigenvalue spectrum based on the LNSE data for $\gamma =0$ . Mean flow obtained from the DNS without FST. (b) The DMD eigenmode corresponding to $St=0.67$ as obtained based on LNSE for $\gamma =0$ ; results based on the mean flow field from the DNS without FST.

3.4. Direct numerical simulations with free-stream turbulence

For meaningful comparisons of numerical simulations and wind-tunnel experiments the influence of FST on the LSB has to be addressed. In the DNS results presented and discussed in this paper, low levels of FST were introduced at the inflow boundary by superimposing the velocity and vorticity components of the steady base flow with velocity and vorticity fluctuations. The methodology adopted here to generate free-stream disturbances at the inflow boundary is based on a superposition of eigenmodes from the continuous spectra of the OrrSommerfeld (OS) and homogeneous Squire (SQ) operators (Grosch & Salwen Reference Grosch and Salwen1978; Jacobs & Durbin Reference Jacobs and Durbin2001; Brandt, Schlatter & Henningson Reference Brandt, Schlatter and Henningson2004; Hosseinverdi et al. Reference Hosseinverdi, Balzer and Fasel2012; Hosseinverdi & Fasel Reference Hosseinverdi and Fasel2018, Reference Hosseinverdi and Fasel2019)

(3.7) \begin{equation} \boldsymbol {u}^\prime (x_i, y, z, t)=\sum _\omega \sum _{\gamma } A(|\boldsymbol {k}|)\,\boldsymbol {\phi }(y; \omega , \beta , \gamma )\,e^{i(\gamma z-\omega t)}\,. \end{equation}

Here, $\omega$ is the angular disturbance frequency, $\beta$ and $\gamma$ are the wall-normal and spanwise wavenumbers, respectively, with $|\boldsymbol {k}|= (\omega ^2+\gamma ^2+\beta ^2 )^{0.5}$ . The dispersion relation $\alpha =\omega /U_\infty$ was used in (3.7) to express the streamwise wavenumber ( $\alpha$ ) in terms of the angular frequency. The inflow disturbance vorticity field is calculated from the disturbance velocity field as $\boldsymbol {\omega }^\prime =-\nabla \times \boldsymbol {u}^\prime$ and the eigenfunction $\boldsymbol {\phi }$ is a normalized, randomized and weighted superposition of OS and SQ eigenmodes. The normalization ensures that the energy of each disturbance mode is unity. The coefficients $A(|\boldsymbol {k}|)$ are used to determine the contribution of the eigenfunctions to the total turbulent kinetic energy. For the simulations presented in this work, the amplitudes of the individual inlet perturbations are assigned such that a von Kármán energy spectrum is obtained. A detailed description of the implementation and validation results are provided in Hosseinverdi et al. (Reference Hosseinverdi, Balzer and Fasel2012) and Hosseinverdi & Fasel (Reference Hosseinverdi and Fasel2018, Reference Hosseinverdi and Fasel2019). Several levels of FST have been considered, where the best agreement with the wind-tunnel experiments was obtained for $Tu=0.02\, \%$ , which is comparable to the measured levels of FST in the experiments (§ 2.2). Addition of small levels of FST in the DNS did not change the velocity profiles upstream of the LSB (figure 8 a), thus maintaining the good agreement for the time-averaged flow field upstream of the LSB.

When isotropic FST with an intensity of $Tu=0.02\, \%$ is introduced at the inflow boundary, owing to the large growth rates associated with the inviscid instability mechanism, the vortical disturbances due to the FST are rapidly amplified by several orders of magnitude (Balzer & Fasel Reference Balzer and Fasel2016; Hosseinverdi & Fasel Reference Hosseinverdi and Fasel2019). This leads to an upstream movement of the transition onset location, which in turn reduces the mean bubble size (both in the $x$ - and $y$ -directions) compared with the zero FST case. Analogous to figure 8(b), instantaneous contours of spanwise vorticity (averaged in the $z$ -direction) from the 3-D DNS with FST are shown in figure 14 together with the time- and spanwise-averaged dividing streamline.

Figure 14.

Instantaneous flow visualizations for $Tu=0.02\, \%$ . Shown are contours of spanwise vorticity (averaged in spanwise direction) together with the mean dividing streamline identified by black dashed lines.

As for the case without FST (see figure 11), the LSB obtained from the DNS with FST was found to be globally stable/unstable with respect to the 2-D/3-D disturbances, respectively. Both LSBs (with and without FST) are asymptotically stable with regard to the 2-D mode and become globally unstable for a bounded range of spanwise wavenumbers ( $\gamma$ ), with an increase in maximum temporal growth rates in the DNS with FST compared with the case without FST. Dynamic mode decomposition results for the case with FST were found to be very similar to those presented in § 3.3 for the case without FST and therefore not repeated here. The DMD eigenvalues associated with the least stable mode ( $St=0.67{-}0.68$ ) agree well with those for the case without FST (see figure 13) which confirms that the least stable DMD mode in the presence of FST still corresponds to the dominant shear-layer mode captured in the DNS.

The structure, location and frequency of the most dominant mode in the stability calculations and DNS provide strong evidence of an inviscid shear-layer instability mechanism based on the inflectional velocity profiles. To confirm this conjecture, LST analysis based on the mean flow obtained from the DNS with FST was performed. The downstream evolution of a 2-D $u^{\prime }$ -velocity disturbance wave with $St=0.67$ obtained from the DNS with FST is compared with the LST results in figure 15(a). Based on the stability diagram in figure 10 the location corresponding to the critical Reynolds number for $St=0.67$ is located between the onset of the APG and the separation location ( $x_{cr}\approx 20.7$ ). While the onset of the exponential amplification of the dominant 2-D mode in the DNS is observed further downstream, results for $x\geqslant 24$ inside the separated region are in excellent agreement with the LST.

In figure 15(b) the wall-normal Fourier amplitude distributions of the $u^{\prime }$ -velocity between DNS and LST are compared at several streamwise locations within the region of exponential growth. The amplitude distributions from both LST and DNS are normalized by their respective maximum values within the boundary layer. The amplitude distributions from the DNS and LST agree remarkably well within this region, supporting the fact that the same flow features are captured. The maximum of the $u^{\prime }$ -velocity occurs close to the displacement thickness ( $y/\delta _1=1$ ) of the local boundary-layer profile, which in turn is very close to the local inflection point in the velocity profile. This further supports that the exponential growth of this most dominant mode is due to the inviscid shear-layer instability.

Figure 15.

Comparison of disturbance evolution and eigenfunction distributions between LST and DNS with FST intensity of $Tu=0.02\, \%$ . (a) Downstream development of the Fourier amplitude (maximum inside the boundary layer) for the $u^{\prime }$ -velocity of the 2-D disturbances with $St=0.67$ ( $f={185}\,\mathrm {Hz}$ ). (b) Wall-normal distribution of $u^{\prime }$ -velocity Fourier amplitude taken at the several downstream locations. Disturbance amplitudes are scaled by their respective maximum value and $y$ -coordinate is normalized with the local displacement thickness.