3.1. Cyclostationary analysis

In this section, we briefly introduce the theory of cyclostationary analysis (Gardner Reference Gardner1988; Gardner et al. Reference Gardner, Napolitano and Paura2006; Antoni Reference Antoni2007; Napolitano Reference Napolitano2019). Let $\boldsymbol{q}({\boldsymbol{x}}, t)$ be a possibly complex, vector-valued random process at time $t$ over the set of independent variables $\boldsymbol{x}$ . Let the mean and correlation function be denoted as

(3.1) \begin{align} \widetilde {\boldsymbol{q}}({\boldsymbol{x}}, t) & \equiv E\{\boldsymbol{q}({\boldsymbol{x}}, t)\}, \end{align}

and

(3.2) \begin{align} {\boldsymbol{R}}\left({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, t, \tau \right) &\equiv E\left\{\boldsymbol{q}^{\prime \prime }({\boldsymbol{x}}, t+\tau /2)\boldsymbol{q}^{\prime \prime *}\left({\boldsymbol{x}}^\prime , t - \tau /2\right)\right\}\!, \end{align}

respectively, where $\tau$ represents the time shift and $(\boldsymbol{\cdot })^*$ is the Hermitian transpose, and where $\boldsymbol{q}^{\prime \prime }({\boldsymbol{x}}, t) \equiv \boldsymbol{q}({\boldsymbol{x}}, t) - \widetilde {\boldsymbol{q}}({\boldsymbol{x}}, t)$ are the fluctuations. The process is (wide-sense) cyclostationary if the mean and correlation function are periodic with period $T_0$ , i.e. $E\{\boldsymbol{q}({\boldsymbol{x}}, t)\} = E\{\boldsymbol{q}({\boldsymbol{x}}, t + T_0)\}$ and $ {\boldsymbol{R}}({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, t, \tau ) = {\boldsymbol{R}}({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, t + T_0, \tau )$ . Assuming cycloergodicity, defined as the equality of ensemble and phase-locked time averages for cyclostationary processes, the expectation can be replaced by the phase average

(3.3) \begin{align} \widetilde {\boldsymbol{q}}({\boldsymbol{x}}, t) = E\{\boldsymbol{q}({\boldsymbol{x}}, t)\} = \lim _{K \rightarrow \infty } \frac {1}{K} \sum _{k = 0}^K \boldsymbol{q}({\boldsymbol{x}}, t + kT_0). \end{align}

The above decomposition into cyclic mean and fluctuations is analogous to the triple decomposition (Hussain & Reynolds Reference Hussain and Reynolds1970, Reference Hussain and Reynolds1972), except that in that work the mean is further decomposed into its zero-frequency and non-zero-frequency components.

Since the correlation function is periodic in time, it can be expressed as a Fourier series

(3.4a) \begin{align} {\boldsymbol{R}}\left({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, t, \tau \right) &= \sum _{k_\alpha =-\infty }^\infty \hat {{\boldsymbol{R}}}_{k_\alpha \alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, \tau \right) e^{i2\pi k_\alpha \alpha _0 t}, \end{align}

where $T_0$ is the period, $\alpha _0$ is the fundamental frequency present in the system (typically the forcing frequency) and $\alpha = k_\alpha \alpha _0$ is the cycle frequency that represents how the statistics are modulated in time. The Fourier series coefficients $\hat {{\boldsymbol{R}}}_{k_\alpha \alpha _0}({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, \tau )$ are known as the cyclic correlation functions at cycle frequency $k_\alpha \alpha _0$ . The cyclostationary equivalent of the cross-spectral density (CSD) is the cyclic CSD (CCSD), which is equal to the Fourier transform of the cyclic correlation

(3.5) \begin{align} \boldsymbol{{S}}_{k_\alpha \alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, f\right) = \int _{-\infty }^\infty \hat {{\boldsymbol{R}}}_{k_\alpha \alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, \tau \right) e^{-i2\pi f\tau } \text{d}\tau , \end{align}

where $f$ is the spectral frequency. We note that the CSD tensor is the $k_\alpha =0$ component of the CCSD, $ \boldsymbol{{S}}_{0}({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, f)$ . Combining (3.5) and (3.4a ), and then taking the inverse Fourier series with respect to the cycle frequencies gives the time-periodic Wigner–Ville (WV) spectrum for each frequency $f$ (Martin Reference Martin1982; Martin & Flandrin Reference Martin and Flandrin1985; Antoni Reference Antoni2007)

(3.6) \begin{align} \boldsymbol{WV}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , t, f\right) = \sum _{k_\alpha =-\infty }^\infty \boldsymbol{{S}}_{k_\alpha \alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , f\right) e^{i2\pi k_\alpha \alpha _0 t}. \end{align}

For ${\boldsymbol{x}} = {\boldsymbol{x}}^\prime$ , the WV spectrum describes how the power-spectral density (PSD) varies as a function of time (or phase over a period). Numerically, the infinite summation is limited to $K_\alpha$ harmonics (i.e. $\sum _{k_\alpha = -{K_\alpha }}^{K_\alpha }$ ), where $K_\alpha$ is chosen such that all cycle frequencies of interest and sufficient energy are considered. The parameter $\boldsymbol{{S}}_{k_\alpha \alpha _0}({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, f)$ shows which cycle frequencies, $\alpha$ , are active for each spectral frequency, $f$ , and $\boldsymbol{WV}({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , t, f)$ shows us the evolution of the CSD (or PSD) of the signal as a function of the phase of the forcing. A simplified example of these tools applied to the harmonically forced Ginzburg–Landau system can be found in Heidt & Colonius (Reference Heidt and Colonius2024b ).

Integrating the CCSD and WV spectrum with respect to frequency (Gardner Reference Gardner1994; Randall et al. Reference Randall, Antoni and Chobsaard2001) results in the instantaneous variance and the cyclic distribution of the instantaneous variance, respectively,

(3.7a) \begin{align} {\boldsymbol{m}}(\boldsymbol{x}, t) &= E\left\{\boldsymbol{{q}}(\boldsymbol{x}, t)\boldsymbol{{q}}^*(\boldsymbol{x}, t)\right\} = \int _{-\infty }^\infty \boldsymbol{{WV}}(\boldsymbol{x}, t, f) \text{d}f, \end{align}

(3.7b) \begin{align} \hat {{\boldsymbol{m}}}_{k_\alpha \alpha _0}(\boldsymbol{x}) &= \int _{-\infty }^\infty \boldsymbol{{S}}_{k_\alpha \alpha _0}(\boldsymbol{x}, f) \text{d}f, \end{align}

where ${{\boldsymbol{m}}}(\boldsymbol{x}, t)$ is the mean variance and $\hat {{\boldsymbol{m}}}_{k_\alpha \alpha _0}(\boldsymbol{x})$ quantifies the mean-variance contribution from each cycle frequency $k_\alpha \alpha _0$ .

To estimate the CCSD in practice, after assuming cycloergodicity, we employ Welch’s method (Welch Reference Welch1967; Antoni Reference Antoni2007) as detailed in Heidt & Colonius (Reference Heidt and Colonius2024b ). We employ a block length of $N_{fft}$ such that $\Delta \textit{St} \approx 0.025$ for each case. This $\Delta St$ is chosen based on the suggestions by Heidt & Colonius (Reference Heidt and Colonius2024a ) where, for the turbulent jets investigated, they showed that $\Delta St \approx 0.025$ is a value that reduces variance while limiting the bias present at low frequencies $\textit{St} \rightarrow 0$ when using $\Delta St \approx 0.05$ as has been employed in previous studies, e.g. Schmidt et al. (Reference Schmidt, Towne, Rigas, Colonius and Brès2018) and Pickering et al. (Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020, Reference Pickering, Rigas, Schmidt, Sipp and Colonius2021a ). An overlap of $67\,\%$ is used to minimise spectral and cyclic leakage (Gardner Reference Gardner1986; Antoni Reference Antoni2007). Generally, the CCSD tensor may exist for any cycle frequency, i.e. $\boldsymbol{{S}}_{\alpha }({\boldsymbol{x}}, {\boldsymbol{x}}^{\prime }, f)$ . To verify that the flow is indeed cyclostationary, we scanned over all possible cycle frequencies, $\alpha$ , and confirmed that, for the forced jets, the only non-zero cycle frequencies occurs only at $\alpha = {z} \alpha _0 = {z} {\textit{St}}_{\!f}, z \in \mathbb{Z}$ (i.e. at the forcing frequency and its harmonics), meaning that the flow is indeed cyclostationary. To compute the WV spectrum (and other cyclostationary statistics throughout this paper) we consider cycle frequencies up to $\alpha = 4.5$ .

3.2. Cyclostationary spectral proper orthogonal decomposition

CS-SPOD is an extension of SPOD to flows with time-periodic statistics such as in harmonically forced turbulent flows (Heidt & Colonius Reference Heidt and Colonius2024b ). Like SPOD, CS-SPOD seeks modes $\boldsymbol{\phi }({\boldsymbol{x}}, t)$ as eigenfunctions of the two-point space-time correlation tensor, $\boldsymbol{{R}}({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , t, t^\prime ) \equiv E\{\boldsymbol{{q}}({\boldsymbol{x}}, t)\boldsymbol{{q}}^*({\boldsymbol{x}}^\prime , t^\prime )\}$ ,

(3.8) \begin{align} \int _{-\infty }^\infty \int _{\varOmega } \boldsymbol{{R}}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , t, t^\prime \right) \boldsymbol{{W}}\left({\boldsymbol{x}}^\prime \right) \boldsymbol{\phi }\left({\boldsymbol{x}}^\prime , t^\prime \right) \mathrm{d}\boldsymbol{x}^\prime \mathrm{d}t^\prime &= \lambda \boldsymbol{\phi }({\boldsymbol{x}}, t), \end{align}

where $\boldsymbol{{W}}({\boldsymbol{x}}^\prime )$ is the desired weight operator. Based on the time periodicity of the correlation tensor, the eigenvectors are assumed to contain a coupled set of spectral components denoted by $\mathcal{\varOmega }_{{\textit{St}}_c} = \{{\cdots} ,{\textit{St}}_c- 2\alpha _0,{\textit{St}}_c- \alpha _0,{\textit{St}}_c,{\textit{St}}_c+ \alpha _0,{\textit{St}}_c + 2\alpha _0,{\cdots} \}$ and oscillate as

(3.9) \begin{align} \boldsymbol{\phi }({\boldsymbol{x}}, t) = \sum _{k_{\!f} = -\infty }^\infty \boldsymbol{\psi }\left({\boldsymbol{x}}, {\textit{St}}_c + k_{\!f}\alpha _0\right) e^{i2\pi ({\textit{St}}_c + k_{\!f}\alpha _0)t}, \end{align}

where ${\textit{St}}_c$ is denoted as the centre frequency, defining which frequency set is being analysed. The associated eigenvalues, $\lambda$ optimally reconstruct the total variance (trace of $\boldsymbol{{R}}$ ). Numerically, $\mathcal{\varOmega }_{{\textit{St}}_c}$ is limited to ${K_{\!f}}$ harmonics such that $\mathcal{\varOmega }_{{\textit{St}}_c} = \{-{{K_{\!f}}}\alpha _0 + {\textit{St}}_c,(-{{K_{\!f}}}+1)\alpha _0 + {\textit{St}}_c,{\cdots} , {\textit{St}}_c, {\cdots} , ({{K_{\!f}}}-1)\alpha _0+ {\textit{St}}_c,{{K_{\!f}}}\alpha _0+ {\textit{St}}_c\}$ . The parameter ${K_{\!f}}$ should be chosen such that $\mathcal{\varOmega }_{{\textit{St}}_c}$ encompasses all frequencies of interest and sufficient energy.

Using the eigenfunction solution ansatz from (3.9), employing the definition of the cyclic cross-correlation functions and CCSD as in (3.5) and (3.4a ), and then performing a harmonic balance procedure as in Heidt & Colonius (Reference Heidt and Colonius2024b ) gives the final CS-SPOD eigenproblem, compactly written as

(3.10) \begin{align} \int _{\varOmega } {\unicode{x1D64E}}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , {\textit{St}}_c\right) {\unicode{x1D652}}\left({\boldsymbol{x}}^\prime \right) \boldsymbol{\varPsi }\left({\boldsymbol{x}}^\prime , {\textit{St}}_c\right) \mathrm{d}\boldsymbol{x}^\prime = \lambda \boldsymbol{\varPsi }({\boldsymbol{x}}, {\textit{St}}_c), \end{align}

where

(3.11a) \begin{align} &{\unicode{x1D64E}}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , {\textit{St}}_c\right) = \\ \nonumber &\begin{bmatrix} \ddots & \ddots & \ddots & \ddots & \ddots \\ \ddots & \boldsymbol{{S}}_0\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , {\textit{St}}_c - \alpha _0\right) & \boldsymbol{{S}}_{-\alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , {\textit{St}}_c - \frac {\alpha _0}{2} \right) & \boldsymbol{{S}}_{-2\alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime ,{\textit{St}}_c\right) & \ddots \\ \ddots & \boldsymbol{{S}}_{\alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime ,{\textit{St}}_c- \frac {\alpha _0}{2}\right) & \boldsymbol{{S}}_0\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime ,{\textit{St}}_c\right) & \boldsymbol{{S}}_{-\alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime ,{\textit{St}}_c+ \frac {\alpha _0}{2} \right) & \ddots \\ \ddots & \boldsymbol{{S}}_{2\alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime ,{\textit{St}}_c\right) & \boldsymbol{{S}}_{\alpha _0}\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime ,{\textit{St}}_c+ \frac {\alpha _0}{2} \right) & \boldsymbol{{S}}_0\left({\boldsymbol{x}}, {\boldsymbol{x}}^\prime ,{\textit{St}}_c+ \alpha _0\right) & \ddots \\ \ddots & \ddots & \ddots & \ddots & \ddots \end{bmatrix}, \end{align}

(3.11b) \begin{align} {\unicode{x1D652}}({\boldsymbol{x}}) = \begin{bmatrix} \ddots & & & & \\ & \hspace {-2mm}\boldsymbol{{W}}({\boldsymbol{x}}) & & & \\ & & \hspace {-2mm}\boldsymbol{{W}}({\boldsymbol{x}}) & & \\ & & & \hspace {-2mm}\boldsymbol{{W}}({\boldsymbol{x}}) & \\ & & & & \hspace {-2mm}\ddots \end{bmatrix}, \end{align}

(3.11c) \begin{align} \boldsymbol{\varPsi }({\boldsymbol{x}}, {\textit{St}}_c) = \begin{bmatrix} \vdots \\ \boldsymbol{\psi }({\boldsymbol{x}},{\textit{St}}_c- \alpha _{0})\\ \boldsymbol{\psi }({\boldsymbol{x}},{\textit{St}}_c)\\ \boldsymbol{\psi }({\boldsymbol{x}},{\textit{St}}_c+ \alpha _{0})\\ \vdots \end{bmatrix}. \end{align}

Here, ${\unicode{x1D64E}}({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , {\textit{St}}_c)$ is the CS-SPOD decomposition tensor, ${\unicode{x1D652}}({\boldsymbol{x}})$ is the concatenated weight tensor and $\boldsymbol{\varPsi }({\boldsymbol{x}}, {\textit{St}}_c)$ are the CS-SPOD eigenvectors in the frequency domain. The time-domain version of the CS-SPOD eigenvectors $\boldsymbol{\phi }({\boldsymbol{x}}, t)$ have Fourier series coefficients, at each $f \in \mathcal{\varOmega }_{{\textit{St}}_c}$ , of $\boldsymbol{\psi }({\boldsymbol{x}}, f)$ . This coupling of frequencies in CS-SPOD occurs because frequency components separated by $\mathbb{Z}\alpha _0$ are correlated to each other. Thus, each frequency component can not be solved individually. This is analogous to the frequency coupling in Floquet analysis or harmonic resolvent analysis (HRA).

Due to this frequency coupling, CS-SPOD performed at $f$ and $f + \alpha _0$ solve the same problem, i.e. giving $\mathcal{\varOmega }_{{\textit{St}}_c} = \mathcal{\varOmega }_{{\textit{St}}_c + z \alpha _0}$ , where $z \in \mathbb{Z}$ . Consequently, CS-SPOD only contains unique problems for the frequency sets corresponding to ${\textit{St}}_c = (-\alpha _0/2,\alpha _0/2]$ .

We note that the diagonal entries of ${\unicode{x1D64E}}({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , {\textit{St}}_c)$ are simply the standard CSD tensors as decomposed in SPOD. Thus, when analysing statistically stationary processes, which do not exhibit correlation between frequencies, all off-diagonal blocks of ${\unicode{x1D64E}}({\boldsymbol{x}}, {\boldsymbol{x}}^\prime , {\textit{St}}_c)$ become zero, and the CS-SPOD modes reduce to harmonic SPOD modes.

The properties of CS-SPOD are analogous to that of SPOD. Solving the eigendecomposition in (3.10) results in a countably infinite set of eigenfunctions $\boldsymbol{\varPsi }_{j}({\boldsymbol{x}}, {\textit{St}}_c)$ at ${\textit{St}}_c$ that are orthogonal to all other modes at ${\textit{St}}_c$ in the spatial inner norm $\langle \boldsymbol{{q}}_1, \boldsymbol{{q}}_2 \rangle _x = \int _{\varOmega } \boldsymbol{{q}}_2^*({\boldsymbol{x}}, t) \boldsymbol{{W}}({\boldsymbol{x}}) \boldsymbol{{q}}_1({\boldsymbol{x}}, t) \mathrm{d}\boldsymbol{x}$ , i.e. $\langle \boldsymbol{\varPsi }_{j}({\boldsymbol{x}}, {\textit{St}}_c), \boldsymbol{\varPsi }_{k}({\boldsymbol{x}}, {\textit{St}}_c) \rangle _x = \delta _{j,k}$ , are optimal in the amount of energy they capture, and form a complete basis of the flow (Heidt & Colonius Reference Heidt and Colonius2024b ). We employ the Chu compressible energy norm (Chu Reference Chu1965) used in previous investigations of jets (Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018), which maximises the sum of kinetic, thermal and potential energies. In addition, like previous investigations of these found jets (Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018; Pickering et al. Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020), we study the coherent structures a single azimuthal mode number, $m$ , at a time.

Heidt & Colonius (Reference Heidt and Colonius2024b ) presented an efficient algorithm for CS-SPOD based on Welch’s method and a frequency-shifting technique. The computational costs (and memory) are similar to standard SPOD. We employ a MATLAB implementation of the algorithms available at https://github.com/CyclostationarySPOD/CSSPOD (Heidt & Colonius Reference Heidt and Colonius2024b ). As before, $N_{fft}$ is selected such that $\Delta St \approx 0.025$ for each case, and an overlap of $67\,\%$ is employed to minimise spectral and cyclic leakage (Gardner Reference Gardner1986; Antoni Reference Antoni2007). To capture all relevant frequencies, we set ${K_{\!f}}$ such that $f_{\textit{max}} = 6$ , where $f_{\textit{max}}$ is the maximum frequency that CS-SPOD captures (i.e. $\text{max}\{|\mathcal{\varOmega }_{{\textit{St}}_c}|\}$ ).

3.3. Cyclostationary Navier–Stokes equations

To determine how the forcing interacts with the turbulence, we perform a cyclostationary decomposition of the governing equations. In particular, we wish to determine how the turbulence is modulated by (is reorganised by) the forcing and to what degree there is energy transfer between the mean flow and turbulent fields.

Thus, we decompose the governing equations in a similar manner to the standard Reynolds decomposition and the triple decomposition (Hussain & Reynolds Reference Hussain and Reynolds1970; Reynolds & Hussain Reference Reynolds and Hussain1972). However, unlike the triple decomposition, we split our state variables into two components; the first-order (deterministic) and second-order (stochastic) components. We do not decompose the first-order component (i.e. the phase average) into the period-averaged mean and phase-locked component since both components belong to the deterministic first-order field, and there is no advantage in separating them.

Let $\boldsymbol{f}(\boldsymbol{{x}}, t)$ be a vector-valued cyclostationary process, which we decompose into its first-order $\boldsymbol{\tilde {f}}(\boldsymbol{{x}}, t)$ (mean) and second-order $\boldsymbol{f^{\prime \prime }}(\boldsymbol{{x}}, t)$ (turbulent) components

(3.12) \begin{align} \boldsymbol{f}(\boldsymbol{{x}}, t) = \boldsymbol{\tilde {f}}(\boldsymbol{{x}}, t) + \boldsymbol{f^{\prime \prime }}(\boldsymbol{{x}}, t). \end{align}

The following three properties then arise from the decomposition:

(3.13a) \begin{align} \langle \boldsymbol{\tilde {f}}(\boldsymbol{{x}}, t) \rangle &= \boldsymbol{\tilde {f}}(\boldsymbol{{x}}, t), \end{align}

(3.13b) \begin{align} \langle \boldsymbol{\tilde {f}}(\boldsymbol{{x}}, t)\boldsymbol{{g}}(\boldsymbol{{x}}, t) \rangle &= \boldsymbol{\tilde {f}}(\boldsymbol{{x}}, t) \boldsymbol{\tilde {g}}(\boldsymbol{{x}}, t), \end{align}

(3.13c) \begin{align} \langle \boldsymbol{f^{\prime \prime }}(\boldsymbol{{x}}, t) \rangle &= \langle \boldsymbol{f}(\boldsymbol{{x}}, t) - \boldsymbol{\tilde {f}}(\boldsymbol{{x}}, t) \rangle = 0, \end{align}

where $\boldsymbol{\tilde {f}}(\boldsymbol{{x}}, t) = \langle \boldsymbol{f}(\boldsymbol{{x}}, t) \rangle$ and $\boldsymbol{g}(\boldsymbol{{x}}, t)$ is another vector-valued cyclostationary process. Combining the second and third properties gives

(3.14) \begin{align} \langle \boldsymbol{\tilde {f}}(\boldsymbol{{x}}, t)\boldsymbol{{f^{\prime \prime }}}(\boldsymbol{{x}}, t)\rangle = 0 , \end{align}

i.e. the mean flow and the turbulence are linearly uncorrelated (but not independent).

We now use this decomposition to understand the interaction between the forcing and the turbulence of periodically forced turbulent flows. For convenience, we use the incompressible form (in the low subsonic jets ( $M_{\!j}=0.4$ ) we consider, the density fluctuations can be expected to be minor in the energy balance), written in Cartesian notation

(3.15a,b) \begin{align} \frac {\partial {u_i}}{\partial {x_i}} = 0,\quad \frac {\partial {u_i}}{\partial {t}} + u_{\!j}\frac {\partial {u_i}}{\partial {x_{\!j}}} = -\frac {1}{\rho }\frac {\partial {p}}{\partial {x_i}}+ \nu \frac {\partial ^2u_i}{\partial x_{\!j} \partial x_{\!j}}, \end{align}

where $u_i, x_i, t, \rho , \nu$ are the velocity component in the $i$ th dimension, the spatial component in the $i$ th dimension, time, density and kinematic viscosity, respectively. Substituting (3.12) into the continuity and momentum equations and taking the mean, we obtain

(3.16a,b) \begin{align} \frac {\partial {\tilde {u}_i}}{\partial {x_i}} = \frac {\partial {u_i^{\prime \prime }}}{\partial {x_i}} = 0, \quad \frac {\partial {\tilde {u}_i }}{\partial {t}} + \frac {\partial {\tilde {u}_i \tilde {u}_{\!j}}}{\partial {x_{\!j}}} + \frac {\partial {\widetilde {u_i^{\prime \prime }u_{\!j}^{\prime \prime }}}}{\partial {x_{\!j}}} = -\frac {1}{\rho }\frac {\partial {\tilde {p}}}{\partial {x_i}} + \nu \frac {\partial ^2\tilde {u}_i }{\partial x_{\!j} \partial x_{\!j}}. \end{align}

To investigate the energy transfer between the mean and turbulent components, we determine the mean kinetic energy equation by multiplying the mean momentum equation by $\tilde {u}_i$ and then taking the mean. After manipulating, we obtain

(3.17) \begin{align} \frac {\tilde {D}\left(\frac {1}{2}\tilde {u}_i \tilde {u}_i \right)}{Dt} &= - \frac {\partial {\tilde {u}_i \widetilde {u_i^{\prime \prime }u_{\!j}^{\prime \prime }}}}{\partial {x_{\!j}}} - \left (- \widetilde {u_i^{\prime \prime }u_{\!j}^{\prime \prime }}\frac {\partial {\tilde {u}_i}}{\partial {x_{\!j}}} \right ) -\frac {1}{\rho } \frac {\partial {\tilde {p} \tilde {u}_i}}{\partial {x_i}} \nonumber \\&\quad + \nu \frac {\partial }{\partial x_{\!j}}\left [\tilde {u}_i \left ( \frac {\partial \tilde {u}_i}{\partial x_{\!j}} + \frac {\partial \tilde {u}_{\!j}}{\partial x_i}\right ) \right ] - \frac {\nu }{2} \left ( \frac {\partial \tilde {u}_i}{\partial x_{\!j}} + \frac {\partial \tilde {u}_{\!j}}{\partial x_i}\right ) \left ( \frac {\partial \tilde {u}_i}{\partial x_{\!j}} + \frac {\partial \tilde {u}_{\!j}}{\partial x_i}\right )\!, \end{align}

where the first term is accumulation + transport, the second, fourth and fifth terms are diffusion, the sixth is dissipation and the third is transfer. The production term

(3.18) \begin{align} - \widetilde {u_i^{\prime \prime }u_{\!j}^{\prime \prime }} \frac {\partial {\tilde {u}_i}}{\partial {x_{\!j}}}, \end{align}

is of particular interest as it represents the drain of energy from the mean to the turbulent field by the action of the mean against the turbulent Reynolds stresses. This term is a source in the turbulent kinetic energy (TKE) equation and a sink in the mean kinetic energy equation. Unlike the Reynolds/triple decomposition, we obtain a phase-dependent energy transfer that allows us to investigate how the energy transfer between the mean and turbulent fields varies under the external forcing.

4.1. Mean flow

We begin with an analysis of the impact of the forcing on the mean flow. It is known that high-Reynolds-number turbulent jets are receptive to a wide range of forcing frequencies and azimuthal mode numbers (Crow & Champagne Reference Crow and Champagne1971; Samimy et al. Reference Samimy, Kim, Kastner, Adamovich and Utkin2007, Reference Samimy, Kim, Kearney-Fischer and Sinha2010, Reference Samimy, Kearney-Fischer, Kim and Sinha2012; Sinha et al. Reference Sinha, Towne, Colonius, Schlinker, Reba, Simonich and Shannon2018). Generally, forcing results in a thickening of the shear layer and a consequent reduction in the period-averaged potential-core length, which can be seen in figures 3 and 4.

Figure 3.

Contours of the (a) mean axial velocity $\widetilde {u}_x(\phi _{\!f})$ at $\phi _{\!f} = 0$ and (b) the $f = 0$ component (period averaged) of the axial velocity $\widetilde {u}_{x, 0}$ . The solid and dashed white lines correspond to lines of $\widetilde {u}_x(\phi _{\!f}) = 0.25, 0.75$ and $\widetilde {u}_{x, 0} = 0.25, 0.75$ for (a) and (b), respectively.

Figure 4.

(a) Period-averaged axial velocity $\widetilde {u}_{x,0}$ , (b) standard deviation of the mean axial velocity $\sigma (\widetilde {u}_x)$ along the centreline and lipline and (c) momentum thickness $\delta _\theta$ .

Figure 4. Long description

Panel A: A line graph shows the period-averaged axial velocity. The horizontal axis is labeled ‘x’ and the vertical axis is labeled ‘ux,0’. Multiple lines represent different different forced jets.The graph shows the values along the centerline and lipline. Panel B: A line graph displays the standard deviation of the mean axial velocity along the centerline and lipline. The horizontal axis is labeled ‘x’ and the vertical axis is labeled ‘σ(ux)’. Different lines correspond to various Strouhal numbers and amplitudes. Panel C: A line graph illustrates the momentum thickness. The horizontal axis is labeled ‘x’ and the vertical axis is labeled ‘δθ’. Multiple lines represent different forced jets.

Figure 3 displays the mean axial velocity. The left column shows the mean flow at a forcing phase $\phi _{\!f} = 0^\circ$ . For ${\textit{St}}_{\!f} = 0.3$ , the mean shows a substantial response to forcing, and is modulated by the roll-up of large-scale vortices in the potential-core region whose strength increases as the amplitude of forcing is increased. In contrast, for ${\textit{St}}_{\!f} = 1.5$ , the effect of the forcing is spatially limited to the near-nozzle region, highlighted in the two zoomed-in contours in figure 3, with the largest phase dependency occurring in the first diameter downstream of the nozzle exit. As we will show in figure 5(a), the ${\textit{St}}_{\!f}=1.5$ flow also features phase-locked roll-up of vortices in the near-nozzle region. These observations, namely the large impact of the ${\textit{St}}_{\!f}=0.3$ forcing on the mean flow, are in accord with past studies of forced jets and the identification of the largest response to forcing (jet preferred mode) at ${\textit{St}}_{\!f} \approx 0.3$ . The largest response of the jet to forcing occurring at ${\textit{St}}_{\!f} \approx 0.3$ is also supported by linear stability analysis (Garnaud et al. Reference Garnaud, Lesshafft, Schmid and Huerre2013). For the ${\textit{St}}_{\!f} = 1.5,\ a_0/U_{\!j} = 10\,\%$ forcing, the high-amplitude acoustic waves have steepened as they propagate downstream and begin to exhibit shock-like structures resulting in the observed structure outside the jet.

Figure 5.

(a) Root-mean-square contours $\sqrt {(\widetilde {u_x^{\prime \prime } u_x^{\prime \prime }})}$ at $\phi _{\!f} = 0$ , (b) period-averaged r.m.s. $\sqrt {(\widetilde {u_x^{\prime \prime } u_x^{\prime \prime }})_0}$ along the centreline and lipline and (c) r.m.s. contours $\sqrt {(\widetilde {u_x^{\prime \prime } u_x^{\prime \prime }})}$ at $\phi _{\!f} = 0,\ \pi /2,\ \pi ,\ 3\pi /2$ for ${\textit{St}}_{\!f}=0.3, a_0/U_{\!j}=10\,\%$ .

Figure 3(b) shows the zero-frequency component of the mean axial velocity (i.e. the period-averaged flow). Despite the jet’s strong harmonic response to the forcing, the impact on the period-averaged field is less obvious. This is further examined in figure 4, where we display the period-averaged axial velocity $\widetilde {u}_{x,0}(x)$ and the standard deviation of the mean axial velocity $\sigma (\widetilde {u}_x(x))$ along the centreline and lipline, along with the momentum thickness $\delta _\theta$ . Similar to previous studies (Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018), $\delta _\theta$ is defined as

(4.1) \begin{align} \delta _\theta (x) = \int _0^{r_{0.05}}\frac {\widetilde {\rho }_0(x, r)\widetilde {u}_{x,0}(x, r)}{\widetilde {\rho }_0(x, 0)\widetilde {u}_{x,0}(x, 0)} \left ( 1 - \frac {\widetilde {u}_{x,0}(x, r)}{\widetilde {u}_{x,0}(x, 0)}\right ) \mathrm{d}r, \end{align}

where $\widetilde {u}_{x,0}(x, r)$ is the time- and azimuthally averaged axial velocity. The integral radial bound $r_{0.05}$ corresponds to the radial location where $\widetilde {u}_{x,0}(x, r_{0.05}) - U_\infty = 0.05 \ \widetilde {u}_{x,0}(x, 0)$ . An increase in the momentum thickness with forcing amplitude and an associated shortening of the period-averaged potential-core length is evident for all cases. The ${\textit{St}}_{\!f} = 0.3,a_0/U_{\!j}=10 \,\%$ forcing is sufficiently strong to result in a non-monotonically decreasing period-averaged centreline velocity. This is also seen close to the nozzle for the ${\textit{St}}_{\!f} = 1.5,a_0/U_{\!j}=10\,\%$ forcing case. The ${\textit{St}}_{\!f}=0.3$ forcing at an amplitude $a_0/U_{\!j} = 1\,\%$ and $10\,\%$ results in maximum standard deviations (as a ratio of the jet velocity) of $15\,\%$ and $33\,\%$ , respectively, along the centreline, and $9.8\,\%$ and $35\,\%$ , respectively, along the lipline. In contrast, the ${\textit{St}}_{\!f} = 1.5$ forcing is more localised to the near-nozzle exit, with substantial lipline phase dependency occurring within $0.15$ diameters of the nozzle exit with a maximum standard deviation of $7.0\,\%$ and $18\,\%$ for the $a_0/U_{\!j} = 1\,\%$ and $10\,\%$ forcing, respectively. The effect along the centreline is reduced, with a maximum phase dependency of just $3.5\,\%$ for the $a_0/U_{\!j} = 10\,\%$ forcing. Curiously, the ${\textit{St}}_{\!f} = 1.5,a_0/U_{\!j}=10\,\%$ forcing results in a slight decrease in mixing far downstream, where the centreline velocity and momentum thickness are slightly larger and smaller, respectively, than the natural jet at $x = 20$ .

4.2. Second-order statistics

In figure 5(a), we display the azimuthally averaged root-mean-square (r.m.s.) turbulent axial velocity $\widetilde {u^{\prime \prime }_x u^{\prime \prime }_x}$ at $\phi _{\!f} = 0$ . The turbulence is modulated and is phase dependent in similar regions as was seen in the mean flow – up to the end of the potential core for ${\textit{St}}_{\!f}=0.3$ and near the nozzle exit (up to $x \approx 2$ ) for ${\textit{St}}_{\!f} = 1.5$ . Both $a_0/U_{\!j}=10\,\%$ forcing cases exhibit clear vortex roll-up. For all cases, the radial location of the maximum phase dependency of the r.m.s. is along the lipline, $r \approx 0.5$ , which is where the centre of the rolled-up vortices pass.

We display the period-averaged axial velocity r.m.s. along the centreline and lipline in figure 5(b). The ${\textit{St}}_{\!f} = 0.3$ forcing results in a minor upstream shift of the centreline r.m.s., which is expected given the observed shortening of the jet. Additionally, an increase in the peak centreline and lipline r.m.s. is observed, indicating an increase in the period-averaged intensity of the turbulence. A more substantial effect on the lipline r.m.s. is seen, with a decrease for the first $5$ diameters and then a substantial increase between $x \in [6, 11]$ diameters downstream. This decrease in the period-averaged r.m.s. is due to the shear layer being modulated towards and away from the centreline due to the forcing, thereby resulting in a lower period-averaged value. The substantial increase in the r.m.s. at $x \in [6, 11]$ is likely caused by the highly energetic vortices generated by the forcing. These vortices advect downstream while breaking down, thereby creating a substantial amount of turbulent energy. For the ${\textit{St}}_{\!f} = 1.5,a_0/U_{\!j}=10\,\%$ forcing, the r.m.s. contours clearly show the presence of vortex roll-up. However, in contrast to the ${\textit{St}}_{\!f} = 0.3$ case, these vortices are smaller in size and are spatially limited to $x \in [0, 2]$ . Consequently, these vortices break down more rapidly and do not increase peak r.m.s. along the centreline. Instead, a slight decrease in the peak centreline r.m.s. is observed (and a slight upstream shift for the $a_0/U_{\!j}=10\,\%$ forcing). Along the lipline, a large increase is observed for $x \in [0,2]$ corresponding to where the ${\textit{St}}_{\!f}=1.5$ forcing is dominant. After $x = 5$ , a slight decrease is observed in the lipline r.m.s. for the forced cases.

In figure 6, we display the azimuthally averaged WV spectrum of the axial velocity at two axial stations, $x=1.25$ and $x=5$ , along the lipline. This shows how the frequency content of $u^{\prime \prime }_x$ varies as a function of the phase of the forcing. These stations were chosen to correspond to the locations in figure 5(a) where the r.m.s. is most phase-dependent for the ${\textit{St}}_{\!f} = 1.5$ and $0.3$ forcing, respectively. The WV spectrum shows that the $a_0/U_{\!j} = 10\,\%$ forcing cases exhibit a phase dependency at both locations, while the $a_0/U_{\!j} = 1\,\%$ cases show only a mild phase dependency that is limited to $x = 1.25$ for the ${\textit{St}}_{\!f} = 1.5$ forcing and $x = 5$ for the ${\textit{St}}_{\!f} = 0.3$ forcing. Looking at the mean axial velocity in figure 3 and the r.m.s. of the turbulent axial velocity as a function of phase in figure 5(c) for the ${\textit{St}}_{\!f} = 0.3,\ a_0/U_{\!j} = 10\,\%$ case at $x = 5, r = 0.5$ , we find that the high-velocity phases correspond to bursts of higher-frequency energy in the WV spectrum associated with the passage of the rolling-up vortex.

Figure 6.

WV spectrum of $u_x^{\prime \prime }$ at (a) $x = 1.25,\ r = 0.5$ and (b) $x = 5,\ r = 0.5$ .

Figure 6. Long description

A contour map of the WV spectrum of the turbulent axial velocity across different conditions. The heat map is organized into a grid layout with two rows labeled (a) and (b), and five columns labeled Natural, Stf = 0.3, a0/Uj = 1 percent, Stf = 0.3, a0/Uj = 10 percent, Stf = 1.5, a0/Uj = 1 percent, and Stf = 1.5, a0/Uj = 10 percent. The x-axis is labeled with phi_f ranging from 0 to 3π/2, and the y-axis is labeled with St ranging from 0 to 1.5. The color scale on the right indicates the magnitude of the data values, with colors ranging from blue (low values) to red (high values). Each panel shows variations in data distribution under different conditions, with notable patterns and clusters visible in specific regions.

For reference, we also display in figure 6 the WV spectrum of the natural jet computed by assuming a forcing frequency of ${\textit{St}}_{\!f}=0.3$ . As the natural jet is stationary, the minor phase dependency of the spectrum with the phase of forcing is an artefact of non-convergence of the statistics. Comparing these spectra with the forced cases allows one to discern a satisfactory statistical convergence – where a phase dependence is absent from the natural jet.

To further isolate where the turbulence is most strongly modulated by forcing, we display contours of the cyclic variance of the turbulent axial velocity as a function of $\alpha /\alpha _0$ given by (3.7b ) in figure 7. This is also equal to the Fourier series expansion of the phase-dependent mean-squared variance. The period-averaged variance, $\alpha /\alpha _0 = 0$ , shows the change in the r.m.s. Next, $\alpha /\alpha _0=1$ and $2$ shows where the turbulence is phase-dependent with a period equal to half, and a quarter, of the forcing frequency, respectively. This shows where the turbulence is being modified by the forcing (i.e. where nonlinearity has resulted in phase-dependent turbulence). Thus, non-zero values of $\hat {{\boldsymbol{m}}}_{\alpha }(x, r) = 0$ for $\alpha /\alpha _0 \ne 0$ identify phase-dependent second-order statistics, whereas their absence indicates that this variance measure is statistically stationary. Thus, for example, we see that despite the extremely strong ${\textit{St}}_{\!f}=1.5,\ a_0/U_{\!j}=10\,\%$ forcing, the flow is statistically stationary for $x \gt 2$ . The period-averaged variance, $\alpha /\alpha _0 = 0$ , when compared with the natural jet, is elevated towards the end of the potential core for the ${\textit{St}}_{\!f}=0.3$ forcing and in the near-nozzle shear layer region for the ${\textit{St}}_{\!f}=1.5$ forcing. We also see that the phase dependency of the variance $\alpha /\alpha _0=1,\ 2$ is more global for the ${\textit{St}}_{\!f}=0.3$ forcing, located $x \in [0, 10]$ , than the ${\textit{St}}_{\!f}=1.5$ forcing which is much more localised to the near-nozzle region $x \in [0, 2]$ .

Figure 7.

Cyclic variance $\hat {{\boldsymbol{m}}}_{\alpha }(x, r)$ contours of the turbulent axial velocity $u^{\prime \prime }_x(x, r)$ for (a) $\alpha /\alpha _0 = 0$ , (b) $\alpha /\alpha _0 = 1$ and (c) $\alpha /\alpha _0 = 2$ .

4.3. Energy transfer

We quantify the interaction between the mean flow and the turbulence by computing the energy transfer (production of TKE)

(4.2) \begin{align} \widetilde {E}(x, r, \theta , \phi _{\!f}) = -\widetilde {u_{i}^{\prime \prime } u_{j}^{\prime \prime }} \frac {\partial \widetilde {u}_{i}}{\partial x_{\!j}}, \end{align}

where, for ease of notation, we use Cartesian coordinates on the right-hand side and cylindrical on the left. We note that any phase-averaged quantities are axisymmetric due to the homogeneity in $\theta$ of the flow and the forcing. In figure 8 we display the energy transfer at three phases of the forcing $\phi _{\!f} = 0,\pi /3$ and $2\pi /3$ . In line with previous results, the ${\textit{St}}_{\!f} = 0.3$ forcing results in a maximum phase dependency of the energy transfer towards the end of the potential core while the ${\textit{St}}_{\!f} = 1.5$ cases are localised to $x \in [0,2]$ . Similar to the r.m.s. contours, the phase dependency of the energy transfer for the $a_0/U_{\!j}=1\,\%$ forcing cases is relatively minor. The phase dependency is substantial for the $a_0/U_{\!j}=10\,\%$ forcing cases, and the localisation of the energy transfer to the vortex cores is clearly observed.

Figure 8.

Energy transfer contours ${\widetilde {E}(x, r, \theta , \phi _{\!f}})$ at (a) $\phi _{\!f} = 0,\ {(b) }\, \phi _{\!f} =2\pi /3$ and $\ {(c) }\, \phi _{\!f} =4\pi /3$ .

Figure 8. Long description

Panel A: A heat map showing energy transfer contours for natural conditions and various St_f and a_0/U_{j} values. The x-axis represents distance (x) and the y-axis represents radius (r). The color scale ranges from 0 to 0.05, with higher values indicated by lighter colors. The heat map shows distinct patterns of energy transfer, with notable variations in intensity across different conditions. Panels b and c show the evolution of the energy transfer across the phase of the forcing.

Next, we analyse the global impact that the forcing has on the energy transfer. In figure 9, we display cumulative energy transfers by first integrating $\widetilde {E}$ in $r$ and $\theta$

(4.3) \begin{align} \widetilde {E}_r(x, \phi _{\!f})=\int _{r} \int _{\theta } \widetilde {E}\big(x, r, \theta , \phi _{\!f}\big) r \partial r \partial \theta , \end{align}

and then, additionally, in $x$

(4.4) \begin{align} \widetilde {E}_{rx}(x, \phi _{\!f}) =\int _{x^\prime = 0}^x \int _{r} \int _{\theta } \widetilde {E}\left(x^\prime , r, \theta , \phi _{\!f}\right) r \partial r \partial \theta \partial x^\prime , \end{align}

and we refer to these as the local and global energy transfer, respectively. In the local measure, we observe a strongly phase-dependent energy transfer for $x \in [0, 10]$ for the ${\textit{St}}_{\!f} = 0.3$ cases, where bands of low and high energy transfer correspond to the roll-up and passage of the large-scale vortices. The width and energy of the bands increase as the vortices grow in size and strength and then, after $x \approx 8$ , begin to decay. A similar phenomenon occurs for the ${\textit{St}}_{\!f} = 1.5$ cases but is spatially restricted to $x = [0, 2]$ , in line with our mean-flow observations. After this, another phase-dependent energy transfer phenomenon becomes apparent with a wavelength of $\approx 1.66$ , which corresponds to the acoustic wavelength associated with the forcing. This is an acoustic consequence of the strong high-frequency forcing and steepening of the plane acoustic waves outside the shear layer (and also evident in figure 3(a)), which apparently interacts with the natural stochastic perturbations present in the jet.

Figure 9.

Contours of (a) $\widetilde {E}_r(x, \phi _{\!f})$ and (b) $\widetilde {E}_{rx}(x, \phi _{\!f})$ as a function of the phase of the forcing $\phi _{\!f}$ and axial location $x$ . The white lines overlaid on the $\widetilde {E}_{rx}(x, \phi _{\!f})$ contours represent lines of $[2,\ 10,\ 50,\ 90]\,\%\ \times 0.25$ energy transfer.

Despite the strong local phase dependency, the global energy transfer is similar across all cases, as seen in the $\widetilde {E}_{rx}(x, \phi _{\!f})$ contours. The white lines represent locations of fixed energy transfer which are phase-independent for the natural case and nearly so with $a_0/U_{\!j}=1\,\%$ forcing. To obtain a modest phase-dependent integrated energy transfer, $a_0/U_{\!j}=10\,\%$ forcing must be used. In figure 10, we quantify this effect by displaying the period-average and minimum/maximum values of $\widetilde {E}_{rx}(x, \phi _{\!f})$ across a forcing period. A similar result is seen: only the $10\,\%$ forcing results in an appreciable phase dependency of the total energy transfer. Even with the $10\,\%$ forcing, the variation as a function of phase is minor after $x \approx 5$ and 1 for ${\textit{St}}_{\!f}=0.3$ and $1.5$ , respectively. This absence of a substantial phase dependence indicates that at a given phase, the higher energy transfer linked to the high-velocity regions is counteracted by the lower energy transfer linked to the low-velocity regions resulting in relatively constant integrated energy transfer despite the large local variation. After $x= 2,0.7$ for the ${\textit{St}}_{\!f}=0.3,1.5$ forcing, the change in the period-averaged energy transfer for the forced jets is large in comparison with the maximum phase-dependent variation. This suggests that the global energy transfer in the jet is not a strong function of the phase of the forcing. Combined with the relatively modest phase-dependent variation of the second-order statistics, this is consistent with the interpretation that forcing influences the turbulence primarily through modifications to the period-averaged mean flow. This appears consistent with the hypothesis (Kœnig et al. Reference Kœnig, Sasaki, Cavalieri, Jordan and Gervais2016; Sinha et al. Reference Sinha, Towne, Colonius, Schlinker, Reba, Simonich and Shannon2018) that the primary change in jet noise is driven by a corresponding change to the period-averaged mean.

Figure 10.

Comparison of $\widetilde {E}_{rx, 0}(x)$ and the maximum fluctuation of $\widetilde {E}_{rx}(x, \phi _{\!f})$ .

Figure 10. Long description

A line graph displays integrated energy across different cases with various parameters. The horizontal axis is labeled ‘x’ and ranges from 10^−1 to 10^1. The vertical axis is labeled with a logarithmic scale ranging from 10^−3 to 10^0. The graph includes multiple lines representing different cases: Natural, St_f = 0.3 with a_0/U_{j} = 1 percent, St_f = 0.3 with a_0/U_{j} = 10 percent, St_f = 1.5 with a_0/U_{j} = 1 percent, and St_f = 1.5 with a_0/U_{j} = 10 percent. Each line shows the trend of energy transfer for these cases. The black line represents the natural case, the blue line represents St_f = 0.3 with a_0/U_{j} = 1 percent, the red line represents St_f = 0.3 with a_0/U_{j} = 10 percent, the purple line represents St_f = 1.5 with a_0/U_{j} = 1 percent, and the green line represents St_f = 1.5 with a_0/U_{j} = 10 percent. The graph also includes dashed lines indicating the minimum and maximum fluctuation of energy transfer.

5.1. CS-SPOD

The results of § 4 showed that forcing generates a strong deterministic response and substantial local phase-dependent modulation of the turbulence, while producing comparatively modest changes in the integrated energy transfer. We now address the central question of this study: How does forcing modify the coherent structures of the non-phase-locked turbulent field? To answer this question, we employ CS-SPOD, which isolates the phase-dependent second-order statistics and enables changes to the turbulent coherent structures to be quantified independently of the deterministic forced response. Similar to previous studies on round jets, we decompose the jets by azimuthal mode number $m$ and study the coherent structures at each $m$ .

Interpreting coherent structures in cyclostationary flows is not as straightforward as in the stationary case, where each frequency is independent. Instead, each CS-SPOD mode oscillates with a set of frequencies $\mathcal{\varOmega }_{{\textit{St}}_c} = \{{\cdots} ,{\textit{St}}_c- 2\alpha _0,{\textit{St}}_c- \alpha _0,{\textit{St}}_c,{\textit{St}}_c+ \alpha _0,{\textit{St}}_c + 2\alpha _0,{\cdots} \}$ , where ${\textit{St}}_c$ is termed the centre frequency. For every ${\textit{St}}_c$ , the amplitude (and energy) of the mode varies periodically in time with period $1/\alpha _0$ . The CS-SPOD eigenvalues are periodic in ${\textit{St}}_c$ with unique values existing in $(-\alpha _0/2,\alpha _0/2]$ . The SPOD modes can be thought of as the limit of CS-SPOD modes as $\alpha _0 \rightarrow 0$ , in which case the amplitude (energy) becomes constant and the set of frequencies at which the real and imaginary modes oscillate collapse to the centre frequency, $\textit{St} = {\textit{St}}_c$ , but the range of ${\textit{St}}_c$ now extends to $[-\infty ,\infty ]$ . Statistically (and strictly for $m=0$ ), the SPOD and CS-SPOD eigenvalues are the same at $\pm {\textit{St}}_c$ , so we only need to examine the range ${\textit{St}}_c \in (0,\alpha _0/2]$ .

Figure 11(a) shows the first two eigenvalues of the CS-SPOD spectrum for the forced jets at $m=0$ and ${\textit{St}}_{\!f}=0.3$ (with $\alpha _0=0.3$ ). The overall trend for ${\textit{St}}_c \ne 0$ is a decrease in energy with increasing ${\textit{St}}_c$ , which can be associated with the typical behaviour of jet turbulence with monotonically decreasing energy with increasing $\textit{St}$ . The energy of the modes increases with forcing amplitude across all frequencies, with the largest increase as ${\textit{St}}_c$ approaches $\pm \alpha _0/2$ . Increasing the forcing amplitude also increases the energy separation between the dominant and subdominant modes, and this effect is also stronger as $\rightarrow \pm \alpha _0/2$ .

Figure 11.

Comparison of the $m = 0$ CS-SPOD eigenspectrum of the natural and (a) ${\textit{St}}_{\!f}=0.3$ and (b) ${\textit{St}}_{\!f}=1.5$ forced jets.

Figure 11(b) shows the results for the ${\textit{St}}_{\!f}=1.5$ forcing. A different behaviour is observed, in that the dominant mode is relatively unaffected by the forcing level for lower ${\textit{St}}_c$ . Then at higher ${\textit{St}}_c$ the energy decreases with increasing forcing level until ticking up again as ${\textit{St}}_c$ approaches $\alpha _0/2$ . The first subdominant modes are unaffected by the forcing for all ${\textit{St}}_c$ . This reveals that the physical mechanism associated with the dominant coherent structures is strengthened (thereby enhancing the low-rank behaviour of the jet) by the forcing at some ${\textit{St}}_c$ , and decreases the low-rank behaviour of the jet at other ${\textit{St}}_c$ .

In figure 11, we also superpose the CS-SPOD spectra of the natural jets. Using CS-SPOD for the natural jet is permissible as stationarity is a subset of cyclostationarity. In fact, up to statistical convergence, we show in Appendix A that SPOD and CS-SPOD produce the same spectrum for the natural jet, regardless of the chosen value of $\alpha _0$ . However, the ordering of the modes is different for each choice of $\alpha _0$ (and for regular SPOD), and this makes a side-by-side comparison with the forced jet challenging. We find that the simplest interpretation of the eigenvalues follows by comparing the natural jet with the forced jets using CS-SPOD with $\alpha _0$ chosen to be the corresponding forcing frequency. We are essentially interpreting the natural jet as forcing with vanishing strength. We see that for the most part, the natural jet fits the trend-wise descriptions for each case given above.

Kœnig et al. (Reference Kœnig, Sasaki, Cavalieri, Jordan and Gervais2016) and Sinha et al. (Reference Sinha, Towne, Colonius, Schlinker, Reba, Simonich and Shannon2018) showed that high-frequency forcing results in a more stable mean flow that reduces the amplification of coherent structures, while the opposite is true for low-frequency forcing. The CS-SPOD results agree with this interpretation, and provide some further clarity. For both forcing frequencies, we find that the spectrum is more impacted as ${\textit{St}}_c \rightarrow \alpha _0/2$ and that the spectrum is affected less as ${\textit{St}}_c \rightarrow 0$ . As we will show shortly, we find that coherent structures at $\textit{St} \lt \lt {\textit{St}}_{\!f}$ (i.e. the components present for ${\textit{St}}_c \rightarrow 0$ ) evolved independently to the phase-dependent mean field. This occurs because the dominant coherent structures at $\textit{St}$ have a wavelength much greater than the wavelength of the forcing. Instead, these structures are modified by corresponding changes to the period-averaged mean flow. Since the period-averaged mean flow is impacted substantially less by the forcing than the phase-averaged mean field, as seen in figures 3 and 4, these coherent structures are impacted less. For structures $\textit{St} \approx {\textit{St}}_{\!f}/2$ (i.e. the low-frequency components present for ${\textit{St}}_c \rightarrow \alpha _0/2$ ), the structures respond more to the phase-dependent mean. Since the phase-dependent variation of the mean is far greater than the period-averaged mean flow, as shown in figures 3 and 4, the impact on the energy/structure of the coherent structures is also far greater.

Samimy et al. (Reference Samimy, Kim, Kastner, Adamovich and Utkin2007) showed that low-frequency forcing ( ${\textit{St}}_{\!f} \approx 0.3$ ) results in a louder jet while high-frequency forcing ( ${\textit{St}}_{\!f} \gt 1$ ) results in a small noise reduction. While the stability of the mean flow is not directly computed in their study, the increase/decrease in the noise is consistent with the increase/decrease in the stability of the period-averaged mean flow found by Kœnig et al. (Reference Kœnig, Sasaki, Cavalieri, Jordan and Gervais2016) and Sinha et al. (Reference Sinha, Towne, Colonius, Schlinker, Reba, Simonich and Shannon2018). Here, we provide further clarity. In Samimy et al. (Reference Samimy, Kim, Kastner, Adamovich and Utkin2007), low-frequency forcing ( ${\textit{St}}_{\!f} = 0.36, m = 0$ ) at $\theta = 30^\circ$ results in a noise increase from $\textit{St} \approx 0.1$ onward with the noise increase being larger at higher $\textit{St}$ . This increase in the high-frequency noise can be explained using Lighthill acoustic analogy (Lighthill Reference Lighthill1952), which states that total sound power is proportional to the eighth power of the jet velocity, $\propto U_{\!j}^8$ . Due to a large difference in the frequency between the low-frequency forcing and the high-frequency wavepackets, the high-frequency wavepackets are likely to respond to the change in the mean flow in a quasi-steady manner. This was seen in the WV spectrum in figure 6, where it was seen that the high-frequency components responded in a quasi-steady manner to the change in the mean. Furthermore, this is similar to the results by Franceschini et al. (Reference Franceschini, Sipp, Marquet, Moulin and Dandois2022), who showed that the high-frequency Kelvin–Helmholtz structures present in the turbulent flow over a square cylinder evolved in a quasi-steady manner, and the structures were phase locked to the vortex shedding. Since low-frequency forcing results in a substantial modulation of the jet velocity (with alternating regions of high and low velocity), the high-frequency wavepackets are more energetic (louder) during the high-velocity phases and less energetic (quieter) during the low-velocity phases. Since acoustic power is highly nonlinear $\propto U_{\!j}^8$ , the increase in acoustic power during the high-velocity phases overwhelms the decrease in acoustic power during low-velocity phases, thereby resulting in a louder jet at high frequencies on average. The same trend is seen for high-frequency ${\textit{St}}_{\!f} = 1.07$ forcing in Samimy et al. (Reference Samimy, Kim, Kastner, Adamovich and Utkin2007), where the forcing still results in an increase in higher-frequency noise $\textit{St} \gt 1$ (where for $\textit{St} \gt \gt {\textit{St}}_{\!f}$ the structures can again be considered to respond to the periodic mean flow in a quasi-steady manner). This is consistent with our quasi-steady structure hypothesis.

To investigate how the forcing impacts the coherent structures, we display the dominant CS-SPOD mode for $m=0$ and ${\textit{St}}_{\!f}=0.3$ jet at ${\textit{St}}_c = 0.1$ in figure 12. As described previously, CS-SPOD modes have a periodic amplitude in time, so we display the modes at several time instances. We also display the corresponding phase of the forcing and overlay contours of $\widetilde {u}_x(\phi _{\!f}) = 0.25$ and $0.75$ to demonstrate how the coherent structures are linked to/modulated by the phase-dependent mean flow. By contrast with the eigenvalues, we compare CS-SPOD eigenvectors with the SPOD (rather than CS-SPOD) eigenvector of the natural jet, which has constant amplitude in time. We select the corresponding SPOD mode by carefully sorting the spectra as discussed in Appendix C. We find that this provides a cleaner comparison, as the CS-SPOD eigenvectors for the natural jet show statistical artefacts that are easy to misinterpret as cyclostationary features.

Figure 12.

Pressure of the dominant $m=0$ mode of (a) the natural jet at $\textit{St}=0.1$ and (b, c) the ${\textit{St}}_{\!f}=0.3$ forced jets (1 % and 10 % forcing, respectively) at ${\textit{St}}_c=0.1$ . The contour limits are $\pm 0.75 |\boldsymbol{\phi }(t)|_\infty$ . The solid and dashed black lines correspond to lines of $\widetilde {u}_x(\phi _{\!f}) = 0.25, 0.75$ , respectively.

Figure 12. Long description

Panel A: A heatmap shows the pressure of the dominant mode in a natural jet at a specific condition. The x-axis is labeled ‘x’ and the y-axis is labeled ‘r’. The contour limits are set to specific values. Solid and dashed black lines represent contours lines of 25% and 75% of the mean axial velocity. Panel B: Three heatmaps display the pressure of the dominant mode in forced jets with 1 percent forcing at different times. The x-axis is labeled ‘x’ and the y-axis is labeled ‘r’. The contour limits are set to specific values. Solid and dashed black lines represent contours lines of 25% and 75% of the mean axial velocity. Panel C: Three heatmaps show the pressure of the dominant mode in forced jets with 10 percent forcing at different times. The x-axis is labeled ‘x’ and the y-axis is labeled ‘r’. The contour limits are set to specific values. Solid and dashed black lines represent contours lines of 25% and 75% of the mean axial velocity.

Looking at the CS-SPOD modes, we find a few salient features. First, all modes are Orr-type modes with a dominant frequency component at $\textit{St} = 0.1$ , which corresponds to the most energetic frequency component in $\mathcal{\varOmega }_{{\textit{St}}_c}$ . The Orr-type mode structure is consistent with the literature for this $m-St$ pair (Pickering et al. Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020); forcing has not changed the dominant mechanism. Second, we see an upstream shift in the location of the Orr mode corresponding to the shift in the period-averaged field. Third, as the forcing amplitude is increased, we observed a strong phase dependency of the dominant Orr-type mode where the highest-amplitude region of the mode follows the region of the greatest velocity/shear present in the jet. This phase dependency of the dominant mode is most clearly observed in figure 12.

In figure 13, we decompose the CS-SPOD modes into their Fourier components, $\boldsymbol{\psi }({\boldsymbol{x}}, f)$ of the CS-SPOD modes $\boldsymbol{\phi }({\boldsymbol{x}}, t)$ displayed in figure 12. We can clearly observe the features seen in the time-domain modes, such as the contraction of the Orr-type mode as the forcing amplitude increases. Furthermore, we note that for the high-amplitude forcing, there is a substantial change in the structure of the $\textit{St} = 0.1$ component as well as significant contributions from the $\textit{St} = -0.5, -0.2, 0.4$ and $0.7$ components. These frequency components are all coupled together to create the time-domain modes seen in figure 12.We note that these $\textit{St} = -0.5, 0.4$ and $0.7$ components do not have a Kelvin–Helmholtz-type structure as would be expected of the dominant SPOD mode at these frequencies, as seen in Schmidt et al. (Reference Schmidt, Towne, Rigas, Colonius and Brès2018) and Pickering et al. (Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020). Thus, the dominant CS-SPOD mode is not simply the summation of multiple dominant SPOD modes. Instead, it forms an entirely new mode.

Figure 13.

Fourier coefficients of the CS-SPOD modes depicted in figure 12; layout the same as figure 12. Contour limits are $\pm 0.75 |\boldsymbol{\psi }({\boldsymbol{x}}, f)|_\infty$ .

In figure 14, we display the dominant CS-SPOD mode for the natural jet and the ${\textit{St}}_{\!f}=0.3$ forced cases at two other centre frequencies, ${\textit{St}}_c = 0.05$ and $0.15$ . For the forced jets, the time-dependent modes are shown at a single instance in time corresponding to the phase $\phi _{\!f} = 0$ . Similar to the previous results at ${\textit{St}}_c = 0.1$ , both the natural and forced jets show a dominant structure consistent with that of an Orr mode consistent with previous results in the literature (Pickering et al. Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020). We find that for both ${\textit{St}}_c$ , as the forcing level is increased there is a reduction in the length of the dominant frequency component and stronger phase-dependent structures.

Figure 14.

Dominant $m=0$ SPOD/CS-SPOD mode of the natural jet and the ${\textit{St}}_{\!f}=0.3$ forced jets at (a) ${\textit{St}}_c=0.05$ and (b) ${\textit{St}}_c=0.15$ . The pressure is shown at $\phi _{\!f} = 0$ for the forced jets. Contour limits are $\pm 0.75 |\boldsymbol{\phi }(\phi _{\!f} = 0)|_\infty$ . The solid and dashed black lines correspond to lines of $\widetilde {u}_x(\phi _{\!f} = 0) = 0.25, 0.75$ , respectively. Note the varying axis limits between different sub-figures.

Figure 14. Long description

Six contour maps in two panels of three. Panel (a) is at a centre frequency of 0.05 and panel (b) at 0.15. Within each panel the maps are stacked vertically: the natural jet on top, then the jet forced at a Strouhal number of 0.3 at 1 percent amplitude, then the same forcing at 10 percent. All six show the pressure of the dominant axisymmetric SPOD mode for the natural jet and CS-SPOD mode for the forced jets, the latter at a forcing phase of zero. The horizontal axis is the axial coordinate and the vertical axis the radial coordinate; axis limits differ between panels, the axial extent reaching about 28 in (a) and 16 in (b). Pressure uses a diverging red and blue scale with symmetric limits at plus and minus 0.75 of the peak mode amplitude, so colour shows the shape of the mode rather than absolute magnitude. Black lines are contours of the mean axial velocity at the same phase: solid for 0.25, the outer edge of the shear layer, and dashed for 0.75, the inner edge.

In figure 15, we compare the dominant modes for the natural jet and the ${\textit{St}}_{\!f}=1.5$ forced jets at ${\textit{St}}_c = 0.05,0.15,0.475$ and $0.75$ . Based on the separation between the centre frequency and $\alpha _0={\textit{St}}_{\!f} = 1.5$ , we expect the dominant mechanisms to be those associated with the natural jet at the centre frequency. Based on Pickering et al. (Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020), we expect ${\textit{St}}_c=0.05$ and ${\textit{St}}_c=0.15$ to correspond to Orr mechanism and ${\textit{St}}_c = 0.475$ and $0.75$ to correspond to Kelvin–Helmholtz mechanism.

Figure 15.

Dominant $m=0$ SPOD/CS-SPOD mode of the natural jet and the ${\textit{St}}_{\!f}=1.5$ forced jets at (a) ${\textit{St}}_c=0.05$ , (b) ${\textit{St}}_c=0.15$ , (c) ${\textit{St}}_c=0.475$ and (d) ${\textit{St}}_c=0.75$ . The pressure is shown at $\phi _{\!f} = 0$ for the forced jets. Contour limits are $\pm 0.75 |\boldsymbol{\phi }(\phi _{\!f} = 0)|_\infty$ . The solid and dashed black lines correspond to lines of $\widetilde {u}_x(\phi _{\!f} = 0) = 0.25, 0.75$ , respectively. Note the varying axis limits between different sub-figures.

The figure supports these expectations, although the 10 % forcing cases display the previously discussed phenomenon (see §§ 4.1 and 4.3) whereby the strong acoustic forcing leads to shocks and direct nonlinear interaction between the waves and the coherent structures. Apart from this effect, forcing has relatively minor impact on the coherent structures. The primary small difference is a slight reduction in the length of the Kelvin–Helmholtz mode at ${\textit{St}}_c = 0.75$ . For ${\textit{St}}_c = 0.05, 0.15$ and $0.475$ , no significant phase dependency in the dominant modes is observed even for the $10\,\%$ forcing (apart from the shock interaction). By contrast, for ${\textit{St}}_c = 0.75$ , some phase dependency is observed for the $10\,\%$ forcing in the near-nozzle region.

The lack of phase dependency for the lower frequencies is due to several reasons. First, the dominant wavelength present in the ${\textit{St}}_c = 0.05, 0.15$ and $0.475$ modes is much greater than the wavelength of the mean-flow deformation caused by the ${\textit{St}}_{\!f}=1.5$ forcing. This wavelength mismatch suggests that the low-frequency coherent structures evolve independently to the phase-dependent mean-flow variation, instead, they respond to changes in the period-averaged mean. In contrast, for ${\textit{St}}_c = 0.75$ , the primary wavelength is closer to that of the forcing, resulting in phase dependency. This trend was also seen for the ${\textit{St}}_{\!f}=0.3$ cases. Second, the spatial support of these dominant low-frequency/longer-wavelength structures is much further downstream than the spatial support of the mean-flow phase dependency, resulting in a limited phase dependency of the low-frequency coherent structures. This is in support of our observations made in figure 7 which shows that the flow associated with ${\textit{St}}_{\!f}=1.5$ forcing reaches statistical stationarity downstream of $x = 2$ . Lastly, Heidt & Colonius (Reference Heidt and Colonius2024b ) showed that CS-SPOD modes are approximated by SPOD modes when there is a large energy difference between the different frequencies in $\mathcal{\varOmega }_{{\textit{St}}_c}$ , which is precisely the case for the ${\textit{St}}_{\!f}=1.5$ forcing for low ${\textit{St}}_c$ .

For the high-amplitude ${\textit{St}}_{\!f}=1.5$ forcing case across all ${\textit{St}}_c$ , but most clearly for ${\textit{St}}_c \in [0.05, 0.55]$ , we observe upstream propagating waves in the dominant CS-SPOD modes in the near-nozzle centreline region $x \in [0, 1], r \in [0, 0.4]$ . These acoustic waves are not phase locked (since the mean is removed before computing CS-SPOD). We hypothesise that these waves are a result of Kelvin–Helmholtz waves propagating downstream and interacting with the periodic mean. Similar waves are produced in screeching jets (Raman Reference Raman1999), impinging jets (Ho & Nosseir Reference Ho and Nosseir1981) and trapped modes in transonic jets (Schmidt et al. Reference Schmidt, Towne, Colonius, Cavalieri, Jordan and Brès2017; Towne et al. Reference Towne, Cavalieri, Jordan, Colonius, Schmidt, Jaunet and Brès2017b ). However, in this case, we do not observe resonance or tonal peaks; the waves occur over a range of ${\textit{St}}_c$ .

In figure 16, we display the CS-SPOD eigenspectrum for $m = 1, 2$ . Since the instantaneous data at $m \gt 0$ are complex, the sampled statistics are not equal for positive and negative frequencies. However, the true statistics are identical (Schmidt & Colonius Reference Schmidt and Colonius2020). Thus, as noted in Schmidt & Colonius (Reference Schmidt and Colonius2020), we average positive and negative ${\textit{St}}_c$ eigenvalues together to obtain a better statistical representation of the flow. Since the spectrum is symmetric about ${\textit{St}}_c = 0$ , we display the spectrum for ${\textit{St}}_c \in [0, \alpha _0/2]$ only.

Figure 16.

The CS-SPOD eigenspectrum of the natural, (a) ${\textit{St}}_{\!f}=0.3$ and (b) ${\textit{St}}_{\!f}=1.5$ forced jets for $m = 1, 2$ .

Figure 16. Long description

Four line graphs in two rows. Row (a) compares the natural jet with jets forced at a Strouhal number of 0.3; row (b) compares it with jets forced at 1.5. Each row holds two graphs, azimuthal wavenumber m = 1 on the left and m = 2 on the right.

Similar to the SPOD spectrum at higher azimuthal mode numbers (Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018), for all jets, the energy separation between the dominant and subdominant modes becomes smaller with increasing azimuthal mode numbers. In addition, as the azimuthal mode number increases, the impact of the forcing on the energy of the dominant mode decreases for all ${\textit{St}}_c \ne 0$ . The effect on the $m = 1$ spectrum is approximately similar to the $m = 0$ for all cases. For $m = 2$ , the ${\textit{St}}_{\!f}=0.3$ forcing has a limited impact with a slight increase in the energy of the dominant mode for the $a_0/U_{\!j}=1\,\%$ forcing as ${\textit{St}}_c \rightarrow 0.15$ . In contrast, the ${\textit{St}}_{\!f}=1.5$ forcing at both amplitudes results in a larger increase in the energy of the dominant mode as ${\textit{St}}_c \rightarrow 0.75$ for $m = 2$ . Interestingly, the ${\textit{St}}_{\!f}=1.5,\ a_0/U_{\!j}=1\,\%$ forcing results in a larger amplification of the dominant mode for approximately ${\textit{St}}_c \in [0.55, 0.75]$ at $m = 1, 2$ than the $a_0/U_{\!j}=10\,\%$ forcing. In contrast, the $a_0/U_{\!j}=10\,\%$ forcing results in a larger suppression of the dominant mode for ${\textit{St}}_c \in [0.25, 0.55]$ for $m = 1, 2$ but a smaller amplification for ${\textit{St}}_c \in [0.6, 0.75]$ .

In Appendix B (figures 20 and 21), we display the dominant CS-SPOD modes for $m = 1, 2$ . The trends in the CS-SPOD modes for $m = 1, 2$ are similar to $m = 0$ . In addition, the $m = 2$ CS-SPOD modes do not appear as well converged as at $m = 0, 1$ due to the decrease in the low-rank behaviour with increasing $m$ as seen in figure 16.

Overall, we find that low- and high-frequency forcing increases and decreases the energy of the axisymmetric dominant coherent structures, respectively. The dominant modes of the forced jets are found to generally share the same primary mechanism as the natural jet with strong forcing required to introduce a substantial phase dependency with a stronger phase dependency occurring as ${\textit{St}}_c \rightarrow \alpha _0/2$ . For increasing $m$ , the impact of the forcing on the energy of the dominant modes decreases while the changes to the modes are generally consistent to that of the $m = 0$ modes.

5.2. Comparison of CS-SPOD and HRA

For unforced turbulent jets, resolvent analysis (RA) has proven to be a powerful tool for understanding the dominant coherent structures and amplification mechanisms, including Kelvin–Helmholtz, Orr and lift-up mechanisms (Garnaud et al. Reference Garnaud, Lesshafft, Schmid and Huerre2013; Schmidt et al. Reference Schmidt, Towne, Rigas, Colonius and Brès2018; Lesshafft et al. Reference Lesshafft, Semeraro, Jaunet, Cavalieri and Jordan2019; Nogueira et al. Reference Nogueira, Cavalieri, Jordan and Jaunet2019; Pickering et al. Reference Pickering, Rigas, Nogueira, Cavalieri, Schmidt and Colonius2020, Reference Pickering, Rigas, Schmidt, Sipp and Colonius2021a ), and for developing reduced-order models of jet noise (Towne et al. Reference Towne, Brès and Lele2017a ; Cavalieri et al. Reference Cavalieri, Jordan and Lesshafft2019; Pickering et al. Reference Pickering, Towne, Jordan and Colonius2021b ; Karban et al. Reference Karban, Bugeat, Towne, Lesshafft, Agarwal and Jordan2023). HRA extends this framework to time-periodic mean flows (Padovan et al. Reference Padovan, Otto and Rowley2020) and has been applied to a variety of periodically forced flows, including vortex pairing in low-Reynolds-number jets (Padovan & Rowley Reference Padovan and Rowley2022) and the Rossiter mechanism in cavity flows (Islam & Sun Reference Islam and Sun2024). HRA bears the same relationship to CS-SPOD that RA bears to SPOD (Heidt & Colonius Reference Heidt and Colonius2024b ). In particular, HRA and CS-SPOD modes are identical when the harmonic resolvent expansion coefficients are uncorrelated, and their modal rankings coincide when the nonlinear forcing is white.

As discussed in Appendix C, the nonlinear forcing in turbulent jets is not white and the dominant HRA mode does not necessarily correspond to the dominant CS-SPOD mode. To account for this effect, we employ a modified HRA operator analogous to the modified resolvent operator used by Towne et al. (Reference Towne, Brès and Lele2017a ), in which the nonlinear forcing is directly modelled. The present forcing model is ad hoc and represents a limitation of the current analysis. Nevertheless, it is constructed exclusively from unsteady data of the natural jet; no unsteady information from the forced jets is employed. Consequently, the comparison between HRA and CS-SPOD provides an independent assessment of the extent to which the observed modifications of the coherent structures can be explained by the modified mean flow.

For selected values of the centre frequency, the dominant CS-SPOD and HRA modes are compared in figures 17 and 18 for the ${\textit{St}}_{\!f}=0.3$ and ${\textit{St}}_{\!f}=1.5$ forcing cases, respectively. For the ${\textit{St}}_{\!f}=0.3$ cases, the agreement between the CS-SPOD and HRA modes is generally good. At low centre frequencies, increasing the forcing amplitude shortens the spatial extent of the low-frequency wavepacket and increases the phase-dependent modulation concentrated within the high-velocity and high-shear regions of the flow. These trends are reproduced by HRA, although the amplitude of the phase-dependent component is somewhat over-predicted. For the ${\textit{St}}_{\!f}=1.5$ cases, the agreement is even stronger. Across all centre frequencies considered, the dominant HRA modes closely resemble the corresponding CS-SPOD modes, both in the large-scale wavepacket structure and in the higher-frequency components localised near the nozzle. The agreement is particularly notable at higher centre frequencies, where both analyses capture the same phase-dependent modulation of the coherent structures.

Figure 17.

Comparison between the dominant $m = 0$ (a, c) SPOD/CS-SPOD and (b, d) RA/HRA modes of the natural and the ${\textit{St}}_{\!f}=0.3$ forced jets at (a, b) ${\textit{St}}_c = 0.05$ and (c, d) ${\textit{St}}_c = 0.15$ . The pressure is shown at $\phi _{\!f} = 0$ . Contour limits are $\pm 0.75 |\boldsymbol{\phi }(\phi _{\!f} = 0)|_\infty$ . The solid and dashed black lines correspond to lines of $\widetilde {u}_x(\phi _{\!f} = 0) = 0.25, 0.75$ , respectively. Note the varying axis limits between different sub-figures.

Figure 18.

Comparison between the dominant $m = 0$ (a, c, e, g) SPOD/CS-SPOD and (b, d, f, h) RA/HRA modes of the natural and the ${\textit{St}}_{\!f}=1.5$ forced jets at (a, b) ${\textit{St}}_c = 0.05$ , (c, d) ${\textit{St}}_c = 0.15$ , (e, f) ${\textit{St}}_c = 0.475$ and (g, h) ${\textit{St}}_c = 0.75$ . The pressure is shown at $\phi _{\!f} = 0$ . Contour limits are $\pm 0.75 |\boldsymbol{\phi }(\phi _{\!f} = 0)|_\infty$ . The solid and dashed black lines correspond to lines of $\widetilde {u}_x(\phi _{\!f} = 0) = 0.25, 0.75$ , respectively. Note the varying axis limits between different sub-figures.

Figure 18. Long description

Twenty-four contour maps in eight panels arranged as four rows of two, each panel holding three maps stacked vertically. The left column, panels (a), (c), (e) and (g), shows dominant axisymmetric SPOD modes for the natural jet and CS-SPOD modes for the forced jets; the right column, panels (b), (d), (f) and (h), shows the corresponding resolvent and harmonic resolvent modes. The four rows are at centre frequencies of 0.05, 0.15, 0.475 and 0.75. Within every panel the three maps are the natural jet, then forcing at a Strouhal number of 1.5 at 1 percent amplitude, then the same forcing at 10 percent. All maps show pressure on a diverging red and blue scale with symmetric limits at plus and minus 0.75 of each mode’s peak amplitude. The horizontal axis is the axial coordinate and the vertical axis the radial coordinate. Axis limits shrink down the figure as centre frequency rises. Solid black lines mark the mean axial velocity contour 0.25 and dashed lines 0.75.

For the ${\textit{St}}_{\!f}=1.5$ , $a_0/U_{\!j}=10\,\%$ case, the dominant CS-SPOD modes exhibit broadband upstream-propagating acoustic waves in the near-nozzle region, as discussed in § 5. Figure 18 shows that these features are also present in the corresponding HRA modes. The ability of a linearised analysis about the periodically varying mean flow to reproduce these structures suggests that they are closely connected to the modified mean state. This observation is consistent with the hypothesis that the upstream-propagating acoustic waves arise from a linear mechanism associated with the interaction of instability waves with the periodic mean flow, although additional analysis would be required to establish the precise generation mechanism.

Overall, the close agreement between the dominant CS-SPOD and HRA modes indicates that a substantial fraction of the observed modifications to the coherent structures can be explained by the corresponding modifications of the mean flow. This conclusion is particularly significant given that the HRA model is constructed without using unsteady information from the forced jets. While the present implementation relies on an ad hoc model for the nonlinear forcing, the results suggest that mean-flow modification plays a central role in determining how harmonic forcing alters the dominant coherent structures of turbulent jets. Future work should focus on replacing the present forcing model with a physics-based representation of the nonlinear forcing and on assessing the generality of these findings across a broader range of forcing conditions.