1. Introduction
1.1. Axisymmetric three-dimensional bluff bodies
The flow past a stationary isolated axisymmetric body may be considered as a simplified case of a more general family of immersed three-dimensional (3-D) bluff bodies, which are ubiquitous in human life and engineering applications. In spite of the symmetries of the body, instabilities are known to generate asymmetric and possibly unsteady flows; different regimes are indeed possible depending on the Reynolds number
${{Re}} = U_\infty H / \nu$
, based on the body cross-stream dimension
$H$
, the free-stream velocity
$U_\infty$
and the fluid kinematic viscosity
$\nu$
. In this work, we focus on the primary symmetry-breaking instability.
The flow past 3-D axisymmetric bodies becomes asymmetric before transitioning to an unsteady regime (Magarvey & Bishop Reference Magarvey and Bishop1961a
,
Reference Magarvey and Bishopb
; Magarvey & MacLatchy Reference Magarvey and MacLatchy1965). The steady axisymmetric flow past a sphere, for example, is known to transition to a steady asymmetric state at
${Re} \approx 211$
(Johnson & Patel Reference Johnson and Patel1999) through the regular bifurcation of an eigenmode of azimuthal wavenumber
$m=1$
(Tomboulides & Orszag Reference Tomboulides and Orszag2000). This bifurcation gives origin to a pair of counter-rotating steady streamwise vortices in the wake, which are not aligned with the flow but possess a reflectional symmetry about a longitudinal plane of arbitrary azimuthal orientation, and exhibit an eccentricity that increases with the distance from the body. The same bifurcation has been observed for the flow past other axisymmetric bodies, e.g. disks and bullet-shaped bodies (slender cylindrical bodies with a smooth leading edge and a blunt trailing edge). The flow past a disk placed perpendicular to the incoming flow exhibits the primary regular bifurcation towards an asymmetric state at Reynolds number
${Re} \approx 115$
(Natarajan & Acrivos Reference Natarajan and Acrivos1993; Fabre, Auguste & Magnaudet Reference Fabre, Auguste and Magnaudet2008; Meliga, Chomaz & Sipp Reference Meliga, Chomaz and Sipp2009a
). The flow past a bullet-shaped body undergoes the same bifurcation, but at a larger Reynolds number (Brucker Reference Brucker2001; Bohorquez et al. Reference Bohorquez, Sanmiguel-Rojas, Sevilla, Jiménez-González and Martínez-Bazán2011) that increases with the length-to-diameter aspect ratio 
of the body, e.g.
${Re} \approx 216$
for 
, and
${Re} \approx 435$
for 
. A similar dependence on 
has been observed by Chrust, Bouchet & Dusek (Reference Chrust, Bouchet and Dusek2010) for the flow past flat cylinders whose axis is parallel to the free stream. They observed that the flow bifurcates towards an asymmetric state at
${Re} \approx 115{-}120$
for 
, and at
${Re} \approx 270$
for 
.
The instability mechanism of the primary bifurcation of the flow past axisymmetric bodies has been investigated extensively over the years. Monkewitz (Reference Monkewitz1988) investigated the linear stability of an analytic two-parameter family of axisymmetric, parallel and incompressible wake profiles. He observed that the first helical mode with
$m=1$
displays the largest growth rate for all cases, and that it is the only mode to become absolutely unstable for velocity profiles approximating those found in the near wake. Pier (Reference Pier2008) studied the local absolute instability features of the flow past a sphere, with the aim of linking the local flow properties with the fundamental mechanism driving the global flow bifurcation. He neglected the strong non-parallelism of the near-wake region, and highlighted the local properties by freezing the flow at different streamwise coordinates, and studying the equivalent axially parallel shear flows. In doing so, he demonstrated the existence of absolutely unstable regions in the near wake, and found that the strength and spatial extent of these regions increase with the Reynolds number. The results of Monkewitz (Reference Monkewitz1988) and Pier (Reference Pier2008) clearly hint at the existence of a wavemaker in the wake (an absolutely unstable region where the fluctuations self-sustain), which feeds the convectively unstable region downstream. However, these results cannot fully explain the symmetry breaking corresponding to the first bifurcation of the flow.
1.2. The Magnaudet and Mougin model for free-slip bodies
In the context of free-slip oblate spheroidal bodies, Magnaudet & Mougin (Reference Magnaudet and Mougin2007) (hereafter referred to as MM) proposed an instability mechanism that accounts for the strongly non-parallel nature of the flow. The theoretical arguments that they put forward are based on the idea that the regular bifurcation of the flow is driven by the azimuthal vorticity
$\omega _\theta = \partial u_r/\partial z - \partial u_z/\partial r$
generated at the body surface and then transported into the wake (where
$(u_r,u_z)$
are the velocity components in the radial and axial directions
$r$
and
$z$
). For free-slip ellipsoids, in fact, they found that the wider range of
${Re}$
for which the flow becomes unstable as 
increases (see figure 7 of their paper) well correlates with the increase in maximum surface vorticity
$\omega _{\theta ,max }$
and normal diffusive vorticity flux at the body surface (see § 3 of their paper). They suggested that the instability originates in a thin region of the flow downstream the body. By inspecting the vorticity field, they observed that two distinct regions may be identified (see figure 21 of their paper and figure 1): (i) very close to the body, isocontours of
$\omega _\theta$
are almost parallel to the free-slip surface; (ii) farther downstream in the wake recirculation region, instead, isocontours of
$\omega _\theta$
are almost parallel to the symmetry axis, because of the tendency of
$\omega _\theta /r$
to approach a constant value as
${Re}$
increases. At large Reynolds numbers, indeed, the spatial distribution of the vorticity in the recirculation region of the steady viscous flow past an axisymmetric body resembles that of a Hill spherical vortex (Hill Reference Hill1894; Batchelor Reference Batchelor1967); see e.g. figure 10 in Fornberg (Reference Fornberg1988) and the related discussion. Thus MM conjectured that the instability arises in the transition region where isocontours of
$\omega _\theta$
have to turn to satisfy the conditions both at the body surface and in the wake recirculation region. As
${Re}$
increases, this transition region shrinks, and isocontours of
$\omega _\theta$
must turn more and more sharply, resulting in stronger and stronger vorticity gradients.
Figure 1. Near-wake distribution of the azimuthal vorticity around an ellipsoid with free-slip surface and 
for
$(a)$
${Re}=200$
and
$(b)$
${Re}=270$
. Black lines are isocontours of the azimuthal vorticity
$\omega _\theta$
, and coloured contours show
$\partial \omega _\theta /\partial r$
. The grey dashed line shows
$u_{z} = 0$
, and the red dashed line shows
$\partial \omega _\theta /\partial r = 0$
. Blue/red circles indicate the negative/positive maxima of
$\partial \omega _\theta /\partial r$
in the near-wake region.
In an attempt to provide a rational instability criterion, MM related the flow bifurcation to the appearance of points where
$\partial \omega _\theta / \partial r = 0$
in the near wake. This comes directly from the balance equation for the steady axisymmetric base-flow vorticity, which in cylindrical coordinates reads
\begin{equation} \underbrace { u_r \dfrac {\partial \omega _\theta }{\partial r } }_{\mathrm{I}} \underbrace {{} - u_r \dfrac {\omega _\theta }{r} }_{\mathrm{II}} + \underbrace { u_z \dfrac { \partial \omega _\theta }{\partial z} }_{\mathrm{III}} = \underbrace { \dfrac {1}{{Re}} \dfrac {\partial }{\partial r} \left ( \dfrac {1}{r} \dfrac { \partial ( r \omega _\theta ) }{\partial r} \right )}_{\mathrm{IV}} + \underbrace { \dfrac {1}{{Re}} \dfrac {\partial ^2 \omega _\theta }{\partial z^2} }_{\mathrm{V}}. \end{equation}
We focus on the near-wake region where
$\omega _\theta \lt 0$
,
$u_z\lt 0$
and
$u_r\gt 0$
(see figure 1). Here, terms
$\mathrm{II}$
and
$\mathrm{III}$
on the left-hand side are positive everywhere. As long as
$\omega _\theta$
becomes more negative with
$r$
(i.e.
$\partial \omega _\theta /\partial r \lt 0$
), term
$\mathrm{I}$
on the left-hand side is negative and may, at least partially, balance terms
$\mathrm{II}$
and
$\mathrm{III}$
. When instead an isocontour of
$\omega _\theta$
becomes perpendicular to the symmetry axis (i.e.
$\partial \omega _\theta /\partial r=0$
), term
$\mathrm{I}$
vanishes and the positive left-hand side can only be balanced by the positive viscous terms
$\mathrm{IV}$
and
$\mathrm{V}$
. Note that where
$\partial \omega _\theta /\partial r = 0$
, term
$\mathrm{IV}$
reduces to
$- \omega _\theta /r^2/{Re} \sim {Re}^{-1}$
. Therefore, as the Reynolds number increases, term
$\mathrm{V}$
is observed to dominate term
$\mathrm{IV}$
, and the balance equation reduces to
\begin{equation} - u_r \frac { \omega _\theta }{r} + u_z \frac { \partial \omega _\theta }{\partial z} \sim \frac {1}{{Re}} \frac { \partial ^2 \omega _\theta }{\partial z^2}. \end{equation}
This implies that to balance the left-hand side, the streamwise gradient of
$\omega _\theta$
in term
$\mathrm{V}$
on the right-hand side has to vary more and more sharply (over a shorter distance along
$z$
) as
${Re}$
increases. The MM idea is that this tendency of the vorticity to become discontinuous in the points where
$\partial \omega _\theta /\partial r = 0$
underlies the instability mechanism.
Notably, the numerical simulations of MM corroborate their arguments, showing a good correlation between the first appearance of points with
$\partial \omega _\theta /\partial r = 0$
and the onset of the bifurcation in the context of free-slip axisymmetric bluff bodies. However, despite some hints provided by the same authors (see the concluding discussion in their paper), it is not clear whether axisymmetric bodies with free-slip and no-slip surfaces share this same instability mechanism driven by the vorticity generated at the surface. It is well known, indeed, that the mechanism of vorticity generation on a surface changes with the boundary condition; see the seminal works of Truesdell (Reference Truesdell1954), Lighthill (Reference Lighthill1963), Morton (Reference Morton1984), Leal (Reference Leal1989), Wu & Wu (Reference Wu and Wu1993), Wu (Reference Wu1995) and Lundgren & Koumoutsakos (Reference Lundgren and Koumoutsakos1999), and the more recent works of Brøns et al. (Reference Brøns, Thompson, Leweke and Hourigan2014) and Terrington, Hourigan & Thompson (Reference Terrington, Hourigan and Thompson2020). For no-slip surfaces, the normal diffusion vorticity flux depends on the tangential shear stress and on the pressure gradient, and it is non-null also in the limit case of a flat wall (Wu & Wu Reference Wu and Wu1993). On a free-slip surface, instead, the vorticity appears as a consequence of the continuity of the tangent stresses, and is non-null only in the case of curved surfaces (Wu Reference Wu1995). Thus the generality of the MM mechanism and the influence of the different vorticity generation mechanisms at no-slip and free-slip surfaces on the primary flow bifurcation are still open questions that deserve further investigation.
1.3. Focus of the present work
To address these questions, we study the flow past no-slip axisymmetric bluff bodies with different shapes and aspect ratios, and focus on the primary symmetry-breaking bifurcation. We consider four geometries – ellipsoids, bullets, cones and bicones – that feature different combinations of smooth/sharp leading edge and smooth/blunt trailing edge (see figure 2). The aspect ratio of the ellipsoid and the bicone is varied over 
, and that of the bullet and the cone over 
.
Figure 2. Sketch of the considered geometries. Here, 
.
We give evidence that the MM instability mechanism captures rather well the primary instability of the flow past axisymmetric bodies with different shapes and surface types. Indeed, we find that the onset of the bifurcation correlates well with the first appearance of points where
$\partial \omega _\theta /\partial r = 0$
in the near wake, not only for free-slip surfaces but also for no-slip surfaces, despite the different vorticity generation mechanism at the surface. As first hinted by MM, this shows that although the vorticity is at the root of the instability, the way it is produced at the wall does not play a major role. We then propose a simple yet effective scaling law that may be used to predict the onset of the instability without the need for an additional stability analysis. Indeed, the large variability of the critical Reynolds number with the geometry (e.g. as mentioned above,
${Re}_c \approx 115$
for the disk,
${Re}_c \approx 211$
for the sphere, and
$ 126 \leqslant {Re}_c \leqslant 435$
for bullet-shaped bodies with 
) clearly shows that the free-stream velocity
$U_\infty$
and the cross-stream dimension
$H$
of the body are not the appropriate velocity and length scales for describing this bifurcation. In the same spirit as Chiarini, Quadrio & Auteri (Reference Chiarini, Quadrio and Auteri2022) in the context of two-dimensional (2-D) symmetric bodies, the new scaling is based on quantities related to the physics of the problem, i.e. measures of the vorticity at the body surface and of the radial extent of the wake recirculation region. We show that the resulting Reynolds number, when evaluated at the onset of the bifurcation, collapses rather well to a constant value for all the considered bodies, and thus provides a simple criterion to assess the flow stability.
2. Methods
2.1. Problem formulation
We investigate the primary symmetry-breaking bifurcation of the incompressible flow past 3-D axisymmetric bluff bodies with different geometries and aspect ratios. The bodies have length
$L$
and maximum diameter
$H$
, and are placed in a uniform flow with velocity
$U_\infty$
aligned with their symmetry axis (see figure 2). Standard cylindrical coordinates are used (with
$r$
,
$\theta$
and
$z$
as the radial, azimuthal and axial directions), and the origin is placed at the leading edge of the bodies. The flow is governed by the incompressible Navier–Stokes (NS) equations for the velocity and pressure fields
$\{\boldsymbol{u},p\}$
:
\begin{equation} \frac {\partial \boldsymbol{u}}{\partial t} + \boldsymbol{u \cdot \nabla } \boldsymbol{u} = - \nabla p + \frac {1}{{Re}} \nabla^2 \boldsymbol{u}, \quad \boldsymbol{\nabla \cdot u} = 0, \end{equation}
with
${Re} = U_\infty H /\nu$
. Unless otherwise stated, all quantities are made dimensionless using
$U_\infty$
and
$H$
.
We consider four different geometries that yield a base flow with distinct features: ellipsoids, bullets, cones and bicones (see figure 2). The aspect ratio of the bodies 
is varied over 
for ellipsoids and bicones, and over 
for bullets and cones. Ellipsoids and bullets possess a smooth leading edge, where the flow separation is driven by the adverse pressure gradient and the position of the separation point changes with the Reynolds number; conversely, cones and bicones feature a zero-thickness sharp leading edge. Bullets and cones have a blunt trailing edge (TE), where the separation point is enforced by the geometry, and the cross-stream size of the wake recirculation region does not depend on
${Re}$
and 
; by contrast, ellipsoids and bicones have a zero-thickness TE, where the flow separation is not set by the geometry, and the size of the wake recirculation region varies with
${Re}$
and 
.
2.2. Linear stability
The onset of the primary instability is studied using linear theory (Theofilis Reference Theofilis2003, Reference Theofilis2011). The field
$\{\boldsymbol{u},p\}$
of velocity and pressure is divided into a time-independent axisymmetric base flow
$\{\boldsymbol{u}_0,p_0\}$
and an unsteady perturbation
$\{\boldsymbol{u}_1,p_1\}$
of small amplitude
$0\lt \epsilon \ll 1$
:
\begin{equation} \boldsymbol{u}(\boldsymbol{x},t) = \boldsymbol{u}_0(\boldsymbol{x}) + \epsilon\, \boldsymbol{u}_1(\boldsymbol{x},t) \quad \text{and} \quad p(\boldsymbol{x},t) = p_0(\boldsymbol{x}) + \epsilon\, p_1(\boldsymbol{x},t). \end{equation}
Introducing this decomposition in the NS equations, one obtains at order
$\epsilon ^0$
the steady nonlinear NS equations for the base flow
$\{\boldsymbol{u}_0,p_0\}$
. At order
$\epsilon ^1$
, one obtains the linearised NS equations for the perturbation field
$\{\boldsymbol{u}_1,p_1\}$
. By using a normal mode ansatz and a Fourier transform in the azimuthal direction, for each mode the perturbation field takes the form
\begin{equation} \{{\boldsymbol{u}}_{1},{p}_{1}\} (\boldsymbol{x},t)= \{ \hat {\boldsymbol{u}}_{1}, \hat {p}_{1} \}(r,z)\, \mathrm{e}^{\lambda t + \mathrm{i} m \theta } + \text{c.c.}, \end{equation}
where
$m$
is the azimuthal wavenumber, and
$\text{c.c.}$
designates the complex conjugate terms. Introducing this expansion in the linearised NS equations yields, for each
$m$
, an eigenvalue problem for the complex eigenvalue
$\lambda = \lambda _r+ i \lambda _i$
and the complex eigenvector
$\{\hat {\boldsymbol{u}}_1,\hat {p}_1\}$
:
\begin{equation} \lambda \hat {\boldsymbol{u}}_1 + \mathcal{L}_m\{\boldsymbol{u}_0,{Re}\}\, \hat {\boldsymbol{u}}_1 + \nabla _m \hat {p}_1 = 0, \quad \nabla _m \boldsymbol{\cdot} \hat {\boldsymbol{u}}_1=0. \end{equation}
Here,
$\nabla _m$
is the gradient operator relative to the azimuthal wavenumber
$m$
, and
$\mathcal{L}_m$
stands for the Fourier-transformed linearised NS operator
\begin{equation} \mathcal{L}_m\{\boldsymbol{u}_0,{Re}\}\,\hat {\boldsymbol{u}}_1 = \mathcal{C}_m ( \hat {\boldsymbol{u}}_1, \boldsymbol{u}_0 ) - \frac {1}{{Re}} \nabla_m^2 \hat {\boldsymbol{u}}_1, \end{equation}
where
\begin{equation} \mathcal{C}_m( \boldsymbol{u}_A, \boldsymbol{u}_B ) = \boldsymbol{u}_A \boldsymbol{\cdot} \nabla _m \boldsymbol{u}_B + \boldsymbol{u}_B \boldsymbol{\cdot} \nabla _m \boldsymbol{u}_A. \end{equation}
The flow stability is ascertained by the solution of the generalised eigenvalue problem (2.4) for the complex frequency
$\lambda$
. When
$\lambda _r\lt 0$
, the flow is linearly stable, while when
$\lambda _r\gt 0$
, the associated global mode is linearly unstable and grows exponentially in time. When
$\lambda _i \neq 0$
, the unstable mode is time-dependent. When
$m \neq 0$
, the unstable mode is modulated in the azimuthal direction. Since the focus of the present work is on the primary symmetry-breaking regular bifurcation of the flow past axisymmetric bodies, we restrict our analysis to
$m=1$
.
2.3. Numerical method
The analysis is based on finite elements, and the simulations have been carried out using the COMSOL Multiphysics software (Multiphysics Reference Multiphysics1998).
Given the axisymmetric nature of the bodies, we consider only the
$(r,z)$
plane, and build the mesh in the numerical domain
$\{r,z \ | \ 0 \leqslant r \leqslant 50; -50 \leqslant z \leqslant 150\}$
, with nodes strongly clustered near the body. We employ a quadrilateral mesh with eight layers and growth factor
$1.2$
in the vicinity of the body, and a triangular mesh in the remaining part of the domain. Several sub-domains are used to control the mesh size. The element size decreases from
$1$
in the far field to
$1/120$
on the body surface; at the corners, the element size is approximately
$10^{-5}$
. The number of elements changes with the geometry and the aspect ratio. For the shortest and longest bodies, the numbers of elements are
$79 \times 10^3$
and
$126 \times 10^3$
for ellipsoids,
$83 \times 10^3$
and
$146 \times 10^3$
for bicones,
$74 \times 10^3$
and
$122 \times 10^3$
for bullets, and
$77 \times 10^3$
and
$143 \times 10^3$
for cones. See Appendix B for details about the mesh convergence.
The low-
${Re}$
steady axisymmetric base flow
$\{\boldsymbol{u}_0,p_0\}$
is obtained by solving the axisymmetric steady version of the NS equations (2.1) using the Newton’s iteration method. The NS equations are completed with the following boundary conditions: uniform velocity field
$\boldsymbol{u}_0 = U_\infty \boldsymbol{e}_z$
at the inlet, stress-free
$p_0\boldsymbol{n} - {Re}^{-1}\, \nabla \boldsymbol{u}_0 \boldsymbol{\cdot} \boldsymbol{n} = \boldsymbol{0}$
at the outlet and far field (where
$\boldsymbol{n}$
is the unit normal vector), no-slip and no-penetration
$\boldsymbol{u}_0=\boldsymbol{0}$
on the body surface, and axisymmetry conditions
$u_{0r} = \partial _r u_{0z} = \partial _r p_0 = 0$
on the axis
$r=0$
. The generalised eigenvalue problem (2.4) for the onset of the primary instability is then solved using the Arnoldi algorithm, with standard boundary conditions:
$\boldsymbol{u}_1=\boldsymbol{0}$
at the inlet and on the body surface, stress-free
$p_1\boldsymbol{n} - {Re}^{-1}\, \nabla \boldsymbol{u}_1 \boldsymbol{\cdot} \boldsymbol{n} = \boldsymbol{0}$
at the outlet and far field, and
$m=1$
conditions
$\partial _r u_{1r} = \partial _r u_{1\theta } = u_{1z} = p_1 = 0$
on the axis. To avoid singularity on the axis, the equations are multiplied by
$r^2$
before taking the variational form. For the computation of the base flow, the finite elements formulation employs the higher-order Lagrange P3 and P2 elements for velocity and pressure, respectively. For the linear stability analysis, instead, P2 and P1 elements are used.
3. The symmetry-breaking bifurcation
3.1. Base flow
The axisymmetric low-
${Re}$
base flow for the four considered geometries is shown in figures 3 and 4 at
${Re} = {Re}_c$
, corresponding to the first onset of the symmetry-breaking instability, for 
and
$4$
. Only a brief characterisation is provided here, as the effect of the different placement of the corners resembles what was found by Chiarini et al. (Reference Chiarini, Quadrio and Auteri2022) for 2-D symmetric bodies.
Figure 3. Base flow for 
at
${Re} \approx {Re}_c$
: streamlines and contours of
$(a)$
azimuthal vorticity,
$(b)$
pressure. Ellipsoid,
${Re}=210$
; bicone,
${Re}=135$
; bullet,
${Re}=220$
; cone,
${Re}=160$
. Dashed line:
$u_{0z}=0$
.
Figure 4. Same as figure 3 for 
. Ellipsoid,
${Re}=1000$
; bicone,
${Re}=775$
; bullet,
${Re}=390$
; cone,
${Re}=260$
.
A shear layer with negative azimuthal vorticity
$\omega _{0\theta }$
separates from the rear part of the body and delimits the axisymmetric wake recirculation region. The vorticity is maximum in the vicinity of the body surface: for cones and bicones,
$\omega _{0\theta }$
is maximum close to the corners, while for ellipsoids and bullets, the region with intense vorticity is more spread. The pressure distribution changes accordingly (see the bottom right panels in figures 3 and 4). For bodies with a blunt TE (cones and bullets), the minimum of the pressure is placed at the TE corners, and the flow streamlines along the lateral side of the body face a favourable pressure gradient. For bodies with zero-thickness TE (ellipsoids and bicones), instead, the pressure is minimum where the cross-stream size of the body is maximum. In this case, the flow streamlines along the lateral side of the body face first a favourable pressure gradient and then an adverse one that promotes the flow separation. Notably, the adverse pressure gradient becomes milder as 
increases. For cones and bullets, the flow separation point is set by the geometry at the TE corner. The cross-stream dimension of the wake recirculation region is thus determined by the body width, and its extent only slightly decreases with 
. On the contrary, for ellipsoids and bicones, the flow separation is driven by the pressure distribution (see figures 3 and 4). In agreement with the milder pressure gradient, the separation point moves downstream for larger 
, yielding a thinner and shorter wake recirculation region (see also figure 19 in Appendix A).
3.2. Neutral curves
We now move to the results of the linear stability analysis. Figure 5 shows the neutral curves for the primary instability that, for all cases, consists in the regular (
$\lambda _i=0$
) symmetry-breaking bifurcation of an eigenmode of azimuthal wavenumber
$m=1$
. Table 1 compares the results of our computations with those of Meliga et al. (Reference Meliga, Chomaz and Sipp2009b
) and Bohorquez et al. (Reference Bohorquez, Sanmiguel-Rojas, Sevilla, Jiménez-González and Martínez-Bazán2011) for validation purposes. We find very good agreement with the results of Meliga et al. (Reference Meliga, Chomaz and Sipp2009b
), while some discrepancies are observed when comparing with Bohorquez et al. (Reference Bohorquez, Sanmiguel-Rojas, Sevilla, Jiménez-González and Martínez-Bazán2011) (we measure a relative difference in the value of
${Re}_c$
of approximately
$ 4\, \%$
for 
). We conjecture that this discrepancy is due to the different numerical method and to the different size of the computational domain.
Figure 5. Critical Reynolds number as a function of the aspect ratio.
Table 1. Comparison of the critical Reynolds number
${Re}_c$
with results from the literature, for some geometries.
Looking at figure 5, a first observation is that although the primary bifurcation is the same for the different bodies, the value of
${Re}_c$
shows large variability, with
${Re}_c \approx 150$
for the 
bicone, and
${Re}_c \approx 1400$
for the 
ellipsoid. This highlights that
$U_\infty$
and
$H$
are not the appropriate velocity and length scales for the characterisation and prediction of this bifurcation. Similarly to the flow past 2-D symmetric cylinders (Chiarini et al. Reference Chiarini, Quadrio and Auteri2022), an increase in the aspect ratio leads to a stabilisation of the base flow regardless of the body geometry; see in figure 5 that
${Re}_c$
increases monotonically with 
. This effect holds also for non-axisymmetric bodies, as shown e.g. by Zampogna & Boujo (Reference Zampogna and Boujo2023) and Chiarini & Boujo (Reference Chiarini and Boujo2025) in the context of 3-D rectangular prisms. Notably, the way
${Re}_c$
varies with 
depends on the shape of the TE of the body. For bodies with a blunt TE (cones and bullets), the critical Reynolds number almost flattens as 
increases. For bodies with zero-thickness TE (ellipsoids and bicones), instead, the critical Reynolds number increases faster than linearly, and very large values of
${Re}_c$
are observed already at intermediate 
. An explanation for the steep increase of
${Re}_c$
with 
for bodies with zero-thickness TE is provided in § 4.
For completeness, we show in figures 6 and 7 the axial velocity of the unstable mode for different geometries. Qualitatively, the structure of the mode does not change among the considered geometries, and agrees with the results of other authors; see e.g. figure 8 of Natarajan & Acrivos (Reference Natarajan and Acrivos1993) and figure 4 of Meliga et al. (Reference Meliga, Chomaz and Sipp2009a ). The perturbation field is confined in the wake, with the largest value found close to the TE of the body within the wake recirculation region. This is consistent with the appearance of a pair of streamwise counter-rotating vortices in the wake (see the discussion in § 1), and with the instability mechanism discussed in § 4. For bodies with zero-thickness TE, the region where the mode is intense extends farther downstream.
Figure 6. Eigenmode for 
at
${Re} \approx {Re}_c$
: streamwise velocity. Ellipsoid,
${Re}=210$
; bicone,
${Re}=135$
; bullet,
${Re}=220$
; cone,
${Re}=160$
. Dashed line:
$u_{0z}=0$
.
Figure 7. Same as figure 6 for 
. Ellipsoid,
${Re}=1000$
; bicone,
${Re}=775$
; bullet,
${Re}=390$
; cone,
${Re}=260$
.
Figure 8. (a,b) Structural sensitivity
$S$
and (c,d) sensitivity of the growth rate to base-flow modifications
$\nabla _{\boldsymbol{u}_0} \lambda _r$
, for 
at
${Re} \approx {Re}_c$
, for (a,c) ellipsoid,
${Re}=210$
, (b,d) bullet,
${Re}=220$
. Dashed line:
$u_{0z} = 0$
.
3.3. Sensitivities
To further characterise the instability, we consider in figures 8 and 9 the structural sensitivity (Giannetti & Luchini Reference Giannetti and Luchini2007) and the sensitivity to base-flow modifications (Marquet, Sipp & Jacquin Reference Marquet, Sipp and Jacquin2008; Meliga, Sipp & Chomaz Reference Meliga, Sipp and Chomaz2010). For brevity, here we show only ellipsoids and bullets, being representative of bodies with zero-thickness TE and with a blunt base.
Figure 9. Same as figure 8 for 
. Ellipsoid,
${Re}=1000$
; bullet,
${Re}=390$
.
The structural sensitivity
$S(\boldsymbol{x})$
is based on the interplay between the direct and adjoint modes, and identifies the region where a generic structural modification of the stability problem provides the largest drift of the leading eigenvalue. It is an upper bound for the eigenvalue variation
$|\delta \lambda |$
induced by a local body force actuation proportional to the signal of a velocity sensor located at the exact same station. The region of the flow where
$S$
is large identifies the wavemaker (Monkewitz, Huerre & Chomaz Reference Monkewitz, Huerre and Chomaz1993). For the primary instability, the structural sensitivity is defined as
\begin{equation} S(\boldsymbol{x}) = \frac { \| \hat {\boldsymbol{u}}_1^\dagger \|\, \| \hat {\boldsymbol{u}}_1 \| }{ \langle \hat {\boldsymbol{u}}_1^\dagger , \hat {\boldsymbol{u}}_1 \rangle }, \end{equation}
where
$\| \boldsymbol{\cdot} \|$
represents the
$R^2$
vector norm,
$\hat {\boldsymbol{u}}_1^\dagger$
is the adjoint mode, and
$ \langle \boldsymbol{u}_A,\boldsymbol{u}_B \rangle = \int _D \boldsymbol{u}_A^* \boldsymbol{\cdot} \boldsymbol{u}_B\, \text{d}\Omega$
is the inner product of
$L^2(D)$
, with the superscript
$*$
denoting the complex conjugate. For all cases,
$S(\boldsymbol{x})$
is substantially close to zero everywhere except close to the body where the product of the direct and adjoint modes is large. Large values are observed along the separation line
$u_{0z}=0$
and within the recirculating region, where the maximum is observed (Meliga et al. Reference Meliga, Chomaz and Sipp2009b
). The similar distribution of
$S(\boldsymbol{x})$
suggests that the different geometries and 
values considered share the same wavemaker and the same instability mechanism, which are spatially located within the recirculating region. Note that for bodies with zero-thickness TE, the spatial extent of the wavemaker decreases as 
increases, in agreement with the smaller recirculating region and with the base-flow stabilisation.
Next, the sensitivity to base-flow modifications quantifies the variation of the complex eigenpair
$(\lambda , \hat {\boldsymbol{u}}_1, \hat {p}_1)$
induced by a small variation of the base flow
$\delta \boldsymbol{u}_0$
. Specifically, the variation of the eigenvalue
$\delta \lambda$
is linked with
$\delta \boldsymbol{u}_0$
by the inner product
$\delta \lambda = \langle \nabla _{\boldsymbol{u}_0} \lambda , \delta \boldsymbol{u}_0 \rangle$
, where
\begin{equation} \nabla _{\boldsymbol{u}_0} \lambda = \frac { - ( \nabla \hat {\boldsymbol{u}}_1 )^H \boldsymbol{\cdot} \hat {\boldsymbol{u}}_1^\dagger + \nabla \hat {\boldsymbol{u}}_1^\dagger \boldsymbol{\cdot} \hat {\boldsymbol{u}}^* } { \langle \hat {\boldsymbol{u}}_1^\dagger , \hat {\boldsymbol{u}}_1 \rangle } \end{equation}
is indeed the sensitivity of
$\lambda$
to base-flow modifications; here, the superscript
$H$
indicates the trans-conjugate. The variation of the growth rate induced by
$\delta \boldsymbol{u}_0$
is thus expressed as
$\delta \lambda _r = ( \nabla _{\boldsymbol{u}_0} \lambda _r, \delta \boldsymbol{u}_0)$
, where the corresponding sensitivity is
$\nabla _{\boldsymbol{u}_0} \lambda _r = \text{Re}( \nabla _{\boldsymbol{u}_0} \lambda )$
. Unlike the structural sensitivity
$S$
, which is a scalar field, the sensitivity
$\nabla _{\boldsymbol{u}_0} \lambda _r$
to base-flow modification is a vector field: the field lines provide the local orientation of the sensitivity field, while the magnitude provides the intensity. As expected, far from the bodies,
$\nabla _{\boldsymbol{u}_0} \lambda _r$
decays to zero due to the spatial separation of the direct and adjoint modes. Large values are instead observed close to the
$u_{0z} = 0$
line and within the recirculating region. For bodies with zero-thickness TE,
$\nabla _{\boldsymbol{u}_0} \lambda _r$
is large within the entire recirculating region, while for bodies with a blunt base, the sensitivity is maximum at the downstream end of the recirculating region, close to the base and in correspondence with the corners (this resembles what is observed for 2-D rectangular cylinders; see Chiarini, Quadrio & Auteri Reference Chiarini, Quadrio and Auteri2021). For all cases, an increase of the backflow within the recirculating region (
$ \delta u_{0z}\lt 0$
) largely destabilises the flow. Again, for bodies with zero-thickness TE, the spatial extent where the sensitivity is large decreases with 
.
3.4. Effect of the Reynolds number
We now investigate the effect of a small variation of
${Re}$
on the growth rate
$\lambda _r$
. We start by looking for an expression for
$\partial \lambda _r/\partial {Re}$
. We consider a small departure from criticality such that
\begin{equation} \frac {1}{{Re}_c} - \frac {1}{{Re}} = \varepsilon, \end{equation}
where
$|\varepsilon | \ll 1$
(
$\varepsilon \lt 0$
for
${Re}\lt {Re}_c$
, and
$\varepsilon \gt 0$
for
${Re}\gt {Re}_c$
). A small-amplitude change in
${Re}$
induces small changes in the base flow, in the eigenmode and in the eigenvalues, i.e.
\begin{equation} \{ \boldsymbol{u}_0, p_0 \} \rightarrow \{ \boldsymbol{u}_0, p_0 \} + \{\delta \boldsymbol{u}_0, \delta p_0 \}, \quad \{\hat {\boldsymbol{u}}_1,\hat {p}_1\} \rightarrow \{\hat {\boldsymbol{u}}_1,\hat {p}_1\} + \{\delta \hat {\boldsymbol{u}}_1, \delta \hat {p}_1 \}, \quad \lambda \rightarrow \lambda + \delta \lambda . \end{equation}
We now set
${Re}={Re}_c$
as reference, and inject these changes into the steady nonlinear NS equations and into the linearised NS equations. By keeping the first-order terms only, we then obtain the equations for the base-flow modification
$\{\delta \boldsymbol{u}_0,\delta p_0\}$
and for the eigenmode modification
$\{\delta \hat {\boldsymbol{u}}_1, \delta \hat {p}_1\}$
, i.e.
\begin{equation} \begin{pmatrix} \mathcal{L}_0\{\boldsymbol{u}_0,{Re}\}(\boldsymbol{\cdot} ) & - \nabla _0(\boldsymbol{\cdot} ) \\[1.5pt] \nabla _0 \boldsymbol{\cdot} (\boldsymbol{\cdot} ) & 0 \end{pmatrix} \begin{pmatrix} \delta \boldsymbol{u}_0 \\[1.5pt] \delta p_0 \end{pmatrix} = \begin{pmatrix} - \varepsilon\, \nabla _0^2 \boldsymbol{u}_0 \\[1.5pt] 0 \end{pmatrix} \end{equation}
and
\begin{equation} \lambda \begin{pmatrix} \mathcal{I} & \boldsymbol{0} \\ \boldsymbol{0}^{\rm T} & 0 \end{pmatrix} \begin{pmatrix} \delta \hat {\boldsymbol{u}}_1 \\ \delta \hat {p}_1 \end{pmatrix} + \begin{pmatrix} \mathcal{L}_m\{\boldsymbol{u}_0,{Re}\}(\boldsymbol{\cdot} ) & - \nabla _m(\boldsymbol{\cdot} ) \\ \nabla _m \boldsymbol{\cdot} ( \boldsymbol{\cdot} ) & 0 \end{pmatrix} \begin{pmatrix} \delta \hat {\boldsymbol{u}}_1 \\ \delta \hat {p}_1 \end{pmatrix} = \begin{pmatrix} \mathcal{F}_m(\delta \boldsymbol{u}_0, \hat {\boldsymbol{u}}_1) \\ 0 \end{pmatrix}, \end{equation}
accompanied by homogeneous Dirichlet boundary conditions at the inlet and at the surface of the body, and stress-free conditions at the outlet and at the far field; here,
\begin{equation} \mathcal{F}_m(\delta \boldsymbol{u}_0,\hat {\boldsymbol{u}}_1) = - \mathcal{C}_m( \hat {\boldsymbol{u}}_1, \delta \boldsymbol{u}_0 ) - \,\varepsilon \nabla ^2_m \hat {\boldsymbol{u}}_1 - \delta \lambda \hat {\boldsymbol{u}}_1. \end{equation}
Equation (3.5) is divided by
$\varepsilon$
and solved for
$\boldsymbol{u}_0^\varepsilon = \delta \boldsymbol{u}_0/\varepsilon$
. Equation (3.6) is projected on the adjoint mode
$\{\hat {\boldsymbol{u}}_1^\dagger ,\hat {p}_1^\dagger \}$
to eliminate
$\delta \hat {\boldsymbol{u}}_1$
and obtain an expression for
$\delta \lambda$
, i.e.
\begin{equation} \delta \lambda = -\varepsilon \frac { \langle \hat {\boldsymbol{u}}_1^\dagger , \mathcal{C}_m( \hat {\boldsymbol{u}}_1, \boldsymbol{u}_0^\varepsilon ) + \nabla ^2_m \hat {\boldsymbol{u}}_1 \rangle } { \langle \hat {\boldsymbol{u}}_1^\dagger , \hat {\boldsymbol{u}}_1 \rangle }. \end{equation}
Interestingly, this expression matches exactly the one for the linear coefficient of the amplitude equation obtained with a standard weakly nonlinear stability analysis (Sipp & Lebedev Reference Sipp and Lebedev2007; Zampogna & Boujo Reference Zampogna and Boujo2023). At this point, we can write the
$\partial \lambda _r/\partial {Re}$
derivative as
\begin{equation} \frac {\partial \lambda _r}{\partial {Re}} = \frac {\partial \lambda _r}{\partial \varepsilon } \frac {\partial \varepsilon }{\partial \rm Re} = \frac {1}{{Re}^2} \frac {\partial \lambda _r}{\partial \varepsilon } = - \frac {1}{{Re}^2}\, \text{Re} \left ( \frac { \langle \hat {\boldsymbol{u}}_1^\dagger , \mathcal{C}_m( \hat {\boldsymbol{u}}_1, \boldsymbol{u}_0^\varepsilon ) + \nabla ^2_m \hat {\boldsymbol{u}}_1 \rangle } { \langle \hat {\boldsymbol{u}}_1^\dagger , \hat {\boldsymbol{u}}_1 \rangle } \right ). \end{equation}
A small change in the Reynolds number thus leads to a modification of the growth rate that depends on the effect of the base-flow modification on the eigenmode, i.e.
$\langle \hat {\boldsymbol{u}}_1^\dagger , \mathcal{C}_m( \hat {\boldsymbol{u}}_1, \boldsymbol{u}_0^\varepsilon ) \rangle$
, and on the direct effect of
${Re}$
on the viscous term of the eigenmode, i.e.
$ \langle \hat {\boldsymbol{u}}_1^\dagger , \nabla ^2_m \hat {\boldsymbol{u}}_1 \rangle$
. The spatial distribution of the integrand in the inner product of (3.9) provides information regarding the region that contributes most to the increase of the growth rate when
${Re}$
increases (recall that we set
${Re}={Re}_c$
as reference such that here
$\partial \lambda _r/\partial {Re}\gt 0$
).
Figure 10 shows the results for the ellipsoid and the bullet with 
. Figures 10(a,b) show the modification of the base-flow azimuthal vorticity
$\omega _{0\theta }^\varepsilon$
(to be compared with figure 3). For all considered geometries, the
$\varepsilon$
Reynolds increase leads to a more negative vorticity in the layer along the separating streamline, and to a more positive vorticity in the base region and in a tiny elongated region along the
$u_{0z} = 0$
line. Interestingly, the latter region closely matches the region where
$\partial \omega _{0\theta }/\partial r$
first changes sign (see the discussion in § 4). Figures 10(c,d) show the integrand of the inner product in (3.9), and reveal that the destabilisation due to the increase of
${Re}$
mostly comes from the upstream region of the recirculating region, which closely matches the wavemaker identified by the structural sensitivity. Moreover, a close inspection of the terms in (3.9) shows that almost the complete destabilisation is due to the convective-like term resulting from the base-flow modification rather than due to the direct effect of
${Re}$
on the eigenmode, as
$|\hat {\boldsymbol{u}}_1^{\dagger ,*} \boldsymbol{\cdot} \mathcal{C}_m( \hat {\boldsymbol{u}}_1, \boldsymbol{u}_0^\varepsilon ) | \gg |\hat {\boldsymbol{u}}_1^{\dagger ,*} \boldsymbol{\cdot} \nabla ^2_m \hat {\boldsymbol{u}}_1|$
at all positions (not shown).
Figure 10. Effect of a small increase in
${Re}$
on the base flow and on the growth rate for 
at
${Re} \approx {Re}_c$
: (a,b) spatial distribution of
$\omega _{\theta ,0}^{\varepsilon }$
, (c,d) spatial distribution of the integrand of
$\delta \lambda _r$
(see text and (3.9)), for (a,c) ellipsoid,
${Re}=210$
, (b,d) bullet,
${Re}=220$
. Dashed line:
$u_{0z} = 0$
.
4. Azimuthal vorticity and mechanism of the primary bifurcation
In this section, we show that the instability mechanism proposed by MM for oblate spheroidal objects with a free-slip surface extends to axisymmetric bodies with a no-slip surface.
4.1. Maximum surface vorticity
As discussed in § 1, the MM arguments are based on the idea that the bifurcation is driven by the vorticity generated at the body surface and transported into the wake. Therefore, we start by assessing the relation between the base-flow surface vorticity and the onset of the bifurcation. To simplify the notation, in this and in the following subsection (§ 4.2), we drop the
$0$
subscript. Introducing
$\omega _{\theta ,max}({Re})$
as the maximum vorticity measured on the surface of a body of given geometry and aspect ratio, figure 11 shows
$\omega _{\theta ,c}({Re})$
, i.e. the maximum surface vorticity
$\omega _{\theta ,max}$
at the critical Reynolds number
${Re}={Re}_c$
. We report values for the ellipsoids and bullets, which have a smooth geometry. For cones and bicones, we note that the maximum surface vorticity diverges at the sharp edges, which leads to
$\omega _{\theta ,max}$
growing unboundedly when refining the numerical mesh, while the critical Reynolds number is well converged. It is therefore unlikely that a simple relationship between
$\omega _{\theta ,max}$
and
${Re}_c$
holds for bodies with a sharp geometry. See Appendix C for more details.
Figure 11. Maximum surface azimuthal vorticity at
${Re}={Re}_c$
.
Similarly to what was observed by MM for free-slip surfaces, figure 11 shows that
$\omega _{\theta ,c}$
has a linear dependence on
${Re}$
, and that the results for the two considered geometries collapse rather well onto the same curve:
\begin{equation} \omega _{\theta ,c}({Re}) \approx a + b{Re}, \end{equation}
where
$a$
and
$b$
are two constants. The amount of vorticity that has to be produced at the body surface to promote the flow instability increases linearly with the Reynolds number, in a way that does not depend on the geometry of the body. This provides a simple criterion for predicting the onset of the instability, which requires only the knowledge of the maximum azimuthal vorticity at the surface (which is, however, not always easily available, as discussed later, in § 5): at a given Reynolds number, the flow is stable when
$\omega _{\theta ,max }({Re}) \lt \omega _{\theta ,c}({Re})$
, and unstable when
$\omega _{\theta ,max }({Re}) \gt \omega _{\theta ,c}({Re})$
. The linear fit of our data with (4.1) gives
$a = 14$
and
$b = 4.4 \times 10^{-2}$
. These values differ from the constants
$a = 12.5$
and
$b=4.3 \times 10^{-3}$
found by MM for ellipsoids with free-slip surface. This is expected, as the maximum surface vorticity at criticality depends on the slip at the body surface; see e.g. figure 4 of Legendre, Lauga & Magnaudet (Reference Legendre, Lauga and Magnaudet2009) for the flow past a circular cylinder at different Knudsen numbers. Notably, the different values of
$b$
show that the critical vorticity at the body surface exhibits a faster variation with
${Re}$
for axisymmetric bodies with a no-slip surface than for those with a free-slip surface.
4.2. Vorticity gradient
We now focus on the near-wake region, and investigate the correlation between the flow bifurcation and the appearance of points where
$\partial \omega _\theta /\partial r = 0$
. Figure 12 shows the vorticity distribution in the near wake for the 
ellipsoid (with the same representation as in figure 1). Due to the no-slip boundary condition, a thin boundary layer arises on the rear side of the body, where the azimuthal vorticity is positive,
$\omega _\theta \gt 0$
, i.e. has the opposite sign as in the rest of the flow (the black thick line in figure 12 denotes
$\omega _\theta =0$
). Outside this boundary layer, the vorticity distribution closely resembles that of free-slip bodies: isocontours of
$\omega _\theta$
are essentially parallel to the solid surface in the vicinity of the body, and farther downstream, they align with the symmetry axis. In the transition region in between, the turning of the
$\omega _\theta$
isocontours becomes sharper and sharper as
${Re}$
gets larger. Two values of the Reynolds number are considered in figure 12, one below (figure 12
$a$
) and one above (figure 12
$b$
) the critical Reynolds number
${Re}_c$
. One can immediately notice that
$\partial \omega _\theta /\partial r \lt 0$
(blue) everywhere in the transition region when
${Re} \lt {Re}_c$
, while a region where
$\partial \omega _\theta /\partial r \geqslant 0$
(red) appears when
${Re} \gt {Re}_c$
, in agreement with the MM arguments. Figures 13 and 14 are as figure 12
$(b)$
, for different geometries and aspect ratios. Notably, the location where the
$\partial \omega _\theta /\partial r \geqslant 0$
region appears depends on the shape of the TE. For a blunt TE, it is placed along the
$\omega _\theta = 0$
isocontour delimiting the base boundary layer due to the no-slip boundary condition. For bodies with zero-thickness TE, instead, it is located farther downstream within the wake recirculation region, and close to the contour of zero streamwise velocity
$u_z=0$
.
Figure 12. Near-wake distribution of the azimuthal vorticity around an ellipsoid with no-slip surface and 
, for
$(a)$
${Re}=250 \lt \widetilde {{Re}}$
,
$(b)$
${Re}=290 \gt \widetilde {{Re}}$
, where
$\widetilde {{Re}}$
is the Reynolds number corresponding to the first appearance of a region with
$\partial \omega _\theta /\partial r \geqslant 0$
with
$\omega _\theta \lt 0$
. Thin black lines are isocontours of the azimuthal vorticity
$\omega _\theta$
. The thick black line delimits the ‘boundary layer’ where
$\omega _\theta \gt 0$
. Coloured contours show
$\partial \omega _\theta /\partial r$
. The grey dashed line shows
$u_z = 0$
, and the red dashed line shows
$\partial \omega _\theta /\partial r = 0$
. Blue/red circles indicate the negative/positive maxima of
$\partial \omega _\theta /\partial r$
in the near-wake region.
Figure 13. As figure 12 for the four considered geometries with 
. The Reynolds number is slightly above
$\widetilde {{Re}}$
, i.e. after the first appearance of a region with
$\partial \omega _\theta /\partial r \geqslant 0$
where
$\omega _\theta \lt 0$
. Here, (a) ellipsoid at
${Re}=270$
, (b) bicone at
${Re}=295$
, (c) bullet at
${Re}=120$
, (d) cone at
${Re}=80$
.
Figure 14. As figure 13 for 
, for (a) ellipsoid at
${Re}=1250$
, (b) bicone at
${Re}=950$
, (c) bullet at
${Re}=175$
, (d) cone at
${Re}=150$
.
It is worth noting that the vorticity distribution for bodies with zero-thickness TE explains the steep increase of
${Re}_c$
with 
observed in figure 5. Due to the specific shapes of these bodies, the
$\omega _\theta$
isocontours close to the body surface become more aligned with the symmetry axis
$r=0$
as 
increases. For these bodies, therefore, the condition
$\partial \omega _\theta / \partial r = 0$
(
$\omega _\theta$
isocontours being locally perpendicular to the symmetry axis) for the instability to occur requires more vorticity to be produced at the body surface, i.e. a larger Reynolds number.
Figure 15 quantifies the correlation between the onset of the instability and the first appearance of points where
$\partial \omega _\theta /\partial r = 0$
. We compare the critical Reynolds number
${Re}_c$
obtained from our stability analysis with the Reynolds number
$\widetilde {{Re}}_c$
corresponding to the first appearance of a point where
$\partial \omega _\theta /\partial r = 0$
in the near wake (we consider only the near wake outside the boundary layer that develops along the TE). Although the collapse with
${Re}_c$
is not perfect,
$\widetilde {{Re}}_c$
captures the trend of the critical Reynolds number rather nicely for all the considered geometries, increasing sharply with 
for bodies with a zero-thickness TE, and much more slowly for bodies with a blunt TE.
Figure 15. Critical Reynolds number
${Re}_c$
(lines) and lowest Reynolds number
$\widetilde {{Re}}$
(open symbols) such that
$\partial \omega _\theta /\partial r \geqslant 0$
where
$\omega _\theta \lt 0$
in the near wake.
To summarise, our results show that the instability mechanism proposed by MM in the context of ellipsoids with free-slip surfaces also extends to no-slip surfaces, and describes fairly well the primary symmetry-breaking bifurcation of the flow past axisymmetric bodies with different geometries and aspect ratios. This clearly hints to the fact that the vorticity generated at the body surface is at the root of the instability mechanism, while the way it is produced, i.e. the nature of the surface, does not play a major role.
4.3. Free-slip versus no-slip
It is worth spending few words on the fact that, as shown in the previous subsections, axisymmetric wakes past free-slip and no-slip bodies share the same instability mechanism. It is well known that the mechanism of vorticity generation on a surface changes with the boundary condition; see e.g. Truesdell (Reference Truesdell1954), Lighthill (Reference Lighthill1963), Morton (Reference Morton1984), Leal (Reference Leal1989), Wu & Wu (Reference Wu and Wu1993), Wu (Reference Wu1995), Lundgren & Koumoutsakos (Reference Lundgren and Koumoutsakos1999), Brøns et al. (Reference Brøns, Thompson, Leweke and Hourigan2014) and Terrington et al. (Reference Terrington, Hourigan and Thompson2020). This is clear when looking at the ‘boundary vorticity flux’ (normal diffusion flux of vorticity) at the surface,
$\boldsymbol{\sigma } = \nu \boldsymbol{n} \boldsymbol{\cdot} \nabla \boldsymbol{\omega }$
, which is a measure of the vorticity creation at the wall. For stationary rigid walls, it reads (see Wu & Wu Reference Wu and Wu1993)
\begin{equation} \boldsymbol{\sigma } = \underbrace { \boldsymbol{n} \times \left ( \frac {\nabla p}{\rho } \right ) }_{\boldsymbol{\sigma }_p} \underbrace {{} - \left ( \boldsymbol{n} \times \nabla \right ) \boldsymbol{\cdot} \left ( \boldsymbol{\tau } \boldsymbol{n} \right ) }_{\boldsymbol{\sigma }_{\tau } }, \end{equation}
where
$\boldsymbol{\tau } = \nu \boldsymbol{\omega } \times \boldsymbol{n}$
is the tangential shear stress. Here,
$\boldsymbol{\sigma }$
depends on the non-uniform distribution of the normal (
$\boldsymbol{\sigma }_p$
) and shear (
$\boldsymbol{\sigma }_\tau$
) stresses. In the limit case of a 2-D flow
$\boldsymbol{u}(x,y)=(u_x,u_y,0)$
near a stationary, no-slip flat wall of normal
$\boldsymbol{n}=\boldsymbol{e}_y$
, the boundary vorticity flux reduces to
$\boldsymbol{\sigma } = \boldsymbol{\sigma }_p = -(1/ \rho ) (\partial p /\partial x) \boldsymbol{e}_z$
, while
$\boldsymbol{\sigma }_\tau =\boldsymbol{0}$
(Lighthill Reference Lighthill1963). By contrast, on a free-slip surface, the vorticity appears as a consequence of the continuity of the tangent stresses and is non-null only in the case of curved surfaces; in the limit case of a 2-D stationary surface,
$\omega = 2 U \kappa$
, where
$U$
is the tangential velocity and
$\kappa$
is the local curvature; accordingly, the boundary vorticity flux is non-null only for curved surfaces (see Wu Reference Wu1995).
The fact that axisymmetric bodies with no-slip and free-slip surfaces exhibit the same instability mechanism means that although the vorticity is at the root of the instability, the way it is produced at the wall does not play a dominant role in the overall triggering mechanism (see also the discussion in MM). In other words, the flow becomes unstable once the amount of vorticity brought into the near wake is large enough, irrespective of the production mechanism. Clearly, the different vorticity creation mechanism has an impact on the critical Reynolds number and on the actual onset of the instability, as it influences the size of the wake recirculation region and the amount of vorticity brought into the wake. The different vorticity creation mechanism, indeed, leads to some substantial differences between no-slip and free-slip surfaces in terms of boundary vorticity dynamics (Moore Reference Moore1963; Wu Reference Wu1995). For large Reynolds numbers, dimensional analysis shows that: (i) in the no-slip case, a boundary vorticity flux
$\sigma \sim O(1)$
is generated, which implies that the surface vorticity increases as
$\omega \sim {Re}^{1/2}$
; (ii) in the free-slip case, the boundary vorticity flux decreases as
$\sigma \sim {Re}^{-1/2}$
, and the amount of surface vorticity is independent of the Reynolds number,
$\omega \sim O(1)$
.
4.4. On the amplification mechanism
We now investigate the contribution of various physical mechanisms to the exponential growth of the mode in the linear regime (
${Re} \approx {Re}_c$
). We start by looking at the energy equation for the perturbation
$\boldsymbol{u}_1$
(see Lanzerstorfer & Kuhlmann Reference Lanzerstorfer and Kuhlmann2012). We sum (2.4) multiplied by
$\hat {\boldsymbol{u}}_1^*$
with the complex conjugate of the same equation multiplied by
$\hat {\boldsymbol{u}}_1$
, and after some manipulation and using the incompressibility constraint, we obtain the equation for
$\hat {\boldsymbol{u}}_1 \boldsymbol{\cdot} \hat {\boldsymbol{u}}_1^*$
, i.e.
\begin{align} 2 \lambda _r \hat {\boldsymbol{u}}_1 \boldsymbol{\cdot} \hat {\boldsymbol{u}}_1^* ={}& \underbrace { - \lbrace( \hat {\boldsymbol{u}}_1 \otimes \hat {\boldsymbol{u}}_1^* + \hat {\boldsymbol{u}}_1^* \otimes \hat {\boldsymbol{u}}_1 ) : \nabla _0 \boldsymbol{u}_0 \rbrace}_{ \mathcal{P} } \underbrace {{} - \frac {1}{{Re}} \nabla _0^2 ( \hat {\boldsymbol{u}}_1 \boldsymbol{\cdot} \hat {\boldsymbol{u}}_1^* ) }_{ \mathcal{D} } \nonumber\\& \underbrace {{} - \nabla _0 \boldsymbol{\cdot} ( \hat {p}_1 \hat {\boldsymbol{u}}_1^* + \hat {p}_1^* \hat {\boldsymbol{u}}_1 ) }_{ \mathcal{T}_p } \underbrace { - \nabla _0 \boldsymbol{\cdot} ( \boldsymbol{u}_0 \hat {\boldsymbol{u}}_1 \boldsymbol{\cdot} \hat {\boldsymbol{u}}_1^* ) }_{ \mathcal{A} } - \underbrace { \frac {2}{{Re}} \left ( \nabla _m \hat {\boldsymbol{u}}_1 : \nabla _m^* \hat {\boldsymbol{u}}_1^* \right ) }_{ \mathcal{E} }. \end{align}
On the right-hand side, in order, we find the production term
$\mathcal{P}$
that accounts for the exchange of energy between the base flow and the perturbation, the viscous diffusion term
$\mathcal{D}$
, the pressure transport term
$\mathcal{T}_p$
, the advection associated with the base flow
$\mathcal{A}$
, and the viscous dissipation
$\mathcal{E}$
. All the terms contribute to the total budget, but differently depending on the flow region.
Figure 16 considers the energy equation for the ellipsoid and the bullet with 
at
${Re}={Re}_c$
. In this case,
$\lambda _r = 0$
and the Reynolds–Orr equation reduces to
$\int _D \mathcal{P}\, \text{d}\Omega = \int _D \epsilon\, \text{d}\Omega$
, with
$\int _D \mathcal{A} + \mathcal{T}_p + \mathcal{D} = 0$
. We focus on the spatial distribution of
$\mathcal{P}$
and
$\mathcal{A}$
to locate the flow regions where the mode amplification is driven by the production and/or the advection. We observe that
$\mathcal{P}\gt 0$
downstream of the bodies, with the maximum placed along the separating streamline and within the recirculating region, where the base-flow velocity gradient is indeed maximum. At the base, the small
$\mathcal{P}\lt 0$
region indicates that energy is moved from the perturbation to the base flow there. For bodies with a blunt base, we observe a more intense production activity, in agreement with the stronger base-flow velocity gradients. The advection term
$\mathcal{A}$
, instead, behaves differently depending on the direction of
$\boldsymbol{u}_0$
. It is a sink (
$\mathcal{A}\lt 0$
) and tends to stabilise the flow along the separating shear layer where
$u_{0z}\gt 0$
, while it is a source (
$\mathcal{A}\gt 0$
) and tends to destabilise the flow within the recirculating region where
$u_{0z}\lt 0$
. In the recirculating region, the base flow transports the perturbations backwards in the near wake so that they can grow more compared to the situation where they are advected downstream. Note that this agrees with § 3.4, where the effect of a small increase in
${Re}$
on
$\lambda _r$
is discussed (see in figure 10 the large positive
$\partial \lambda _r/\partial {Re}\gt 0$
region in the core of the recirculating region).
Figure 16. Energy budget for 
at
${Re} \approx {Re}_c$
: (a,b) production term
$\mathcal{P}$
, (c,d) advection term
$\mathcal{A}$
, for (a,c) ellipsoid,
${Re}=210$
, (b,d) bullet,
${Re}=220$
. Dashed line:
$u_{0z} = 0$
.
Similar conclusions can be drawn by using the concept of endogeneity introduced by Marquet & Lesshafft (Reference Marquet and Lesshafft2015). The endogeneity
$E(\boldsymbol{x})$
characterises the contributions of localised flow regions to the global dynamics of the eigenmode, and allows us to separate the contributions from individual mechanisms such as production, base-flow advection, pressure forces and viscous diffusion. The endogeneity is defined as
$E(\boldsymbol{x}) = \hat {\boldsymbol{u}}_1^{\dagger *} \boldsymbol{\cdot} ( - \mathcal{L}_m \hat {\boldsymbol{u}}_1 - \nabla _m \hat {p}_1 )$
, and has the property that its integral equals the eigenvalue, i.e
$\int _D E(\boldsymbol{x})\, \text{d}\Omega = \lambda$
. We now replace the definition of
$\mathcal{L}_m$
in the above relation to obtain
\begin{equation} E(\boldsymbol{x}) = - \hat {\boldsymbol{u}}_1^{\dagger *} \boldsymbol{\cdot} ( (\boldsymbol{u}_0 \boldsymbol{\cdot} \nabla _m ) \hat {\boldsymbol{u}}_1 ) - \hat {\boldsymbol{u}}_1^{\dagger *} \boldsymbol{\cdot} ( (\hat {\boldsymbol{u}}_1 \boldsymbol{\cdot} \nabla _m ) \boldsymbol{u}_0 ) - \hat {\boldsymbol{u}}_1^{\dagger *} \boldsymbol{\cdot} \nabla \hat {p}_1 + \frac {1}{{Re}}\, \hat {\boldsymbol{u}}_1^{\dagger *} \boldsymbol{\cdot} \nabla _m^2 \hat {\boldsymbol{u}}_1, \end{equation}
where, in order, the different terms account for the contribution of the base-flow advection, the production due to the base-flow velocity gradients, pressure forces and viscous diffusion. The real part of the endogeneity and of (4.4) isolates the contributions to the growth rate
$\lambda _r$
, while the imaginary part isolates the contributions to the eigenfrequency
$\lambda _i$
. Figure 17 shows the production and advection contribution to
$\text{Re} (E(\boldsymbol{x}))$
for the ellipsoid and the bullet with 
at
${Re}={Re}_c$
. The endogeneity is concentrated close to the body around the shear layers and within the recirculating region. The production has a stabilising/destabilising effect depending on the region, with the positive contribution to
$\text{Re} (E(\boldsymbol{x}))$
being maximum at the centre of the recirculating region. The advection, instead, has a strongly destabilising effect within the recirculating region where
$u_{0z}\lt 0$
, but a stabilising effect in the first portion of the shear layers separating from the bodies where
$u_{0z}\gt 0$
. This partially agrees with the energy budget analysis.
Figure 17. Spatial distribution of the (a,b) production and (c,d) advection contributions to the endogeneity for 
at
${Re} \approx {Re}_c$
, for (a,c) ellipsoid,
${Re}=210$
, (b,d) bullet,
${Re}=220$
. Dashed line:
$u_{0z} = 0$
.
5. A new scaling
In the previous section, we have shown that the base-flow azimuthal vorticity drives the primary instability of the flow past 3-D axisymmetric bodies. We now use this information to define a new Reynolds number that is more suitable for predicting the onset of the instability than the standard Reynolds number based on
$U_\infty$
and
$H$
. This is done in the same way as by Chiarini et al. (Reference Chiarini, Quadrio and Auteri2022), who proposed a new scaling for the description of the primary instability of the flow past 2-D symmetric bluff bodies. They used a measure of the spatial extent of the separation bubble as a length scale, and the largest reverse flow speed within it as the velocity scale, and found that the ensuing Reynolds number evaluated at criticality only marginally changed with the geometry and aspect ratio of the body. This choice of scales was inspired by their link to the local amplification of the unstable wave packets and to the extension of the absolute instability region (Hammond & Redekopp Reference Hammond and Redekopp1997; Chomaz Reference Chomaz2005). Here, we introduce a new scaling that can be used to estimate whether the low-
${Re}$
steady flow past axisymmetric bodies is unstable, without the need for a computationally expensive stability analysis. This new scaling is based on the inspection of the base flow only (here we again drop the
$0$
subscript), and relies on quantities that are related to the physics of the problem and are easily accessible in both experiments and numerical simulations.
The onset of the primary symmetry-breaking bifurcation depends on the thickness of the transition region where the instability mechanism takes place (see § 4). A thinner transition region, indeed, means that the
$\omega _\theta$
isocontours exhibit sharper turnings, and points with
$\partial \omega _\theta /\partial r = 0$
are more likely to arise. The thickness of this region depends on (i) the maximum vorticity
$\omega _{\theta ,max }$
generated at the body surface (see figure 11), and (ii) the radial extent of the wake recirculation region (see MM). Thus one may be tempted to use
$\omega _{\theta ,max }$
and a measure of
$r_0$
to define a new Reynolds number. However, the measure of
$\omega _{\theta ,max }$
at the body surface is not always available. First, vorticity diverges at a geometric singularity such as a sharp leading edge, as observed for cones and bicones. Second, in experiments, the flow region close to the body surface is not easily accessible. A different estimation of the surface vorticity is thus needed. Here, we propose
$\widetilde {\omega }_\theta = \max _r (\omega _\theta (r,z_{max}))$
, where
$z_{max}$
is the streamwise location of maximum reverse flow
$U_{rev}$
within the wake recirculation region. We use
$\widetilde {\omega }_\theta$
as a vorticity-related scale, and the radius of the
$u_z=0$
isocontour at the same
$z_{max}$
streamwise location as a measure of
$r_0$
, and define the new Reynolds number as
$\widehat {{Re}} = \widetilde {\omega }_\theta r_0^2/\nu$
. Figure 18
$(a)$
shows that at criticality, this new Reynolds number
$\widehat {{Re}}$
collapses approximately to the same value for all the considered bodies. The collapse is not perfect, yet the relative variation is one order of magnitude smaller than that observed with the standard Reynolds number based on
$U_\infty$
and
$H$
. To be quantitative, we measure e.g.
$\widehat {{Re}}_c \approx 133$
for the 
cone, and
$\widehat {{Re}}_c \approx 147$
for the 
ellipsoid, to be compared with
${Re}_c \approx 167$
and
${Re}_c \approx 1360$
, respectively. More generally, the new proposed Reynolds number provides a simple and effective global criterion for the prediction of the instability onset: considering a variation in
$\widehat {{Re}}_c$
of one standard deviation around its mean value (both computed with all the cases considered in this study and with weights such that each of the four geometries contributes equally), we obtain
$125 \pm 20$
and can therefore conclude that the flow is likely to be stable when
$\widehat {{Re}} \lesssim 105$
and likely to be unstable when
$\widehat {{Re}} \gtrsim 145$
. For comparison, the original critical Reynolds number has a much larger relative spread,
${Re}_c = 480 \pm 340$
.
Figure 18. Critical Reynolds number
${Re}_c$
(dashed lines) and newly defined Reynolds numbers (solid lines):
$(a)$
$\widehat {{Re}}_c = \widetilde {\omega }_\theta r_0^2/ \nu$
,
$(b)$
$\overline {{Re}}_c = U_{rev} r_0/\nu$
.
It is worth noting that the new scaling works well also when other vorticity-related quantities are used instead of
$\widetilde {\omega }_\theta$
, provided that they are a measure of the vorticity generated at the body surface. As an example, in figure 18
$(b)$
, we keep
$r_0$
as the length scale and use
$\widehat {\omega }_\theta = U_{rev}/r_0$
as a vorticity scale, which may indeed be seen as a rough estimate of
$\widetilde {\omega }_\theta$
. Again, at criticality, the new proposed Reynolds number
$\overline {{Re}}_c = \widehat {\omega }_\theta r_0^2/\nu = U_{rev} r_0/\nu$
collapses rather well for all the considered aspect ratios and geometries to
$ 34 \pm 6$
. The relevance of
$U_{rev}$
on the physics of the instability agrees with the map of the sensitivity to base-flow modifications shown in figures 8 and 9.
6. Conclusion
In this study, we have investigated the primary symmetry-breaking bifurcation of the flow past axisymmetric bodies with a no-slip surface. We have computed the neutral curves 
and investigated the instability mechanism using various quantities: structural sensitivity, sensitivity to a base-flow modification or to an increase in
${Re}$
, dynamics of the perturbation kinetic energy, and endogeneity. We have also proposed a new scaling that is suitable for predicting the onset of the bifurcation. The generality of our conclusions has been assessed by considering bodies with different geometries (ellipsoids, bullets, cones and bicones) and different aspect ratios (
).
The instability is driven by the azimuthal vorticity generated at the body surface. We have shown that the mechanism proposed by MM in the context of free-slip spheroidal bubbles well describes the onset of the instability for bodies with a no-slip surface too. In both cases, indeed, the instability arises in a thin region in the near wake, and its onset is strongly related to the occurrence of points where
$\partial \omega _\theta /\partial r = 0$
, i.e. points where isocontours of
$\omega _\theta$
align with the radial direction and are perpendicular to the symmetry axis. At these points, the streamwise vorticity gradient
$\partial \omega _\theta /\partial z$
largely increases with the Reynolds number, favouring the instability.
Having characterised the physical mechanism, we have then proposed a new scaling for the prediction of the instability, in the spirit of the study of Chiarini et al. (Reference Chiarini, Quadrio and Auteri2022) in the context of 2-D bluff bodies. The new scaling is based on measures of the near-wake azimuthal vorticity,
$\widetilde {\omega }_\theta$
, and the radial extent of the wake recirculation region,
$r_0$
. When computed at criticality, the ensuing Reynolds number
$\widehat {{Re}} = \widetilde {\omega }_\theta r_0^2/\nu$
has been shown to collapse approximately to the same value across all geometries and aspect ratios. This observation can be used to readily estimate whether the steady base flow past axisymmetric bodies is unstable, avoiding a computationally expensive stability analysis, and relying only on an inspection of the base flow.
While the present study focuses on the first instability, it is worth concluding with a word about turbulent wakes. In the turbulent regime, while the mean (time-averaged) flow past an axisymmetric bluff body is itself axisymmetric, the most probable instantaneous flow state is a static
$m=1$
symmetry-breaking deflection (Grandemange, Gohlke & Cadot Reference Grandemange, Gohlke and Cadot2014; Rigas et al. Reference Rigas, Oxlade, Morgans and Morrison2014). The wake visits different azimuthal orientations
$\theta$
through sudden switches that are triggered by turbulent fluctuations, and that are both random and rare events: the residence time
$\Delta t$
between successive switches follows an exponentially decreasing distribution with no preferred frequency, and the mean residence time
$\langle \Delta t\rangle$
is orders of magnitude larger than the convective time
$L/U_\infty$
. Qualitatively similar dynamics is observed in non-axisymmetric bodies with a planar symmetry, such as the Ahmed body (Grandemange, Gohlke & Cadot Reference Grandemange, Gohlke and Cadot2013; Barros et al. Reference Barros, Borèe, Cadot, Spohn and Noack2017). For both types of bluff bodies, static symmetry-breaking deflections of the wake are reminiscent of the primary, laminar pitchfork bifurcation. The understanding of the mechanism of the low-
${Re}$
instabilities may thus play an important role for the understanding of the flow physics at larger
${Re}$
and, for example, to locate the most sensitive regions in order to develop new flow control strategies. While the notion of critical Reynolds number may not be relevant in the turbulent regime, it would be interesting to analyse mean turbulent wakes through the lens of the present study. For example, one could evaluate the newly defined Reynolds numbers
$\widehat {{Re}}$
and
$\overline {{Re}}$
, and investigate how they depend on the geometry, aspect ratio and
${Re}$
. It may also be worth determining whether regions with
$\partial \omega _\theta /\partial r\lt 0$
can be sustained, or if the competition between the production and advection of vorticity is such that
$\partial \omega _\theta /\partial r \gt 0$
everywhere in the near wake.
Funding
This research received no specific grant from any funding agency, commercial or not-for-profit sectors.
Declaration of interests
The authors report no conflict of interest.
Appendix A. Characteristic base-flow quantities at
$\boldsymbol{Re_{c}}$
In § 5, we introduced two new Reynolds numbers, i.e.
$\widehat {{Re}} = \widetilde {\omega }_\theta r_0^2/\nu$
and
$\overline {{Re}} = \widehat {\omega }_\theta r_0^2/\nu = U_{rev} r_0/\nu$
. Figure 19 shows the base-flow quantities involved in the definitions of
$\widehat {{Re}}$
and
$\overline {{Re}}$
, and measured at
${Re}={Re}_c$
. They are the maximum reverse flow
$U_{rev}$
(figure 19
$a$
), the streamwise location
$z_{max}$
of maximum reverse flow (figure 19
$b$
), the maximum azimuthal vorticity
$\widetilde {\omega }_\theta$
(and, for reference, the maximum shear
$\widetilde {\tau }_{zr}$
) at
$z=z_{max}$
(figure 19
$c$
), and the radius
$r_0$
of zero streamwise velocity (and, for reference, the radius
$r_\omega$
of maximum azimuthal vorticity) at
$z=z_{max}$
(figure 19
$d$
). We note in particular that for bodies with a zero-thickness TE,
$\widetilde {\omega }_\theta$
increases with 
, and
$r_0$
decreases, while for bodies with a blunt TE,
$\widetilde {\omega }_\theta$
decreases and
$r_0$
remains approximately constant.
Figure 19. Base-flow quantities at the critical Reynolds number
${Re}_c$
.
$(a)$
Maximum reverse flow
$U_{rev}$
.
$(b)$
Streamwise location
$z_{max}$
of maximum reverse flow.
$(c)$
Maximum azimuthal vorticity
$\widetilde {\omega }_\theta$
(solid lines) and maximum shear
$\widetilde {\tau }_{zr}$
(dashed lines) at
$z=z_{max}$
.
$(d)$
Radius
$r_0$
of zero streamwise velocity (solid lines) and radius
$r_\omega$
of maximum azimuthal vorticity (dashed lines) at
$z=z_{max}$
.
Appendix B. Mesh convergence
Table 2 reports the critical Reynolds numbers obtained on three different numerical meshes, for all the geometries considered in the present study and for 
and
$4$
. From one mesh to the next, all mesh sizes are divided by
$\sqrt {2}$
, such that the number of elements
$N_{elmts}$
increases by a factor of approximately 2. Additionally, table 3 reports the base-flow quantities involved in the definition of the newly introduced Reynolds numbers
$\widehat {{Re}}$
and
$\overline {{Re}}$
. For illustration purposes, we arbitrarily chose the flow past the 
cone at
${Re}=160$
. In the present study, we used mesh M2, for which all results are very well converged.
Table 2. Values of the critical Reynolds numbers obtained on three different meshes M1–M3 for the four different geometries, and two aspect ratios each.
Table 3. Values of the base-flow quantities of figure 19 obtained on three different meshes M1–M3 for the 
cone at
${Re}=160$
.
Figure 20. Azimuthal vorticity
$\omega _\theta (z)$
on the lateral surface of a cone for 
at
${Re}=160$
.
$(a)$
Sharp leading edge, varying local mesh size
$h$
.
$(b)$
Rounded leading edge, varying fillet radius
$R$
.
Appendix C. Rounding sharp leading edges
In § 4.1, we reported the maximum vorticity at the wall for ellipsoids and bullets, i.e. bodies with a smooth leading edge. Here, we give more details about bodies with a sharp leading edge.
Figure 20
$(a)$
shows the azimuthal vorticity
$\omega _\theta (z)$
obtained numerically on the surface of an 
cone with a sharp leading edge, at
${Re}=160$
. Different curves correspond to different values of the local mesh size
$h$
at the leading edge. The vorticity computed at the sharp leading edge keeps increasing when the numerical mesh is refined locally.
Next, we investigate the effect of smoothing the leading edge with a fillet of radius
$R$
. We use a local mesh size
$h=R/10$
, and verify convergence with
$h=R/30$
. Figure 20
$(b)$
shows that the maximum vorticity is obtained slightly downstream of the fillet, and increases as
$R$
decreases. However, the critical Reynolds number remains essentially constant (not shown), which is consistent with all the curves
$\omega _\theta (z)$
collapsing downstream enough (for
$z \gtrsim 0.4$
in the present case). Therefore,
${Re}_c$
and the maximum vorticity are not related.
Open access
