3.1 Test of base flow structure computations

The flow past a 180-degree sharp bend is a deceptively difficult problem to fully resolve, especially at large Reynolds number, due to the sensitivity of the solution to the mesh structure mainly near the sharp bend. This section describes the tests used to validate the numerical algorithm and to select appropriate meshes and element order. For the spatial resolution study, we varied element polynomial degree from $N=4$ to $N=8$ for a mesh based on domain length parameters from the mesh domain. For consistency with the domain size study, the mesh employed in this study models a 180-degree sharp bend with opening ratio $\unicode[STIX]{x1D6FD}=1$ and $\mathit{Re}=500$ . In this regime, the flow is steady, with two recirculation bubbles. Figure 2 shows the detail of the mesh with $N=3$ . The mesh is structured and refined in the vicinity of the sharp bend as well as in the downstream channel in order to capture the detailed structure of the flow that passes around the bend.

Figure 2. Details of the mesh around the turning area with polynomial order $N=3$ .

Figure 3. Sketch of separation and reattachment points defining the locations of all recirculations. The unbroken arrow represents the direction of the main bulk flow as it navigates the bend.

To demonstrate the accuracy of computing recirculation length in the base flow, table 1 shows the relative error of several measured quantities as a function of polynomial order. From computations of error on the recirculation length (table 1), it was found that the polynomial order $N=5$ provides good accuracy for running the base flow computations. To examine the effect of downstream channel length on the solutions, a convergence study was conducted on the lengths of the first recirculation bubble on the bottom wall of the downstream channel and the secondary recirculation bubble on the top wall. The results shown in table 2 demonstrate that an outlet length of 10 results in an error in the determination of the primary recirculation bubble length of approximately 0.007 %, and the error is approximately 4.5 % for the secondary bubble. The relatively large error in the size of the secondary bubble is caused by its proximity to the outlet. Outlet lengths of 20–100 are required to achieve at least 5 significant figures of accuracy. Hence, an outlet length of 30 is considered adequately long to be used throughout this study. Barton (Reference Barton1997) and Cruchaga (Reference Cruchaga1998) studied the entrance effect for backward-facing step flow with an expansion ratio of 2 and found that inlet lengths of $10h$ and $2h$ (where $h$ is the step height) return slightly different numerical solutions compared to that with zero inlet length. In this study, the inlet length is 15, which is adequately long for the velocity flow to be fully developed before reaching the bend.

3.2 Test of eigenvalue computations

The precision of the eigenvalue $\unicode[STIX]{x1D707}$ and eigenmode $\hat{\boldsymbol{u}}$ produced by subspace iteration was quantified by the residual

where $\Vert \cdot \Vert$ is the standard vector norm and eigenmodes were normalised ( $\Vert \hat{\boldsymbol{u}}\Vert =1$ ). The linear stability analysis technique relied on an iterative process to obtain the leading eigenvalues and eigenmodes of the system. The process ceased when $r<10^{-7}$ was achieved. Nevertheless, the eigenmodes are also resolution dependent. Table 3 reports the accuracy of the eigenvalue computations as a function of element polynomial degree  $N$ . The leading eigenvalue for $\mathit{Re}=500$ , $\unicode[STIX]{x1D6FD}=1$ and $k=6.4$ is real and linearly stable. It is found that the eigenvalue converged to an error of merely 0.0003180 % at $N=5$ , which is employed hereafter.

Figure 4. Length of the first recirculation bubble ( $L_{R_{1}}$ ) against Reynolds number ( $\mathit{Re}$ ) for $\unicode[STIX]{x1D6FD}=1$ , comparing the present results to those of Chung, Tucker & Roychowdhury (Reference Chung, Tucker and Roychowdhury2003) and Zhang & Pothérat (Reference Zhang and Pothérat2013).

3.3 Code validation

Finally, the present model was validated against published two-dimensional simulations providing the position and length of the recirculation bubble. In a viscous flow (featuring boundary layers adjacent to no-slip surfaces), the point of reattachment can be precisely measured by finding the location where the wall shear stress $\unicode[STIX]{x1D70F}_{wall}=-\unicode[STIX]{x1D707}\unicode[STIX]{x2202}u/\unicode[STIX]{x2202}y$ is zero.

Figure 4 shows a comparison between the recirculation length of the first bubble ( $L_{R_{1}}$ ) in a flow with $\unicode[STIX]{x1D6FD}=1$ as a function of Reynolds number between the present and previously reported results digitised from figures in Chung et al. (Reference Chung, Tucker and Roychowdhury2003) and Zhang & Pothérat (Reference Zhang and Pothérat2013). The coefficient of determination, $R^{2}$ , between the present data and those of these previous studies differ by just 2.2 % and 0.2 %, respectively. The comparison displays a strong agreement between the studies, with each curve increasing rapidly and linearly at Reynolds numbers $\mathit{Re}\lesssim 200$ , before transitioning to a more gradual linear regime of further bubble elongation beyond $\mathit{Re}\approx 300$ . This regime terminates with the onset of unsteady flow. The sudden drop in the data from Zhang & Pothérat (Reference Zhang and Pothérat2013) at $\mathit{Re}\approx 600$ coincides with the onset of unsteady flow in that study. The present computations return a steady state (in agreement with Chung et al. Reference Chung, Tucker and Roychowdhury2003) up to $\mathit{Re}\approx 700$ .

4.1 Flow regimes

In this section, we focus on the behaviour of the two-dimensional base flow, especially around the sharp bend and along the downstream channel. For the range of opening ratios $\unicode[STIX]{x1D6FD}$ studied, four regimes are identified (figure 5). The first regime exhibits only a single recirculation bubble immediately behind the sharp bend. The second regime sees the emergence of a bubble at the opposite wall slightly downstream of the first bubble. The third regime reveals the appearance of a small counter-rotating recirculation bubble between the primary recirculation bubble and the bottom wall. Finally, the fourth regime marks the development of an unsteady two-dimensional flow. The Reynolds numbers at the onset of each regime are denoted by $\mathit{Re}_{R_{1}}$ , $\mathit{Re}_{R_{2}}$ , $\mathit{Re}_{in}$ and $\mathit{Re}_{c}$ , respectively, as shown in figure 6. The results of the current study agree with those of Zhang & Pothérat (Reference Zhang and Pothérat2013), but small discrepancies are found for $\mathit{Re}_{in}$ at $\unicode[STIX]{x1D6FD}\geqslant 1$ , which are attributed to the different node density in the meshes at high $\mathit{Re}$ in both studies. The present study also extends the lower end of the range of $\unicode[STIX]{x1D6FD}$ from $\unicode[STIX]{x1D6FD}=0.1$ (Zhang & Pothérat Reference Zhang and Pothérat2013) to $\unicode[STIX]{x1D6FD}=0.0125$ . The constriction at the bend at smaller $\unicode[STIX]{x1D6FD}$ leads to high velocities and shear in that region.

Figure 5. Streamlines of steady two-dimensional base flow at (a)  $\mathit{Re}=10$ , (b)  $\mathit{Re}=200$ and (c)  $\mathit{Re}=600$ ; (d) flooded contours of vorticity magnitude at $\mathit{Re}=800$ for $\unicode[STIX]{x1D6FD}=0.5$ . In (d), white represents zero vorticity (no rotation) and darker shading denotes arbitrarily higher vorticity magnitude levels.

Figure 6. Reynolds number where the primary recirculation bubble ( $\mathit{Re}_{R_{1}}$ ), secondary recirculation bubble ( $\mathit{Re}_{R_{2}}$ ) and inside recirculation bubble ( $\mathit{Re}_{in}$ ) start to appear in the two-dimensional flow; $\mathit{Re}_{c}$ is the transition from steady to unsteady. The lines are from the current study and the symbols are results from Zhang & Pothérat (Reference Zhang and Pothérat2013): ○, ▵, ▫ and ▿ represent $\mathit{Re}_{R_{1}}$ , $\mathit{Re}_{R_{2}}$ , $\mathit{Re}_{in}$ and $\mathit{Re}_{c}$ , respectively.

For clarity, we shall highlight the main features of the regimes shown in figure 6, but a more detailed description can be found in Zhang & Pothérat (Reference Zhang and Pothérat2013). The onset of all regimes is delayed consistently to larger $\mathit{Re}$ at larger $\unicode[STIX]{x1D6FD}$ over $\unicode[STIX]{x1D6FD}\lesssim 1$ , and it remains almost constant when $\unicode[STIX]{x1D6FD}\gtrsim 1$ . It is worth noting that the behaviour of the flow is different between $\unicode[STIX]{x1D6FD}<0.2$ , $0.2<\unicode[STIX]{x1D6FD}<1$ and $\unicode[STIX]{x1D6FD}\gtrsim 1$ . For $\unicode[STIX]{x1D6FD}\lesssim 0.2$ , a flow resembling jet flow is created through the narrow bend orifice, causing the flow to accelerate transversely and hit the top wall before deflecting towards the streamwise direction. This causes the onset of each regime to occur at low $\mathit{Re}$ . By contrast, at $0.2\lesssim \unicode[STIX]{x1D6FD}\lesssim 1$ , the width of the opening of the bend is sufficient for the flow to turn smoothly to the streamwise direction around the bend; hence, $\unicode[STIX]{x1D6FD}$ influences the onset of each regime. However, at $\unicode[STIX]{x1D6FD}\gtrsim 1$ , the onsets do not vary significantly because a recirculation bubble develops at the outer wall of the bend, which confines the turning flow, such that the true breadth of the bend opening is no longer apparent. Thus the flow behaves similarly to that of $\unicode[STIX]{x1D6FD}=1$ and the onset $\mathit{Re}$ values for the flow regimes remain almost constant for $\unicode[STIX]{x1D6FD}\gtrsim 1$ .

At very low $\mathit{Re}$ , the flow in the downstream and upstream channel is almost symmetrical with respect to $y=0$ . However, as $\mathit{Re}$ increases, the flow streamlines at the bottom wall move upward until an inflexion point becomes visible behind the edge of the sharp bend. Consequently, flow separation occurs and a recirculation bubble is formed. In the range of $\unicode[STIX]{x1D6FD}$ studied, the primary recirculation bubble appears at very small $\mathit{Re}$ because of a strong adverse pressure gradient behind the bend, as can be seen in figure 5(a). Theoretically a sharp edge always causes separation to a flow, even at very low $\mathit{Re}\rightarrow 0$ (Taneda Reference Taneda1979; Moffatt Reference Moffatt1985). The finite $\mathit{Re}_{R_{1}}$ captured in this study is a numerical artefact of finite spatial resolution at the sharp bend corner. The finite discretisation obscures the bubble at very low $\mathit{Re}$ .

As $\mathit{Re}$ increases, the secondary recirculation bubble appears when the flow streamlines in the bulk flow above the primary recirculation bubble move away from the top wall. This results in a strong adverse pressure gradient at the wall, which causes another separation to occur at $\mathit{Re}_{R_{2}}$ (figure 5 b). Under the same circumstances, an inner counter-rotating recirculation bubble, as seen in figure 5(c), is formed between the primary recirculation bubble and the bottom wall when the backflow in the primary recirculation bubble moves away from the wall.

Zhang & Pothérat (Reference Zhang and Pothérat2013) found an additional regime between $\mathit{Re}_{in}$ and $\mathit{Re}_{c}$ , where they reported a small-scale vortices structure far downstream of the channel; this state is not observed in the current study. The $\mathit{Re}_{c}$ in figure 6 indicates an unsteady flow where large eddies or vortices are shed downstream from the sharp corner of the bend as depicted in figure 5(d); $\mathit{Re}_{c}$ is found to be in agreement with the values in the Zhang & Pothérat (Reference Zhang and Pothérat2013) study, and it is important to note that both studies started the simulations from rest to obtain this $\mathit{Re}_{c}$ . In a further analysis, hysteretic behaviour has been observed when the simulation is started from different initial conditions. Considering the unsteady flow as a departure from the base steady flow at the same Reynolds number, its amplitude $|A|$ is measured as the ${\mathcal{L}}^{2}$ norm (the integral of the magnitude of velocity over the computational domain) of the difference between velocity fields in these states. Figure 7 shows the variations of $|A|$ , and therefore regimes of unsteady flow where $|A|>0$ , when the Reynolds number is incrementally varied for $\unicode[STIX]{x1D6FD}=1$ . Two distinct onsets of two-dimensional unsteadiness are found by initiating the simulation from three different initial conditions: (i) flow at rest, (ii) a snapshot of unsteady flow computed at a slightly lower Reynolds number and (iii) the steady flow solution obtained at a slightly lower Reynolds number. It is evident from the figure that the simulations starting from the flow at rest and from an unsteady flow yield the same value of $\mathit{Re}_{c}=742$ , whereas the simulations starting from a steady flow become unsteady at a higher $\mathit{Re}$ , $\mathit{Re}_{c}\approx 1150$ .

Figure 7. Hysteretic behaviour described by the fluctuation of the integral of velocity magnitude throughout the domain as a function of $\mathit{Re}$ at $\unicode[STIX]{x1D6FD}=1$ : reducing $\mathit{Re}$ from an unsteady flow (◂) and increasing $\mathit{Re}$ from a steady flow (▸) give different $\mathit{Re}_{c}$ ; $|A|=0$ indicates the steady flow solution.

The location where shedding initiates at the onset depends on the initial conditions too. When a simulation is initiated from a flow at rest, the vortex shedding can be seen to emerge from the sharp corner of the bend, as illustrated in figure 8(a). It is likely that large-amplitude perturbations caused by the impulsive initiation of the flow are sufficient to induce a shedding from the bend that bypasses the downstream destabilisation, and beyond $\mathit{Re}_{c}\approx 742$ (at least in these simulations) this is self-sustaining. A similar observation was made when reducing $\mathit{Re}$ from an unsteady flow: a large sudden decrease could revert the unsteady flow to steady state at a Reynolds number at which the unsteady state could be preserved via a gradual decrease in Reynolds number.

Figure 8. Flooded contours of vorticity magnitude demonstrating unsteady saturated flows for (a)  $\mathit{Re}=800$ at $\unicode[STIX]{x1D6FD}=1$ (which was initiated from rest) and (b–d)  $\mathit{Re}=1152$ at $\unicode[STIX]{x1D6FD}=1$ (which was initiated from a saturated steady-state flow solution at a lower Reynolds number). Contour levels are as per figure 5(d). Panels (b–d) show domain segments $-45\leqslant x\leqslant -29$ , $-30\leqslant x\leqslant -14$ and $-15\leqslant x\leqslant 1$ , respectively.

When $\mathit{Re}$ is increased gradually along the steady flow branch, the flow remains steady up to $\mathit{Re}\simeq 1150$ . Figure 8(b–d) depicts the saturated flow at $\mathit{Re}=1152$ . This final state exhibits a wavy disturbance extending 25 diameters downstream of the bend. In this case, unsteadiness first manifested in the shear layers behind the secondary recirculation bubble in the form of small eddies, but over time the flow in this region more proximate to the bend restabilised and the unsteady region retreated to its ultimate position further downstream. This has some resemblance to the regime described by Zhang & Pothérat (Reference Zhang and Pothérat2013) before the flow becomes unsteady in their study. It is likely that they found this regime at lower $\mathit{Re}$ due to the high sensitivity to mesh resolution of this feature. Since the flow was already unstable far downstream of the bend, noise tended to be amplified and created small vortical structures. In testing this hypothesis, we found that by increasing the resolution of the mesh in our study, the onset of unsteadiness could be delayed significantly by increasing $\mathit{Re}$ gradually. It is also plausible that the length of the outlet channel would influence the upper limit of the steady-state branch, though this was not tested. Finally, it is noted that when comparing the flow pattern and the region of the flow producing the unsteady flow features between figures 8(a) and 8(b–d), the two depicted unsteady branches are different. It will be shown in § 5 that the two-dimensional flows are unstable to three-dimensional perturbations at Reynolds numbers below those of this hysteretic zone, so it will not be characterised further here.

4.2 Bubble separation points (steady flows)

Figure 9 shows the Reynolds number dependence of the separation and reattachment points (expressed by their distance downstream of the bend,  $x_{s}$ ) on both the bottom and top walls of the channel for $\unicode[STIX]{x1D6FD}=1$ . The empty circle symbol indicates the limit of the primary recirculation bubble behind the sharp corner. At $\mathit{Re}=240$ , a secondary recirculation bubble appears on the top wall at $x_{s}=3.77$ . At $\mathit{Re}\approx 700$ , two additional recirculation bubbles appear at the bottom wall; one is a small bubble in the primary recirculation bubble and the other is near the reattachment point of the secondary recirculation bubble.

From figure 9, we notice that there are regions where the locations of the separation and reattachment points are almost linear functions of $\mathit{Re}$ . It can clearly be seen in figure 9 that when the bifurcation regions around $\mathit{Re}=300$ and 750 are excluded, the locations of both points for all bubbles behave almost linearly with $\mathit{Re}$ . When a new recirculation bubble appears, the growth of the recirculation bubble just upstream of it is affected and so is the variation with Reynolds number of its reattachment point. The same behaviour was also found in the backward-facing step flow by Erturk (Reference Erturk2008).

Figure 9. A plot showing the Reynolds number dependence of the locations of separation points (measured by their distances from the bend – left of the origin), $x_{s}$ , for the base flows. Open circles represent the stagnation points for the primary recirculation bubbles. Solid circles denote the stagnation points for the secondary recirculation bubbles which form at $\mathit{Re}\approx 240$ and $x_{s}\approx 3.77$ . Solid and open squares represent the stagnation points for the third bubble and inner recirculation bubble, respectively, which are formed at higher  $\mathit{Re}$ . The top panel shows the streamlines of the base flow and the separation points of primary, secondary, inner and third recirculation bubbles at $\mathit{Re}=1100$ in the outlet channel.

The effect of $\unicode[STIX]{x1D6FD}$ on the size of the primary recirculation bubble is illustrated in figure 10. At small values of $\unicode[STIX]{x1D6FD}$ , the main bulk flow is accelerated by the small jet opening, making the primary recirculation bubble elongated in the streamwise direction. This causes the bubble to be bigger than that at larger  $\unicode[STIX]{x1D6FD}$ ; $L_{R_{1}}$ decreases as $\unicode[STIX]{x1D6FD}$ increases, but as $\unicode[STIX]{x1D6FD}>1$ and $\mathit{Re}\gtrsim 200$ , the size of the bubble increases due to the effect of the recirculation bubble at the far end of the bend wall.

Figure 10. A plot of the length of the primary recirculation bubble as a function of $\unicode[STIX]{x1D6FD}$ for Reynolds numbers as labelled.

Since for $\unicode[STIX]{x1D6FD}>1$ part of the flow in the bend is trapped in a closed eddy outside of the through flow taking the bend, Zhang & Pothérat (Reference Zhang and Pothérat2013) defined an effective opening ratio $\unicode[STIX]{x1D6FD}_{eff}$ based on the horizontal thickness of the flow effectively turning from inlet to outlet at $y=0$ ; that is, the distance from the inner vertical wall at the bend to the dividing streamline separating the turning flow from the closed recirculation. The same definition is used in the present study. For $\unicode[STIX]{x1D6FD}\leqslant 1$ , $\unicode[STIX]{x1D6FD}_{eff}=\unicode[STIX]{x1D6FD}$ (e.g. see figure 11 a), but for $\unicode[STIX]{x1D6FD}>1$ , starting from $\mathit{Re}\approx 30$ , $\unicode[STIX]{x1D6FD}_{eff}$ is smaller than $\unicode[STIX]{x1D6FD}$ (e.g. see figure 11 b). At $\mathit{Re}>200$ , $\unicode[STIX]{x1D6FD}_{eff}$ saturated close to 0.7. This degradation of $\unicode[STIX]{x1D6FD}_{eff}$ to values below $\unicode[STIX]{x1D6FD}$ , and its saturation behaviour at large $\mathit{Re}$ , can be observed in figure 11(c). For $\unicode[STIX]{x1D6FD}=2$ in figure 11(b), the effective $\unicode[STIX]{x1D6FD}$ is found to be $\unicode[STIX]{x1D6FD}_{eff}\approx 0.8$ , which is slightly higher than what was found by Zhang & Pothérat (Reference Zhang and Pothérat2013).

Figure 11. Streamlines of steady two-dimensional flow for $\mathit{Re}=600$ : (a)  $\unicode[STIX]{x1D6FD}=1$ and (b)  $\unicode[STIX]{x1D6FD}=2$ . A big recirculation bubble appears at the far end of the wall when $\unicode[STIX]{x1D6FD}>1$ , causing the effective bend opening ratio to be less than the actual opening; $\unicode[STIX]{x1D6FD}_{eff}(\mathit{Re})$ is plotted in (c) for several values of  $\unicode[STIX]{x1D6FD}$ .

5.1 Growth rates and marginal stability

This subsection analyses the dependence of the perturbation growth on Reynolds number $\mathit{Re}$ , spanwise wavenumber $k$ and opening bend ratio  $\unicode[STIX]{x1D6FD}$ . Figure 12 shows the predicted growth rates as a function of the Reynolds number and spanwise wavenumber $k$ for $\unicode[STIX]{x1D6FD}=0.2$ , 0.5, 1 and 2. The primary linear spanwise instability is obtained via polynomial interpolation to determine the lowest Reynolds number that first produces $\unicode[STIX]{x1D70E}=0$ , and the wavenumber at which this occurs. Both are reported in table 4. For $\unicode[STIX]{x1D6FD}>0.2$ , at very low wavenumbers $k\lesssim 0.3$ , a local maximum is observed, but the corresponding mode is always stable. Between this local maximum and the primary maximum, a small range of $k$ produces leading non-real eigenvalues. Larger wavenumbers that those shown in figure 12 were also studied and found to decay increasingly fast at large  $k$ . This trend of the growth rate as a function of $\mathit{Re}$ and $k$ shows a good resemblance with that of backward-facing step flow (Barkley et al. Reference Barkley, Gomes and Henderson2002). However, for all $\unicode[STIX]{x1D6FD}$ , the critical $\mathit{Re}$ values for the flow to become three-dimensional are found to be much lower compared to the flow in backward-facing step (Armaly et al. Reference Armaly, Durst, Pereira and Schonung1983; Barkley et al. Reference Barkley, Gomes and Henderson2002) and partially blocked channel (Griffith et al. Reference Griffith, Thompson, Leweke, Hourigan and Anderson2007). This finding is not surprising, as the two-dimensional flow around the sharp bend becomes unsteady at much lower $\mathit{Re}$ compared to those geometries.

Figure 12. Growth rates of leading eigenmodes as a function of spanwise wavenumber $k$ for (a)  $\unicode[STIX]{x1D6FD}=0.2$ and $\mathit{Re}\leqslant 200$ , (b)  $\unicode[STIX]{x1D6FD}=0.5$ and $\mathit{Re}\leqslant 600$ , (c)  $\unicode[STIX]{x1D6FD}=1$ and $\mathit{Re}\leqslant 700$ and (d)  $\unicode[STIX]{x1D6FD}=2$ and $\mathit{Re}\leqslant 700$ . Solid symbols represent real leading eigenvalues, while hollow symbols represent complex-conjugate pairs of non-real leading eigenvalues. Solid lines connect all dominant leading eigenvalues from several branches of the same Reynolds number.

Figure 12 shows that $\unicode[STIX]{x1D6FD}=0.2$ always has an oscillatory leading mode for the range of $\mathit{Re}$ investigated. Meanwhile, $\unicode[STIX]{x1D6FD}=0.5$ has a synchronous leading mode at $\mathit{Re}\leqslant 400$ , before an oscillatory mode becomes dominant at higher $\mathit{Re}$ . On the other hand, $\unicode[STIX]{x1D6FD}>1$ always has a synchronous leading mode. The transition from oscillatory to synchronous behaviour as $\unicode[STIX]{x1D6FD}$ increases is similar to what was observed by Lanzerstorfer & Kuhlmann (Reference Lanzerstorfer and Kuhlmann2012) in backward-facing step flow, where they found that the flow with very large step height is destabilised by an oscillatory mode, and this changes from oscillatory to synchronous if the step height is further decreased. To clarify the shift of the leading mode from synchronous to oscillatory, several of the leading eigenvalues have been computed for $\unicode[STIX]{x1D6FD}=0.5$ at each wavenumber and three different Reynolds numbers.

The results are shown in figure 13. The curves closely resemble those for the flow over a backward-facing step (Barkley et al. Reference Barkley, Gomes and Henderson2002), which consists of two branches of real eigenvalues at low wavenumber coalescing into a single branch of non-real eigenvalues as $k$ increases. In Barkley et al.’s (Reference Barkley, Gomes and Henderson2002) case, however, the primary leading eigenvalues appear at higher wavenumbers, and (as with the present study for $\unicode[STIX]{x1D6FD}>0.5$ ) the real branch at high wavenumber is the first to become unstable. All branches shift to higher $\unicode[STIX]{x1D70E}$ as $\mathit{Re}$ increases. However, as $\mathit{Re}$ increases further, it can be seen that the leading oscillatory mode is more stable than the real one (figure 13 b) before it becomes dominant at $\mathit{Re}\approx 600$ as shown in figure 13(c). A similar observation was made by Natarajan & Acrivos (Reference Natarajan and Acrivos1993) in flow past spheres and disks, by Tomboulides & Orszag (Reference Tomboulides and Orszag2000) in the weak turbulent flow past a sphere, and by Johnson & Patel (Reference Johnson and Patel1999) in a numerical and experimental study on flow past a sphere up to $\mathit{Re}_{d}=300$ .

Figure 13. Growth rates of leading eigenmodes plotted against spanwise wavenumber at $\unicode[STIX]{x1D6FD}=0.5$ for (a)  $\mathit{Re}=200$ , (b)  $\mathit{Re}=400$ and (c)  $\mathit{Re}=600$ . Solid symbols represent real leading eigenvalues, while hollow symbols represent complex-conjugate pairs of leading non-real eigenvalues.

5.2 Structure of the eigenvalue spectra

Figure 14(a,d,g) show the eigenvalue spectra for three different cases. Each has a different type of leading eigenmode, and the dashed line indicates the onset of instability ( $\unicode[STIX]{x1D70E}=0$ ). The real parts of the eigenmodes are visualised via plots of spanwise vorticity on the plane $z=0$ in figure 14(b,c,e,f,h,i). Figure 14(a) depicts the eigenvalue spectrum for $\unicode[STIX]{x1D6FD}=0.2$ , $\mathit{Re}=120$ and $k=2.8$ , which is near to the onset of instability. Two complex-conjugate pairs of non-real eigenvalues are the fastest growing eigenvalues in the spectrum. The leading pair (e.g. figure 14 b) exhibits a strong growth rate, with strong perturbation structure in the primary recirculation bubble, while the second pair (e.g. figure 14 c) has perturbation structure mainly localised in the bulk flow near the secondary recirculation bubble.

Figure 14. Eigenvalue spectra for (a)  $\unicode[STIX]{x1D6FD}=0.2$ , $\mathit{Re}=120$ and $k=2.8$ , (d)  $\unicode[STIX]{x1D6FD}=0.5$ , $\mathit{Re}=400$ and $k=3.5$ and (g)  $\unicode[STIX]{x1D6FD}=0.5$ , $\mathit{Re}=600$ and $k=4.5$ . (b,c,e,f,h,i) Visualisations of the real part of the complex eigenvector fields via spanwise vorticity on the plane at $z=0$ : (b,e,h) show the respective leading eigenvalues, while (c,f,i) show the corresponding field for the second most dominant eigenmode. In these vorticity plots, zero vorticity is represented by the mid-level shading, while darker and lighter shading respectively show negative and positive vorticity.

In contrast to $\unicode[STIX]{x1D6FD}=0.2$ , larger values of $\unicode[STIX]{x1D6FD}$ tend to favour synchronous leading modes; $\unicode[STIX]{x1D6FD}=0.5$ is particularly interesting because this is the point at which a transition from synchronous ( $\mathit{Re}=400$ in figure 14 d) to oscillatory ( $\mathit{Re}=600$ in figure 14 g) leading mode is seen. The leading eigenvalue for $\mathit{Re}=400$ is synchronous, while the second largest eigenmode is oscillatory. However, at higher $\mathit{Re}$ , the oscillatory mode has higher growth rate compared to the synchronous mode, as can be seen in figure 14(g). The perturbation fields associated with these modes are strongest in the same vicinity, which is in the primary recirculation bubble near to the reattachment point. From the contours represented in figure 14(h,i), we can see that the contours in figure 14(b,e,i) have a similar synchronous eigenmode structure; meanwhile figure 14(c,f,h) have a consistent oscillatory eigenmode structure. The oscillatory mode structure is distinguished from the synchronous mode structure by the presence of an array of chevron-shaped vorticity structures following the core flow downstream from the aft end of the primary recirculation bubble.

5.3 Dependence on $\unicode[STIX]{x1D6FD}$ and analogy with related flows

From the known influence of the expansion ratio on the three-dimensional characteristics of the backward-facing step flow (Barkley et al. Reference Barkley, Gomes and Henderson2002; Lanzerstorfer & Kuhlmann Reference Lanzerstorfer and Kuhlmann2012), the opening ratio in a partially blocked channel flow (Griffith et al. Reference Griffith, Thompson, Leweke, Hourigan and Anderson2007) and the two-dimensional characteristics of 180-degree sharp bend flow (Zhang & Pothérat Reference Zhang and Pothérat2013), we expect that the three-dimensional flow in a 180-degree sharp bend is also dependent on the same physical parameter (here  $\unicode[STIX]{x1D6FD}$ ).

Figure 15. Critical (a) Reynolds number and (b) wavenumber as a function of the opening ratio. In (a) the critical Reynolds number is also plotted against the effective opening ratio $\unicode[STIX]{x1D6FD}_{eff}$ . For $\unicode[STIX]{x1D6FD}$ data, synchronous modes and oscillatory modes are represented respectively by ‘▫’ and ‘ $+$ ’. For $\unicode[STIX]{x1D6FD}_{eff}$ data, synchronous and oscillatory modes are respectively represented by ‘♢’ and ‘ $\times$ ’.

Figure 15 illustrates the effect of $\unicode[STIX]{x1D6FD}$ on $\mathit{Re}$ and $k$ at the onset of three-dimensional instability. In the range of $\unicode[STIX]{x1D6FD}$ studied, only the case $\unicode[STIX]{x1D6FD}=0.2$ becomes three-dimensional through the onset of an oscillatory mode, while all other cases transition through the onset of a non-oscillatory mode. Apparently, $\mathit{Re}_{3D}$ increases steadily as the bend opening becomes larger until $\unicode[STIX]{x1D6FD}\approx 1$ ; $\mathit{Re}_{3D}$ at $\unicode[STIX]{x1D6FD}=2$ is lower than at $\unicode[STIX]{x1D6FD}=1$ because of the appearance of the recirculation bubble at the far end of the bend wall that limits the width of the bulk flow in the bend, causing a reduction in the effective value of $\unicode[STIX]{x1D6FD}$ . If we look closely, the value of $\mathit{Re}_{3D}$ at $\unicode[STIX]{x1D6FD}=2$ is almost the same as that at $\unicode[STIX]{x1D6FD}=0.8$ , which testifies to the importance of $\unicode[STIX]{x1D6FD}_{eff}$ in determining the stability of the flow. In figure 15(a), $\mathit{Re}_{3D}$ is also plotted against $\unicode[STIX]{x1D6FD}_{eff}$ to illustrate this behaviour.

The eigenvector fields at the onset of instability for all synchronous modes ( $\unicode[STIX]{x1D6FD}\gtrsim 0.3$ ) have perturbation structure located in the recirculation bubble, similar to what was found in the backward-facing step flow (Alam & Sandham Reference Alam and Sandham2000; Barkley et al. Reference Barkley, Gomes and Henderson2002), partially blocked channel flow (Griffith et al. Reference Griffith, Thompson, Leweke, Hourigan and Anderson2007) and separation flow (Hammond & Redekopp Reference Hammond and Redekopp1998). Barkley et al. (Reference Barkley, Gomes and Henderson2002) concluded that the size and shape of the bubble directly affected the instability. Conversely, Hammond & Redekopp (Reference Hammond and Redekopp1998) and Alam & Sandham (Reference Alam and Sandham2000) confirmed in their studies that the onset of local absolute instability depended on the backflow of the bubble. Interestingly, $\unicode[STIX]{x1D6FD}=0.8$ and $\unicode[STIX]{x1D6FD}=2$ have comparable bubble size, peak backflow velocity and location of the peak backflow velocity (with discrepancies of only 4.3 %, 1.0 % and 0.4 %, respectively). This supports the view that for $\unicode[STIX]{x1D6FD}>1$ , the stability of the flow is characterised by $\unicode[STIX]{x1D6FD}_{eff}$ ; in both of these cases $\unicode[STIX]{x1D6FD}_{eff}=0.8$ .

The dependence of $k_{c}$ on $\unicode[STIX]{x1D6FD}$ is also illustrated in figure 15. As $\unicode[STIX]{x1D6FD}$ increases from 0.3 to 0.5, the dominant wavelength of the instability increases from $\unicode[STIX]{x1D706}\approx 2\unicode[STIX]{x03C0}/3$ to ${\approx}\unicode[STIX]{x03C0}$ . This is perhaps due to the transition from jet-like flow around the bend to a broader turning flow. The behaviour of $k_{c}$ at $\unicode[STIX]{x1D6FD}=0.2$ does not follow the trend observed between $\unicode[STIX]{x1D6FD}=0.3$ and 0.5, due to the different mechanism. For $\unicode[STIX]{x1D6FD}\geqslant 0.5$ , the critical wavenumber exhibits little dependence on  $\unicode[STIX]{x1D6FD}$ .

Figure 16. Marginal stability curves for sharp 180-degree bend flow with $\unicode[STIX]{x1D6FD}=0.2$ , 0.5, 1 and 2. Regions on the right of the curves represent flow conditions that are linearly unstable to three-dimensional perturbations for that particular  $\unicode[STIX]{x1D6FD}$ .

Figure 17. Length of primary recirculation bubble as a function of (a) critical Reynolds number and (b) critical wavelength.

Figure 16 shows the marginal stability curves at several values of $\unicode[STIX]{x1D6FD}$ . The marginal curves are obtained by interpolating $\unicode[STIX]{x1D70E}$ ( $k,\mathit{Re}$ ) to zero growth rate for each $\mathit{Re}$ from figure 12. For $\unicode[STIX]{x1D6FD}\leqslant 1$ , with increasing $\unicode[STIX]{x1D6FD}$ the neutral stability curve shifts to the right as expected, as the flow with wider bend opening ratio is more stable than those with smaller opening ratios. Beyond $\unicode[STIX]{x1D6FD}\approx 1$ , the stability curve recesses slightly towards lower Reynolds numbers, occupying the region between the marginal stability curves for $\unicode[STIX]{x1D6FD}=0.5$ and 1. This is explained by the decrease in $\unicode[STIX]{x1D6FD}_{eff}$ with $\unicode[STIX]{x1D6FD}$ for $\unicode[STIX]{x1D6FD}>1$ (see § 4). The stability curves also show a decrease in dominant wavenumber with increasing  $\unicode[STIX]{x1D6FD}$ . This is because the instability is scaled with bubble size as discussed by Barkley & Henderson (Reference Barkley and Henderson1996). The length of the primary recirculation bubble as a function of critical Reynolds number and critical wavelength is depicted in figure 17. It can be seen from figure 17(a) that the flow becomes unstable at higher $\mathit{Re}$ at bigger  $\unicode[STIX]{x1D6FD}$ . For the synchronous modes ( $0.3\lesssim \unicode[STIX]{x1D6FD}\lesssim 1$ ), the critical wavelength increases as the size of the primary recirculation bubble increases.

Figure 18. Parameter space for flow regimes and leading peak eigenvalues of three-dimensional stability. The onset of three-dimensional instability is represented by ▵ symbols. Lines are included for guidance and the shaded region shows the parameter values exhibiting unsteady two-dimensional solutions.

In order to consider the three-dimensional stability of these flows in the context of the underlying two-dimensional flows, figure 18 summarises the ( $\mathit{Re},\unicode[STIX]{x1D6FD}$ ) parameter space explored. Within this range, instability to three-dimensional perturbations always occurs in the regime where both primary and secondary recirculation bubbles exist. This agrees with Armaly et al. (Reference Armaly, Durst, Pereira and Schonung1983), who found for the backward-facing step that three-dimensionality appeared in the flow after the secondary recirculation bubble had formed. In our study, $\mathit{Re}_{3D}$ increases monotonically with $\unicode[STIX]{x1D6FD}$ before slightly reducing and becoming independent of $\unicode[STIX]{x1D6FD}$ as $\unicode[STIX]{x1D6FD}$ exceeds unity. The critical eigenvalue for $\unicode[STIX]{x1D6FD}=0.2$ is found to be non-real, while it is real for $\unicode[STIX]{x1D6FD}=0.5$ , 1 and 2. Across all considered $\mathit{Re}$ values studied at $\unicode[STIX]{x1D6FD}=0.2$ , the dominant eigenvalues are non-real. On the other hand, for $\unicode[STIX]{x1D6FD}\gtrsim 1$ , the dominant eigenvalues are consistently real. A transition from real to non-real (from solid to hollow symbol) can be seen at $\unicode[STIX]{x1D6FD}=0.5$ in figure 18, which occurs near the regime where the inside recirculation appears in the primary recirculation bubble.

5.4 Instability mode structure

In this section, the mechanisms by which the three-dimensional infinitesimal perturbations are amplified are addressed. The obvious mechanism seen in the two-dimensional flow is Kelvin–Helmholtz instability. Zhang & Pothérat (Reference Zhang and Pothérat2013) found that in regime IV, the shear layers around both bubbles are subject to it, leading to unsteadiness. However, as this study demonstrates, absolute three-dimensional instabilities are found at $\mathit{Re}\ll \mathit{Re}_{c}$ involving different mechanisms.

As mentioned earlier, the structure of instability affecting the primary recirculation bubble bears a strong similarity to both the flow over a backward-facing step and the flow in a partially blocked channel. Ghia et al. (Reference Ghia, Osswald and Ghia1989) suggested that the appearance of the secondary bubble introduced a concave curvature in the streamlines of the bulk flow, thus inducing Taylor–Görtler instability. However, this scenario has been ruled out (Barkley et al. Reference Barkley, Gomes and Henderson2002; Griffith et al. Reference Griffith, Thompson, Leweke, Hourigan and Anderson2007) because instability arises neither in the secondary recirculation bubble nor in the main bulk flow between the primary and secondary recirculation bubble zones. The leading instability mode from our analysis is also localised in a different location, although the second leading, lower-wavenumber eigenmode is found to be located in these regions.

Figure 19. Visualisation of the real part of leading eigenmodes at (a)  $\unicode[STIX]{x1D6FD}=0.2$ , (b)  $\unicode[STIX]{x1D6FD}=0.5$ , (c)  $\unicode[STIX]{x1D6FD}=1$ and (d)  $\unicode[STIX]{x1D6FD}=2$ , consisting of (i) a three-dimensional visualisation of the streamwise ( $x$ -component of) vorticity, (ii) spanwise ( $z$ -component of) vorticity and (iii) spanwise velocity contours overlaid with the base flow streamlines. Dark, mid and light shading represent negative, zero and positive levels, respectively, while (ii) and (iii) are plotted on the plane $z=0$ .

Figure 19 visualises the real part of the leading eigenmode at gap ratios spanning $0.2\leqslant \unicode[STIX]{x1D6FD}\leqslant 2$ . Three-dimensionality appears at the reattachment and separating points of the primary recirculation bubble as shown by the spanwise velocity component in figure 19(b (iii),c (iii),d (iii)). The perturbation is strongest within the closed streamlines of the primary recirculation bubble (see the isosurface plots in figure 19), with further perturbation structure also present in the secondary bubble and propagating downstream into the core jet flow.

Spanwise perturbation vorticity contour plots shown in figure 19 exhibit perturbation vorticity structures that resemble those arising from an elliptic instability, namely a pair of counter-rotating vortices inside the recirculation bubble. A similar interaction of counter-rotating vortices in two-dimensional elliptic streamlines was seen in Thompson, Leweke & Williamson (Reference Thompson, Leweke and Williamson2001). Leweke & Williamson (Reference Leweke and Williamson1998) suggested that when two counter-rotating vortices balance each other, the radial component of the strain field leads the disturbance to grow exponentially.

Lanzerstorfer & Kuhlmann (Reference Lanzerstorfer and Kuhlmann2012) observed that the combination of the flow deceleration near the reattachment point, a lift-up process on both sides of the bulk flow between the primary and the secondary recirculation bubbles, and an amplification due to streamline convergence near and in the separated flow regions is the cause of the flow instability in the flow over a backward-facing step with an expansion ratio of 0.5. The same observation was made by Wee et al. (Reference Wee, Yi, Annaswamy and Ghoniem2004), who found that the backward-facing step flow was locally absolutely unstable near the middle of the primary recirculation bubble. The backflow was found to be high and the shear layer was sufficiently thick to support an absolutely unstable mode; hence, an absolute mode was more likely to originate in the middle of the bubble. By contrast, Marquillie & Ehrenstein (Reference Marquillie and Ehrenstein2003) studied a flow behind a bump and observed that the structural changes near the reattachment point of the primary recirculation bubble behind the bump triggered an abrupt local transition from convective to absolute instability. It is observed that for the flow around a 180-degree sharp bend, as $\mathit{Re}$ increases, the location of the peak backflow in the primary recirculation bubble shifts towards the reattachment point. The close gap between the peak backflow and the reattachment point means that the flow is strongly decelerated upon approaching the reattachment point. Interestingly, the peak backflow at the onset of instability in this study for $0.3\lesssim \unicode[STIX]{x1D6FD}\lesssim 2$ is consistently located approximately $0.22L_{R}$ from the reattachment point, which is also the same location of the peak perturbation spanwise velocity from the linear stability analysis (figure 19 b (iii),c (iii),d (iii)). This suggests that as in the other geometries mentioned, the instability in the 180-degree sharp bend for $\unicode[STIX]{x1D6FD}\geqslant 0.3$ is localised near the peak of backflow intensity.

The oscillatory critical mode at $\unicode[STIX]{x1D6FD}=0.2$ exhibits strong spanwise velocities at the upstream end of the primary recirculation bubble and quite strong values around the intense vortex in the bubble. This is highly consistent with a mechanism involving a centrifugal instability around the intense vortex. A similar mode was seen by Lanzerstorfer & Kuhlmann (Reference Lanzerstorfer and Kuhlmann2012) in the flow over a backward-facing step with a very small opening. There, the perturbations were found to be a spanwise travelling wave that displaces the jet and the intense vortex periodically.

The structure of the mode that destabilises the flow at $\unicode[STIX]{x1D6FD}=0.5$ , $\mathit{Re}=278$ and $k=2$ is shown in figure 19(b) in the isosurface plots of streamwise vorticity, spanwise vorticity and $w$ velocity contours. The isosurface consists of positive (light) and negative (dark) vorticity contours located almost entirely in the primary recirculation bubble, near the separation and reattachment points. The spanwise vorticity contour plot shows that there is a pair of counter-rotating vortices in the primary recirculation bubble which resemble the flow in a partially blocked channel (Griffith et al. Reference Griffith, Thompson, Leweke, Hourigan and Anderson2007). The spanwise velocity contours are also qualitatively similar to those of the unstable mode in the flow over a backward-facing step (Barkley et al. Reference Barkley, Gomes and Henderson2002), where a ‘flat roll’ (i.e. a roll in a horizontal plane about a vertical axis) mode structure exists in the bubble near the reattachment point. Across the opening ratios studied, the same mode structure has been found for all real primary leading eigenmodes indicated by the solid symbols in figure 18.

The same type of plots describing the dominant eigenmodes for $\unicode[STIX]{x1D6FD}=0.2$ , $\mathit{Re}=123$ and $k=2.8$ are shown in figure 19(a). The mode appears to grow in the primary recirculation bubble near the separation point and upstream of the intense vortex close to the reattachment point. Both of these modes have a strong $x$ -component of perturbation vorticity in the upstream part of the primary recirculation bubble, which is where Zhang & Pothérat (Reference Zhang and Pothérat2013) found secondary instability in their three-dimensional simulation of an unsteady flow at $\mathit{Re}=2000$ and $\unicode[STIX]{x1D6FD}=1$ with a spanwise-periodic domain of length 2 units (this is equivalent to the case of wavenumber $k=\unicode[STIX]{x03C0}$ in the present notation).

Figure 20. Structure of the real part of the leading eigenmode at $\unicode[STIX]{x1D6FD}=0.5$ , $\mathit{Re}=400$ and $k=0.2$ depicted on the plane $z=0$ : flooded contours of (a) streamwise vorticity and (b)  $u$ -, (c)  $v$ - and (d)  $w$ -velocity overlaid with the base flow streamlines. For clarity, only the vicinity of the bend is shown. Dark, mid and light shading represent negative, zero and positive levels, respectively.

The structure of the most unstable eigenmode at very small wavenumber is shown in figure 20 as a plot of spanwise vorticity and $(u,v,w)$ velocity contours. The structure consists of spanwise vortices in the main bulk flow located between the primary and secondary recirculation bubbles. The spanwise velocity contour in figure 20(d) clearly shows that the bifurcating mode is located in the main bulk flow near the closed streamlines of both bubbles.

5.5 Two-dimensional instability

Section 5.1 described how for $\mathit{Re}<\mathit{Re}_{c}$ at any $\unicode[STIX]{x1D6FD}$ , the flow is stable to two-dimensional infinitesimal perturbations ( $k=0$ ). The linear stability analysis found that all eigenvalues for $k=0$ are real and have negative growth rate. The two-dimensional flow around a sharp 180-degree bend becomes unsteady at small $\mathit{Re}$ depending on $\unicode[STIX]{x1D6FD}$ , which is almost twice the critical value for the onset of three-dimensional instability. From two-dimensional simulations, Zhang & Pothérat (Reference Zhang and Pothérat2013) observed that two-dimensional instability starts in the shear layer between the two steady recirculation bubbles and sheds throughout the whole width of the channel. This agrees with the mechanism found for the two-dimensional leading eigenmode. As $k\rightarrow 0$ , a leading real eigenmode splits into two complex-conjugate pairs. Figure 21 is plotted similarly to figure 13 in Barkley et al. (Reference Barkley, Gomes and Henderson2002); $\log (-\unicode[STIX]{x1D70E})$ will approach $-\infty$ as $\unicode[STIX]{x1D70E}\rightarrow 0$ from below. However, the data here show that the eigenmode growth rate remains negative as the critical Reynolds number is approached, lacking any evidence of a departure towards $-\infty$ . While it was not possible due to computing time limitations to obtain stability data closer to $\mathit{Re}_{c,upper}$ , this nevertheless suggests that the transition to unsteadiness of the steady-state two-dimensional solution branch is not due to a global linear instability. This is supported by our earlier observation that $\mathit{Re}_{c,upper}$ was resolution sensitive. Similar observations to those shown in figure 21 were made for flow behind a backward-facing step (Barkley et al. Reference Barkley, Gomes and Henderson2002) and through a partially blocked channel (Griffith et al. Reference Griffith, Thompson, Leweke, Hourigan and Anderson2007). The data exhibit a kink at $\mathit{Re}\approx 700$ , where a sudden steepening in gradient is observed, before the data revert to a nearly horizontal trend to higher Reynolds numbers. Barkley et al. (Reference Barkley, Gomes and Henderson2002) observed a similar behaviour for the backward-facing step flow (in that case occurring at $\mathit{Re}\approx 1250$ ), and found by inspection of secondary eigenvalues that the kink occurred due to an ‘avoided crossing’ of the two eigenmode branches.

Figure 21. Leading two-dimensional ( $k=0$ ) eigenvalues for $\unicode[STIX]{x1D6FD}=1$ . Due to the hysteresis, there are two $\mathit{Re}_{c}$ for this flow: $\mathit{Re}_{c,lower}$ is the lowest $\mathit{Re}$ at which the flow can become unsteady, while $\mathit{Re}_{c,upper}$ is the highest $\mathit{Re}$ for which the flow remains steady. The leading eigenvalues remain to have finite value of $\log (-\unicode[STIX]{x1D70E})$ as $\mathit{Re}\rightarrow \mathit{Re}_{c}$ . Note that the flow will become linearly neutrally stable when $\unicode[STIX]{x1D70E}\rightarrow 0^{-}$ and $\log (-\unicode[STIX]{x1D70E})\rightarrow -\infty$ .

6.1 Nonlinear evolution of the unstable modes

Figure 23. Visualisation of the three-dimensional disturbances via isosurface plots of the (streamwise) $x$ -component of vorticity for (a)  $\unicode[STIX]{x1D6FD}=0.2$ , $\mathit{Re}=160$ , $k=4$ and (b)  $\unicode[STIX]{x1D6FD}=1$ , $\mathit{Re}=600$ , $k=4.5$ . In each case (i) shows the leading eigenmode predicted by the linear stability analysis and (ii) shows the saturated state of a three-dimensional DNS solution. The saturated solution in panel (a (ii)) is oscillatory, and that in (b (ii)) is steady state.

Three-dimensional simulations were subsequently performed at selected $\unicode[STIX]{x1D6FD}$ and $\mathit{Re}$ combinations. The spanwise wavenumber in each case is deliberately set to match the corresponding linear instability eigenmode. It is acknowledged that this choice excludes long-wavelength features that may or may not arise. However, it facilitates an isolation of the instability mode under scrutiny. Figure 22 shows the time history of the spanwise velocity in these three-dimensional simulations for (a)  $\unicode[STIX]{x1D6FD}=0.2$ , $\mathit{Re}=160$ , $k=4$ and (b)  $\unicode[STIX]{x1D6FD}=1$ , $\mathit{Re}=600$ , $k=4.5$ . The oscillatory and the synchronous behaviour in figures 22(a) and 22(b), respectively, agree well with the behaviour of the leading eigenmodes predicted using linear stability analysis (figure 12). From the time history of the $w$ -velocity, the growth rate of the perturbation can be calculated. Table 6 shows strong agreement between the growth rates of the perturbation field obtained from the linear stability analysis and from the three-dimensional DNS for five different cases. The combination $\unicode[STIX]{x1D6FD}=2$ and $k=4.5$ is chosen to demonstrate the accuracy of the predictions on flow with fast-growing perturbations. Figure 23 shows three-dimensional isosurface plots of streamwise vorticity for the same parameters as in figure 22, comparing the predicted three-dimensional eigenmode in each case with the actual three-dimensional state produced once the flow saturates following instability growth. The streamwise vorticity from the linear stability analysis (figure 23 a,b (i)) has a strong resemblance to that obtained from the three-dimensional DNS (figure 23 a,b (ii)). The strong agreement seen between the predicted eigenmode structure and the resulting saturated three-dimensional structure verifies that the linear stability analysis provides meaningful predictions of the three-dimensional nature of the flow. The Reynolds numbers in figures 23(a) and 23(b) are 28 % and 51 % higher than the critical Reynolds numbers for $\unicode[STIX]{x1D6FD}=0.2$ and 1, respectively. Both cases produce non-zero Fourier mode energy at saturation of magnitude $10^{-2}$ relative to the base flow energy, which is very small. The smaller the disturbance energy compared to the base flow energy, the closer the saturated state will be to the predicted infinitesimal eigenmode because the contribution of nonlinear terms is weaker.

From the three-dimensional DNS, we found that the unsteady saturated state for $\unicode[STIX]{x1D6FD}=0.2$ persists to at least $\mathit{Re}=160$ ; three-dimensional simulations were not conducted beyond $\mathit{Re}=160$ , though it might be anticipated based on the close agreement between the saturated three-dimensional flow structure and the corresponding predicted linear instability mode that a similar behaviour would extend to higher Reynolds numbers: linear stability analysis was performed up to $\mathit{Re}=200$ , continuing to capture this mode. Meanwhile, for $\unicode[STIX]{x1D6FD}=0.5$ , as $\mathit{Re}$ increases, the saturated state changed from a steady state at $\mathit{Re}=300$ and 400 to an unsteady saturated state at $\mathit{Re}=500$ with $k=4.5$ . The same is predicted in figure 12, demonstrating the applicability of linear stability analysis to this flow and the rich tapestry of flow regimes across the $\mathit{Re}$ – $\unicode[STIX]{x1D6FD}$ parameter space.

6.2 Stuart–Landau model analysis

In this subsection, analysis of the nonlinear features of the instability mode evolution is performed using a truncated Stuart–Landau equation. The Stuart–Landau model is valid in the vicinity of the transition Reynolds number, and has found wide application for classification of the nonlinear characteristics of bifurcations in fluid flows. Examples include analysis of the Hopf bifurcation from steady-state flow past a circular cylinder producing the classical Kármán vortex street (Provansal, Mathis & Boyer Reference Provansal, Mathis and Boyer1987; Duŝek, Le Gal & Fraunié Reference Duŝek, Le Gal and Fraunié1994; Schumm, Berger & Monkewitz Reference Schumm, Berger and Monkewitz1994; Albarède & Provansal Reference Albarède and Provansal1995; Thompson & Le Gal Reference Thompson and Le Gal2004), the regular (steady-to-steady) bifurcation breaking axisymmetry in the flow behind a sphere (Thompson et al. Reference Thompson, Leweke and Williamson2001) and three-dimensional transition behind a cylinder (Henderson & Barkley Reference Henderson and Barkley1996; Sheard, Thompson & Hourigan Reference Sheard, Thompson and Hourigan2003), staggered cylinders (Carmo et al. Reference Carmo, Sherwin, Bearman and Willden2008) and rings (Sheard, Thompson & Hourigan Reference Sheard, Thompson and Hourigan2004a ,Reference Sheard, Thompson and Hourigan b ). The model describes the growth and saturation of perturbation as (Landau & Lifshitz Reference Landau and Lifshitz1976)

where $A$ is the complex amplitude of the evolving instability as a function of time, and the right-hand side of the equation represents the first two terms of a series expansion. The growth rate and angular frequency of the mode in the linear regime ( $|A|\rightarrow 0$ ) are respectively denoted by $\unicode[STIX]{x1D70E}$ and $\unicode[STIX]{x1D714}$ , while weakly nonlinear properties are determined by the second term on the right-hand side. The sign of $l$ dictates whether the mode evolution is via a supercritical ( $l>0$ ) or subcritical ( $l<0$ ) bifurcation, and any frequency shift is described by the Landau constant $c$ (Duŝek et al. Reference Duŝek, Le Gal and Fraunié1994; Le Gal, Nadim & Thompson Reference Le Gal, Nadim and Thompson2001; Thompson et al. Reference Thompson, Leweke and Williamson2001). A common treatment (Le Gal et al. Reference Le Gal, Nadim and Thompson2001; Thompson et al. Reference Thompson, Leweke and Williamson2001) is to decompose $A(t)$ into magnitude and phase components as

where $\unicode[STIX]{x1D70C}(t)=|A(t)|$ and $\unicode[STIX]{x1D719}(t)=\arg (A(t))$ . Substitution of (6.2) into (6.1), separation into real and imaginary parts, and simplification yields separate equations for amplitude and phase, i.e.

It is convenient to then manipulate (6.3) as

Hence a positive slope ( $-l$ ) in a plot of $\text{d}(\log |A|)/\text{d}t$ against $|A|^{2}$ will indicate a subcritical bifurcation, while a negative slope corresponds to a supercritical bifurcation. Whether a supercritical bifurcation is of pitchfork or Hopf type is dependent on whether the growing mode is synchronous or oscillatory.

Figure 24. The time derivative of the mode amplitude logarithm plotted against the square of the amplitude for (a)  $\unicode[STIX]{x1D6FD}=0.2$ , $\mathit{Re}=160$ , $k=4$ and (b)  $\unicode[STIX]{x1D6FD}=0.5$ , $\mathit{Re}=400$ , $k=3.5$ , demonstrating supercritical behaviour. The solid circle represents the growth rate predicted by linear stability analysis.

The envelope of energy in the non-zero spanwise Fourier modes of the three-dimensional simulations was taken as a global amplitude measure ( $|A|^{2}$ ) of the growing three-dimensional instabilities. Figure 24 plots the time derivative of amplitude logarithm against the square of the amplitude for two different cases: (a)  $\unicode[STIX]{x1D6FD}=0.2$ , $\mathit{Re}=160$ and $k=4$ , which grows from an oscillatory mode, and (b)  $\unicode[STIX]{x1D6FD}=0.5$ , $\mathit{Re}=400$ and $k=3.5$ , which grows from a synchronous mode. The behaviour shown for $\unicode[STIX]{x1D6FD}=0.5$ is consistent with that found at larger  $\unicode[STIX]{x1D6FD}$ . The nearly linear variation with a negative gradient towards small $|A|^{2}$ shown in both plots indicates that the transition in both cases occurs through a supercritical bifurcation. This is consistent with other types of confined flow featuring recirculation bubbles, such as the flow past a backward-facing step (Kaiktsis, Karniadakis & Orszag Reference Kaiktsis, Karniadakis and Orszag1991), through a sudden expansion in a circular pipe (Mullin et al. Reference Mullin, Seddon, Mantle and Sederman2009) and past a sphere (Tomboulides & Orszag Reference Tomboulides and Orszag2000). The oscillatory and synchronous behaviour of the modes in (a) and (b) reveal them to occur through supercritical Hopf and pitchfork bifurcations, respectively.

In addition to the aforementioned sub- and supercritical bifurcation scenarios, another possibility is that of a transcritical bifurcation. A dynamical system producing transcritical behaviour takes the form

which differs from the amplitude part of the Stuart–Landau model (6.3) by the replacement of $\unicode[STIX]{x1D70C}^{3}$ with $\unicode[STIX]{x1D70C}^{2}$ in the nonlinear term on the right-hand side. Under an analogous manipulation,

The $\sqrt{\unicode[STIX]{x1D70C}^{2}}$ term indicates that an infinite gradient would occur at $|A|^{2}=0$ in a plot of $\text{d}(\log |A|)/\text{d}t$ against $|A|^{2}$ . The absence of such behaviour in figure 24 rules out a transcritical bifurcation. Similarly, a transcritical dynamical system expressed in terms of the complex amplitude  $A$ ,

features a real component that simplifies to

In the limit $|A|=\unicode[STIX]{x1D70C}\rightarrow 0$ , the oscillation described by the trigonometric term produces gradients in $\text{d}(\log |A|)/\text{d}t$ as a function of $|A|^{2}$ that approach infinity; the absence of this behaviour in figure 24 again supports the present classification of these bifurcations as being supercritical.