1. Introduction
Triadic interactions are fundamental to the transfer of energy in turbulent flow. These interactions arise from the quadratic nonlinearity of the convective term in the Navier–Stokes equation. A triad is a combination of three wave vectors that form a triangle, represented as
where
$\boldsymbol{k}_{\!j},\boldsymbol{k}_k,\boldsymbol{k}_l$
denote wavenumber vectors, and
$f_{\!j}, f_k, f_l$
denote frequencies. Triadic interactions and their role in energy transfer have been investigated in homogeneous turbulence (Domaradzki & Rogallo Reference Domaradzki and Rogallo1990; Waleffe Reference Waleffe1992; Moffatt Reference Moffatt2014). These studies categorise triadic interactions as either local or non-local. Local interactions refer to energy transfer between wavelengths of similar scales, while non-local interactions involve energy transfer between length scales that differ. For homogeneous turbulence, Domaradzki & Rogallo (Reference Domaradzki and Rogallo1990) demonstrated that non-local triadic interactions are primarily responsible for energy transfer. Triadic interactions have also been examined in other types of flows, such as oceanographic flows (Phillips Reference Phillips1960; Garrett & Munk Reference Garrett and Munk1975; McComas & Bretherton Reference McComas and Bretherton1977) and transitional flows (Craik Reference Craik1971; Smith & Stewart Reference Smith and Stewart1987; Mankbadi, Wu & Lee Reference Mankbadi, Wu and Lee1993; Wu Reference Wu1995, Reference Wu2019). However, triadic interactions and their associated energy transfer in turbulent jets have yet to be explored in detail. Triadic interactions in oceanographic and transitional flows typically involve waves that all correspond to eigenmodes. In contrast, this work examines interactions among three waves, which may or may not correspond to eigenmodes.
Coherent structures in turbulent jets are the primary sources of aft-angle noise (Jordan & Colonius Reference Jordan and Colonius2013) and have been investigated for many years. Pioneering work by Crow & Champagne (Reference Crow and Champagne1971) identified the presence of large-scale coherent structures in turbulent jets, and subsequent researchers (Crighton & Gaster Reference Crighton and Gaster1976; Michalke Reference Michalke1977; Tam & Morris Reference Tam and Morris1980; Tam & Burton Reference Tam and Burton1984; Crighton & Huerre Reference Crighton and Huerre1990) modelled these coherent structures, with varying degrees of success, as linear disturbances to the turbulent mean flow. Batchelor & Gill (Reference Batchelor and Gill1962), using temporal linear stability analysis for high Reynolds number flows, showed that in the fully developed region, only helical modes with with azimuthal wavenumbers
$m = \pm 1$
are unstable. In more recent work, the resolvent approach, pioneered by McKeon & Sharma (Reference McKeon and Sharma2010) for wall-bounded flows, has clarified the fundamental linear mechanisms leading to amplification of disturbances in turbulent jets (Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018; Pickering et al. Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020). These mechanisms include Kelvin–Helmholtz (KH) wavepackets (Crighton & Gaster Reference Crighton and Gaster1976), Orr wavepackets (Tissot et al. Reference Tissot, Lajús, Cavalieri and Jordan2017a
; Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018), and streaks (Nogueira et al. Reference Nogueira, Cavalieri, Jordan and Jaunet2019; Pickering et al. Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020). The resolvent framework also shows that the optimal nonlinear forcing for the Orr wavepackets and streaks occurs extensively throughout the jet shear layers and fully developed region, whereas the KH wavepackets are convectively amplified from spatially compact regions much closer to the nozzle exit. These interpretations are also consistent with local (spatial) stability analysis, where KH instability is modal, whereas streaks and Orr wavepackets are non-modal. These three wavepackets are all associated with amplification behaviour: streaks and Orr wavepackets exhibit transient growth, while KH wavepackets exhibit long-term exponential growth. Note that the term ‘Orr wavepackets’ is used loosely, as their structures resemble the Orr mechanism in boundary layers.
In jets, especially those with initially turbulent boundary layers, much less is understood about nonlinear mechanisms and the presence of self-sustaining processes involving the dominant amplification mechanisms. In boundary layers, nonlinear interactions are crucial to the self-sustaining near-wall cycle (Waleffe Reference Waleffe1997; Jiménez Reference Jiménez2013; Bae, Lozano-Duran & McKeon Reference Bae, Lozano-Duran and McKeon2021), supplying energy to rolls and extracting energy from streaks (Waleffe & Kim Reference Waleffe and Kim1998; Dessup et al. Reference Dessup, Tuckerman, Wesfreid, Barkley and Willis2018). A few studies have explored how nonlinearities specifically contribute to coherent structures in jets (Sandham & Salgado Reference Sandham and Salgado2008; Wu & Huerre Reference Wu and Huerre2009; Suponitsky, Sandham & Morfey Reference Suponitsky, Sandham and Morfey2010; Wu Reference Wu2019; Zhang & Wu Reference Zhang and Wu2023). For instance, Zhang & Wu (Reference Zhang and Wu2023) demonstrated that streaks are tied to helical KH modes in a weakly nonlinear expansion of a weakly non-parallel jet flow. Nonlinear dynamics of coherent structures have also been explored by acoustic (and other) forcing near the nozzle. These flows exhibit tonal responses indicative of strong nonlinear resonant interactions. Forcing the jet at the azimuthal wavenumbers
$\pm m$
triggers nonlinear interactions that result in mean flow distortion (Cohen & Wygnanski Reference Cohen and Wygnanski1987; Long & Petersen Reference Long and Petersen1992; Corke & Kusek Reference Corke and Kusek1993). In particular, the mutual interaction between helical modes
$\pm m$
at the same frequency generates streamwise rolls, which subsequently produce streaks through the lift-up mechanism (Hall & Sherwin Reference Hall and Sherwin2010). Samimy et al. (Reference Samimy, Kim, Kastner, Adamovich and Utkin2007) found that forcing helical modes significantly shortens the jet potential core, and substantially enhances entrainment. Other studies (Husain & Hussain Reference Husain and Hussain1995; Paschereit, Wygnanski & Fiedler Reference Paschereit, Wygnanski and Fiedler1995; Broze & Hussain Reference Broze and Hussain1996; Shaabani-Ardali, Sipp & Lesshafft Reference Shaabani-Ardali, Sipp and Lesshafft2019; Nekkanti, Colonius & Schmidt Reference Nekkanti, Colonius and Schmidt2025) have examined the response of initially laminar jets. Broze & Hussain (Reference Broze and Hussain1994, Reference Broze and Hussain1996) studied how forcing an initially laminar jet at different frequencies and amplitudes leads to periodic and chaotic patterns as the flow transitions to turbulence.
An outstanding question is whether, in a natural (unforced) initially turbulent jet, there are important nonlinear interactions amongst the dominant coherent structures, and if so, which structures, and at which frequencies and azimuthal mode numbers? Do interactions differ in the shear layer and fully developed regions? To answer these questions, we use bispectral mode decomposition (BMD), a technique recently proposed by Schmidt (Reference Schmidt2020), which extends classical bispectral analysis to multi-dimensional and multivariate data. The classical bispectrum, or its normalised form, the bicoherence, measures quadratic phase coupling at a single spatial point, as detailed by Hasselmann, Munk & MacDonald (Reference Hasselmann, Munk and MacDonald1963), Brillinger (Reference Brillinger1965) and Jenkins (Reference Jenkins1965). The bispectrum has been used to analyse homogeneous turbulence flows (Yeh & Van Atta Reference Yeh and Van Atta1973; Lii, Rosenblatt & Van Atta Reference Lii, Rosenblatt and Van Atta1976; Herring Reference Herring1980) and laminar–turbulent transitions (Corke & Mangano Reference Corke and Mangano1989; Hajj et al. Reference Hajj, Miksad and Powers1992, Reference Hajj, Miksad and Powers1993). The BMD facilitates the extraction of spatial flow structures associated with triadic interactions by maximising the spatially integrated bispectrum. It has previously been applied to various flow configurations, including laminar–turbulent transition on a flat plate (Goparaju & Gaitonde Reference Goparaju and Gaitonde2022), swirling flows (Moczarski et al. Reference Moczarski, Treleaven, Oberleithner, Schmidt, Fischer and Kaiser2022), bluff body wakes (Nekkanti et al. Reference Nekkanti and Schmidt2023a ), forced jets (Maia et al. Reference Maia, Jordan, Heidt, Colonius, Nekkanti and Schmidt2021; Nekkanti et al. Reference Nekkanti, Maia, Jordan, Heidt, Colonius and Schmidt2022, Reference Nekkanti, Schmidt, Maia, Jordan, Heidt and Colonius2023b ), and twin rectangular jets (Yeung & Schmidt Reference Yeung and Schmidt2025).
This paper is organised as follows. In § 2, the numerical methodology for BMD is discussed, and the symmetries of the bispectrum are examined. Results from BMD analysis are presented in § 3. Section 4 details the spatial localisation of energy transfer for different triads. The results are summarised in § 5.
2. Cross-BMD
We employ the cross-BMD, recently developed by Schmidt (Reference Schmidt2020), to identify the spatially coherent structures arising from triadic interactions in a turbulent flow. Here, we provide a brief overview of the method. The reader is referred to Schmidt (Reference Schmidt2020) for further details of the derivation and mathematical properties of the method.
Given a statistically stationary flow field, let
$\boldsymbol{q}_i =\boldsymbol{q}(t_i)$
denote the mean subtracted snapshots, where
$i =1,2,\ldots, n_t$
are
$n_t$
number of snapshots. In this study,
$\boldsymbol{q}$
comprises the three components of the velocity vector, but in theory, BMD can be applied to
$\boldsymbol{q}$
containing any set of variables. For spectral estimation, using Welch’s approach, the dataset is first segmented into
$n_{{blk}}$
overlapping (by 50 %) blocks with
$n_{\textit{fft}}$
snapshots in each block. Each block is then Fourier transformed in time, and all Fourier realisations of the
$l$
th frequency,
$\hat {\boldsymbol{q}}^{(j)}_l$
, are arranged in a matrix:
The auto-bispectral matrix is then computed as
where
$\hat {\boldsymbol{Q}}_{k\circ l}^{H} = \hat {\boldsymbol{Q}}_{k}^{*}\circ \hat {\boldsymbol{Q}}_{l}^{*}$
,
$\circ$
denotes the Hadamard (or element-wise) product,
$({\cdot })^*$
denotes the complex conjugate, and
$\boldsymbol{W}$
is the diagonal matrix containing the spatial quadrature weights. The auto-bispectral density matrix measures the interactions between different frequencies at the same azimuthal wavenumber. To estimate the interactions between the azimuthal wavenumber triad
$[m_1,m_2,m_3]$
, where
$m_1 + m_2=m_3$
, we construct the cross-bispectral matrix
Here,
$\hat {\boldsymbol{Q}}_{k}$
,
$\hat {\boldsymbol{R}}_{l}$
,
$\hat {\boldsymbol{S}}_{k+l}$
comprise all the Fourier realisations of the velocity field at the
$k$
th frequency of the azimuthal wavenumber
$m_1$
, the
$l$
th frequency of the azimuthal wavenumber
$m_2$
, and the
$(k+l)$
th frequency of the azimuthal wavenumber
$m_3$
, respectively. The velocity field obtained from large-eddy simulations (LES) is Fourier transformed first in the azimuthal direction and then in time to construct the matrices
$\hat {\boldsymbol{Q}}_{k}$
,
$\hat {\boldsymbol{R}}_{l}$
and
$\hat {\boldsymbol{S}}_{k+l}$
. Owing to the non-Hermitian nature of the bispectral matrix, the optimal expansion coefficients
$\boldsymbol{a}_1$
are obtained by maximising the absolute value of the Rayleigh quotient of
$\boldsymbol{B}_c$
:
The complex mode bispectrum is then obtained as
Finally, the leading-order bispectral modes and the cross-frequency fields are recovered as
respectively. By construction, the bispectral modes and cross-frequency fields have the same set of expansion coefficients. This explicitly ensures phase coupling between the resonant frequency triad
$(f_k,f_l,f_k + f_l)$
, where
$\hat {\boldsymbol{Q}}_{k} \circ \hat {\boldsymbol{R}}_{l}$
is the cause, and
$\hat {\boldsymbol{S}}_{k+l}$
is the effect. The complex mode bispectrum
$\lambda _1$
measures the intensity of the triadic interaction, and the bispectral mode
${\phi }_{k+l}$
represents the structures that result from the nonlinear triadic interaction. Another useful quantity is the interaction map, defined as
This interaction map identifies the spatial regions of activity of the triadic interaction for the frequency triad
$ (f_k,f_l,f_k+f_l )$
and wavenumber
$ [m_1,m_2,m_3 ]$
. In what follows, we differentiate between the frequency and azimuthal wavenumber triads using the notation
$({\cdot },{\cdot },{\cdot })$
and
$[{\cdot },{\cdot },{\cdot }]$
, respectively.
In §§ 3 and 4, we use the terms local and global analysis. Local analysis refers to evaluating the above equations at each
$(x,r)$
location, followed by radial integration. Global analysis refers to considering the spatial correlation, where the above equations are applied to a vector
$\boldsymbol{q}$
that contains the velocity field at all spatial points in the
$x{-}r$
domain.
2.1. Symmetries of the bispectrum
Consider the pointwise cross-bispectrum for a frequency triad
$(f_1,f_2,f_3)$
and azimuthal wavenumber triad
$[m_1,m_2,m_3]$
as
where
$f_3 = f_1 + f_2$
and
$m_3 = m_1 + m_2$
, and
$\mathbb{E}[{\cdot }]$
denotes the expectation operator. The bispectral modes maximise the spatial integral of the pointwise bispectrum, and share the same symmetries. For simplicity, we suppress the azimuthal wavenumber triad notation, and consider the cross-bispectrum for three different signals
$\boldsymbol{q},\boldsymbol{r},\boldsymbol{s}$
as
There exist different possibilities for
$\boldsymbol{q},\boldsymbol{r},\boldsymbol{s}$
in a turbulent flow that is homogeneous in one or more directions, and their symmetries are discussed below.
2.1.1. Case 1:
$\boldsymbol{q}=\boldsymbol{r}=\boldsymbol{s} \in \mathbb{R}$
The bispectrum is symmetric about the line
$f_1 = f_2$
, as
$B_{\boldsymbol{q}\boldsymbol{q}\boldsymbol{q}}(f_1,f_2) = B_{\boldsymbol{q}\boldsymbol{q}\boldsymbol{q}}(f_2,f_1)$
. For a real signal, the Fourier transform satisfies
$\hat {\boldsymbol{q}}(f)=\hat {\boldsymbol{q}}^{*}(-f)$
, which results in the bispectrum being symmetric about the line
$f_1=-f_2$
:
\begin{align} B_{\boldsymbol{q}\boldsymbol{q}\boldsymbol{q}}(-f_2,-f_1) &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(-f_2) \circ \hat {\boldsymbol{q}}^{*}(-f_1) \circ \hat {\boldsymbol{q}}(-f_3)\right ]\big | \nonumber \\ &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}(f_2) \circ \hat {\boldsymbol{q}}(f_1) \circ \hat {\boldsymbol{q}}^{*}(f_3)\right ] \big | \nonumber \\ & = B_{\boldsymbol{q}\boldsymbol{q}\boldsymbol{q}}^{*}(f_1,f_2). \end{align}
Additionally, it can be shown that
Due to these symmetries, the triangle enclosed by the vertices
$(f_1,f_2)=(0,0)$
,
$(0,f_N)$
,
$(f_N/2,f_N/2)$
for
$f_1 \geqslant 0$
,
$f_2 \geqslant 0$
, suffices to give us the total information. Here,
$f_N$
is the Nyquist frequency. Physically, this case corresponds to the azimuthal wavenumber triad
$[m_1,m_2,m_3]=[0,0,0]$
.
2.1.2. Case 2:
$\boldsymbol{q}=\boldsymbol{r}=\boldsymbol{s} \in \mathbb{C}$
Similar to the case above, the bispectrum is symmetric about the line
$f_1 = f_2$
as
$B_{\boldsymbol{q}\boldsymbol{q}\boldsymbol{q}}(f_1,f_2) = B_{\boldsymbol{q}\boldsymbol{q}\boldsymbol{q}}(f_2,f_1)$
. For a complex-valued signal,
$\hat {\boldsymbol{q}}(f) \neq \hat {\boldsymbol{q}}^{*}(-f)$
, therefore the bispectrum is not symmetric about the line
$f_1=-f_2$
:
\begin{align} B_{\boldsymbol{q}\boldsymbol{q}\boldsymbol{q}}(-f_2,-f_1) & = \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(-f_2) \circ \hat {\boldsymbol{q}}^{*}(-f_1) \circ \hat {\boldsymbol{q}}(-f_3)\right ] \big | \nonumber \\ & \neq \big |\mathbb{E}\left [\hat {\boldsymbol{q}}(f_2) \circ \hat {\boldsymbol{q}}(f_1) \circ \hat {\boldsymbol{q}}^{*}(f_3)\right ] \big | \nonumber \\ & \neq B_{\boldsymbol{q}\boldsymbol{q}\boldsymbol{q}}^*(f_1,f_2). \end{align}
This case is generally not encountered in turbulent flows with different wavenumber triads, as
$m_3 = m_1 + m_2$
. A commonly encountered triad is
$[m,m,2m]$
, i.e.
$\boldsymbol{q} = \boldsymbol{r} \neq \boldsymbol{s} \in \mathbb{C}$
, which is also symmetric only about the line
$f_1 = f_2$
.
2.1.3. Case 3: three different real signals
For three different real signals, i.e.
$\boldsymbol{q} \neq \boldsymbol{r} \neq \boldsymbol{s} \in \mathbb{R}$
, the cross-bispectrum exhibits
$180^\circ$
rotational symmetry. Consequently, the first and third quadrants, as well as the second and fourth, contain the same information, making it sufficient to analyse only one half-plane:
\begin{align} B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{s}}(f_1,f_2) &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(f_1) \circ \hat {\boldsymbol{r}}^*(f_2) \circ \hat {\boldsymbol{s}}(f_3)\right ]\big |, \nonumber \\[3pt] B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{s}}(-f_1,-f_2) &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(-f_1) \circ \hat {\boldsymbol{r}}^*(-f_2) \circ \hat {\boldsymbol{s}}(-f_3)\right ]\big | \nonumber \\ &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}(f_1) \circ \hat {\boldsymbol{r}}(f_2) \circ \hat {\boldsymbol{s}}^*(f_3)\right ]\big | \nonumber \\ &= B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{s}}^*(f_1,f_2). \end{align}
2.1.4. Case 4:
$\boldsymbol{q} = \boldsymbol{r}^* \in \mathbb{C},\ \boldsymbol{s} \in \mathbb{R}$
Now let us consider the case where
$\boldsymbol{r} = \boldsymbol{q}^{*} \in \mathbb{C}$
and
$\boldsymbol{s} \in \mathbb{R}$
. This forms a special triad
$[m,-m,0]$
, which contributes to the mean flow deformation. Here,
$\hat {\boldsymbol{r}}(f) = \hat {\boldsymbol{q}}^{*}(-f)$
, and as a result, the cross-bispectrum is symmetric about the line
$f_1 = -f_2$
:
\begin{align} B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{s}}(f_1,f_2) &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(f_1) \circ \hat {\boldsymbol{r}}^*(f_2) \circ \hat {\boldsymbol{s}}(f_3)\right ]\big |, \nonumber \\[3pt] B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{s}}(-f_2,-f_1) &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(-f_2) \circ \hat {\boldsymbol{r}}^*(-f_1) \circ \hat {\boldsymbol{s}}(-f_3)\right ] \big | \nonumber \\ &= \big |\mathbb{E}\left [\hat {\boldsymbol{r}}(f_2) \circ \hat {\boldsymbol{q}}(f_1) \circ \hat {\boldsymbol{s}}_3^*(f_3)\right ]\big | \nonumber \\ &= B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{s}}^*(f_1,f_2). \end{align}
An additional symmetry arises if
$\boldsymbol{q}$
and
$\boldsymbol{s}$
are jointly strict-sense stationary, in which case the cross-bispectrum is symmetric about
$f_1=f_2$
:
\begin{align} B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{s}}(f_1,f_2) &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(f_1) \circ \hat {\boldsymbol{r}}^*(f_2) \circ \hat {\boldsymbol{s}}(f_3)\right ]\big | \nonumber \\ &= \big | \mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(f_1) \circ \hat {\boldsymbol{q}}(-f_2) \circ \hat {\boldsymbol{s}}(f_3)\right ]\big |, \end{align}
\begin{align} B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{s}}(f_2,f_1) &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(f_2) \circ \hat {\boldsymbol{r}}^*(f_1) \circ \hat {\boldsymbol{s}}(f_3)\right ]\big | \nonumber \\ &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}(-f_1) \circ \hat {\boldsymbol{q}}^*(f_2) \circ \hat {\boldsymbol{s}}(f_3)\right ]\big | \end{align}
The joint strict-sense stationary property ensures that the right-hand sides of (2.17) and (2.18) are equivalent, confirming the equality in (2.19).
2.1.5. Case 5:
$\boldsymbol{q} \in \mathbb{R},\ \boldsymbol{r}=\boldsymbol{s} \in \mathbb{C}$
A
$[0,m,m]$
triad is representative of the case with
$\boldsymbol{q} \in \mathbb{R}$
and
$\boldsymbol{r}= \boldsymbol{s} \in \mathbb{C}$
. The cross-bispectrum exhibits the symmetry
\begin{align} B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{r}}(f_1,f_2) &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(f_1) \circ \hat {\boldsymbol{r}}^*(f_2) \circ \hat {\boldsymbol{r}}(f_3)\right ]\big |, \nonumber \\[3pt] B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{r}}(-f_1,f_3) &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}^{*}(-f_1) \circ \hat {\boldsymbol{r}}^*(f_3) \circ \hat {\boldsymbol{r}}(f_2)\right ]\big | \nonumber \\ &= \big |\mathbb{E}\left [\hat {\boldsymbol{q}}(f_1) \circ \hat {\boldsymbol{r}}(f_2) \circ \hat {\boldsymbol{r}}^*(f_3)\right ]\big | \nonumber \\ &= B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{r}}^*(f_1,f_2) . \end{align}
This symmetry can be understood as translating the ordinate by the value of the abscissa and reflecting it across the
$y$
-axis.
2.1.6. Case 6: three different complex signals
For three different complex signals,
$\boldsymbol{q}_1 \neq \boldsymbol{q}_2 \neq \boldsymbol{q}_3 \neq \boldsymbol{q} \in \mathbb{C}$
, the cross-bispectrum does not exhibit any symmetries.
2.1.7. Additional symmetries
There are additional symmetries when considering the relationship between bispectra with different azimuthal triads. The most obvious one is
which indicates that the bispectrum of an azimuthal wavenumber triad
$[m_1,m_2,m_3]$
can be derived from that of
$[m_2,m_1,m_3]$
by reflecting it across the line
$f_1=f_2$
. Second,
implying that the bispectrum of the triad
$[-m_1,-m_2,-m_3]$
is a
$180^\circ$
rotated version of the bispectrum corresponding to the triad
$[m_1,m_2,m_3]$
. Additionally,
\begin{align} B_{\boldsymbol{q}\boldsymbol{r}\boldsymbol{s}}(f_1,f_2) &= B_{\boldsymbol{s}\boldsymbol{r^*}\boldsymbol{q}}^*(f_1+f_2,-f_2) \nonumber \\ &= B_{\boldsymbol{s}\boldsymbol{q^*}\boldsymbol{r}}^*(f_1+f_2,-f_1), \end{align}
revealing that the bispectra of the triads
$[m_3,-m_2,m_1]$
and
$[m_3,-m_1,m_2]$
can be deduced from the bispectrum of
$[m_1,m_2,m_3]$
. In summary, if the bispectrum of the azimuthal wavenumber triad
$[1,2,3]$
is known, then bispectra of the triads
$[2,1,3]$
,
$[-1,-2,-3]$
,
$[-2,-1,-3]$
,
$[3,-1,2]$
,
$[3,-2,1]$
,
$[1,-3,2]$
and
$[2,-3,-1]$
can be readily determined.
3. Globally dominant triads and their flow structures
3.1. The LES data
We utilise LES data for an isothermal subsonic turbulent jet, at Mach number
$M_{\!j}=U_{\!j}/c_{\infty }=0.4$
and Reynolds number
$Re=\rho _{\!j} U_{\!j} D/\mu _{\!j} = 450\,000$
computed by Brès & Lele (Reference Brès and Lele2019). Here,
$\rho$
is density,
$U$
is velocity,
$D$
is nozzle diameter,
$\mu$
is dynamic viscosity, and
$c$
is speed of sound. The simulations were carried out using the compressible flow solver Charles on an unstructured grid using a finite-volume method with subgrid turbulence generation in the nozzle region (Brès et al. Reference Brès, Ham, Nichols and Lele2017, Reference Brès, Jordan, Jaunet, Le Rallic, Cavalieri, Towne, Lele, Colonius and Schmidt2018). The LES database consists of 20 000 snapshots sampled at an interval
$\Delta tc_{\infty }/D = 0.2$
acoustic time units. Data interpolated on a cylindrical grid spanning
$x,r \in [0, 30] \times [0, 6]$
were used in this analysis. The flow is non-dimensionalised by the nozzle exit values, namely velocity by
$U_{\!j}$
, pressure by
$\rho _{\!j} U_{\!j}^2$
, length by the nozzle diameter
$D$
, and time by
$D/U_{\!j}$
. Frequencies are reported in terms of the Strouhal number
$St = \textit{fD}{\kern-0.5pt}/{\kern-0.5pt}U_{\!j}$
. We focus our analysis on the most energetic wavenumbers
$m = 0, \pm 1, \pm 2, \pm 3, \pm 4$
and frequencies
${St}_3 \in [-2,2]$
within the domain
$x/D \in [0, 20]$
and
$r/D \in [0, 3]$
. The higher wavenumbers and frequencies are important near the nozzle’s exit and will be considered wherever relevant.
To provide context for the triadic interactions and the choice of wavenumbers and frequencies discussed in the following sections, we briefly summarise the key flow physics of the turbulent jet. The jet exhibits three dominant instability mechanisms: KH, Orr and lift-up. Among these mechanisms, KH modes exhibit the highest spatial growth rates, whereas Orr and lift-up mechanisms exhibit the largest transient growth rates for
$m=0$
and
$m\gt 0$
wavenumbers, respectively. In the near-nozzle region, within the first few jet diameters, KH instability dominates as the shear layer develops. Farther downstream, beyond the end of the potential core, the lift-up mechanism becomes more significant. The Orr mechanism is also observed beyond the potential core, although it generally leads to less energetic structures compared to KH and lift-up. Coherent structures associated with the KH instability are observed over a broad frequency range, approximately
$St = 0.3{-}0.8$
, across azimuthal wavenumbers
$m = 0$
,
$\pm 1$
,
$\pm 2$
and
$\pm 3$
, with the most amplified mode occurring near
$St \approx 0.4$
and
$m = 1$
. Although KH structures persist at higher azimuthal wavenumbers, their growth rates diminish. Structures linked to the lift-up mechanism appear at zero frequency (
$St = 0$
) for non-axisymmetric wavenumbers (
$m \gt 0$
) manifesting as streaks. The Orr mechanism, on the other hand, is most prominent for
$m = 0$
and at lower frequencies,
$St \lesssim 0.2$
.
The production of a single frequency and azimuthal wavenumber is defined as
where
$\bar{\strut{\cdot}}$
denotes the temporal mean, and
$\mathcal{R}$
denotes the real part. This quantity represents the energy transfer between the fluctuations and the mean flow at specific
$(St,m)$
. A positive value of
$\mathcal{P}$
indicates that the fluctuation corresponding to this frequency–wavenumber combination extracts energy from the mean flow, while a negative value indicates that it transfers energy back to the mean. In § 3.4, we will discuss production in more detail in the context of the energy budget of streaks.
Figure 1 shows the production term for azimuthal wavenumbers
$m = 0$
, 1 and 2 at three frequencies,
$St = 0.05$
, 0.4 and 0.8. It is observed that fluctuations at all these scales extract energy from the mean flow. However, the location of energy extraction varies with frequencies: lower frequencies tend to extract energy further downstream, while higher frequencies draw energy from the mean flow closer to the nozzle. The production associated with the Orr mechanism (
$m = 0$
,
$St = 0.05$
) occurs downstream of the potential core, and is concentrated near the centreline. In contrast, the KH wavepackets (figure 1
$b$
,
$e$
,
$h$
) extract energy from the mean flow in the initial shear layer, over the range
$x \approx 1$
to
$x \approx 9$
.
Production fields at different frequencies and azimuthal wavenumbers: (
$a$
)
$m=0,St= 0.05$
; (
$b$
)
$m=0,St= 0.4$
; (
$c$
)
$m=0,St= 0.8$
;(
$d$
)
$m=1,St= 0.05$
; (
$e$
)
$m=1,St= 0.4$
; (
$f$
)
$m=1,St= 0.8$
; (
$g$
)
$m=2,St= 0.05$
; (
$h$
)
$m=2,St= 0.4$
; (i)
$m=2,St= 0.8$
.

3.2. Nonlinear transfer and BMD
Quadratic nonlinearities arise from the convective term in the Navier–Stokes equation. Correspondingly, they contribute to the turbulent kinetic energy (TKE) through the nonlinear transport term (Pope Reference Pope2001) given by
Owing to the low Mach number, we have used the incompressible form of nonlinear transport. Using the Fourier transform, the transfer for an azimuthal wavenumber triad
$[m_1,m_2,m_3]$
is given as
Figure 2 shows the spatial field for the total nonlinear transfer and the wavenumber triad
$[0,0,0]$
. Positive (red) and negative (blue) values represent energy gained and lost by the fluctuations, respectively. The total nonlinear energy transfer
$\mathcal{T}_{nl}$
reveals that all nonlinear interactions occur in the shear layer and are concentrated around the lip line (
$r=0.5$
) within the first five jet diameters. Along
$r=0.5$
, the fluctuating velocities lose energy, while they gain energy in the radially outward part of the shear layer. On the other hand, for the self-interaction of the axisymmetric component
$[0,0,0]$
, the nonlinearity is concentrated near the centreline beyond the end of the potential core. Later, in § 4, we will demonstrate that the nonlinear activity near the nozzle’s exit is due to higher wavenumber interactions, and the nonlinear activity at the end of the potential core results from interactions between lower azimuthal wavenumbers.
Spatial field of the nonlinear energy transfer: (
$a$
) total; (
$b$
) azimuthal wavenumber triad
$[m_1,m_2,m_3]=[0,0,0]$
. The white dashed lines indicate the envelope of the jet within which the TKE is greater than 10 % of the maximum TKE.

Figure 3 shows the scaling of the TKE spectrum and bispectrum in the fully developed region of the jet, i.e. we specifically choose the point
$x=20,\ r=0.5$
. The TKE spectrum exhibits an
$St^{-5/3}$
scaling, which is characteristic of the inertial subrange of the turbulent flow. The summed bispectrum in figure 3(
$b$
) is computed as
$ S(f_3) = \sum \nolimits _{f_1+f_2=f_3} E [\hat {\boldsymbol{u}}^{*}_i(x,f_1) \circ \hat {\boldsymbol{u}}^{*}_i(x,f_2) \circ \hat {\boldsymbol{u}}_i(x,f_3) ]$
. The bispectrum demonstrates a
$-3$
scaling. This scaling is typical for the inertial range, as shown by Van Atta (Reference Van Atta1979) and Herring (Reference Herring1980) for isotropic turbulence in the inertial range. Later, in figure 17, we will show how the pattern of the bispectrum changes from the initial shear layer to the fully developed region of the turbulent jet.
The (a) TKE spectrum and (b) bispectrum at
$x=20,\ r=0.5$
, with
$-5/3$
and
$-3$
scalings shown in (
$a$
) and (
$b$
), respectively.

We now perform BMD to identify the dominant triadic interactions. It is computed for blocks containing
$n_{\textit {fft}}=256$
snapshots with
$50\,\%$
overlap, resulting in a total number
$n_{{blk}} = 158$
blocks, which represents a reasonable trade-off between spectral resolution and converged statistics (Schmidt & Colonius Reference Schmidt and Colonius2020; Nekkanti & Schmidt Reference Nekkanti and Schmidt2023). The cross-bispectra for the twelve most dominant azimuthal triads are shown in figure 4. The high-intensity regions signify the dominant triads that arise from the interactions of the two displayed frequencies. All twelve cases exhibit broadband behaviour, with the most intense interactions occurring at low frequencies. Similarly, in the wavenumber space, the intensity decreases for higher azimuthal wavenumber interactions, as seen in
$[0,3,3]$
,
$[1,3,4]$
and
$[2,2,4]$
. A band is observed along
$St_3=0$
for
$[0,0,0]$
,
$[1,-1,0]$
and
$[2,-2,0]$
triads, which is due to spectral leakage, and caution is advised while making interpretations. Spectral leakage was identified by varying the number of frequency bins. As the number of bins increases, the intensity of the observed bands decreases, confirming that these features arise from spectral leakage. While some leakage is present across all triads, it is most pronounced for triads of the form
$[m, -m, 0]$
, since
$m_3 = 0$
and
$St_3 = 0$
correspond to the azimuthally uniform mean flow distortion. The cross-bispectra for the azimuthal triads
$[0,0,0]$
,
$[1,-1,0]$
and
$[2,-2,0]$
exhibit a fourfold symmetry about the lines
$St_1=St_2$
and
$St_1=-St_2$
, whereas the azimuthal triads
$[1,1,2]$
and
$[2,2,4]$
are symmetric about the line
$St_1=St_2$
. Following the discussion in § 2.1, the bispectra of the triads
$[0,1,1]$
,
$[0,2,2]$
,
$[0,3,3]$
exhibit the reflection–translation symmetric (2.20). The bispectra of other triads do not exhibit any symmetries.
Cross-bispectra of twelve azimuthal wavenumber triads
$ [m_1, m_2, m_3]$
. Here, the magnitude of the complex cross-bispectral measure
$|\lambda _1|$
is presented.

Figure 5 shows the cross-bispectra summed over all frequency triads and summed over all diagonals of slope
$-1$
such that
$St_1 + St_2 =St_3$
. The former represents the total intensity of triadic interactions for a single azimuthal wavenumber triad, while the latter denotes the intensity of the resulting frequency for a single wavenumber triad. Figure 5(
$a$
) reveals that the four dominant azimuthal wavenumber triads are
$[1,1,2]$
,
$[0,0,0]$
,
$[2,-1,1]$
and
$[-1,1,0]$
. The self-interaction of
$m=1$
that results in
$m=2$
is the most significant triad. This is expected since
$m=1$
and
$m=2$
are the most energetic components (Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018). The total intensity of the triadic interaction
$[1,-1,0]$
equals
$[0,1,1]$
,
$[2,-2,0]$
equals
$[0,2,2]$
, and
$[1,2,3]$
equals
$[3,-2,1]$
, confirming the bispectrum symmetries discussed in § 2.1. The cross-bispectra, summed along diagonals, are shown in figure 5(
$b$
) for the four dominant wavenumber triads. The four curves decrease monotonically. The triads
$[0,0,0]$
,
$[1,-1,0]$
,
$[1,1,2]$
peak at
${St}_3 \rightarrow 0$
, whereas
$[2,-1,1]$
peaks at the first non-zero frequency,
${St}_3 \approx 0.05$
, as shown in the inset of figure 5(
$b$
). For
${St}_3 \geqslant 1$
, the triad
$[1,1,2]$
exhibits the largest values, indicating that the strongest triadic interactions lead to the azimuthal wavenumber
$m=2$
within this frequency range. The summed cross-bispectrum is symmetric about
${St}_3 = 0$
for
$[0,0,0]$
and
$[1,-1,0]$
triads, but is not symmetric for the
$[1,1,2]$
and
$[2,-1,1]$
triads. This follows from the symmetry of the bispectrum, i.e. if the bispectrum is symmetric about the line
${St}_1=-{St}_2$
, then the summed bispectrum is symmetric about
${St}_3 = 0$
.
Cross-bispectra summed by magnitude: (
$a$
) over all frequency triads; (
$b$
) over diagonals of slope
$-1$
such that
$St_1 + St_2 =St_3$
.

Next, we visualise the spatial flow structures associated with the dominant triads. Figure 6 shows the bispectral modes due to the interaction between
$m=\pm 1$
,
${St}= \pm 0.3$
. These interactions lead to four possible modes: (
$i$
)
$[1,1,2]$
,
$(0.3,0.3,0.6)$
; (
$ii$
)
$[1,1,2]$
,
$(0.3,-0.3,0.0)$
; (
$iii$
)
$[1,-1,0]$
,
$(0.3,0.3,0.6)$
; (
$iv$
)
$[1,-1,0]$
,
$(0.3,-0.3,0.0)$
. We specifically chose
$m=\pm 1$
and
${St}= \pm 0.3$
because these exhibit the largest linear resolvent gain (Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018). In both azimuthal wavenumber triads, sum-interaction of
${St}=0.3$
generates structures reminiscent of KH wavepackets, while difference-interaction of
${St}=0.3$
produces elongated structures localised beyond the end of the potential core. Notably, the structure in figure 6(
$b$
) resembles a large-scale streaky structure, suggesting that the KH mode of
$m=1,\ {St}=0.3$
interacts with its counter-rotating counterpart
$m=1,\ {St}=-0.3$
, giving rise to streaks. These observations are in line with a number of previous studies in free shear flow and jets that have shown the creation of vortices that ultimately give rise to streaks via the interaction of KH modes (Pierrehumbert & Widnall Reference Pierrehumbert and Widnall1982; Wu & Huerre Reference Wu and Huerre2009).
Real component of the fluctuating streamwise velocity of the cross-bispectral mode for two frequency triads and two azimuthal wavenumber triads: (
$a$
,
$b$
)
$[1,1,2]$
and (c,d)
$[1,-1,0]$
. We depict the frequency triads (a,c)
$(0.3,0.3,0.6)$
and (b,d)
$(0.3,-0.3,0.0)$
, respectively. Contours (
) are given by
$\pm \|\boldsymbol{\phi }_1: u_x \|_\infty$
of each mode. The white dashed lines indicate the envelope of the jet within which the TKE is greater than 10 % of the maximum TKE.

To further confirm that the structure in figure 6(
$b$
) is a streak, the magnitude of its three velocity components are examined in figure 7. Figure 7(a–c) reveal that the magnitude of the streamwise velocity is significantly greater (by a factor of 3) than the other two velocity components, in accordance with previous studies (Wu & Huerre Reference Wu and Huerre2009; Boronin, Healey & Sazhin Reference Boronin, Healey and Sazhin2013; Pickering et al. Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020; Nekkanti et al. Reference Nekkanti and Schmidt2023a
). For comparison, figure 7(d–f) show the magnitude of the three velocity components for the frequency triad
$(0.3,0.3,0.6)$
, which corresponds to a KH-type wavepacket. In this case, the magnitudes of the three velocity components are comparable, clearly distinguishing the characteristics of KH wavepackets from streaky structures. Note that only some, not all, BMD modes correspond to KH wavepackets. Specifically, the modes that exhibit spatial patterns similar to figure 6(
$a$
) can be classified as KH wavepackets.
Magnitudes of three velocity components of the cross-bispectral mode for the frequency triads (a–c)
$(0.3,-0.3,0.0)$
and (d–f) (
$0.3,0.3,0.6$
) at the wavenumber triad
$[1,1,2]$
: (a,d)
$u_x$
, (b,e)
$u_r$
, (c, f)
$u_{\theta }$
. Contours (
) range from zero to maximum amplitude of the velocity component
$u_i$
for the BMD mode , i.e.
$ \|\boldsymbol{\phi }_1: u_i \|_\infty$
. The white dashed lines indicate the envelope of the jet within which the TKE is greater than 10 % of the maximum TKE.

3.3. Streaks
We further investigate the role of triadic interactions in streak dynamics. Figure 8 shows the bispectral modes for three azimuthal wavenumber triads,
$[1,1,2]$
,
$[1,2,3]$
and
$[2,2,4]$
, and three frequency triads,
$(0.05,-0.05,0.0)$
,
$(0.4,-0.4,0.0)$
and
$(0.8,-0.8,0.0)$
. Streaks are azimuthally non-uniform structures and are not present in the axisymmetric component (
$m = 0$
). Pickering et al. (Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020) demonstrated that
$m=2$
, 3 and 4 are the most significant contributors to streaks. Therefore, we choose the dominant azimuthal triads that result in these wavenumbers. All nine modes display an elongated region of pronounced streamwise velocity. The lowest azimuthal wavenumber streaks are located beyond the end of the potential core at lower frequencies, shifting upstream at higher frequencies and azimuthal wavenumbers. This upstream shift of the streaks at higher wavenumbers is much more evident in figure 8(c, f,i).
Real components of the fluctuating streamwise velocity of the cross-bispectral mode for three azimuthal wavenumber triads (a–c)
$[1,1,2]$
, (d–f)
$[1,2,3]$
and (g–i)
$[2,2,4]$
, and three frequency triads (a,d,g)
$(0.05,-0.05,0.0)$
, (b,e,h)
$(0.4,-0.4,0.0)$
and (c, f,i)
$(0.8,-0.8,0.0)$
. Contours (
) are given by
$\pm 0.5\,\|\boldsymbol{\phi }_1: u_x \|_\infty$
of each mode. The white dashed lines indicate the envelope of the jet within which the TKE is greater than 10 % of the maximum TKE.

The structure of the BMD streaks is qualitatively similar to those identified using spectral proper orthogonal decomposition (SPOD). To quantify the correspondence, we define an alignment metric between BMD modes and the leading SPOD mode as a normalised inner product:
\begin{equation} \beta _{k,l} = \frac {\big |\left ({\phi }_{k+l}^{\textit{BMD}}\right )^*\boldsymbol{W} {\phi }_{k+l}^{\textit{SPOD}}\big |}{\sqrt {\left ({\phi }_{k+l}^{\textit{BMD}}\right )^*\boldsymbol{W} {\phi }_{k+l}^{\textit{BMD}}}\sqrt {\left ({\phi }_{k+l}^{\textit{SPOD}}\right )^*\boldsymbol{W} {\phi }_{k+l}^{\textit{SPOD}}}}. \end{equation}
The similarity index
$0 \leqslant \beta \leqslant 1$
is unity if the structures are the same. The similarity of BMD modes of the
$[1,1,2]$
triad with the leading SPOD mode of
$m=2$
is shown in figure 9. The SPOD modes are computed via the eigenvalue decomposition of the cross-spectral density matrix
$ \boldsymbol{Q}_l \boldsymbol{Q}_l^* \boldsymbol{\varPhi }_l = \boldsymbol{\varPhi }_l \boldsymbol{\varLambda }_l$
. The columns of
$\boldsymbol{\varPhi }_l = [\boldsymbol{\phi }_l^{(1)}, \boldsymbol{\phi }_l^{(2)}, \ldots , \boldsymbol{\phi }_l^{(n_{{blk}})} ]$
represent the SPOD modes, with
$\boldsymbol{\phi }_l^{(1)}$
being the most energetic or leading mode. The corresponding modal energies are given by the diagonal matrix of eigenvalues
$\boldsymbol{\varLambda }_l = \mathrm{diag}(\lambda _l^{(1)}, \lambda _l^{(2)}, \ldots , \lambda _l^{(n_{{blk}})})$
, ordered such that
$\lambda _l^{(1)} \geqslant \lambda _l^{(2)} \geqslant \cdots \geqslant \lambda _l^{(n_{{blk}})}$
. For further details on the SPOD algorithm and mode computation, see Towne, Schmidt & Colonius (Reference Towne, Schmidt and Colonius2018), Schmidt et al. (Reference Schmidt, Towne, Rigas, Colonius and Brès2018), Schmidt & Colonius (Reference Schmidt and Colonius2020) and Nekkanti & Schmidt (Reference Nekkanti and Schmidt2021). Figure 9(
$a$
) shows that the BMD modes associated with difference interactions, i.e.
$(\pm {St}_1, \mp {St}_2)$
, denoted in black, exhibit high similarity with the leading SPOD mode. This suggests that the difference interactions, rather than sum interactions, exhibit flow structures corresponding to the most energetic coherent structures, which is consistent with the theory of Wu & Huerre (Reference Wu and Huerre2009). Additionally, BMD modes generated due to self-interaction along the line
${St}_1={St}_2$
for
$|{St}_1|=|{St}_2| \leqslant 0.45$
exhibit high similarity, implying that self-interactions in this frequency range also result in the most energetic coherent structures. Since streaky structures are present along the line
${St}_1 + {St}_2=0$
, the similarity index
$\beta$
along this line is as shown in figure 9(
$b$
). The triadic interaction
$[1,1,2]$
,
$(0.05,-0.05,0.0)$
has the highest similarity with the leading SPOD mode. This is expected as the spatial support of the frequency
$St=0.05$
is present downstream where the most energetic streaky structure is also present. This suggests that although any frequency and its counter-rotating wave can interact to produce streaky structures, the most energetic streak occurs due to the interaction between counter-rotating helical modes at the lowest frequency.
Similarity index
$\beta$
in (3.4) of BMD modes and the leading SPOD mode interaction amongst low-frequency modes (leading SPOD mode) that result in a streak (BMD mode).

Here, the spatial domain
$x, r \in [0,30] \times [0,6]$
is used to compute both the BMD and SPOD modes. Within this domain, SPOD identifies the most energetic streaks in the downstream region, which align well with BMD modes associated with low-frequency triadic interactions (
$St_1 \leqslant 0.1$
). Reducing the spatial extent to a smaller domain, such as
$0 \leqslant x \leqslant 15$
, shifts the dominant SPOD structures upstream and enhances their alignment with BMD modes corresponding to higher-frequency triads, such as
$(0.4, -0.4, 0.0)$
. In general, all triads of the form
$(St, -St, 0)$
,
$[1, 1, 2]$
, contribute to streak formation. However, their spatial structure and location are frequency-dependent: lower-frequency interactions produce streaks that dominate downstream, while higher-frequency interactions generate streaks that are dominant upstream.
Figure 10 shows the schematic representation of the frequency triads
$(0.05, -0.05, 0.0)$
and
$(0.4, -0.4, 0.0)$
for the azimuthal wavenumber triad
$[1, 1, 2]$
. In figure 10(
$a$
), the structure corresponding to
$m=1$
,
$St=\pm 0.05$
is the leading SPOD mode, while the structure corresponding to
$m=2$
,
$St=0.0$
is the BMD mode. This schematic illustrates how two helical modes at the same frequency but with opposite sense of rotation are coupled to streaks. Note that the mode corresponding to the low frequency
$St = 0.05$
with azimuthal wavenumber
$m = 1$
lies on the borderline between modal and non-modal instability, as shown by Pickering et al. (Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020, figure 13), and may contain contributions from both mechanisms. However, upon visual inspection, we find that the spatial structure of the mode is reminiscent of a KH-type wavepacket. A similar triad, but at a higher frequency, unambiguously belonging to the KH range, is shown in figure 10(
$b$
). The structure associated with
$m=1$
,
$St=\pm 0.4$
occurs upstream, causing the streaks associated with this triadic interaction to be located further upstream compared to figure 10(
$a$
). A previous study by Wu & Huerre (Reference Wu and Huerre2009) found that perturbing the flow at
$m = \pm 1$
and the same frequency led to large, streamwise elongated structures that modulated the mean flow. Similar findings were reported by Kusek, Corke & Reisenthel (Reference Kusek, Corke and Reisenthel1990), Corke & Kusek (Reference Corke and Kusek1993) and Samimy et al. (Reference Samimy, Kim, Kastner, Adamovich and Utkin2007). These studies, along with the current findings, suggest that the difference interactions between the helical structures of the same frequency can lead to streaks, which are captured by cross-bispectral modes.
Schematic of the triadic interaction between KH wavepacket, its conjugate, and streak for the azimuthal wavenumber triad
$[m_1,m_2,m_3]=[1,1,2]$
for two frequency triads: (
$a$
)
$(0.05,-0.05, 0.0)$
, (
$b$
)
$(0.4,-0.4, 0.0)$
. In (
$a$
) and (
$b$
), the left-hand subplots display SPOD modes for
$m=1,\ St=\pm 0.05$
, and
$m=1,\ St=\pm 0.4$
, while the right-hand plots display BMD modes for
$m=2,\ St=0.0$
. Modes present the streamwise velocity reconstructed in physical space with red–blue isocontours representing
$\pm 0.5\,\|\boldsymbol{\psi }:u_x\|$
for SPOD and
$\pm 0.5\,\|\boldsymbol{\phi }:u_x\|$
for cross-BMD.

3.4. Spectral energy budget of streaks
Now, to quantify the contribution of the linear and nonlinear terms to streaks, we use the spectral TKE equation
\begin{equation} \def \?{\vphantom {\displaystyle \sum _{i=1}}} \frac {\partial \hat {k}}{\partial t} = \mathcal{R} \!\left [ \!\underbrace {\?-\bar {u}_i \frac {\partial \hat {k}}{\partial x_i}}_{\mathcal{A}} \underbrace {\?- \overline {\hat {u}_{\!j}^*\,\widehat { u_i\frac {\partial u_{\!j}}{\partial x_i}}}}_{\mathcal{T}_{nl}} \underbrace {\?-\overline {\hat {u}_{\!j}^*\hat {u}_i\frac {\partial \bar {u}_{\!j}}{\partial x_i}}}_{\mathcal{P}} \underbrace {\? - \frac {2}{\textit{Re}} \overline {\hat {s}_{\textit{ij}}^* \hat {s}_{\textit{ij}}}}_{ \mathcal{D}} \underbrace {\? -\frac {\partial }{\partial x_{\!j}}\!\left (\overline {\hat {u}_{\!j}^*\hat {p}}\right )}_{\mathcal{T}_{p}} \underbrace {\?+\frac {2}{\textit{Re}} \frac {\partial }{\partial x_i} \!\left (\overline {\hat {u}_{\!j}^*\hat {s}_{\textit{ij}}}\right )}_{\mathcal{T}_{\nu }}\right ]\!, \end{equation}
where
$\hat {s}_{\textit{ij}} = 1/2 (\partial \hat {u}_i/\partial x_{\!j} + \partial \hat {u}_{\!j}/\partial x_i )$
is the spectral strain rate,
$\mathcal{R}$
denotes the real part, and
$\hat {k}(f,m)= \overline {\hat {u}_i^{*}(f,m)\,\hat {u}_i(f,m)}/2$
. For statistically stationary flows, the left-hand side is zero. On the right-hand side, the terms
$\mathcal{A}$
,
$\mathcal{T}_{nl}$
,
$\mathcal{P}$
,
$\mathcal{D}$
,
$\mathcal{T}_{p}$
and
$\mathcal{T}_{\nu }$
denote the advection, nonlinear transfer, production, dissipation, pressure transport and viscous transport, respectively. Readers are referred to Bolotnov et al. (Reference Bolotnov, Lahey, Drew, Jansen and Oberai2010) and Nekkanti et al. (Reference Nekkanti, Colonius and Schmidt2025) for the derivation of the spectral TKE equation. The linear contribution to streaks, i.e. the lift-up mechanism, is through the production term
where
$\hat {u}_i (f,m )$
is obtained by collecting all Fourier realisations at frequency
$f$
for the
$m$
th azimuthal wavenumber, similar to (2.1). The triadic contribution is through the nonlinear transfer term
\begin{equation} \mathcal{T}_{nl}\left (f, m\right ) = - \mathcal{R}\overline {\left [\hat {u}_{\!j}^*\left (f , m\right )\,\widehat {u_i\frac {\partial u_{\!j}}{\partial x_i}}\left (f, m\right )\right ]}, \end{equation}
where
$\widehat {u_i\, {\partial u_{\!j}}/{\partial x_i}}$
contains the interactions between different frequencies and azimuthal wavenumbers. Following Cho, Hwang & Choi (Reference Cho, Hwang and Choi2018), we isolate the energy transfer of individual triads, i.e. frequency and wavenumber triplets related by the resonance conditions,
$m_1 \pm m_2 \pm m_3=0$
and
$f_1 \pm f_2 \pm f_3=0$
, by splitting
$\widehat {u_i\, {\partial u_{\!j}}/{\partial x_i}}$
using the discrete convolution to obtain
\begin{equation} \mathcal{T}_{nl}\left (f_3 , m_3\right ) = - \mathcal{R}\overline {\left [\hat {u}_{\!j}^*\left (f_3,m_3\right )\sum \limits _{\substack {f_1 +f_2=f_3 \\ m_1 +m_2=m_3}}\hat {u}_i\left (f_1,m_1\right )\frac {\partial \hat {u}_{\!j}}{\partial x_i}\left (f_2,m_2\right )\right ]}, \end{equation}
\begin{equation} \mathcal{T}_{nl}\left (St=0 , m=2\right ) = - \mathcal{R}{\left [\overline {\sum \limits _{St_1,m_1}\hat {u}_{\!j}^*\left (St=0,m=2\right ) \hat {u}_i\left (St_1,m_1\right )\frac {\partial \hat {u}_{\!j}}{\partial x_i}\left (-St_1,2-m_1\right )}\right ]}. \end{equation}
The incompressible TKE budget is performed for the azimuthal wavenumber
$m=2$
and frequency
${St}=0$
. Only advection, production and nonlinear transfer are significant, and we confirmed that the sum of all terms on the right-hand side of (3.5) is likewise negligible. Figure 11 shows the three relevant terms, of which
$\mathcal{P}$
and
$\mathcal{T}_{nl}$
are more dominant. Interestingly, these two fields are spatially similar but have opposite signs. This indicates that production supplies energy to streaks, while nonlinear transfer extracts energy from them. The production term represents the lift-up mechanism, in which streamwise vortices (or rolls, as identified from the forcing modes in resolvent analysis; Pickering et al. Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020) lift up high-speed fluid from the jet centre and bring down low-speed fluid from the outer jet, inducing streaks. The cumulative effect of nonlinear interactions is the backscatter of energy from the streaks, which ultimately leads to their attenuation. Our results corroborate the concept of the self-sustaining cycle (Waleffe Reference Waleffe1992; Hamilton, Kim & Waleffe Reference Hamilton, Kim and Waleffe1995; Brandt Reference Brandt2014), where the linear lift-up process supplies energy to the streaks, and the nonlinear interactions draw energy from them. Specifically, it has been shown that the nonlinear interactions between streaks and waves feed energy into the waves, whose interactions in turn drive the rolls (Waleffe Reference Waleffe1997; Waleffe & Kim Reference Waleffe and Kim1998; Dessup et al. Reference Dessup, Tuckerman, Wesfreid, Barkley and Willis2018). Finally, the advection term reveals that it extracts energy in the region
$3\leqslant x \leqslant 10$
and supplies energy in the region
$x\geqslant 16$
. However, this term is of much lower significance than production and nonlinear transfer, serving primarily to balance the energy budget. As shown in figure 10(
$b$
), the KH wavepackets are located upstream, whereas the streaks resulting from the counter-rotating KH interaction (i.e. the
$(0.4, -0.4, 0.0)$
mode) appear farther downstream. This spatial separation suggests that the difference interaction of KH waves gives rise to streaks, while the net nonlinear effects lead to their eventual breakdown, as illustrated in figure 11(
$b$
). The findings in this section collectively suggest the following process: (i) nonlinear interactions of KH wavepackets generate streamwise vortices; (ii) these vortices transfer energy to the streaks through the linear lift-up mechanism; and (iii) the streaks ultimately attenuate due to further nonlinear interactions that drain their energy.
Spectral TKE budget (3.5) for
$m=2$
,
${St}=0$
: (
$a$
) advection,
$\mathcal{A}$
; (
$b$
) production
$\mathcal{P}$
; (
$c$
) nonlinear transfer
$\mathcal{T}_{nl}$
.

4. Spatial regions of nonlinear energy transfer and energy cascade
We next focus on identifying the spatial regions where the dominant triadic interactions occur, and investigate the direction of energy transfer in azimuthal wavenumber space. We will use the terms inverse and forward energy cascade to imply the energy transfer from higher azimuthal wavenumbers (or higher frequencies) to lower azimuthal wavenumbers (or lower frequencies), and vice versa.
We first identify the regions where different azimuthal wavenumber triads are most active. To achieve this, the weighted interaction maps are summed over all frequency triads as
$\sum \limits _{k,l} \lambda _1(St_k,St_l)\, | \phi _{k \circ l} \circ \phi _{k+l}|$
. Figure 12 shows the weighted interaction maps for six azimuthal wavenumber triads. Figure 12(a–c) highlight three out of the four dominant wavenumber triads, while figure 12(d–f) illustrate interactions resulting in higher azimuthal wavenumbers. For the dominant triads, the interaction is present in the shear layer and beyond the end of the potential core. The region of nonlinear activity extends farthest downstream for the
$[0,0,0]$
triad. Our findings agree with the works of Tissot et al. (Reference Tissot, Zhang, Lajús, Cavalieri and Jordan2017b
) and Zhang & Wu (Reference Zhang and Wu2023), demonstrating that nonlinearity becomes pronounced downstream of the potential core for the axisymmetric component. The interaction region shifts closer to the nozzle exit for interactions involving higher azimuthal wavenumbers, further supporting our observation from figure 2 that interactions of higher wavenumbers play a more significant role in the nonlinearities closer to the nozzle exit.
Integral interaction maps
$\sum _{k,l} \lambda _1(St_k,St_l)\, | \phi _{k \circ l} \circ \phi _{k+l}|$
for azimuthal triads
$[m_1,m_2,m_3]$
: (
$a$
)
$[0,0,0]$
, (
$b$
)
$[1,-1,0]$
, (
$c$
)
$[1,1,2]$
, (
$d$
)
$[2,2,4]$
, (
$e$
)
$[3,3,6]$
, (
$f$
)
$[5,5,10]$
. Contours (
) range from zero to maximum amplitude of the integral interaction map, i.e.
$ \|\sum _{k,l} \lambda _1(St_k,St_l)\, \phi _{k \circ l} \circ \phi _{k+l}\|_\infty$
. The white dashed lines indicate the envelope of the jet within which the TKE is greater than 10 % of the maximum TKE.

Figure 13 displays the interaction maps for six representative frequency triads of the most dominant azimuthal wavenumber triad
$[1,1,2]$
. This figure highlights the regions where nonlinear interactions are most active in frequency space, complementary to the wavenumber interactions shown in figure 12. The interaction maps indicate that significant nonlinear interactions between lower frequencies occur downstream, whereas for interactions involving higher frequencies, this region shifts towards the nozzle exit. Overall, figures 12 and 13 reveal that small-scale structures are more actively involved in nonlinear interactions in the initial shear layer, whereas large-scale structures become more active downstream. This observation is not surprising and is consistent with predictions from linear stability analysis that as the shear layer spreads downstream, it supports structures with larger wavelengths farther downstream than those with smaller wavelengths.
Interaction maps
$| \phi _{k \circ l} \circ \phi _{k+l}|$
of the azimuthal wavenumber triad
$[1,1,2]$
for different frequency triads
$(St_1,St_2,St_3)$
: (
$a$
)
$(0.15,-0.15,0.0)$
, (
$b$
)
$(0.15,0.15,0.3)$
, (
$c$
)
$(0.3,0.3,0.6)$
, (
$d$
)
$(0.5,0.5,1.0)$
, (
$e$
)
$(0.75,0.75,1.5)$
, (
$f$
)
$(1.0,1.0,2.0)$
. Contours (
) range from zero to maximum amplitude of the interaction map, i.e.
$ \|\phi _{k \circ l} \circ \phi _{k+l}\|_\infty$
. The white dashed lines indicate the envelope of the jet within which the TKE is greater than 10 % of the maximum TKE.

Figure 14 shows the normalised TKE and bispectrum as functions of streamwise location. The TKE is computed at each spatial location and then integrated radially for each azimuthal wavenumber. The bispectral measure
$b(x,m,{\textit{St}})$
is computed as
Finally, at each streamwise location, the TKE and bispectral measure are normalised by total TKE and total bispectral measure
$\sum _{m,St} b (x,m,St )$
. A similar trend is observed in both the TKE and bispectrum. Near the nozzle exit, higher azimuthal wavenumbers are more energetic and show greater triadic activity. As the flow progresses downstream, the energy and bispectral measure shift to lower azimuthal wavenumbers. The most energetic and triadically dominant azimuthal wavenumber scale as
$1/\sqrt {x}$
. The consistent scaling and trends between figures 14(
$a$
) and 14(
$b$
) suggest that the azimuthal wavenumbers involved in the strongest nonlinear interactions are also the most energetic. Previous studies (Arndt, Long & Glauser Reference Arndt, Long and Glauser1997; Citriniti & George Reference Citriniti and George2000; Sasaki et al. Reference Sasaki, Cavalieri, Jordan, Schmidt, Colonius and Brès2017) have shown that in the region
$x\leqslant 3$
, modes up to
$m \leqslant 10$
contain significant energy. Downstream of this region, modes with
$m \gt 10$
lose energy, while those with
$m \lt 10$
gain energy. Eventually, in the self-similar region, most of the energy is found for
$m=0$
, 1 and 2. This explains why higher azimuthal wavenumbers are more involved in nonlinear interactions upstream, while lower azimuthal wavenumbers are more active downstream.
Normalised (
$a$
) TKE and (
$b$
) bispectrum as functions of streamwise location and azimuthal wavenumber.

Dominant azimuthal wavenumber triads at different streamwise locations.

To illustrate how the bispectral measure varies across different frequencies, figure 15(a,b) show the bispectral measure at
$x=1$
and 10 for
$m=2$
and 12. At
$x=1$
, the bispectral measure for
$m=2$
and
$m=12$
exhibits peaks, with global maxima at
${St}=1$
and
${St}=0.4$
, respectively. Further downstream at
$x=10$
, the bispectral measure for both
$m=2$
and
$m=12$
becomes broadband. Since the flow is developing in the initial shear layer, the bispectral measure exhibits peaks at relatively higher frequencies. Downstream of the potential core, the flow becomes fully developed, and the bispectrum is broadband. This behaviour resembles the inertial subrange energy spectrum, indicating a forward energy cascade in this region.
Bispectral measure (4.1), for azimuthal wavenumbers
$m=2$
and
$m=12$
, at two streamwise location: (
$a$
)
$x=1$
, (
$b$
)
$x=10$
.

Table 1 tracks the dominant azimuthal wavenumber triads as the flow progresses downstream. To identify these dominant triads, we first compute the spatial triple correlation of velocity components for various azimuthal wavenumber triplets. Next, we average these correlations over time. Finally, the triads with the highest values are selected. Note that for this spatial triple correlation, the intensity of the triad
$[m_1, m_2, m_3]$
is equal to the intensities of
$[m_3, -m_2, m_1]$
,
$[m_3, -m_1, m_2]$
and
$[-m_1, -m_2, -m_3]$
. Therefore, only the positive triads
$[m_1, m_2, m_3]$
where
$m_1, m_2, m_3 \geqslant 0$
are reported in table 1. Additionally, the
$[0, 0, 0]$
triad is omitted from the analysis, as it does not provide information on the direction of energy transfer in azimuthal wavenumber space. In the initial shear layer at
$x=1$
, the most dominant triad is
$[6,6,12]$
. At
$x=1$
,
$m=12$
has the highest TKE, suggesting that the most energetic wavenumber is also involved in the dominant nonlinear interactions. Moving downstream, the dominant azimuthal triad transitions to
$[3,3,6]$
at
$x=2$
,
$[2,2,4]$
at
$x=3$
, and
$[1,1,2]$
at
$x=4$
, remaining the same until
$x=8$
. The transition of the dominant triad from higher to lower azimuthal wavenumbers is evident up to
$x \approx 4$
, suggesting that the inverse energy cascade is active in the developing region of the jet. Downstream of the potential core, the dominant triad increases to
$[2,2,4]$
, indicating that the forward cascade is active in the fully developed region. This is not surprising, as the TKE spectrum in this region is broadband and exhibits inertial subrange scaling (see figure 3), which is typical of a forward energy cascade.
The radially integrated TKE at different streamwise locations and azimuthal wavenumbers in the developing region of the jet is shown in figure 16. Figure 16(
$a$
) shows that the spectrum initially exhibits a distinct global maximum at
$St = 0.74$
for
$x = 1$
, which gradually transitions into a broadband spectrum by
$x = 10$
. As the flow evolves downstream, the frequency corresponding to the peak TKE shifts towards lower values, further indicating the presence of an inverse cascade in frequency space. Furthermore, at downstream locations, the radially integrated TKE spectrum is elevated, suggesting that other frequencies are also gaining energy from the mean flow, indicating that production remains active in this initial shear layer region. The spectra at downstream locations,
$x \leqslant 8$
, exhibit a sharp drop-off at higher frequencies due to the modelling of unresolved small scales. Figure 16(
$b$
) presents the radially integrated TKE at
$x = 2$
for different azimuthal wavenumbers
$m$
, showing that the peak frequency shifts to lower values as
$m$
increases. This trend is indicative of an inverse cascade in wavenumber space. These observations, together with the production fields in figure 1 and the data in table 1, reveal that in the developing region of the jet, low frequencies gain energy from both the mean flow and smaller scales via inverse cascade.
Radially integrated TKE spectra at different streamwise locations and azimuthal wavenumbers: (a) spectra at
$x=1,2,3,\ldots,10$
for
$m = 1$
; (
$b$
) spectra at
$x = 2$
for azimuthal wavenumbers
$m=0,1,2,3,4$
.

Mode bispectra for the azimuthal wavenumber triads (a–d)
$[1,1,2]$
, (e–h)
$[2,2,4]$
, at different streamwise locations: (
$a$
,
$e$
)
$x=1$
, (
$b$
,
$f$
)
$x=5$
, (
$c$
,
$g$
)
$x=15$
, (
$d$
,
$h$
)
$x=30$
.

Next, we examine the patterns of the bispectrum in different spatial regions. Figure 17 shows the mode bispectra of two azimuthal wavenumber triads,
$[1,1,2]$
and
$[2,2,4]$
, computed at four different streamwise locations,
$x=1$
, 5, 15 and 30, i.e. ranging from the initial shear layer to the self-similar region. These mode bispectra are computed by performing BMD on different streamwise planes for the
$[1,1,2]$
and
$[2,2,4]$
triads. At
$x=1$
, the bispectra in figure 17(
$a$
,
$e$
) exhibit a six-lobed pattern, with lobes along the lines
${St}_1={St}_2$
,
${St}_1=-2\, {St}_2$
and
${St}_1=- {St}_2/2$
. Previous work by Herring & Métais (Reference Herring and Métais1992) on decaying turbulence also identified a similar bispectrum pattern, which they referred to as a ‘six-leaf rose’. As the flow evolves downstream, the six lobes coalesce into one lobe about
$({St}_1,{St}_2) = (0,0)$
. This is expected, as the
$m=1$
mode in the initial shear layer exhibits spectral peaks. The sum interaction of the frequencies corresponding to these peaks results in lobes along the line
${St}_1={St}_2$
, while the difference interactions generate lobes along the lines
${St}_1=-2\,{St}_2$
and
${St}_1=- {St}_2/2$
. Downstream of the potential core, the spectrum of
$m=1$
becomes broadband, and its self-interaction results in a single lobe at the origin. Similar trends are observed for the
$[2,2,4]$
triad. This figure implies that for a developing flow, the bispectrum pattern is similar to figure 17(
$a$
), while for a fully developed broadband flow, it is similar to figure 17(
$f$
).
Next, the direction of energy transfer for the triads reported in table 1 is investigated. Figure 18 shows the nonlinear energy transfer (3.3) for eight triads,
$[12,-6,6]$
,
$[6,6,12]$
,
$[6,-3,3]$
,
$[3,3,6]$
,
$[4,-2,2]$
,
$[2,2,4]$
,
$[2,-1,1]$
and
$[1,1,2]$
. For a triad
$[m_1,m_2,m_3]$
, the red colour denotes the energy transfer from
$m_1$
or
$m_2$
to
$m_3$
, and the blue colour denotes the energy extracted from
$m_3$
by
$m_1$
or
$m_2$
. The spatial fields in figure 18(a,c,e,g) display a region of positive nonlinear transfer, followed by a region of negative energy transfer. In contrast, the energy transfer patterns in figure 18(b,d, f,h) are opposite to those in figure 18(a,c,e,g). However, both sets convey the same information. For example, the spatial field corresponding to the
$[12,-6,6]$
triad demonstrates energy transfer from
$m= 12$
to
$m=6$
until
$x \approx 2.8$
, after which there is a backscatter of energy from
$m=6$
to
$m=12$
. Similarly, the
$[6,6,12]$
triad indicates that initially,
$m=6$
extracts energy from
$m=12$
, but later, the energy transfers from
$m=6$
to
$m=12$
. These observations clearly show that the inverse energy cascade is active in the initial shear layer, while the forward energy cascade becomes prominent at locations further downstream.
Nonlinear transfer fields for different azimuthal triads
$[m_1,m_2,m_3]$
: (
$a$
)
$[12,-6,6]$
, (
$b$
)
$[6,6,12]$
, (
$c$
)
$[6,-3,3]$
, (
$d$
)
$[3,3,6]$
, (
$e$
)
$[4,-2,2]$
, (
$f$
)
$[2,2,4]$
, (
$g$
)
$[2,-1,1]s$
, (
$h$
)
$[1,1,2]$
.

Interestingly, the spatial fields in figure 18(
$g$
,
$h$
) demonstrate that the energy transfer is unidirectional, with
$m=2$
supplying energy to
$m=1$
. It is noted that within the current streamwise extent (
$x \leqslant 30$
), there is no observed energy transfer from
$m=1$
to
$m=2$
; however, this does not rule out the possibility of such transfer occurring farther downstream. Additionally, the triad
$[1,-1,0]$
reveals that
$m=0$
extracts energy from
$m=1$
, not shown here for brevity. To summarise, for the spatial domain under investigation,
$m\geqslant 2$
are involved in both forward and inverse energy cascades, while the azimuthal wavenumbers
$m=0$
and 1 act as energy sinks. Specifically,
$m=1$
gains energy from
$m\geqslant 2$
, and
$m=0$
gains energy from all
$m$
.
Figure 19 shows the net nonlinear energy transfer into the azimuthal wavenumber
$m$
. This is computed as
\begin{equation} \mathcal{T}_{nl}(m) = \sum_{m_1} -\mathcal{R}\Bigg [\overline {\hat {\boldsymbol{u}}^*_{\!j}(m)\, \hat {\boldsymbol{u}}_i(m_1)\frac {\partial \hat {\boldsymbol{u}}_{\!j}}{\partial x_i}(m-m_1)} \Bigg ]. \end{equation}
All the triadic interactions supply energy to the axisymmetric component
$m=0$
. For
$m=1$
and
$m=2$
, triadic interactions extract energy in the shear layer, but supply energy near the centreline towards the end of the potential core. This trend is also observed for
$m=2$
. For higher azimuthal wavenumbers, the loss of energy is more prominent than the gain. This is more clearly reflected in figure 19(c–f). Although individual triads both supply and extract energy to each azimuthal wavenumber, as shown in figure 18, the net triadic interactions result in an energy gain for lower azimuthal wavenumbers, and an energy loss for higher azimuthal wavenumbers.
Net nonlinear energy transfer fields (4.2) for different azimuthal wavenumbers: (
$a$
)
$m=0$
, (
$b$
)
$m=1$
, (
$c$
)
$m=2$
, (
$d$
)
$m=3$
, (
$e$
)
$m=4$
, (
$f$
)
$m=5$
.

5. Discussions and conclusions
Triadic interactions transfer energy between scales in turbulent flows. Using bispectral mode decomposition (BMD), we characterise these interactions in a Mach 0.4 turbulent jet. Our findings show that BMD can effectively identify the dominant triadic interactions and extract related flow structures, even in broadband turbulence.
The Kelvin–Helmholtz (KH) mechanism, the most amplified linear mechanism, is the primary component of triadic interactions. The BMD reveals that its self-interaction leads to higher helical structures, while its difference self-interaction leads to elongated structures. Specifically, we identify a strong triadic interaction between KH wavepackets corresponding to counter-rotating helical modes and streaks, spanning various azimuthal wavenumbers. The most dominant interactions occur between
$(m=1,{St})$
,
$(m=1, -{St})$
and
$(m=2, {St}=0)$
. While this interaction is present across a broad frequency range, the most energetic streaks arise from lower-frequency interactions, such as the triad
$(0.05,-0.05,0.0)$
. These findings provide the first direct data verification of the theoretical result that the mutual interaction of KH waves associated with wavenumbers
$\pm m$
forms vortices (or rolls), which subsequently give rise to streaks (Benney & Lin Reference Benney and Lin1960; Pierrehumbert & Widnall Reference Pierrehumbert and Widnall1982). This behaviour is further supported by Wu & Huerre (Reference Wu and Huerre2009) and Zhang & Wu (Reference Zhang and Wu2023), who demonstrated that the quadratic interaction of helical modes generates streaks of large amplitude. Furthermore, a spectral energy budget analysis of streaks reveals that the contribution of nonlinear interactions is comparable in magnitude to the linear production term. While the linear mechanism supplies energy to streaks, the nonlinear interactions serve to extract energy from them. This energy loss due to the triadic interactions is analogous to the self-sustaining mechanism in wall-bounded flows, where the nonlinear interactions drain energy from streaks, feed it into the waves, and subsequently drive the rolls (Waleffe Reference Waleffe1997; Waleffe & Kim Reference Waleffe and Kim1998). Note that this is not inconsistent with the observation that interactions among KH wavepackets lead to streak formation. Essentially, two distinct types of nonlinear interactions are present. First, triadic interactions among KH wavepackets generate streamwise vortices or rolls that act as a forcing for the linear lift-up mechanism, which subsequently produces streaks. Second, once the streaks have formed, nonlinear interactions become significant again, specifically between the streaks and all other frequencies, resulting in energy loss and eventual attenuation of the streaks.
Next, a global and local analysis of triadic interactions across different scales is presented. Globally, across the spatial domain considered, the triadic interactions are most prominent among lower azimuthal wavenumbers
$m \leqslant 2$
, and lower frequencies
${St} \leqslant 0.5$
. This finding suggests that in a turbulent jet with a tripped boundary layer, the larger scales, which exhibit the highest energy, also play a dominant role in triadic interactions. In contrast, the local analysis reveals a distinct spatial variation: triadic interactions among smaller scales are more pronounced in the initial shear layer, whereas interactions among larger scales dominate downstream, beyond the potential core. More specifically, in the shear layer region
$0 \leqslant x \leqslant 6$
, the dominant azimuthal wavenumber triad transitions from
$[6,6,12]$
to
$[3,3,6]$
, then to
$[2,2,4]$
, and finally to
$[1,1,2]$
. However, this trend reverses in the fully developed region
$x \geqslant 6$
, with the dominant triad evolving from
$[1,1,2]$
to
$[2,2,4]$
. By computing the nonlinear transfer term, it is found that energy is transferred from higher azimuthal wavenumbers to lower wavenumbers in the shear layer, while beyond the end of the potential core, energy transfer occurs from lower to higher wavenumbers. These results suggest that as the flow exits the nozzle, it initially experiences an inverse energy cascade, which then transitions to a forward energy cascade as the flow progresses downstream.
Our results leave some open questions for future research. Previous studies by Wu & Huerre (Reference Wu and Huerre2009) and Suponitsky et al. (Reference Suponitsky, Sandham and Morfey2010) show that the quadratic interaction of
$m=\pm 1$
radiates low-frequency sound at angles
$60{-}70^\circ$
relative to the jet axis. Therefore, one possible direction is to investigate the effect of nonlinear interactions on the far-field jet noise. Another promising avenue is to examine the impact of nonlinear interactions on mean-flow distortion (Long & Petersen Reference Long and Petersen1992; Zhang & Wu Reference Zhang and Wu2022). This distortion is particularly pronounced in non-axisymmetrically forced jets. Investigating how such forcing mechanisms (Corke, Shakib & Nagib Reference Corke, Shakib and Nagib1991; Long & Petersen Reference Long and Petersen1992; Nekkanti et al. Reference Nekkanti and Schmidt2023b
) contribute to mean flow distortion could yield valuable insights. It is also worth investigating whether an interaction similar to KH–streak interaction occurs in supersonic jets, and if nonlinear interactions exist between screech tones and streaks.
Funding
The authors are grateful for support from the US Office of Naval Research under contract N00014-23-12650, with Dr S. Martens as programme manager. O.T.S. acknowledges support from ONR grant N00014-20-1-2311.
Declaration of interests
The authors report no conflict of interest.
















































































































































































