3.2.1 Formulation

Let the output response, $\unicode[STIX]{x1D643}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D714},m}$ , where $\unicode[STIX]{x1D643}$ is an observation operator, be rewritten in terms of the forcing as follows:

where $R_{\bar{\boldsymbol{q}},\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D714},m}=\unicode[STIX]{x1D643}{\mathcal{L}}_{\bar{\boldsymbol{q}},\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D714},m}^{-1}\unicode[STIX]{x1D63D}$ is the resolvent operator (Schmid Reference Schmid2007). Once the problem is discretised, the operators become matrices, and will be considered as such in what follows. The goal of the analysis is to determine the forcing such that the norm of the associated response is maximised,

where $\ast$ applied to matrices implies conjugate transpose, the weighted norm, $\Vert \cdot \Vert _{\unicode[STIX]{x1D652}}$ , is associated with the inner product, $(\cdot ,\unicode[STIX]{x1D652}\cdot )$ , with $\unicode[STIX]{x1D652}$ a definite positive Hermitian matrix, $(\cdot \,,\cdot )$ the Euclidean inner product and $\Vert \cdot \Vert$ the Euclidean norm. For simplicity, it is considered that the weighted norm is the same in the input and output spaces. The maximisation of the Rayleigh quotient (3.4) can be achieved by means of the following singular value decomposition:

where $\unicode[STIX]{x1D652}^{1/2}$ is defined by the Cholesky decomposition, $\unicode[STIX]{x1D652}=\unicode[STIX]{x1D652}^{1/2}(\unicode[STIX]{x1D652}^{1/2})^{\ast }$ , $\unicode[STIX]{x1D650}=(\unicode[STIX]{x1D650}_{1},\ldots ,\unicode[STIX]{x1D650}_{N})$ and $\unicode[STIX]{x1D651}=(\unicode[STIX]{x1D651}_{1},\ldots ,\unicode[STIX]{x1D651}_{N})$ are orthonormal matrices and $\unicode[STIX]{x1D72E}$ is a diagonal matrix. We define the optimal forcing modes (OF) as $\unicode[STIX]{x1D731}_{i}=(\unicode[STIX]{x1D652}^{-1/2})^{\ast }\unicode[STIX]{x1D651}_{i}$ and the associated optimal response modes (OR) as $\unicode[STIX]{x1D733}_{i}=(\unicode[STIX]{x1D652}^{-1/2})^{\ast }\unicode[STIX]{x1D650}_{i}$ . The Hankel singular values, $\unicode[STIX]{x1D70E}_{i}$ , diagonal elements of $\unicode[STIX]{x1D72E}$ sorted in descending order, indicate the gains associated with forcing–response mode pairs. Indeed, we obtain the relationship

which shows that the optimal response is that of the optimal forcing amplified by the corresponding Hankel singular value $\unicode[STIX]{x1D70E}_{i}$ .

3.2.2 Resolvent analysis in an axisymmetric compressible jet

Resolvent analysis is applied in a locally parallel jet configuration and results are compared with measurements. As discussed previously, linear wavepackets computed using the homogeneous linearised system closely match measurements performed in an $M=0.4$ jet up to the end of the potential core, downstream of which the agreement deteriorates. We use resolvent analysis to explore whether this discrepancy might be remedied by considering optimally forced wavepackets.

In performing the resolvent analysis, we assume the jet flow to be viscous, compressible, isothermal and axisymmetric in a locally parallel setting. The base flows are taken from the measurements of Cavalieri et al. (Reference Cavalieri, Rodriguez, Jordan, Colonius and Gervais2013). In order to numerically converge the resolvent analysis, a lower Reynolds number than that of the experiment must be used (Meseguer & Trefethen Reference Meseguer and Trefethen2003); this is chosen as $Re=10\,000$ . As discussed in § 3.1, this viscosity regularises the solution, avoiding infinite responses. Due to the locally parallel flow assumption, only the axial component of the base-flow velocity is used, and the thermodynamic variables are determined using the Crocco–Busemann relation and the perfect gas law, thus defining the base flow as $\bar{\boldsymbol{q}}=(\bar{U},0,\bar{T},\bar{\unicode[STIX]{x1D70C}})^{\text{T}}$ . The viscosity $\bar{\unicode[STIX]{x1D707}}$ is determined using Sutherland’s law. The bulk viscosity is taken as $\unicode[STIX]{x1D707}_{B}=0.6\unicode[STIX]{x1D707}$ for air in accordance with Freund (Reference Freund1997). The equations are supplemented with homogeneous Dirichlet boundary conditions in the far field, as done in several previous studies (Hanifi, Schmid & Henningson Reference Hanifi, Schmid and Henningson1996; Gloor, Obrist & Kleiser Reference Gloor, Obrist and Kleiser2013). Application of a non-reflecting boundary condition is more accurate, using asymptotic (Gudmundsson Reference Gudmundsson2010) or characteristic methods (Poinsot & Lele Reference Poinsot and Lele1992). Nonetheless, the exponentially decaying modes in a sufficiently large computation domain will generate only negligible reflected waves. The full set of equations in cylindrical coordinates under such assumptions is provided in appendix A.

The equations are discretised and solved on a Chebyshev grid using a spectral collocation method (Weideman & Reddy Reference Weideman and Reddy2000) in the radial direction. To compute the singular value decomposition of the resolvent operator, one must define the matrices $\unicode[STIX]{x1D63D}$ and $\unicode[STIX]{x1D643}$ , and the norm $\Vert \cdot \Vert _{\unicode[STIX]{x1D652}}$ . We make the choice of restricting the input–output spaces to clarify physical interpretation of the results. Given the relatively low Mach number, nonlinearities of the continuity equation can be neglected (see the low-Mach-number limit in Alazard Reference Alazard2006). Moreover, since we are studying a cold jet, we can also neglect the nonlinear terms arising in the energy equation, such as the transport of the fluctuating temperature by the velocity fluctuations. We thus constrain the nonlinearities to act solely on the momentum equation. Moreover, in order to simplify interpretation, we focus on outputs in the momentum equations also, which leads to the matrices

Numerical tests show that the results remain qualitatively the same if the output matrix $\unicode[STIX]{x1D643}$ is the identity matrix. We recall that, despite the simplifications in the input/output behaviour, the full linearised equations for compressible flow are used.

Figure 5.

Characterisation of the optimal forcing–response modes. (a) Optimal response energy gain $\unicode[STIX]{x1D70E}_{i}^{2}$ for all axial positions for $St=0.6$ . (b) Optimal forcing/response at $x/D=8$ . For clarity, the forcing modes $\unicode[STIX]{x1D731}_{i}$ , in solid lines, have been multiplied by $10^{-3}$ . The associated responses $\unicode[STIX]{x1D70E}_{i}\unicode[STIX]{x1D733}_{i}$ are in dashed lines. Vertical solid: critical layer, vertical dashed: inflection point.

For compressible flows, as discussed in Hanifi et al. (Reference Hanifi, Schmid and Henningson1996) and Rowley, Colonius & Murray (Reference Rowley, Colonius and Murray2004), there exists no unique choice for the norm definition. Since we have restricted the inputs and outputs to the momentum equations, these norms simplify to the standard incompressible energy norm; $\unicode[STIX]{x1D652}$ is then defined such that $\Vert \unicode[STIX]{x1D643}\widetilde{\boldsymbol{q}}\Vert _{\unicode[STIX]{x1D652}}^{2}=\int _{\unicode[STIX]{x1D6FA}}(\bar{\unicode[STIX]{x1D70C}}/2)(\widetilde{u}^{\ast }\widetilde{u}+\widetilde{v}^{\ast }\widetilde{v})r\,\text{d}r$ .

In what follows, the phase velocity, $c$ , is that predicted by a PSE computation for the same jet, detailed in § 4.1.1. We have taken a real-valued $c$ , obtained as $\unicode[STIX]{x1D714}/\text{Re}(\unicode[STIX]{x1D6FC})$ , $\unicode[STIX]{x1D6FC}$ being the local wavenumber in the PSE computation. The justification of this choice is given in appendix B. The axisymmetric mode, $m=0$ , is considered, as it is the most acoustically efficient (Cavalieri et al. Reference Cavalieri, Jordan, Colonius and Gervais2012, Reference Cavalieri, Rodriguez, Jordan, Colonius and Gervais2013), and $St=0.6$ is arbitrarily chosen to illustrate the salient features of forced-wavepacket behaviour.

The resolvent analysis identifies a basis of forcing modes spanning $\text{Ker}(\unicode[STIX]{x1D643}{\mathcal{L}}_{\bar{\boldsymbol{q}},\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D714},m}^{-1}\unicode[STIX]{x1D63D})^{\text{T}}$ and a basis of response modes spanning $\text{Im}(\unicode[STIX]{x1D643}{\mathcal{L}}_{\bar{\boldsymbol{q}},\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D714},m}^{-1}\unicode[STIX]{x1D63D})$ . Figure 5(a) shows the optimal energy gain $\unicode[STIX]{x1D70E}_{i}^{2}$ of the five first modes as a function of axial position. This figure shows that for all axial positions, the first optimal forcing–response mode pair clearly dominates. This suggests that the optimal forcing–response mode pair is the most relevant for the dynamics of wavepackets forced by background turbulence. This idea is supported by the recent work of Beneddine et al. (Reference Beneddine, Sipp, Arnault, Dandois and Lesshafft2016), which shows that when the first optimal energy gain is clearly higher than the others and if the turbulent forcing does not project preferentially onto a suboptimal mode, the response is dominated by the first optimal response mode. Figure 5(b) shows the radial structure of the first four optimal forcing modes at $x/D=8$ and of the associated responses. Beyond the fact that the first optimal response mode dominates by one order of magnitude, the shape of the forcing and the associated response are typical of the critical-layer mechanism highlighted in the direct forcing study in Tissot et al. (Reference Tissot, Zhang, Lajús, Cavalieri, Jordan and Colonius2015). While the forcing has broad radial support, its maximum is at the critical layer and the response has here one dip and a phase jump of approximately $\unicode[STIX]{x03C0}$ (phases are shown in appendix C, figure 22), consistent with the shape of critical-layer modes, and indicative of a vortex-like structure, where axial velocity fluctuations at either side of the vortex are in phase opposition, and zero amplitude is found at the vortex centre (see the discussion in Cavalieri et al. Reference Cavalieri, Rodriguez, Jordan, Colonius and Gervais2013). The higher modes show several dips and can be interpreted as functions that span the response to less relevant turbulent structures. The dominance of the first mode pair is encouraging vis-à-vis an eventual reduced-order model, and so in what follows we consider only this first mode pair.

Figure 6.

Comparison of optimal forcing modes (red dashed lines) and optimal response modes (red solid lines) of axial velocity fluctuation at $x/D=3$ (a) and $x/D=8$ (b) with experimental data (blue squares) and the solution of the homogeneous linear problem (LEE, green line with solid dots; PSE, black dashed line). Here, $St=0.6$ , $Re=10\,000$ , $M=0.4$ . The phase velocity $c$ predicted by the PSE is used. The vertical dashed line is the critical-layer position and the vertical dash-dotted line is the inflectional point position.

In figure 6, we show the power spectral density (PSD) of the axial perturbation velocity obtained by the resolvent analysis compared with time-resolved particle image velocimetry (PIV) measurements (TRPIV), and also with two linear wavepacket models, solutions of either the linearised Euler equations (LEE) or the PSE. The LEE solutions are taken from Zhang et al. (Reference Zhang, Jordan, Lehnasch, Cavalieri and Agarwal2014). We focus the comparison on two emblematic positions: $x/D=3$ , which typifies results upstream of the end of the potential core, where the homogeneous linear framework matches the data closely, and $x/D=8$ , typifying results obtained downstream of the potential core, where linear theory fails. Results for other axial stations are shown in appendix C, figure 21.

In the experimental results, upstream of $x/D\sim 5$ , the first proper orthogonal decomposition (POD) mode of the TRPIV data for $m=0$ shows that the turbulent jet flow comprises a vortex-like signature in this region (Cavalieri et al. Reference Cavalieri, Rodriguez, Jordan, Colonius and Gervais2013). We can also see that, in this region, the linear wavepacket models (PSE and LEE) correctly capture this structure, as these curves collapse with the TRPIV data. This indicates that the homogeneous linear wavepacket model is able to predict the behaviour of the K–H instability in the upstream part of the jet, which undergoes exponential growth in $x$ . The resolvent optimal response modes, although corresponding to a forced-wavepacket scenario, clearly dynamically different from the homogeneous linear wavepackets, also match the measurement upstream of $x/D\sim 5$ . The similarities between the results of the homogeneous and the forced problems suggest that in the upstream region resolvent analysis leads to a most amplified response corresponding to a pseudo-resonance of the Kelvin–Helmholtz mode; the flow response is approximately the same for both upstream excitation (in the homogeneous problem) and forcing by nonlinearities (in the forced problem).

Downstream of the end of the potential core, on the other hand, the homogeneous linear models (PSE or LEE) predict a quite different radial structure, as can be seen at $x/D=8$ , figure 6, or in figure 21. The measurement shows that the flow field in these regions retains a vortex-like structure possessing a local dip in amplitude at approximately $r/D\sim 0.4$ (this dip in amplitude is a signature of phase opposition of disturbances on either side of the lipline, as previously discussed) which is not captured by the linear homogeneous wavepacket model. When forcing by nonlinearities is considered in the resolvent analysis, the optimal response displays a radial structure similar to that of the experimentally educed wavepacket. This is an encouraging indication that the notion of forced linear wavepackets may correspond to what is happening in the flow downstream of the potential core.

We note here that the position of the dip does not always correspond exactly with experimental results. This position is sensitive to the perturbation phase velocity $c$ , which fixes the critical-layer position. The forced wavepacket, which we assume to be present in experiments, may have a slightly different phase velocity from the unforced wavepacket, explaining the observed discrepancies. Indeed, as shown in appendix B, the agreement at $x/D=8$ can be slightly improved with a modification of $c$ . Selection of $c$ from linear PSE results remains preferable, as this is a choice based on first principles and not a posteriori.

The comparison can be extended to other Strouhal numbers and axial stations by means of an agreement metric between any pair taken from TRPIV, PSE and OR; we thus define the metric $\unicode[STIX]{x1D6FD}$ (shown here between PSE and TRPIV, for example),

where the ${\mathcal{Y}}$ are the PSD of velocity components from experimental measurements and

$r_{b}$ and $r_{t}$ defining the integration range for each variable. For the axisymmetric mode considered, only $\widetilde{u}$ and $\widetilde{v}$ are involved in the calculation, since these are the quantities available from the experimental results. The results for the absolute values of $\unicode[STIX]{x1D6FD}$ ( $|\unicode[STIX]{x1D6FD}|$ ) between TRPIV and PSE or OR are shown in figure 7. When two quantities have similar radial structure, the value of $|\unicode[STIX]{x1D6FD}|$ is close to unity.

In figure 7(a), we see that PSE matches the experimental data well upstream of the end of the potential core, especially for $St=0.3{-}0.9$ , and the agreement deteriorates in the downstream region. This result has been noted in Cavalieri et al. (Reference Cavalieri, Rodriguez, Jordan, Colonius and Gervais2013), for example. For the comparison between TRPIV and OR, we see in figure 7(b) that the values of the $\unicode[STIX]{x1D6FD}$ metric are globally higher than those obtained for TRPIV and PSE. Moreover, $|\unicode[STIX]{x1D6FD}|$ does not decay in the downstream region as significantly as in figure 7(a). This result has already been inferred and discussed from figure 6 above, but figure 7 is a useful means by which to compress the comparisons for all frequencies and axial positions into a single figure.

Figure 7.

Values of $|\unicode[STIX]{x1D6FD}|$ . Here, $c$ is determined from the PSE calculation.

Figure 8(a) shows a map of the optimal forcing mode $\unicode[STIX]{x1D731}_{1}$ , whose norm is fixed to 1. With this arbitrary normalisation we cannot evaluate the downstream evolution of forced wavepackets; the purpose of the plot is to examine the radial structure of optimal forcing and response for different axial stations. It is clear that for all $x/D$ the optimal forcing is concentrated around the critical layer. The associated response, $\unicode[STIX]{x1D70E}_{1}\unicode[STIX]{x1D733}_{1}$ , is shown in figure 8(b). The response also follows the critical layer, with a dip in amplitude consistently close to the critical-layer position, and a clear signature of a vortex-like structure. A maximum response appears at $x/D=3$ , in keeping with the peak at that position of the optimal energy gain in figure 5(a).

Figure 8.

The PSD of the first optimal forcing and response modes along the axial direction at $St=0.6$ , with $c$ determined by PSE. The solid line is the critical-layer position and the dashed line is the inflection point position.

Finally, spatial reconstructions of the optimal forcings and responses (figure 9) at $x/D=3$ and $x/D=8$ show that the optimal forcing is inclined against the shear, suggesting an exploitation of non-normality by the Orr mechanism.

Figure 9.

Spatial reconstruction of the streamwise component of the optimal response (colour) and of the optimal forcing (isolines, solid for positive, dashed for negative).

To summarise, the first optimal forcing mode dominates the higher-order ones and shows a maximum at the critical layer for all axial positions. The response is consistent with the observed flow behaviour downstream of the end of the potential core, supporting the idea that a forced-wavepacket model may be what is needed to describe features present in this region of the jet. Finally, a spot of high sensitivity to forcing is found at $x/D=3$ for $St=0.6$ , which is where the K–H mode becomes neutral and thus close, in the complex plane, to the real wavenumber associated with the forcing. Because of the locally parallel assumption of the resolvent analysis, conclusions about relationships between different axial positions can be difficult to draw. In § 4.1, we therefore complement the study with a non-parallel approach, PSE-4D-Var, which accommodates the true structure of the mean flow and permits a more complete consideration of the spatially convective nature of the flow.

4.1.1 Formulation

The PSE (Herbert Reference Herbert1997) allow instability wave evolution over a slowly varying base flow to be modelled. The formulation used here is similar to that used in previous works (Gudmundsson Reference Gudmundsson2010; Cavalieri et al. Reference Cavalieri, Rodriguez, Jordan, Colonius and Gervais2013; Sasaki Reference Sasaki2015). The assumption of a slowly varying base flow allows (Crighton & Gaster Reference Crighton and Gaster1976) decomposition of the perturbation associated with mode $(\unicode[STIX]{x1D714},m)$ into slowly and rapidly varying (wave-like) parts (Herbert Reference Herbert1997),

Here, $\boldsymbol{q}(x,r)$ is the slowly varying part and $\text{e}^{\text{i}\int _{0}^{x}\unicode[STIX]{x1D6FC}(\unicode[STIX]{x1D709})\,\text{d}\unicode[STIX]{x1D709}}$ is the wave-like part. The $(\unicode[STIX]{x1D714},m)$ index is dropped for compactness. The decomposition (4.1) is introduced into the compressible Navier–Stokes equations. The nonlinear terms are then neglected, assuming small perturbations, and the first axial derivatives of $\unicode[STIX]{x1D6FC}$ and second axial derivatives of $\boldsymbol{q}$ are neglected, assuming a slow variation of these variables in the streamwise direction. We thus obtain an equation of the form

where the matrix $\unicode[STIX]{x1D63C}$ contains discrete radial derivatives. The domain, $x/D=[0,L]$ , $r/D=[0,\infty ]$ , is discretised using $N_{x}$ points in the axial direction and $N_{r}$ Chebyshev collocation points in the radial direction. Details of derivation of the PSE for compressible flow in cylindrical coordinates are given in Sasaki (Reference Sasaki2015).

The decomposition (4.1) is a priori not unique because of the presence of the $x$ variable in $\boldsymbol{q}(x,r)$ as well as in $\unicode[STIX]{x1D6FC}(x)$ , and because no evolution equation is given for $\unicode[STIX]{x1D6FC}$ . In order to overcome this ambiguity, a normalisation constraint is added (see Herbert Reference Herbert1997) such that the exponential dependence (real and imaginary) is absorbed by the wave-like term $\text{e}^{\text{i}\int _{0}^{x}\unicode[STIX]{x1D6FC}(\unicode[STIX]{x1D709})\,\text{d}\unicode[STIX]{x1D709}}$ :

Starting with $(\boldsymbol{q}_{0},\unicode[STIX]{x1D6FC}_{0})$ , often taken as the K–H mode, solution of the locally parallel problem at the initial position, $\boldsymbol{q}(x,r)$ , can be obtained by integrating (4.2) in space using an implicit Euler scheme, with iterations in $\unicode[STIX]{x1D6FC}$ such that the constraint (4.3) is respected (Gudmundsson Reference Gudmundsson2010). The linear PSE procedure allows wavepackets for each frequency–azimuthal mode pair to be propagated independently.

In the 4D-Var under weak dynamical constraint, we would like the solution, $\boldsymbol{q}(x,r)$ , to follow the observations, and we consider that this solution respects the model – here the linear PSE – only in a weak sense. In order to relax this constraint, we introduce a forcing term, $\boldsymbol{f}$ , interpreted as the effect of nonlinearities, and an initial condition perturbation, $\boldsymbol{\unicode[STIX]{x1D702}}$ . Combining the nonlinear forcing with (4.2) and (4.3), the forced parabolised stability equations can be written as

where $\unicode[STIX]{x1D652}_{r}$ is the discrete weight matrix associated with the inner product $(\cdot \,,\cdot )_{r}$ . It should be noted that we implicitly enforce the shape of the forcing as

and we can define a forcing in physical space as

This means that because a harmonic framework is used, the optimal forcing is constrained to have a wavepacket shape. This shape is indeed able to act efficiently on the instability wave because a forcing whose wavenumber matches that of the wavepacket should lead to the most significant dynamic effects. In principle, $\boldsymbol{f}$ could come, for instance, from a nonlinear PSE computation taking account of a reduced number of wave–wave interactions. In the present study, all possible nonlinear interactions are potentially considered via the inverse problem.

The objective of 4D-Var is to seek $(\boldsymbol{f},\boldsymbol{\unicode[STIX]{x1D702}})$ minimising

where $\mathbb{H}(\boldsymbol{q},\unicode[STIX]{x1D6FC})$ is an observation operator allowing passage from the state space to the observation space. The observation operator $\mathbb{H}(\boldsymbol{q},\unicode[STIX]{x1D6FC})$ is related to the output matrix $\unicode[STIX]{x1D643}$ , introduced in § 2, since both define the output space. Here, it is simply generalised to a nonlinear operator. The ${\mathcal{Y}}$ are associated observations, i.e. experimental measurements. A standard sensitivity problem is achieved choosing

For the 4D-Var problem, we choose

and ${\mathcal{Y}}$ as the PSD of streamwise and radial velocity measurements. Here, $\unicode[STIX]{x1D64C}$ is a rectangular matrix of zeros and ones. Similarly to the output matrix $\unicode[STIX]{x1D643}$ in (3.7), rows of $\unicode[STIX]{x1D64C}$ contain a ‘one’ in order to select the components and the radial positions that are observed, here the axial and radial velocities for the positions where there are available TRPIV results. The terms $\boldsymbol{q}_{L}$ , $\unicode[STIX]{x1D6FC}_{L}$ and ${\mathcal{Y}}_{L}$ refer to the values of the variables or parameters at $x/D=L$ , where $L$ is the downstream end of the computational domain; we penalise errors in the terminal condition in order to avoid transient behaviour at the end of the domain of the adjoint equation, presented later.

The norms $\Vert \cdot \Vert _{\unicode[STIX]{x1D652}_{o}}$ , $\Vert \cdot \Vert _{\unicode[STIX]{x1D652}_{T}}$ , $\Vert \cdot \Vert _{\unicode[STIX]{x1D652}_{f}}$ and $\Vert \cdot \Vert _{\unicode[STIX]{x1D652}_{\unicode[STIX]{x1D702}}}$ are respectively associated with the inner products, $(\unicode[STIX]{x1D64C}^{\ast }\cdot ,\unicode[STIX]{x1D64C}^{\ast }\unicode[STIX]{x1D652}_{o}\cdot )_{r}$ , $(\unicode[STIX]{x1D64C}^{\ast }\cdot ,\unicode[STIX]{x1D64C}^{\ast }\unicode[STIX]{x1D652}_{T}\cdot )_{r}$ , $(\cdot ,\unicode[STIX]{x1D652}_{f}\cdot )_{r}$ and $(\cdot ,\unicode[STIX]{x1D652}_{\unicode[STIX]{x1D702}}\cdot )_{r}$ . The weight matrices, $\unicode[STIX]{x1D652}_{o}$ , $\unicode[STIX]{x1D652}_{T}$ , $\unicode[STIX]{x1D652}_{f}$ and $\unicode[STIX]{x1D652}_{\unicode[STIX]{x1D702}}$ , allow different degrees of confidence to be specified in the observations or in the PSE model. Finally, for later use, we define, for any $\unicode[STIX]{x1D648}\in \mathbb{R}^{N_{c}N_{r}\times N_{c}N_{r}}$ , the adjoint matrix $\unicode[STIX]{x1D648}^{+}=\unicode[STIX]{x1D652}_{r}^{-1}\unicode[STIX]{x1D648}^{\ast }\unicode[STIX]{x1D652}_{r}$ , such that $(\unicode[STIX]{x1D648}\boldsymbol{x},\boldsymbol{y})_{r}=(\boldsymbol{x},\unicode[STIX]{x1D648}^{+}\boldsymbol{y})_{r}$ .

4.1.2 Optimality system

The minimisation problem exposed in § 4.1.1 is a constrained optimisation problem. It consists in minimising (4.7) under the constraint (4.4). This problem is usually solved using a quasi-Newton algorithm, which requires the gradient of the cost functional, ${\mathcal{J}}$ , with respect to the control parameters, $(\boldsymbol{f},\boldsymbol{\unicode[STIX]{x1D702}})$ . Because the control parameter space is of large dimension, we use an adjoint method, where the gradient of the cost functional ${\mathcal{J}}$ , also called sensitivity, is determined by the integration of an adjoint equation (Gunzburger Reference Gunzburger2003). The direct and adjoint systems are then solved iteratively until convergence to the optimal solution.

The continuous adjoint PSE can be determined in different ways. The first consists in parabolising the adjoint LEE, as done by Dobrinsky (Reference Dobrinsky2002) for a receptivity problem. The second is to directly derive the adjoint PSE as in Pralits et al. (Reference Pralits, Airiau, Hanifi and Henningson2000) and Ansaldi & Airiau (Reference Ansaldi and Airiau2015) for a sensitivity analysis or as in Pralits, Hanifi & Henningson (Reference Pralits, Hanifi and Henningson2002) and Airiau et al. (Reference Airiau, Bottaro, Walther and Legendre2003) for flow control. For reasons of consistency between PSE and adjoint PSE, the latter method has been chosen.

In order to derive the adjoint equation, we introduce additional variables, $\unicode[STIX]{x1D740}$ , $\unicode[STIX]{x1D701}$ and $\unicode[STIX]{x1D743}$ , called Lagrange multipliers, whose role is to enforce the constraint (4.4). We then transform the constrained optimisation problem to an unconstrained one, defining the Lagrangian

The derivation of the Lagrangian with respect to each variable is detailed in appendix D and leads to the following optimality system.

Direct system:

Adjoint system:

The terms $RHS_{1}$ , $RHS_{2}$ and $RHS_{3}$ are expressed in appendix D and depend on the definition of the cost functional (4.7), especially on the observation operator.

Optimality condition:

The minimisation problem can be solved iteratively (Bewley, Moin & Temam Reference Bewley, Moin and Temam2001) in a computationally efficient manner using this optimality system. Once the direct system (4.11) is solved, $(\boldsymbol{q},\unicode[STIX]{x1D6FC})$ are known and one can integrate the adjoint system (4.12) backwards in space to obtain $(\unicode[STIX]{x1D740},\unicode[STIX]{x1D701},\unicode[STIX]{x1D743})$ , and finally obtain the gradients $(\unicode[STIX]{x2202}J/\unicode[STIX]{x2202}\boldsymbol{f},\unicode[STIX]{x2202}J/\unicode[STIX]{x2202}\boldsymbol{\unicode[STIX]{x1D702}})$ , using the optimality condition (4.13). Using a quasi-Newton algorithm, $(\boldsymbol{f},\boldsymbol{\unicode[STIX]{x1D702}})$ can be updated and a new iteration can begin. The procedure stops when convergence is reached, giving access to the optimal forcing $\boldsymbol{f}$ and the optimal initial condition perturbation $\boldsymbol{\unicode[STIX]{x1D702}}$ . As a verification of the procedure described above, the gradient in (4.13) was also computed using finite differences; several tests demonstrated the consistency of the method.

4.2.1 Standard sensitivity

In the sensitivity study, we explore the effect of a small variation of the control parameters capable of improving the cost functional. For the sensitivity, we will mainly comment on the sensitivity with respect to the forcing, since it is related to the effect of nonlinearity. The sensitivity analysis in this section highlights regions where the wavepacket is most sensitive to nonlinear effects. Sensitivity with respect to changes in the initial condition will be exploited in § 4.2.3.

The gradient of the cost functional, ${\mathcal{J}}$ , with respect to the wavepacket forcing, $\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}$ , can be ambiguous to define and to interpret because of the scale separation assumption (4.1). We choose therefore to define the forcing that will be used at the second iteration step of the optimisation problem. Considering a Newton-like gradient algorithm, the $(n+1)$ th forcing term is found using the following rule:

where the relaxation parameter $\unicode[STIX]{x1D6FE}^{(n)}$ is found using a line-search algorithm (Nocedal & Wright Reference Nocedal and Wright1999). Because we consider the first iteration, $\boldsymbol{f}^{(1)}=0$ , and a step sufficiently small such that $\unicode[STIX]{x1D6FC}$ varies sufficiently slightly, we have the forcing in physical space expressed as

We then define

which represents the direction of the small forcing of the linear PSE, expressed in physical space, and which will improve the cost functional. This infinitesimal forcing, defined up to a constant, will be referred to henceforth as the forcing sensitivity of the cost functional, ${\mathcal{J}}$ . An infinitesimal response, $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}$ , can be defined in a similar manner. Let $\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{h}$ be the solution of the homogeneous PSE and $\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{f}$ the solution of the PSE forced by $\unicode[STIX]{x1D6FE}\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}$ with $\unicode[STIX]{x1D6FE}=10^{-8}$ ; we define $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}=\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{f}-\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{h}$ .

Sensitivity of PSE with respect to an external forcing is here explored for the $M=0.4$ round jet presented in Cavalieri et al. (Reference Cavalieri, Rodriguez, Jordan, Colonius and Gervais2013), the same as used for the resolvent analysis in § 3.2. The axisymmetric wavepacket ( $m=0$ ) at $St=0.6$ is again used as an example to explore the salient traits of the problem. Complementary results for other Strouhal numbers are shown in appendix E. In this standard sensitivity, we choose $\unicode[STIX]{x1D652}_{T}=0$ to be as parameter-free as possible.

Figure 10.

Standard sensitivity study. Superposition of the infinitesimal forcing $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}$ (isocontours) with the associated response $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}$ (colourmap). The critical-layer position (thick solid line) and the inflection point position (thick dashed line) are also reported on these graphs. Panel (e) respects the physical aspect ratio. (a) Real part of axial velocity. Positive contours for forcing $\text{Re}(\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u})$ with solid lines, negative with dashed lines. (b) Real part of radial velocity. Positive contours for forcing $\text{Re}(\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u})$ with solid lines, negative with dashed lines. (c) PSD of axial velocity. Solid lines are isocontours of forcing $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u}|^{2}$ . Dotted lines are isocontours of the response $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{u}|^{2}$ . (d) PSD of radial velocity. Solid lines are isocontours of forcing $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u}|^{2}$ . Dotted lines are isocontours of the response $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{u}|^{2}$ . (e) Vector field of $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u}$ , coloured by the vorticity of $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{u}$ . Positive contours for forcing $\text{Re}(\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u})$ with solid lines, negative with dashed lines.

Results for the forcing sensitivity $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}$ and associated response $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}$ for $\mathbb{H}(\boldsymbol{q},\unicode[STIX]{x1D6FC})=\boldsymbol{q}\text{e}^{\text{i}\int _{0}^{x}\unicode[STIX]{x1D6FC}(\unicode[STIX]{x1D709})\,\text{d}\unicode[STIX]{x1D709}}$ and ${\mathcal{Y}}=0$ are presented in figure 10. Since a local wavenumber of the perturbation is defined by PSE, we can draw the critical-layer position, as we did in the resolvent analysis, on $x$ – $r$ maps of forcing sensitivity and response. It can be seen in figure 10(c,d) that the highest sensitivity follows the critical layer and has a maximum value at approximately $x/D=3$ ; we note that the locally parallel analysis also identified a maximal sensitivity at this axial station. This axial position is where the wavepacket growth rate $-\text{Im}(\unicode[STIX]{x1D6FC})$ becomes small (see figure 13 b), such that the phase speed is nearly real and a critical layer is present. Similar computations have been performed for different Strouhal numbers (see appendix E), and the axial position of the maximum value of the forcing sensitivity is reported in figure 11, compared with the axial position at which the wavepacket is neutrally stable, i.e. where $-\text{Im}(\unicode[STIX]{x1D6FC})=0$ . The maximum sensitivity appears systematically just upstream of the neutral position and follows the same trend. In the locally parallel analysis, the maximum sensitivity is expected for a neutral wave. Here, the mean velocity profile changes with axial position due to the jet divergence. For lower $x/D$ , the profile is steeper and the critical and inflection points are closer; proximity of these two points also leads to higher sensitivity to forcing, as studied in Tissot et al. (Reference Tissot, Zhang, Lajús, Cavalieri, Jordan and Colonius2015). This may explain why the maximum sensitivity appears slightly upstream of the neutral position.

Figure 11.

Axial position $x_{\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{max}}$ of the maximum forcing sensitivity $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}|^{2}$ (standard and 4D-Var) and axial position $x_{\unicode[STIX]{x1D6FC}=0}$ where the wavepacket is neutrally stable, i.e. where $-\text{Im}(\unicode[STIX]{x1D6FC})=0$ , function of the Strouhal number.

From $x/D=1$ to $x/D=5$ , it can be seen in figure 10(a) that the infinitesimal forcing and response are in phase; both follow the critical layer and are inclined against the shear in a manner similar to that of the OF modes of the parallel flow (figure 9). Then, from $x/D=5$ to $x/D=9$ , the forcing becomes small and the response is convected downstream. Figure 10(a,e) shows how advecting vortices are centred around the critical layer; an observation consistent with the minimum of axial velocity magnitude that follows the critical layer in figure 10(c). These inclined perturbations appear to be tilted by the shear around the critical layer as well. Moreover, figure 10(d) shows, around $r/D=0.5$ , a growth and decay from $x/D\sim 6$ to $x/D\sim 9$ of the PSD of the radial component of the response, while the axial component remains strong and follows the critical layer. This behaviour suggests the presence of an Orr mechanism.

Figure 12.

Radial profiles of the PSD. Forcing sensitivity $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}|^{2}$ (green bullet), the first optimal forcing mode $\unicode[STIX]{x1D731}_{1}$ (black triangle) and the associated response $\unicode[STIX]{x1D70E}_{1}\unicode[STIX]{x1D733}_{1}$ (red circle). The critical-layer position (thick solid line) and the inflection point position (thick dashed line) are also reported on these graphs.

A comparison between forcing sensitivity and response with the locally parallel resolvent analysis of § 3.2.2 can be performed. Figure 12 compares radial profiles of sensitivity $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}$ with the optimal forcing mode $\unicode[STIX]{x1D731}_{1}$ and associated response $\unicode[STIX]{x1D70E}_{1}\unicode[STIX]{x1D733}_{1}$ at $x/D=3$ and $x/D=8$ . The first location is close to the region of peak sensitivity to forcing. The second is more characteristic of structures downstream of the potential core. A similar trend of the forcing between resolvent and PSE is shown, with maxima located near the critical layer, especially for $x/D=3$ . From these results we can conclude that the sensitivity analysis based on locally parallel base flows leads to general trends that are representative of the behaviour for a slowly diverging jet. At $x/D=8$ , the forcing is three orders of magnitude weaker. The vertical component of the forcing closely resembles the optimal forcing mode, while the axial component does not, but both peak close to the critical layer. Differences from locally parallel resolvent analysis may be attributed to cumulative effects in the streamwise direction, not accounted for in the resolvent analysis but included in PSE, which considers the entire wavepacket evolution. Disturbances thus have memory of forcing applied upstream.

The sensitivity of wavepackets predicted by PSE corroborates the results of resolvent analysis which have shown high sensitivity around the critical layer with features suggestive of an Orr mechanism. This non-parallel computation shows how the response is active downstream of the forcing locations. Small differences in comparison with the locally parallel analysis have also been noticed due to flow divergence and memory of the upstream perturbations. In order to further probe the effect of nonlinearity, we will drive this sensitivity study using the experimental data, and then push the study to convergence.

4.2.2 Sensitivity of 4D-Var

In order to explain the mismatch between linear models and experimental measurements, the sensitivity is central, but it has to be weighted by the nonlinear term actually present. With 4D-Var, by solving an inverse problem, we use the observation error to get the minimal nonlinearity sufficient to explain the mismatch.

The observations ${\mathcal{Y}}$ considered are the PSD of the axial and radial components of the TRPIV measurements radially filtered by POD (Cavalieri et al. Reference Cavalieri, Rodriguez, Jordan, Colonius and Gervais2013). The observation of the axial component on the centreline, $\widetilde{u}_{r=0}$ , is compared in figure 13(a) with the values $\mathbb{H}(\boldsymbol{q},\unicode[STIX]{x1D6FC})$ predicted by PSE. It shows the good agreement of the growth of the instability wave, until this wave becomes neutrally stable at approximately $x/D\sim 5$ and then begins to decay, while experiments show an increase after the end of the potential core. Figure 13(b) shows the real and imaginary parts of the local wavenumber, $\unicode[STIX]{x1D6FC}$ , predicted by PSE.

Figure 13.

Characterisation of the observations and of the PSE results. (a) Centerline axial fluctuating kinetic energy. Comparision between the observations ${\mathcal{Y}}$ (line with circles) with the PSE (dashed line) and the associated observation operator $\mathbb{H}(\boldsymbol{q},\unicode[STIX]{x1D6FC})$ (squares). (b) Local wavenumber predicted by PSE. Real part and imaginary part.

Figure 14(a,b) shows $x$ – $r$ maps of the observed (experimental) PSD of the axial and radial components, while figure 14(c,d) shows the axial and radial components of the homogeneous PSE solution. These plots show again how linear PSE and experiments match until approximately $x/D\sim 4$ . Then, looking at the PSD of the axial velocities of the observations (figure 14 a) and PSE (figure 14 c), we note that PSE has higher values around ( $x/D=5,r/D=0.25$ ) than the experiments. The experiments, on the other hand, have a significant increase of axial fluctuating kinetic energy around the centreline at $x/D=8$ , while PSE decays. These differences are clearly visible in the error map figure 15(a). Moreover, this error map shows that, downstream of the end of the potential core, an error follows the critical layer until the end of the domain.

Figure 14.

Comparison between observations (first POD mode of TRPIV measurements) and PSE solution. (a) ${\mathcal{Y}}$ axial component. (b) ${\mathcal{Y}}$ radial component. (c) PSE solution, axial component of $|\tilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}|^{2}$ . (d) PSE solution, radial component of $|\tilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}|^{2}$ .

Figure 15.

Superposition of the observation error $(\mathbb{H}(\boldsymbol{q},\unicode[STIX]{x1D6FC})-{\mathcal{Y}})^{2}$ (grey level and thin dotted isocontours) and the forcing sensitivity $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}|^{2}$ (thin solid isocontours). The dashed lines correspond also with the forcing sensitivity $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}|^{2}$ , but the isocontours are 10 times closer to each other in order see the smaller ‘spot’ ${\mathcal{S}}_{2}$ . The axial component is considered. The critical-layer position (thick solid line) and the inflection point position (thick dashed line) are also reported on these graphs. In figure 15(b), the final axial position $L$ is chosen such that the errors near $x/D=8$ around the centreline are not observed. As a consequence, the near-centreline forcing ${\mathcal{S}}_{2}$ disappears.

In order to explore these errors, the sensitivity has been computed using the new observation space, $\unicode[STIX]{x1D652}_{o}=\mathbb{I}$ and $\unicode[STIX]{x1D652}_{T}=\unicode[STIX]{x1D70E}_{T}\mathbb{I}$ , with $\unicode[STIX]{x1D70E}_{T}=1$ . The terms $\unicode[STIX]{x1D652}_{f}$ and $\unicode[STIX]{x1D652}_{\unicode[STIX]{x1D702}}$ are unimportant for the sensitivity since $\boldsymbol{f}=0$ and $\boldsymbol{\unicode[STIX]{x1D702}}=0$ . Two spots of high sensitivity to forcing can be observed. Both are upstream of the locations of disagreement between linear PSE and experiments. The first, denoted ${\mathcal{S}}_{1}$ , peaks at approximately ( $x/D=3,r/D=0.45$ ), and is similar to the standard sensitivity presented in § 4.2.1. The second, denoted ${\mathcal{S}}_{2}$ , is located near the centreline at approximately $x/D=5$ , where the wavepacket amplitude is highest. To further explore what underlies these sensitivities, another computation is performed (results shown in figure 15 b), where the domain is truncated to $x_{\max }/D=7$ , such that the increase of experimental values near the centreline at $x/D=8$ is not observed; ${\mathcal{S}}_{2}$ disappears, and ${\mathcal{S}}_{1}$ is slightly deformed but remains qualitatively the same, showing that the small additional centreline forcing identified by 4D-Var appears to arise in order to reduce the downstream error. This suggests that even if the wavepacket is not particularly sensitive in that region, high centreline turbulent fluctuating energy near the end of the potential core contributes to reducing the downstream error. This location has been identified in the literature (Jordan & Gervais Reference Jordan and Gervais2008, and references therein) to comprise fluctuations with significant correlation with the radiated sound, but not as the source itself, as this is more consistent with a spatially extended wavepacket (Jordan & Colonius Reference Jordan and Colonius2013).

The high sensitivity at ${\mathcal{S}}_{1}$ is also predicted by the standard PSE sensitivity of § 4.2.1 and by the locally parallel resolvent analysis of § 3.2.2. It is interesting to see that 4D-Var also identifies this region, now considering observations, as sensitive to forcing, such that inclusion of nonlinearities would lead to improvements of model results. Figure 16 is the counterpart of figure 10 for the 4D-Var sensitivity. The similarity between the results indicates that the sensitive regions are actually excited by the background turbulence. Moreover, the position of maximum sensitivity versus $St$ is reported in figure 11, and the radial profiles are plotted in figure 17, showing again the same trend.

Figure 16.

Superposition of the infinitesimal forcing $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}$ (isocontours) with the associated response $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}$ (colourmap). The critical-layer position (thick solid line) and the inflection point position (thick dashed line) are also reported on these graphs. Panel (e) respects the physical aspect ratio. (a) Real part of axial velocity. Positive contours for forcing $\text{Re}(\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u})$ with solid lines, negative with dashed lines. (b) Real part of radial velocity. Positive contours for forcing $\text{Re}(\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u})$ with solid lines, negative with dashed lines. (c) PSD of axial velocity. Solid lines are isocontours of forcing $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u}|^{2}$ . Dotted lines are isocontours of the response $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{u}|^{2}$ . (d) PSD of radial velocity. Solid lines are isocontours of forcing $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{v}|^{2}$ . Dotted lines are isocontours of the response $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{v}|^{2}$ . (e) Vector field of $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}$ , colored by the vorticity of $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}$ . Positive contours for forcing $\text{Re}(\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u})$ with solid lines, negative with dashed lines.

Figure 17.

Radial profiles of the PSD of the axial component: observations ${\mathcal{Y}}$ (blue squares), PSE (black dashed line), forcing sensitivity $|\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}|^{2}$ (green bullet), the first optimal forcing mode $\unicode[STIX]{x1D731}_{1}$ (black triangle) and the associated response $\unicode[STIX]{x1D70E}_{1}\unicode[STIX]{x1D733}_{1}$ (red circle). The critical-layer position (thick solid line) and the inflection point position (thick dashed line) are also reported on these graphs.

Figure 18.

Comparison between observations (first POD mode of TRPIV measurements), homogeneous PSE solution and PSE solution forced by 4D-Var for $St=0.6$ . (a) Homogeneous PSE solution, axial component of $|\tilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}|^{2}$ . (b) Homogeneous PSE solution, radial component of $|\tilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}|^{2}$ . (c) Forced PSE solution, axial component of $|\tilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}|^{2}$ . (d) Forced PSE solution, radial component of $|\tilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}|^{2}$ . (e) ${\mathcal{Y}}$ axial component. (f) ${\mathcal{Y}}$ radial component.

The response (figure 16 c) and the error (figure 15 a) can be directly compared; their correspondence shows that the response actually decreases the error. The observed mechanism can be described as follows: nonlinearities induce disturbances inclined against the shear in the critical-layer region; advection of such disturbances leads to a tilting of structures with an amplitude growth characteristic of the Orr mechanism; this growth increases amplitudes in the downstream region, and wavepacket amplitudes get closer to experimental results.

All of the results presented in this section are for $St=0.6$ , but other Strouhal numbers display the same trends. Some results for other Strouhal numbers are shown in appendix E, and are consistent with the phenomena described here, with some minor differences.

To summarise, the PSE-4D-Var reveals that a high sensitivity to forcing at approximately $x/D=3$ promotes the generation of convected perturbations that decrease the error between PSE and experiment at the end of the potential core. This mechanism is consistent with the locally parallel results and is believed to be associated with the critical layer. Downstream of the potential core, results suggest that an Orr mechanism may be active around the critical layer, tilting the structures and generating a growth and decay of the PSD of the radial velocity. The results suggest that nonlinearities may act on wavepackets via this mechanism, with significant effects in the downstream region. Another important region was found on the jet centreline, close to the end of the potential core. This region did not show high sensitivity to forcing in the previous section. Since PSE-4D-Var determines forcings that minimise the error between model and reference measurements, the results suggest that high-amplitude forcings (i.e. significant nonlinear effects) are present in the identified region.

4.2.3 Convergence to experimental results

The solution of the minimisation of (4.7) is the minimal nonlinear term necessary to best approach the measurements. We now solve this problem by iterating the direct/adjoint procedure until convergence. A steepest descent algorithm is used with a backtracking Armijo-type linsearch method. The covariance matrices are chosen as $\unicode[STIX]{x1D652}_{o}=\mathbb{I}$ , $\unicode[STIX]{x1D652}_{\unicode[STIX]{x1D702}}=\unicode[STIX]{x1D70E}_{\unicode[STIX]{x1D702}}\mathbb{I}$ , $\unicode[STIX]{x1D652}_{f}=\unicode[STIX]{x1D70E}_{f}\mathbb{I}$ and $\unicode[STIX]{x1D652}_{T}=\unicode[STIX]{x1D70E}_{T}\mathbb{I}$ , with $\unicode[STIX]{x1D70E}_{f}=\unicode[STIX]{x1D70E}_{\unicode[STIX]{x1D702}}=1$ and $\unicode[STIX]{x1D70E}_{T}=10^{-2}$ . These values have been determined by the L-curve method (Hansen Reference Hansen1992; Cordier, El Majd & Favier Reference Cordier, El Majd and Favier2010).

Figure 18 compares homogeneous PSE, forced PSE and measurements for $St=0.6$ , showing that the forcing is successful in producing waves with PSD that match those of the experiments.

The forcing $\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}$ and the response $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}=\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{f}-\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{h}$ are plotted in figure 19. The forcing magnitude is approximately $10\,\%$ of the wavepacket amplitude, which is significant but remains small. From a reduced-order-modelling perspective, the relatively low amplitude of the required forcing suggests that linear wavepacket models might be substantially improved, in terms of agreement with experiment, by slight modifications via the inclusion of some relevant nonlinear interactions.

The forcing shape is more complex than what was found in the sensitivity analyses, but it possesses the same features. A high-amplitude forcing follows the critical layer, and another important forcing region is concentrated around the centreline. There is, moreover, a third region of relevant forcing that appears to follow the inflection point but vanishes quickly. The two major differences compared with the 4D-Var sensitivity of § 4.2.2 are that the forcing continues along the critical layer much further downstream, and the centreline forcing is as large as the critical-layer forcing, while it is one order of magnitude less for the 4D-Var sensitivity. The downstream continuation of the forcing can be explained by the fact that as long as the forcing excites the wavepacket, nearly neutral growth of amplitudes can be sustained. A critical layer is then still present and the wavepacket continues to be sensitive to turbulence. This scenario remains plausible since high turbulent levels are still present at these downstream stations. Since wavy disturbances are continuously excited while propagating downstream, the tilting of structures by an Orr-like mechanism is a little less clear in figure 19 than in the sensitivity study, but the overall trend can still be perceived. Concerning the strong centreline forcing, the present result confirms the analysis that even though the wavepacket is not particularly sensitive near the centreline, it has to be forced there in order to lead to the observed experimental amplitudes.

Figure 19.

Superposition of the converged 4D-Var forcing $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}$ (isocontours) with the associated response $\unicode[STIX]{x1D6FF}\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}=\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{f}-\widetilde{\boldsymbol{q}}_{\unicode[STIX]{x1D714},m}^{h}$ (colourmap); real part of the axial velocity; positive contours for forcing $\text{Re}(\widetilde{\boldsymbol{f}}_{\unicode[STIX]{x1D714},m}^{u})$ with solid lines, negative with dashed lines. The critical-layer position (thick solid line) and the inflection point position (thick dashed line) are also reported. The physical aspect ratio is respected.

Finally, the observation errors seem to be insensitive to the inflow initial condition (sensitivity not plotted), which is probably a consequence of the dominance of the K–H mode in the upstream region and a fast decay of all of the other modes.

We have learnt by converging 4D-Var that the minimal forcing excites the wavepacket along the critical layer throughout its evolution, which, in the downstream region, maintains nearly constant disturbance amplitudes as these are advected. Moreover, detection of a relevant forcing near the centreline, despite the relatively low sensitivity of this region to external forcing, suggests the presence of strong active turbulence in that region, capable of considerably influencing wavepacket evolution.