3.1. Externally forced compressible NS equations

To account for the rate of change of perturbation density, momentum and total energy, we model unsteady external disturbances $\boldsymbol {{d}}(\boldsymbol{x}, t)$ to the compressible NS equations in the conservative form (2.1) as volumetric sources of excitation,

(3.1)\begin{equation} \frac{\partial {\boldsymbol{\varPsi}} (\boldsymbol{x}, t)}{\partial t} = \mathcal{F} ({\boldsymbol{\varPsi}} (\boldsymbol{x}, t)) + \boldsymbol{{d}}(\boldsymbol{x}, t), \end{equation}

where $\boldsymbol{x} \mathrel {\mathop :}= (x,y,z)$ is the vector of streamwise, normal and spanwise spatial coordinates. By decomposing the flow field ${\boldsymbol {\varPsi }}$ into the sum of the base $\bar {{\boldsymbol {\varPsi }}}$ and fluctuating ${\boldsymbol {\psi }}$ parts,

(3.2)\begin{equation} {\boldsymbol{\varPsi}} (\boldsymbol{x}, t) = \bar{{\boldsymbol{\varPsi}}} (\boldsymbol{x}) +{\boldsymbol{\psi}} (\boldsymbol{x}, t), \end{equation}

we obtain the equation that governs the dynamics of flow fluctuations around $\bar {{\boldsymbol {\varPsi }}} (\boldsymbol{x})$,

(3.3)\begin{equation} \frac{\partial {\boldsymbol{\psi}} (\boldsymbol{x}, t)}{\partial t} = \mathcal{F} (\bar{{\boldsymbol{\varPsi}}} (\boldsymbol{x}) + {\boldsymbol{\psi}} (\boldsymbol{x}, t)) + \boldsymbol{{d}}(\boldsymbol{x}, t). \end{equation}

For disturbances with small amplitude $\epsilon$,

(3.4a)\begin{equation} \boldsymbol{{d}}(\boldsymbol{x}, t) = \epsilon \boldsymbol{{d}}^{(1)}(\boldsymbol{x}, t), \end{equation}

a weakly nonlinear analysis can be utilized to examine externally forced compressible NS equations (3.3) and determine the corrections to the steady laminar 2-D base flow $\bar {{\boldsymbol {\varPsi }}} (\boldsymbol{x})$. Up to a second order in $\epsilon$, the vector of flow fluctuations ${\boldsymbol {\psi }}$ can be represented as

(3.4b)\begin{equation} {\boldsymbol{\psi}} (\boldsymbol{x}, t) = \epsilon {\boldsymbol{\psi}}^{(1)}(\boldsymbol{x}, t) + \epsilon^2 {\boldsymbol{\psi}}^{(2)}(\boldsymbol{x}, t) + {O} (\epsilon^3), \end{equation}

where ${\boldsymbol {\psi }}^{(1)}(\boldsymbol{x}, t)$ satisfies the linearized flow equations

(3.5a)\begin{equation} \left[\frac{\partial}{\partial t} - \mathcal{A} \! \left( \bar{{\boldsymbol{\varPsi}}} \right)\right] {\boldsymbol{\psi}}^{(1)} = \boldsymbol{{d}}^{(1)}, \end{equation}

and ${\boldsymbol {\psi }}^{(2)}(\boldsymbol{x}, t)$ satisfies

(3.5b)\begin{equation} \left[\frac{\partial}{\partial t} - \mathcal{A} \! \left( \bar{{\boldsymbol{\varPsi}}} \right)\right] {\boldsymbol{\psi}}^{(2)} = \mathcal{N}^{(2)} \! \left( {\boldsymbol{\psi}}^{(1)} \right). \end{equation}

Equations (3.5a) and (3.5b) are respectively obtained by neglecting $O (\epsilon ^2)$ and $O (\epsilon ^3)$ terms upon substitution of (3.4a) and (3.4b) to the compressible NS equations (3.3). The compressible NS operator resulting from linearization around the base flow $\bar {{\boldsymbol {\varPsi }}}$ is determined by $\mathcal {A} (\bar {{\boldsymbol {\varPsi }}})$ (see Candler, Subbareddy & Nompelis (Reference Candler, Subbareddy and Nompelis2015b, equation (23)) and Sidharth et al. (Reference Sidharth, Dwivedi, Candler and Nichols2018, equations (A1)–(A2))) and $\mathcal {N}^{(2)} ( {\boldsymbol {\psi }}^{(1)} )$ accounts for quadratic interactions at $O (\epsilon ^2)$; see Appendix A for details.

Several recent studies demonstrated the utility of the compressible energy norm (Chu Reference Chu1965; Hanifi et al. Reference Hanifi, Schmid and Henningson1996),

(3.6a)\begin{equation} { E \; \mathrel{\mathop:}= \; \dfrac{1}{2} \, \int_{\varOmega} \left( \dfrac{\bar{p}}{\bar{\rho}^{2}} \, \rho^{\prime 2} + \bar{\rho} \left|\boldsymbol{u}^{\prime}\right|^{2} + \dfrac{C_{v} \bar{\rho}}{\bar{T}} \, T^{\prime 2} \right) \,\mathrm{d}\varOmega,} \end{equation}

in quantifying the growth of fluctuations in high-speed boundary layer flows (Franko & Lele Reference Franko and Lele2013; Sidharth et al. Reference Sidharth, Dwivedi, Candler and Nichols2018; Quintanilha et al. Reference Quintanilha, Paredes, Hanifi and Theofilis2022). This quantity is determined by the weighted $L_2$ norm of the vector of flow fluctuations ${\boldsymbol {\phi }} \mathrel {\mathop :}= (\phi _{1}, {\boldsymbol {\phi }}_{2}, \phi _{3}) = (\rho ^\prime, \boldsymbol {u}^\prime, T^\prime )$ in primitive variables,

(3.6b)\begin{equation} { E = \| {\boldsymbol{\phi}} \|_{E}^{2} = \langle {\boldsymbol{\phi}} , {\boldsymbol{\phi}} \rangle_E = \langle {\boldsymbol{\phi}} , \mathcal{W} {\boldsymbol{\phi}} \rangle_2,} \end{equation}

where $\langle\, {\cdot }\, , \,{\cdot }\, \rangle _2$ is the standard $L_2$ inner product over the domain $\varOmega$, $C_{v}$ is the specific heat at constant volume in $\varOmega$ and

(3.6c)\begin{equation} {\mathcal{W} \; \mathrel{\mathop:}= \; \dfrac{1}{2} \left[ \begin{array}{ccc} \bar{p} / \bar{\rho}^{2} & 0 & 0 \\ 0 & \bar{\rho} & 0 \\ 0 & 0 & C_{v} \bar{\rho}/\bar{T} \end{array} \right]} \end{equation}

is the multiplication operator determined by the pressure $\bar {p}$, density $\bar {\rho }$ and temperature $\bar {T}$ of the base flow $\bar {{\boldsymbol {\varPsi }}}$. For small amplitude disturbances, we can represent ${\boldsymbol {\phi }}$ as

(3.7a)\begin{equation} {{\boldsymbol{\phi}} (\boldsymbol{x}, t) = \epsilon {\boldsymbol{\phi}}^{(1)}(\boldsymbol{x}, t) + \epsilon^2 {\boldsymbol{\phi}}^{(2)}(\boldsymbol{x}, t) + {O} (\epsilon^3),} \end{equation}

where

(3.7b)\begin{equation} \left. \begin{gathered} {\boldsymbol{\phi}}^{(1)}(\boldsymbol{x}, t) = \mathcal{C} \, {\boldsymbol{\psi}}^{(1)}(\boldsymbol{x}, t),\\ {\boldsymbol{\phi}}^{(2)}(\boldsymbol{x}, t) = \mathcal{C} \, {\boldsymbol{\psi}}^{(2)}(\boldsymbol{x}, t) + \mathcal{D} \left[ \begin{array}{c} -{\phi^{(1)}_{1}} {{\boldsymbol{\phi}}^{(1)}_{2}} \\ -2 C_{v} {\phi^{(1)}_{1}} {\phi^{(1)}_{3}} - \bar{\varPhi}_{1} |{{\boldsymbol{\phi}}^{(1)}_{2}}|^{2} \end{array} \right]. \end{gathered} \right\} \end{equation}

As shown in Appendix B, $\mathcal {C}$ and $\mathcal {D}$ are multiplication operators parameterized by the laminar 2-D base flow $\bar {{\boldsymbol {\varPsi }}}$; see (B4) for their definition.

In the remainder of this section, we identify oblique waves as the most energetic responses of the linearized flow equations (3.5a) in the presence of unsteady harmonic disturbances $\boldsymbol {{d}}^{(1)}$. In § 4 we utilize (3.5b) to demonstrate that steady streaks can arise from quadratic interactions of unsteady oblique waves.

3.2. Amplification of exogenous disturbances to the linearized flow equations

The linearized equations (3.5a) describe the evolution of the fluctuation vector ${\boldsymbol {\psi }}^{(1)}$ in the presence of external disturbances $\boldsymbol {d}^{(1)}$ and they can be written using the state-space formulation (Schmid & Henningson Reference Schmid and Henningson2001)

(3.8)\begin{equation} \left. \begin{aligned} \frac{\partial {\boldsymbol{\psi}}^{(1)}}{\partial t} & = \mathcal{A} {\boldsymbol{\psi}}^{(1)} + \mathcal{B} \boldsymbol{d}^{(1)}, \\ {\boldsymbol{\phi}}^{(1)} & = \mathcal{C} {\boldsymbol{\psi}}^{(1)}. \end{aligned} \right\} \end{equation}

Here, $\boldsymbol {d}^{(1)}$ is a spatially distributed and temporally varying disturbance (input), ${\boldsymbol {\psi }}^{(1)}$ is the state of the linearized system (which is determined by the vector of flow fluctuations in conserved variables), ${\boldsymbol {\phi }}^{(1)}$ is the quantity of interest (output) whose weighted $L_2$ norm determines the energy of flow fluctuations (3.6) and $\mathcal {A}$ is the generator of the linearized compressible NS dynamics. The input operator $\mathcal {B}$ in (3.8) allows us to specify spatial support of body forcing inputs, and the output operator $\mathcal {C}$ relates the state of the linearized system ${\boldsymbol {\psi }}^{(1)}$ to the vector of flow fluctuations in primitive variables ${\boldsymbol {\phi }}^{(1)}$.

For the parameters shown in figure 2, the linearized system is globally stable and for a time-periodic input with frequency $\omega$, $\boldsymbol {d}^{(1)} (\boldsymbol {x}, t) = \hat {\boldsymbol{d}}^{(1)} (\boldsymbol {x}, \omega ) \mathrm {e}^{\mathrm {i} \omega t}$, the steady-state output of (3.8) is determined by ${\boldsymbol {\phi }}^{(1)} (\boldsymbol {x},t) = \hat {{\boldsymbol {\phi }}}^{(1)} (\boldsymbol {x}, \omega ) \mathrm {e}^{\mathrm {i} \omega t}$, where

(3.9)\begin{equation} \hat{{\boldsymbol{\phi}}}^{(1)} (\boldsymbol{x}, \omega) = \left[ \mathcal{H} (\omega) \hat{\boldsymbol{d}}^{(1)} (\, {\cdot}\, , \omega) \right] (\boldsymbol{x}), \end{equation}

$\mathcal {H} (\omega )$ is the frequency response operator,

(3.10)\begin{equation} \mathcal{H} (\omega) = \mathcal{C}(\mathrm{i} \omega \mathcal{I} \, - \, \mathcal{A})^{{-}1} \mathcal{B}, \end{equation}

and $\mathcal {R} (\omega ) = (\mathrm {i} \omega \mathcal {I} - \mathcal {A})^{-1}$ is the resolvent operator associated with the linearized model (3.8). At any $\omega$, the SVD of $\mathcal {H} (\omega )$ can be used to quantify amplification of time-periodic inputs (Jovanović Reference Jovanović2004; Schmid Reference Schmid2007; Jovanović Reference Jovanović2021),

(3.11)\begin{equation} \hat{{\boldsymbol{\phi}}}^{(1)} (\boldsymbol{x}, \omega) = \left[ \mathcal{H} (\omega) \hat{\boldsymbol{d}}^{(1)} ( \,{\cdot}\, , \omega) \right] (\boldsymbol{x}) = \displaystyle{\sum_i} \, \sigma_i (\omega) {\boldsymbol{\phi}}_{i} (\boldsymbol{x}, \omega) \langle \boldsymbol{d}_{i} ( \,{\cdot}\, , \omega), \hat{\boldsymbol{d}} ( \,{\cdot}\, , \omega) \rangle_E, \end{equation}

where $\sigma _i (\omega )$ denotes the $i$th singular value of $\mathcal {H} (\omega )$, $\langle \,{\cdot }\, , \,{\cdot }\, \rangle _E$ is the inner product in (3.6) that induces the compressible energy norm, and $\boldsymbol {d}_{i} (\boldsymbol {x}, \omega )$ and ${\boldsymbol {\phi }}_{i} (\boldsymbol {x}, \omega )$ are the left and right singular functions of $\mathcal {H} (\omega )$ that provide orthonormal bases of the corresponding input and output spaces (with respect to $\langle \,{\cdot }\, , \,{\cdot }\, \rangle _E$).

The frequency response operator $\mathcal {H} (\omega )$ maps the $i$th input mode $\boldsymbol {d}_{i} (\boldsymbol {x}, \omega )$ into the response whose spatial profile is specified by the $i$th output mode ${\boldsymbol {\phi }}_{i} (\boldsymbol {x}, \omega )$ and the amplification is determined by the corresponding singular value $\sigma _{i} (\omega )$. In other words, for

(3.12)\begin{equation} \hat{\boldsymbol{d}}^{(1)} (\boldsymbol{x}, \omega) = \boldsymbol{d}_i (\boldsymbol{x},\omega) \ \Rightarrow \ \hat{{\boldsymbol{\phi}}}^{(1)} (\boldsymbol{x}, \omega) = \left[ \mathcal{H} (\omega) \boldsymbol{d}_i ( \,{\cdot}\, , \omega) \right] (\boldsymbol{x}) = \sigma_i (\omega) {\boldsymbol{\phi}}_{i} (\boldsymbol{x}, \omega), \end{equation}

and $\| \hat {{\boldsymbol {\phi }}}^{(1)} ( \,{\cdot }\, , \omega ) \|_{E} = \sigma _i (\omega )$. Note that, at any $\omega$,

(3.13)\begin{equation} G (\omega) \; \mathrel{\mathop:}= \; \sigma_{1} (\omega) = \dfrac{\| \mathcal{H} (\omega) \boldsymbol{d}_{1} ( \,{\cdot}\, , \omega)\|_{E}}{\| \boldsymbol{d}_{1} ( \,{\cdot}\, ,\omega)\|_{E}} = \dfrac{\| \sigma_{1} (\omega) {\boldsymbol{\phi}}_{1} ( \,{\cdot}\, , \omega) \|_{E}} {\| \boldsymbol{d}_{1} ( \,{\cdot}\, ,\omega)\|_{E}} \end{equation}

determines the largest induced gain with respect to a compressible energy norm, where ($\boldsymbol {d}_1 (\boldsymbol {x}, \omega ),{\boldsymbol {\phi }}_1 (\boldsymbol {x}, \omega )$) identify the spatial structure of the dominant input–output pair. We use a second-order central finite-volume discretization (Sidharth et al. Reference Sidharth, Dwivedi, Candler and Nichols2018) to obtain a finite-dimensional approximation of the linearized model (3.8) and employ matrix-free Arnoldi iterations (Jeun et al. Reference Jeun, Nichols and Jovanović2016; Dwivedi Reference Dwivedi2020) to compute the singular values $\sigma _i (\omega )$ of $\mathcal {H} (\omega )$.

3.3. Frequency response analysis

We utilize the resolvent analysis to study amplification of harmonic disturbances with frequency $\omega$ to the linearized flow equations. In double-wedge geometry the laminar 2-D base flow $\bar {{\boldsymbol {\varPsi }}}$ is a function of streamwise normal coordinates, $\bar {{\boldsymbol {\varPsi }}} (x,y)$, and owing to homogeneity in the spanwise direction, the 3-D fluctuations in (3.5a) take the form

(3.14)\begin{equation} {\boldsymbol{\psi}}^{(1)} (x,y,z,t) = \hat{{\boldsymbol{\psi}}}^{(1)} (x,y; \beta,\omega) \mathrm{e}^{\mathrm{i}(\beta z + \omega t)}, \end{equation}

where $\beta =2{\rm \pi} /\lambda _z$ is the spanwise wavenumber. Thus, in addition to $\omega$, the frequency response operator is also parameterized by $\beta$,

(3.15)\begin{equation} {\mathcal{H}_\beta (\omega) = \mathcal{C} (\mathrm{i} \omega \mathcal{I} - \mathcal{A}_\beta)^{{-}1} \mathcal{B},} \end{equation}

where $\mathcal {A}_\beta$ denotes the Fourier symbol of the operator $\mathcal {A}$ in (3.8) obtained by replacing the spanwise differential operator $\partial _z$ with $\mathrm {i} \beta$. At any pair $(\omega, \beta )$, $\mathcal {H}_\beta (\omega )$ maps the input function $\hat {\boldsymbol {d}}^{(1)}$ of $x$ and $y$ into the output function $\hat {{\boldsymbol {\phi }}}^{(1)}$ of $x$ and $y$,

(3.16)\begin{equation} \hat{{\boldsymbol{\phi}}}^{(1)} (x, y; \beta, \omega) = \left[ \mathcal{H}_\beta (\omega) \hat{\boldsymbol{d}}^{(1)} ( \,{\cdot}\, , \,{\cdot}\, ; \beta, \omega) \right] (x,y), \end{equation}

and SVD of $\mathcal {H}_\beta (\omega )$ can be used to study amplification across spatio-temporal frequencies.

We first set ${\mathcal {B}} = {\mathcal {I}}$, i.e. we introduce body forcing inputs to excite flow at every spatial location in the computational domain $\varOmega$ and we choose the output operator ${\mathcal {C}}$ to examine the impact of forcing on the compressible energy norm of ${\boldsymbol {\phi }}^{(1)}$ in the entire $\varOmega$. The resolvent analysis is done using a resolution that yields grid-independent outputs with $545$ cells in the streamwise direction, $249$ cells in the normal direction, and numerical sponge boundary conditions near the leading edge ($x=1$) with the outflow boundaries are utilized.

Figure 4 shows the dependence of the input–output gain $G(\omega,\lambda _z)$ on the frequency $\omega$ and the wavelength $\lambda _z$. There are two major amplification regions with the respective peaks at $(\omega =0,\lambda _{z}=1.5)$ and $(\omega =0.4,\lambda _{z}=3)$. The first peak in $G$ identifies the largest amplification and the corresponding output is determined by reattachment streaks that result from steady vortical disturbances upstream of the recirculation zone. We observe selective amplification of disturbances with $\lambda _z \approx 1.5$ and low-pass filtering features over $\omega$. The gain $G$ experiences rapid decay beyond the roll-off frequency $\omega \approx 0.4$ and it attains its largest value at $\omega =0$ for $\lambda _z$ that scales with the reattaching shear layer thickness (Dwivedi et al. Reference Dwivedi, Sidharth, Nichols, Candler and Jovanović2019). In contrast to Dwivedi et al. (Reference Dwivedi, Sidharth, Nichols, Candler and Jovanović2019), which focused on disturbances with $\omega = 0$, we examine unsteady disturbances that trigger oblique waves in the reattaching shear layer, as identified by the second peak in $G$. This amplification region takes place in a narrow band of temporal frequencies $\omega$ over a fairly broad range of spanwise wavelengths $\lambda _{z}$.

Figure 4. (a) Input–output gain $G (\omega,\lambda _z)$ associated with the resolvent operator across temporal frequency $\omega$ and spanwise wavelength $\lambda _z$. (b) Isosurfaces of streamwise vorticity and density fluctuations corresponding to the input–output modes $\boldsymbol {d}_{1}$ and ${\boldsymbol {\phi }}_1$.

As demonstrated in figure 4, even when we allow disturbances to enter through the entire computational domain the largest amplification is caused by inputs that are localized upstream of the corner and the resulting response is localized downstream of the corner. The upstream disturbances are the most effective way to excite the flow because of the large convection velocity of the laminar 2-D base flow (Chomaz Reference Chomaz2005; Schmid Reference Schmid2007) and the dominant output emerges in the separated and reattached regions.

Experimental studies of oblique transition in channel and boundary layer flows (Elofsson & Alfredsson Reference Elofsson and Alfredsson1998; Berlin & Henningson Reference Berlin and Henningson1999) often utilize streamwise localized disturbances and a common criticism of the resolvent analysis is that the identified global input modes represent excitation sources that are not easy to realize experimentally. In contrast, traditional approach to the analysis of boundary layers utilizes spatially localized fluctuation sources and evaluates the streamwise growth of fluctuations as they convect downstream (Herbert Reference Herbert1997). However, in the presence of flow separation a parabolized approximation of the NS equations cannot be made. To evaluate amplification in different spatial regions, we restrict inputs and outputs to belong to a plane but still account for the global nature of the separated flow through the resolvent operator $(\mathrm {i}\omega \mathcal {I} - \mathcal {A}_{\beta })^{-1}$. As illustrated in figure 5(a), this is accomplished via a proper selection of the operators $\mathcal {B}$ and $\mathcal {C}$ in (3.8) by fixing the input location before flow separation, at $x_{in} = 25$, and by evaluating the output at different locations downstream of the separation, $x_{{out}}$. In this set-up, $G_{x_{{out}}}$ quantifies the largest amplification at $x_{{out}}$ of disturbances that are introduced at $x_{in} = 25$.

Figure 5. Spatial input–output analysis: (a) input is introduced at a streamwise location $x_{in} = 25$ before separation and output is evaluated at $x_{{out}}$; (b) dependence of the input–output gain $G_{x_{{out}}}$ on the streamwise location $x_{{out}}$ for streaks and oblique waves. Unsteady oblique waves with $(\omega =0.4,\lambda _{z}=3)$ are strongly amplified throughout the separation zone.

Figure 5(b) shows the dependence of $G_{x_{{out}}}$ on $x_{{out}}$ for upstream disturbances with $(\omega =0,\lambda _{z}=1.5)$ and $(\omega =0.4,\lambda _{z}=3)$. The gain associated with the steady fluctuations begins to grow in the latter half of the recirculation zone, especially near the reattachment location. In contrast, unsteady perturbations with $\omega = 0.4$ experience significant amplification throughout the separation zone. This observation suggests that the separation zone plays a critical role in the amplification of unsteady fluctuations.

3.4. Amplification of oblique waves: physical mechanism

We now analyse physical mechanisms responsible for amplification of flow fluctuations within the separation zone in the presence of upstream unsteady disturbances. In particular, we examine the global response of the linearized equations to the input with ($\omega =0.4,\lambda _{z}=3.0$) introduced prior to separation (i.e. at $x_{in}=25$) that triggers the largest amplification in the entire domain $\varOmega$. The spatial structure of flow fluctuations is studied in the $(s,n,z)$ coordinate system which is locally aligned with the streamlines of the laminar 2-D base flow $\bar {{\boldsymbol {\varPsi }}}$; see figure 3 for an illustration. In this frame of reference, $\bar {\boldsymbol {u}} = (\bar {u}_s,0,0)$ with $\bar {u}_{s} \geq 0$, and, as discussed in Finnigan (Reference Finnigan1983), Patel & Sotiropoulos (Reference Patel and Sotiropoulos1997), Dwivedi et al. (Reference Dwivedi, Sidharth, Nichols, Candler and Jovanović2019), this coordinate system is convenient for the analysis of separated boundary layers especially within the separation zone.

The streamwise specific kinetic energy $\mathcal {E}_{s} \mathrel {\mathop :}= u_s^{\prime } u_s^{\prime }$ obeys the transport equation,

(3.17)\begin{equation} \frac{\partial \mathcal{E}_{s}}{\partial t} + \bar{u}_{s}\frac{\partial \mathcal{E}_{s}}{\partial s} = \mathcal{P} + \mathcal{S} + \mathcal{V} + {\mathcal{K}} + \mathcal{F}, \end{equation}

where $\mathcal {P}$, $\mathcal {S}$, $\mathcal {V}$ and $\mathcal {K}$ are the production, source, viscous and curvature terms (see Appendix C), and $\mathcal {F}$ is the work done by external disturbances (Dwivedi et al. Reference Dwivedi, Sidharth, Nichols, Candler and Jovanović2019, appendix C). The production term $\mathcal {P}$ quantifies interactions of fluctuation stresses with the base flow gradients, the source term $\mathcal {S}$ corresponds to the perturbation component of the inviscid material derivative, the viscous term $\mathcal {V}$ determines dissipation of kinetic energy by viscous stresses and $\mathcal {K}$ accounts for the curvature that arises from a coordinate transformation. Our computations indicate that while the dissipative viscous term $\mathcal {V}$ is negative throughout the domain, the production term $\mathcal {P}$ is sign indefinite, $\mathcal {S}$ and $\mathcal {K}$ are negligible, and $\mathcal {F}$ is zero downstream of the forcing plane.

Insight into physical mechanisms can be gained by the analysis of dominant production terms in (3.17) associated with the global linearized response to upstream oblique disturbances. Averaging over the time period $T = {2{\rm \pi} }/{\omega }$ and the spanwise wavelength $\lambda _{z} = 2 {\rm \pi}/ \beta$, $\langle \,{\cdot }\,\rangle \mathrel {\mathop :}= ( T \lambda _{z} )^{-1}\int ^T_{0} \int ^{\lambda _{z}}_{0} (\,{\cdot }\,) \, \mathrm {d} z \,\mathrm {d} t$, and neglecting the terms that do not contribute significantly to the production of the averaged streamwise specific kinetic energy $E_{s} \mathrel {\mathop :}= \langle \mathcal {E}_s \rangle$ yields the following approximation to the transport equation (3.17):

(3.18)\begin{equation} \bar{u}_{s} \frac{\partial E_{s}}{\partial s} + 2 ( \partial_s \bar{u}_{s} ) E_{s} \approx{-} 2 ( \partial_n \bar{u}_{s} ) {R}_{sn}. \end{equation}

Here ${R}_{sn} \mathrel {\mathop :}= \langle {u_s^{\prime } u_n^{\prime }}\rangle $ denotes the averaged shear stress of the streamwise velocity fluctuations. The second term on the left-hand side represents the production of fluctuations’ energy that arises from the streamwise gradient of the base flow $\partial _{s}\bar {u}_{s}$, and the term on the right-hand side determines the production term that originates from interactions of the base flow shear $\partial _n \bar {u}_{s}$ with the fluctuation shear stress ${R}_{sn}$.

To understand the mechanism that facilitate the growth of $E_{s}$, we now investigate the streamwise transport of ${R}_{sn}$. In contrast to the transport equation for $E_{s}$, both the production $\mathcal {P}$ and curvature $\mathcal {K}$ terms contribute significantly to the streamwise transport of ${R}_{sn}$ for fluctuations with $\omega =0.4$ and $\lambda _{z}=3.0$. As demonstrated in Appendix D, omitting negligible terms leads to the following approximate transport equation for ${R}_{sn}$:

(3.19a)\begin{equation} \bar{u}_{s} \frac{\partial {R}_{s n}}{\partial s} + ( \partial_s \bar{u}_{s} + K_{s} ) {R}_{s n} \approx 2 K_{c} E_{s}. \end{equation}

Here $K_{c}$ and $K_{s}$ denote contributions that arise from the curvature normal to the streamlines and from deceleration along the streamline direction, respectively,

(3.19b)\begin{equation} {K_c ={-}(\bar{\varOmega} + \partial_{n}{\bar{u}_{s}}), \quad K_{s} ={-}\partial_{s} \bar{u}_{s}.} \end{equation}

In the $(s,n,z)$ coordinate system, $\bar {\varOmega } = \partial _{x}\bar {v} - \partial _{y}\bar {u}$ denotes the spanwise vorticity of the base flow in Cartesian coordinates (Finnigan Reference Finnigan1983) and using the definition of $K_{s}$, (3.19a) simplifies to

(3.19c)\begin{equation} \bar{u}_{s} \frac{\partial {R}_{s n}}{\partial s} \approx 2 K_{c} E_{s}. \end{equation}

In summary, (3.18) and (3.19c) determine a coupled system of linear equations that governs the streamwise transport of ${E_{s}}$ and ${R}_{s n}$ in the separation zone for oblique fluctuations with ($\omega =0.4,\lambda _{z}=3.0$),

(3.20)\begin{equation} \left[ \begin{array}{ccc} {\bar{u}_s}&{0}\\{0}&{\bar{u}_s}\end{array} \right]\frac{\partial}{\partial s}\left[ \begin{array}{ccc} {E_{s}}\\{{R}_{s n}}\end{array} \right]\approx \left[ \begin{array}{ccc}{-2 \partial_s \bar{u}_{s}}&{-2 \partial_n \bar{u}_{s}}\\{2 K_c}&{0}\end{array} \right]\left[ \begin{array}{ccc}{E_{s}}\\{{R}_{s n}} \end{array} \right]. \end{equation}

Oblique waves experience the largest amplification in the separated shear layer above the recirculation bubble, i.e. in the region where the presence of flow separation leads to concave flow curvature $K_c<0$. Figure 6(a) shows this negative curvature along the separation streamline and figure 6(b) illustrates the physical mechanism which is absent in the attached boundary layers because of negligibly small positive streamwise curvature.

Figure 6. (a) Curvature ($-K_{c}$) of the laminar base flow along the separation streamline; (b) illustration of a physical mechanism that facilitates growth of the averaged streamwise specific kinetic energy $E_{s}$ of oblique fluctuations in the separated shear layer.

Concave base flow curvature (i.e. $K_{c} < 0$) in the shear layer provides the destabilizing effect in system (3.20) that can be understood by simplifying (3.20) for oblique waves. In figure 7(a) we compare $\bar {u}_{s} \partial _s E_{s} \mathrel {\mathop :}= \bar {u}_{s} \partial E_{s} / \partial s$ with the dominant production term in (3.18) to illustrate that $( \partial _n \bar {u}_{s} ){R}_{sn}$ dictates the streamwise growth of $E_{s}$. Furthermore, since the base shear $\partial _n \bar {u}_{s}$ remains almost constant throughout the shear layer (cf. figure 7b), its streamwise derivative can be neglected, thereby leading to

(3.21)\begin{equation} {\bar{u}^2_{s} \, \frac{\partial^2 {E}_{s}}{\partial s^2} + \bar{u}_{s} (\partial_s \bar{u}_s) \, \frac{\partial {E}_{s}}{\partial s} + 4 K_c ( \partial_n \bar{u}_s ) {E}_{s} \approx 0.} \end{equation}

This second-order differential equation for $E_s$ is obtained by taking the derivative of the equation for $E_s$ in (3.20), keeping the dominant terms and substituting the equation for $R_{sn}$ from (3.20) into the resulting expression. Figure 8(a) shows the streamwise evolution of $E_s$ and figure 8(b) compares the coefficients in (3.21). Since the effect of $\partial _{s} \bar {u}_{s}$ is negligible, (3.21) can further be simplified to obtain

(3.22)\begin{equation} { \bar{u}^2_{s} \, \frac{\partial^2 {E}_{s}}{\partial s^2} + 4 K_c ( \partial_n \bar{u}_s ) {E}_{s} \approx 0.} \end{equation}

Figure 7. Streamwise variation of (a) $\bar {u}_s \partial _{s} E_{s}$ along with the dominant production term in (3.18); (b) average fluctuation shear stress ${R}_{sn}$ and base flow shear $\partial _{n} \bar {u}_{s}$.

Figure 8. (a) Streamwise evolution of the wall-normal integral of $E_{s} \mathrel {\mathop :}= \langle u^\prime _s u^\prime _s \rangle$ for the primary output resolvent mode of oblique fluctuations with ($\omega =0.4,\lambda _{z}=3.0$) along with contours of $E_{s}$ in the separated shear layer near reattachment (inset); (b) the coefficients in (3.21) along the separation shear layer.

As shown in figure 8(b), the concave base flow curvature (i.e. $K_c < 0$) provides the destabilizing influence throughout the separated shear layer and a simple mechanical analogy can be used to explain amplification of oblique waves. In the regions where $K_c < 0$ the ‘spring constant’ $4 K_c ( \partial _n \bar {u}_s)$ in (3.22) is negative and this system behaves as an inverted pendulum, which enables the spatial growth of $E_s$.

In summary, we have utilized the resolvent analysis to identify the spatial structure of oblique fluctuations that amplify rapidly in the separation zone. Furthermore, by conducting transport analysis of the most energetic fluctuations, we have demonstrated that the resulting amplification arises from concave curvature of the laminar 2-D base flow.

4.1. Streaks generated by oblique waves: a weakly nonlinear analysis

In the presence of a small external disturbance,

(4.1)\begin{equation} \boldsymbol{d}(x,y,z,t) = \epsilon \, ({\hat{\boldsymbol{d}}^{(1)}_{+}}(x,y) \mathrm{e}^{\mathrm{i} \omega t} + { \hat{\boldsymbol{d}}^{(1)}_{-}}(x,y)\mathrm{e}^{-\mathrm{i} \omega t})\,\mathrm{e}^{\mathrm{i} \beta z}, \end{equation}

a weakly nonlinear analysis of § 3.1 can be utilized to represent the flow state components in compressible NS equations (3.3) as

(4.2)\begin{equation} \left. \begin{aligned} {{O} (1):} & \quad \bar{{\boldsymbol{\varPsi}}} ({\boldsymbol{x}}, t) = \bar{{\boldsymbol{\varPsi}}} (x,y),\\ {{O} (\epsilon):} & \quad {\boldsymbol{\psi}}^{(1)}({\boldsymbol{x}}, t) = \left( \hat{{\boldsymbol{\psi}}}^{(1)}_{+}(x,y)\mathrm{e}^{\mathrm{i} \omega t} + \hat{{\boldsymbol{\psi}}}^{(1)}_{-}(x,y)\mathrm{e}^{-\mathrm{i} \omega t} \right) \mathrm{e}^{\mathrm{i} \beta z},\\ {{O} (\epsilon^2):} & \quad {\boldsymbol{\psi}}^{(2)}({\boldsymbol{x}}, t) = \left( \hat{{\boldsymbol{\psi}}}^{(2)}_{0} (x,y) + \hat{{\boldsymbol{\psi}}}^{(2)}_{+} (x,y) \mathrm{e}^{2 \mathrm{i} \omega t} + \hat{{\boldsymbol{\psi}}}^{(2)}_{-} (x,y) \mathrm{e}^{{-}2 \mathrm{i} \omega t} \right) \mathrm{e}^{2 \mathrm{i} \beta z}. \end{aligned} \right\} \end{equation}

Here, $\bar {{\boldsymbol {\varPsi }}} (x,y)$ is the 2-D laminar base flow, $\hat {\boldsymbol {d}}^{(1)}_{\pm }$ and $\hat {{\boldsymbol {\psi }}}^{(1)}_{\pm }$ are the principal oblique input and state modes resulting from the linearized analysis of § 3, whereas $\hat {{\boldsymbol {\psi }}}^{(2)}_{0}$ and $\hat {{\boldsymbol {\psi }}}^{(2)}_{\pm }$ are the steady and harmonic components of the state at ${O}(\epsilon ^2)$. At ${O}(\epsilon ^2)$, the fluctuation's dynamics is governed by (3.5b), where the steady component $\hat {{\boldsymbol {\psi }}}^{(2)}_{0}$ satisfies

(4.3)\begin{equation} { \left[ \mathcal{A}_{2 \beta} \, \hat{{\boldsymbol{\psi}}}^{(2)}_{0} ( \,{\cdot}\, , \,{\cdot}\, ) \right] (x,y) ={-} \hat{\boldsymbol{d}}^{(2)}_0 (x,y).} \end{equation}

Here, $\mathcal {A}_{2 \beta }$ is the Fourier symbol of the dynamical generator in the linearized state-space model (3.8) and $\hat {\boldsymbol {d}}^{(2)}_0 \mathrel {\mathop :}= \mathcal {N}^{(2)}_0 ( \hat {{\boldsymbol {\psi }}}^{(1)}_{\pm } )$ is the forcing term that arises from quadratic interactions of ${O} (\epsilon )$ oblique waves with the spanwise wavenumber $\beta$; see Appendix A for details. Thus, the resolvent operator associated with (3.8) evaluated at ($\omega = 0$, $2 \beta$) maps the nonlinear modulation $\hat {\boldsymbol {d}}^{(2)}_0$ of ${O}(\epsilon )$ oblique waves to ${O}(\epsilon ^{2})$ steady streamwise streaks,

(4.4)\begin{equation} {\hat{{\boldsymbol{\psi}}}^{(2)}_{0} (x,y) = \left[ \mathcal{R}_{2 \beta} (0) \hat{\boldsymbol{d}}^{(2)}_{0} ( \,{\cdot}\, , \,{\cdot}\, ) \right] (x,y) ={-} \left[ \mathcal{A}_{2 \beta}^{{-}1} \, \hat{\boldsymbol{d}}^{(2)}_{0} ( \,{\cdot}\, , \,{\cdot}\, ) \right] (x,y).} \end{equation}

To investigate the emergence of streaks from unsteady disturbances, we introduce forcing inputs with ($\omega =\pm 0.4, \lambda _{z}=3$) and examine a weakly nonlinear evolution of the resulting oblique waves. These forcing inputs are introduced at the upstream plane $x_{{in}} = 25$ and their spatial structure is identified using the resolvent analysis of § 3 to generate the most energetic response at the reattachment (i.e. at $x_{out}=60$).

Figure 9(a) illustrates the set-up in which disturbances corresponding to a pair of oblique input modes $\hat {\boldsymbol {d}}^{(1)}_{\pm }$ with ($\omega =\pm 0.4,\ \lambda _{z}=3$) are introduced at $x_{in}= 25$. The resulting response of the linearized dynamics consists of oblique waves with opposite phase velocities, leading to a checkerboard wave pattern in the spanwise direction. Figure 9(b) shows the streamwise velocity component of the steady response $\hat {{\boldsymbol {\phi }}}^{(2)}_{0} (x,y)$ at ${O} (\epsilon ^2)$ that arises from weakly nonlinear interactions of ${O} (\epsilon )$ oblique waves. The steady response is given by streamwise streaks with half the spanwise wavelength $\lambda ^{\textit {streaks}}_{z} = \lambda ^{\textit {oblique}}_{z}/2 = 1.5$ of the oblique input.

Figure 9. (a) Set-up for weakly nonlinear analysis: a pair of dominant input modes with ($\omega =\pm 0.4, \lambda _{z}=3$) resulting from resolvent analysis is introduced at $x_{{in}} = 25$ and the corresponding streamwise velocity fluctuations arise as the output of the linearized dynamics. Panel (b) shows that ${O} (\epsilon ^2)$ steady streamwise streaks with $\lambda _z = 1.5$ are triggered by weakly nonlinear interactions of ${O} (\epsilon )$ oblique waves with $\lambda _{z}=3$.

A weakly nonlinear analysis allows us to demonstrate that steady streaks at ${O}(\epsilon ^2)$ arise from quadratic interactions of ${O}(\epsilon )$ oblique waves. Figure 10(a) utilizes a wall-aligned ($\xi,\eta$) coordinate system to illustrate the forcing term $\hat {\boldsymbol {d}}^{(2)}_0 \mathrel {\mathop :}= \mathcal {N}^{(2)}_0 ( \hat {{\boldsymbol {\psi }}}^{(1)}_{\pm } )$ in (4.3), where $\xi$ and $\eta$ denote the directions parallel and normal to the wall, respectively. Large amplification of oblique waves that result from linearized analysis in the reattachment region triggers strongest forcing $\hat {\boldsymbol {d}}^{(2)}_{0}$ in that region. Figure 10(b) shows the wall-normal profiles of the forcing term to the mass, momentum and temperature equations in (4.3) before reattachment at $x = 58$. We observe the strongest contribution of the forcing to the wall-normal and spanwise components of the momentum equations, thereby demonstrating its vortical nature. Figure 10(c) illustrates the spatial structure of the forcing near reattachment in the ($z,\eta$) plane. The forcing term $\hat {\boldsymbol {d}}^{(2)}_{0}$ which forms counter-rotating vortices in the separated shear layer is $90^{\circ }$ out of phase relative to the induced streak response $u^\prime _s$. In contrast to the dominant vortical forcing resulting from the linearized analysis, the vortical source term that arises from weakly nonlinear interactions of oblique waves primarily lies downstream of the recirculation zone.

Figure 10. (a) Real part of the streamwise vorticity forcing component that arises from weakly nonlinear interactions of oblique waves $\hat {\boldsymbol {d}}^{(2)}_0$ in the ($x,y$) plane along with the base flow streamlines. (b) Wall-normal profiles of the forcing terms to the streamwise vorticity, density, temperature and streamwise velocity equations that originate from $\hat {\boldsymbol {d}}^{(2)}_0$ before reattachment, at $x = 58$. (c) Spatial structure of the forcing $\hat {\boldsymbol {d}}^{(2)}_0$ near reattachment in the ($z,\eta$) plane (colour plots) along with the resulting streak response $\hat {{\boldsymbol {\phi }}}^{(2)}_0$ (contour lines) at $x = 58$.

We utilize DNS to verify the predictions of weakly nonlinear analysis. In DNS the input oblique modes $\hat {\boldsymbol {d}}^{(1)}_{\pm }$ (with $\omega =\pm 0.4$, $\lambda _z = 3$ and reference amplitude $A_{0} = 1.0$) are introduced in the plane $x_{in} = 25$. Additional details about the grid resolution and numerical implementation in our DNS study are provided in § 5. Figure 11 shows the spatial evolution of steady streamwise velocity fluctuations $u^\prime _s$ along the base flow separation streamline that are triggered by unsteady oblique disturbances of different amplitudes. As shown in figure 11(a), steady streaks undergo similar spatial growth across the range of forcing amplitudes. Figure 11(b) illustrates that the magnitude of streamwise velocity fluctuations collapses when scaled with $A^2_0$. This demonstrates an excellent agreement between streak profiles resulting from DNS and a weakly nonlinear analysis (that leads to (4.3)). Deviations from predictions of the weakly nonlinear analysis are only observed for the largest amplitude considered and they are manifested by the saturation of streaks in the post-reattachment region.

Figure 11. (a) Streamwise velocity fluctuations $u_s^{\prime }$ along the separation streamline associated with ${O}(\epsilon ^{2})$ steady streaks at $\lambda _z = 1.5$; DNS results for various amplitudes of oblique disturbances with ($\omega =\pm 0.4, \lambda _{z}=3$) are shown. (b) Direct numerical simulation results normalized with the square of the amplitude $A_{0}$ of oblique disturbances are compared with the results of weakly nonlinear analysis.

Collapse of spatial profiles that characterize the amplification of streaks irrespective of the amplitude of oblique disturbances, which differ by ${O}(10^4)$, shows that streaks generated via quadratic interactions of oblique waves undergo linear amplification in the separation zone. This demonstrates the predictive power of the weakly nonlinear analysis across the range of forcing amplitudes. In what follows, we utilize the input–output modes obtained from the resolvent analysis to characterize evolution of ${O} (\epsilon ^2)$ steady streamwise streaks and uncover corresponding amplification mechanisms.

4.2. Resolvent mode representation of ${O}(\epsilon ^{2})$ streaks

As described in § 3.2, the left and right singular function of the frequency response operator provide orthonormal bases of input and output spaces that can be used to study responses of the double-wedge flow to external excitations. In particular, ${O} (\epsilon ^{2})$ steady streaks resulting from weakly nonlinear interactions of ${O} (\epsilon )$ unsteady oblique waves (cf. (4.4)) can be represented using SVD of the frequency response operator associated with the linearized system (3.8) at $\omega = 0$ and the spanwise wavenumber $2 \beta$,

(4.5)\begin{equation} { \hat{{\boldsymbol{\phi}}}^{(2)}_{0} (x,y) = \left[ \mathcal{H}_{2 \beta} (0) \hat{\boldsymbol{d}}^{(2)}_{0} ( \,{\cdot}\, , \,{\cdot}\, ), \right] (x,y) = \displaystyle{\sum_i} \, a_i {\boldsymbol{\phi}}_{i} (x, y).} \end{equation}

Here, $a_{i} \mathrel {\mathop :}= \sigma _i \langle \boldsymbol {d}_{i}, \hat {\boldsymbol {d}}^{(2)}_{0} \rangle _E$, with $\hat {\boldsymbol {d}}^{(2)}_{0} \mathrel {\mathop :}= \mathcal {N}^{(2)}_{0} (\hat {{\boldsymbol {\psi }}}_{\pm }^{(1)})$, quantifies the contribution of the $i$th output mode ${\boldsymbol {\phi }}_{i}$ of $\mathcal {H}_{2 \beta } (0)$ to ${O} (\epsilon ^{2})$ steady streaks $\hat {{\boldsymbol {\phi }}}^{(2)}_{0}$, $\sigma _{i}$ is the $i$th singular value, ($\boldsymbol {d}_{i},{\boldsymbol {\phi }} _{i}$) are the corresponding input–output modes of $\mathcal {H}_{2 \beta } (0)$ and, for $\lambda _z = 1.5$, the inner product $\langle \,{\cdot }\, , \,{\cdot }\, \rangle _{E}$ is carried over the entire flow domain in ($x,y$).

Figure 12(a) shows the $15$ largest singular values of the resolvent for the linearized system (3.8) with ($\omega = 0, \lambda _{z} = 1.5$). Even though the principal singular value $\sigma _1$ is an order of magnitude larger than $\sigma _2$, figure 12(b) demonstrates that the second output mode ${\boldsymbol {\phi }}_2$ contributes most to $\hat {{\boldsymbol {\phi }}}^{(2)}_{0}$. Figure 13(a) shows the wall-normal profiles of the streamwise velocity component $u^{\prime }$ associated with $\hat {{\boldsymbol {\phi }}}^{(2)}_{0}$ and the first two output modes (${\boldsymbol {\phi }}_1,{\boldsymbol {\phi }}_2$) of the resolvent. We observe a striking similarity between ${O} (\epsilon ^{2})$ steady streaks and ${\boldsymbol {\phi }}_{2}$ in the post-reattachment region, at $x = 65$. Similarly, figure 13(b) compares the wall-normal shapes of the corresponding input modes $\boldsymbol {d}_1$ and $\boldsymbol {d}_2$ with the forcing $\hat {\boldsymbol {d}}^{(2)}_0$ that arises from quadratic interactions. The input modes are visualized in the reattaching shear layer, at $x = 57$, and the streamwise vorticity component of $\boldsymbol {d}_2$ provides a good approximation to the vortical forcing that captures interactions of unsteady oblique fluctuations.

Figure 12. (a) Singular values of the resolvent operator associated with the linearized dynamics around the 2-D laminar flow evaluated at $(\omega,\lambda _z) = (0,1.5)$. (b) Contribution of the $i$th output mode ${\boldsymbol {\phi }}_{i}$ to the energy of ${O}(\epsilon ^{2})$ steady streaks $\hat {{\boldsymbol {\phi }}}^{(2)}_{0}$ that are triggered by weakly nonlinear interactions of ${O}(\epsilon )$ oblique waves.

Figure 13. The wall-normal profiles of the real part of streamwise velocity fluctuations resulting from weakly nonlinear interactions of oblique waves and corresponding to (a) ${O}(\epsilon ^{2})$ steady streaks $\hat {{\boldsymbol {\phi }}}^{(2)}_{0}$ and the first two output resolvent modes after reattachment, at $x=65$; (b) steady vortical forcing $\hat {\boldsymbol{d}}^{(2)}_{0}$ and the first two input resolvent modes in the reattaching shear layer, at $x=57$.

As demonstrated above, near reattachment, ${O}(\epsilon ^{2})$ streaks are well approximated by the second output mode of the resolvent associated with the linearized equations at ($\omega = 0, \lambda _z = 1.5$). To gain insight into the amplification mechanism that generates ${O} (\epsilon ^{2})$ streaks, we examine the dominant terms in the energy transport equation for the output mode ${\boldsymbol {\phi }}_{2}$. Similar to the analysis in § 3.4, we utilize the ($s,n,z$) coordinate system which is locally aligned with the streamlines of the base flow $(\bar {u}_{s},0,0)$. The transport of the spanwise-averaged specific kinetic energy of streamwise velocity fluctuations $E_{s} = \langle {u^\prime _{s} u^{\prime }_s} \rangle$ and fluctuation shear stress ${R}_{sn} = \langle u^\prime _{s} u^{\prime }_n \rangle$ associated with the second output mode is approximately governed by (3.20) in § 3.4.

Figure 14(a) shows the streamwise evolution of the wall-normal integral of $E_{s}$ for the mode ${\boldsymbol {\phi }}_{2}$ associated with ${O} (\epsilon ^{2})$ streaks with ($\omega = 0,\lambda _z = 1.5$). In contrast to the oblique waves (cf. figure 8a), which experience monotonic amplification throughout the separation shear layer above the recirculation bubble, we observe a non-monotonic $x$-dependence of $E_s$ for the streaks within the bubble. To gain physical insight, we evaluate the terms of the coefficient matrix in (3.20) for steady streaks. Figure 14(b) shows that both streamwise deceleration (i.e. $\partial _{s} \bar {u}_{s} < 0$) and shear $\partial _{n} \bar {u}$ contribute to amplification of ${O} (\epsilon ^{2})$ streaks in different regions of the separated laminar flow. Towards reattachment, the base flow curvature $K_c$ is positive inside the recirculation bubble (cf. figure 14c) and the streamwise deceleration term (i.e. $\partial _{s} \bar {u}_{s} < 0$) is primarily responsible for the energy $E_{s}$ amplification. In this region, the transport equation for $E_{s}$ can be approximated by

(4.6)\begin{align} { \bar{u}_{s} \frac{\partial E_{s}}{\partial s} \approx{-} 2 ( \partial_s \bar{u}_{s} ) {E}_{s}.} \end{align}

This is in concert with the analysis of the contribution of the first (i.e. most amplified) output mode ${\boldsymbol {\phi }}_{1}$ to the amplification of streaks in a compression ramp flow (Dwivedi et al. Reference Dwivedi, Sidharth, Nichols, Candler and Jovanović2019), which also revealed dominance of streamwise deceleration near reattachment. On the other hand, in the post-reattachment region the shear term $\partial _{n} \bar {u}_{s} < 0$ dominates and the coupled system of (3.20) for $E_s$ and $R_{sn}$ simplifies to the following second-order equation (3.22) for $E_s$:

(4.7)\begin{align} { \bar{u}^2_{s} \, \frac{\partial^2 {E}_{s}}{\partial s^2} + 4 K_c ( \partial_n \bar{u}_s ) {E}_{s} \approx 0.} \end{align}

Thus, in the post-reattachment region the concave streamline curvature of the laminar 2-D base flow (i.e. $K_c < 0$) is primarily responsible for amplification of ${O} (\epsilon ^{2})$ streaks.

Figure 14. (a) Streamwise evolution of the wall-normal integral of $E_{s} = \langle {u^\prime _{s} u^{\prime }_s} \rangle$ for the second output mode associated with streaks at ($\omega = 0, \lambda _z = 1.5$) along with contours of $E_{s}$ near reattachment (inset); (b) streamwise variation of deceleration and shear components of the spanwise-averaged production term $\langle \mathcal {P} \rangle$ in the transport equation for $E_{s}$; and (c) base flow curvature $K_{c}$ near reattachment.

5.1. Numerical set-up

To study the onset of transition, we extend the computational domain in the streamwise direction downstream of reattachment. We place the inflow boundary at $x=20$ and, at this location, we interpolate the inflow profile that results from 2-D base flow computations. To avoid spurious numerical errors, we utilize a non-reflecting numerical sponge near the outflow. The wall is assumed to be adiabatic and periodic boundary conditions in the spanwise direction are applied.

A computational domain of size $80 \times 13 \times 9$ in the streamwise, wall-normal and spanwise directions is discretized using $900 \times 249 \times 384$ grid points (i.e. $86$ million cells). The grid is constructed to ensure uniform spacing in the streamwise and spanwise directions and, in the normal direction, the mesh is stretched to ensure $y^{+} < 0.22$ at the wall. In Appendix E we provide a grid convergence study and compare the numerical resolution that we use with recent simulations of supersonic and hypersonic transitional and turbulent flows. Our simulations utilize low-dissipation sixth-order spatially accurate kinetic-energy-consistent (KEC) fluxes (Subbareddy & Candler Reference Subbareddy and Candler2009) for spatial discretization. The low-dissipation KEC fluxes were previously employed in high-fidelity simulations of transitional and turbulent hypersonic boundary layers (Subbareddy & Candler Reference Subbareddy and Candler2012; Subbareddy, Bartkowicz & Candler Reference Subbareddy, Bartkowicz and Candler2014). The time integration is carried out using the explicit third-order Runge–Kutta scheme with a Courant–Friedrichs–Lewy number of $0.5$.

5.2. Secondary instability and breakdown

To simulate the breakdown of streaks, we excite the laminar 2-D flow using the oblique input modes with ($\omega = \pm 0.4, \lambda _z=3$). The disturbance amplitude is set to $A_{ob} = 2.50 \times 10^{2} {A}_{0}$ and is determined based on the results reported in figure 11. With this amplitude, fluctuations grow linearly in the recirculation region but saturate nonlinearly post-reattachment (i.e. beyond $x=60$).

In figure 15 we show the time evolution of streaks in the plane close to the wall, at $\eta = 0.25$, for three different values of $x$. These plots demonstrate that, in the presence of unsteady oblique disturbances, ${O}(\epsilon ^{2})$ steady streaks undergo spanwise motion close to reattachment. The identified spanwise oscillations become stronger as we progress downstream and their period corresponds to the time period of oblique wave inputs. We note the presence of a ‘sinuous-subharmonic’ motion, where two adjacent low-speed streaks oscillate out of phase, and observe that the amplitude of spanwise oscillations almost doubles as we move from $x=55$ to $x = 65$.

Figure 15. Streamwise streaks in the plane close to the wall, $\eta = 0.25$, at three different streamwise locations: (a) before reattachment, at $x = 55$; (b) at reattachment, at $x = 60$; and (c) after reattachment, at $x=65$. The low-speed streaks are marked in blue (solid lines) and high-speed streaks are marked in red (dashed lines).

In figure 16 we illustrate the effect of streak oscillations on the mean flow by visualizing the root mean square (r.m.s.) of temporal fluctuations in the streamwise velocity. At different streamwise locations, fluctuations are plotted against $\partial \bar {u}/\partial z$, which characterizes steady spanwise variations of the mean flow. Initially, unsteadiness is restricted to the oblique wave pair, which is located further away from the wall relative to the streaks. However, at $x=65$, unsteadiness is observed in the region of the largest spanwise shear because of the presence of streaks. This ‘locking-in’ of $u_{rms}$ with streak oscillations identifies late stages of the streak evolution just before smaller scales (associated with higher frequencies) set in.

Figure 16. Colour plots of the spanwise gradient of time-averaged streamwise velocity ${\partial \bar {u}}/{\partial z}$ and contour lines of $u^{\prime }_{rms}$ at the same streamwise locations as in figure 15.

The amplification of streaks and their unsteadiness induce rapid steepening of spanwise and wall-normal mean flow gradients, thereby leading to the emergence of inflection points in the mean (time- and spanwise-averaged) flow. Figure 17(a) shows the resulting inflectional mean flow profile and figure 17(b) shows the location of unsteady fluctuations with respect to the streaks. As we move from $x=65$ to $x=70$, we note the appearance of higher spanwise wavenumbers in the streaks as well as in the unsteady fluctuations. As discussed by Hall & Horseman (Reference Hall and Horseman1991) and Yu & Liu (Reference Yu and Liu1994), inflectional points in the mean velocity serve as an indicator of its susceptibility to the growth of broadband fluctuations. Amplification of high-frequency harmonics is also observed in the temporal spectra of the fluctuations. Therefore, at $x=70$, flow is well within fully nonlinear stages of transition. We also note the strong spatial correlation of unsteady fluctuations with the spanwise shear associated with the streaks, even at this nonlinear stage.

Figure 17. (a) Wall-normal profiles of laminar and mean streamwise velocity components at $x = 70$. (b) Colour plots of $u^{\prime }_{rms}$ along with contour lines of time-averaged streamwise velocity.

In figure 18 we plot $Q$-isosurfaces of the instantaneous flow field. These visualize vortical structures that arise from the sinuous-subharmonic motion of the streaks, prior to breakdown of the flow. We note the appearance of staggered lambda vortices, similar to those identified in transition induced by oblique waves in incompressible boundary layers (Berlin et al. Reference Berlin, Wiegel and Henningson1999). These vortical structures are also sometimes referred to as $X$-vortices, for reasons illustrated in the ($\eta,z$) slice array plot in figure 18. Initially, at $x \approx 62$, we observe that the interaction with the oblique waves causes the roll-up of the low-speed streak as they come together due to sinuous-subharmonic motion. Further downstream, as the flow structures associated with the $X$-vortices spread apart, we see that wall-normal jets of low-speed flow form mushroom structures. As these low-speed regions oscillate further, they interact to generate fluctuations with higher spanwise wavenumbers. Finally, at $x\approx 70$, the shear introduced by the upward jets of low-speed streaks becomes significant enough to cause a local Kelvin–Helmholtz instability that induces streamwise rollers (Reddy et al. Reference Reddy, Schmid, Baggett and Henningson1998). At this stage, the laminar boundary layer flow has started to break down to turbulence.

Figure 18. (a) Instantaneous isosurface of the $Q$-criterion coloured by contours of the streamwise velocity; (b) contour plot of the wall-shear stress show the formation of streaks before transition and $(\eta,z)$ slices of instantaneous streamwise velocity contours along with contour lines of the $Q$-criteria.

6.1. Wall statistics: skin friction

Boundary layer transition is characterized by a rapid increase in wall friction. Figure 19(a) utilizes instantaneous wall-shear distribution to identify transitional and turbulent regions. The wall-shear stress experiences sinuous-subharmonic modulation downstream of the reattachment (at $x=60$, i.e. $Re_{x} = 8.2 \times 10^{5}$) and finer spanwise scales emerge after $x=80$ (i.e. $Re_{x} \approx 11 \times 10^{5}$). As shown in figure 19(b), this location coincides with the highest value of time-averaged wall shear. Even though significant attenuation of spanwise variations of the wall-shear stress associated with the initial streaks takes place by $x = 80$, nonlinear interactions within the transition zone lead to the appearance of new streaks further downstream.

Figure 19. Streamwise variations of (a) instantaneous and (b) time-averaged streamwise shear stresses as well as (c) wall skin-friction coefficient.

The spatial extent of the transition zone is visualized in figure 19(c) by showing streamwise development of the skin-friction coefficient,

(6.1)\begin{equation} C_{f} = \frac{1}{T}\frac{1}{L_{z}}\int_{0}^{T} \int_{0}^{L_{z}} \frac{\tau^{*}_{w}}{\frac{1}{2} \rho^*_{e} (U_{e}^{*})^2}\,\mathrm{d} z \,\mathrm{d} t. \end{equation}

Here, $\tau ^*_{w}$ is the dimensional shear stress at the wall, $\rho ^*_{e}$ and $u^*_{e}$ are the dimensional density and streamwise velocity at the boundary layer edge, $L_{z}$ is the spanwise extent of the computational domain and $T = 4 L_{x}/u^*_{\infty }$. The values of $C_{f}$ in laminar flows over a $12$ degree wedge and over the double wedge are also plotted for comparison. The skin-friction coefficient drops because of flow separation but it grows again after reattachment. Near and downstream of reattachment, we observe a significant difference between laminar and turbulent values of $C_f$. After $Re_{x} = 11.5 \times 10^{5}$, when $C_f$ starts to decay with $x$, the wall-friction coefficient is approximately seven times larger than its laminar counterpart. A comparison with the Van-Driest turbulent correlation for standard compressible boundary layers demonstrates that the flow on the second wedge is approaching fully developed turbulent values towards the end of the computational domain. In addition to skin friction, transition also has a significant impact on wall temperature. Additional details about the mean temperature and the wall statistics are included in Appendix F.

6.2. Towards a turbulent compressible boundary layer

To illustrate the onset of turbulence, we also examine first- and second-order statistics in the latter stages of transition. At a given streamwise location, we report statistics in terms of the inner coordinate $\eta ^{+}$ that is obtained by non-dimensionalizing the wall-normal coordinate with the local viscous length $\delta _{\nu }$. Figure 20(a) shows the mean streamwise velocity $u^{+}$ non-dimensionalized by the local friction velocity $u_{\tau }$. The mean profile is obtained by averaging in time over $2 L_{x}/u_{\infty }$ and in the spanwise direction over the extent of the computational domain $L_{z}$. In the presence of density variation along the compressible boundary layer, we utilize the Van-Driest transformation,

(6.2)\begin{equation} u_{c}(\eta) = \int_{0}^{\langle{u}\rangle} \sqrt{\dfrac{\langle{\rho}\rangle}{\langle{\rho}\rangle_{w}}} \,\mathrm{d} \hat{u} = \int_{0}^{\eta}\sqrt{\dfrac{\langle{\rho}\rangle}{\langle{\rho}\rangle_{w}}} \, \frac{\partial \! \langle{u}\rangle}{\partial \hat{\eta}} \, \,\mathrm{d}\hat{\eta}, \end{equation}

to compare the transitional velocity profiles with the ‘log law’ observed in incompressible turbulent boundary layers (White & Majdalani Reference White and Majdalani2006), where $\langle \,{\cdot }\,\rangle$ denotes averaging over time and spanwise direction, and $\langle {\rho }\rangle_{w}$ is the mean density at the wall. A fully developed turbulent boundary layer is characterized by the presence of the viscous sublayer close to the wall ($\eta ^{+} < 10$), where ${u^+_{c}} = \eta ^+$, and further away from the wall the mean velocity obeys the log law,

(6.3)\begin{equation} {u^+_{c}} =(1/\kappa) \operatorname{ln}(\eta^{+}) + C, \end{equation}

where $\kappa = 0.41$ and $C = 5.2$.

Figure 20. The wall-normal profiles of mean (a) streamwise velocity and (b) temperature.

Figure 20(a) shows that the mean velocity in the boundary layer approaches the fully turbulent profile at $Re_{x} = 1.2\times 10^6$. Upstream of this location, within the transition zone, the boundary layer has a significantly greater momentum and most of it lies away from the wall. Furthermore, closer to the wall, the slope of the streamwise velocity profile is larger than the slope obtained using the viscous sublayer profile of a fully turbulent flow. This observation is consistent with the overshoot in skin-friction coefficient within the transition zone; cf. figure 19(b).

In addition to the mean velocity, the mean temperature profile plays an important role in heat transfer and material response analysis of hypersonic flows. Walz's modified Crocco-Busemann relation,

(6.4)\begin{equation} \frac{\langle{T}\rangle}{\langle{T_{e}}\rangle} = \frac{\langle{T_{w}}\rangle}{\langle{T_{e}}\rangle} + \frac{\langle{T_{r}\rangle}-\langle{T_{w}}\rangle}{\langle{T_{e}}\rangle} \frac{\langle{u}\rangle}{\langle{u_{e}}\rangle} + \frac{\langle{T_{e}\rangle}-\langle{T_{r}}\rangle}{\langle{T_{e}}\rangle} \left( \frac{\langle{u}\rangle}{\langle{u_{e}}\rangle} \right)^{2}, \end{equation}

is commonly used to describe the relation between temperature and velocity in a zero pressure gradient turbulent boundary layer. Here, $\langle {u_{e}}\rangle$ and $\langle {T_{e}}\rangle$ denote the mean boundary layer edge velocity and temperature, respectively, $\langle {T_{r}}\rangle$ is the mean recovery temperature and, since the wall is adiabatic, the mean wall temperature is determined by $\langle {T_{w}}\rangle = \langle {T_{r}}\rangle$. In spite of pressure variations post-reattachment, figure 20(b) demonstrates that the quadratic relation in (6.4) holds throughout the transition zone. As the flow approaches a fully turbulent profile, we see a slight deviation from this relation in the outer region away from the wall; this observation is in agreement with prior studies of fully turbulent compressible boundary layers (Duan, Beekman & Martín Reference Duan, Beekman and Martín2011; Franko & Lele Reference Franko and Lele2013).

We next examine spatial development of flow fluctuations by evaluating second-order statistics and comparing them with those observed in a Mach 5 fully developed turbulent flat-plate boundary layer (Duan et al. Reference Duan, Beekman and Martín2011). Figure 21(a–c) shows the streamwise variation of the density-weighted r.m.s. values of the streamwise, wall-normal and spanwise velocity fluctuation components. In the transition zone we observe large values of fluctuations away from the wall in all three plots. Further downstream, at $Re_{x} = 1.2 \times 10^6$, the r.m.s. values closer to the wall are in agreement with fully turbulent values. However, away from the wall, the r.m.s. values of $u'$ and $w'$ deviate from the flat-plate profiles. We attribute this deviation to the presence of 3-D oblique waves and streaks that persist downstream because of continuous upstream forcing in our set-up.

Figure 21. The wall-normal profiles of density weighted r.m.s. values of (a) streamwise, (b) wall-normal and (c) spanwise velocity fluctuations at different Reynolds numbers. (d) Energy spectrum as a function of spanwise wavenumber $\beta$ at $Re_{x} = 1.2 \times 10^6$.

To illustrate the broadband nature of velocity fluctuations in the later stages of transition, in figure 21(d) we evaluate the spanwise wavenumber dependence of the individual contributions of velocity fluctuations to the one-dimensional energy spectrum at $Re_{x} = 1.2 \times 10^{6}$ and $\eta ^{+}=40$. At this location, the energy spectrum of streamwise velocity fluctuations clearly exhibits inertial and dissipative subrange features, indicating that the flow is approaching a fully turbulent stage.

Our DNS also provides data for evaluating non-dimensional parameters which can be utilized for low-complexity modelling of high-speed compressible turbulent flows using time-averaged NS equations (Hirsch Reference Hirsch2007). In particular, we examine the spatial evolution of the fluctuation Mach number, $M_{t}$, the fluctuation Prandtl number, $Pr_t$, and Huang's modified strong Reynolds analogy parameter (Huang, Coleman & Bradshaw Reference Huang, Coleman and Bradshaw1995), $\psi _{sra}$,

(6.5a–c)\begin{align} M_t = \sqrt{\frac{{\langle{u_i^\prime u_i^\prime}}\rangle}{\gamma \,R \langle{T}\rangle}}, \quad Pr_t = \frac{\langle{ u^\prime v^\prime}\rangle} {\langle{ v^\prime T^\prime}\rangle} \frac{\partial_{\eta} \! \langle{T}\rangle}{\partial_{\eta} \! \langle{u}\rangle},\quad \psi_{{sra}} = \frac{\left(1-\dfrac{\partial \langle{T_{0}}\rangle}{\partial \langle{T}\rangle}\right) \, Pr_t}{(\gamma-1)\langle{M}\rangle^2}\,\dfrac{\langle{u}\rangle}{\langle{T}\rangle}\,\frac{T'_{rms}}{u'_{rms}}. \end{align}

Here, $\langle {u_i^\prime u_i^\prime }\rangle = \langle {u^\prime u^\prime }\rangle + \langle {v^\prime v^\prime }\rangle + \langle {w^\prime w^\prime }\rangle$, $\gamma$ and $R$ denote the specific heat ratio and the gas constant, respectively, and $T_{0}$ is the stagnation temperature. In a fully developed turbulent boundary layer, $\psi _{sra}$ provides a measure of correlation between velocity and temperature fluctuations (Huang et al. Reference Huang, Coleman and Bradshaw1995).

Figure 22 shows that towards the end of the computational domain, at $Re_x= 1.2\times 10^{6}$, these parameters are close to those observed in canonical hypersonic boundary layers, i.e. $M_t \approx 0.2-0.3$, $Pr_{t} \approx 0.9$ and $\psi _{{sra}} \approx 1$ (Pirozzoli & Bernardini Reference Pirozzoli and Bernardini2011). However, upstream of the breakdown region, at $Re_x= 1.0\times 10^{6}$, there is a significant deviation compared with these canonical values. Here, $M_t$ becomes as high as $0.45$ which suggests that compressibility effects on flow fluctuations cannot be neglected in the transition zone. Similarly, $Pr_t$ can reach values close to $1.5$ that correspond to decreased fluctuation temperature fluxes in the breakdown region. Furthermore, $Pr_t$ exhibits large variations away from the wall, which is in contrast to the observations in flat-plate turbulent boundary layers, where $Pr_t$ has a value of $0.9$ throughout the boundary layer (Saffman & Wilcox Reference Saffman and Wilcox1974; Smith & Smits Reference Smith and Smits1993; Pirozzoli & Bernardini Reference Pirozzoli and Bernardini2011). Similar deviations from canonical values in the outer region of the boundary layer are also observed in $\psi _{sra}$. We conjecture that, in the presence of persistent upstream excitations, the resulting unsteady fluctuations are primarily responsible for these discrepancies.

Figure 22. The wall-normal profiles of (a) the fluctuation Mach number, $M_t$; (b) the fluctuation Prandtl number, $Pr_t$; and (c) Huang's modified strong Reynolds analogy parameter, $\psi _{sra}$.