%20-%20https://sketch.com%20--%3e%3ctitle%3e3ECBBF50-FD08-412A-943D-2CFC5E025A16%3c/title%3e%3cdesc%3eCreated%20with%20sketchtool.%3c/desc%3e%3cg%20id='Stylesheet'%20stroke='none'%20stroke-width='1'%20fill='none'%20fill-rule='evenodd'%3e%3cg%20transform='translate(-926.000000,%20-668.000000)'%20id='BUTTONS'%3e%3cg%20transform='translate(100.000000,%20597.000000)'%3e%3cg%20id='BTN-/-RIGHTS'%20transform='translate(816.000000,%2065.000000)'%3e%3cg%20id='ICON-/-Rights'%20transform='translate(10.000000,%206.000000)'%3e%3cpath%20d='M13.5752841,0%20L1.34943182,0%20C0.59375,0%200,0.586666667%200,1.33333333%20L0,22.6666667%20C0,23.4133333%200.620738636,24%201.34943182,24%20L17.6505682,24%20C18.40625,24%2019,23.3866667%2019,22.6666667%20L19,5.36%20L13.5752841,0%20Z%20M13.5752841,1.89333333%20L17.0838068,5.36%20L13.5752841,5.36%20L13.5752841,1.89333333%20Z%20M1.34943182,22.6666667%20L1.34943182,1.33333333%20L12.1988636,1.33333333%20L12.1988636,5.36%20C12.1988636,6.10666667%2012.7926136,6.69333333%2013.5482955,6.69333333%20L17.6235795,6.69333333%20L17.6505682,22.6666667%20L1.34943182,22.6666667%20Z'%20id='Shape'%20fill='%23595959'%20fill-rule='nonzero'%3e%3c/path%3e%3ccircle%20id='Oval'%20stroke='%23595959'%20stroke-width='1.5'%20fill='%23595959'%20cx='9.5'%20cy='13.5355339'%20r='5.5'%3e%3c/circle%3e%3cpolyline%20id='Path-2'%20stroke='%23FFFFFF'%20stroke-width='1.5'%20transform='translate(9.500000,%2013.035534)%20rotate(45.000000)%20translate(-9.500000,%20-13.035534)%20'%20points='8%2016.0355339%2011%2016.0355339%2011%2010.0355339'%3e%3c/polyline%3e%3c/g%3e%3c/g%3e%3c/g%3e%3c/g%3e%3c/g%3e%3c/svg%3e)
1. Introduction
1.1. Overview
The response of shear layers to different spatio-temporal perturbations is a crucial step in better understanding their stability properties and acoustic responses. Various linear and nonlinear techniques have been developed to ascertain such responses. Most employ the steady laminar or the time-mean turbulent flow as the basic state, as obtained numerically from the Navier–Stokes equations. Experimental data do not provide all variables of interest and have thus been used only sparingly. The goal of the present work is to extract acoustic-related perturbation responses from time-resolved schlieren video data. For this, we combine a physics-based method with a data-driven reduced model, and demonstrate its effectiveness on a canonical high-speed shear layer for which spatio-temporally resolved experimental schlieren data are available.
Outline of the proposed perturbation method using schlieren data. Raw schlieren images,
$q$
, are first subjected to MPT irrotational filtering to extract the irrotational schlieren field,
$\Theta$
. Both raw and filtered schlieren data are then used in a DMD-ROM forcing-response analysis, under impulse and harmonic forcing conditions
$F$
.
Figure 1 provides a schematic overview of the two steps in the proposed method to process schlieren data to obtain the forced acoustic response of a turbulent mixing layer. In the first step, discussed in § 2.1, momentum potential theory (MPT), a physics-informed Kovasznay-like decomposition (Doak Reference Doak1998), is adapted to filter the density gradients in raw schlieren data (Prasad & Gaitonde Reference Prasad and Gaitonde2022), thus yielding the irrotational (in this case mainly acoustic) component. The second step, discussed in § 2.2, applies the dynamic mode decomposition reduced-order model (DMD-ROM) data-driven procedure of Stahl & Gaitonde (Reference Stahl and Gaitonde2024), to efficiently obtain the transient or asymptotic response, including energy amplification, of any defined forcing. The benefits of applying the DMD-ROM to the MPT filtered schlieren rather than the raw data clarifies the acoustic response for impulse and harmonic forcing. This novel combination of post-processing techniques isolates acoustic stability responses by perturbing schlieren video, without the need for actuated experiments, thereby offering a new and efficient method for studying the receptivity of noise emissions, providing valuable insights to inform potential actuated control strategies.
To demonstrate the new method, we use as our test case, a planar compressible shear layer formed by two parallel streams initially separated by a splitter plate for which spatio-temporally resolved schlieren data are available (Kim et al. Reference Kim, Elliott and Craig Dutton2020b
), as summarized in § 2.3. The Mach numbers of the upper (subscript
$1$
) and lower (
$2$
) streams are
$M_1=U_1/a_1=2.46$
and
$M_2/a_2=0.175$
, respectively, where
$M$
is the Mach number,
$U$
is the streamwise velocity and
$a$
is the speed of sound. The parameter that principally characterizes the dynamics of the shear layer is the convective Mach number,
$M_c = (U_1 - U_2)/(a_1 + a_2)$
. Among the different cases of Kim et al. (Reference Kim, Elliott and Craig Dutton2020b
), the
$M_c = 0.88$
case is considered because it displays an interesting combination of coherent structures, shocks and attendant compressibility effects, which influence the shear-layer instabilities. The objectives are thus to obtain transient and asymptotic stability results for disturbance sensitivity and acoustic amplification, which are not apparent in the turbulent schlieren video. Of particular interest are the most receptive forcing frequencies and locations in the upper and lower streams that generate shock related sound, which have not been examined for these benchmark experiments.
1.2. Context of the proposed approach
A review of stability approaches provides context for the proposed method. Stability analysis examines the growth of small perturbations on a basic state. Depending on the flow features, the initial linear perturbation growth may develop into intricate, nonlinear, multi-scale, chaotic phenomena. Stability analysis has evolved substantially for wall-bounded (Mack Reference Mack1984; Herbert Reference Herbert1997) as well as free shear-layer (Helmholtz Reference Helmholtz1868; Kelvin Reference Kelvin1871) flows. Foundational developments covering a range of sophisticated theoretical and computational methodologies have been extensively reviewed in the literature (Schmid Reference Schmid2007b ; Theofilis Reference Theofilis2011). Among the most prominent techniques are linear stability theory (LST) and the parabolized stability equations (PSE), early uses of which may be found in Mack (Reference Mack1984) and Herbert (Reference Herbert1997), respectively. These methods are subject to inherent assumptions and constraints, such as parallel flow (for LST), or slow downstream change (PSE). Nonlinear extensions are also becoming popular as for example nonlinear PSE. An alternative, complementary linear approach is resolvent analysis (Trefethen Reference Trefethen1999; Jovanovic & Bamieh Reference Jovanovic and Bamieh2005; McKeon & Sharma Reference McKeon and Sharma2010), which focuses on the identification of inputoutput forcing-response mechanisms and sensitivity measures to determine optimal growth pathways for instabilities (Taira et al. Reference Taira, Yeh and Munday2017). Assumptions and exclusion of nonlinearity in many of these approaches can result in poor prediction of experimental observations, particularly for high Reynolds and Mach number flows with complex geometry (Stetson Reference Stetson1988; Lachowicz et al. Reference Lachowicz, Chokani and Wilkinson1996; Lei & Zhong Reference Lei and Zhong2012). Extensions to two-dimensional(2-D) or 3-D analyses, necessary when the requirements of LST or PSE become untenable, substantially increase the computational burden (Juniper et al. Reference Juniper, Hanifi and Theofilis2014; Johnson et al. Reference Johnson, Candler and Alba2010; Paredes Reference Paredes2014) especially for more accurate methods such as bi-global or tri-global methods Theofilis (Reference Theofilis2011).
This work seeks to model stability responses directly from experimental data. The use of such data is relatively rare, since experimental measurements are generally incomplete, i.e. do not provide all variables necessary for stability analyses. For example, Particle image velocimetry (PIV) measures velocity (Willert & Gharib Reference Willert and Gharib1991; Adrian & Westerweel Reference Adrian and Westerweel2011), but not density or pressure, and schlieren images measure only density gradients in specific flow directions. Nonetheless, an interesting example of the use of PSE with mean velocity profiles from PIV and pressure probe data to understand jet noise may be found in Sinha et al. (Reference Sinha, Gudmundsson, Xia and Colonius2016).
The proposed method seeks to extract stability information, but directly and solely from time-resolved schlieren data, without augmentation with other measurements or Navier–Stokes simulations. This goal is advantageous as well as scope limiting as we now discuss. There have been continual improvements in obtaining higher-resolution time-resolved schlieren and shadowgraph video and currently such measurements are relatively common even at very high speeds. (Settles & Hargather Reference Settles and Hargather2017; Kennedy et al. Reference Kennedy, Jewell and Laurence2024; Paquin et al. Reference Paquin, Skinner and Laurence2023) However, their primary use has been to discern features in the flow to guide a better understanding of the dynamics; although pixel intensity scales may be referenced, post-processing schlieren data is often focused on qualitative detection of phenomena, such as shocks and shear- ayers, and coherent structures associated with dominant mechanisms. Our approach considerably enhances the value of these high-resolution data by facilitating their use in perturbation growth analysis.
However, schlieren data contain limited information, specifically direction-specific density gradients, which poses a challenge. To address this, after converting the schlieren data based on pixel intensity, we introduce a processing step to filter an irrotational component of acoustic nature (§ 2.1). When all variables are available, as from a scale-resolved simulation, MPT exactly decomposes ‘momentum-density’ fluctuations,
$(\rho \boldsymbol {u})'$
, where
$\rho$
is density, and
$\boldsymbol{u}$
is the velocity vector, into acoustic (irrotational-isentropic), thermal (irrotational-isobaric) and vortical (rotational) components, regardless of flow nonlinearity or fluctuation magnitude. This decomposition has been used to study instabilities in jets (Unnikrishnan & Gaitonde Reference Unnikrishnan and Gaitonde2019) and hypersonic transition (Unnikrishnan & Gaitonde Reference Unnikrishnan and Gaitonde2021). For schlieren data, an MPT extension by Prasad & Gaitonde (Reference Prasad and Gaitonde2022) extracts the irrotational component; in the unheated flows of interest here, where thermal fluctuations are not significant, the resulting field represents the acoustic dynamics. Filtering the schlieren data successfully extracts various aspects of shear layers not apparent from processing the raw schlieren data, including Mach wave radiation, feedback mode shapes, screech signatures and coupling between adjacent twin jets (Prasad et al. Reference Prasad, Gaitonde, Esfahani, Webb and Samimy2022).
Although MPT-based filtering effectively isolates schlieren content representative of acoustic phenomena, it does not by itself provide perturbation growth insights. For this, the second component of the proposed method in figure 1 uses a reduced-order model (ROM); when perturbations are imposed on the model, the transient and asymptotic linear forced response is obtained. Data-driven modal decomposition of the unsteady fluctuations is a natural choice upon which to build the desired model. The most common approach, proper orthogonal decomposition (POD), yields modes of correlated spatial structures ranked by energy (Berkooz et al. Reference Berkooz, Holmes and Lumley1993; Sirovich Reference Sirovich1987). Spectral POD (Lumley Reference Lumley1970; Towne et al. Reference Towne, Schmidt and Colonius2018; Schmidt Reference Schmidt2020) and conditional spacetime variants (Schmidt & Schmid Reference Schmidt and Schmid2019) (SPOD and CST-POD, respectively) can also generate insights on instabilities, including signatures of absolute (tonal) and convective instabilities and have been used for various problems including shear layers (Stack et al. Reference Stack, Adler and Gaitonde2019; Lu et al. Reference Lu, Huang, Wang and Hong2019) and related jets (Berkooz et al. Reference Berkooz, Holmes and Lumley1993; Boree Reference Boree2003; Karami & Soria Reference Karami and Soria2018; Towne et al. Reference Towne, Schmidt and Colonius2018; Weightman et al. Reference Weightman, Amili, Honnery, Soria and Edgington-Mitchell2018). Recent examples for processing schlieren data to extract coherent structures in such flows may be found in Weightman et al. (Reference Weightman, Amili, Honnery, Soria and Edgington-Mitchell2017), Liu et al. (Reference Liu, Li, Hao, Zhang and He2021) and Berry et al. (Reference Berry, Magstadt and Glauser2017).
For the present goal to discern perturbation growth, we seek a method that captures the time-accurate nonlinear behaviour encapsulated in the data through a linear approximation. Dynamic mode decomposition, or DMD, (Schmid Reference Schmid2007a ; Schmid et al. Reference Schmid, Li, Juniper and Pust2011) is an appropriate choice. Dynamic mode decomposition takes a dynamical systems approach to form a ‘best fit’, reduced linear operator that approximates the nonlinear dynamics based on Koopman theory. As a modal decomposition technique, time-dependent datasets are reduced by POD coordinates, enabling DMD eigen-modes analogous to those of the linearized Navier–Stokes equations (LNSE) operator. Connections between DMD modes obtained from statistically stationary data, SPOD and resolvent modes have been clearly established (Towne et al. Reference Towne, Schmidt and Colonius2018). Several stability related results have been obtained when the full scale-resolved simulation data are available. For example, CST-POD may be used to isolate causal forcing mechanisms and discern shear-layer instabilities (Stahl et al. Reference Stahl, Prasad, Goparaju and Gaitonde2023). Other approaches have leveraged the connection between DMD and stability modes in combination with LNSE-generated data. The mean flow perturbation (MFP) method obtains unstable tri-global modes that, unlike local LST, propagate linear disturbances on a Navier–Stokes-based 3-D basic state before subjecting them to DMD (Ranjan et al. Reference Ranjan, Unnikrishnan and Gaitonde2020).
In the present forced-response application, DMD is not used for modal analysis, but rather as a reduced-order model (DMD-ROM) where the linear operator is exploited to evolve the dynamics forward in a reduced subspace, as discussed in § 2.2. The stability analysis is formulated by applying a forcing term (see figure 1), representing the imposed perturbation, to the DMD-ROM, and analysing the temporal response of the disturbed flow. As demonstrated by Stahl & Gaitonde (Reference Stahl and Gaitonde2024) and Stahl et al. (Reference Stahl, Schwartz, Benton and Gaitonde2024) for scale-resolved numerical data, the approach constitutes a versatile and efficient framework, suitable for parametric studies of the instability dynamics, the receptivity of arbitrary forcing parameters and to guide flow control. Here, we considerably extend the approach by combining the DMD-ROM with MPT filtered data to induce and isolate acoustic-related instabilities directly from schlieren video alone.
1.3. Compressible shear-layer instabilities
The pertinent behaviour of a compressible shear layer is reviewed to place the acoustic response sought from the combined MPT and DMD-ROM method in context. Shear-layer instabilities are ubiquitous and associated with inflection points in velocity gradients across adjacent fluid layers (Browand Reference Browand1966; Michalke Reference Michalke1972). Small perturbations from the inflow environment or near the splitter-plate trailing edge excite instability waves that grow downstream. At low speeds, these waves roll up into vortical structures (Nitsche Reference Nitsche, Françoise, Naber and Tsun2006); subsequently, 3-D effects deform the structures through spanwise undulations, eventually resulting in a fully developed turbulent region characterized by self-similar time-averaged velocity and Reynolds stress profiles (Launder et al. Reference Launder, Reece and Rodi1975; Rogers & Moser Reference Rogers and Moser1994). The convective Mach number,
$M_c$
, is related to the speeds of the largest structures in the mixing-layer region, and also identifies compressibility effects related to the growth rates of instability waves (Papamoschou & Roshko Reference Papamoschou and Roshko1988a
). Compressibility effects typically become significant when
$M_c\gt 0.5$
and increase rapidly with the Mach number of the supersonic stream, as shocks of increasing strength augment the nonlinear development of large-scale structures; these effects asymptote beyond
$M_c\gt 1$
(Papamoschou & Roshko Reference Papamoschou and Roshko1988b
). Under such conditions, the generalized inflection point, which also factors density or temperature profile variations, becomes important (Leib Reference Leib1991). The chosen
$M_c = 0.88$
experiments of Kim et al. (Reference Kim, Elliott and Craig Dutton2020b
,Reference Kim, Elliott and Duttona)and Dutton et al. (Reference Dutton, Elliott and Kim2019) provide a rich compressible dataset to test the new analysis method proposed in this work.
Flow features of free shear layers for (a)
$M_c = 0.19$
(b)
$M_c = 0.88$
. Large Brown–Roshko roller structures at lower speeds are replaced by fine-scale turbulent mixing and the presence of shocks at higher convective Mach numbers. Schlieren images from Kim et al. (Reference Kim, Elliott and Craig Dutton2020b
) with labels added to highlight the features.
Schlieren images provide excellent visualization of compressible shear-layer instabilities. For illustration, instantaneous snapshots from the data of Kim et al. (Reference Kim, Elliott and Craig Dutton2020b
) are referenced in figure 2 at two convective Mach numbers, (a)
$M_c=0.19$
and (b)
$M_c=0.88$
. At lower
$M_c$
, the coherent structure roll-up is apparent in figure 2(a). Brown & Roshko (Reference Brown and Roshko1974) provide a comprehensive discussion of the shape, organization and dynamics of such large-scale Brown–Roshko roller structures found in incompressible planar mixing layers. The initially round structures transition to flattened and elongated waves with increasing
$M_c$
. Elliott et al. (Reference Elliott, Samimy and Arnette1995) similarly observed stretched shear-layer structures undergoing a pairing process at
$M_c=0.51$
, resembling Kelvin–Helmholtz-type structures. However, as the convective Mach number increases to
$M_c=0.85$
, no large-scale structures are observed as the mixing layer quickly transitions to turbulence. This behaviour is also evident in the
$M_c=0.88$
experiments of Kim et al. (Reference Kim, Elliott and Craig Dutton2020b
), shown in figure 2(b), where the mixing layer behind the shock rapidly transitions to fine-scale turbulence. The increased compressibility enhances the anisotropic behaviour of streamwise normal stresses, consistent with the observed stretching and breakdown of smaller structures. Although the turbulent mixing layer appears noisy in the schlieren image of figure 2(b), obscuring the source of acoustic waves emanating into the lower stream, the MPT filtering approach successfully extracts the details as demonstrated in this work.
Additional complexities are encountered at higher convective Mach number flows with shocks. Returning to figure 2(b), the shock system is initiated at the splitter-plate trailing edge and undergoes a series of expansions and reflections in the upper stream between the top wall and upper boundary of the mixing layer. The strength of the compression depends on the pressure difference across the splitter plate, which is equalized through the shocks, weakening downstream. Kim et al. (Reference Kim, Elliott and Craig Dutton2020b ) examined the effect of the shock, reporting that the turbulent statistics of the mixing layer were augmented by increased Reynolds stresses. An important byproduct of the shock interaction is emission of acoustic energy, which was not characterized by the schlieren or PIV results, but is relevant for shear-layer instabilities, for example, in the context of jet noise. Elsewhere, in simple and complex shear-layer configurations (Berry et al. Reference Berry, Magstadt and Glauser2017), unsteady shock interactions arise, with prominent acoustic waves from shedding at the splitter plate again suggesting the underlying importance of acoustic processes.
The background disturbance environment plays an important role in influencing shear-layer instabilities, and is therefore relevant to the present forced-response analysis. Tunnel noise is unavoidable, and ‘natural’ instability mechanisms from stochastic disturbances are well documented (Hussain Reference Hussain1986). Samimy et al. (Reference Samimy, Webb and Crawley2018) provide a recent overview of flow structures in high-speed shear layers, noting that naturally generated perturbations mostly belong to a broadband spectrum. The proximity of the excitation frequency to the most-amplified instability frequency can lead to a lock-on control situation, suppressing the natural pairing process of large-scale structures, motivating the impulse and harmonic forcing analyses shown in figure 1. For mixing layers, a numerical study indicated that 2-D and weakly oblique structures form downstream of the splitter plate (Sandham & Sandberg Reference Sandham, Sandberg, Armenio, Geurts and Fröhlich2010). However, these structures lack the high coherence observed in experiments (Papamoschou Reference Papamoschou1995), pointing to the existence of forcing noise components that are not accounted for in the simulations.
Other experiments investigate shear-layer stability by intentionally forcing the flow. For example, McLaughlin et al. (Reference McLaughlin, Martens and Kinzie1992) used a glow discharge forcing technique to excite the shear layer. At higher Mach numbers, additional oblique instability waves dominate, which are absent in the unforced experiments. Martens et al. (Reference Martens, Kinzie and McLaughlin1996) similarly conducted detailed measurements of high-speed Mach 3 and 4 streams with a low-speed Mach 1.2 stream, using a glow discharge to excite planar and oblique waves. Planar excitation remained planar, but oblique excitation changed the observations significantly, with the emergence of different wave propagation angles than initially used for excitation. Increasing the convective Mach number yielded slightly more oblique instability waves, consistent with the results of McLaughlin et al. (Reference McLaughlin, Martens and Kinzie1992). Forcing insights such as these have lead to promising control techniques for supersonic shear layers (Doshi et al. Reference Doshi, Gaitonde, Gist and Glauser2022; Oka et al. Reference Oka, Ozawa, Handa and Nonomura2024).
The aforementioned research motivates the desire to understand shear-layer receptivity and its response to different forcing inputs such as frequency, location and orientation. For experiments that intentionally force the shear layer (McLaughlin et al. Reference McLaughlin, Martens and Kinzie1992; Martens et al. Reference Martens, Kinzie and McLaughlin1996; Papamoschou & Roshko Reference Papamoschou and Roshko1988a ; Samimy et al. Reference Samimy, Webb and Crawley2018), these questions can be directly pursued through parametric studies. However, many experiments are limited in where the flow can be actuated or only examine the unforced system and thereby only elicit the natural instabilities excited by wind-tunnel noise, leaving uncertain the response to specific forcing. While PIV measurements may be used in a post facto analysis to test forced stability or noise generation mechanisms under these circumstances (Sinha et al. Reference Sinha, Alkandry, Kearney-Fischer, Samimy and Colonius2012; Ghasemi et al. Reference Ghasemi, Tuna and Li2018), no equivalent analysis using schlieren data has yet been conducted, which we address in this paper.
The paper is organized as follows. Each of the two steps in the procedure is described in § 2, including the filtering process in § 2.1 and the DMD-ROM forcing framework § 2.2, using figure 1 as a reference. The features of the raw schlieren data are outlined in § 2.3. An ancillary SPOD analysis in § 3 demonstrates the advantages of the filtered versus the raw schlieren data. In § 4, the DMD-ROM is implemented to examine the forced response of the system. An impulse forcing is first applied to isolate the shear-layer convective instability, followed by a harmonic forcing analysis to examine the asymptotic dynamics. Parametric receptivity studies test the frequencies and locations in the upper and lower streams that most amplify the acoustic response. Finally, conclusions are summarized in § 5.
2. Methodology
2.1. Irrotational filtering
The MPT procedure used to extract the irrotational component is first summarized in its original form (Doak Reference Doak1989), suitable for use when complete data are available as from scale-resolved simulations (e.g. Unnikrishnan & Gaitonde (Reference Unnikrishnan and Gaitonde2016)) before describing its adaptation to unsteady schlieren data (Prasad & Gaitonde Reference Prasad and Gaitonde2022). The MPT applies the Helmholtz decomposition to the ‘momentum-density’ vector field,
$\rho \boldsymbol{u}$
, to obtain mean and fluctuating rotational (solenoidal) components,
$\overline {\boldsymbol{B}}$
and
$\boldsymbol{B}'$
, and the irrotational fluctuating component,
$\nabla \psi '$
:
\begin{equation} \rho \boldsymbol{u}=\overline {\boldsymbol{B}} + \boldsymbol{B}^{\prime} - \nabla \psi ^{\prime}; \nabla \cdot \overline {\boldsymbol{B}}=0,\,\,\,\,\,\nabla \cdot \boldsymbol{B}^{\prime}=0,\,\,\,\,\,\nabla \times \psi '=0. \end{equation}
The scalar potential
$\psi$
comprises mean,
$\overline {\psi }$
, and fluctuating components,
$\psi ^{\prime}$
. Given a statistically stationary flow, the mean scalar potential gradient (
$\nabla \overline {\psi }$
) is zero. The scalar potential is calculated by substituting (2.1) into the conservation of mass equation
\begin{equation} \frac {\partial \rho }{\partial t} + \nabla \cdot ( \rho \boldsymbol{u})=0, \end{equation}
resulting in the Poisson equation
\begin{equation} \nabla ^{2}\psi ^{\prime}=\frac {\partial \rho ^{\prime}}{\partial t}. \end{equation}
The fluctuating scalar potential may be split into acoustic (isentropic)
$\psi ^{\prime}_{A}$
, and thermal (isobaric)
$\psi ^{\prime}_{T}$
components:
\begin{equation} \psi ^{\prime} = \psi ^{\prime}_{A} + \psi ^{\prime}_{T}, \end{equation}
which are given by
\begin{equation} \nabla ^{2}\psi ^{\prime}_{A}=\frac {1}{c^{2}}\frac {\partial P^{\prime}}{\partial t}\,,\,\,\,\,\,\,\,\,\,\,\,\,\, \nabla ^{2}\psi ^{\prime}_{T}=\frac {\partial \rho }{\partial S}\frac {\partial S^{\prime}}{\partial t}. \end{equation}
Here,
$c$
is the instantaneous local speed of sound,
$P$
is the thermodynamic pressure and
$S$
is the entropy.
Schlieren data provide density gradients; the form of (2.3), where the time derivative of density appears as the forcing term, can thus be leveraged to isolate the irrotational component. As discussed in Prasad & Gaitonde (Reference Prasad and Gaitonde2022), since the raw schlieren pixel intensity,
$q$
, is proportional to the streamwise density gradient,
\begin{equation} q \propto \frac {\partial \rho }{\partial x}, \end{equation}
the MPT may be used to solve for the appropriate derivative (streamwise in this case) of the irrotational scalar potential
$\psi '$
in the Poisson relation of (2.3):
\begin{equation} \nabla ^{2}\left(\frac {\partial \psi '}{\partial x}\right)=\frac {\partial }{\partial t}\left(\frac {\partial \rho }{\partial x}\right)^{\prime}. \end{equation}
A quantity proportional to the irrotational momentum in the streamwise direction,
$\Theta '$
, can be directly solved for by substituting the time derivative of pixel intensity fluctuations
\begin{equation} \nabla ^{2}\Theta ^{\prime}=\frac {\partial q^{\prime}}{\partial t}. \end{equation}
The Poisson equation is limited to a 2-D field of the irrotational vector component in the streamwise direction due to the inherently 2-D line-of-sight schlieren integration of the measured streamwise density gradients; thus, other vector components cannot be accounted for, but are generally less important compared with the streamwise direction. The irrotational filtered quantity,
$\Theta '$
, has contributions from both acoustic and thermal components. Previous numerical studies (Unnikrishnan & Gaitonde Reference Unnikrishnan and Gaitonde2018; Prasad & Morris Reference Prasad and Morris2020, Reference Prasad and Morris2021) have shown the thermal component exhibits coherent structures similar to its acoustic counterpart but is non-radiating, even for heated flows. For unheated flows, changes in the irrotational component are dominated by those in the acoustic component. Outside of the shear-layer region, radiated density fluctuations are purely acoustic and thus captured in
$\Theta '$
(Prasad & Gaitonde Reference Prasad and Gaitonde2022).
The solution of a Poisson equation requires suitable boundary conditions. Several prior efforts (Unnikrishnan & Gaitonde Reference Unnikrishnan and Gaitonde2018; Prasad & Morris Reference Prasad and Morris2020; Stahl & Gaitonde Reference Stahl and Gaitonde2021; Stahl & Gaitonde Reference Stahl and Gaitonde2021) have explored the sensitivity of the results to boundary condition specifics. When the boundaries are far from regions of rapid gradients, suitable conditions may be obtained by integrating the momentum vector to obtain the irrotational boundary conditions. The solution is relatively insensitive to boundary conditions when applied far from the regions of interest, since the potential function is dominated by source term fluctuations. For schlieren data on the other hand, the
$\rho \boldsymbol {u}$
is not available and the boundaries may not lie in regions of low vorticity. Following Prasad & Gaitonde (Reference Prasad and Gaitonde2022), an effective procedure is to use a Dirichlet boundary condition of
$\Theta '=0$
with a hyperbolic tangent sponge zone extending eight pixels from the border. The thin splitter plate is not treated as a boundary condition due to the limited pixel resolution and sponge zone which would damp the boundary-layer fluctuations. While this decreases spatial resolution across the splitter-plate region, the resulting wavepacket signatures are nonetheless dominated by the Poisson equation source term (i.e. pixel fluctuations) and represent the dynamics on the upper and lower sides, with no residual effect further downstream. An alternative approach to process schlieren data using MPT in the frequency domain has been discussed by Padilla-Montero et al. (Reference Padilla-Montero, Rodríguez, Jaunet and Jordan2024).
2.2. The DMD-ROM forcing and response
The DMD-ROM forcing methodology. (a) Input flow snapshots provide the DMD linear operator
$\boldsymbol {\tilde {A}}$
and POD modes
$\textbf{U}$
for snapshot reduction. (b) User defined forcing snapshots. (c) The DMD-ROM time marching scheme featuring the reduction of the initial condition, time advancement with forcing and the final expansion of the forced solution.
The shear-layer stability and receptivity study is enabled by applying a forcing term to a DMD-ROM, as outlined in the framework of Stahl & Gaitonde (Reference Stahl and Gaitonde2024), where it was employed with the unsteady pressure field obtained from a scale-resolving simulation. An overview of the method is shown in figure 3. The DMD-ROM is derived from the (a) snapshot data matrix of either the raw or MPT-filtered schlieren data,
$\unicode{x1D64C}=[q_1,\ldots ,q_N]$
, where there are
$N$
individual snapshots,
$q_n$
, of size
$M$
pixels. The DMD ultimately seeks a linear operator,
$\unicode{x1D63C}$
, that marches the system from one snapshot,
$q_n$
, to the next,
$q_{n+1}$
:
\begin{equation} q_{n+1}=\unicode{x1D63C} q_{n}. \end{equation}
However,
$\unicode{x1D64C}$
and
$\unicode{x1D63C}$
are the full snapshot size (
$M \times N$
and
$M \times M$
, respectively). For computational economy, the snapshots and dynamic operator are projected onto POD modes, obtained by first splitting and shifting the snapshot matrix into
$ \unicode{x1D64C}_1= [q_1,\ldots ,q_{N-1}]$
and
$\unicode{x1D64C}_2=[q_2,\ldots ,q_{N}]$
, and then performing the singular value decomposition
\begin{equation} \unicode{x1D64C}_1=\unicode{x1D650}\boldsymbol {\Sigma } \unicode{x1D651}^{*} \end{equation}
to produce POD modes
$\textbf{U}$
, which are used to calculate the reduced DMD operator
$\tilde {\kern-2pt\unicode{x1D63C}}$
:
\begin{equation} \tilde {\kern-2pt\unicode{x1D63C}}=\unicode{x1D650}^{*}\unicode{x1D63C}\unicode{x1D650}=\unicode{x1D650}^{*}\unicode{x1D64C}_2\unicode{x1D651}\boldsymbol {\Sigma ^{-1}}. \end{equation}
The
$\tilde {\kern-2pt\unicode{x1D63C}}$
matrix is then of size
$N\times N$
. A POD mode truncation may be used to further reduce or stabilize
$\tilde {\kern-2pt\unicode{x1D63C}}$
as becomes necessary for the filtered schlieren data (Appendix A).
Reduced snapshots,
$\tilde {\kern-2pt\unicode{x1D64C}}=[\tilde {q}_1,\ldots ,\tilde {q}_N]$
, are projected onto POD modes
\begin{equation} \tilde {\kern-2pt\unicode{x1D64C}}=\unicode{x1D650}^{-1} \unicode{x1D64C}, \end{equation}
and consistently evolved forward in time, reproducing the input solution. The stability and receptivity analysis is enabled by the addition of a forcing term
\begin{equation} \tilde {q}^{\prime}_{n+1}=\tilde {\kern-2pt\unicode{x1D63C}} \tilde {q}_{n} + \tilde {f}_{n}. \end{equation}
As shown in figure 3(b), user-defined forcing snapshots,
$\unicode{x1D641}=[f_1,\ldots ,f_N]$
, provide the forcing term and are of consistent size with the schlieren snapshots because they are projected onto the reduced subspace using the POD modes of the flow data
\begin{equation} \tilde {\kern-2pt\unicode{x1D641}}=[\tilde {f}_1,\ldots ,\tilde {f}_N]=\unicode{x1D650}^{-1} \unicode{x1D641}. \end{equation}
After the DMD operator, schlieren, and forcing snapshots have been reduced, the system is iterated forward in time under the influence of forcing, according to (2.13) and shown in figure 3(c). The response snapshots,
$\tilde {\kern-2pt\unicode{x1D64C}}'=[\tilde {q}^{\prime}_1,\ldots ,\tilde {q}^{\prime}_N]$
, are collected and expanded back to the full sized domain for observation by
\begin{equation} \unicode{x1D64C}'=\unicode{x1D650}\tilde {\unicode{x1D64C}}'. \end{equation}
Apart from directly observing the forced-response snapshots, the energy amplification is measured using gain metrics that compare the L2 energy norm of the response with the forcing. For parametric studies, the gain is more efficiently calculated in the reduced space to avoid the expansion step (2.15). The total gain over all snapshots is defined as
\begin{equation} \sigma =\frac {||\,\tilde {\kern-2pt\unicode{x1D64C}}'\,||_2}{||\,\tilde {\kern-2pt\unicode{x1D641}}\,||_2}, \end{equation}
and is used for the harmonic forcing analysis. A transient gain is also defined for individual snapshots
\begin{equation} \sigma (t)=\frac {||{\tilde {q}(t)}'||_2}{||\,\tilde {\kern-2pt\unicode{x1D641}}\,||_2}. \end{equation}
For impulse forcing, this gain is defined with respect to the initial forcing snapshot
$f_o$
. The behaviour of decaying amplitude perturbations may be exaggerated by using the transient gain to normalize the impulse response
$q(t)'/\sigma (t)$
, as discussed in § 4.
In previous studies using the DMD-ROM (Stahl & Gaitonde Reference Stahl and Gaitonde2024), the forcing was applied to the unsteady flow, emulating a simpler analogue of the synchronized large eddy simulations approach of Adler & Gaitonde (Reference Adler and Gaitonde2018). From an alternative perspective of a perturbed basic or equilibrium state, the forcing would be applied to the mean flow to better isolate the transient response, similar to MFP analysis (Ranjan et al. Reference Ranjan, Unnikrishnan and Gaitonde2020). This latter perspective is reproduced here with the DMD-ROM by replacing the initial condition snapshot with the mean flow. However, because DMD is able to distinguish the mean flow as a mode without a growth rate or frequency, the mean flow implicitly does not evolve with time. As such, the dynamic operator only evolves the perturbations in time without contributions of the mean flow, which would be subtracted anyway in the case of MFP. Due to its redundancy, here the mean-flow initial condition is replaced by a zero-vector snapshot base flow.
The forcing snapshots are populated with impulse or harmonic forcing at locations just upstream of the splitter-plate trailing edge. The impulse is applied only at the first time step, while the harmonic forcing lasts the entire duration of the schlieren video,
$f(t)=A\cos (\omega t)$
. The impulse forcing thus targets the growth over a shorter time horizon, while the harmonic forcing investigates the asymptotic, statistically stationary response. Note that, in contrast to typical stability methods that necessarily consider small perturbations, the DMD-ROM forcing amplitude tests conducted by Stahl & Gaitonde (Reference Stahl and Gaitonde2024) revealed that larger amplitudes were more effective in eliciting the relative linear gain rates across forcing parameters (i.e. receptivity). This is primarily dependent on the size of the perturbations relative to the fluctuations of the unsteady flow. Here, the perturbation response is calculated without the comparatively larger unsteady fluctuations or mean flow. Under these conditions, the output response is simply scaled by the input amplitude due to the linearity of the DMD-ROM; this scaling is accounted for by the normalization of gain in (2.16). Therefore, the results do not change as a function of forcing amplitude, which is arbitrarily set to unity for all cases.
(a) Wind-tunnel configuration and test section. (b) Splitter plate, incoming boundary layers and flow conditions. Not to scale.
Flow conditions of upper and lower streams.
2.3. Experiments and schlieren data
Turbulent shear-layer measurements are used from Kim et al. (Reference Kim, Elliott and Craig Dutton2020b
), who obtained results at several test conditions reflecting different levels of compressibility effects; these are documented at https://wiki.illinois.edu/wiki/display/NCSLF. Case 5, corresponding to
$M_c=0.88$
(
$M_1=2.461$
and
$M_2=0.175$
) is chosen for reasons discussed in § 1.3. A schematic of the test apparatus is introduced in figure 4(a) with the splitter plate, incoming boundary layers and shear-layer configuration isolated in (b). For reference, the flow conditions are reproduced in table 1 with boundary-layer characteristics listed in table 2. The splitter-plate thickness relative to the boundary-layer height is an influential parameter for the mixing-layer dynamics, since it controls the relative importance of shedding versus Kelvin–Helmholtz instabilities; a recent effort by Stack et al. (Reference Stack, Adler and Gaitonde2019) systematically highlights the principal effects as the plate thickness is changed. Here, the splitter plate tapers to a thin trailing edge with a thickness of
$T=0.5\,\textrm{mm}$
, and marks the coordinate system origin for the streamwise (
$X$
), wall-normal (
$Y$
), and spanwise (
$Z$
) directions. The ratio of the splitter-plate thickness to boundary-layer height is approximately
$T/\delta _1=1/6$
and
$T/\delta _2=1/8$
for the upper and lower streams respectively; this suggests immediate mixing of small-scale structures, as evidenced in figure 2. The spanwise extent of the configuration is
$127\,\textrm{mm}$
, which, compared with the thin sidewall boundary layers, establishes spanwise homogeneous flow in the centre of the duct where schlieren line-of-sight integration is focused.
Boundary-layer conditions on splitter plate.
(a) Time-mean velocity field from PIV data. (b) Time-mean and (c) instantaneous schlieren images.
Optical measurements were obtained across the entire
$127 \times 127 \times 762\,\textrm{mm}$
test section, as depicted in figure 4(a). Due to spatio-temporal resolution considerations, the test section is partitioned into multiple smaller subsections. For reference, PIV measurements of the time-mean velocity field in the subdomain immediately downstream of the splitter plate are presented in figure 5(a), showing the mixing-layer region and velocity profiles at multiple streamwise locations. Details of the PIV set-up with in-depth analysis of the velocity field and Reynolds stresses may be found in Kim et al. (Reference Kim, Elliott and Craig Dutton2020b
). More pertinent here are schlieren measurements, which use the angular deflection of collated light passing through the spanwise direction of the flow to record pixel intensities modulated by streamwise density gradients.
Time-mean and instantaneous schlieren images are shown in figures 5(b) and 5(c), respectively, better illustrating the shock structures in the upper stream. In this work, the subdomain closest to the splitter plate is used, which has slightly different dimensions for PIV and schlieren windows. The schlieren field of view extends
$126.5 \times 39.5$
mm around the splitter plate and is resolved with
$384 \times 128$
pixels. The smaller region affords a higher sampling rate of
$F_s=124\,\textrm{kHz}$
, which is sufficient to calculate the time derivatives for the irrotational filtering and DMD-ROM analysis. A total of
$2{,}500$
images, spanning
$20.2\,\textrm{ms}$
, are used to form the models. Since our interest is in the transient response of instabilities, time is reported in non-dimensional form throughout the paper by
$t^*=tF_s$
.
3. Spectral analysis
Prior to investigating the transient forcing and response, a conventional spectral analysis is conducted to demonstrate the advantages of the filtered schlieren data over the raw values. Spectral POD is used with the block Welch formulation of Schmidt (Reference Schmidt2020). The
$2{,}500$
images are partitioned into three blocks with 50 % overlap, resolving a frequency range of
$150$
to
$62{,}000$
Hz. The resulting spectra of the first three modes are shown in figure 6 for the (a) raw and (b) filtered data. For the purposes of this work, the relative amplitude of frequencies within each spectrum is more important than the amplitude comparison across (a) and (b); the quantities are of different magnitudes, but have been normalized for convenience. The (a) raw schlieren is dominated by ultra low frequencies; this bias is likely not physical and is attributed to camera effects. In the dynamically significant range beyond
$10^3$
Hz, the constant broadband spectra yield little insight into which modes are important. In contrast, the (b) filtered schlieren shows the irrotational energy is biased towards a higher frequency band, and is clearly dominated by several tones. The largest tones (
$12.23$
and
$15.25\,\textrm{kHz})$
are selected for further analysis.
Spectral POD spectra for the (a) raw and (b) filtered schlieren data.
Spectral POD modes for the (a) raw and (b) filtered schlieren data at the peak acoustic tones (12.23 and 15.25 kHz). The dashed lines represent the approximate locations of the shock and the shear layer extended from the splitter plate.
The corresponding SPOD modes of the identified tones are presented in figure 7, again comparing the (a) raw and (b) filtered schlieren. The black dotted lines reference the shock and shear-layer line extending between supersonic and subsonic streams, as obtained from figure 5(b), of which the supersonic stream contains stronger fluctuations. Here, the differences between the data sets are clear; the (a) raw schlieren modes display low coherence associated with the mixing of smaller turbulent structures, primarily in the upper stream, while the (b) filtered schlieren modes feature larger and well-defined wavepacket structures. For the most energetic
$12.23\,\textrm{kHz}$
mode, the largest irrotational structures occur in the supersonic stream. The higher frequency
$15.25\,\textrm{kHz}$
mode displays smaller waves, with a distinct change in wavepacket structure before and after the shock. Also of significance are the acoustic waves propagating in the lower stream that are not evident in the raw modes. The filtered schlieren SPOD modes are very revealing in terms of isolating the general structure of upper-stream wavepackets (
$W$
) and lower-stream acoustics (
$A$
), and of note for later, the acoustic wavelength
$\lambda$
. The wavepackets comprise acoustic and thermal/entropic components from the shear layer. The former propagate into the free stream, while the latter are prominent in the shear region and decay as they convect However, the convective phenomenon that relates these different wavepacket structures across all frequencies is not entirely clear in isolated SPOD modes, motivating the transient impulse analysis to further distinguish the evolution of wavepackets and acoustic radiation.
4. Forced-response analysis
The forced response is first obtained for an impulse applied simultaneously across both boundary layers, which yields insight into the the overall stability properties of the shear layer. This is followed by two more precise receptivity studies, where specific locations are independently perturbed with impulse or harmonic forcing. For the DMD-ROM method, both forcing types (see § 2.2) are implemented by constructing empty snapshots and populating individual pixels in the wall-normal direction of the splitter plate with the prescribed forcing function over the duration of the video.
4.1. Impulse forcing across both boundary layers
Impulse response comparison of the raw (left) and filtered (right) schlieren, progressing in time from top to bottom. The impulse location in shown in the first frame. Wavepacket (
$W$
) and acoustic (
$A$
) structures are identified convecting downstream.
Figure 8 displays the impulse response for the situation where the disturbance is applied across a single vertical column of pixels at
$X=-12\,\textrm{mm}$
(outlined at time
$t^*=2$
) such that it perturbs both boundary layers simultaneously. Results from the raw and irrotational filtered schlieren data are presented in the left and right columns, respectively, and scaled by the transient gain, as discussed in § 2.2, to better visualize the structures in each frame. The signature of the initial disturbance is observed around this location and induces a relatively quiescent background of small-scale noise; an artefact of the DMD model which is quickly damped out. The response progresses in time from the top frame to the bottom frame, with the disturbance convecting downstream and growing into a large wavepacket structure (w) that is noticeably more coherent in the filtered data. By time
$t^*=20$
, it is clear that the wavepacket in the upper stream disperses when crossing the shock (w1 and w2) and two different speeds become apparent. In contrast, the initial impulse disturbance in the lower stream is much weaker compared with the larger waves in the upper stream. However, an acoustic wave (
$A$
) generated later at the trailing edge of the splitter plate (
$t^*=10$
) propagates through this subsonic region, and is followed by a series of weaker waves emitted from the shock interaction. Such acoustic waves that propagate through the subsonic region are evidently related to the corresponding supersonic convection of hydrodynamic structures in the upper stream (Hall et al. Reference Hall, Dimotakis and Rosemann1993). Nonetheless, as the raw schlieren response progresses in time, these disturbed structures are barely discernible as they become indistinguishable in the turbulent mixing layer.
(a) Data extraction location (dashed line) in the upper stream. Corresponding space–time plots from the (b) raw and (c) filtered schlieren impulse responses. The initial wavepacket
$w$
, bifurcates across the shock intersection location (arrow in (a)), maintaining speed and structure with
$w1$
and producing a slower wavepacket
$w2$
.
The impulse wavepackets are further investigated with space–time plots of the response in figure 9. The data are extracted from pixels located on the dashed probe line in figure 9(a) to capture the wavepacketshock interaction region. Space–time plots are compared for the figure 9(b) raw and (c) filtered schlieren results, again demonstrating a clearer representation of the wave propagation pattern in the irrotational data. As time progresses, downstream propagating waves are interpreted as streaks travelling down and to the right, with a speed corresponding to the slope of the line. Examining figure 9(c), the incipient wavepacket (
$w$
) has a rightward running path that maintains a constant slope after the shock intersection location (
$w1$
), corresponding to an approximate speed of
$535\, \textrm{m/s}$
. The second wavepacket (
$w2$
) is generated after the shock and has a steeper slope, which corresponds to a slower
$306\, \textrm{m/s}$
.
Schematic of the shear-layer instability and dispersion of wavepackets across the shock. The initial wavepacket
$w1$
emits acoustic wave
$A1$
at the splitter-plate trailing edge. The wavepacket–shock interaction produces the slow acoustic wavepacket
$w2$
and acoustic emission
$A2$
.
Figure 10 presents a generalized illustration of the acoustic wavepacket dynamics produced by the shear-layer instability. The horizontal dashed line represents the location where space–time plots were extracted to compute the wavepacket velocities shown in the diagram. Starting at
$t^*=5$
, the incoming disturbance along the upper boundary layer reaches the trailing edge of the splitter plate, emitting the acoustic wave
$A_1$
in the lower stream. Further downstream, at approximately
$t^*=15$
, the disturbance interacts with the shock, bifurcating the upstream wavepacket
$w1$
and
$w2$
.
It is well established that such a shock interaction results in a complicated dispersion of vortical, acoustic and entropic waves, with varying phase speeds and orientations based on the shock angle and strength (Duck et al. Reference Duck, Lasseigne and Hussaini1995). The filtered schlieren data isolate the acoustic components of this process. Given the measured velocities and nominal flow conditions in the upper stream (
$M_1=2.46, U_1=564\, \textrm{m/s}$
and
$a_1=229\, \textrm{m/s}$
, before the shock), wavepacket
$w1$
is limited by the supersonic flow speed
$u_1=535\, \textrm{m/s}$
, while
$w2$
matches the slow acoustic wave speed
$u_1-a_1=306\, \textrm{m/s}$
. This interaction also generates an acoustic wavefront
$A2$
, distinguished from
$w1$
for its propagation out of the supersonic stream and into the subsonic flow. As shown at
$t^*=25$
, dispersion effects from the variable Mach number across the shear layer cause refracted portions of the acoustic wave (
$A_1'$
and
$A_2'$
) to travel at slower speeds when passing through supersonic flow, although, this is more difficult to measure with limited pixel resolution in the mixing region. Finally, the time delay between the production of the trailing edge acoustic wave
$A1$
and the slower shock-acoustic wave
$A2$
, defines the distance
$\lambda$
, observed earlier, for instance in the SPOD mode wavelength of figure 7(b).
4.2. Impulse location study
(a) Parametric impulse forcing locations and acoustic probe location where the response is measured. All analyses are conducted on the irrotational filtered variable
$\Theta$
. The most (
$Y=1.8$
) and least (
$Y=-17.5$
) amplified responses are compared. (b) Transient gain,
$\sigma$
, measures the energy growth and decay across the entire domain. (c) Acoustic probe fluctuations
$\Theta$
are used to measure the time lag
$\Delta T$
between successive acoustic waves. (d) The acoustic signal is normalized by the transient gain,
$\Theta /\sigma$
, to better compare the asymptotic energy at the probe location, as measured by the variance of the fluctuations.
The time delay between waves in the acoustic field may be explicitly measured by using isolated disturbances to further analyse the response, including the value of
$\lambda$
. Individual pixels are perturbed along the same vertical column to find the most sensitive locations in each stream. The acoustic response is measured by a probe in the lower stream at
$X=45\,\textrm{mm}$
, identified by the yellow circle in figure 11(a). The instantaneous raw schlieren image in figure 11(a) is only shown for reference, while the analysis is conducted on the irrotational flow,
$\Theta '$
, to measure the acoustic output produced by each disturbance location. To illustrate the disparate effects of forcing location, the most and least amplified responses from the upper (red,
$Y=1.8$
) and lower (blue,
$Y=-17.5$
) streams are compared.
The impulse response is measured in three ways in figure 11: the (b) transient gain encompassing the energy of the entire flow field, the (c) acoustic probe signal and the (d) probe signal normalized by the transient gain. Examining the transient gain of (b), the initial energy amplification of the convective instability is followed by asymptotic decay. Both forcing locations experience their peak growth at
$t^*=3$
but have different decay rates. The probe signal in (c) illustrates the initial response of acoustic fluctuations, representing the passage of waves, which similarly decay over time. For the
$Y=1.8$
input forcing, the time between the two largest peaks is
$\Delta T_A=15=0.121\, \textrm{ms}$
and corresponds to the lag between acoustic waves
$A1$
and
$A2$
. Given the speed of sound measured in the experiments, this yields an instantaneous acoustic wave distance
$\lambda =a_2 \Delta T_A = 41\,\textrm{mm}$
. This is comparable to the characteristic wavelength observed in the
$15.25\,\textrm{kHz}$
SPOD mode of figure 7(b) in this probe region. The acoustic signal in (d), normalized by the transient gain, provides a clearer comparison of fluctuation size over time, resulting in a periodic signal as the transients decay. The asymptotic frequency is
$31\,\textrm{kHz}$
, which is a harmonic of the peak SPOD tone; this frequency represents an artefact of the slowest decaying mode in the DMD model. Regardless, the transient normalization is useful because, for the same forcing amplitude, the relative energy output is preserved over time, as measured by the variance of the normalized signal (d), which is utilized in the upcoming parametric study.
(a) Variance of the normalized acoustic signals (measured at the probe in figure 11 a) as a function of wall-normal impulse location. The (b) schlieren image and streamwise velocity profile provide reference for these locations. The variance peaks are identified and normalized relative to the average of the lower stream, showing significant multiples of acoustic output in the upper stream relative to free-stream disturbances.
The amplitude response of all individual forcing locations is compared using the probe signal variance in figure 12(a). This metric effectively determines the relative acoustic energy output of disturbances placed along the upper (red) and lower (blue) streams. These input locations are mapped to the flow field in figure 12(b), depicting a schlieren image overlaid with time-mean streamwise velocity profiles, taken from PIV data immediately after the splitter plate. The overall trend in variance shows several peaks in the upper stream that are much larger than the lower stream, which is nearly constant apart from the boundary-layer region. The average lower-stream variance (vertical dashed line in (a)) of acoustic disturbances in the subsonic free stream is not produced from the splitter plate and is more representative of ambient wind-tunnel noise.
This background noise is used to normalize the relative energy amplification of other peak locations, marked with dashed lines in figure 12. The peaks correspond to the lower-stream boundary layer (
$Y=-1.8\,\textrm{mm}$
), and upper-stream boundary layer (
$Y=0.6,1.8$
and
$4.0$
mm). The obvious conclusion is that the maximum acoustic energy is generated from perturbations in the upper-stream, outer boundary-layer region and is approximately nine times larger than subsonic free-stream disturbances. In comparison, the maximum energy in the lower stream is from forcing the turbulent boundary layer (Y = −1.8 mm) and is only three times greater than disturbances in the free stream, where the acoustic density fluctuations are even weaker under the splitter plate and have little impact on the shear-layer instability growth. These results illustrate how the combined post-processing techniques can yield relative amplification data from qualitative schlieren data, and are consistent with receptivity studies showing higher-speed flow is more receptive at the boundary-layer edge to free-stream disturbances which can sustain growth downstream (Schmid & Henningson Reference Schmid and Henningson2001).
4.3. Harmonic forcing sweep
Harmonic forcing gain
$\sigma$
of the (a) raw and (b) acoustic filtered schlieren data as a function of forcing pixel location in the splitter-plate boundary layer.
The receptivity and stability properties of the shear layer are now investigated with parametric harmonic forcing over a range of forcing frequencies (
$0$
to
$60$
kHz) at the same wall-normal pixel locations as the previous analysis. The DMD operator, forcing amplitude and initial conditions are the same as used for the impulse study. The total gain (2.16) across all response snapshots is first presented in figure 13, comparing (a) raw and (b) irrotational schlieren data. Both gain plots affirm that the supersonic stream is more receptive to amplifying disturbances than the lower stream, for nearly all frequencies. The raw schlieren, (a), is generally noisy with flat broadband spectra. As with the SPOD results, this observation does not provide much insight into the most important frequencies. In contrast, the filtered schlieren analysis, (b), shows several prominent tones in the supersonic stream in the 20–40
$\textrm{kHz}$
range.
The (a) general impulse response predicts the (b) most receptive harmonic forcing conditions, labelled A through D, bounded by the sub- and super-harmonics of the general impulse response frequency
$f_I=21\,\textrm{kHz}$
. Individual impulse locations provide more precises responses, for example with the most amplified harmonic forcing
$f_B=17.5 Hz$
at
$Y=3.3\textrm{mm}$
.
The most receptive conditions from the filtered results are further analysed in figure 14 with connection to the earlier (a) general impulse response and (b) prominent forcing frequencies, identified and labelled A–D. To make this connection, consider that for a linear time-invariant system, an impulse effectively triggers all frequencies in the system, with the response indicative of those that display highest receptivity. figure 14(a) captures that general impulse response with the acoustic signal extracted just downstream of the wall-normal forcing locations AD. A characteristic impulse time scale,
$\Delta T_I$
, is measured to estimate such a naturally preferred frequency, which corresponds to
$f_I=1/\Delta T_I=21\,\textrm{kHz}$
.
By mapping the general impulse response frequency onto the gain results of figure 14(b), the most amplified forcing frequencies are identified and bound the overall sensitive frequency range with its first sub- and super-harmonic,
$1/2 f_I$
and
$2f_I$
, respectively; the reasons for these bounds are discussed later. The most receptive conditions (B) occur at frequency of
$\sim 17.5\,\textrm{kHz}$
at
$Y=3.3\textrm{mm}$
near the outer edge of the upper-stream boundary layer, where the velocity gradients change rapidly. Repeating the impulse at only this location (cutout in (b)) yields a longer response time scale
$\Delta T_B$
which matches this peak frequency. Similar results follow for (A) and (C) with nearby frequencies located in the upper stream, but are not shown for brevity. Although the lower-stream gain is much weaker, a distinct frequency of
$16\,\textrm{kHz}$
is observed for (D), and is also selected for further analysis.
Instantaneous snapshots of the harmonic forcing response showing the raw (left) and filtered (right) schlieren results for the forcing conditions at (a)
$Y=9.9\,\textrm{mm}, 22.5\,\textrm{kHz}$
, (b)
$Y=3.3\,\textrm{mm}, 17.5\,\textrm{kHz}$
, (c)
$Y=2.5\,\textrm{mm}, 29\,\textrm{kHz}$
and (d)
$Y=-1.8\,\textrm{mm}, 16\,\textrm{kHz}$
.
The flow-field response at the identified harmonic forcing conditions is investigated by capturing instantaneous snapshots during the asymptotic stage of the forcing. The instances selected in figure 15 represent the salient dynamics throughout the periodic cycle corresponding to the forcing frequency. Figure 15 shows both the raw schlieren (left column) and filtered schlieren (right column) under forcing parameters (A) through (D). The raw schlieren generally shows small-scale turbulent structures in the shear layer, while the filtered images reveal larger wavepacket structures and acoustic radiation that emphasize longer wavelengths. Also noticeable in the filtered results are the incipient disturbances near the forcing locations, marked by arrows. In the case of (A) and (C), these are accompanied by larger wavepackets above the splitter plate, while (B) and (D) display structures primarily downstream of the shock. The acoustic radiation wavelength for case (A) is consistent in both the upper and lower stream. Case (B), which is the most energetic response, has a more complicated structure, particularly around the shock region. Case (C) demonstrates that this forcing is receptive, but not prone, to generating acoustic radiation, which is lacking in the lower stream. This suggests potential forcing conditions that the flow could ‘lock-on’ to without amplifying emitted noise. Finally, case (D) shows the forcing in the lower stream has a unique angle of acoustic propagation compared with the other cases. This case also demonstrates that, under harmonic forcing conditions, the DMD-ROM response is closely related to spectral modes, as compared with similar SPOD modes of figure 7(b) at the nearby
$15.25\,\textrm{kHz}$
frequency. In this case, however, the DMD-ROM can target the precise locations that are sensitive to these frequencies, suggesting the integral role of the lower stream in noise generation.
4.4. Receptivity discussion
Finally, we discuss the unification of the impulse and harmonic results to better characterize the mechanisms contributing to the most receptive locations. The instabilities of a nominally unforced compressible shear layer are well known in the classical literature. Receptive frequencies correspond to natural modes associated with the Kelvin--Helmholtz (KH) instability (Ray Reference Ray1982). At various upstream locations in the compressible boundary layer, the interaction between acoustic waves and vortical structures plays a critical role in influencing the characteristics of KH instabilities, including the frequency, wavelength, phase speed and pairing of coherent structure observed in the shear layer (Martens et al. Reference Martens, Kinzie and McLaughlin1994). For instance, at supersonic conditions, as disturbances propagate through subsonic regions of a turbulent boundary layer, upstream and downstream acoustic waves may interact, leading to increased wavepacket growth over the splitter plate, modulated here by the nearby expansion-shock system reflected off the trailing edge. This incipient point can augment nonlinear phenomenon (such as the convective shock-acoustic dispersion observed in the transient impulse response of figure 10), distorting competing instability modes and their harmonics, resulting in the complex frequency spectrum found in the turbulent mixing layer and a strong radiating source of acoustic emissions.
The role of the excitation (forcing) frequency in supersonic shear layers has been discussed by Samimy et al. (Reference Samimy, Webb and Crawley2018), who used arc plasma actuators of relatively high bandwidth. They note the importance of the
$1/2$
excitation frequency sub-harmonic, which can interact with the natural shear-layer modes, or their harmonics, leading to a special ‘lock-on case’, where input disturbances grow to saturate the response; these being the most receptive forcing conditions. In particular, they compared different excitation frequencies with an impulse forcing, finding the harmonic response to be approximated by a linear superposition of impulse responses, as measured in the irrotational acoustic field, with Mach wave radiation similar to that observed in our study. This is consistent with the analysis of figure 14, where the general impulse response is related to the lower bound of the receptive harmonic forcing frequency range at the
$1/2f_I$
sub-harmonic, representing its naturally preferred instability. The role of the excitation frequency sub-harmonic suggests that the harmonic forcing for the present conditions is strongest when resonating with
$1/2f_I$
, manifesting in the most receptive forcing conditions A through D, all of which peak near
$f_I$
. Regarding the
$2f_I$
upper bound of receptive frequencies in figure 14(b), we hypothesize that any faster forcing time scales do not support resonance interactions with the larger wavelength structures present in the boundary layer. They would also be ineffective in the turbulent mixing layer where the vortical structures break into smaller, less organized structures (Kim et al. Reference Kim, Elliott and Dutton2020a
), and thus higher frequency dynamics, are not receptive in this limit either.
5. Conclusions
A novel combination of two post-processing methods to analyse time-resolved schlieren data has been advanced and applied to obtain acoustic receptivity and stability insights. A transonic shear layer characterized by
$M=2.46$
and
$M=0.175$
upper and lower streams, respectively, and convective Mach number
$M_c=0.88$
is used as a testbed. The first step is a schlieren-specific application of MPT, which extracts the irrotational component of the flow, thus isolating acoustic wavepackets and propagation within the turbulence. Spectral POD demonstrates the advantages of the filtering; whereas the raw schlieren shows only broadband spectra, the filtered results very clearly display dominant tones. The second component is a DMD-ROM used to test the forced response of the shear layer. This linear framework is efficient, flexible and enables parametric studies of forced responses and energy gains from high-spatio-temporal schlieren data covering long time periods compared with those from scale-resolved simulations.
Results with impulse forcing along a vertical line using both raw and filtered schlieren data to trigger a prototypical shear-layer instability indicate that the perturbation evolves into a wavepacket structure that emits an acoustic disturbance from the splitter-plate trailing edge. The downstream interaction of this wavepacket with a shock produces a bifurcation of structures travelling at two different speeds, inferred to be fast and slow acoustic waves, in addition to acoustic wavefronts which refract into the subsonic stream. In a local receptivity study, the impulse was applied to each pixel in the upper and lower streams to identify the location of the optimal acoustic response. The highest amplification factor is obtained with perturbations in the upper-stream outer boundary layer, where acoustic amplification is almost an order of magnitude larger than that of background acoustic disturbances measured in the lower stream. On the other hand, instabilities triggered in the lower-stream boundary layer only yield threefold acoustic gains. A harmonic parametric forcing analysis was also conducted across all pixels in the upper and lower streams. These results reveal multiple tones from the upper stream that are significantly larger than the single tone produced in the lower stream. Connections are made to the impulse forcing, yielding that the natural impulse response
$1/2$
sub-harmonic and second harmonic bound the most receptive range of forcing frequencies. Individual responses of these largest harmonic forcing values demonstrate spatial structures that adopt the unique disturbance signatures of each forcing location. While most of the selected conditions trigger the shear-layer instability such that the acoustic response is amplified, one forcing location and frequency was found to be receptive but produced little acoustic radiation. Results of this type indicate ways by which the DMD-ROM may be used to augment traditional analysis of qualitative schlieren data.
Acknowledgments
The authors are grateful for comments from J.C. Dutton and G.S. Elliott. They also acknowledge support from the Collaborative Center for Aeronautical Sciences and the Office of Naval Research. The first author is also currently supported by an NRC fellowship award at AFRL.
Declaration of interests
The authors declare no conflict of interest.
Appendix A. The DMD-ROM stabilization
The numerical stability of data-driven reduced-order models sometimes poses a challenge, particularly when employed with measurements containing background noise. The problem manifests in an unbounded solution, as is the case for the filtered schlieren data in this work. In the context of DMD, stability characteristics are dependent on the nature of the dynamic operator,
$A$
, or more precisely, the reduced operator
$\tilde {A}$
.
To stabilize the DMD-ROM, a POD rank truncation of the first
$r$
modes is applied to the derivation of
$\tilde {A}$
. This effectively removes spurious low energy modes which can cause unstable eigenvalues. The stabilization is implemented by truncating the output of the singular value decomposition in (2.10), producing
$\unicode{x1D650}_r$
,
$\boldsymbol{\Sigma}_r$
and
$\unicode{x1D651}_r$
, and resulting in the new dynamic operator
\begin{equation} \tilde {\kern-2pt\unicode{x1D63C}}_r=\unicode{x1D650}_r^* \unicode{x1D64C}_2 \unicode{x1D651}_r \boldsymbol{\Sigma}_r^{-1}. \end{equation}
Transient gain of the stabilized DMD model for various POD mode rank truncation levels
$r$
, for the (a) raw and (b) irrotational filtered data.
The behaviour of the truncated ROM is a function of the number of modes adopted, yielding a trade-off between energy retained versus stability. A brief study of the effect of truncation rank
$r$
is presented to determine effects appropriate for the receptivity analysis. Specifically, the first impulse forcing study of § 4 is repeated for values of
$r=10, 100, 1{,}000$
and
$2{,}000$
, with the transient gain used to measure energy behaviour. Results are presented in figure 16 for the (a) raw and (b) filtered schlierens. First, considering the full ROM (
$r=N=2{,}500$
), the (a) raw schlieren is stable while (b) has unbounded growth after only a few timesteps. Figure 16(b) demonstrates how the truncation dampens the solution over time, where after a certain value, the decay rate is approximately inversely proportional to
$r$
. Although the gain decay of the truncated ROM drives the solution towards zero over time, the relative trends in amplitude for different forcing parameters such as location and frequency are preserved. Furthermore, the asymptotic behaviour of the solution may be recovered by normalizing the forced snapshots by the transient gain at every instant in time
$q'(t)/\sigma (t)$
, as shown in figure 11.
(a) Transient gain of the harmonically forced responses identified in § 4.3. (b) Beat frequency associated with the
$29\,\textrm{kHz}$
forced response is independent of the truncation level
$r$
, and is a product of the second forcing harmonic and a residual
$50\,\textrm{kHz}$
frequency, as observed in the transient gain spectra of (c). The spectrum from the acoustic probe demonstrates that the resulting response is not effected by these artefacts.
Finally, the behaviour of the truncated DMD-ROM is considered for the harmonically forced conditions. The transient gain of the most receptive forcing frequencies found in § 4.3 are reported in figure 17(a), they initially grow, then asymptote to constant levels. Under some conditions, such as the
$29\,\textrm{kHz}$
forcing for example, a unique beat frequency appears. This is isolated in figure 17(b) and verified using different DMD mode truncation levels (r = 100, 1000 and 1500). For all cases, the beat frequency persists, indicating this is not dependent on the truncation of the DMD model. More insight is obtained from the forced-response spectra of figure 17(c), calculated using the same transient gain signal (
$r=1000$
) over the entire available time trace, as well as the acoustic probe spectra taken at the same probe location in figure 11. The beating phenomenon at frequency
$f_b=8\,\textrm{kHz}$
appears in the transient gain spectra, and can be calculated from the other tones
$f_b=f_1-f_2$
, where in this case,
$f_1$
happens to be twice the frequency of the harmonic forcing input, and
$f_2$
is a residual
$50$
kHz tone in the signal. Although this aliasing is a product of the second harmonic and the anomalous
$50$
kHz tone, its genesis is unclear. Nonetheless, the predicted acoustic field is unaffected, as demonstrated by the expected
$29$
kHz tone measured from the harmonic forcing response.
