1. Introduction
Jet screech has been extensively studied due to its potential to cause structural damage to aircraft (Powell Reference Powell1953; Panda Reference Panda1999; Edgington-Mitchell Reference Edgington-Mitchell2019). In this aeroacoustic resonance loop high-amplitude tones can interact with the natural frequencies of aircraft structures, increasing the risk of structural failure (Berndt Reference Berndt1984). Powell (Reference Powell1953) proposed a model to describe the aeroacoustic resonance loop. This model outlines four essential stages: (i) the downstream propagation of energy, (ii) the transformation of this energy into an upstream-travelling wave, (iii) the upstream propagation of energy and (iv) the subsequent conversion back into a downstream-travelling wave, thereby sustaining the resonance cycle. This fundamental mechanism remains central to any resonance process, as it requires the presence of two wave-like structures capable of facilitating energy transport in upstream and downstream directions. The wave propagating downstream has long been identified as the Kelvin–Helmholtz (KH) wavepacket (Mollo-Christensen Reference Mollo-Christensen1967; Crow & Champagne Reference Crow and Champagne1971; Cavalieri et al. Reference Cavalieri, Rodríguez, Jordan, Colonius and Gervais2013; Jordan & Colonius Reference Jordan and Colonius2013), while the upstream-travelling wave has more recently been recognised as a guided jet mode (GJM) (Tam & Hu Reference Tam and Hu1989; Bogey & Gojon Reference Bogey and Gojon2017; Gojon, Bogey & Mihaescu Reference Gojon, Bogey and Mihaescu2018; Edgington-Mitchell et al. Reference Edgington-Mitchell, Jaunet, Jordan, Towne, Soria and Honnery2018; Nogueira et al. Reference Nogueira, Cavalieri, Martini, Towne, Jordan and Edgington-Mitchell2024).
Recent research has significantly advanced our understanding of the screech mechanism (Edgington-Mitchell et al. Reference Edgington-Mitchell, Wang, Nogueira, Schmidt, Jaunet, Duke, Jordan and Towne2021; Nogueira et al. Reference Nogueira, Jaunet, Mancinelli, Jordan and Edgington-Mitchell2022a ,Reference Nogueira, Jordan, Jaunet, Cavalieri, Towne and Edgington-Mitchell b ; Edgington-Mitchell et al. Reference Edgington-Mitchell, Li, Liu, He, Wong, Mackenzie and Nogueira2022). Nogueira et al. (Reference Nogueira, Jordan, Jaunet, Cavalieri, Towne and Edgington-Mitchell2022b ) demonstrated that screech can be predicted solely by considering the periodicity of the mean flow, without explicitly considering this periodicity as arising from shocks. Nogueira et al. (Reference Nogueira, Jaunet, Mancinelli, Jordan and Edgington-Mitchell2022a ) demonstrated that different regions of quasi-periodicity could be responsible for the screech closure; these different regions are represented in the spectral domain as sub-optimal peaks in the streamwise-wavenumber spectrum calculated from the time-averaged flow. Edgington-Mitchell et al. (Reference Edgington-Mitchell, Li, Liu, He, Wong, Mackenzie and Nogueira2022) extended this analysis across a broad parameter space, showing that screech arises from triadic interactions between the KH wavepacket and locally quasi-periodic regions of the flow associated with the shock structure. These interactions energise a range of wavenumbers, including the GJM, which propagates upstream to complete the resonance loop. However, these studies, and to the best of our knowledge all existing models, only consider the mean flow once screech has already reached a limit-cycle state. By the time limit-cycle state is established, the mean flow has already been altered by screech; thus it is difficult to ascertain what role the resonance loop has played in modifying the effective resonance-free mean flow. Screech-prediction models that utilise a modelled mean flow typically only assume a one-way coupling between the waves involved, neglecting how the resonance loop may have distorted the mean flow in order for the jet to reach its limit cycle. While this effect may not be as important for jets with a single resonance loop, it is a crucial part of the jet dynamics when multiple, mutually exclusive or intermittent tones are present in the flow.
As screech is almost always studied in its established limit-cycle state, it is not straightforward to interrogate the relationship between a given screech mode and the turbulent mean atop which it sits. To delineate between the long-time mean with and without screech requires a means by which screech itself can be switched ‘on’ and ‘off’ in the flow. While techniques to cancel out screech without otherwise modifying the shear layer do exist in the form of carefully placed acoustic reflectors (Nagel, Denham & Papathanasiou Reference Nagel, Denham and Papathanasiou1983), these cannot be activated and de-activated on time scales relevant to the flow dynamics. An alternative is to find cases where the jet exhibits intermittent behaviour, with the resonance switching itself on or off due to some other mechanism. One possible mechanism is multimodality; some screech processes appear to be mutually exclusive, and thus a competition between them produces intermittency (Mancinelli et al. Reference Mancinelli, Jaunet, Jordan and Towne2019; Wong et al. Reference Wong, Stavropoulos, Beekman, Towne, Nogueira, Weightman and Edgington-Mitchell2023). Thus to examine mean-flow distortion by the screech loop, we seek cases where jet screech is multimodal, or intermittent for some other reason. Multimodality has been reported for screeching jets arising from both axisymmetric (Gao & Li Reference Gao and Li2010) and non-axisymmetric (Karnam, Saleem & Gutmark Reference Karnam, Saleem and Gutmark2023) nozzles. In axisymmetric jets, the occurrence of multiple screech tones has historically been observed near transitions between different azimuthal modes (Powell, Umeda & Ishii Reference Powell, Umeda and Ishii1992). Recent investigations into axisymmetric jets indicate that, under certain conditions, multiple screech tones can exhibit similar azimuthal modes. In such cases, the KH wavepacket interacts with distinct peaks in the shock-cell spectra corresponding to each screech tone (Edgington-Mitchell et al. Reference Edgington-Mitchell, Li, Liu, He, Wong, Mackenzie and Nogueira2022). For non-axisymmetric jets, shifts in equivalent symmetry and variations in the shock-cell spectral peaks have also been reported in the presence of multiple screech modes (Raman Reference Raman1997; Mazharmanesh et al. Reference Mazharmanesh, Nogueira, Weightman and Edgington-Mitchell2025). Additionally, research has explored whether multiple screech modes occur simultaneously or exclusively. While some studies have concluded that multiple screech modes do not coexist (Raman Reference Raman1997; Mancinelli et al. Reference Mancinelli, Jaunet, Jordan and Towne2019), others have provided evidence that screech tones can be present simultaneously across a range of operating conditions (Walker, Gordeyev & Thomas Reference Walker, Gordeyev and Thomas1997; Shen & Tam Reference Shen and Tam2002). In this work, we will exploit cases with mutually exclusive tones to isolate short-time mean flows associated with each resonance loop, and study them separately to the long-time mean.
The link between mean-flow distortion and triadic interactions has been highlighted in the recent literature (Schmidt Reference Schmidt2020; Kinjangi & Foti Reference Kinjangi and Foti2023). Triadic interactions arise due to the quadratic nonlinearity present in the Navier–Stokes equations. They serve as the fundamental mechanism for energy transfer within fluid dynamics, and are represented in Fourier space as a triplet of wavenumbers or frequencies that sum to zero (Kraichnan Reference Kraichnan1971; Waleffe Reference Waleffe1992). This zero-sum condition is given by
where
$\textit{St}_{\!j}=\textit{St}_k+\textit{St}_l$
represents sum interactions, and
$\textit{St}_{\!j}=\textit{St}_k-\textit{St}_l$
corresponds to difference interactions represented as function of Strouhal number
$\textit{St}$
. A particular case of difference interaction, known as difference self-interaction, occurs when the resultant frequency,
$\textit{St}_{\!j}$
, is zero. The mode associated with this zero-frequency bin characterises mean-flow distortion; it represents how a flow structure at a given frequency, interacting with itself, modifies the temporal mean (zero frequency). Figure 1 presents a schematic of the sum and difference interactions, based on the illustration reported by Schmidt (Reference Schmidt2020), providing a conceptual visualisation of these processes. With the advent of recent techniques such as bispectral mode decomposition (Schmidt Reference Schmidt2020), we now have suitable tools to study mean-flow distortion in a quantitative sense. To date, however, these tools have not been applied to screeching jets, where the mean-flow distortion is of both particular significance and unusual complexity. The selection of the screech mode is dictated by the shock-cell spacing, but the mean-flow distortion produced by the resonance loop in turn modifies the shock-cell spacing.
Illustration of typical frequency triads: (a) sum interaction; (b) difference interaction; (c) harmonic generated by sum self-interaction; (d) mean-flow deformation generated by difference self-interaction.

The primary objectives of this study are twofold: first, to examine how screech acts to modify the mean flow, and second, to explicitly link the wavenumber of the quasi-periodic flow responsible for closing screech to a localised region in the spatial domain. To that end, we systematically investigate screech modes and their sum and difference interactions in jets exhibiting intermittency and multimodality. High-speed schlieren visualisations are analysed via spectral proper orthogonal decomposition (SPOD), SPOD-based frequency–time analysis and bispectral mode decomposition (BMD) to examine the closure mechanisms of screech modes and the triadic interactions between screech modes themselves. The paper is structured as follows: the experimental set-up and analysis techniques are detailed in § 2. Section 3 identifies mean-flow distortion in jets exhibiting intermittent behaviour. For this part of the study, two different nozzle geometries are considered: an axisymmetric jet and an elliptical jet. Section 3 then examines the influence of screech on mean-flow distortion. Section 4 presents the triadic interactions responsible for mean-flow distortion. Section 5 explores the two-way coupling between screech and the mean flow, with a focus on the elliptical jet. Finally, the conclusions are summarised in § 6.
2. Method and theory
2.1. Experimental database
The experimental database is a series of underexpanded screeching jets measured in the Supersonic Gas-Jet Facility at the Shock Lab, Monash University. A schematic of the supersonic-jet set-up is shown in figure 2(a). The system is supplied with a continuous flow of high-pressure air with nozzle-pressure ratios (
$\mathrm{NPR}$
s) ranging from 2.0 to 5.0, wherein
$\mathrm{NPR}$
represents the ratio of plenum to ambient pressure. To condition the flow prior to its entry into the nozzle, a series of screens and honeycombs are employed within the plenum chamber. The two nozzle geometries considered are presented in table 1: a purely converging axisymmetric nozzle and a purely converging elliptical nozzle (Mazharmanesh et al. Reference Mazharmanesh, Nogueira, Weightman and Edgington-Mitchell2025).
Observations of the supersonic jet are conducted using a Z-type Toepler schlieren system comprising two mirrors with a focal length of 2032 mm to establish a collimated light path through the test section. Images are captured using a Photron Fastcam SA-Z 2100k operating at
$150\,000$
frames per second with a resolution
$512\times 200$
pixels. Illumination is provided by a pulsed light-emitting diode (LED) (Willert, Mitchell & Soria Reference Willert, Mitchell and Soria2012) with a pulse width of
$0.5\unicode{x03BC} s$
to capture
$400\,000$
images containing light intensities proportional to path-integrated density gradient in the streamwise direction,
$(\partial \rho _l/\partial x)$
, per operating condition. Figures 2(b) and 2(c) indicate which semi-axis is aligned with the line of sight of the schlieren measurements. These figures present exemplar instantaneous schlieren images of the axisymmetric jet at
$\mathrm{NPR}=2.20$
and the elliptical jet operating at
$\mathrm{NPR}=3.60$
, viewed in the minor-axis plane.
Nozzle geometries used in this study. Here,
$D_e$
is the equivalent circular diameter at the nozzle exit, and
$l/D_e$
is the mean non-dimensional lip thickness.

Table 1. Long description
A table comparing nozzle geometries with their equivalent circular diameters and mean non-dimensional lip thicknesses. The table has two rows and four columns. The columns are labeled Exit geometry, Internal geometry, D subscript e (mm), and l over D subscript e. The first row shows Axisymmetric, Converging, 16, and 0.31. The second row shows Elliptical, Converging, 10, and 0.05.
(a) Schematic of the supersonic-jet experimental facility. (b) Instantaneous schlieren image of an axisymmetric jet operating at
$\mathrm{NPR}=2.20$
. (c) Instantaneous schlieren image of an elliptical jet operating at
$\mathrm{NPR}=3.60$
, viewed in the minor-axis plane.

2.2. Data decomposition
2.2.1. Spectral proper orthogonal decomposition
In this work, we utilise two modal-decomposition techniques, SPOD and BMD, to analyse the schlieren data. The SPOD (Lumley Reference Lumley1967, Reference Lumley1970; Towne, Schmidt & Colonius Reference Towne, Schmidt and Colonius2018) is a method that isolates energy-ranked structures at each frequency. It has demonstrated efficacy in delineating coherent structures as function of frequency within turbulent flows (Gudmundsson & Colonius Reference Gudmundsson and Colonius2011; Sinha et al. Reference Sinha, Rodríguez, Brès and Colonius2014; Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Bres2018; Sano et al. Reference Sano, Abreu, Cavalieri and Wolf2019; Abreu et al. Reference Abreu, Tanarro, Cavalieri, Schlatter, Vinuesa, Hanifi and Henningson2021). This technique will be employed in the present study to identify the screech modes
$\xi$
and corresponding SPOD eigenvalues
$\lambda$
at the screech frequency. Here, we apply SPOD on the light-intensity data derived from schlieren. It is important to note that the magnitude of the SPOD eigenvalues reflects the extent of light-intensity fluctuations rather than any intrinsic flow variable. A brief overview of the method is provided here, while a detailed discussion of the mathematical derivations and algorithmic implementation of SPOD can be found in Towne et al. (Reference Towne, Schmidt and Colonius2018).
In statistically stationary flows, consider
${ \unicode{x1D666}}_{\!j} = { \unicode{x1D666}}(t_{\!j})$
as the mean-subtracted snapshots, where
$j = 1, 2, \ldots , n_t$
indexes the number of snapshots. For spectral estimation, the dataset is segmented into
$n_{\textit{blk}}$
overlapping blocks, each containing
$n_{\textit{fft}}$
snapshots. Adjacent blocks overlap by
$n_{\textit{ovlp}}$
snapshots, with
$n_{\textit{ovlp}} = n_{\textit{fft}}/2$
for the present case (
$50\,\%$
overlap). Each of the
$n_{\textit{blk}}$
blocks undergoes a Fourier transform in time, and all Fourier realisations at the
$l$
th frequency, denoted as
${ \unicode{x1D666}}^{(i)}_l$
, are organised into a matrix
$\hat { \unicode{x1D64C}_l} = [\begin{smallmatrix} \hat {{ \unicode{x1D666}}}^{(1)}_l,&\hat {{ \unicode{x1D666}}}^{(2)}_l,& \ldots ,&\hat {{ \unicode{x1D666}}}^{(n_{\textit{blk}})}_l \end{smallmatrix}]$
. The SPOD eigenvalues
$\boldsymbol{\varLambda }_l$
are determined by solving the eigenvalue problem
$({1}/{n_{\textit{blk}}}) \hat {\unicode{x1D64C}}_l^* \unicode{x1D652} \hat {\unicode{x1D64C}}_l \boldsymbol{\varPsi }_l = \boldsymbol{\varPsi }_l \boldsymbol{\varLambda }_l$
, where
$ \unicode{x1D652}$
is a positive–definite Hermitian matrix accounting for the component-wise and numerical quadrature weights, and
$( \boldsymbol{\cdot })^*$
indicates the complex conjugate. The SPOD modes for the
$l$
th frequency are then obtained as:
$\boldsymbol{\varXi }_l = ({1}/{\sqrt {n_{\textit{blk}})}} \hat {\unicode{x1D64C}}_l \boldsymbol{\varPsi }_l \boldsymbol{\varLambda }_l^{-1/2}$
. The eigenvalues
$\boldsymbol{\varLambda }_l = \text{diag}(\lambda _l^{(1)}, \lambda _l^{(2)}, \ldots , \lambda _l^{(n_{\textit{blk}})})$
represent the energies of the corresponding SPOD modes, with
$\lambda _l^{(1)} \geqslant \lambda _l^{(2)} \geqslant {\cdots} \geqslant \lambda _l^{(n_{\textit{blk}})}$
. The SPOD modes, denoted by the columns of the matrix
$\boldsymbol{\varXi }_l = [\boldsymbol{\xi }_l^{(1)}, \boldsymbol{\xi }_l^{(2)}, \ldots , \boldsymbol{\xi }_l^{(n_{\textit{blk}})}]$
, represent the coherent structures in the flow at the
$l$
-h frequency, with
$\boldsymbol{\xi }_l^{(i)}$
being the
$i$
th mode.
A more recent extension of SPOD is in its time–frequency form, producing a spectrogram that captures the temporal evolution of spatio-temporally coherent flow structures (Nekkanti & Schmidt Reference Nekkanti and Schmidt2021). Specifically, we leverage the ability of the time-continuous SPOD expansion coefficients to provide insights into the global evolution of modal structures. For computational efficiency, we determine the expansion coefficients
for the time index,
$j$
, with the Fourier transform
$\hat {\unicode{x1D64C}}_l$
taken one snapshot at a time (Nekkanti & Schmidt Reference Nekkanti and Schmidt2021). The magnitudes of the complex-valued expansion coefficients,
$|\mathrm{\sum _{i=1}^{n_{\textit{mode}}}}a_i(\textit{St},t)|$
, are the time-dependent SPOD amplitudes. As the SPOD expansion coefficients are complex quantities for non-zero frequencies, we focus on their absolute values in § 3. The SPOD expansion coefficients at each
$\mathrm{NPR}$
are individually normalised by their respective maxima over all frequencies, such that the range of values at each condition lies between 0 and 1. For
$\textit{St} = 0.0$
,
$a_i(\textit{St},t)$
is a real value, and its normalised value varies between –1 and 1. To increase the time resolution the spectral-estimation parameters used for this study are
$n_{\textit{fft}}=1024$
and
$n_{\textit{ovlp}}=512$
, resulting in
$n_{\textit{blk}}=780$
SPOD modes for each frequency. The dataset used herein has
$n_t=400\,000$
snapshots.
2.2.2. Bispectral mode decomposition
The BMD method, introduced by Schmidt (Reference Schmidt2020), isolates the flow structures resulting from triadic interactions, which are the fundamental mechanism for energy transfer in fluid dynamics. This technique detects the primary triads in the flow by maximising the spatially integrated bispectrum. Additionally, it identifies the spatial areas where nonlinear interactions among the coherent structures occur. In this study, we will utilise the BMD method to analyse the extent and intensity of nonlinear interactions among the screech modes and their impact on the distortion of the mean flow. Comprehensive mathematical derivations and algorithmic details of original implementation of BMD can be found in Schmidt (Reference Schmidt2020). In this study, we employ the BMD method with bicoherence normalisation, as recently developed by Yeung & Schmidt (Reference Yeung and Schmidt2025). This enhanced implementation utilises the iterative algorithm of Mengi & Overton (Reference Mengi and Overton2005), ensuring convergence to the global maximum. A brief overview of the method is presented here, while readers are referred to Schmidt (Reference Schmidt2020) and Yeung & Schmidt (Reference Yeung and Schmidt2025) for detailed derivations and a discussion of the mathematical properties of the approach.
The BMD technique extends traditional bispectral analysis to multidimensional data. In classical bispectral analysis, the bispectrum is defined as the double Fourier transform of the third moment of a time signal. For a time series
$y(t)$
with zero mean, the bispectrum is
where
$R_{yyy}(\tau _1, \tau _2) = \mathrm{E} [ y(t) y(t - \tau _1) y(t - \tau _2) ]$
is the third moment of
$y(t)$
, and
$\mathrm{E}[\boldsymbol{\cdot }]$
is the expectation operator. The bispectrum is a signal-processing technique for analysing one-dimensional time series, focusing on quadratic phase coupling at a single spatial point. In comparison, BMD evaluates the strength of triadic interactions across a given spatial domain and identifies the associated spatially coherent structures. Additionally, BMD generalises the classical bispectrum to flow fields,
$q(\boldsymbol{x},t)$
, by defining the spatially integrated bispectrum
Here,
$\hat {{ \unicode{x1D666}}}$
is the temporal Fourier transform of flow data
${ \unicode{x1D666}}$
,
$( \boldsymbol{\cdot })^*$
indicates the complex conjugate and
$\circ$
denotes the Hadamard (or elementwise) product such that
$\hat {{ \unicode{x1D666}}}_{k \circ l} = \hat {{ \unicode{x1D666}}}_k \circ \hat {{ \unicode{x1D666}}}_l$
. The weighted inner product is defined as
$\langle { \unicode{x1D666}}_1, { \unicode{x1D666}}_2\rangle _{\boldsymbol{x}}={ \unicode{x1D666}}^*_{2} \unicode{x1D652}{ \unicode{x1D666}}_1$
, where the diagonal matrix
$ \unicode{x1D652}$
contains the numerical quadrature weights. In the refined formulation of the BMD method, the objective is to evaluate the expected quadratic phase coupling between the components
$\hat {{ \unicode{x1D666}}}_{k+l}$
and
$\hat {{ \unicode{x1D666}}}_{k \circ l}$
, independently of the energy of each component. To achieve this, each component is normalised by the square root of its average integral energy, namely
$\sqrt {\mathrm{E} \{{||\hat {{ \unicode{x1D666}}}_{k+l}||_{\boldsymbol{x}}^2} \}}$
and
$\sqrt {\mathrm{E} \{{||\hat {{ \unicode{x1D666}}}_{k \circ l}||_{\boldsymbol{x}}^2}\}}$
, respectively, where the norm is defined as
${||{{ \unicode{x1D666}}}||_{\boldsymbol{x}}^2}=\langle { \unicode{x1D666}}, { \unicode{x1D666}} \rangle _{\boldsymbol{x}}$
.
For each triad, the estimation of the BMD mode and corresponding mode bispectrum is performed using
$n_{\textit{blk}}$
statistically independent realisations of the Fourier modes. These realisations are assembled into weighted data matrices as follows:
\begin{equation} \hat { \kern-2pt\unicode{x1D64C}}_{k \circ l} = \frac { \unicode{x1D652}^{1/2}\begin{bmatrix} \hat {{ \unicode{x1D666}}}^{(1)}_{k \circ l},&\hat {{ \unicode{x1D666}}}^{(2)}_{k \circ l},& \ldots ,&\hat {{ \unicode{x1D666}}}^{(n_{\textit{blk}})}_{k \circ l} \end{bmatrix} }{\left \lVert \unicode{x1D652}^{1/2}\begin{bmatrix} \hat {{ \unicode{x1D666}}}^{(1)}_{k \circ l},&\hat {{ \unicode{x1D666}}}^{(2)}_{k \circ l},& \ldots ,&\hat {{ \unicode{x1D666}}}^{(n_{\textit{blk}})}_{k \circ l} \end{bmatrix}\right \rVert _{F}}, \end{equation}
and
\begin{equation} \hat { \kern-2pt\unicode{x1D64C}}_{k + l} = \frac { \unicode{x1D652}^{1/2}\begin{bmatrix} \hat {{ \unicode{x1D666}}}^{(1)}_{k + l},&\hat {{ \unicode{x1D666}}}^{(2)}_{k + l},& \ldots ,&\hat {{ \unicode{x1D666}}}^{(n_{\textit{blk}})}_{k + l} \end{bmatrix} }{\left \lVert \unicode{x1D652}^{1/2}\begin{bmatrix} \hat {{ \unicode{x1D666}}}^{(1)}_{k + l},&\hat {{ \unicode{x1D666}}}^{(2)}_{k + l},& \ldots ,&\hat {{ \unicode{x1D666}}}^{(n_{\textit{blk}})}_{k + l} \end{bmatrix}\right \rVert _{F}}, \end{equation}
where each matrix is normalised such that
$||\hat { \kern-2pt\unicode{x1D64C}}_{k \circ l}||_{F}=||\hat { \kern-2pt\unicode{x1D64C}}_{k + l}||_{F}=1$
, with
$||\boldsymbol{\cdot }||_{F}$
denoting the Frobenius norm. We define the cross-frequency field,
${\boldsymbol{\phi }}_{k \circ l} = \unicode{x1D652}^{-1/2}\hat { \kern-2pt\unicode{x1D64C}}_{k \circ l} {\boldsymbol{\alpha }}_{k,l}$
, and bispectral mode,
${\boldsymbol{\phi }}_{k + l} = \boldsymbol{W}^{-1/2}\hat {\boldsymbol{Q}}_{k + l} {\boldsymbol{\alpha }}_{k,l}$
. They are linked through the expansion coefficients,
$\boldsymbol{\alpha }_{k,l}$
, that maximise the magnitude of the mode bispectrum
where
$ \unicode{x1D63D}_{k,l}=\hat { \kern-2pt\unicode{x1D64C}}^*_{k \circ l}\hat { \kern-2pt\unicode{x1D64C}}_{k+l}$
is the bispectral density matrix. The complex mode bispectrum is recovered as
$\beta _{k,l}=\boldsymbol{\alpha }^*_{k,l} \unicode{x1D63D}_{k,l}\boldsymbol{\alpha }_{k,l}$
and satisfies
$|\beta _{k,l}|\leqslant 1$
.
In the BMD technique, three key outputs are of particular importance: the complex mode bispectrum,
$\beta$
, which could be interpreted as strength of the interactions; the bispectral mode,
${\boldsymbol{\phi }}_{k+l}$
, which represents the flow structures arising from nonlinear triadic interactions; and the field
$\boldsymbol{\psi }_{k,l} = |\boldsymbol{\phi }_{k \circ l} \circ \boldsymbol{\phi }_{k +l}|$
, which serves as an interaction map. This interaction map highlights the spatial regions where nonlinear interactions between coherent structures occur.
The SPOD eigenvalue spectra of the leading mode for the axisymmetric jet operating at (a)
$\mathrm{NPR}=2.19$
where the A1 mode dominates; (b)
$\mathrm{NPR}=2.20$
where a transition from the A1 to the A2 mode occurs.
$\textit{St}_{1,2}$
indicates screech frequency; (c)
$\mathrm{NPR}=2.25$
where the A2 mode dominates.

Figure 3. Long description
Three line graphs depict eigenvalue spectra for different jet modes. Panel A: A line graph shows the eigenvalue spectra for an axisymmetric jet where the A1 mode dominates. The x-axis represents the Strouhal number (St) ranging from 0 to 1.5, and the y-axis represents the eigenvalue magnitude (λ) ranging from 10^6 to 10^9. Key peaks are labeled as St1 and 2St1. The Noise Power Ratio (NPR) is 2.19. Panel B: A line graph shows the eigenvalue spectra for an axisymmetric jet where a transition from the A1 to the A2 mode occurs. The x-axis represents the Strouhal number (St) ranging from 10^-2 to 10^1, and the y-axis represents the eigenvalue magnitude (λ) ranging from 10^6 to 10^8. Key peaks are labeled as St2 and St1. The Noise Power Ratio (NPR) is 2.20. Panel C: A line graph shows the eigenvalue spectra for an axisymmetric jet where the A2 mode dominates. The x-axis represents the Strouhal number (St) ranging from 0 to 1.5, and the y-axis represents the eigenvalue magnitude (λ) ranging from 10^6 to 10^9. Key peaks are labeled as St1 and 2St1. The Noise Power Ratio (NPR) is 2.25.
3. Mean-flow distortion
3.1. Intermittency from mutual exclusivity
We begin the analysis by identifying the frequencies associated with screech using the SPOD spectra. In this study, the dimensionless frequency, or Strouhal number, is defined as
$\textit{St}=(( {\textit{fD}_e})/{U_{\!j}})$
, where
$f$
is dimensional frequency,
$D_e$
represents the equivalent diameter (based on matched area) and
$U_{\!j}$
denotes the ideally expanded jet exit velocity. Figure 3 shows the SPOD eigenvalue spectra of the leading mode for the axisymmetric jet at three operating conditions. A common scenario in which jet screech exhibits multimodal behaviour occurs near the mode staging location. Mode staging refers to the tendency of screech to undergo abrupt frequency shifts in response to small variations in jet supply pressure (Powell Reference Powell1953; Merle Reference Merle1957; Powell et al. Reference Powell, Umeda and Ishii1992; Umeda & Ishii Reference Umeda and Ishii2001; Gao & Li Reference Gao and Li2010; Mancinelli et al. Reference Mancinelli, Jaunet, Jordan and Towne2019; Edgington-Mitchell et al. Reference Edgington-Mitchell, Li, Liu, He, Wong, Mackenzie and Nogueira2022; Nogueira et al. Reference Nogueira, Jaunet, Mancinelli, Jordan and Edgington-Mitchell2022a
).
In the present study, the axisymmetric jet is examined over the range
$\mathrm{NPR}=2.19$
to
$\mathrm{NPR}=2.25$
, capturing the multimodal behaviour during the transition from A1 to A2. Three regimes are considered: before staging (mode A1, figure 3
a), during staging (modes A1 and A2, figure 3
b) and after staging (mode A2, figure 3
c). For
$\mathrm{NPR}=2.19$
the distinct peaks in figure 3(a) correspond to the fundamental screech frequency (
$\textit{St}_1=0.607$
) and its first harmonic (
$2\textit{St}_{1}=1.221$
). The spatial structure of the mode is shown in figures 4(a) and 4(b). The screech tone is axisymmetric and associated with the
$A1$
mode prior to mode staging.
Screech modes (
$\textit{St}_1$
and
$\textit{St}_2$
) and the mean-flow distortion mode (
$\textit{St}=0$
) for the axisymmetric jet at
$\mathrm{NPR}=2.19$
, 2.20 and 2.25. (a,c,e,g,i) Real part of SPOD spatial mode. (b,d, f,h, j) Absolute value of SPOD spatial mode.

The SPOD-based frequency–time diagrams derived from the leading mode at each frequency for the axisymmetric jet operating at (a)
$\mathrm{NPR}=2.25$
where the A2 mode dominates. The horizontal axis is normalised by
$T$
, the period of the screech tone at
$\textit{St} = 0.665$
. Panel (b) shows
$\mathrm{NPR}=2.20$
where a transition from the A1 to the A2 mode occurs. The horizontal axis is normalised by
$T$
, the period of the screech tone at
$\textit{St} = 0.694$
. Only the absolute values of
$a_i(\textit{St},t)$
are shown in (a) and (b). (c) The absolute values of the SPOD expansion coefficients at
$\textit{St} = 0.606$
and
$\textit{St} = 0.694$
, along with the real part of the SPOD expansion coefficient at
$\textit{St}=0$
for
$\mathrm{NPR}=2.20$
.

Figure 5. Long description
Panel A: A heatmap displays the absolute values of a variable over time normalized by the period of the screech tone at a specific Strouhal number. The vertical axis is labeled with Strouhal number, and the horizontal axis is normalized time. The color scale ranges from 0 to 1. Panel B: Another heatmap shows the absolute values of the same variable over time normalized by the period of the screech tone at a different Strouhal number. The vertical axis is labeled with Strouhal number, and the horizontal axis is normalized time. The color scale ranges from 0 to 1. Panel C: A line graph shows the absolute values of the SPOD expansion coefficients at two different Strouhal numbers, along with the real part of the SPOD expansion coefficient at another Strouhal number. The horizontal axis is normalized time, and the vertical axis represents the normalized amplitude. Three lines are plotted in different colors representing different Strouhal numbers.
At
$\mathrm{NPR}=2.20$
, figure 3(b) reveals two distinct peaks (
$\textit{St}_1$
,
$\textit{St}_2$
) that are not harmonically related, suggesting multimodal behaviour. The SPOD spatial modes corresponding to the screech frequencies
$\textit{St}_1=0.694$
and
$\textit{St}_2=0.606$
are presented in figure 4(c–h). Both screech modes are symmetric about the centreline and exhibit non-zero amplitude along it, confirming their axisymmetric nature. As demonstrated in the studies of Mancinelli et al. (Reference Mancinelli, Jaunet, Jordan and Towne2019), Edgington-Mitchell et al. (Reference Edgington-Mitchell, Li, Liu, He, Wong, Mackenzie and Nogueira2022), and Nogueira et al. (Reference Nogueira, Jaunet, Mancinelli, Jordan and Edgington-Mitchell2022a
), during A1–A2 mode staging the higher-frequency tone, here
$\textit{St}_1=0.694$
, is associated with the A2 mode of jet screech, while the lower-frequency tone, here
$\textit{St}_2=0.606$
, corresponds to the A1 mode. Note that the A1–A2 mode staging occurs without any change in azimuthal symmetry. Consequently, a very careful spectral analysis is required to distinguish between these two modes. Nogueira et al. (Reference Nogueira, Jaunet, Mancinelli, Jordan and Edgington-Mitchell2022a
) attributed the frequency shift observed between the A1 and A2 modes to differences in the resonance-closing mechanism: the A1 mode is closed through the interaction between the KH wavepacket and the primary shock-cell wavenumber
$(k_{s_{1}})$
, whereas the A2 mode results from interaction with a secondary wavenumber peak
$(k_{s_{2}})$
arising from spatial variations in shock-cell spacing. The frequency–time analysis presented in figures 5(b) and 5(c) further confirms the existence of these two distinct modes. We discuss figure 5 in greater detail later in this section, where the distinction between the A1 and A2 modes becomes particularly evident. The SPOD spatial mode corresponding to the zero-frequency bin is presented in figures 4(g) and 4(h). Since the temporal mean corresponds to zero frequency, the SPOD spatial mode at
$\textit{St} = 0$
is associated with the mean flow. However, because the long-time mean is subtracted in SPOD analysis, any remaining features in the
$\textit{St} = 0$
spatial mode indicate distortions in the flow slower than those represented by the first non-zero frequency, i.e. deviations from the long-time mean. As shown in figures 4(g) and 4 (h), the SPOD spatial mode at
$\textit{St} = 0$
exhibits shock-cell-like structures, with non-zero values along the centreline in both the real part and the absolute value.
With a slight increase in jet supply pressure, the transition from A1 to A2 is completed, and the screech enters the second stage (A2 mode). For
$\mathrm{NPR}=2.25$
the distinct peaks in figure 3(c) correspond to the fundamental screech frequency (
$\textit{St}_1=0.665$
) and its first harmonic (
$2\textit{St}_{1}=1.324$
). The spatial structure of the mode is shown in figures 4(i) and 4(j). The screech tone remains axisymmetric and is associated with the A2 mode after mode staging. Once again, a very careful spectral analysis is required to distinguish between modes A1 and A2. Despite the similarity of their SPOD spatial mode shapes and SPOD eigenvalue spectra, the two modes exhibit distinct wavenumber spectra. The A1 screech mode is associated with the interaction between the KH wavepacket and the fundamental peak in the shock-cell wavenumber spectrum, whereas the A2 mode arises from interaction with the first suboptimal peak. Further details of the spectral domain and its connection to the spatial domain are discussed in § 5.
The presence of two tones in the SPOD spectra could be indicative of either contemporaneous or mutually exclusive processes; we now examine the temporal relationship between these tones by way of SPOD-based frequency–time analysis, via the SPOD expansion coefficients,
$a_i(\textit{St},t)$
. Figures 5(a) and 5(b) presents frequency–time diagrams, where the contour represents the absolute value of the expansion coefficients
$|{a_1(\textit{St},t)|}$
for the leading mode at each frequency. Figure 5(a) indicates that when a single resonance loop is active, the screech tone is strong and remains stable over time. Conversely, in the multimodal case shown in figure 5(b), the screech tones are highly intermittent; there is switching between the A1 and A2 modes, as was observed in Mancinelli et al. (Reference Mancinelli, Jaunet, Jordan and Towne2019). While mutually exclusive screech modes in competition must be associated with intermittency, the converse is not true; intermittency does not automatically imply mutual exclusivity. To determine whether multiple screech modes coexist or switch exclusively, joint-probability-density functions (JPDF) of the absolute value of the SPOD expansion coefficients at the screech modes,
$|{a_1(\textit{St}_s)}|$
, are computed. The temporal evolution of these SPOD expansion coefficients is shown in figure 5(c) in cyan and magenta.
Figure 6(a) presents the JPDF of the SPOD expansion coefficients corresponding to the two distinct screech modes. The figure reveals that the two modes are nearly perfectly mutually exclusive: when the amplitude of
$a_1(\textit{St}_1)$
(or
$a_1(\textit{St}_2)$
) is close to one, the amplitude of
$a_1(\textit{St}_2)$
(or
$a_1(\textit{St}_1)$
) approaches zero. This indicates that when the screech tone at
$\textit{St}_1$
is ‘on’, the tone at
$\textit{St}_2$
is ‘off’, and vice versa. Note that due to the finite temporal resolution of the frequency–time analysis, some windows will represent a transition between the two tones, which is likely why figure 6(a) has a weak region in between the two points of mutual exclusivity in the JPDF. To quantify the correlation between the two screech modes, the Pearson correlation coefficient
$\rho$
is computed using
$|{a_1(\textit{St}_1)}|$
and
$|{a_1(\textit{St}_2)}|$
. The parameter
$\rho$
ranges from −1 to 1, with 1 indicating perfect correlation, −1 indicating perfect anti-correlation and 0 indicating that the two modes are uncorrelated. The calculated correlation coefficient,
$\rho (|{a_1(\textit{St}_1)}|, |{a_1(\textit{St}_2)}|) = -0.82$
, indicates a significant degree of anti-correlation. This suggests that the system switches between two mutually exclusive states while remaining statistically stationary.
(a) The JPDF of two distinct screech modes for the axisymmetric jet operating at
$\mathrm{NPR}=2.20$
. The correlation coefficient is
$\rho (|{a_1(\textit{St}_1)}|, |{a_1(\textit{St}_2)}|) = -0.82$
. (b,c) The JPDFs between the screech modes and the zero-frequency component at
$\mathrm{NPR}=2.20$
. The correlation coefficients are
$\rho (\text{Re} ({a_1(\textit{St}_0)}), |{a_1(\textit{St}_1)}|) = 0.96$
and
$\rho (\text{Re} ({a_1(\textit{St}_0)}), |{a_1(\textit{St}_2)}|) = -0.82$
, where
$\textit{St}_0$
denotes the zero-frequency bin,
$\text{Re} (\boldsymbol{\cdot })$
represents the real part of the SPOD expansion coefficient, and
$a_1(\boldsymbol{\cdot })$
denotes the SPOD expansion coefficient from the first leading mode.

Figure 6. Long description
Panel A: A contour plot shows the joint probability density function (JPDF) of two distinct screech modes for an axisymmetric jet operating at a specific condition. The x-axis represents the magnitude of the SPOD expansion coefficient for the first screech mode, and the y-axis represents the magnitude for the second screech mode. The correlation coefficient is -0.82. Panel B: Another contour plot displays the JPDF between the screech mode and the zero-frequency component. The x-axis shows the real part of the SPOD expansion coefficient for the zero-frequency bin, and the y-axis shows the SPOD expansion coefficient from the first leading mode. The correlation coefficient is 0.96. Panel C: The third contour plot illustrates the JPDF between the screech mode and the zero-frequency component at a different condition. The x-axis represents the real part of the SPOD expansion coefficient for the zero-frequency bin, and the y-axis represents the magnitude of the SPOD expansion coefficient for the screech mode. The correlation coefficient is -0.82.
Comparison between the long-time mean and short-time mean for the axisymmetric jet at
$\mathrm{NPR} = 2.20$
, with the two regions separated by the red dashed line. The upper section above the dashed line represents the base flow (
$\bar {I}$
), while the lower section shows the short-time mean flow modified by (a)
$\textit{St}_1 = 0.694$
and (b)
$\textit{St}_2 = 0.606$
. Here,
$a_{\mathit{max}}$
and
$a_{\mathit{min}}$
denote the maximum and minimum values of
${a_1(\textit{St} = 0)}$
corresponding to
${a_1(\textit{St} = 0)} \gt 0$
and
${a_1(\textit{St} = 0)} \lt 0$
, respectively. The directory including the reconstructed flow field video and the accompanying Jupyter notebook can be accessed at https://www.cambridge.org/S0022112026117960/JFM-Notebooks/files/Figure7/Figure7.ipynb.

Another noteworthy feature in figure 5(b) is the high SPOD expansion coefficient associated with
$\textit{St}=0$
for
$\mathrm{NPR}=2.20$
. Since a mode in the zero-frequency bin represents mean-flow distortion, it can be concluded that the SPOD expansion coefficient at
$\textit{St}=0$
in figure 5(b) represents the detectable mean-flow distortion generated by the difference self-interaction of fundamental modes (Schmidt Reference Schmidt2020), as will be discussed in § 4. The presence of intermittency, whether due to mutual exclusivity, or simply each modes own intermittency, is required for the zero-frequency bin to contain information resulting from the distortion of the mean flow by the screech process. When the amplitude of the screech mode remains stable with time, its difference self-interaction is entirely contained within the long-time mean, and thus does not appear in the zero-frequency bin. This is reflected in the low SPOD expansion coefficient at
$\textit{St}=0$
in figure 5(a). To assess the relationship between the two screech modes and the mean-flow distortion, the JPDF and the Pearson correlation coefficient are computed between the absolute value of the SPOD expansion coefficient at the screech modes (
$|{a_1(\textit{St}_s)}|$
) and the real part of the SPOD expansion coefficient at
$\textit{St} = 0$
(
$\text{Re} ({a_1(\textit{St}_0)})$
). The temporal evolution of these SPOD expansion coefficients is shown in figure 5(c). As discussed in § 2, the normalised value of
$\text{Re} ({a_1(\textit{St}_0)})$
ranges between –1 and 1. Figures 6(b) and 6(c) illustrate the correlation between mean-flow distortion and both screech modes. In interpreting these results, it is important to be clear that the long-time mean one might calculate for this flow is itself an average of two different mean-flow states; the long-time mean itself does not correspond to any mean-flow structure present in the flow. There is a short-time mean associated with
$|{a_1(\textit{St}_1)}|$
, and another associated with
$|{a_1(\textit{St}_2)}|$
; there is no short-time mean associated with both modes being active, as this state does not occur in this system, nor is there a state in this dataset associated with a complete cessation of both tones. As shown in figures 5(c) and 6(b) and 6(c),
$|\text{Re} ({a_1(\textit{St}_0)})| \geqslant 0.6$
for most points in time, which indicates that the short-time mean almost never resembles the long-time mean. The decomposition produced by SPOD could represent these two possible short-time mean-flow states with a single eigenmode, or with the leading mode and several suboptimal modes. In this case, it appears that the leading spatial mode at
$\textit{St} = 0$
is sufficient to reconstruct a virtually unmodified base flow, when
${a_1(\textit{St}=0)}\approx 0$
, and two deviations from it associated with the screech loops:
${a_1(\textit{St}=0)} \gt 0$
for modification by the
$\textit{St}_1=0.694$
loop or
${a_1(\textit{St}=0)} \lt 0$
for modification by the
$\textit{St}_2=0.606$
loop (see figure 7). Comparison of figures 7(a) and 7(b) indicates that the modification of the shock-cell structure relative to the long-time mean is stronger when the screech tone
$\textit{St}_2=0.606$
is active. This is due to the larger value of
$|{a_1(\textit{St}=0)|}$
when
$\textit{St}_2=0.606$
is active, as shown in figure 5(c). A video constructed from the reconstructed flow field has been linked to figure 7 via a Jupyter notebook to illustrate the temporal evolution of the mean-flow distortion. The link between the different states of mean-flow distortion and the screech modes is clear from the JPDFs in figures 6(b) and 6(c). When
$|{a_1(\textit{St}_1)}|$
is high the mean-flow distortion is strong, as indicated by the peak in the upper-right corner of figure 6(b), and is associated with a positive value of
$\text{Re} ({a_1(\textit{St}_0)})$
. Since the two screech tones are anti-correlated, this peak also appears in the lower-right corner of figure 6(c), corresponding to instances when
$|{a_1(\textit{St}_2)}\approx 0|$
. Furthermore, when
$|{a_1(\textit{St}_2)}|$
is high the mean-flow distortion is also strong, as indicated by the peak in the upper-left corner of figure 6(c), but it is associated with negative values of
$\text{Re} ({a_1(\textit{St}_0)})$
. Again, due to the anti-correlated relationship between the two screech tones, the mean-flow distortion induced by the second tone
$(\textit{St}_2)$
appears in the lower-left corner of figure 6(b), corresponding to instances when
$|{a_1(\textit{St}_1)}\approx 0|$
. The high values of correlation coefficient
$|\rho | \geqslant 0.8$
indicate a very strong correlation between the screech tones and the modification of the mean flow. Figure 6 represents, to the best of our knowledge, the first clear demonstration that different screech loops at the same condition are associated with different short-time means.
(a) The SPOD-based time–frequency diagrams derived from the leading mode at each frequency for the elliptical jet operating at
$\mathrm{NPR} = 3.60$
. (b) Normalised SPOD expansion coefficients corresponding to two screech modes (
$\textit{St} = 0.314$
and
$\textit{St} = 0.269$
) and the mean-flow distortion mode (
$\textit{St} = 0$
), shown for modes 1 and 2. The horizontal axis is normalised by
$T$
, the period of the screech tone at
$\textit{St} = 0.269$
. The elliptical jet is viewed in the minor-axis plane.

Figure 8. Long description
Panel A: A heatmap displays the SPOD-based time-frequency diagrams derived from the leading mode at each frequency for the elliptical jet operating at a specific condition. The horizontal axis is labeled t/T, representing the normalized time, and the vertical axis is labeled St, representing the Strouhal number. The color scale on the right indicates the magnitude of the normalized SPOD expansion coefficients, ranging from 0 to 1. Panel B: A line graph shows the normalized SPOD expansion coefficients corresponding to two screech modes and the mean-flow distortion mode. The horizontal axis is labeled t/T, normalized by the period of the screech tone at a specific Strouhal number. The vertical axis is labeled a/a_max, representing the normalized amplitude. Different lines represent different modes: St = 0 mode 1, St = 0 mode 2, St = 0.269, and St = 0.314. The graph illustrates the temporal evolution of these modes.
3.2. Intermittency from mode-specific dynamics
The axisymmetric jet at
$\mathrm{NPR}=2.20$
illustrates mean-flow distortion characterised by two distinct mean states, each associated with the activation of the two screech tones. In this case, the leading spatial mode at
$\textit{St} = 0$
adequately captures the unmodified base flow, as well as the two deviations from it. We now turn to a case where there is a more complex temporal dynamics, and the leading spatial mode at
$\textit{St} = 0$
is insufficient to represent the mean-flow modification: for this example we select an elliptical jet. Figure 8 shows a frequency–time analysis for an elliptical jet operating at
$\mathrm{NPR} = 3.60$
, a condition at which two distinct screech tones are evident in the SPOD spectra. As shown in this figure the dominant screech tone at
$\textit{St} = 0.269$
is replaced by a dominant tone at
$\textit{St} = 0.314$
at approximately
$\mathrm{t/T} = 1.45 \times 10^4$
. Before proceeding with further analysis of the elliptical jet, it is important to note that, unlike the axisymmetric case, this system is not statistically stationary, as it undergoes a single transition from one state to another over the duration of the time series. However, comparison of results obtained using the full dataset with those from the pre- and post-switch subsets shows that the modes identified in the full dataset have direct counterparts in the respective subsets. In particular, the SPOD expansion coefficient at
$\textit{St}=0.269$
is identical for the full dataset and the pre-switch dataset, while the coefficient at
$\textit{St}=0.314$
is identical for the full dataset and the post-switch dataset. This confirms that SPOD does not introduce artefacts associated with the transition between modes. Evidence from comparisons between SPOD results obtained using the full dataset and the subsets is presented in Appendix A.
(a) The JPDF of two distinct screech modes for the elliptical jet operating at
$\mathrm{NPR}=3.60$
. The correlation coefficient is
$\rho (|{a_1(\textit{St}_1)}|, |{a_1(\textit{St}_2)|}) = -0.92$
. (b, c) The JPDFs involving screech modes and the zero-frequency bin at
$\mathrm{NPR}=3.60$
. The correlation coefficients are
$\rho (\text{Re} ({a_1(\textit{St}_0))}, |{a_1(\textit{St}_1)}|) = 0.97$
and
$\rho (\text{Re} ({a_1(\textit{St}_0)}), |{a_1(\textit{St}_2)}|) = -0.91$
. (d) The JPDF involving screech mode and the zero-frequency bin at mode 2 at
$\mathrm{NPR}=3.60$
. The correlation coefficient is
$\rho (\text{Re} ({a_2(\textit{St}_0)}), |{a_1(\textit{St}_2)}|) = 0.83$
. Here,
$\textit{St}_0$
denotes the zero-frequency bin,
$\text{Re} (\boldsymbol{\cdot })$
represents the real part of the SPOD expansion coefficient and
$a_{1,2}(\boldsymbol{\cdot })$
denotes the SPOD expansion coefficient from the first and second leading modes. The elliptical jet is viewed in the minor-axis plane.

The SPOD expansion coefficients for the elliptical jet differ from those presented in figure 5 in that here the screech mode associated with
$\textit{St} = 0.314$
exhibits a degree of intermittency. As expected from the time–frequency diagram, the two screech modes are nearly perfectly mutually exclusive, exhibiting a correlation coefficient of
$\rho (|{a_1 (\textit{St}_1)}|, |{a_1(\textit{St}_2)}|) = -0.92$
(see figure 9
a), and two prominent peaks in the JPDF. Unlike the axisymmetric jet, capturing the mean-flow modification in this case requires the two leading spatial modes at
$\textit{St} = 0$
. Mode 1, as in the axisymmetric case, reflects the difference between the short-time mean associated with each screech loop and the long-time mean; once again the long-time mean is an average of the two short-time means, and does not correspond to a flow state that is ever observed. Mode 2, which has no equivalent in the axisymmetric case, exists to account for the intermittent nature of
$\textit{St}_2$
. For
$t/T \leqslant (t/T)_{\textit{switch}}$
, when only
$\textit{St}_1$
is active,
${a_1(\textit{St}=0)} \gt 0$
, representing a contraction of the shock cells compared with the long-time mean. For
$t/T \geqslant (t/T)_{\textit{switch}}$
, when
$\textit{St}_2$
is active,
${a_1(\textit{St}=0)} \lt 0$
, representing an extension of the shock cells relative to the long-time mean. For
$t/T \leqslant (t/T)_{\textit{switch}}$
,
${a_2(\textit{St}=0)}\approx 0$
, but once
$\textit{St}_2$
is active,
${a_2(\textit{St}=0)}\neq 0$
, rapidly oscillating between positive and negative values. For the relationship between the two screech tones, and between the tones and the leading mode at
$\textit{St} = 0$
, the correlations are incredibly strong, as shown in figure 9(a–c). These cross-frequency correlations are indicative of triadic interactions, motivating the use of the BMD and the analysis of the BMD interaction maps, as will be discussed in §§ 4 and 5. The JPDF between
$\textit{St}_2$
and
$a_2(\textit{St}=0)$
demonstrates a strong correlation between the strength of
$\textit{St}_2$
and the degree to which the mean flow is modified by
$a_2(\textit{St}=0)$
; when
$\textit{St}_2$
is weak,
${a_2(\textit{St}=0)}\leqslant 0$
, when it is strong,
${a_2(\textit{St}=0)}\geqslant 0$
, and at roughly its median value,
${a_2(\textit{St}=0)}\approx 0$
. It is worth noting that the coefficient
$a_2(\textit{St}=0)$
obtained from the full dataset is equivalent to the coefficient
$a_1(\textit{St}=0)$
obtained from the post-switch dataset
$(t/T \geqslant (t/T)_{\textit{switch}})$
. Unsurprisingly, there is no correlation between
$\textit{St}_1$
and
$a_2(\textit{St}=0)$
, with a coefficient of
$\rho \lt 0.05$
.
Three screech modes (
$\textit{St}=0.314$
,
$\textit{St}=0.273$
and
$\textit{St}=0.269$
) and the mean-flow distortion mode (
$\textit{St}=0$
) for the elliptical jet at
$\mathrm{NPR}=3.60$
. (a,c,e,g,i) Real part of SPOD spatial mode. (b,d, f,h, j) Absolute value of SPOD spatial mode.

Figure 10. Long description
Panel A: A heatmap shows the real part of the SPOD spatial mode for screech mode 1. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2. The Strouhal number is 0.314. Panel B: A heatmap shows the absolute value of the SPOD spatial mode for screech mode 1. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2. Panel C: A heatmap shows the real part of the SPOD spatial mode for screech mode 2. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2. The Strouhal number is 0.273. Panel D: A heatmap shows the absolute value of the SPOD spatial mode for screech mode 2. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2. Panel E: A heatmap shows the real part of the SPOD spatial mode for screech mode 3. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2. The Strouhal number is 0.269. Panel F: A heatmap shows the absolute value of the SPOD spatial mode for screech mode 3. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2. Panel G: A heatmap shows the real part of the SPOD spatial mode for the mean-flow distortion mode 1. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2. The Strouhal number is 0. Panel H: A heatmap shows the absolute value of the SPOD spatial mode for the mean-flow distortion mode 1. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2. Panel I: A heatmap shows the real part of the SPOD spatial mode for the mean-flow distortion mode 2. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2. The Strouhal number is 0. Panel J: A heatmap shows the absolute value of the SPOD spatial mode for the mean-flow distortion mode 2. The x-axis is labeled x/D_e and ranges from 0 to 12. The y-axis is labeled y/D_e and ranges from -2 to 2.
Another noteworthy feature in figure 8(a) is for
$t/T \geqslant (t/T)_{\textit{switch}}$
, in addition to the higher-frequency screech mode at
$\textit{St} = 0.314$
, a lower-frequency screech mode at
$\textit{St} = 0.273$
also appears. This lower-frequency mode is relatively weaker than the higher-frequency one and exhibits intermittent behaviour. The nature of two screech modes appearing in an elliptical jet operating at
$\mathrm{NPR}=3.60$
were studied in our previous work (Mazharmanesh et al. Reference Mazharmanesh, Nogueira, Weightman and Edgington-Mitchell2025), but here with a significantly longer time sample, three distinct modes are apparent. The SPOD spatial modes corresponding to the screech frequencies
$\textit{St}=0.269$
,
$\textit{St}=0.273$
and
$\textit{St}=0.314$
, are reproduced in figure 10(a–f). The SPOD analysis was performed using the full dataset,
$0\leqslant t/T \leqslant 3 \times 10^4$
, as well as on two subsets:
$t/T \leqslant (t/T)_{\textit{switch}}$
and
$t/T \geqslant (t/T)_{\textit{switch}}$
. The spatial mode corresponding to
$\textit{St}=0.269$
is identical for the full dataset and the first subset. Similarly, the spatial modes associated with
$\textit{St}=0.273$
and
$\textit{St}=0.314$
are the same when comparing the full dataset with the second subset,
$t/T \geqslant (t/T)_{\textit{switch}}$
. For all screech modes, the real part of the spatial mode is antisymmetric about the centreline, while the absolute value is zero at the centreline when the jet is viewed in the minor-axis plane. This pattern is consistent with either a helical or flapping motion. Our previous work demonstrated that the screech mode,
$\textit{St}=0.273$
, exhibits characteristics of a quasi-helical mode, whereas the screech modes,
$\textit{St}=0.269$
and
$\textit{St}=0.314$
, closely resemble a flapping mode. The SPOD spatial modes at the zero-frequency bin for the first two leading modes are also shown in figure 10(g–j). Here, the SPOD is performed only on the full dataset, as mean-flow modification is not observable when using the first subset. As illustrated in this figure, the SPOD spatial mode at
$\textit{St}=0$
, whether in mode 1 or mode 2, features shock-cell-like structures with non-zero values along the centreline in both the real and absolute parts. Further examples of the relationship between screech and mean-flow distortion are presented in Appendix B, including an example where the
$\textit{St}_0$
mode appears to be linked to slow variation in supply pressure, as opposed to the rapid switching observed here.
4. Triadic interactions
Having established clear correlations between different screech modes and different short-time mean-flow states, we now use BMD to interrogate the specific triadic interactions involved. For brevity, we will consider only the aforementioned elliptical jet at
$\mathrm{NPR}=3.60$
in the remainder of the paper. Before discussing the BMD results in detail, we must re-emphasise that the decomposition is being performed on schlieren data, thus the bicoherence mode bispectrum (
$\beta$
) is only representative of the relative strength of path-integrated density fluctuations. As a corollary to this, BMD in this form cannot identify the direction of energy transfers, only that transfers occur between specific frequencies and symmetries; a newer form of BMD can identify the direction of energy transfer in triadic interactions, but requires volumetric velocity and density information (Nekkanti et al. Reference Nekkanti, Pickering, Schmidt and Colonius2025). Furthermore, as with the SPOD analysis, BMD was performed on the full dataset (
$0 \leqslant t/T \leqslant 3 \times 10^4$
) as well as on two subsets:
$t/T \leqslant (t/T)_{\textit{switch}}$
and
$t/T \geqslant (t/T)_{\textit{switch}}$
. It is found that the bispectral modes and interaction maps computed from the full dataset are consistent with those from the first subset (
$t/T \leqslant (t/T)_{\textit{switch}}$
) for
$\textit{St}=0.269$
and, with those from the second subset (
$t/T \geqslant (t/T)_{\textit{switch}}$
) for
$\textit{St}=0.273$
and
$\textit{St}=0.314$
. Evidence from comparisons between BMD results obtained using the full dataset and the subsets is presented in Appendix A.
The relative strength of triadic interactions are depicted using the mode bispectrum
$\beta$
, shown in figure 11. The upper region, denoted by
$\textit{St}_l\gt 0$
and
$\textit{St}_k\gt 0$
, signifies the sum-interaction region (
$\textit{St}_{\!j}=\textit{St}_k+\textit{St}_l$
), while the lower region, denoted by
$\textit{St}_l\lt 0$
and
$\textit{St}_k\gt 0$
, indicates the difference-interaction region (
$\textit{St}_{\!j}=\textit{St}_k-\textit{St}_l$
). Figure 11 compares the BMD mode bispectra interactions when the mean is removed or included. The bispectrum with the mean included shows higher magnitudes along
$\textit{St}_l=0$
and
$\textit{St}_{\!j}=0$
than the bispectrum without the mean. The mean-included bispectrum also displays local maxima along
$\textit{St}_l=0$
. The remaining regions of the bispectrum are independent of mean inclusion or removal. Difference self-interaction with
$\textit{St}_{\!j}=0$
represents the mean-flow distortion when the long-time mean is removed from data, and the long-time mean itself when it is included in the data. The two peaks (black dots) along the line
$\textit{St}_l=0$
in figure 11 align with the screech frequencies
$\textit{St}_1=0.269$
and
$\textit{St}_2=0.314$
shown in figure 10, respectively. Although the triadic interaction
$\textit{St}_{\!j}=0.273-0.0=0.273$
is present in figure 11 it does not appear as a separate peak (black dot) because
$\textit{St}=0.273$
is very close to
$\textit{St}_1=0.269$
. In what follows, we perform BMD exclusively with the mean included.
Mode bispectrum of the elliptical jet operating at
$\mathrm{NPR}=3.60$
: (a) the long-time mean is removed from the data; (b) the long-time mean is included. Higher values of
$\beta$
indicate the most dominant triadic interactions, defined by
$\textit{St}_{\!j} = \textit{St}_k \pm \textit{St}_l$
.

Bispectral modes of the elliptical jet operating at
$\mathrm{NPR}=3.6$
, corresponding to the triads in figure 11: (a)
$({0.314 - 0.0 = 0.314})$
; (b)
$({0.273 - 0.0 = 0.273})$
; (c)
$({0.269 - 0.0 = 0.269})$
; (d)
$({0.269 - 0.269 = 0.0})$
. The long-time mean is removed from the data. For panel (e)
$({0.269{ -} 0.269 = 0.0})$
the long-time mean is included in the data. Only the absolute value of each mode is shown.

To identify the spatially coherent structures associated with the triadic interaction, bispectral modes
${\boldsymbol{\phi }}_{k+l}$
are calculated. These modes can be interpreted as observable physical structures (Schmidt Reference Schmidt2020), and the bispectral modes of
$\textit{St}_{\!j}$
are qualitatively equivalent to the SPOD spatial modes (see figures 10 and 12). Figure 12(a–c) shows the absolute value of several relevant bispectral modes associated with the dominant interactions along the line
$\textit{St}_l=0$
shown in figure 11. The screech modes, denoted as black dots, have been discussed in the context of figure 10 and are associated with a flapping mode (
$\textit{St}=0.314$
and
$0.269$
) and a quasi-helical mode (
$\textit{St}=0.273$
) (Mazharmanesh et al. Reference Mazharmanesh, Nogueira, Weightman and Edgington-Mitchell2025). The difference self-interactions of the three screech tones resemble the
$\textit{St}=0$
SPOD mode when the long-time mean is removed from the data, with most fluctuations occurring in regions further downstream. For compactness, only the difference self-interaction of one screech tone (
$\textit{St}_{\!j}=0.269-0.269=0.0$
) is shown in figure 12(d), as they are identical and qualitatively similar to the SPOD spatial mode shown in figures 10(h) and 10(j). Note that figure 12(e) shows the long-time mean itself, as the long-time mean is included in the data.
The BMD modes themselves provide little new insight, given their similarity to SPOD modes. Of significantly more interest are the interaction maps
$\boldsymbol{\psi }_{k,l}$
, which spatially localise the regions of the flow where triadic interactions are taking place. Figure 13 presents the interaction maps associated with
$0.269-0.0=0.269$
,
$0.273-0.0=0.273$
and
$0.314-0.0=0.314$
, which can be interpreted as the regions where interaction between
$\textit{St}_{\textit{screech}}$
and the mean flow occurs. Despite the three screech tones having relatively similar mode shapes, their interaction maps are significantly different: the map associated with
$\textit{St}=0.273$
peaks relatively far downstream at
$x/D_e \approx 8$
, where the interaction peaks for
$\textit{St}=0.269$
and
$\textit{St}=0.314$
occur closer to
$x/D_e \approx 6$
and
$x/D_e \approx 5$
, respectively. These regions of activity are consistent with the position of the peak wavepacket amplitude observed in figure 10. Note that, since the jet is invariant under reflection about the
$y=0$
axis, a BMD workflow was employed that leverages the jet’s geometric symmetry to reduce memory requirements and accelerate statistical convergence.
Interaction map of the elliptical jet operating at
$\mathrm{NPR}=3.60$
, corresponding to the triads in figure 11: (a)
$({0.314{ -} 0.0 = 0.314})$
; (b)
$({0.273 - 0.0 = 0.273})$
; (c)
$({0.269 - 0.0 = 0.269})$
. The long-time mean is included in the data. Only the absolute value of each mode is shown.

5. The two-way relationship between screech and the mean flow
In the previous section, BMD was used to identify the triadic interactions that underpin screech. In the present section, we address the second research objective of this study by further exploring the screech-loop closure mechanism. This is done by explicitly linking the wavenumber of the quasi-periodic flow responsible for closing the screech loop to a localised region in the spatial domain. Screech can only occur when the wavenumbers from the triadic interaction between the KH wavepacket and the spatially varying shock structures align with the wavenumbers at which the flow can support the upstream-travelling wave known as the GJM (Nogueira et al. Reference Nogueira, Cavalieri, Martini, Towne, Jordan and Edgington-Mitchell2024)
It has recently become apparent that the GJM need only be supported over a limited streamwise domain where the triadic interaction is taking place (Ivelja et al. Reference Ivelja, Edgington-Mitchell, Maigler and Nogueira2024); identifying this region has so far generally been done via qualitative arguments regarding the streamwise peak of the wavepacket-amplitude envelope. In the physical domain, different downstream regions of the flow are associated with different shock-cell spacings; in the spectral domain, the streamwise variation in shock structure appears as multiple peaks in the wavenumber spectrum (Nogueira et al. Reference Nogueira, Jaunet, Mancinelli, Jordan and Edgington-Mitchell2022a
). The KH wavepacket can interact not only with the primary wavenumber peak
$(k_{s_{1}})$
but also with the sub-optimal peaks
$(k_{s_{2,3}})$
; in the physical domain, this is equivalent to the triadic interaction occurring at different distances downstream from the nozzle. Efforts to either model or describe screech to date using (5.1) have not considered that
$k_{\textit{GJM}}$
and
$k_s$
are intrinsically linked; the relevant peak of
$k_s$
is determined by where the triadic interactions take place, and the streamwise location also determines the range of frequencies at which the GJM is supported. In this section, we will demonstrate a method by which the appropriate
$k_s$
can be extracted, if time-resolved data are available. In work to date, the peaks
$k_{s_{n}}$
have been obtained from temporally averaged particle image velocimetry (Edgington-Mitchell et al. Reference Edgington-Mitchell, Wang, Nogueira, Schmidt, Jaunet, Duke, Jordan and Towne2021), schlieren (Edgington-Mitchell et al. Reference Edgington-Mitchell, Li, Liu, He, Wong, Mackenzie and Nogueira2022) or numerical data (Wang et al. Reference Wang, Li, Sheng, He, Hao and Zhang2025). Specifically, an axial Fourier transform is performed along the centreline of the jet to extract the wavenumber spectrum, as illustrated in figure 14. In this method,
$k_{s}$
is measured after the screech has already reached its limit-cycle state; and the flow has already been altered by the screech. Additionally, this representation includes the entire flow, with no a priori indication as to which regions might be relevant for the screech process.
(a) Temporally averaged schlieren images for the elliptical jet operating at
$\mathrm{NPR}=3.60$
. (b) Results of axial Fourier transforms performed at the centreline of the time-averaged shock structures image at
$\mathrm{NPR}=3.60$
. Red dots represent
$k_{s_1}$
and
$k_{s_2}$
.

We now propose a novel approach to determine the relevant wavenumbers of the shock-cell structure. The BMD interaction maps for
$({\textit{St}_{\textit{screech}}-0=\textit{St}_{\textit{screech}}})$
indicate the spatial locations at which triadic interactions between the screech tones and the mean flow are occurring. As the mean flow in this case includes the shock cells, rather than extracting a spectral representation of the time-averaged shock cells and attempting a post hoc identification of the relevant peak, we can directly calculate the
$k_s$
relevant to each screech tone by identifying the dominant wavenumber in the interaction map. An example of this method is illustrated in figure 15 for
$\mathrm{NPR}=3.6$
. An axial Fourier transform is applied to the absolute value of the interaction maps to obtain a wavenumber spectra, from which the wavenumber associated with only the shock cells most relevant to screech is extracted, denoted
$k_S$
. In this approach, the interaction maps are integrated along the y-axis to minimise the sensitivity of the spectral peak to the specific radial position. The peak near
$k_xD_e \approx 0$
with a value close to unity in figure 15(a–c) arises from the wavepacket envelope of the interaction map. Specifically, the structure grows in the initial region of the jet, reaches a saturation point, and subsequently decays further downstream. As a result, the signal to which the Fourier transforms is applied is not centred around zero, leading to a peak near
$k_xD_e \approx 0$
.
Wavenumber spectra associated with the interaction maps of the screech modes for the elliptical jet operating at
$\mathrm{NPR}=3.60$
: (a)
$({0.314 - 0.0 = 0.314})$
; (b)
$({0.273 - 0.0 = 0.273})$
; (c)
$({0.269 - 0.0 = 0.269})$
. The interaction maps are integrated along the
$y$
-axis to minimise the sensitivity of peak location in the interaction map, followed by an axial Fourier transform.

If the BMD interaction maps have indeed correctly identified the shock cells involved in each screech resonance loop, then setting
$k_s = k_S$
in (5.1) should result in the relation being satisfied. To test whether this is true, it is necessary to determine the wavenumbers of the KH wavepacket
$(k_{\textit{kh}})$
and the GJM
$(k_{\textit{GJM}})$
. These wavenumbers are educed through the application of an axial Fourier transform to the SPOD spatial modes associated with
$\textit{St}=0.269$
,
$\textit{St}=0.273$
and
$\textit{St}=0.314$
. In the spectral domain, wavenumbers
$k_xD_e\gt 0$
correspond to waves with a positive phase velocity, representing downstream-propagating waves, whereas wavenumbers
$k_xD_e\lt 0$
correspond to waves with a negative phase velocity, usually indicative of upstream-propagating waves. Figure 16(a–c) presents the wavenumber spectra as a function of radial position associated with screech modes
$\textit{St}=0.269$
,
$\textit{St}=0.273$
and
$\textit{St}=0.314$
at
$\mathrm{NPR}=3.6$
. Overlaid on these images are lines indicating the wavenumber of the KH wavepacket and the GJM, calculated at the jet lipline (
$y/D_e=\pm 0.35$
). The corresponding values,
$k_{\textit{kh}}$
and
$k_{\textit{GJM}}$
, are highlighted by the red and green dashed lines, respectively. Now that all components of (5.1) are known, the next step is to evaluate whether it is satisfied. The difference wavenumber
$(k_{\textit{kh}}-k_{\textit{GJM}})$
is compared with
$k_{s_{1}}$
,
$k_{s_{2}}$
and
$k_S$
for the dominant screech mode at
$\textit{St}=0.269$
in figure 16(d). The results show that the BMD-derived
$k_S$
predict screech more accurately than
$k_{s_{1}}$
. It is also worth noting that the wavenumber peaks
$(k_{s_{1,2}})$
are estimated from temporally averaged schlieren images for
$t/T \leqslant (t/T)_{\textit{switch}}$
. Using the long-time mean instead results in changes of less than 1
$\,\%$
in
$k{s_{1,2}}$
.
Wavenumber spectra associated with SPOD mode for the elliptical jet operating at
$\mathrm{NPR}=3.60$
with (a)
$\textit{St}=0.314$
, (b)
$\textit{St}=0.273$
, (c)
$\textit{St}=0.269$
. The red and green dashed lines mark
$k_{\textit{kh}}$
and
$k_{GJM}$
, respectively. Horizontal orange dotted lines indicate the jet lipline at
$y/D=\pm 0.35$
. (d) Comparison of wavenumbers associated with the (
$k_{\textit{kh}}-k_{\textit{GJM}}$
) and the primary (
$k_{s_1}$
) and sub-optimal (
$k_{s_2}$
) peaks of the shock-cell structure and BMD-derived
$k_S$
for the dominant screech mode at
$\textit{St}=0.269$
.

Frequency spectra as a function of NPR for the elliptical jet with
$\mathrm{AR} = 2.0$
, viewed in the minor-axis plane. The dominant screech tone at each
$\mathrm{NPR}$
is highlighted using white markers. White crosses indicate the lower-frequency screech modes, while white asterisks represent the higher-frequency screech modes.

To assess how accurately the newly identified
$k_S$
values satisfy the triadic-interaction model, as compared with the conventional approach using temporally averaged schlieren images
$(k_{s_{1}}, k_{s_{2}})$
, we conduct a parametric study with
$\mathrm{NPR}$
varying from 2.30 to 4.90. The dataset used in this part of the study contains
$n_t = 50{\,}000$
snapshots. The SPOD-eigenvalue spectra for the leading mode for the elliptical jet across all considered
$\mathrm{NPR}$
values are presented as a contour in figure 17. For
$2.3 \leqslant \mathrm{NPR} \leqslant 3.3$
, the elliptical jet exhibits a single flapping mode, whereas for
$3.4 \leqslant \mathrm{NPR} \leqslant 4.9$
, a multimodal behaviour with highly intermittent screech tones occurs, similar to the condition observed at
$\mathrm{NPR}=3.6$
for
$t/T \geqslant (t/T)_{\textit{switch}}$
in figure 8(a). The dominant screech tone at each
$\mathrm{NPR}$
, as identified from the SPOD eigenvalue spectra, is overlaid in figure 17 using white markers. A comparison of the SPOD eigenvalues associated with screech (
$\mathrm{\lambda _s}$
) at each
$\mathrm{NPR}$
for the elliptical jet in the multimodal region shows that the higher-frequency screech mode (flapping mode) is dominant for
$3.4\leqslant \mathrm{NPR}\leqslant 4.0$
, whereas the lower-frequency screech mode (quasi-helical mode) dominates at
$4.1\leqslant \mathrm{NPR}\leqslant 4.9$
.
To extract the value of the sum of wavenumbers
$k_{\textit{kh}}$
and
$|k_{\textit{GJM}}|$
, we in fact obtain the value from the standing wave wavelength
$k_{sw}$
from the SPOD modes, according to the relation
$k_{sw}=k_{\textit{kh}}-k_{\textit{GJM}}$
established in Edgington-Mitchell et al. (Reference Edgington-Mitchell, Li, Liu, He, Wong, Mackenzie and Nogueira2022). This approach reduces scatter in the data, and permits a rapid calculation over the full range of
$\mathrm{NPR}$
values. Similar to the analysis presented in figure 16 (d), the standing wave wavenumbers
$(k_{sw}=k_{\textit{kh}}-k_{\textit{GJM}})$
is compared with
$k_{s_{1}}$
,
$k_{s_{2}}$
and
$k_S$
. For the lower-frequency screech modes, indicated by white crosses in figure 17,
$(k_{\textit{kh}}-k_{\textit{GJM}})$
is represented by black crosses (
$\times$
) in figure 18. For the higher-frequency screech modes, marked by white asterisks (*) in figure 17,
$(k_{\textit{kh}}-k_{\textit{GJM}})$
is shown as black asterisks (*) in figure 18. Additionally,
$k_{s_{1}}$
,
$k_{s_{2}}$
and
$k_S$
are represented by red circles, blue squares and green triangles, respectively. The results reveal that when a single screech loop is active
$(2.3 \leqslant \mathrm{NPR} \leqslant 3.3)$
, screech is generated by the interaction between the KH wavepacket and the primary wavenumber peak
$(k_{s_{1}})$
(see figure 18). In the multimodal region
$(3.4 \leqslant \mathrm{NPR} \leqslant 4.9)$
, both screech loops are closed by waves originating from the sub-optimal peak
$(k_{s_{2}})$
. Noted, while
$k_{s_{2}}$
closes multiple screech modes, the spatial structures of the resonant modes differ (Mazharmanesh et al. Reference Mazharmanesh, Nogueira, Weightman and Edgington-Mitchell2025). A comparison of the wavenumbers obtained from BMD with those of the standing waves demonstrates that the
$k_S$
values from BMD predict screech more accurately than
$k_{s_{1}}$
and
$k_{s_{2}}$
in both the single-mode and multimodal regions. Specifically, in the multimodal region, figure 18 shows that
$k_S$
aligns equally well or even better with
$k_{\textit{kh}}-k_{\textit{GJM}}$
than
$k_{s_2}$
for
$4.1 \leqslant \mathrm{NPR} \leqslant 4.9$
, where the lower-frequency screech mode dominates. Similarly, for
$3.4 \leqslant \mathrm{NPR} \leqslant 4.0$
, where the higher-frequency screech mode dominates,
$k_S$
matches
$k_{\textit{kh}}-k_{\textit{GJM}}$
more closely than
$k_{s_2}$
. Overall, these results indicate that
$k_S$
values obtained from BMD provide a more precise determination of shock-cell spacing compared with those derived from temporally averaged schlieren images.
Elliptical jet. Comparison of the wavenumbers associated with
$(k_{sw}=k_{\textit{kh}} - k_{\textit{GJM}})$
, and the primary
$(k_{s_1})$
and sub-optimal
$(k_{s_2})$
peaks of the shock-cell structure, along with the BMD-derived
$k_S$
. Black crosses indicate the lower-frequency screech modes
$(\textit{St}_L)$
, while black asterisks represent the higher-frequency screech modes
$(\textit{St}_H)$
, as identified in figure 17.

Having used the wavenumber spectra to demonstrate that
$k_S$
correctly identifies the region of quasi-periodicity involved in screech, we can now return to the spatial domain to localise these regions. The interaction maps from which
$k_S$
were calculated can also be directly compared with the time-averaged flow structures; we overlay the interaction maps of
$({\textit{St}_{\textit{screech}}-0=\textit{St}_{\textit{screech}}})$
onto temporally averaged schlieren images as demonstrated for
$\mathrm{NPR} = 3.6$
in figure 19(a–c). The interaction maps shown in figure 13(d) are integrated along the
$y$
-axis and then superimposed on the temporally averaged schlieren images in figure 19(d). For
$\mathrm{NPR} = 3.6$
, it is evident that each of the screech modes is associated with a different streamwise region of the flow; while the interactions are spread across extended regions, each peaks at a different location. The screech mode associated with
$(\textit{St}_{\!j}=0.314)$
peaks around the fifth shock cell of the jet, but with measurable interaction occurring over a span of eight shock cells, with a roughly symmetrical envelope (see yellow line in figure 19
d). The interactions for
$(\textit{St}_{\!j}=0.273)$
peak around the seventh shock cell, with a much more skewed envelope; a slow increase from the first shock cell to the seventh, then a rapid decay (see cyan line in figure 19
d). The interactions associated with the highest-amplitude screech mode at
$(\textit{St}_{\!j}=0.269)$
peak at the fourth shock cell, though with non-trivial amplitudes from the second to the sixth; outside this region the amplitude of interaction decays more rapidly than the other two tones (see red line in figure 19
d).
Interaction maps overlaid on the short-time mean shock structures at
$\mathrm{NPR}=3.60$
for the triad interactions: (a)
$({0.314-0.0=0.314})$
; (b)
$({0.273-0.0=0.273})$
; (c)
$({0.269-0.0=0.269})$
. (d) Integrated interaction maps along the
$y$
-axis overlaid on the long-time shock structures for the triad interactions
$({0.314-0.0=0.314})$
(yellow solid line),
$({0.273-0.0=0.273})$
(cyan solid line) and
$({0.269-0.0=0.269})$
(red solid line). For visualisation purposes, the
$y$
-axis of the interaction maps is not scaled consistently with that of the mean flow.

6. Conclusions
This work aims, firstly, to explore the manner in which screech alters the mean flow, and secondly, to identify the spatial region associated with the wavenumber of the quasi-periodic flow that reinforces the screech feedback loop. By applying SPOD to high-speed schlieren images, we identify distinct mean-flow distortions. Frequency–time analysis reveals that intermittency, whether arising from mutual exclusivity or from the intrinsic intermittency of individual screech modes, leads to the emergence of multiple short-time ‘mean-flow’ states, none of which align with the long-time average. Each of these short-time mean flows is associated with a specific screech mode. Analysis of the JPDFs and the Pearson correlation coefficients derived from the SPOD expansion coefficients reveals a strong correlation between screech modes and mean-flow distortion, with each screech mode inducing a unique modification of the mean flow.
Bispectral mode decomposition of schlieren images is employed to investigate the mechanisms responsible for generating mean-flow distortion. The BMD analysis reveals that all dominant triad interactions contributing to mean-flow distortion arise directly from the difference self-interaction of screech modes.
Interaction maps derived from the BMD technique are utilised to investigate the regions of the flow where triadic interaction between the screech mode and the shock occurs. Recent studies have emphasised that a comprehensive understanding of the screech phenomenon requires analysing the triadic interaction between the KH wavepacket and the spatially varying shock structures within the flow. Building upon this, the present study introduced a novel method to determine the peak wavenumbers of the shock-cell structure. This approach enabled the association of suboptimal peaks, arising from streamwise variations in the shock pattern, with specific regions of the jet, allowing
$k_S$
to be determined independently for each screech loop. By applying an axial Fourier transform to the interaction maps of the screech modes, the corresponding wavenumber spectra were obtained, and the peak in each spectrum was identified as
$k_S$
. Overall, our results indicate that the
$k_S$
values obtained from BMD predict screech either equally well or more accurately than
$k_{s_{1,2}}$
, which is derived from the overall mean flow, particularly for the dominant screech mode in the multimodal region. Finally, we establish a connection between the spectral results and the spatial domain, enabling the identification of specific shock cells responsible for closing each screech loop.
Supplementary material
Computational Notebook files are available as supplementary material at https://doi.org/10.1017/jfm.2026.11796 and online at https://www.cambridge.org/S0022112026117960/JFM-Notebooks.
Funding
This work was supported by the Australian Research Council through the Discovery Project scheme: DP220103873. O.T.S. and B.Y. gratefully acknowledge support from Office of Naval Research award N00014-23-1-2457, under the supervision of Dr S. Martens.
Declaration of interests
The authors report no conflict of interest.
Appendix A. Comparison of SPOD and BMD results for the full dataset and subsets
In §§ 3.2 and 4, we stated that the SPOD and BMD results obtained from the full dataset are consistent with those obtained from the pre- and post-switch subsets for the elliptical jet operating at
$\mathrm{NPR}=3.60$
. Therefore, SPOD and BMD remain applicable in this specific condition despite the flow being non-statistically stationary. In this section, we provide evidence supporting this argument. First, we compare the SPOD spatial modes associated with screech tones at
$\textit{St}=0.269$
,
$\textit{St}=0.273$
and
$\textit{St}=0.314$
, obtained using the full dataset and the subsets. To quantitatively assess the similarity between the spatial structures, a normalised inner product,
$\eta _{i,k}$
, is calculated between each mode
$\xi _{i,k}$
obtained from subset
$i$
and the corresponding mode
$\xi _k$
obtained from the complete dataset
where
$\langle \boldsymbol{\cdot }, \boldsymbol{\cdot }\rangle$
denotes the inner product and the double vertical bars denote the Frobenius norm. The index
$i=(1,2)$
represents the two subsets: pre-switch
$(t/T \leqslant (t/T)\textit{switch})$
and post-switch
$(t/T \geqslant (t/T)\textit{switch})$
, while
$k$
denotes the SPOD mode number. The parameter
$\eta$
ranges from
$-1$
to
$1$
. The screech mode is considered highly similar when the correlation coefficient
$\eta$
approaches unity, indicating that the SPOD spatial mode obtained from the subset closely matches that obtained from the full dataset. Conversely,
$\eta =-1$
represents anti-correlated structures, which are not observed in the present case, while values of
$\eta$
close to zero indicate uncorrelated screech modes.
Three screech modes (
$\textit{St}=0.314$
,
$\textit{St}=0.273$
and
$\textit{St}=0.269$
) for the elliptical jet at
$\mathrm{NPR}=3.60$
: (a,c,e) using either the first subset (
$t/T \leqslant (t/T)_{\textit{switch}}$
) or the second subset (
$t/T \geqslant (t/T)_{\textit{switch}}$
); (b,d, f) using full dataset (
$0\leqslant t/T \leqslant 3.05 \times 10^4$
). Only the absolute value of each mode is shown.

Figure 20. Long description
Panel A: A heat map showing screech mode for the elliptical jet using the first subset. The x-axis is labeled x/D_e and the y-axis is labeled y/D_e. The color intensity represents the magnitude of the screech mode, with higher values shown in lighter colors. The range of values is from 0 to approximately 12 on the x-axis and -2 to 2 on the y-axis. The heat map shows distinct clusters of higher intensity regions. Panel B: A heat map showing screech mode for the elliptical jet using the full dataset. The x-axis is labeled x/D_e and the y-axis is labeled y/D_e. The color intensity represents the magnitude of the screech mode, with higher values shown in lighter colors. The range of values is from 0 to approximately 12 on the x-axis and -2 to 2 on the y-axis. The heat map shows distinct clusters of higher intensity regions. Panel C: A heat map showing screech mode for the elliptical jet using the second subset. The x-axis is labeled x/D_e and the y-axis is labeled y/D_e. The color intensity represents the magnitude of the screech mode, with higher values shown in lighter colors. The range of values is from 0 to approximately 12 on the x-axis and -2 to 2 on the y-axis. The heat map shows distinct clusters of higher intensity regions. Panel D: A heat map showing screech mode for the elliptical jet using the full dataset. The x-axis is labeled x/D_e and the y-axis is labeled y/D_e. The color intensity represents the magnitude of the screech mode, with higher values shown in lighter colors. The range of values is from 0 to approximately 12 on the x-axis and -2 to 2 on the y-axis. The heat map shows distinct clusters of higher intensity regions. Panel E: A heat map showing screech mode for the elliptical jet using the first subset. The x-axis is labeled x/D_e and the y-axis is labeled y/D_e. The color intensity represents the magnitude of the screech mode, with higher values shown in lighter colors. The range of values is from 0 to approximately 12 on the x-axis and -2 to 2 on the y-axis. The heat map shows distinct clusters of higher intensity regions. Panel F: A heat map showing screech mode for the elliptical jet using the full dataset. The x-axis is labeled x/D_e and the y-axis is labeled y/D_e. The color intensity represents the magnitude of the screech mode, with higher values shown in lighter colors. The range of values is from 0 to approximately 12 on the x-axis and -2 to 2 on the y-axis. The heat map shows distinct clusters of higher intensity regions.
Figure 20(a–f) shows the SPOD spatial modes and the corresponding
$\eta$
values for the leading SPOD mode (
$k=1$
) of the elliptical jet operating at
$\mathrm{NPR}=3.6$
. It is evident from this figure that the subsets and the full dataset yield nearly identical SPOD spatial modes, with correlation coefficients greater than
$0.99$
. This confirms that, under the present conditions, SPOD applied to the full dataset produces consistent results.
We also examined the sensitivity of the SPOD expansion coefficients to the use of subsets versus the full dataset. The SPOD expansion coefficient associated with
$\textit{St}=0.269$
was computed using the first
$n_t=100{\,}000$
snapshots prior to
$(t/T)\textit{switch}$
, while the coefficient associated with
$\textit{St}=0.314$
was computed using the last
$n_t=100{\,}000$
snapshots after
$(t/T)\textit{switch}$
. These were then compared with the coefficients obtained from the full dataset. Figure 21 clearly shows that the SPOD expansion coefficient at
$\textit{St}=0.269$
is identical for the full dataset and the pre-switch subset, while the coefficient at
$\textit{St}=0.314$
is identical for the full dataset and the post-switch subset.
Normalised SPOD expansion coefficients corresponding to the two screech modes. Solid lines represent results obtained using the full dataset, while dashed lines correspond to the first
$n_t=100\,000$
snapshots prior to
$(t/T)_{\textit{switch}}$
for
$\textit{St} = 0.269$
and the last
$n_t=100\,000$
snapshots after
$(t/T)_{\textit{switch}}$
for
$\textit{St} = 0.314$
.

Bispectral modes of three screech modes (
${\textit{St}_{\!j}}=0.314$
,
${\textit{St}_{\!j}}=0.273$
and
${\textit{St}_{\!j}}=0.269$
) for the elliptical jet at
$\mathrm{NPR}=3.60$
: (a,c,e) using either the first subset (
$t/T \leqslant (t/T)_{\textit{switch}}$
) or the second subset (
$t/T \geqslant (t/T)_{\textit{switch}}$
); (b,d, f) using full dataset (
$0\leqslant t/T \leqslant 3.05 \times 10^4$
). The long-time mean is included in data. Only the absolute value of each mode is shown.

Similar to the SPOD spatial modes, we use (A1) to evaluate the similarity of the BMD results obtained from the full dataset and the subsets. Figure 22 shows the bispectral modes of the three screech tones obtained using the different datasets. Here, the normalised inner product is computed using
$\boldsymbol{\phi }$
, which represents the flow structures arising from nonlinear triadic interactions. It is evident from figure 22 that the subsets and the full dataset yield nearly identical bispectral modes, with correlation coefficients greater than
$0.99$
.
We also examine the similarity of the interaction maps obtained from the full dataset and the subsets. In this case, the normalised inner product is computed using
$\boldsymbol{\psi }$
. Figure 23 demonstrates that the subsets and the full dataset produce nearly identical interaction maps, again with correlation coefficients greater than
$0.99$
. This figure, together with figure 22, confirms that under the present conditions BMD applied to the full dataset yields consistent results.
Interaction map of three screech modes (
${\textit{St}_{\!j}}=0.314$
,
${\textit{St}_{\!j}}=0.273$
and
${\textit{St}_{\!j}}=0.269$
) for the elliptical jet at
$\mathrm{NPR}=3.60$
: (a,c,e) using either the first subset (
$t/T \leqslant (t/T)_{\textit{switch}}$
) or the second subset (
$t/T \geqslant (t/T)_{\textit{switch}}$
); (b,d, f) using full dataset (
$0\leqslant t/T \leqslant 3.05 \times 10^4$
). The long-time mean is included in data. Only the absolute value of each mode is shown.

Figure 23. Long description
The image contains six panels of heatmap graphs showing the interaction map of three screech modes for an elliptical jet at a specific condition. Each panel represents different subsets of data and conditions. Panel A, C, and E use either the first subset or the second subset of data, while Panel B, D, and F use the full dataset. The x-axis is labeled x/D_e and the y-axis is labeled y/D_e. Each panel shows the absolute value of each mode. Panel A shows data for t/T greater than or equal to (t/T)_switch with an efficiency of 0.9892 and a Strouhal number of 0.314. Panel B shows data for 0 less than or equal to t/T less than or equal to 3.05 times 10^4. Panel C shows data for t/T greater than or equal to 3.05 times 10^4 with an efficiency of 0.9902 and a Strouhal number of 0.273. Panel D shows data for 0 less than or equal to t/T less than or equal to 3.05 times 10^4. Panel E shows data for t/T less than or equal to (t/T)_switch with an efficiency of 0.9910 and a Strouhal number of 0.269. Panel F shows data for 0 less than or equal to t/T less than or equal to 3.05 times 10^4. The long-time mean is included in the data.
(a) The SPOD-based frequency–time diagrams derived from the leading mode at each frequency for the axisymmetric jet, operating at
$\mathrm{NPR}=2.19$
where the A1 mode dominates. (b) Normalised SPOD expansion coefficients corresponding to the screech mode (
$\textit{St} = 0.607$
) and the mean-flow distortion mode (
$\textit{St} = 0$
), shown for mode 1. (c) The JPDFs involving screech mode and the zero-frequency bin at
$\mathrm{NPR}=2.19$
. (d) Real part of SPOD spatial mode associated with the screech mode (
$\textit{St}=0.607$
) and the mean-flow distortion mode (
$\textit{St}=0$
) for the axisymmetric jet at
$\mathrm{NPR}=2.19$
.

Figure 24. Long description
Panel A: A heatmap displays the SPOD-based frequency-time diagrams for an axisymmetric jet, with the vertical axis labeled St and the horizontal axis labeled t in seconds. The color scale indicates the magnitude of the leading mode at each frequency. Panel B: A line graph shows the normalized SPOD expansion coefficients over time for two modes, with the vertical axis labeled |a/a_max| and the horizontal axis labeled t in seconds. Two lines, one red and one black, represent different modes. Panel C: A contour plot illustrates the joint probability density function (JPDF) involving the screech mode and the zero-frequency bin, with the vertical axis labeled |a1(St1 = 0.607)| and the horizontal axis labeled Re(a1(St0)). The color scale indicates the density values. Panel D: Two heatmaps show the real part of the SPOD spatial mode for the screech mode and the mean-flow distortion mode. The vertical axis is labeled y/De and the horizontal axis is labeled x/De. The color scale indicates the magnitude of the spatial mode.
Appendix B. Mean-flow distortion in the absence of intermittency
The presence of intermittency can lead to mean-flow distortion; however, the converse is not necessarily true, mean-flow distortion may arise from other sources. As shown in figures 4 and 10, mean-flow distortion appears as shock-cell-like structures in the SPOD spatial mode at
$\textit{St}=0$
. However, the presence of such structures at the zero-frequency bin does not necessarily indicate that the distortion is caused by intermittent screech modes. Figure 24 illustrates this with an example of an axisymmetric jet operating at
$\mathrm{NPR}=2.19$
, where the A1 mode is dominant. Frequency–time analysis reveals a strong and stable screech tone over time. Nonetheless, a slow and gradual variation is observed in the screech tone at
$\textit{St}=0.607$
(see figure 24
a). This variation is accompanied by a corresponding gradual change in the SPOD expansion coefficient at
$\textit{St}=0$
, as shown in figure 24(b). This long-time-scale variation is attributed to fluctuations in the air pressure supplied during the experiment. Therefore, although the JPDF between the screech mode and the zero-frequency bin shows a strong correlation, this relationship does not necessarily stem from screech intermittency or self-interaction. Rather, it may simply reflect changes in the experimental pressure supply. This finding highlights the importance of carefully examining frequency–time analysis results alongside SPOD spatial modes in order to correctly identify the origin of mean-flow distortion.


De
l/De
NPR=2.20
NPR=3.60
NPR=2.19
NPR=2.20
St1,2
NPR=2.25
St1
St2
St=0
NPR=2.19
NPR=2.25
T
St=0.665
NPR=2.20
T
St=0.694
ai(St,t)
St=0.606
St=0.694
St=0
NPR=2.20
NPR=2.20
ρ(|a1(St1)|,|a1(St2)|)=−0.82
NPR=2.20
ρ(Re(a1(St0)),|a1(St1)|)=0.96
ρ(Re(a1(St0)),|a1(St2)|)=−0.82
St0
Re(⋅)
a1(⋅)
NPR=2.20
I¯
St1=0.694
St2=0.606
amax
amin
a1(St=0)
a1(St=0)>0
a1(St=0)<0
NPR=3.60
St=0.314
St=0.269
St=0
T
St=0.269
NPR=3.60
ρ(|a1(St1)|,|a1(St2)|)=−0.92
NPR=3.60
ρ(Re(a1(St0)),|a1(St1)|)=0.97
ρ(Re(a1(St0)),|a1(St2)|)=−0.91
NPR=3.60
ρ(Re(a2(St0)),|a1(St2)|)=0.83
St0
Re(⋅)
a1,2(⋅)
St=0.314
St=0.273
St=0.269
St=0
NPR=3.60
NPR=3.60
β
Stj=Stk±Stl
NPR=3.6
(0.314−0.0=0.314)
(0.273−0.0=0.273)
(0.269−0.0=0.269)
(0.269−0.269=0.0)
(0.269−0.269=0.0)
NPR=3.60
(0.314−0.0=0.314)
(0.273−0.0=0.273)
(0.269−0.0=0.269)
NPR=3.60
NPR=3.60
ks1
ks2
NPR=3.60
(0.314−0.0=0.314)
(0.273−0.0=0.273)
(0.269−0.0=0.269)
y
NPR=3.60
St=0.314
St=0.273
St=0.269
kkh
kGJM
y/D=±0.35
kkh−kGJM
ks1
ks2
kS
St=0.269
AR=2.0
NPR
(ksw=kkh−kGJM)
(ks1)
(ks2)
kS
(StL)
(StH)
NPR=3.60
(0.314−0.0=0.314)
(0.273−0.0=0.273)
(0.269−0.0=0.269)
y
(0.314−0.0=0.314)
(0.273−0.0=0.273)
(0.269−0.0=0.269)
y
St=0.314
St=0.273
St=0.269
NPR=3.60
t/T⩽(t/T)switch
t/T⩾(t/T)switch
0⩽t/T⩽3.05×104
nt=100000
(t/T)switch
St=0.269
nt=100000
(t/T)switch
St=0.314
Stj=0.314
Stj=0.273
Stj=0.269
NPR=3.60
t/T⩽(t/T)switch
t/T⩾(t/T)switch
0⩽t/T⩽3.05×104
Stj=0.314
Stj=0.273
Stj=0.269
NPR=3.60
t/T⩽(t/T)switch
t/T⩾(t/T)switch
0⩽t/T⩽3.05×104
NPR=2.19
St=0.607
St=0
NPR=2.19
St=0.607
St=0
NPR=2.19