2.1. Flow configuration

In the present study, an oblique shock wave/laminar-boundary-layer interaction is considered. The flow configuration is displayed in figure 1. To enable a comprehensive analysis of the unsteadiness of the structures developing in such a flow, a set of flow cases are run with parameters referenced in table 1. Hereafter, the coordinate system is the following: $x$ is oriented in the streamwise direction and $y$ is associated with the direction normal to the wall. Also $L_{sep}$ and $L_{int}$ the separation and interaction lengths, respectively.

Figure 1.

The computational domain and a schematic representation of oblique shock wave/boundary-layer interaction.

Table 1.

Flow parameters for laminar OSWBLIs.

2.2. Governing equations

The two-dimensional compressible Navier–Stokes equations for a compressible perfect gas are considered. These equations govern the evolution of the system state $\boldsymbol{q}=[{\it\rho},{\it\rho}\boldsymbol{u},{\it\rho}E]^{\text{T}}$ in the conservative form, where ${\it\rho}$ , $\boldsymbol{u}$ and $E$ are the fluid density, the velocity vector and the total energy, respectively. Written in non-dimensional form, these equations are

For an ideal and Newtonian fluid, the non-dimensional pressure $p$ and energy $E$ are related to the temperature through the equation of state for an ideal gas.

${\bf\tau}$ is the viscous stress tensor and is written as

The dynamics viscosity ${\it\mu}(T)$ is computed using Sutherland’s law, i.e.,

where $T_{s}=110.4~\text{K}$ . The coefficient of heat conductivity ${\it\kappa}(T)$ is given in terms of the Prandtl number ${\it\kappa}(T)={\it\mu}(T)/Pr$ . Here, $M_{\infty }$ is the free-stream Mach number at the inflow and $Re_{{\it\delta}^{\star }}$ is the Reynolds number based on the free-stream velocity downstream of the separated shock, $U_{e}$ , and the boundary layer displacement thickness at the theoretical position where the incident shock impinges the boundary-layer, ${\it\delta}^{\star }$ . The time scales are normalized with ${\it\delta}^{\star }/U_{e}$ and pressure with ${\it\rho}_{\infty }U_{e}^{2}$ . Finally, the specific heat capacity ratio is equal to ${\it\gamma}=1.4$ and the Prandtl number is taken as $Pr=0.72$ .

In the following, the system (2.1) can be recast in the formal conservative form:

where ${\mathcal{R}}$ is the differential nonlinear Navier–Stokes operator.

2.3. Compressible global stability analysis

Linear stability analysis assumes the existence of an equilibrium solution to the system (2.5) referred to as the base flow and defined by ${\mathcal{R}}(\boldsymbol{q}_{b})=0$ . In the following, the base flow is assumed to be two-dimensional $\boldsymbol{q}_{b}(x,y)$ . Using the standard small-perturbation technique, the instantaneous flow is decomposed into a base flow and a small disturbance:

The resulting equations are further simplified by considering that the perturbation is infinitesimal, i.e. the nonlinear fluctuating terms are neglected. Finally, the compressible Navier–Stokes equations are transformed into linear partial differential equations

where the vector $\boldsymbol{q}^{\prime }=({\it\rho}^{\prime },{\it\rho}^{\prime }\boldsymbol{u}_{b}+{\it\rho}_{b}\boldsymbol{u}^{\prime },{\it\rho}^{\prime }E_{b}+{\it\rho}_{b}E^{\prime })^{\text{T}}$ represents the conservative variables for the perturbation and ${\mathcal{J}}=\partial {\mathcal{R}}/\partial \boldsymbol{q}|_{\boldsymbol{q}_{b}}$ is the Jacobian operator which corresponds to the linearization of the Navier–Stokes operator ${\mathcal{R}}$ around the base flow $\boldsymbol{q}_{b}$ .

where the perturbed stress tensor is given by

The perturbed pressure, total energy, dynamic viscosity and heat conductivity coefficient are written

3.1. Navier–Stokes solver and boundary conditions

All numerical simulations in this paper were run with an in-house CFD solver named Phoenix, both to compute the base flow and to solve the linearized and the full Navier–Stokes equations. The numerical method implemented in our solver is based on the finite-volume approach and on a cell-centred discretization. The code solves the compressible Navier–Stokes equations on multi-block structured grids and these equations are discretized in space using an upwind scheme. Roe’s flux difference splitting scheme (Roe Reference Roe1981) is employed to obtain advective fluxes at the cell interface for all N–S equations. The MUSCL approach extends the spatial accuracy to third order. All viscous terms are centrally differentiated. For unsteady computations, the dual time stepping method, proposed by Jameson (Reference Jameson1991), was used to tackle the lack of numerical efficiency of the global time stepping approach. The derivative with respect to the physical time is discretized by a second-order formula. The use of an implicit scheme with respect to the dual time provides fast convergence to the time-accurate solution. Between each time step, the solution is advanced in a dual time and acceleration strategies developed for steady problems can be used to speed up the convergence in fictitious time. The implicit time-integration procedure leads to a system which can be solved directly or iteratively. The direct inversion can be memory intensive and computationally expensive. Therefore, an implicit relaxation procedure is preferred and the point Jacobi relaxation algorithm was chosen. The initialization of the derivative with respect to the physical time was performed with a first-order formula.

The boundary conditions applied to the computational domain are no-slip velocity, adiabatic temperature and pressure extrapolation at the wall. The steady two-dimensional Navier–Stokes solution is imposed at the inflow. To minimize the reflection of waves into the domain at the outflow and at the upper boundary, a characteristic method is used for the conservative variables. The oblique shock is imposed by the Rankine–Hugoniot jump relations at the top boundary.

3.2. Linearization of discrete Navier–Stokes equations

Equations (2.8) are discretized with the same type of scheme as the nonlinear N–S equations (2.1). However, the spatial schemes as well as boundary conditions have to be adapted to the new set of equations (2.8). To achieve such a task, the different schemes and boundary conditions are built by considering the flux function associated with the linearized equations (2.8). The Roe scheme is adopted in this study. As suggested by Crouch, Garbaruk & Magidov (Reference Crouch, Garbaruk and Magidov2007), the Roe scheme adapted for the linearized compressible N–S variables is based on the Jacobian matrix of the new flux function associated with the linearized equations, as discussed above. The latter matrix being only a function of the base flow, there is no need for a Roe average. In addition, the boundary conditions are also only driven by the base flow, yielding similar modifications of our CFD solver. Finally the semi-discrete system can be written as

where $\unicode[STIX]{x1D645}$ is $(n\times n)$ Jacobian matrix and $n=\text{dim}(\boldsymbol{q}^{\prime })\times n_{x}\times n_{y}$ .

3.3. Stability algorithm

To determine the eigenmodes of discrete linearized Navier–Stokes equations (3.1) the solution is sought under the form of a normal mode $\boldsymbol{q}^{\prime }(\boldsymbol{x},t)=\widehat{\boldsymbol{q}}(\boldsymbol{x})\exp ({\it\lambda}t)$ , the problem (3.1) becomes an eigenproblem $\unicode[STIX]{x1D645}\widehat{\boldsymbol{q}}={\it\lambda}\widehat{\boldsymbol{q}}$ where ${\it\lambda}={\it\sigma}+\text{i}{\it\omega}$ is the eigenvalue and $\widehat{\boldsymbol{q}}$ is the eigenvector. The sign of the leading eigenvalue’s real part ${\it\sigma}$ then determines whether the fixed point $\boldsymbol{q}_{b}$ is linearly stable or unstable, whereas its imaginary part ${\it\omega}$ characterizes the stationary or oscillatory nature of the associated eigenvector. To solve this eigenproblem, two types of methods can be used: (i) methods where the Jacobian matrix is explicitly formed or (ii) methods wherein the matrix is not directly and fully calculated, these methods are known as matrix-free methods. In the present configuration, the size of our problem ( $\unicode[STIX]{x1D645}$ would be approximately a $10^{5}\times 10^{5}$ matrix) is not limiting and both methods are possible. The matrix-free method is selected where snapshots of the solution of (3.1) are used to approximate the Jacobian matrix. In the compressible regime, this method has a number of advantages because it is written in a fully-discrete formalism allowing to better take into account the complex boundary conditions such as the presence or the emission of acoustic waves while maintaining an accurate description of the hydrodynamic instability. Using a snapshot method based on the discrete equation (3.1), we solve the following eigenvalue problem $\unicode[STIX]{x1D648}({\rm\Delta}t_{s})\boldsymbol{q}^{\prime }=\widehat{{\it\lambda}}\boldsymbol{q}^{\prime }$ , where $\unicode[STIX]{x1D648}({\rm\Delta}t_{s})=\exp (\unicode[STIX]{x1D645}{\rm\Delta}t_{s})$ is the exponential propagator of the system (3.1). It appears at first sight that this so-called time stepper approach (Edwards et al. Reference Edwards, Tuckerman, Friesner and Sorensen1994; Bagheri et al. Reference Bagheri, Akervik, Brandt and Henningson2009) does not really simplify the initial problem of memory footprint. Indeed, not only is $\exp (\unicode[STIX]{x1D645}{\rm\Delta}t_{s})$ a matrix of very large dimensions, but it is well known that computing a matrix exponential can be quite a challenging computational task. The most interesting feature of such a time-stepper approach however is that, though $\exp (\unicode[STIX]{x1D645}{\rm\Delta}t_{s})$ is as computationally expensive to construct explicitly as the Jacobian $\unicode[STIX]{x1D645}$ , the action of this exponential propagator on a vector $\boldsymbol{q}_{0}^{\prime }$ can easily be approximated by simply time marching the linearized Navier–Stokes equations (3.1) with $\boldsymbol{q}_{0}^{\prime }$ as the initial condition. The eigenmodes of $\unicode[STIX]{x1D648}$ are the same as those of the matrix $\unicode[STIX]{x1D645}$ if ${\rm\Delta}t_{s}$ is chosen in accordance with the Nyquist criterion. The dataset length $N_{s}{\rm\Delta}t_{s}$ is sufficient to capture the slowest frequency of the system by a factor 20 (see table 2). The iterative technique is based on the orthogonal projection of the large matrix $\unicode[STIX]{x1D648}$ onto a lower-dimensional Krylov subspace which is spanned by snapshots taken from the flow field $\boldsymbol{q}^{\prime }$ separated by a constant time interval ${\rm\Delta}t_{s}$ . The resulting system is a significantly smaller system and can be solved using classical direct or iterative methods:

To increase the conditioning of the Krylov basis, a Gram–Schmidt orthogonalisation process is used. The exponential propagator matrix $\exp (\unicode[STIX]{x1D645}{\rm\Delta}t_{s})$ is then projected on this unitary basis $\boldsymbol{U}$ :

where $\boldsymbol{U}_{k}$ is an orthonormal set of vectors, $\unicode[STIX]{x1D643}_{k}$ a $k\times k$ upper Hessenberg matrix and $\boldsymbol{r}_{k}\boldsymbol{e}_{k}^{\text{T}}$ is the residual vector indicating how far from an invariant subspace of $\unicode[STIX]{x1D648}$ is $\boldsymbol{U}_{k}$ . Because of its relatively small dimension, the eigenpairs ( $\widehat{{\it\lambda}}_{H},\boldsymbol{X}$ ) of the Hessenberg matrix can be directly computed which are a good approximation to those of $\unicode[STIX]{x1D648}$ . Since this Hessenberg matrix is a low-dimensional approximation of the exponential propagator $\unicode[STIX]{x1D648}=\exp (\unicode[STIX]{x1D645}{\rm\Delta}t_{s})$ , these Ritz pairs are linked to the eigenpairs of the Jacobian matrix $\unicode[STIX]{x1D645}$ by:

The easiest way to solve the eigenvalue problem for a matrix-free method is (i) to force the system (3.1) by an impulse-shaped parietal forcing; (ii) to generate the set of snapshots; and (iii) to solve the eigenproblem by an Arnoldi method. Although this method is used to calculate with good accuracy the unstable modes, the stable part of the spectrum is more difficult to obtain. When the flow is globally stable (selective noise amplifier behaviour) this problem becomes critical. The origin of the problem is two-fold. For the open flows, the linearized evolution operator is often highly non-normal, the eigenvalues are then highly sensitive to perturbations (Trefethen & Embree Reference Trefethen and Embree2005; Schmid Reference Schmid2007). In addition, a snapshot method, which is similar to a signal processing method, is ill conditioned for temporally stable flows because over time the modes fade, making it more difficult to extract their frequency content. To solve this problem, the computation of the matrix $\unicode[STIX]{x1D643}$ from numerical simulations is accomplished by an in-house Arnoldi method (Arnoldi Reference Arnoldi1951; Lehoucq, Sorensen & Yang Reference Lehoucq, Sorensen and Yang1997; Barkley, Blackburn & Sherwin Reference Barkley, Blackburn and Sherwin2008) directly integrated into the temporal loop of our CFD solver (see Loiseau et al. Reference Loiseau, Robinet, Cherubini and Leriche2014).

4.1. Base flows

Seven flow cases will be considered hereafter, referenced A–G in table 3. For all flow cases, a preliminary laminar boundary-layer flow computation without an impinging shock is carried out. Then we denote by ${\it\delta}^{\star }$ the displacement thickness at the position of impact of the incident shock. For flow case A, the flow conditions are matched to the experimental and numerical investigation of Degrez et al. (Reference Degrez, Boccadoro and Wendt1987). The free-stream Mach number is $M_{\infty }=2.15$ and the angle of the incident shock measured clockwise from the horizontal axis is ${\it\theta}=30.8^{\circ }$ (corresponding to a flow deflection angle of $3.81^{\circ }$ ). The Reynolds number based on ${\it\delta}^{\star }$ is fixed to $Re_{{\it\delta}^{\star }}=1050$ . The computational domain ${\mathcal{D}}$ has dimensions $[x_{0}/{\it\delta}^{\star };x_{n}/{\it\delta}^{\star }]\times [y_{0}/{\it\delta}^{\star };y_{n}/{\it\delta}^{\star }]=[0;400]\times [0;80]$ in the streamwise and wall-normal directions, respectively. The total number of points that are clustered in the normal-to-wall direction, as well as in the interaction region, is $(N_{x}\times N_{y})=(560\times 100)$ . Numerical validations (not shown here for the sake of conciseness) have proven that this guarantees grid independency. The dual time step is switched off and the CFL (Courant–Friedrichs–Lewy) number is held to 10.

Table 2.

Numerical set-up for the different stability computations. ${\rm\Delta}t$ is the time stepping used in linearized Navier–Stokes solver, ${\rm\Delta}t_{s}$ is the time between two consecutive snapshots, $N_{s}$ is the number of snapshots, ${\it\varepsilon}_{ev}$ is the minimal eigenvalues convergence.

Table 3.

Flow parameters for the OSWBLI.

After having reached a steady state for the laminar boundary-layer flow, the incident shock is introduced at the inflow condition through the Rankine–Hugoniot relations. The latter shock impinges the boundary-layer, which creates a separated flow. It thus leads to compression waves behind the bubble which coalesce to form the reflected or separated shock. The converged base flow is shown in figure 2(a) in terms of streamwise velocity fields and streamlines. Comparisons with both experimental and numerical results given by Degrez et al. (Reference Degrez, Boccadoro and Wendt1987) are displayed in figure 3. The latter figure shows that our simulation is in good agreement with the mean-flow characteristics reported by Degrez et al. (Reference Degrez, Boccadoro and Wendt1987), both in terms of friction coefficient and wall pressure measurements.

Figure 2.

Streamlines and streamwise velocity component of the computed steady base flow for different values for Reynolds and Mach numbers and incident angles. The incident shock (1), the reflected shock (2), the expansion fan (3) and compression waves (4) are also represented. (a) $Re_{{\it\delta}^{\star }}=1050$ , $M_{\infty }=2.15$ and ${\it\theta}=30.8^{\circ }$ , (b) $Re_{{\it\delta}^{\star }}=1050$ , $M_{\infty }=2.15$ and ${\it\theta}=32^{\circ }$ , (c) $Re_{{\it\delta}^{\star }}=1090$ , $M_{\infty }=2.20$ and ${\it\theta}=30.8^{\circ }$ , (d)  $Re_{{\it\delta}^{\star }}=1950$ , $M_{\infty }=2.15$ and ${\it\theta}=30.8^{\circ }$ .

Figure 3.

Flow case A. (a) Friction coefficient and (b) wall pressure distribution normalized by $P_{\infty }$ shown as a function of the streamwise position scaled by ${\it\delta}^{\star }$ . Full line: our computation. Full circle, the results given by Degrez et al. (Reference Degrez, Boccadoro and Wendt1987). (a) The dashed-line is associated with the friction coefficient of the laminar boundary layer without the impinging incident shock.

In order to investigate the influence of various physical parameters (bubble shapes, $M_{\infty }$ and $Re_{{\it\delta}^{\star }}$ ) on stability properties, equilibrium states associated with B–G are also computed (see figure 2 b–d).

Figure 4(a,b) show the distribution of skin-friction coefficients and dimensionless pressure, respectively, as a function of the streamwise position $x$ for ${\it\theta}$ varying from $30.8^{\circ }$ to $32.5^{\circ }$ when the Reynolds number and the Mach number are fixed to 1050 and 2.15, respectively. It is observed that when the incident angle increases, the interaction zone increases almost linearly. Figure 4(a) also shows that for an angle close to $32.5^{\circ }$ (computed for the same flow conditions in terms of Mach and Reynolds numbers), there is a critical point (around $x/{\it\delta}^{\star }\simeq 99.2$ ) where the skin-friction coefficient becomes positive. In an incompressible regime, Theofilis et al. (Reference Theofilis, Hein and Dallmann2000) and Rodriguez & Theofilis (Reference Rodriguez and Theofilis2010) have shown that the latter behaviour is associated with a bifurcation of the flow towards a three-dimensional state. This behaviour has also been observed for supersonic separated flows by Boin et al. (Reference Boin, Robinet, Corre and Deniau2006) and Robinet (Reference Robinet2007). Hence, the figure 4(a) suggests that ${\it\theta}=32.5^{\circ }$ is a relevant upper limit when considering a two-dimensional equilibrium state. Figure 4(b) shows that the pressure ratio reaches the value provided by the Rankine–Hugoniot relationship downstream of the reflected shock for the different incident angles ${\it\theta}$ studied. Hereafter, the stability analysis will be limited to the incident angles ${\it\theta}\leqslant 32^{\circ }$ . Similar behaviour is observed when the Reynolds number and/or Mach number are increased.

Figure 4.

$M_{\infty }=2.15$ , $Re_{{\it\delta}^{\star }}=1050$ . Distributions of (a) skin-friction coefficient and (b) dimensionless pressure with the streamwise position. The values provided by the Rankine–Hugoniot relationship are also provided as horizontal full lines labelled R–H. Full line: ${\it\theta}=30.8^{\circ }$ , long dashed line: ${\it\theta}=31.4^{\circ }$ , dash-dot line: ${\it\theta}=32^{\circ }$ and dashed line: ${\it\theta}=32.5^{\circ }$ .

4.2. Global mode analysis and linear dynamics

4.3. Linear regime: noise amplifier dynamics and receptivity

In the previous sections, we have shown that OSWBLI is globally stable for a wide range of angles of the incident shock, and of Reynolds and Mach numbers. The flow behaviour will thus depend on the receptivity (i.e. associated with both the level and the shape of the external noise; and the way it affects the flow) and a frequency-selective noise amplification of the external perturbation due to non-orthogonal interactions between the modes (see Schmid (Reference Schmid2007) for a review).

To give insight into the selection frequency process of such a flow, we investigate the linear flow response to a white noise excitation (i.e. the linear version of the CFD code is used). For this purpose, a spatially localized forcing near the wall is chosen as actuator. In particular, the latter has a Gaussian envelope:

where $x_{f}^{\star }=x_{f}/{\it\delta}^{\star }$ and $y_{f}^{\star }=x_{f}/{\it\delta}^{\star }$ fix the forcing location and $A_{f}$ accounts for its magnitude. Its spatial extent along the flat plate is controlled with ${\it\xi}_{x}$ and ${\it\xi}_{y}$ . $W_{v}(t)$ is a white noise signal varying between $[-0.5;0.5]$ . The actuator is placed at two locations: one starting upstream of the separation point and a second one localized at the centre of the recirculation zone (referenced as $F1$ and $F2$ respectively). Four pressure probes are placed at several abscissa to measure the response of the flow. The various parameters are specified in table 4. The flow response is investigated through the analysis of frequency-weighted power spectral densities (WPSDs) at various pressure measurement points. The power spectral density (PSD) is estimated using the overlapping segments technique. Three segments with 50 % overlap are used over pressure time series data with a total record length 1000. The WPSDs associated with the flow response at the separation point ( $P_{1}$ ), far from the wall close to the separated shock ( $P_{2}$ ) and inside the bubble ( $P_{3}$ ) are shown in figures 15 and 16 (for $F1$ and $F2$ , respectively) as a function of the non-dimensional frequency $St$ . The choice of the position of the probe $P_{2}$ allows us to evaluate the influence of the separation shock with respect to the interaction zone. For the upstream broadband white noise forcing, $F1$ , the WPSDs exhibit a broadband amplification frequency centred at a non-dimensional frequency in the range of $St\in [0.05;0.5]$ for $P_{1}$ and $P_{3}$ . This selectivity is mainly due to the shear layer. Similar results are obtained for $P_{4}$ (not shown here). In addition, figure 15(b) shows that the measurement associated with probe $P_{2}$ exhibits a low-frequency broadband energy content. In this position, the linearized dynamics presents no particular selectivity and beyond $St\approx 0.1$ , all frequencies are attenuated. Different grids has been tested, showing that these features are independent of the numerical choices.

Table 4.

Linear and nonlinear regimes: forcing parameters and measurement coordinates.

Figure 15.

$M_{\infty }=2.15$ , ${\it\theta}=30.8^{\circ }$ and $Re_{{\it\delta}^{\star }}=1050$ (Case A). Linear regime. Frequency response to upstream broadband white noise forcing ( $F1$ ) at the (a) separating point; (b) separated shock; (c) downstream position. The WPSDs normalized by the forcing amplitude are represented as a function of the dimensionless frequency ( $St$ ). (a) $P_{1}$ – separating point, (b) $P_{2}$ – separated shock, (c) $P_{3}$ – separated zone.

Figure 16.

$M_{\infty }=2.15$ , ${\it\theta}=30.8^{\circ }$ and $Re_{{\it\delta}^{\star }}=1050$ (Case A). Linear regime. Frequency response to internal broadband white noise forcing ( $F2$ ) at the (a) separating point; (b) separated shock; (c) downstream position. The WPSDs normalized by the forcing amplitude are represented as a function of the dimensionless frequency ( $St$ ). (a) $P_{1}$ – separating point, (b) $P_{2}$ – separated shock, (c) $P_{4}$ – downstream point.

For the internal broadband white noise forcing, $F2$ , the WPSDs have a low-frequency broadband content with no particular selectivity for $P_{1}$ and $P_{2}$ . $P_{3}$ and $P_{4}$ display a selective-frequency response similar to the $F1$ case. In addition, it is interesting to note that when the forcing is localized near the separation point, the amplitude of the response is larger than the one associated with forcing $F2$ (i.e. inside the bubble). These results are in agreement with those of Sansica et al. (Reference Sansica, Sandham and Hu2014), only the frequency range involved is different and can be largely explained by an effect of the Mach number (1.5 in their case). To characterize the spatial shape of structures for a single frequency associated with the response to a white noise forcing, 10 000 every ${\rm\Delta}t{\it\delta}^{\star }/U_{e}=0.148$ of streamwise velocity and pressure fluctuation snapshots are stored and Fourier transformed. In figure 17(a,b), Fourier modes corresponding to $St\approx 0.145$ and $St\approx 0.056$ are represented for case A. On the one hand, figure 17(a) shows that most amplified frequencies are associated with not radiating Kelvin–Helmholtz modes reminiscent of subsonic Kelvin–Helmholtz global modes. On the other hand, low-frequency Fourier modes exhibit a spatial amplification of vortices along the separated zone with emission of an acoustic wave into the expansion fan, similar to supersonic Kelvin–Helmholtz global modes. To further illustrate this, the Fourier modes are projected into a set of global modes using an orthogonal projection that relies on Gram–Schmidt orthonormalization. The procedure is detailed in Song et al. (Reference Song, Alizard, Robinet and Gloerfelt2013). Hence, a single Fourier mode may be decomposed as

with $W_{i}$ being the coefficient of projection onto the global mode $\hat{\boldsymbol{q}}_{i}$ . Hereafter, $N=100$ global modes are considered. The weighted spectra are shown in figure 18(a,b). The figures show that, while the Fourier mode associated with $St\approx 0.145$ involves mainly global modes associated with regions II and III, the global modes of region I concentrate more energy for the Fourier mode corresponding to $St\approx 0.056$ , consistent with the acoustic radiation observed in figure 17(a). Consequently, it is clear from figure 17(a,b) that Fourier modes may be considered to be a train of wavepackets beating to a single frequency given by the Fourier transform. In particular, the latter wavepackets are composed of many frequency components, each of them being associated with a global mode. In addition, figure 17(a,b) also show the consequence of dispersive effects associated with boundary-layer modes and Kelvin–Helmholtz modes (region I, II and III, respectively) that travel with different phase velocities. In particular, one may observe that the wavepacket spreads out for low frequencies in comparison to medium frequencies (as depicted in figure 17 a,b). From the discussion above, it seems clear that OSWBLI filters and amplifies a specific range of frequencies, characteristic of a selective noise amplifier (Sipp et al. Reference Sipp, Marquet, Méliga and Barbagallo2010) and also give a physical explanation of individual global modes as a wavepacket component. The Fourier analysis shows that actuation in the vicinity of the separation provides a more efficient response for the most amplified frequencies than an actuation inside the separation zone. Hence, one may conclude that subsonic Kelvin–Helmholtz global modes have greater receptivity near the separation point whereas supersonic Kelvin–Helmholtz global modes are more receptive to a forcing localized inside the bubble. In addition, the predominance of supersonic Kelvin–Helmholtz global modes for Fourier modes with $St\approx 0.056$ provide physical insight useful to the interpretation of the amplification of lower frequencies measured with $P_{2}$ as shown in figures 15(b) and 16(b). Indeed, in the region near the probe $P_{2}$ , the flow is mainly driven by acoustic disturbances which are associated with supersonic Kelvin–Helmholtz modes.

Figure 17.

$M_{\infty }=2.15$ , ${\it\theta}=30.8^{\circ }$ and $Re_{{\it\delta}^{\star }}=1050$ (Case A). Linear regime. Fourier analysis: contours of the pressure fluctuations for (a) mode with $St\approx 0.145$ (i.e. $St\times L_{int}/{\it\delta}^{\star }=5.52$ ) and (b) mode with $St\approx 0.056$ (i.e. $St\times L_{int}/{\it\delta}^{\star }=2.05$ ). In the subfigures, the corresponding spanwise vorticity is represented.

Figure 18.

$M_{\infty }=2.15$ , ${\it\theta}=30.8^{\circ }$ and $Re_{{\it\delta}^{\star }}=1050$ (Case A). Linear regime. Contributions of individual global modes to (a) the single Fourier mode associated with $St\approx 0.145$ (see figure 17 a) and (b) the single Fourier mode associated with $St\approx 0.056$ (see figure 17 b). The circles are proportional to the corresponding expansion coefficients associated with the projection of the Fourier modes (referenced as $W_{i}$ in (4.4)).

The frequency selection observed in figures 15(a,c) and 16(a,c), can be explained in a global framework through the analysis of the global resolvent, ${\mathcal{R}}({\it\omega})$ (see Sipp et al. (Reference Sipp, Marquet, Méliga and Barbagallo2010) for the technical details). The latter quantity represents the maximum response of the flow due to a harmonic forcing at a frequency ${\it\omega}$ and is strongly associated with the non-orthogonality of the global modes. ${\mathcal{R}}({\it\omega})$ can be approximated by using the Hessenberg matrix from the Arnoldi algorithm (see Toh & Trefethen Reference Toh and Trefethen1996). The relevancy of our approach is discussed in § A.2. In figure 19, the distribution of ${\mathcal{R}}$ is shown as a function of the non-dimensional frequency $St$ for flow cases A and C. Apart from the amplification of a wide range of frequencies in between, the latter curve exhibits a maximum for $St\approx 0.145$ consistent with the linear simulation. In addition, the most amplified Strouhal number is seen to be not affected by an increase of ${\it\theta}$ from 30.8° to 32°. This result is in agreement with the scaling for the global spectra analysed in the previous section.

Figure 19.

Distribution of the global resolvent as a function of the dimensionless frequency $St$ for case A (full line) and case C (dashed lines). The maximum amplitude of ${\mathcal{R}}$ is defined up to a multiplicative constant.

4.4. Influence of nonlinearities

Figure 20.

$M_{\infty }=2.15$ , ${\it\theta}=30.8^{\circ }$ and $Re_{{\it\delta}^{\star }}=1050$ (Case A). Frequency response to upstream broadband white noise forcing for nonlinear Navier–Stokes equations. (a) $P_{1}$ – separating point, (b) $P_{2}$ – separated shock.

Figure 21.

$M_{\infty }=2.15$ , ${\it\theta}=30.8^{\circ }$ and $Re_{{\it\delta}^{\star }}=1050$ (Case A). Frequency response to internal broadband white noise forcing for nonlinear Navier–Stokes equations. (a) $P_{2}$ – separated shock, (b) $P_{4}$ – downstream.

To give insight into the influence of nonlinearities, we investigate the nonlinear response to a white noise forcing for case A. For that purpose, the same forcing as the one described in the previous section is used. Besides the linear regime, different simulations with various amplitudes are carried out. In particular, $A_{f}=10^{-8}$ and $A_{f}=10^{-4}$ are considered for $F1$ and $A_{f}=10^{-6}$ is considered for $F2$ . In figure 20, the WPSDs associated with $F1$ are reported for the probes $P_{1}$ and $P_{2}$ . It is interesting to note that the frequency selection observed in the linear regime is not strongly affected by nonlinearities. In particular, WPSDs associated with the probe $P_{1}$ exhibit an energy peak close to $St\approx 0.14$ and a wide range of amplified frequencies between $St\approx 0.05$ and $St\approx 0.5$ for both $A_{f}=10^{-8}$ and $A_{f}=10^{-4}$ , in agreement with the linear frequency response. However, the figure 20 also shows that low frequencies are energetically more significant when the forcing amplitude is increased for both $P_{1}$ and $P_{2}$ . This suggests a stronger selectivity of the lower frequency when nonlinearities have a significant contribution. These results are compatible with the results obtained by Sansica et al. (Reference Sansica, Sandham and Hu2014). However, in their work, the influence of nonlinearity is more pronounced because the interaction zone is larger which yields intense nonlinear dynamics of the shear layer. The WPSDs associated with $F_{2}$ are displayed in figure 21. In (a), the figure shows that the frequency response associated with pressure measurements in $P_{2}$ are not affected by nonlinearities. Similar results are obtained for $P_{1}$ . In (b), one may observe an increase of the energy content over the low frequencies when increasing the forcing amplitude inside the recirculation zone. On the one hand, this suggests that the nonlinear receptivity to frequencies associated with regions I and II (see spectrum 8) exhibit a similar behaviour as the one discussed in the section devoted to the linear regime. In particular, while the nonlinear receptivity of frequencies for region II are localized near the separation point, frequencies for region I have a greater receptivity to a forcing inside the recirculation zone. On the other hand, it is clear from the figure 21(b) that nonlinearities tend to increase this behaviour. Finally, in figure 22, the Fourier mode associated with $St=0.02$ ( $St\times L_{int}/{\it\delta}^{\star }\approx 0.73$ ) is represented. As for the linear regime, we may observe that its spatial support is mainly localized in the separated region. This indicates that nonlinearities tend to favour the emergence of low-frequency modes strongly associated with the recirculation zone. We can also remark that no acoustic radiation is observed in figure 22(a). This thus provides evidence that the mechanism associated with the nonlinear Fourier mode shown in figure 21 is distinct from the one linked to the supersonic Kelvin–Helmholtz modes observed in the linear regime.

Figure 22.

Case A: nonlinear regime. Fourier mode for $St=0.02$ . (a) The pressure component. (b) The spanwise vorticity component. The black solid line denotes a streamline above the separation line.