2.2.1 Two different behaviours

Axisymmetric harmonic forcing at the nozzle of a laminar round jet excites, over a wide range of frequencies, a linear shear instability of the steady flow state, leading to exponential growth of the perturbation amplitude along the axial direction. As the amplitude reaches nonlinear levels, the shear layer rolls up into a regular street of vortex rings, which form and convect at the frequency of the applied forcing. Depending on flow parameters and forcing frequency, these vortices may undergo subsequent pairing (Hussain & Zaman Reference Hussain and Zaman1980; Zaman & Hussain Reference Zaman and Hussain1980; Ho & Huerre Reference Ho and Huerre1984; Shaabani-Ardali et al. Reference Shaabani-Ardali, Sipp and Lesshafft2019), and if the ambient flow is sufficiently quiet and the harmonic forcing is well-controlled, this pairing takes place in a perfectly regular fashion.

In cases where pairing occurs, two neighbouring vortices merge into one, such that the passage frequency of vortices downstream of the pairing location is exactly half that of the imposed forcing. If the forcing is characterised by the time period $T_{f}$, such that $\unicode[STIX]{x1D714}_{f}=2\unicode[STIX]{x03C0}/T_{f}$, the ‘paired state’ is globally $2T_{f}$-periodic ($T_{f}$-periodic upstream of the pairing and $2T_{f}$-periodic downstream). An example of this behaviour, obtained by DNS, is shown in figure 5(a).

Figure 5.

Vorticity snapshots of the periodic paired and unpaired states, obtained naturally for two different parameter settings. Simulation details, forcing amplitude, Reynolds and Strouhal numbers are defined in § 2.2.2. (a) Paired state ($2T_{f}$-periodic) at $St=0.6$, $Re=2000$ and $A=5\,\%$ at $t=3T_{f}/2$. (b) Unpaired state ($T_{f}$-periodic) at $St=0.6$, $Re=1300$ and $A=5\,\%$ at $t=T_{f}/2$. (c) Unpaired stabilised state ($T_{f}$-periodic) at $St=0.6$, $Re=2000$ and $A=5\,\%$ at $t=T_{f}/2$.

Another case at different parameter settings, where no pairing is found to occur, is shown in figure 5(b). Vortices roll up close to the nozzle and advect downstream, until they are dissipated by viscosity. Such a flow state is (globally) $T_{f}$-periodic and will be called hereafter an ‘unpaired state’.

It has been shown (Shaabani-Ardali et al. Reference Shaabani-Ardali, Sipp and Lesshafft2019) that the unpaired state exists for all parameter values, but that this state may be Floquet-unstable (Floquet Reference Floquet1883). In particular, it was found that, when unstable, the most unstable mode was always subharmonic at $\unicode[STIX]{x1D714}_{f}/2$, leading after destabilisation to a paired state.

In the case of a bifurcating jet, none of the vortices undergo pairing, as shown in figures 1 and 2. In addition, the physical mechanism of Parekh et al. (Reference Parekh, Leonard and Reynolds1988) explains bifurcation by an alternating shifting of the ring vortices, which is achieved through helical forcing at half the vortex frequency. In the presence of pairing, paired vortices convect at the same frequency as the helical forcing, cancelling the alternating behaviour. Therefore, for the optimisation of the helical forcing, it is natural to seek an unpaired base flow, even if this flow may be unstable with respect to axisymmetric perturbations (Shaabani-Ardali et al. Reference Shaabani-Ardali, Sipp and Lesshafft2019). To do so, a simple stabilisation tool (Shaabani-Ardali, Sipp & Lesshafft Reference Shaabani-Ardali, Sipp and Lesshafft2017) is used to compute these unpaired states, which consists in adding a forcing term of the form

(2.4)$$\begin{eqnarray}\boldsymbol{f}=-\unicode[STIX]{x1D706}(\boldsymbol{U}(t)-\boldsymbol{U}(t-T_{f}))\end{eqnarray}$$

to the right-hand side of the Navier–Stokes equations (2.5). This allows to suppress $2T_{f}$-periodic fluctuations on $T_{f}$-periodic dynamics. In this framework, $\unicode[STIX]{x1D706}$ is a forcing parameter that needs to be prescribed; an optimal value of $0.044\unicode[STIX]{x1D714}_{f}$ has been identified by Shaabani-Ardali et al. (Reference Shaabani-Ardali, Sipp and Lesshafft2017). In the same way as the selective frequency damping technique for steady flows (Åkervik et al. Reference Åkervik, Brandt, Henningson, Hœpffner, Marxen and Schlatter2006), the forcing term vanishes as the system converges towards a $T_{f}$-periodic unpaired state, such that the recovered state is a consistent solution of the unforced Navier–Stokes equations. The results of this technique are shown in figure 5(c), where the stabilised flow is plotted.

2.2.2 Configuration, governing equations and numerical discretisation

DNS were carried out using Nek5000, an incompressible spectral element code (Fischer et al. Reference Fischer, Kruse, Mullen, Tufo, Lottes and Kerkemeier2008). An axisymmetric laminar jet is described in cylindrical coordinates $(z,r)$, $z$ being the main flow direction and $r$ being the radial distance from the jet axis. The flow is assumed to be governed by the incompressible Navier–Stokes equations with zero azimuthal velocity, written in dimensionless form as

(2.5)$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\boldsymbol{U}}{\unicode[STIX]{x2202}t}+(\boldsymbol{U}\boldsymbol{\cdot }\unicode[STIX]{x1D735})\boldsymbol{U}=-\unicode[STIX]{x1D735}P+\frac{1}{Re}\unicode[STIX]{x0394}\boldsymbol{U},\quad \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{U}=0.\end{eqnarray}$$

The velocity $\boldsymbol{U}$ has axial and radial components $U$ and $V$, and $P$ denotes pressure. The jet diameter $D$ and the inlet centreline velocity $\overline{U}$ are used to render the problem non-dimensional, defining the Reynolds number as $Re=\overline{U}D/\unicode[STIX]{x1D708}$, with $\unicode[STIX]{x1D708}$ the kinematic viscosity. The computational domain extends over $15\times 5$ diameters in the axial and radial directions, and it is discretised with 6600 spectral elements, each containing 64 mesh points. Boundary conditions are specified as follows.

1. (i) In the inlet plane, $z=0$, a hyperbolic-tangent velocity profile is imposed. In dimensionless form, its amplitude is modulated in time as

   (2.6)$$\begin{eqnarray}\displaystyle \boldsymbol{U}(r,t) & = & \displaystyle \frac{1}{2}\left\{1-\tanh \left[\frac{1}{4\unicode[STIX]{x1D703}_{0}}\left(r-\frac{1}{4r}\right)\right]\right\}(1+A\cos (\unicode[STIX]{x1D714}_{f}t))\boldsymbol{e}_{\boldsymbol{z}}\nonumber\\ \displaystyle & = & \displaystyle U_{0}(r)(1+A\cos (\unicode[STIX]{x1D714}_{f}t))\boldsymbol{e}_{\boldsymbol{z}},\end{eqnarray}$$
   where $A$ is the forcing amplitude of the jet, $\unicode[STIX]{x1D703}_{0}$ is the initial dimensionless mixing layer thickness and $\unicode[STIX]{x1D714}_{f}$ is the dimensional forcing frequency. The forcing period is given by $T_{f}=2\unicode[STIX]{x03C0}/\unicode[STIX]{x1D714}_{f}$, and the Strouhal number is defined as $St=\unicode[STIX]{x1D714}_{f}D/(2\unicode[STIX]{x03C0}\overline{U})$.

2. (ii) On the centreline of the jet, $r=0$, axisymmetric boundary conditions are imposed,

   (2.7)$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}U}{\unicode[STIX]{x2202}r}=V=\frac{\unicode[STIX]{x2202}P}{\unicode[STIX]{x2202}r}=0.\end{eqnarray}$$

3. (iii) In the outlet plane, $z=Z_{m}=15$, and on the outer radial boundary, $r=R_{m}=5$, a stress-free outflow condition is applied,

   (2.8)$$\begin{eqnarray}-P\boldsymbol{n}+\frac{1}{Re}(\unicode[STIX]{x1D735}\boldsymbol{U})\boldsymbol{\cdot }\boldsymbol{n}=0,\end{eqnarray}$$
   with $\boldsymbol{n}$ the normal vector at the boundary.

The flow configuration is thus characterised by the Reynolds number $Re$, the Strouhal number $St$, the dimensionless mixing layer thickness $\unicode[STIX]{x1D703}_{0}$ and the forcing amplitude $A$. In this study, we keep the jet parameters fixed by considering a round jet at $Re=2000$ with a mixing layer thickness of $\unicode[STIX]{x1D703}_{0}=0.025$. The amplitude of the axisymmetric forcing is also kept fixed at $A=0.05$, whereas the forcing Strouhal number is varied. For all considered Strouhal numbers, convergence towards a periodic unpaired base state has been achieved. Shaabani-Ardali et al. (Reference Shaabani-Ardali, Sipp and Lesshafft2019) found that the band $0.5\leqslant St\leqslant 0.8$ is linearly unstable to subharmonic perturbations that lead to pairing; however, the time-delay feedback stabilisation in the present calculations successfully yields periodic unpaired base states over the entire considered range of Strouhal numbers.

Numerical mesh convergence in these calculations, for the example case $St=0.6$, is demonstrated in figure 6 by comparing the evolution of residuals in the stabilised computations for different spectral polynomial orders ($n=4$, 6, 8 and 10). A higher order is equivalent to a higher number of collocation points within each spectral element. The plotted quantity $e(t)$ in figure 6 is the norm of the difference between the flow states at time $t$ and at time $t-T_{f}$, as defined in Shaabani-Ardali et al. (Reference Shaabani-Ardali, Sipp and Lesshafft2017). It is seen that all orders yield the same convergence of the stabilisation routine; between resolutions $n\geqslant 6$, the flow states at each instant in the simulations are practically identical. The order $n=8$ is used in all following computations.

Figure 6.

For different mesh resolutions (polynomial order $n$ of individual spectral elements), the norm $e(t)$ of the residual difference between the flow velocity at time $t$ and time $t-T_{f}$ is traced as a function of $t$. These are the stabilised computations (with source term (2.4)) performed to find the periodic base flow at $Re=2000$ and $St=0.6$. The converged setting $n=8$ is used in the following.

2.3.1 Equations governing helical perturbations

Before writing the linearised Navier–Stokes equations, the $m=1$ perturbation field $(\boldsymbol{u}_{\mathbf{1}},p_{1})$ from (2.1) is rewritten as $(\boldsymbol{u}_{\mathbf{2}},p_{2})$ with $\boldsymbol{u}_{\mathbf{2}}=(u_{1,r},iu_{1,\unicode[STIX]{x1D703}},u_{1,z})^{t}$ and $p_{2}=p_{1}$ in order to remove imaginary coefficients. The linearised Navier–Stokes equations for the perturbation $(\boldsymbol{u}_{\mathbf{2}},p_{2})$ around the $T$-periodic base flow $\boldsymbol{U}(r,z,t)$ then read

(2.9)$$\begin{eqnarray}\displaystyle & \displaystyle \frac{1}{r}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}r}(ru_{2,r})+\frac{1}{r}u_{2,\unicode[STIX]{x1D703}}+\frac{\unicode[STIX]{x2202}u_{2,z}}{\unicode[STIX]{x2202}z}=0, & \displaystyle\end{eqnarray}$$

(2.10)$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}u_{2,r}}{\unicode[STIX]{x2202}t}+U_{r}\frac{\unicode[STIX]{x2202}u_{2,r}}{\unicode[STIX]{x2202}r}+u_{2,r}\frac{\unicode[STIX]{x2202}U_{r}}{\unicode[STIX]{x2202}r}+U_{z}\frac{\unicode[STIX]{x2202}u_{2,r}}{\unicode[STIX]{x2202}z}+u_{2,z}\frac{\unicode[STIX]{x2202}U_{r}}{\unicode[STIX]{x2202}z}=-\frac{\unicode[STIX]{x2202}p_{2}}{\unicode[STIX]{x2202}r}+\frac{1}{Re}\unicode[STIX]{x1D6E5}_{r}\boldsymbol{u}_{\mathbf{2}}, & \displaystyle\end{eqnarray}$$

(2.11)$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}u_{2,\unicode[STIX]{x1D703}}}{\unicode[STIX]{x2202}t}+U_{r}\frac{\unicode[STIX]{x2202}u_{2,\unicode[STIX]{x1D703}}}{\unicode[STIX]{x2202}r}+U_{z}\frac{\unicode[STIX]{x2202}u_{2,\unicode[STIX]{x1D703}}}{\unicode[STIX]{x2202}z}+\frac{U_{r}u_{2,\unicode[STIX]{x1D703}}}{r}=\frac{1}{r}p_{2}+\frac{1}{Re}\unicode[STIX]{x1D6E5}_{\unicode[STIX]{x1D703}}\boldsymbol{u}_{\mathbf{2}}, & \displaystyle\end{eqnarray}$$

(2.12)$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}u_{2,z}}{\unicode[STIX]{x2202}t}+U_{r}\frac{\unicode[STIX]{x2202}u_{2,z}}{\unicode[STIX]{x2202}r}+u_{2,r}\frac{\unicode[STIX]{x2202}U_{z}}{\unicode[STIX]{x2202}r}+U_{z}\frac{\unicode[STIX]{x2202}u_{2,z}}{\unicode[STIX]{x2202}z}+u_{2,z}\frac{\unicode[STIX]{x2202}U_{z}}{\unicode[STIX]{x2202}z}=-\frac{\unicode[STIX]{x2202}p_{2}}{\unicode[STIX]{x2202}z}+\frac{1}{Re}\unicode[STIX]{x1D6E5}_{z}\boldsymbol{u}_{\mathbf{2}}, & \displaystyle\end{eqnarray}$$

with the Laplacian in cylindrical coordinates

(2.13)$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6E5}_{r}\boldsymbol{u}_{\mathbf{2}}=\frac{1}{r}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}r}\left(r\frac{\unicode[STIX]{x2202}u_{2,r}}{\unicode[STIX]{x2202}r}\right)+\frac{\unicode[STIX]{x2202}^{2}u_{2,r}}{\unicode[STIX]{x2202}z^{2}}-\frac{2u_{2,r}}{r^{2}}-\frac{2u_{2,\unicode[STIX]{x1D703}}}{r^{2}}, & \displaystyle\end{eqnarray}$$

(2.14)$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6E5}_{\unicode[STIX]{x1D703}}\boldsymbol{u}_{\mathbf{2}}=\frac{1}{r}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}r}\left(r\frac{\unicode[STIX]{x2202}u_{2,\unicode[STIX]{x1D703}}}{\unicode[STIX]{x2202}r}\right)+\frac{\unicode[STIX]{x2202}^{2}u_{2,\unicode[STIX]{x1D703}}}{\unicode[STIX]{x2202}z^{2}}-\frac{2u_{2,\unicode[STIX]{x1D703}}}{r^{2}}-\frac{2u_{2,r}}{r^{2}}, & \displaystyle\end{eqnarray}$$

(2.15)$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6E5}_{z}\boldsymbol{u}_{\mathbf{2}}=\frac{1}{r}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}r}\left(r\frac{\unicode[STIX]{x2202}u_{2,z}}{\unicode[STIX]{x2202}r}\right)-\frac{u_{2,z}}{r^{2}}+\frac{\unicode[STIX]{x2202}^{2}u_{2,z}}{\unicode[STIX]{x2202}z^{2}}. & \displaystyle\end{eqnarray}$$

The following boundary and compatibility conditions are imposed:

(2.16)$$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{u}_{\mathbf{2}}=\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{1}}(r)\cos (\unicode[STIX]{x1D714}_{f}t/2)+\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{2}}(r)\sin (\unicode[STIX]{x1D714}_{f}t/2)\quad \text{at}~z=0, & \displaystyle\end{eqnarray}$$

(2.17)$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}u_{2,r}}{\unicode[STIX]{x2202}r}=\frac{\unicode[STIX]{x2202}u_{2,\unicode[STIX]{x1D703}}}{\unicode[STIX]{x2202}r}=u_{2,z}=0\quad \text{at}~r=0, & \displaystyle\end{eqnarray}$$

(2.18)$$\begin{eqnarray}\displaystyle & \displaystyle u_{2,r}+u_{2,\unicode[STIX]{x1D703}}=0\quad \text{at}~r=0, & \displaystyle\end{eqnarray}$$

(2.19)$$\begin{eqnarray}\displaystyle & \displaystyle \frac{1}{Re}\frac{\unicode[STIX]{x2202}u_{2,r}}{\unicode[STIX]{x2202}r}-p=\frac{\unicode[STIX]{x2202}u_{2,\unicode[STIX]{x1D703}}}{\unicode[STIX]{x2202}r}=\frac{\unicode[STIX]{x2202}u_{2,z}}{\unicode[STIX]{x2202}r}=0\quad \text{at}~r=R_{m}, & \displaystyle\end{eqnarray}$$

(2.20)$$\begin{eqnarray}\displaystyle & \displaystyle \frac{1}{Re}\frac{\unicode[STIX]{x2202}u_{2,z}}{\unicode[STIX]{x2202}z}-p=\frac{\unicode[STIX]{x2202}u_{2,r}}{\unicode[STIX]{x2202}z}=\frac{\unicode[STIX]{x2202}u_{2,\unicode[STIX]{x1D703}}}{\unicode[STIX]{x2202}z}=0\quad \text{at}~z=Z_{m}, & \displaystyle\end{eqnarray}$$

with $\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{1}}$ and $\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{2}}$ the corresponding inlet forcing modified from $\boldsymbol{u}_{\mathbf{1}\boldsymbol{f},\mathbf{1}}$ and $\boldsymbol{u}_{\mathbf{1}\boldsymbol{f},\mathbf{2}}$. In this new formalism, equations for $\boldsymbol{u}_{\mathbf{2}}$ only involve real coefficients, which is much more convenient for our numerical computations. The stress-free conditions (2.19), (2.20), consistent with the Nek5000 boundary treatment (2.8), arise as the ‘natural boundary conditions’ (see § 10.2.4 in Dick Reference Dick and Wendt2009) in this finite-element formulation. Their quality as non-reflecting outflow conditions have been examined by Theofilis (Reference Theofilis2017) and by Lesshafft (Reference Lesshafft2018).

2.3.2 Optimisation procedure

We aim to find the optimal inlet forcing that maximally shifts the vortices off the axis, at a given optimisation time $T_{o}$, per unit input energy. Formally, we want to solve

(2.21)$$\begin{eqnarray}\max _{\Vert \boldsymbol{u}^{\prime }(r,\unicode[STIX]{x1D703},z=0,t)\Vert _{f}=1}\Vert \boldsymbol{u}^{\prime }(r,\unicode[STIX]{x1D703},z,t=T_{o})\Vert ,\end{eqnarray}$$

with $T_{o}\longrightarrow \infty$. Again, the forcing $\boldsymbol{u}^{\prime }(r,\unicode[STIX]{x1D703},z=0,t)$ is assumed periodic with frequency $\unicode[STIX]{x1D714}_{f}/2$ (see (2.3)). The different norms used for forcing and response are detailed in § 2.4.

As will be demonstrated in § 3, this problem is well-posed, because the response converges towards an $\unicode[STIX]{x1D714}_{f}/2$-periodic limit cycle when $T_{o}\longrightarrow \infty$. This shows that the base flow is stable with respect to $m=1$ subharmonic perturbations.

If the base flow were steady, the forced equations for linear perturbations could be written as

(2.22)$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\boldsymbol{u}^{\prime }}{\unicode[STIX]{x2202}t}=L\boldsymbol{u}+\hat{f}\text{e}^{\text{i}\unicode[STIX]{x1D714}t},\end{eqnarray}$$

which would lead to a traditional resolvent analysis (Schmid Reference Schmid2007). As in our case $L$ depends on $t$, this type of analysis is not possible. However, several algorithms are available for the problem at hand.

First, a classical direct-adjoint algorithm, already employed for periodic-flow optimal perturbation analysis by Blackburn, Sherwin & Barkley (Reference Blackburn, Sherwin and Barkley2008), can be used. However, in the present optimisation, the flow perturbation field $\boldsymbol{u}_{\mathbf{2}}$ is two-dimensional, whereas the prescribed inlet forcing distributions $\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{1}}(r)$ and $\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{2}}(r)$ are only one-dimensional. Therefore, if a direct-adjoint calculation is carried out, the adjoint field eventually needs to be projected on the one-dimensional inlet, and extra care needs to be taken to ensure that this is done properly (see Boujo & Gallaire Reference Boujo and Gallaire2015).

Another technique relies on the decomposition of the linearised operator, the forcing and response into Fourier components, and to solve the equations in Fourier space, as it is done in harmonic-balance methods. This has not been attempted here because of the considerable changes needed in the code compared with a standard time-stepping method, and because of the strong stiffness of these Newton-like problems.

In this work, the optimisation is carried out through a classical direct time-marching scheme detailed in § 2.3.3. A basis of forcing functions is chosen and each individual basis function is advanced in time. As the forcing only varies along the single radial dimension, a 300-dimensional basis has been found to be sufficient for an accurate resolution of the forcing space. Once time-stepping is carried out, thanks to linearity, the optimal combination of all 300 forcing vectors is calculated. This technique has several advantages compared with the two mentioned previously.

First, because the time-stepping and the optimisation are uncoupled in this process, many optimisation parameters, such as the final time $T_{o}$, the norm chosen for the optimisation or other norm-related parameters, do not need to be fixed beforehand. These can be adjusted afterwards with the time-stepping results in order to check the dependence on these technical parameters. Changes in these hyper-parameters can be made almost for free.

Second, in our case the time-marching relies on an implicit scheme based on matrix inversions (see § 2.3.3); as the base flow changes in time, the matrix needs to be recalculated at every time step, requiring a significant computational effort. Therefore, time-marching $N\gg 1$ vectors simultaneously, instead of one at a time for direct-adjoint optimisation, saves computational time. In addition, our approach removes the need to compute the adjoint part, thereby saving additional time.

Finally, compared with harmonic-balance techniques, which only capture the infinite time-horizon behaviour, our technique also allows to resolve transient dynamics.

2.3.3 Implementation

The evolution of the perturbation equations defined in § 2.3.1 is calculated using a fully implicit finite-difference time-stepping scheme of second order, implemented in the finite-element solver FreeFem++ (Hecht Reference Hecht2012). The mesh has the same size and resolution as that used in the Nek5000 calculations. We use $P2$ finite elements for the velocity perturbation, whereas $P1$ finite elements are used for the pressure.

Once time-stepping is completed, the calculation of the optimal forcing combination can be done easily with Octave/Matlab routines. The resulting optimal forcing and response fields from this combination can then be obtained in FreeFem++.

This study is carried out for 10 values of the Strouhal number, between 0.35 and 0.8 in increments of 0.05, at a Reynolds number of 2000.

2.4.1 ${\mathcal{L}}_{2}$ norm for the forcing

For the forcing $\boldsymbol{u}^{\prime }$, the standard ${\mathcal{L}}_{2}$ energy norm is used. This norm is derived from the following scalar product, which involves an integration over an entire subharmonic forcing period $2T_{f}$:

(2.23)$$\begin{eqnarray}\displaystyle \langle \boldsymbol{u}^{\prime },\boldsymbol{v}^{\prime }\rangle _{f} & = & \displaystyle \frac{1}{2T_{f}}\frac{1}{2\unicode[STIX]{x03C0}}\int _{0}^{2T_{f}}\!\int _{0}^{2\unicode[STIX]{x03C0}}\!\int _{0}^{R_{m}}r\langle \boldsymbol{u}^{\prime }(r,\unicode[STIX]{x1D703},z=0,t),\boldsymbol{v}^{\prime }(r,\unicode[STIX]{x1D703},z=0,t)\rangle \,\text{d}r\,\text{d}\unicode[STIX]{x1D703}\,\text{d}t,\qquad \quad\end{eqnarray}$$

(2.24)$$\begin{eqnarray}\displaystyle & = & \displaystyle \int _{0}^{R_{m}}r\left(\langle \text{Re}(\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{1}}),\text{Re}(\boldsymbol{v}_{\mathbf{2}\boldsymbol{f},\mathbf{1}})\rangle +\langle \text{Im}(\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{1}}),\text{Im}(\boldsymbol{v}_{\mathbf{2}\boldsymbol{f},\mathbf{1}})\rangle \right.\nonumber\\ \displaystyle & & \displaystyle +\left.\langle \text{Re}(\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{2}}),\text{Re}(\boldsymbol{v}_{\mathbf{2}\boldsymbol{f},\mathbf{2}})\rangle +\langle \text{Im}(\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{2}}),\text{Im}(\boldsymbol{v}_{\mathbf{2}\boldsymbol{f},\mathbf{2}})\rangle \right)\text{d}r,\end{eqnarray}$$

(2.25)$$\begin{eqnarray}\displaystyle & = & \displaystyle \frac{1}{2}\int _{0}^{R_{m}}r(\langle \boldsymbol{u}_{\boldsymbol{ 2}\boldsymbol{f},\mathbf{1}},\boldsymbol{v}_{\mathbf{2}\boldsymbol{f},\mathbf{1}}\rangle +\langle \overline{\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{1}}},\overline{\boldsymbol{v}_{\mathbf{2}\boldsymbol{f},\mathbf{1}}}\rangle +\langle \boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{2}},\boldsymbol{v}_{\mathbf{2}\boldsymbol{f},\mathbf{2}}\rangle +\langle \overline{\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{2}}},\overline{\boldsymbol{v}_{\mathbf{2}\boldsymbol{f},\mathbf{2}}}\rangle )\,\text{d}r.\qquad \quad\end{eqnarray}$$

This scalar product is different from the standard scalar product for complex numbers. While this must be kept in mind during the optimisation process, the norm it induces for the forcing $\boldsymbol{u}^{\prime }$ is the standard ${\mathcal{L}}_{2}$ norm:

(2.26)$$\begin{eqnarray}\Vert \boldsymbol{u}^{\prime }\Vert _{f}^{2}=\int _{0}^{R_{m}}r(\langle \boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{1}},\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{1}}\rangle +\langle \boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{2}},\boldsymbol{u}_{\mathbf{2}\boldsymbol{f},\mathbf{2}}\rangle )\,\text{d}r.\end{eqnarray}$$

2.4.2 ${\mathcal{L}}_{2}$ norm for the flow response

For the response, one possible norm is the energy norm based on the scalar product

(2.27)$$\begin{eqnarray}\displaystyle & & \displaystyle \!\!\!\langle \boldsymbol{u}^{\prime }(T_{o}),\boldsymbol{v}^{\prime }(T_{o})\rangle _{2}\nonumber\\ \displaystyle & & \displaystyle \!\!\!\quad =\frac{1}{2}\frac{1}{2\unicode[STIX]{x03C0}}\int _{0}^{Z_{opt}}\int _{0}^{2\unicode[STIX]{x03C0}}\int _{0}^{R_{m}}r\langle \boldsymbol{u}^{\prime }(r,\unicode[STIX]{x1D703},z,T_{o}),\boldsymbol{v}^{\prime }(r,\unicode[STIX]{x1D703},z,T_{o})\rangle \,\text{d}r\,\text{d}\unicode[STIX]{x1D703}\,\text{d}z\nonumber\\ \displaystyle & & \displaystyle \!\!\!\quad =\int _{0}^{Z_{opt}}\int _{0}^{R_{m}}r(\langle \text{Re}(\boldsymbol{u}_{\mathbf{2}}),\text{Re}(\boldsymbol{v}_{\mathbf{2}})\rangle +\langle \text{Im}(\boldsymbol{u}_{\mathbf{2}}),\text{Im}(\boldsymbol{v}_{\mathbf{2}})\rangle )\,\text{d}r\,\text{d}z,\end{eqnarray}$$

(2.28)$$\begin{eqnarray}\displaystyle & & \displaystyle \!\!\!\quad =\frac{1}{2}\int _{0}^{Z_{opt}}\int _{0}^{R_{m}}r(\langle \boldsymbol{u}_{\boldsymbol{ 2}}(r,z,T_{o}),\boldsymbol{v}_{\mathbf{2}}(r,z,T_{o})\rangle +\langle \overline{\boldsymbol{u}_{\mathbf{2}}}(r,z,T_{o}),\overline{\boldsymbol{v}_{\mathbf{2}}}(r,z,T_{o})\rangle )\,\text{d}r\,\text{d}z.\qquad \quad\end{eqnarray}$$

As we are interested in the evolution of perturbations close to the inlet, this scalar product is calculated on a reduced domain $(r,z)\in [0;R_{m}]\times [0;Z_{opt}]$. Unless specified otherwise, a value of $Z_{opt}=5$ is chosen in the following. The effect of this value is investigated in § 3.4.1. Once again, this scalar product is different from the standard scalar product for complex numbers. However, the corresponding ${\mathcal{L}}_{2}$ norm for the response remains the same,

(2.29)$$\begin{eqnarray}\Vert \boldsymbol{u}^{\prime }(T_{o})\Vert _{2}^{2}=\int _{0}^{Z_{opt}}\int _{0}^{R_{m}}r\langle \boldsymbol{u}_{\mathbf{2}}(r,z,T_{o}),\boldsymbol{u}_{\mathbf{2}}(r,z,T_{o})\rangle \,\text{d}r\,\text{d}z.\end{eqnarray}$$

Note that, for the response, both ${\mathcal{L}}_{2}$ scalar product and norm are not integrated over a full forcing period since this would require to store every time-step of every response field: storage is only done once in every forcing period $T_{f}$. A zero time-shift with respect to the axisymmetric forcing is chosen: $T_{o}$ is taken as $2nT_{f}$, with $n\in \mathbb{N}$. The influence of this shift is investigated in § 3.4.4, and the influence of $n$ is discussed in § 3.4.5.

2.4.3 Displacement norm for response

The bifurcation scenario (Lee & Reynolds Reference Lee and Reynolds1985; Parekh et al. Reference Parekh, Leonard and Reynolds1988) relies on a slight tilt and shift of the vortices near the jet inlet. Therefore, in this section, we define an alternative norm that quantifies the displacement of vortices from their unperturbed position. The location of a given vortex $i$ at $T_{o}$ in a $(r,z)$-plane, cut at an angle $\unicode[STIX]{x1D703}_{c}$, can be defined as the centroid of its vorticity $\unicode[STIX]{x1D714}_{\unicode[STIX]{x1D703}}$ along $\unicode[STIX]{x1D703}$

(2.30)$$\begin{eqnarray}\boldsymbol{x}_{i}(\unicode[STIX]{x1D703}_{c},T_{o})=\frac{\displaystyle \iint _{r,z\in \unicode[STIX]{x1D734}_{i}}r\unicode[STIX]{x1D714}_{\unicode[STIX]{x1D703}}(r,\unicode[STIX]{x1D703}_{c},z,T_{o})\boldsymbol{x}\,\text{d}r\,\text{d}z}{\displaystyle \iint _{r,z\in \unicode[STIX]{x1D734}_{i}}r\unicode[STIX]{x1D714}_{\unicode[STIX]{x1D703}}(r,\unicode[STIX]{x1D703}_{c},z,T_{o})\,\text{d}r\,\text{d}z},\end{eqnarray}$$

with $\unicode[STIX]{x1D734}_{i}$ the domain of the $i\text{th}$ vortex. This domain is defined as the rectangle $(z,r)\in [Z_{i}-\unicode[STIX]{x1D6FF}_{z}/2;Z_{i}+\unicode[STIX]{x1D6FF}_{z}/2]\times [0;1]$, with $Z_{i}$ the axial component of the centroid of the $i\text{th}$ vortex and $\unicode[STIX]{x1D6FF}_{z}$ conveniently fixed at 0.8. Examples of such boxes are represented in figure 7. The effect of $\unicode[STIX]{x1D6FF}_{z}$ on the optimisation is investigated in § 3.4.2.

Figure 7.

Vorticity cut of the $Re=2000$, $St=0.50$ base unpaired flow, with vortex domains represented for the first three vortices. For the second vortex, four $\unicode[STIX]{x1D6FF}_{z}$ sizes are sketched (0.4, 0.6, 0.8, 1.0), whereas for the first and third ones, only the $\unicode[STIX]{x1D6FF}_{z}=0.8$ box is represented.

The vorticity field $\unicode[STIX]{x1D714}_{\unicode[STIX]{x1D703}}$ can be written as the sum of the vorticity $\unicode[STIX]{x1D6FA}_{\unicode[STIX]{x1D703}}$ induced by the base flow $\boldsymbol{U}(r,z,t)$ and an infinitesimal perturbation $\unicode[STIX]{x1D714}_{\unicode[STIX]{x1D703}}^{\prime }$ induced by the $m=1$ perturbation $\boldsymbol{u}^{\prime }$:

(2.31)$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6FA}_{\unicode[STIX]{x1D703}}(r,z,t) & = & \displaystyle \frac{\unicode[STIX]{x2202}U_{r}}{\unicode[STIX]{x2202}z}-\frac{\unicode[STIX]{x2202}U_{z}}{\unicode[STIX]{x2202}r},\end{eqnarray}$$

(2.32)$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D714}_{\unicode[STIX]{x1D703}}^{\prime }(r,\unicode[STIX]{x1D703},z,t) & = & \displaystyle \frac{\unicode[STIX]{x2202}u_{r}^{\prime }}{\unicode[STIX]{x2202}z}-\frac{\unicode[STIX]{x2202}u_{z}^{\prime }}{\unicode[STIX]{x2202}r}\nonumber\\ \displaystyle & = & \displaystyle \left(\frac{\unicode[STIX]{x2202}u_{2,r}}{\unicode[STIX]{x2202}z}-\frac{\unicode[STIX]{x2202}u_{2,z}}{\unicode[STIX]{x2202}r}\right)\text{e}^{\text{i}\unicode[STIX]{x1D703}}+\text{c.c.}=\unicode[STIX]{x1D714}_{2,\unicode[STIX]{x1D703}}(r,z,t)\text{e}^{\text{i}\unicode[STIX]{x1D703}}+\text{c.c.}\end{eqnarray}$$

Specifying a perturbation amplitude $\unicode[STIX]{x1D716}\ll 1$, we obtain

(2.33)$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D714}_{\unicode[STIX]{x1D703}}(r,\unicode[STIX]{x1D703},z,t) & = & \displaystyle \unicode[STIX]{x1D6FA}_{\unicode[STIX]{x1D703}}(r,z,t)+\unicode[STIX]{x1D716}\unicode[STIX]{x1D714}_{\unicode[STIX]{x1D703}}^{\prime }(r,\unicode[STIX]{x1D703},z,t),\end{eqnarray}$$

(2.34)$$\begin{eqnarray}\displaystyle & = & \displaystyle \unicode[STIX]{x1D6FA}_{\unicode[STIX]{x1D703}}(r,z,t)+\unicode[STIX]{x1D716}\unicode[STIX]{x1D714}_{2,\unicode[STIX]{x1D703}}(r,z,t)\text{e}^{\text{i}\unicode[STIX]{x1D703}}+\unicode[STIX]{x1D716}\overline{\unicode[STIX]{x1D714}_{2,\unicode[STIX]{x1D703}}}(r,z,t)\text{e}^{-\text{i}\unicode[STIX]{x1D703}}.\end{eqnarray}$$

The displacement norm with respect to the $i\text{th}$ vortex $\unicode[STIX]{x1D734}_{i}$ can then be defined as

(2.35)$$\begin{eqnarray}\Vert \boldsymbol{u}^{\prime }(T_{o})\Vert _{disp,\unicode[STIX]{x1D734}_{i}}^{2}=\frac{1}{2}\frac{1}{2\unicode[STIX]{x03C0}}\int _{0}^{2\unicode[STIX]{x03C0}}\left\Vert \left.\frac{\text{d}\boldsymbol{x}_{i}(\unicode[STIX]{x1D703}_{c},T_{o})}{\text{d}\unicode[STIX]{x1D716}}\right|_{\unicode[STIX]{x1D716}=0}\right\Vert _{2}^{2}\text{d}\unicode[STIX]{x1D703}_{c}.\end{eqnarray}$$

For each vortex considered, it is possible to define a corresponding displacement norm. In the present study, following the bifurcation scenario (Parekh et al. Reference Parekh, Leonard and Reynolds1988), the displacement norm for the first three vortices is considered. Then, the axial position of the third vortex, depending on the Strouhal number, approximately corresponds to 5 jet diameters, which is consistent with the optimisation domain chosen for the ${\mathcal{L}}_{2}$ norm (2.27). Full details on this norm are given in appendix A.

Note that the displacement norm is not strictly a norm but a semi-norm: $\Vert \boldsymbol{u}^{\prime }(T_{o})\Vert _{disp,\unicode[STIX]{x1D734}_{i}}=0$ does not necessarily imply that $\boldsymbol{u}^{\prime }(T_{o})=\mathbf{0}$. As can be seen in appendix A and in figure 8, the displacement norm being equal to zero only implies that the vorticity centroid of the base flow and of the perturbation are at the same location. Therefore, the alternating pattern of figure 8(b) is a good candidate for a large displacement norm at fixed ${\mathcal{L}}_{2}$ norm. In this figure, one can easily understand that a large ${\mathcal{L}}_{2}$ norm does not necessarily imply a large displacement norm. Conversely, because the ${\mathcal{L}}_{2}$ norm is a true norm, a large displacement norm must imply a large ${\mathcal{L}}_{2}$ norm.

Figure 8.

Schematic of base vortices (contour lines) superposed with a vorticity perturbation (colour field). The points represent the centroid of the superposition of the base flow and a finite-amplitude perturbation. (a) The perturbation does not change the centroid location, hence the displacement norm for this perturbation is zero. However, the perturbation is not identically zero, and neither is its ${\mathcal{L}}_{2}$ norm. (b) The perturbation moves the centroid, and thus both norms are non-zero in this case.

Eventually, because of the definition of this norm that focuses on a given vortex $i$, one displacement norm can be defined for each vortex ring. In our study, we assume the linear behaviour to be only valid in the vicinity of the inlet, where perturbations remain small; therefore, we only consider the displacement norm for each of the first three vortex rings.

3.4.1 Length of the optimisation box for the energy norm

For the ${\mathcal{L}}_{2}$ norm, the optimisation is carried out by calculating the response norm on the domain $0\leqslant r\leqslant R_{m}$ and $0\leqslant z\leqslant Z_{opt}$, with $Z_{opt}=5$. However, changing $Z_{opt}$ could change the results. In this section, we investigate the effect of changing $Z_{opt}$, in the $(Re=2000,~St=0.5)$ regime, by considering values of $Z_{opt}$ from 3 to 7.5 separated by 0.5 increments. When all the optimal forcings are calculated, scalar products between each couple of them are subsequently computed, and the most distant optimal forcings are those for $Z_{opt}=3$ and $Z_{opt}=7.5$, with a scalar product of 0.99961. The optimal forcing is thus independent of the domain length considered.

3.4.2 Length of the boxes around vortices for displacement norm

For the displacement norm, the optimal forcing is calculated by integrating perturbation quantities around the $i\text{th}$ vortex core ($i=1$, 2 or 3). To perform this integration, a rectangular box $[Z_{i}-\unicode[STIX]{x1D6FF}_{z}/2;Z_{i}+\unicode[STIX]{x1D6FF}_{z}/2]\times [0;1]$, with $Z_{i}$ the axial position of the vortex (calculated as a centroid) and $\unicode[STIX]{x1D6FF}_{z}$ fixed to 0.8, is first chosen. However, the choice of $\unicode[STIX]{x1D6FF}_{z}$ could affect the results of the optimisation.

This effect is now studied for the $Re=2000$ and $St=0.5$ case. Four different values of $\unicode[STIX]{x1D6FF}_{z}$ are tested: 0.4, 0.6, 0.8 and 1, as represented in figure 7. Above $\unicode[STIX]{x1D6FF}_{z}=1$, boxes for consecutive vortices would overlap. The first three vortices are considered here.

The results show that the optimal forcing barely depends on the box size: by considering scalar products between optimal forcing structures on the same vortex for different $\unicode[STIX]{x1D6FF}_{z}$, all scalar products are above 0.992.

3.4.3 Effect of the dimension of the forcing space

The convergence of the optimal forcing shape with respect to the number $N$ of basis functions is now evaluated. As these functions are orthonormal, the optimal forcing for a basis of dimension $N^{\prime }<N$ can easily be computed by taking the first $N^{\prime }$ components of the optimal $N$-forcing. Therefore, the evolution of energy component of the optimal forcing over each basis function, as represented in figure 16, shows the convergence of the algorithm as a function of $N$, for $Re=2000$ and $St=0.5$.

Figure 16.

Energy component of the optimal forcing on each of the different basis forcing functions considered, for $Re=2000$ and $St=0.50$. In this case, $N=50$.

It can be stated that the quickest mode to converge are the $k=0$, 1, 3 and 4 modes, corresponding to forcing in the radial and azimuthal directions. However, in all cases, convergence is properly achieved, with energy levels contained in the last modes about 10 000 times lower than the ones in the first modes.

3.4.4 Effect of the time-shift with respect to the axisymmetric forcing

Our optimisation over a periodic flow is based on cost function evaluations only once per subharmonic cycle, $T_{o}=2nT_{f}$. However, a large gain at these synchronized $T_{o}$ does not necessarily imply large gain for different $T_{o}=(2n+\unicode[STIX]{x1D711})T_{f}$ phase-shifted of $2\unicode[STIX]{x03C0}\unicode[STIX]{x1D711}$ with the fundamental forcing; in other words, there is no guarantee that a forcing that shifts well the vortices at a given objective time will not perform poorly at another time.

As a test, we evaluate the influence of a few optimisation time-shifts: $\unicode[STIX]{x1D711}=0$ (phase 0), 1/4 (phase $\unicode[STIX]{x03C0}/2$), 1/2 (phase $\unicode[STIX]{x03C0}$) and 3/4 (phase $3\unicode[STIX]{x03C0}/2$). The different optimal ${\mathcal{L}}_{2}$-norm forcing functions found with these shifts are compared in table 4. Again, the optimal forcing is almost independent of the phase shift. This can be understood using the same arguments of the large gain separation as in § 3.3, but also physically by considering the jet bifurcation as a very receptive physical mechanism: in the interpretation of Parekh et al. (Reference Parekh, Leonard and Reynolds1988), the vortices are shed out of the axis, and they do not oscillate around a fixed equilibrium. Therefore, there is no chance that at a particular time-shift, the forcing has no effect on the vortices. If, at some time-step in the forcing period, the vortices are offset, then this offset continues through the entire cycle.

Table 4.

Scalar product of the ${\mathcal{L}}_{2}$ norm optimal forcings for different phase shifts. Here $T_{o}=(2n+\unicode[STIX]{x1D711})T$ with $2n=18$ for phase 0 and 6 for others, with $Re=2000$ and $St=0.50$.

3.4.5 Final objective time $T_{o}$

The optimal forcing and response can be calculated for different objective times $T_{o}$. The shape of forcing functions and the associated gain evolves with time. At some point, the optimal forcing does not change and the gain reaches a plateau. In this study, the convergence condition we set for the final state is a variation of the gain by less than 0.5 % compared with the gain achieved two forcing periods before. In this section, we discuss the evolution of these quantities for the ${\mathcal{L}}_{2}$ norm and the $St=0.50$ case.

Figure 17(a) shows that, after an initial phase of exponential growth, the gain reaches a plateau. This confirms the stability of the axisymmetric unpaired base flow with respect to $m=1$ perturbations.

Figure 17(b) displays the shape evolution with $T_{o}$ of the optimal forcing and response by comparing each of these with the optimal quantities for $T_{o}=2$ and 18. Regarding forcing structures, little change occurs between $T_{o}=2$ and 18; the initial forcing is almost optimal, its shape is fixed. Therefore, if one is only interested in optimal forcing, there is no need to run the optimisation for long $T_{o}$. Concerning response structures, because of the spatial development of the perturbation field (more than two forcing periods are necessary to fill the five-diameter-long optimisation domain, because the initial perturbation is zero), the initial shape of the response, at $T_{o}=2$, has little resemblance with its final shape. After six forcing periods, they are almost identical.

The gain values obtained with different time shifts over six flow periods, as discussed in § 3.4.4, are shown as red crosses in figure 17(a). These values also give a measure of the moderate gain variations, within 12 % of their mean, that can be expected for time horizons in between integer multiples of the forcing period.

Figure 17.

Evolution of optimal forcing, response and gain as a function of optimisation time $T_{o}$ for the $Re=2000$, $St=0.50$ case. All quantities are optimised with respect to the ${\mathcal{L}}_{2}$ norm. (a) Evolution of the optimal gain. Black circles denote values obtained for the standard cases without phase shift between harmonic and subharmonic forcing, red crosses represent the gain obtained with the phase shifts discussed in § 3.4.4. (b) Evolution of the optimal forcing and response, with scalar product comparisons with $T_{o}=2$ and 18. The scalar products are all carried out over quantities normalized with respect to the relevant ${\mathcal{L}}_{2}$ norm.

4.1.1 Amplitude and phase of the optimal forcing

Adding the optimal helical forcing, with an amplitude of $B\unicode[STIX]{x1D6E4}$, the jet inlet velocity (2.3) is prescribed as

(4.4)$$\begin{eqnarray}\displaystyle \boldsymbol{u}(r,\unicode[STIX]{x1D703},z=0,t) & = & \displaystyle U_{0}(r)(1+A\cos (\unicode[STIX]{x1D714}_{f}t))\boldsymbol{e}_{\boldsymbol{z}}\nonumber\\ \displaystyle & & \displaystyle +\,2B\unicode[STIX]{x1D6E4}\left[\text{Re}(\boldsymbol{u}_{\mathbf{1}\boldsymbol{f},\mathbf{1}})\cos (\unicode[STIX]{x1D703})\cos (\unicode[STIX]{x1D714}_{f}t/2)-\text{Im}(\boldsymbol{u}_{\mathbf{1}\boldsymbol{f},\mathbf{1}})\sin (\unicode[STIX]{x1D703})\cos (\unicode[STIX]{x1D714}_{f}t/2)\right.\nonumber\\ \displaystyle & & \displaystyle +\left.\text{Re}(\boldsymbol{u}_{\mathbf{1}\boldsymbol{f},\mathbf{2}})\cos (\unicode[STIX]{x1D703})\sin (\unicode[STIX]{x1D714}_{f}t/2)-\text{Im}(\boldsymbol{u}_{\mathbf{1}\boldsymbol{f},\mathbf{2}})\sin (\unicode[STIX]{x1D703})\sin (\unicode[STIX]{x1D714}_{f}t/2)\right]\!.\end{eqnarray}$$

A normalisation factor $\unicode[STIX]{x1D6E4}$ is introduced, with a numerical value determined as 0.239 for the present flow case, such that the energy ratio between axisymmetric and helical forcing is given by $A^{2}/B^{2}$.

Without loss of generality, the azimuthal phase of the bifurcation plane is set to $\unicode[STIX]{x1D703}=0$ by imposing all coefficients $u_{2f,\cdot ,\cdot }$ to be real-valued functions of $r$. This is possible because the perturbation equations for $\boldsymbol{u}_{\mathbf{2}}$ only involve real coefficients (see appendix B). Therefore, the total inflow forcing is finally recovered as

(4.5)$$\begin{eqnarray}\displaystyle \boldsymbol{u}(r,\unicode[STIX]{x1D703},0,t) & = & \displaystyle (1+A\cos (\unicode[STIX]{x1D714}_{f}t))U_{0}\boldsymbol{e}_{\boldsymbol{z}}\nonumber\\ \displaystyle & & \displaystyle +\,2B\unicode[STIX]{x1D6E4}\left[\begin{array}{@{}c@{}}(u_{2f,1,r}\cos (\unicode[STIX]{x1D714}_{f}t/2)+u_{2f,2,r}\sin (\unicode[STIX]{x1D714}_{f}t/2))\cos (\unicode[STIX]{x1D703})\\ (u_{2f,1,\unicode[STIX]{x1D703}}\cos (\unicode[STIX]{x1D714}_{f}t/2)+u_{2f,2,\unicode[STIX]{x1D703}}\sin (\unicode[STIX]{x1D714}_{f}t/2))\sin (\unicode[STIX]{x1D703})\\ (u_{2f,1,z}\cos (\unicode[STIX]{x1D714}_{f}t/2)+u_{2f,2,z}\sin (\unicode[STIX]{x1D714}_{f}t/2))\cos (\unicode[STIX]{x1D703})\end{array}\right].\end{eqnarray}$$

In this study, the amplitude $B$ is fixed at 0.01, five times lower than the amplitude of the axisymmetric forcing prescribed; if this helical forcing had a constant shape in the radial direction harmonically modulated with time (as the axisymmetric), this would have corresponded to a maximum velocity amplitude with respect to the infinity norm of 1 % of the mean flow. This amplitude is much lower than what has been used in the literature: 5 % of the mean flow in Gohil et al. (Reference Gohil, Saha and Muralidhar2015), 1 % to 15 % in Tyliszczak & Geurts (Reference Tyliszczak and Geurts2014) or 15 % in Danaila & Boersma (Reference Danaila and Boersma1998, Reference Danaila and Boersma2000).

4.1.2 Amplitude of the ad hoc forcing

The simple, non-optimal helical forcing (3.3) is accordingly cast in normalised form as

(4.6)$$\begin{eqnarray}\displaystyle \boldsymbol{u}(r,\unicode[STIX]{x1D703},z=0,t)=U_{0}(r)\left(1+A\cos (\unicode[STIX]{x1D714}_{f}t)+B^{\prime }\unicode[STIX]{x1D6E4}^{\prime }\cos (\unicode[STIX]{x1D714}_{f}t/2)\cos (\unicode[STIX]{x1D703})\left(\frac{2r}{D}\right)\right)\boldsymbol{e}_{\boldsymbol{z}}, & & \displaystyle\end{eqnarray}$$

with a fixed zero-phase shift. A value $\unicode[STIX]{x1D6E4}^{\prime }=1.850$ gives an energy ratio $A^{2}/B^{\prime 2}$ between axisymmetric and helical components, which is consistent with the above parametrisation of the optimal forcing. To enable fair comparisons, $B^{\prime }$ is fixed at 0.01, the same as $B$ in the previous section, and much lower than forcing levels commonly used in the literature.