2.1. Large-eddy simulations

The simulations of the turbulent supersonic twin-rectangular jet are carried out using Cadence’s unstructured, compressible LES solver ‘Charles’ (Brès et al. Reference Brès, Ham, Nichols and Lele2017; Brès & Lele Reference Brès and Lele2019). The jet is nominally ideally expanded and cold. The nozzle exit conditions are characterised by the jet Mach number, $M_{\!{j}}=u_{\!{j}}/c_{\!{j}}=1.5$ , acoustic Mach number, $M_{{a}}=u_{\!{j}}/c_\infty =1.25$ , pressure, $p_{\!{j}}/p_\infty =1$ , temperature, $T_{\!{j}}/T_\infty =0.69$ , and Reynolds number, $Re_{\!{j}}=\rho _{\!{j}} u_{\!{j}} D_{{e}}/\mu _{\!{j}}=1.07\times 10^6$ , where $u$ is the streamwise velocity, $c$ the speed of sound, $\rho$ the density, $D_{{e}}/h=1.6$ the equivalent nozzle diameter, $h$ the nozzle height, $\mu$ the dynamic viscosity, and $(\boldsymbol{\cdot })_{\!{j}}$ and $(\boldsymbol{\cdot })_\infty$ refer to jet exit and ambient conditions, respectively. The nozzle pressure (NPR) and temperature ratios (NTR) are $p_{{t}}/p_\infty =3.671$ and $T_{{t}}/T_\infty =1$ , respectively. The biconical nozzle geometry, with an aspect ratio of two, is included in the computational domain and closely matches the experimental set-up of Samimy et al. (Reference Samimy, Webb, Esfahani and Leahy2023). The centre-to-centre spacing between the nozzles is $3.6h$ . Each nozzle has a cavity cut into the internal wall just upstream of the nozzle exit. Housed within the cavity are eight pairs of electrodes per nozzle, allowing for up to eight plasma actuators per nozzle. Due to the sharp throat and cavity of the nozzles, shocks are present despite the jet being nominally ideally expanded. Aided by the high Reynolds number, the cavity also facilitates a turbulent boundary layer state at the exit, captured using a wall model. The computational grid contains approximately 77 million cells. The simulation time step, $\textrm{d}{t}\,c_\infty /h$ , is 0.002 for the natural jet and 0.001 for the forced jet (see § 2.2). These parameters are summarised in table 1. For full details of the simulation, including its validation against the experiments of Samimy et al. (Reference Samimy, Webb, Esfahani and Leahy2023), we refer the reader to Brès et al. (Reference Brès, Yeung, Schmidt, Esfahani, Webb, Samimy and Colonius2021, Reference Brès, Bose, Ivey, Emory and Ham2022). The flow variables are non-dimensionalised by the ambient conditions, $\rho _\infty$ , $T_\infty$ and $c_\infty =\sqrt {\gamma p_\infty /\rho _\infty }$ , where $\gamma$ is the ratio of specific heats. Lengths are non-dimensionalised by $h$ . Frequencies are non-dimensionalised by $u_{\!{j}}/D_{{e}}$ , and reported as the Strouhal number, $fD_{{e}}/u_{\!{j}}$ . For notational compactness, we denote frequencies by $f$ , but emphasise that $f$ is dimensionless.

Table 1. LES parameters for the natural and forced jets.

2.2. Plasma actuation model

To include the effect of plasma actuation in the LES, we adapt the model proposed by Kim et al. (Reference Kim, Afshari, Bodony and Freund2009). We insert a cylinder-shaped, volumetric source term into the right-hand side of the energy equation. The energy source is expressed as

(2.1) \begin{equation} S(r,z,t) = k(t) l(r) m(z) P/\big(\pi r_0^2 L\big), \end{equation}

where $P$ is the amplitude with the dimensions of power, $r_0$ is the radius of the cylinder, $L$ is its length, and $r$ and $z$ are local coordinates. For notational simplicity only, in (2.1), we have aligned the local $z$ -axis with the axis of the cylinder, such that $r=\sqrt {x^2+y^2}$ . The exact geometry, including the locations and orientations of the actuators, can be found from Samimy et al. (Reference Samimy, Webb, Esfahani and Leahy2023). In actual implementation, the cylinder – like the plasma arc – is aligned with the nozzle lips. The envelope of the cylinder is given by the shape functions

(2.2) \begin{equation} l(r) = (1/2)[ \tanh (-\sigma (r/r_0-1/2)) +1] \end{equation}

and

(2.3) \begin{equation} m(z) = -(1/2)[\tanh {(-\sigma (z+L/2)/r_0)}+1] + (1/2)[\tanh {(-\sigma (z-L/2)/r_0)}+1], \end{equation}

where $\sigma$ is a smoothness parameter. The energy source is modulated in time by the smoothed rectangular wave

(2.4) \begin{equation} k(t) = (1/2)[ \tanh ((t-n\tau -t_{\textit{on}})/t_{{r}}) - \tanh ((t-n\tau -t_{\textit{off}})/t_{{r}}) ], \end{equation}

where $t_{\textit{on}}$ and $t_{\textit{off}}$ are the time instants within each forcing period at which the signal is switched on and off, respectively, $t_{{r}}$ is the rise time, $\tau$ is the forcing period, and $n=\lfloor t/\tau \rfloor$ .

We assume $r_0$ to be half the depth of the cavity and $L$ to be the distance between each pair of electrodes. The value for $\sigma$ is taken directly from Kim et al. (Reference Kim, Afshari, Bodony and Freund2009). The amplitude, $P$ , and temporal parameters, $t_{\textit{on}}$ , $t_{\textit{off}}$ and $t_{{r}}$ , are adapted from voltage and current measurements by Samimy et al. (Reference Samimy, Webb, Esfahani and Leahy2023) of a typical actuation cycle. Due to the short pulse duration, $t_{\textit{off}}-t_{\textit{on}}$ , the rectangular wave, $k$ , in (2.4) is similar to an impulse train. In general, these parameters can be set independently for each actuator; in this work, however, all actuators behave identically. The parameters are summarised in table 2. Based on our previous plasma modelling effort (Brès et al. Reference Brès, Yeung, Schmidt, Esfahani, Webb, Samimy and Colonius2021), we apply additional grid refinement to the vicinity of the cavity, locally reducing the cell width by a factor of two relative to the natural jet (see Appendix A). The overall grid size is nearly unchanged. However, due to explicit time integration, the reduced minimum grid spacing necessitates a smaller time step (see table 1).

Table 2. Non-dimensionalised plasma actuator modelling parameters. Ambient temperature and pressure are assumed to be 293 K and 1 atm, respectively.

Figure 1. Instantaneous snapshots of the (a,c) natural and (b,d) forced jets: (a,b) $Q$ -criterion isocontours of $Q=5$ ; (c,d) numerical schlieren, $|\boldsymbol{\nabla }\rho |$ , on the major-axis plane, $y=0$ . In panel (a,b), the contours are coloured by temperature fluctuations, $T'$ . The colours saturate at $|T'|=\pm 0.05$ . In panel (c,d), the shading varies from white, $|\boldsymbol{\nabla }\rho |=0$ , to black, $|\boldsymbol{\nabla }\rho |\geqslant 10$ . Dashed lines mark the locations of grid density transitions.

Since it is challenging to obtain experimental measurements of the local flow environment around a plasma arc, rather than attempt to replicate a specific experiment, we seek to demonstrate the control authority of plasma actuation by exciting specific frequency and $D_2$ symmetry components in the twin-rectangular jet, as motivated in § 1.

4.1. SPOD modal energies and modes

Figure 4. Leading SPOD eigenvalue spectra of the natural (solid lines) and forced (faded lines) jets: (a) SS; (b) SA; (c) AS; (d) AA symmetry components. The dotted line in panel (a) marks the SS forcing frequency, $f_0$ . The inset in panel (c) zooms in on the frequency range $f\in [0,0.1]$ and shows $\lambda$ on a linear scale.

The leading SPOD eigenvalue spectra of the natural and forced jets are reported and compared in figure 4. In both cases, the modal energy distribution evidently depends on spatial symmetry. For the SS and SA symmetry components in figure 4(a,b), the energy of the natural jet is predominantly broadband. Over the frequency range $f\gtrsim 0.04$ , the energy decays with frequency, as is expected in a turbulent flow. At very low frequencies, the energy rapidly falls off as $f\to 0$ because the finite domain size limits the resolvability of structures with the longest wavelength (and thus lowest frequency). In the SS spectrum, a small peak can also be discerned at $f=0.28$ . For the AS and AA symmetries in figure 4(c,d), however, the spectra of the natural jet reveal a prominent tone at $f=0.29$ as a result of screech resonance. The screech mechanism is thus well represented by AS and AA modes. For the present nozzle geometry and operating condition, screech modes antisymmetric about the major axis were previously observed by Esfahani, Webb & Samimy (Reference Esfahani, Webb and Samimy2021). The screech tone is more energetic in the AS component – a flapping screech mode – than in AA – a double-flapping screech mode. Screech tones are also present at higher frequencies, in particular at $f\approx 0.34$ and 0.75, which are not harmonics of the dominant screech tone. They have substantially lower energy and stem from shock cells having non-uniform spacing (Edgington-Mitchell et al. Reference Edgington-Mitchell, Li, Liu, He, Wong, Mackenzie and Nogueira2022). Compared with the other three symmetries, the AA spectrum has significantly reduced energy at low frequencies, below $f\approx 0.1$ . Conversely, at low frequencies, the SS symmetry is the most energetic. This is consistent with the strong in-phase coupling between the left and right jets observed at low frequencies by Samimy et al. (Reference Samimy, Webb, Esfahani and Leahy2023). The suboptimal eigenvalues do not display tones and are negligibly dependent on symmetry or influenced by forcing, and thus are not reported.

Whereas in the natural jet, tones are found in the AS and AA symmetries, in the forced jet, they have migrated to the SS and AS symmetries. These tones occur at the fundamental forcing frequency, $f_0=0.29$ , as well as its harmonics. They make up part of the response of the jet to the exogenous forcing. Non-harmonic tones found in the AS spectrum of the natural jet are eliminated in the forced jet. The forced jet also responds by completely suppressing the double-flapping tones in the AA component, producing a broadband spectrum. Nonetheless, the appearance of additional tones in SS and AS means that the present forcing strategy does not achieve overall noise reduction. The SA spectrum remains broadband and is nearly unchanged from the natural jet. Given the SS symmetry of the forcing, the generation of harmonics in the AS component and quelling of tones in AA are clear signs that nonlinear mechanisms are at work.

For the AS symmetry in figure 4(c), a comparison between the natural and forced jets shows that the leading modes of both cases possess approximately equal energy at $f_0$ . Below $f\approx 0.1$ , however, the two display distinct characteristics. As the inset in figure 4(c) makes clear, for $0.01\lesssim f\lesssim 0.07$ , the forced jet has higher energy, suggesting it exhibits a greater prevalence of slowly evolving structures. These changes are unsurprising considering the drastically modified mean flow shown in figure 3. Below $f\approx 0.01$ , the separation between the natural and forced spectra initially shrinks, but diverges again as $f\to 0$ . In the $f=0$ frequency bin, the forced jet has tenfold the energy of the natural jet. By definition, the long-time mean cannot be spatially antisymmetric. As such, we attribute the energy of the AS component of the forced jet in the $f=0$ bin to coherent structures energised by the forcing that oscillate at frequencies below the finite bin width, $\Delta f=1/(n_{\textit{fft}}\Delta t)$ . For the remaining three symmetry components, we observe no significant change to the low-frequency region of each of their spectra. We will revisit the modes at near-zero frequencies in §§ 5.4 and 5.5.

Figure 5. Leading SPOD modes, scaled by their SPOD amplitudes, of the (a–c) natural and (d–f) forced jets at $f=0.29$ : (a,d) SS; (b,e) AS; (c, f) AA symmetry components. For each mode, isocontours of $\sqrt {\lambda _1}{Re}\{\boldsymbol{\phi} _1(\boldsymbol{x})\}=\pm d$ are shown in red and blue. The isovalue, $d$ , is shared in each column: (a,d) $d=0.1$ ; (b,e) $d=0.05$ ; (c, f) $d=0.02$ . Rectangles mark the exits of the twin nozzles.

The leading SPOD modes of the natural and forced jets at frequency $f_0$ are shown in figure 5. By construction, SPOD modes are normalised such that they have unit energy, $\|{\boldsymbol{\phi} _1}\|^2_{\boldsymbol{x}}=1$ . To facilitate comparison of the relative presence of the modes as they occur in the data, we scale each mode by its SPOD amplitude, i.e. the square root of its corresponding mode energy, $\sqrt {\lambda _1}$ . As the eigenvalue spectra in figure 4 reveal, only the SS, AS and AA symmetry components of the natural or forced jet (or both) display tones. For this reason, we focus on the modes of these three symmetries in figure 5.

The three-dimensional (3-D) mode shapes clearly demonstrate the contrasting coupling behaviour of each symmetry component. For the SS symmetry, figure 5(a,d) shows the contours of the amplitude-scaled leading mode, $\sqrt {\lambda _1}{Re}\{\boldsymbol{\phi} _1(\boldsymbol{x})\}=\pm 0.1$ , in red (positive) and blue (negative). At this isovalue, no structure is visible in the natural jet due to the low mode energy. Conversely, in the forced jet, we observe coherent structures phase-locked to the strong SS excitation. The structures exit the twin nozzles independently, then rapidly begin to merge, so that by $x\approx 5$ , they have morphed into a single, large wavepacket that couples the left ( $z\gt 0$ ) and right ( $z\lt 0$ ) jets in phase with each other. For the AS symmetry, contours of $\sqrt {\lambda _1}{Re}\{\boldsymbol{\phi} _1(\boldsymbol{x})\}=\pm 0.05$ are shown in figure 5(b,e). Because the natural and forced jets possess almost equal energy at $f=0.29$ in this symmetry (as the spectra in figure 4 c indicate), the modes of both cases resemble each other. Like the SS modes, the AS modes represent in-phase coupling between the left and right jets. Unlike SS, however, the AS modes are antisymmetric about the major-axis plane, $y=0$ . The AA modes in figure 5(c, f) display stark differences between the natural and forced cases. Whereas in the natural jet the contours of $\sqrt {\lambda _1}{Re}\{\boldsymbol{\phi} _1(\boldsymbol{x})\}=\pm 0.02$ reveal distinct structures, in the forced jet, no discernible structure is present for the same isovalue. This is a direct consequence of the suppression of screech in the AA symmetry by the forcing. The envelope of the AA mode in the natural jet is qualitatively similar to the AS mode. That said, the former is distinguished by its antisymmetry about the minor-axis plane, $z=0$ . The left and right jets are coupled perfectly out of phase, with a phase difference of $\pi$ between each other. For the modes considered in figure 5, while the mode amplitudes are (in the cases of SS and AA) significantly altered by forcing, the mode shapes of the natural and forced jets are similar.

As we will see, SPOD and BMD modes of the twin-rectangular jet at the same frequency appear visually indistinguishable. Both methods successfully capture the most prevalent flow structures. We therefore defer a more detailed discussion of the physical mechanisms revealed by the modes until § 5.

4.2. Intermittency

Figure 6. SPOD-based time-frequency analysis of the (a–d) natural and (e–h) forced jets: (a,e) AS and (b, f) AA symmetry components; the (c,g) left and (d,h) right jets primarily located in the $z\gt 0$ and $z\lt 0$ half-domains, respectively. All panels share the same colour contours.

For the natural jet, the tonal peak at $f=0.29$ in the AS and AA symmetry components (see figure 4) signifies the prevalence of (double)-flapping structures in the flow, but only in a statistical sense. In figure 6, we investigate the temporal dynamics of these structures by conducting a time-frequency analysis of each symmetry. Specifically, we leverage the ability of the time-continuous SPOD expansion coefficients to provide insights into the global (as opposed to pointwise) evolution of modal structures (Schmidt et al. Reference Schmidt, Colonius and Brès2017a ; Towne & Liu Reference Towne and Liu2019). For computational efficiency, we determine the expansion coefficients

(4.4) \begin{equation} {\unicode{x1D63C}}_k^{\{i\}} = \boldsymbol{\varPhi }_k^*\unicode{x1D652}\hat {\unicode{x1D64C}}_k^{\left\{i-\tfrac {n_{\textit{fft}}}{2},\ldots ,i+\tfrac {n_{\textit{fft}}}{2}-1\right\}} \end{equation}

for the time index $i$ , with the Fourier transform $\hat {\unicode{x1D64C}}_k$ taken one snapshot at a time (Nekkanti & Schmidt Reference Nekkanti and Schmidt2021). The magnitudes of the complex-valued expansion coefficients are the time-dependent SPOD amplitudes. The coefficients corresponding to the leading mode are shown in figure 6(a,b) for the AS and AA symmetries, respectively, in the form of a time-frequency diagram. In figure 6(a), the AS symmetry displays a dark band at $f=0.29$ , as expected from the dominant screech tone. The persistence of the band with only small fluctuations suggests the strength of the AS screech tone remains relatively stable over time. The frequency of the tone also remains fixed. However, the AA screech tone shown in figure 6(b) is intermittent. While its frequency is steady, its strength fluctuates significantly such that it vanishes, and later re-emerges, several times in the data record. Compared with the AA screech, the dominance and steadiness of the AS screech suggest the absence of competition between the two screech symmetries that would have otherwise led to symmetry-switching. This is in accordance with observations by Wong et al. (Reference Wong, Stavropoulos, Beekman, Towne, Nogueira, Weightman and Edgington-Mitchell2023) that twin-round jets exhibit steady jet coupling at high NPR. The corresponding time-frequency diagrams for the forced jet are also shown in figure 6(e, f). As expected, the tonal peaks at frequency $f_0$ and its harmonics are clearly visible in the AS component and do not manifest intermittency. No peak is visible in the AA component.

For completeness, we also perform independent time-frequency analyses of the left and right jets. For the left jet, SPOD modes are computed from the data in the $z\gt 0$ half-domain, without recourse to $D_2$ symmetry decomposition. The time-continuous expansion coefficients are obtained from (4.4) as before. For the right jet, modes and coefficients are calculated in an analogous manner, just in the $z\lt 0$ half-plane. For the natural jet, the expansion coefficients of the left and right jets are reported in figure 6(c,d), respectively. The magnitudes of these coefficients are larger than those belonging to the AS and AA symmetries, because the left and right jets include all four symmetry components. Both the left and right jets share the same screech frequency. They appear intermittent and lack a clear phase relationship between each other. At the screech frequency, we have confirmed that the magnitudes of the expansion coefficients of the two jets exhibit low correlation with each other. The corresponding coefficients for the forced jet are displayed in figure 6(g,h). The harmonic peaks are visible and not intermittent, again as expected.

In light of the absence of symmetry-switching – as demonstrated by figure 6(a,b) – and to simplify our analysis, in the remainder of this work, we focus exclusively on results where $D_2$ symmetry is enforced.

5.1. Methodology

5.1.1. BMD with bicoherence normalisation

For a one-dimensional, statistically stationary random signal, $q(t)$ , and its Fourier transform, $\hat q(f)$ , the classical bispectrum is defined as the triple correlation

(5.1) \begin{equation} b(f_k,f_l) = E \{\hat q^*(f_k)\hat q^*(f_l)\hat q(f_k+f_l)\}, \end{equation}

where $ E \{\boldsymbol{\cdot }\}$ is the expectation operator. The frequency triplet, $(f_k,f_l,f_k+f_l)$ , constitutes a triad. BMD generalises the classical bispectrum to flow fields, $q(\boldsymbol{x},t)$ , by defining the spatially integrated bispectrum,

(5.2) \begin{equation} b(f_k,f_l) = E \left \{\int _\varOmega \hat q^*(\boldsymbol{x},f_k)\hat q^*(\boldsymbol{x},f_l)\hat q(\boldsymbol{x},f_k+f_l) \,\textrm{d}{\boldsymbol{x}}\right \} = E \left\{ {\left \langle \hat {\boldsymbol{q}}_{k+l},\hat {\boldsymbol{q}}_{k\circ l}\right \rangle }_{\boldsymbol{x}}\right\}\!, \end{equation}

where $\hat {\boldsymbol{q}}$ are the spatially resolved, temporal Fourier modes computed from the data, $\hat {\boldsymbol{q}}_{k\circ l}:=\hat {\boldsymbol{q}}_k\circ \hat {\boldsymbol{q}}_l$ , $\circ$ denotes a Hadamard product and $(\boldsymbol{\cdot })^*$ denotes the conjugate (transpose). The weighted inner product, $\langle {\boldsymbol{\cdot },\boldsymbol{\cdot }}\rangle _{\boldsymbol{x}}$ , is the same as that used for the SPOD in § 4. Our goal is to measure the expected quadratic phase coupling between the components $\hat {\boldsymbol{q}}_{k+l}$ and $\hat {\boldsymbol{q}}_{k\circ l}$ , independently of the power of each component. To this end, we normalise each component by the square root of its average integral power, $\sqrt { E \{{ \|\hat {\boldsymbol{q}}_{k+l} \|}^2_{\boldsymbol{x}}\}}$ and $\sqrt { E \{{ \|\hat {\boldsymbol{q}}_{k\circ l} \|}^2_{\boldsymbol{x}}\}}$ , respectively, with the norm defined as ${\|\boldsymbol{q}\|}^2_{\boldsymbol{x}}=\langle {\boldsymbol{q},\boldsymbol{q}}\rangle _{\boldsymbol{x}}$ . Thus normalised, the integrated bispectrum becomes equivalent to the integrated bicoherence,

(5.3) \begin{equation} b(f_k,f_l) = \frac { E \left\{ \left \langle {\hat {\boldsymbol{q}}_{k+l},\hat {\boldsymbol{q}}_{k\circ l}}\right \rangle _{\boldsymbol{x}}\right\}}{\sqrt { E \left\{{\left \|\hat {\boldsymbol{q}}_{k+l}\right \|}^2_{\boldsymbol{x}}\right\} E \left\{{\left \|\hat {\boldsymbol{q}}_{k\circ l}\right \|}^2_{\boldsymbol{x}}\right\}}}. \end{equation}

Like the classical bicoherence for one-dimensional signals, the magnitude of the integrated bicoherence is bounded. Specifically,

(5.4a) \begin{align} |b(f_k,f_l)| &= \frac {\left| E \left\{ \left \langle {\hat {\boldsymbol{q}}_{k+l},\hat {\boldsymbol{q}}_{k\circ l}}\right \rangle _{\boldsymbol{x}}\right\}\right|}{\sqrt { E \left\{{\left \|\hat {\boldsymbol{q}}_{k+l}\right \|}^2_{\boldsymbol{x}}\right\} E \left\{{\left \|\hat {\boldsymbol{q}}_{k\circ l}\right \|}^2_{\boldsymbol{x}}\right\}}} \end{align}

(triangle inequality) \begin{align} &\leqslant \frac { E \left\{\left| \left \langle {\hat {\boldsymbol{q}}_{k+l},\hat {\boldsymbol{q}}_{k\circ l}}\right \rangle _{\boldsymbol{x}}\right|\right\}}{\sqrt { E \left\{{\left \|\hat {\boldsymbol{q}}_{k+l}\right \|}^2_{\boldsymbol{x}}\right\} E \left\{{\left \|\hat {\boldsymbol{q}}_{k\circ l}\right \|}^2_{\boldsymbol{x}}\right\}}} \end{align}

(Cauchy-Schwarz inequality) \begin{align} &\leqslant \frac { E \left\{ {\left \|\hat {\boldsymbol{q}}_{k+l}\right \|}_{\boldsymbol{x}}{\left \|\hat {\boldsymbol{q}}_{k\circ l}\right \|}_{\boldsymbol{x}}\right\}}{\sqrt { E \left\{{\left \|\hat {\boldsymbol{q}}_{k+l}\right \|}^2_{\boldsymbol{x}}\right\} E \left\{{\left \|\hat {\boldsymbol{q}}_{k\circ l}\right \|}^2_{\boldsymbol{x}}\right\}}} \leqslant 1. \\[6pt] \nonumber \end{align}

In the following, we will assume the Fourier modes have been normalised to unit expected integral power.

For each triad, to estimate the BMD mode and mode bispectrum, $n_{\textit{blk}}$ independent realisations of the Fourier modes are assembled in the weighted data matrices,

(5.5a) \begin{align} \hat {\unicode{x1D64C}}_{k\circ l} &= \frac {\unicode{x1D652}^{1/2} \left[ \hat {\boldsymbol{q}}^{(1)}_{k\circ l}, \hat {\boldsymbol{q}}^{(2)}_{k\circ l}, \cdots , \hat {\boldsymbol{q}}^{(n_{\textit{blk}})}_{k\circ l} \right ]}{{\left \|\unicode{x1D652}^{1/2} \left[ \hat {\boldsymbol{q}}^{(1)}_{k\circ l}, \hat {\boldsymbol{q}}^{(2)}_{k\circ l}, \cdots , \hat {\boldsymbol{q}}^{(n_{\textit{blk}})}_{k\circ l} \right]\right \|}_F} \end{align}

(5.5b) \begin{align} \textrm {and} \,\,\hat {\unicode{x1D64C}}_{k+l} &= \frac {\unicode{x1D652}^{1/2} \left[ \hat {\boldsymbol{q}}^{(1)}_{k+l}, \hat {\boldsymbol{q}}^{(2)}_{k+l}, \cdots , \hat {\boldsymbol{q}}^{(n_{\textit{blk}})}_{k+l} \right ]}{\left \|{\unicode{x1D652}^{1/2} \left[ \hat {\boldsymbol{q}}^{(1)}_{k+l}, \hat {\boldsymbol{q}}^{(2)}_{k+l}, \cdots , \hat {\boldsymbol{q}}^{(n_{\textit{blk}})}_{k+l} \right ]}\right \|_F}, \end{align}

which have been normalised such that $\|\hat {\unicode{x1D64C}}_{k\circ l}\|_F = \|\hat {\unicode{x1D64C}}_{k+l}\|_F = 1$ , where ${\|\boldsymbol{\cdot }\|}_F$ is the Frobenius norm. In principle, a factor of $1/\sqrt {n_{\textit{blk}}}$ is present in the numerator and denominator of each $\hat {\unicode{x1D64C}}$ , which cancels out when $\hat {\unicode{x1D64C}}$ is normalised. For both $\hat {\unicode{x1D64C}}_{k\circ l}$ and $\hat {\unicode{x1D64C}}_{k+l}$ , by definition,

(5.6) \begin{equation} \big\|\hat {\unicode{x1D64C}}\big\|_F = \sqrt {\textrm{tr}(\hat {\unicode{x1D64C}}^*\hat {\unicode{x1D64C}})} = \frac {\sqrt {\sum _{i=1}^{n_{\textit{blk}}} \big(\hat {\boldsymbol{q}}^{(i)}\big)^* \unicode{x1D652} \hat {\boldsymbol{q}}^{(i)}}}{\big \|{\unicode{x1D652}^{1/2} \big[ \hat {\boldsymbol{q}}^{(1)}, \hat {\boldsymbol{q}}^{(2)}, \cdots , \hat {\boldsymbol{q}}^{(n_{\textit{blk}})}\big]}\big \|_F}, \end{equation}

where $\textrm{tr}(\boldsymbol{\cdot })$ denotes the trace. The Frobenius norm thus expresses the square root of the average power in $ \hat {\boldsymbol{q}}$ , which we normalise to one. The weight and spectral estimation parameters are identical to those used for SPOD (see table 3). We can define a single set of expansion coefficients, $\boldsymbol{a}_{k,l}$ , that relates the cross-frequency field, $\boldsymbol{\phi }_{k\circ l}=\unicode{x1D652}^{\kern1pt-1/2}\hat {\unicode{x1D64C}}_{k\circ l}\boldsymbol{a}_{k,l}$ , and bispectral mode, $\boldsymbol{\phi }_{k+l}=\unicode{x1D652}^{\kern1pt-1/2}\hat {\unicode{x1D64C}}_{k+l}\boldsymbol{a}_{k,l}$ , to the corresponding data matrices in (5.5). BMD then seeks the particular set of expansion coefficients that maximises the magnitude of the mode bispectrum estimated from $\boldsymbol{\phi }_{k\circ l}$ and $\boldsymbol{\phi }_{k+l}$ . That is,

(5.7) \begin{equation} \boldsymbol{a}_{k,l} = \mathop {\arg \max }\limits _{\|{\boldsymbol{a}}\|=1} | \boldsymbol{\phi }^*_{k\circ l}\unicode{x1D652}\boldsymbol{\phi }_{k+l} | = { \mathop {\arg \max }\limits _{\|{\boldsymbol{a}}\|=1} |\boldsymbol{a}^*\unicode{x1D63D}_{k,l}\boldsymbol{a}| }, \end{equation}

where $\unicode{x1D63D}_{k,l} = \hat {\unicode{x1D64C}}^*_{k\circ l}\hat {\unicode{x1D64C}}_{k+l}$ is the bispectral density matrix. The mode bispectrum is recovered as

(5.8) \begin{equation} { \beta _{k,l} = \boldsymbol{a}_{k,l}^*\unicode{x1D63D}_{k,l}\boldsymbol{a}_{k,l} \in \mathbb C. } \end{equation}

Equation (5.7) is the numerical radius problem for $\unicode{x1D63D}_{k,l}$ .

The magnitude of the BMD mode bispectrum is similarly bounded by unity, which can be proven as follows. The numerical radius of a matrix is bounded by its spectral norm (Goldberg & Tadmor Reference Goldberg and Tadmor1982; Horn & Johnson Reference Horn and Johnson1985),

(5.9) \begin{equation} |\beta (\unicode{x1D63D}_{k,l})| \leqslant \|\unicode{x1D63D}_{k,l}\|_2, \end{equation}

with the latter also equal to the leading singular value, $\sigma _1(\unicode{x1D63D}_{k,l})$ . To obtain an upper bound on $|\beta (\unicode{x1D63D}_{k,l})|$ , it thus suffices to consider $\|\unicode{x1D63D}_{k,l}\|_2$ . From the submultiplicativity of the spectral norm (Horn & Johnson Reference Horn and Johnson1985), it follows that

(5.10) \begin{equation} \|\unicode{x1D63D}_{k,l}\|_2 = \big\|\hat {\unicode{x1D64C}}^*_{k\circ l}\hat {\unicode{x1D64C}}_{k+l}\big\|_2 \leqslant \big\|\hat {\unicode{x1D64C}}^*_{k\circ l}\big\|_2 \big\|\hat {\unicode{x1D64C}}_{k+l}\big\|_2 = \big\|\hat {\unicode{x1D64C}}_{k\circ l}\big\|_2 \big\|\hat {\unicode{x1D64C}}_{k+l}\big\|_2. \end{equation}

Since the spectral norm is bounded by the Frobenius norm (Horn & Johnson Reference Horn and Johnson1985),

(5.11) \begin{equation} \big\|\hat {\unicode{x1D64C}}_{k\circ l}\big\|_2 \big\|\hat {\unicode{x1D64C}}_{k+l}\big\|_2 \leqslant \big\|\hat {\unicode{x1D64C}}_{k\circ l}\big\|_F \big\|\hat {\unicode{x1D64C}}_{k+l}\big\|_F = 1. \end{equation}

In other words, the mode bispectrum is bounded by the square root of the product of the average power in $\hat {\boldsymbol{q}}_{k\circ l}$ and $\hat {\boldsymbol{q}}_{k+l}$ , in this case normalised to one.

In the original implementation of BMD by Schmidt (Reference Schmidt2020), the numerical radius is calculated using the iterative algorithm of He & Watson (Reference He and Watson1997). The algorithm is not guaranteed to converge. In practice, non-convergence often leads to a mode bispectrum with a noisy appearance. In this work, we propose an improved implementation of BMD that uses the iterative algorithm of Mengi & Overton (Reference Mengi and Overton2005), which is guaranteed to converge to the global maximum. The clean appearance of the mode bispectra later reported is a direct consequence of the robustness of this new implementation.

In addition to considering frequency triads $(f_k,f_l,f_k+f_l)$ , we will also use BMD to investigate the nonlinear coupling among triadically compatible $D_2$ symmetry components. Using classical bispectral analysis, Walker & Thomas (Reference Walker and Thomas1997) previously uncovered experimental evidence for symmetry triads in supersonic single rectangular jets, with the caveat that symmetry was only assessed about one axis, which precludes the distinction between all four $D_2$ symmetry components. In the present work, $D_2$ symmetry triads are examined by assembling the Fourier modes $\hat {\boldsymbol{q}}_k$ , $\hat {\boldsymbol{q}}_l$ and $\hat {\boldsymbol{q}}_{k+l}$ from different, triadically compatible symmetry components. Ten of these non-redundant symmetry triads exist and are illustrated in figure 7.

Figure 7. Non-redundant $D_2$ symmetry triads, colour-coded by symmetry. The mode bispectra of the SS–SS interaction, (SS,SS,SS), AS–AS interaction, (AS,AS,SS), and SS–AS interaction, (SS,AS,AS), are shown in figure 9. The remaining triads are shown in figure 19.

5.1.2. Recovery of SPOD along the abscissa or ordinate

In the limit of infinitely long data and high frequency resolution, if the long-time mean is removed, the BMD mode bispectrum should vanish along $f_l=0$ , $f_k=0$ and $f_{k+l}=0$ . In practice, for limited data, the bispectrum tends to display finite values along these lines, indicative of unresolved low-frequency structures or slow trends captured by the zero-frequency bin of the DFT. Since the zero-frequency Fourier mode is real, it carries no phase information. For any triad that includes frequency zero, i.e. $(f,0,f)$ , $(0,f,f)$ or $(f,-f,0)$ , it can be shown that the triple correlation in the classical bispectrum defined in (5.1) may be rewritten as

(5.12) \begin{equation} \hat q^*(f)\hat q^*(0)\hat q(f) = \hat q^*(0)\hat q^*(f)\hat q(f) = \hat q^*(f)\hat q^*(-f)\hat q(0) = \hat q(0)|\hat q(f)|^2, \end{equation}

which is a real quantity. Here, we invoked the conjugate symmetry of the Fourier transform of real data. The spatial integral of this triple correlation, equivalent to the inner product $\langle {\hat {\boldsymbol{q}}_{k+l},\hat {\boldsymbol{q}}_{k\circ l}} \rangle _{\boldsymbol{x}}$ , is also real. The expectation, $ E \{ \langle {\hat {\boldsymbol{q}}_{k+l},\hat {\boldsymbol{q}}_{k\circ l}}\rangle _{\boldsymbol{x}}\}$ , which defines the integral bispectrum in (5.2), thus cannot measure the strength of phase coupling for such triads. Along these special lines, the classical and integral bispectra (and bicoherence) become phase-blind and are sensitive only to a weighted (by the zero-frequency Fourier mode) version of power.

This observation motivates an extension of BMD to recover the SPOD spectrum and modes along the abscissa, $f_l=0$ , or ordinate, $f_k=0$ . For flow data with a spatially uniform mean, $\bar q$ , recovery of the SPOD can be achieved simply by not subtracting the mean when solving the BMD problem. Neglecting block-to-block variations of the mean, which are likely small relative to $\bar q$ , the zero-frequency Fourier mode corresponds to the constant scalar mean, $\bar q$ . Along the abscissa, for $f_l=0$ , the normalised bispectral density matrix simplifies to

(5.13) \begin{equation} \unicode{x1D63D}_{k,0} = \hat {\unicode{x1D64C}}^*_{k\circ 0}\hat {\unicode{x1D64C}}_{k} = \hat {\unicode{x1D64C}}^*_k\hat {\unicode{x1D64C}}_k. \end{equation}

Because $\unicode{x1D63D}_{k,0}$ is now Hermitian and, by extension, normal, its numerical radius problem as given by (5.7) coincides with the method-of-snapshots SPOD eigenvalue problem in (4.3) (Goldberg & Tadmor Reference Goldberg and Tadmor1982; Horn & Johnson Reference Horn and Johnson1985). Along $f_l=0$ , BMD hence recovers the leading SPOD eigenvalues and modes of the normalised data. Similarly, along $f_k=0$ ,

(5.14) \begin{equation} \unicode{x1D63D}_{0,l} = \hat {\unicode{x1D64C}}^*_{0\circ l}\hat {\unicode{x1D64C}}_{l} = \hat {\unicode{x1D64C}}^*_l\hat {\unicode{x1D64C}}_l, \end{equation}

which recovers the same SPOD spectrum and modes.

No SPOD is recovered along $f_{k+l}=0$ , even with the uniform mean included. Again noting the conjugate symmetry of the Fourier transform of real data, the bispectral density matrix may be written as

(5.15) \begin{equation} \unicode{x1D63D}_{k,-k} = \hat {\unicode{x1D64C}}^*_{k\circ (-k)}\hat {\unicode{x1D64C}}_0 = \big(\big|\hat {\unicode{x1D64C}}_k\big|^2\big)^{\textrm{T}}\hat {\unicode{x1D64C}}_0, \end{equation}

where $|\boldsymbol{\cdot }|^2$ denotes the element-wise squared absolute value. This is the correlation between the mean and the squared magnitude of the Fourier mode $\hat {\boldsymbol{q}}_k$ , and is clearly not Hermitian, and thus does not correspond to an SPOD problem.

Figure 8. BMD mode bispectra of SS–SS interactions: (a) the long-time mean is removed from the data; (b) the long-time mean is included and SPOD is recovered on the $f_l=0$ axis (see Appendix B).

Whether SPOD is recovered on the abscissa or ordinate of BMD depends also on $D_2$ symmetry considerations. Because the mean flow has SS symmetry, SPOD can only be recovered from the symmetry triads that include SS as the first or second component. Specifically, among the non-redundant triads in figure 7, the BMD for (SS,SS,SS), (SS,SA,SA), (SS,AS,AS) and (SS,AA,AA) recover the SPOD for SS, SA, AS and AA, respectively. In the case of SS–SS interactions, both the abscissa and ordinate yield the SPOD. For SS–SA, SS–AS and SS–AA interactions, due to the (arbitrary) ordering of the symmetries, the ordinate yields the SPOD.

For inhomogeneous flows, the turbulent mean is unlikely to be truly uniform. However, for the present jet data, the mean pressure is nearly uniform (see figure 3). Even close to the nozzle, at $x=5$ , the mean pressure deviates locally from the ambient pressure, $p_\infty$ , by at most four and five percent for the natural and forced jets, respectively. The relationships in (5.13) and (5.14) are therefore satisfied approximately, as we demonstrate in detail in Appendix B. Figure 8 compares the BMD mode bispectra for SS–SS interactions when the mean is removed or included. A full discussion of the bispectrum is deferred until § 5.2. Here, we only remark that the bispectrum with the mean included shows higher magnitudes along $f_l=0$ and $f_{k+l}=0$ than the bispectrum without the mean. The mean-included bispectrum also displays local maxima along $f_l=0$ , consistent with the forced SPOD spectrum for SS in figure 4(a). The remainder of the bispectrum is independent of mean inclusion or removal. In what follows, we perform BMD exclusively with the mean included, combining spectral and bispectral information, that is, energetics and triadic coupling, in a single plot.

5.2. Notes on symmetries of the mode bispectrum under $D_2$ decomposition

Prior to interpreting the physics, we remark on some important symmetries of the mode bispectrum that the reader might observe. We note, however, that these bispectral symmetries are not salient in the physical interpretation. For discussions of the physics, the reader may safely skip to § 5.3.

Figure 9. BMD mode bispectra: (a) SS–SS interactions; (b) AS–AS interactions; (c) SS–AS interactions. All three bispectra share the same contour levels. Each bispectrum displays only its non-redundant region. The redundant regions can be recovered from the non-redundant regions via reflection (solid arrow, ) or reflection and complex conjugation (dotted arrow, ).

Of the ten possible $D_2$ symmetry triads illustrated in figure 7, only the (SS,SS,SS), (AS,AS,SS) and (SS,AS,AS) triads reveal clear signatures of active nonlinear interactions. The (SS,SS,SS) and (AS,AS,SS) triads involve interactions between coherent structures with the same symmetry, which we term symmetry-self interactions. Symmetry-self interactions, like AS–AS, yield a mode with SS symmetry. Conversely, the (SS,AS,AS) triad involves interactions between structures with different symmetries, which we term symmetry-cross interactions. The BMD mode bispectra of the SS–SS, AS–AS and SS–AS interactions are shown in figure 9(a,b,c), respectively. Only the non-redundant regions are shown (Schmidt Reference Schmidt2020). Each bispectrum displays a grid pattern of local maxima that indicate phase coupling between interacting waves. For the symmetry-self interactions in figure 9(a,b), the mode bispectrum is symmetric about the red line, $f_l=f_k$ , because the first two symmetry components can simply be swapped. The white region, $(f_k\gt 0,f_l\gt f_k)$ , can thus be recovered from a reflection of the non-redundant region about the red line, i.e.

(5.16) \begin{equation} \beta (f_k,f_l)=\beta (f_l,f_k)\quad \textrm {for}\,\, f_k\gt 0,\,f_l\gt f_k. \end{equation}

Owing to the conjugate symmetry of the Fourier transform of real data, i.e. $\hat {p}^*(f) = \hat {p}(-f)$ , the mode bispectrum is also symmetric about the blue line, $f_l=-f_k$ . The grey region, $(f_k\gt 0,f_l\lt -f_k)$ , can be recovered from a reflection of the non-redundant region about the blue line, followed by complex conjugation, i.e.

(5.17) \begin{equation} \beta (f_k,f_l)=\beta ^*(-f_l,-f_k) \quad \textrm {for } f_k\gt 0,\,f_l\lt -f_k. \end{equation}

Similarly, the symmetry-cross interactions in figure 9(c) enjoy conjugate symmetry. In this case, all of $f_l\gt -f_k$ is non-redundant. The grey region, $f_l\lt -f_k$ , can be recovered from the former via a reflection about $f_k=0$ , then another reflection about $f_l=0$ , followed by conjugation. Concisely,

(5.18) \begin{equation} \beta (f_k,f_l)=\beta ^*(-f_k,-f_l) \quad \textrm {for} \, f_l\lt -f_k. \end{equation}

Implicit in figure 9(c) is an additional symmetry between SS–AS and AS–SS interactions. It can be inferred from the form of the bispectral density that the mode bispectra of this pair of interactions are identical to each other under reflection about $f_l=f_k$ . In other words,

(5.19) \begin{equation} \beta _{{SS-AS}}(f_k,f_l)=\beta _{{AS-SS}}(f_l,f_k). \end{equation}

The subscript, $(\boldsymbol{\cdot })_{{SS-AS}}$ , denotes the spatial symmetry corresponding to each frequency component: SS and AS symmetries for the $f_k$ and $f_l$ components, respectively. The mode bispectral symmetries illustrated here also generalise to other spatial symmetry triads of the twin jet.

No bispectral symmetry exists for the SS–AS interaction in figure 9(c) about $f_l=f_k$ . That is,

(5.20) \begin{equation} \beta _{{SS-AS}}(f_k,f_l)\neq \beta _{{SS-AS}}(f_l,f_k). \end{equation}

In the non-redundant region, the strength of the $(f_k,f_l,f_{k+l})$ triad for the SS–AS interaction is, on average, higher for $f_l\gt f_k$ than for $f_l\lt f_k$ . Though this is not important to the present study, the likely explanation is that at the same harmonic frequency, the SS component is more energetic than the AS component. Given the identity (5.19), the statistical asymmetry expressed in inequality (5.20) can be restated as

(5.21) \begin{equation} \beta _{{SS-AS}}(f_k,f_l)\neq \beta _{{AS-SS}}(f_k,f_l). \end{equation}

Specifically, for $f_l\gt f_k$ , the $(f_k,f_l,f_{k+l})$ triad for the SS–AS interaction is more strongly coupled than the same frequency triad but for the AS–SS interaction. For $f_l\lt f_k$ , the reverse is true.

The relationship between the signs of $f_k$ and $f_l$ allows us to further classify triadic interactions as sum or difference interactions. Sum interactions are characterised by frequency doublets that satisfy $f_kf_l\gt 0$ . They couple primary waves at frequencies $f_k$ and $f_l$ to a secondary wave at a higher frequency, $f_{k+l}=f_k+f_l$ (Phillips Reference Phillips1960; Kim & Powers Reference Kim and Powers1979). Conversely, difference interactions, which satisfy $f_kf_l\lt 0$ , couple primary waves to a secondary wave at a lower frequency. A special case of difference interactions are mean flow deformations, which obey the condition that $f_k=-f_l$ . Throughout this work, we will sometimes refer to a triad whose secondary wave contributes to the primary wave of another triad as the precursor or progenitor of the latter. When using nomenclature like precursor or progenitor, we refer only to this primary–secondary relationship between the two triads, and do not imply a causal relationship.

5.3. Bispectral analysis

As alluded to in § 2.2.1, all harmonics of $f_0$ are directly excited by the rectangular-wave forcing. However, in a linear system, SS forcing only gives rise to harmonics in SS. The harmonic peaks in the AS SPOD spectrum (figure 4) are direct evidence that the forced jet is in the nonlinear regime. In the case of SS harmonics, direct excitation and nonlinear interactions play complementary roles. In the case of AS harmonics, nonlinear interactions are solely responsible. BMD enables us to confirm the presence, and quantify the significance, of these interactions.

The (SS,SS,SS) triad in figure 9(a) reveals both sum and difference interactions. The strongest nonlinear phase coupling identified is the difference interaction between SS modes at the fundamental and third harmonic frequencies, with a mode bispectrum magnitude of $|\beta _{{SS-SS}}(3f_0,-f_0)|=0.61$ . Overall, however, the sum and difference interactions are of comparable coupling strength. However, the (AS,AS,SS) triad in figure 9(b) is significantly biased towards difference interactions. No sum interaction is detected along $f_k=f_l$ . The most strongly coupled frequency triad is the difference interaction between AS modes at the fundamental and second harmonic frequencies, with a magnitude of $|\beta _{{AS-AS}}(2f_0,-f_0)|=0.29$ . The SS–AS mode bispectrum in figure 9(c) shows strong sum interactions in the $(f_k\gt 0,f_l\gt 0)$ quadrant as well as strong difference interactions in the $(f_k\lt 0,f_l\gt 0)$ quadrant, but weak difference interactions in the $(f_k\gt 0,f_l\lt 0)$ quadrant. The strongest triad is the sum interaction between an SS second harmonic mode and an AS fundamental mode, with a magnitude of $|\beta _{{SS-AS}}(2f_0,f_0)|=0.37$ . In a BMD mode bispectrum, sum interactions are an indicator of the classic frequency-doubling cascade that links the fundamental frequency component to successively higher harmonics. Conversely, difference interactions are often interpreted as the excitation of intermediate and lower frequencies (see e.g. Kim & Powers Reference Kim and Powers1979; Elgar & Guza Reference Elgar and Guza1985), though we reiterate that causal relationships cannot be deduced from BMD. As figure 9 shows, both (SS,SS,SS) and (SS,AS,AS) symmetry triads participate in such frequency-doubling cascades, whereas (AS,AS,SS) does not. This observation supports the hypothesis that while SS–SS and SS–AS interactions nonlinearly generate harmonics in SS and AS, respectively, AS–AS interactions do not generate harmonics in SS.

Figure 10. Dominant triads from figure 9: (a) SS–SS interactions; (b) AS–AS and SS–AS interactions. The SS–SS and AS–AS interactions, which couple to modes with SS symmetry, are represented by red spheres. The SS–AS interactions, which couple to modes with AS symmetry, are represented by green spheres.

In the following, we focus on symmetric mode couplings via the (SS,SS,SS) triad, and antisymmetric-symmetric mode couplings via the (AS,AS,SS) and (SS,AS,AS) triads. A representative set of dominant (SS,SS,SS) triads from figure 9(a) are reproduced in figure 10(a), which we will interpret in the following. The triadic network is formed not from a limited set of triadic cascades, but should rather be understood as the collective interaction of all temporal and spatial scales identified by the BMD bispectrum and modes. As such, we will not exhaustively catalogue all possible pathways. Instead, in § 5.4, we will focus on the most salient interactions that link the dominant triads in figure 10(a).

The (AS,AS,SS) triad cannot form a closed triad network, that is, the network cannot be visualised in a single bispectrum. This is because two primary waves with AS symmetry contribute to a secondary wave with SS symmetry. This secondary wave, however, cannot serve as the progenitor of another (AS,AS,SS) triad. Analogously, in an (SS,AS,AS) triad, primary waves with SS and AS symmetries contribute to a secondary wave with AS symmetry. This secondary wave cannot provide the SS component of another (SS,AS,AS) triad, so the network is again unclosed. The combination of the two symmetry triads, however, does form a closed network. SS–SS and SS–AS interactions contribute to secondary waves with both AS and SS symmetries. The secondary AS waves may in turn serve as the primary waves in another (AS,AS,SS) triad. Similarly, a pair of secondary AS and SS waves may participate in another (SS,AS,AS) triad. To elucidate these inter-triad relationships, we combine the dominant triads from the AS–AS and SS–AS bispectra in figure 9(b,c) into a three-dimensional bispectrum, shown in figure 10(b). The dominant (AS,AS,SS) triads are displayed on one pair of frequency axes, the dominant (SS,AS,AS) triads on another. The pathways between them will be examined in detail in § 5.5.

The remaining seven symmetry triads not included in figure 9 are documented in Appendix E. These triads involve the SA or AA symmetries, and do not manifest significant nonlinear interactions. Using BMD, we have also confirmed that the natural jet does not exhibit dominant triads. While nonlinear interactions do take place in the natural jet, they are significantly strengthened in the presence of periodic forcing, which provides a phase reference allowing for sustained phase coupling between symmetry and frequency components.

5.4. Symmetric mode interactions

Figure 11. Five representative triads from the (SS,SS,SS) mode bispectrum in figure 9(a). Large red spheres highlight the following triads: (a) $(f_0,f_0,2f_0)_{{SS-SS-SS}}$ ; (b) $(2f_0,2f_0,4f_0)_{{SS-SS-SS}}$ ; (c) $(2f_0,f_0, 3f_0)_{{SS-SS-SS}}$ ; (d) $(3f_0,-f_0,2f_0)_{{SS-SS-SS}}$ . Their possible precursor triads are marked by the small red spheres. All panels share the same axes. In each panel, the large sphere indicates the secondary wave at $f_{k+l}$ . Red solid ( ) and dashed ( ) lines distinguish between $f_k$ and $f_l$ , respectively. Dotted lines ( ) connect $f_k$ to $f_l$ . Arrows point towards or away from $f_{k+l}$ in sum or difference interactions, respectively. The translucent sphere denotes a complex conjugate.

Figure 12. Bispectral modes of the SS–SS interactions in figure 11. The corresponding interactions and modes are labelled with the same panel indices, (a–d). The $z=1.8$ plane is displayed. In this and the following figures of BMD modes, the colours saturate at $|\boldsymbol{\phi }|/\max |\boldsymbol{\phi }|=\pm 1$ . See supplementary movie 1 for an animation.

The grid of local maxima in figure 9(a) implies an intricate network of interconnected symmetry-self interactions between modes that possess SS symmetry. In figure 11, we illustrate them using five representative triads. The corresponding bispectral modes are shown in figure 12. For brevity, we display only the $z=1.8$ plane, which passes through the centreline of one nozzle and is parallel to the minor-axis plane. This choice is motivated by the observation that the SS and AS spatial symmetries, which are the only two components to engage in triadic interactions, are both symmetric about the minor-axis plane. They are thus most easily distinguished by their opposing symmetries about the major-axis plane, $y=0$ .

Figure 11(a) considers the $(f_0,f_0,2f_0)_{{SS-SS-SS}}$ triad, where a mode at the fundamental forcing frequency, $f_0$ , is coupled to a mode at the second harmonic, $2f_0$ . To each primary-wave frequency, $f_{k}$ or $f_{l}$ , a number of triads can contribute collectively, in particular, those indicated by the local maxima in the bispectrum along the diagonal of slope $-1$ associated with that $f_{k}$ or $f_{l}$ . In this case, the first harmonic primary mode of the $(f_0,f_0,2f_0)_{{SS-SS-SS}}$ triad coincides with the secondary modes of a number of other triads, e.g. $(f_0,0,f_0)_{{SS-SS-SS}}$ and $(2f_0,-f_0,f_0)_{{SS-SS-SS}}$ . Here, we highlight the pathway from $(f_0,0,f_0)_{{SS-SS-SS}}$ to $(f_0,f_0,2f_0)_{{SS-SS-SS}}$ . The $(f_0,0,f_0)_{{SS-SS-SS}}$ triplet recovers the leading SPOD mode and modal energy at $f_0$ . As shown in figure 12(a), excitations delivered by the plasma actuators to the initial shear layer are spatially amplified by the well-known Kelvin–Helmholtz (KH) instability. The bispectral mode captures symmetric KH-type wavepackets with extended spatial support in the streamwise direction. The mode shape bears qualitative resemblance to the symmetric modes extracted in the experiments of Samimy et al. (Reference Samimy, Webb, Esfahani and Leahy2023) from the same twin-rectangular jet, albeit overexpanded and forced in a different pattern. Trapped in the potential core are waves of a higher wavenumber. We have confirmed that these core modes are also present in the SA component. Due to the low supersonic jet velocity, the core modes are equivalent to the subsonic instability waves identified and extensively investigated by Tam & Hu (Reference Tam and Hu1989) and Towne et al. (Reference Towne, Cavalieri, Jordan, Colonius, Schmidt, Jaunet and Brès2017) in round jets. As they are transversely confined to the supersonic jet core in both the major- and minor-axis planes, the subsonic instability waves cannot propagate upstream (see supplementary movie 1 for an animation). Their maximum amplitude is found much closer to the nozzle than similar waves in supersonic, nearly shock-free round jets (Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018; Towne et al. Reference Towne, Schmidt and Brès2019). For shock-containing jets, these waves are documented in detail by Edgington-Mitchell et al. (Reference Edgington-Mitchell, Wang, Nogueira, Schmidt, Jaunet, Duke, Jordan and Towne2021). Evidence of these modes has previously surfaced in single-rectangular (Zaman et al. Reference Zaman, Fagan, Bridges and Brown2015; Gojon et al. Reference Gojon, Bogey and Marsden2016, Reference Gojon, Gutmark and Mihaescu2019; Semlitsch et al. Reference Semlitsch, Malla, Gutmark and Mihăescu2020; Tam & Chandramouli Reference Tam and Chandramouli2020; Zaman, Fagan & Upadhyay Reference Zaman, Fagan and Upadhyay2022; Ferreira et al. Reference Ferreira, Fiore, Parisot-Dupuis and Gojon2023; Wu, Lele & Jeun Reference Wu, Lele and Jeun2023; Karnam, Saleem & Gutmark Reference Karnam, Saleem and Gutmark2023) as well as twin-rectangular jets (Samimy et al. Reference Samimy, Webb, Esfahani and Leahy2023; Jeun et al. Reference Jeun, Wu and Lele2024). Here, the trapped modes attain an amplitude comparable to the KH-type instability waves, suggesting they too are indirectly energised by the symmetric forcing.

The fundamental forcing mode, $(f_0,0,f_0)_{{SS-SS-SS}}$ , is coupled to the second harmonic mode via the self interaction $(f_0,f_0,2f_0)_{{SS-SS-SS}}$ . The second harmonic mode, at twice the frequency of the fundamental, has approximately double the dominant streamwise wavenumber as well. This follows from the linear dispersion relation of KH-type waves in the initial shear-layer region of jets (Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018), and directly demonstrates that symmetry and frequency triads, which are enforced by our analysis, naturally form triads in wavenumber space. The wavepacket of the second harmonic mode is more compact, with its spatial support confined both axially and transversely, relative to the fundamental mode. Trapped modes are again observed. Compared with the KH-type waves, which are confined to $x\lesssim 5$ , the trapped waves appear to have more extended support in the axial direction. The presence of trapped waves in the bispectral mode $(f_0,f_0,2f_0)_{{SS-SS-SS}}$ suggests these waves are also quadratically phase-coupled with the fundamental mode.

In figures 11(b) and 12(b), the bispectral mode $(2f_0,2f_0,4f_0)_{{SS-SS-SS}}$ couples the second harmonic mode at frequency $2f_0$ to the fourth harmonic mode at frequency $4f_0$ through another frequency-doubling self interaction. Analogous to the relationship between the waveforms of the second harmonic and fundamental modes, the fourth harmonic mode also doubles the dominant streamwise wavenumber of the second harmonic mode. At the same time, the streamwise extent of the former is approximately halved, since only a thinner shear layer can support the higher-wavenumber KH-type wave. Similar trends are observed for the $(2f_0,f_0,3f_0)_{{SS-SS-SS}}$ triad in figures 11(c) and 12(c), where the second harmonic and fundamental modes undergo sum interaction, thereby becoming coupled to the third harmonic mode at frequency $3f_0$ . Together, the $(f_0,0,f_0)_{{SS-SS-SS}}$ , $(f_0,f_0,2f_0)_{{SS-SS-SS}}$ , $(2f_0,f_0,3f_0)_{{SS-SS-SS}}$ and $(2f_0,2f_0,4f_0)_{{SS-SS-SS}}$ triads form a hierarchy of coherent structures at four different frequencies, and exemplify a cascade of successively higher frequency and wavenumber components.

The $(3f_0,-f_0,2f_0)_{{SS-SS-SS}}$ triad in figure 11(d) indicates a difference interaction between the third harmonic and fundamental modes. One candidate for such an interaction is the coupling between the $(2f_0,f_0,3f_0)_{{SS-SS-SS}}$ and $(0,-f_0,-f_0)_{{SS-SS-SS}}$ bispectral modes. The latter is the complex conjugate of the $(f_0,0,f_0)_{{SS-SS-SS}}$ mode in the non-redundant region of the mode bispectrum. These three interacting modes are reported in figure 12(d). They illustrate the triadic coupling of a primary wave with high frequency and wavenumber, $(2f_0,f_0,3f_0)_{{SS-SS-SS}}$ , to a secondary wave with lower frequency and wavenumber, $(3f_0,-f_0,2f_0)_{{SS-SS-SS}}$ .

5.5. Antisymmetric-symmetric mode interactions

Figure 13. Same as figure 11, but for the (AS,AS,SS) (red) and (SS,AS,AS) (green) symmetry triads. The following frequency triads are highlighted by large spheres: (a) $(2f_0,-f_0,f_0)_{{AS-AS-SS}}$ ; (b) $(3f_0,-f_0, 2f_0)_{{AS-AS-SS}}$ ; (c) $(3f_0,-2f_0,f_0)_{{AS-AS-SS}}$ ; (d) $(-f_0,2f_0,f_0)_{{SS-AS-AS}}$ ; (e) $(f_0,f_0,2f_0)_{{SS-AS-AS}}$ ; ( f) $(2f_0,f_0,3f_0)_{{SS-AS-AS}}$ . Their possible precursor triads are exemplified by small red or green spheres. As in figure 11, solid ( ) and dashed () lines distinguish between $f_k$ and $f_l$ , respectively, while dotted lines ( ) connect $f_k$ to $f_l$ . Red ( ) and green ( ) spheres represent SS and AS symmetries, respectively.

Figure 14. Same as figure 12 but for AS–AS (left column) and SS–AS (right column) interactions. The panel indices, (a–f), correspond to the interactions in figure 13. See supplementary movie 2 for an animation.

In this section, we examine the network of spatial symmetry triads that is enabled by the wave coupling between modal structures with SS and AS symmetries. As in § 5.4, we make no attempt to enumerate all such couplings. Rather, we focus on a truncated triad network made up of the most dominant triads, previously shown in figure 10(b), and reproduced in figure 13. For each of these dominant triads, we investigate one pathway that could contribute to the triad.

As we observed in § 5.3, the (AS,AS,SS) symmetry triad mainly undergoes difference interactions. Figure 13(a) considers the coupling between an AS second harmonic mode, at frequency $2f_0$ , an AS fundamental mode, at frequency $-f_0$ , and an SS first harmonic mode, at frequency $f_0$ . This coupling is detected in the AS–AS mode bispectrum as the $(2f_0,-f_0,f_0)_{{AS-AS-SS}}$ frequency triad. One of its possible precursors is the SS–AS bispectral mode of the $(0,-f_0,-f_0)_{{SS-AS-AS}}$ triad, shown in figure 14(a) on the bottom left, which is the complex conjugate of the $(0,f_0,f_0)_{{SS-AS-AS}}$ triad, i.e. the leading AS SPOD mode at $f_0$ . The $(0,-f_0,-f_0)_{{SS-AS-AS}}$ triad has a secondary wave frequency of $f_{k+l}=-f_0$ , matching the primary wave frequency of the $(2f_0,-f_0,f_0)_{{AS-AS-SS}}$ triad, $f_l=-f_0$ . The $(0,-f_0,-f_0)_{{SS-AS-AS}}$ bispectral mode is dominated by a wavepacket structure located in the shear layer. The wavepacket flaps antisymmetrically about the major axis plane, $y=0$ . The envelope of the packet forms a standing wave in the near-jet (see supplementary movie 2 for an animation), and is produced by the interference between downstream-propagating KH-type waves and upstream-propagating guided-jet modes (G-JM) or free stream acoustic waves (Edgington-Mitchell et al. Reference Edgington-Mitchell, Li, Liu, He, Wong, Mackenzie and Nogueira2022). This behaviour is well known in screeching jets. It was first noted in the experiments by Davies & Oldfield (Reference Davies and Oldfield1962), Westley & Woolley (Reference Westley, Woolley and Ribner1968) and others, and investigated in depth by Panda (Reference Panda1999); see also Edgington-Mitchell (Reference Edgington-Mitchell2019) for a recent review. No standing waves are visible in SS or SA bispectral modes, which are symmetric about the major axis. Examples include the SS first and second harmonic modes on the right of figures 14(a) and 14(b), respectively. The absence of symmetric standing waves in twin-rectangular jets is consistent with previous observations (Jeun et al. Reference Jeun, Karnam, Wu, Lele, Baier and Gutmark2022; Samimy et al. Reference Samimy, Webb, Esfahani and Leahy2023).

The role of the complex conjugate in facilitating AS–AS difference interactions is also demonstrated by the $(3f_0,-f_0,2f_0)_{{AS-AS-SS}}$ triad in figure 13(b), with its corresponding modes in figure 14(b). It couples a pair of AS modes at the third and first harmonics, with frequencies $3f_0$ and $-f_0$ , respectively, to an SS mode at the second harmonic. One such example of a pair of AS primary modes consists of the SS–AS bispectral modes of the $(2f_0,f_0,3f_0)_{{SS-AS-AS}}$ and $(0,-f_0,-f_0)_{{SS-AS-AS}}$ triads. The latter is the conjugate of the $(0,f_0,f_0)_{{SS-AS-AS}}$ triad. Similarly, the $(3f_0,-2f_0,f_0)_{{AS-AS-SS}}$ triad in figure 13(c) couples third and second harmonic AS modes, with frequencies $3f_0$ and $-2f_0$ , respectively, to an SS mode at the first harmonic. An example pair of primary waves contributing to the $(3f_0,-2f_0,f_0)_{{AS-AS-SS}}$ triad are the SS–AS bispectral modes resulting from the $(2f_0,f_0,3f_0)_{{SS-AS-AS}}$ and $(-f_0,-f_0,-2f_0)_{{SS-AS-AS}}$ triads. The latter is the conjugate of the $(f_0,f_0,2f_0)_{{SS-AS-AS}}$ triad. The corresponding modes are shown in figure 14(c).

Whereas the triads reported in figure 13(a–c) involve symmetry-self interactions between AS mode pairs, the triads in figure 13(d–f) entail symmetry-cross interactions between SS and AS modes. Unlike AS–AS interactions, the most dominant SS–AS interactions include both sum and difference interactions. The $(-f_0,2f_0,f_0)_{{SS-AS-AS}}$ triad in figure 13(d) indicates a difference interaction between an SS first harmonic mode at frequency $-f_0$ and an AS second harmonic mode at frequency $2f_0$ . These primary waves couple to an AS secondary wave at the first harmonic, $f_0$ . The AS bispectral mode of this triad is displayed on the right of figure 14(d), and shows the KH and G-JM interference pattern discussed previously. The $(f_0,f_0,2f_0)_{{SS-AS-AS}}$ triad in figure 13(e) indicates the sum interaction between two first harmonic modes with SS and AS symmetries and an AS secondary mode at the second harmonic. The bispectral mode at the second harmonic is shown on the right of figure 14(e). It reveals both KH-type instability waves as well as core modes, but not clear signatures of G-JM. The AS core modes at the second harmonic share similar streamwise distribution as the SS core modes at the same frequency on the right of figure 14(b), but differ in their distribution in $y$ . Specifically, the AS core modes straddle the centreline and have double the number of anti-nodes in $y$ , as required by antisymmetry, compared with the SS core modes. In short, AS bispectral modes support both types of subsonic instability waves: those trapped in the core, and the G-JM that propagates in the shear layer and the ambient flow, with the preponderance of each type varying according to frequency. By contrast, SS bispectral modes appear to support only core modes. Sensitivity of the G-JM to $D_2$ symmetry has also been predicted using linear stability models for twin-round jets (Stavropoulos et al. Reference Stavropoulos, Mancinelli, Jordan, Jaunet, Weightman, Edgington-Mitchell and Nogueira2023), and generalises from circular symmetry for single-round jets (Tam & Hu Reference Tam and Hu1989). A subsonic jet mode can propagate upstream against the supersonic jet flow only if it has partial support outside the jet, in the slow ambient flow (Tam & Hu Reference Tam and Hu1989). Because of this, the difference in transverse spatial decay between subsonic modes with AS and SS symmetries, and between such modes at disparate frequencies, impacts the possibility and strength of feedback. This suggests an intimate link between the $D_2$ symmetry of the subsonic mode and its permitted direction of propagation. It forms part of the explanation for the symmetry-dependence of tones in the natural and forced twin-rectangular jets.

Analogous to the cascade of triads exemplified by the SS–SS interactions in figures 11 and 12, SS–AS interactions also facilitate the coupling between harmonic frequencies. In figure 13(e), the SS–AS mode bispectrum detects the phase coupling between SS and AS modes at the first harmonic, and an AS mode at the second harmonic, enabled by the $(f_0,f_0,2f_0)_{{SS-AS-AS}}$ triad. Similarly, figure 13( f) shows the $(2f_0,f_0,3f_0)_{{SS-AS-AS}}$ triad, which couples SS and AS modes at the second and first harmonics, respectively, to an AS mode at the third harmonic. The corresponding bispectral modes are reported in figures 14(e) and 14( f), respectively. They show the familiar doubling and tripling of axial wavenumbers expected from modes that span three successive frequency harmonics.

The AS–AS bispectral mode of the $(3f_0,-f_0,2f_0)_{{AS-AS-SS}}$ triad, shown on the right of figure 14(b), is indistinguishable from the SS–SS mode belonging to the same frequency triad on the right of figure 12(d). The equivalence between the SS–SS and AS–AS bispectral modes indicates that the symmetric and antisymmetric structures at frequencies $3f_0$ and $-f_0$ are phase-locked to the same symmetric structure at $2f_0$ . Although we have investigated the SS–SS bispectrum independently of the AS–AS and SS–AS bispectra, the close agreement between the SS–SS and AS–AS bispectral modes at the same $f_{k+l}$ confirms the notion that all three symmetry triads are in fact part of the same interconnected network of triad interactions.