3.1. Validation

The numerical approach adopted here is validated further by simulating flow around a sphere at Re _D $=\textit{UD}/\nu$ = 10 000 and comparing the results with data available in the literature (Rodriguez et al. Reference Rodriguez, Lehmkuhl, Soria, Gómez, Domínguez-Pumar and Kowalski2019). The results are presented in non-dimensional form, where the length of the minor axis (D) and the free stream velocity (U) have been used to normalize relevant variables. The time-averaged streamwise velocity along the centreline in the wake region ( $\lt u_{1}\gt$ ), the turbulence intensity in terms of the root mean square value of the streamwise velocity ( $u_{\textit{rms}}$ ) and the pressure coefficient ( $C_{p}$ ) distribution along the surface of the sphere are shown in figure 2. Here, $C_{p}$ is defined as follows:

(3.1) \begin{equation}C_{p}=\frac{p-p_{\mathrm{\infty }}}{\dfrac{1}{2}\rho U^{2}},\end{equation}

where, $p_{\infty }$ is the pressure at the centre line at the inlet boundary, and $\rho$ denotes the density of the fluid. Figure 2 indicates good agreement between the present simulations and the reference data. We note that oscillations in $u_{\textit{rms}}$ profiles for sphere wakes have been observed by several studies in the streamwise range $1\leq x/D\leq 3$ , even at significantly lower Reynolds numbers. For instance, $u_{\textit{rms}}$ data from Rodriguez et al. (Reference Rodriguez, Borell, Lehmkuhl, Perez Segarra and Oliva2011) at $Re_{D}=3700$ , and from Tomboulides & Orszag (Reference Tomboulides and Orszag2000) at $Re_{D}=300$ display similar oscillations. The authors of these studies note that the maximum in $u_{\textit{rms}}$ occurs near the tail-end of the recirculation region. This agrees well with the current data shown in figure 2, where the recirculation region ends at $x/D=2.05$ (obtained using the zero-crossing of the velocity from figure 2 a) and the maximum in $u_{\textit{rms}}$ occurs at $x/D=2.045$ (figure 2 b).

Figure 2. Validation for (a) time-averaged streamwise velocity in the wake along the domain centre line, (b) root mean square fluctuation of the streamwise velocity and (c) coefficient of pressure on the sphere surface. In these plots, x denotes the downstream distance from the centre of the ellipsoid normalized with D, and (d) θ is the azimuthal angle starting from the upstream stagnation point.

3.2. Boundary layer separation

Unlike for a sphere, the cross-sections and flow patterns for ellipsoids with non-unity values of AR differ in the meridional (x–y) and equatorial (x–z) planes. To account for these differences introduced by body anisotropy, the flow behaviour is examined and compared for various combinations of AR = H/D = 5 : 1, 5 : 2, 5 : 3, 5 : 4, and in some instances for 1 : 1 (i.e. a sphere), in the meridional and equatorial planes. These orthogonal planes are shown in figure 3 for the 5 : 1 ellipsoid. All the results discussed here are presented in non-dimensional form, where U and D have been used to non-dimensionalize the relevant quantities. Various properties within the boundary layer are examined first, along with the formation of the separated shear layers. The enstrophy transport equation is examined in a later section to study the phenomenon of enstrophy generation and depletion in the separated flow regions and in the wake. The interactions between the vorticity vector and the strain rate tensor are also examined at a later point to better understand the origin of enstrophy generation and depletion.

Figure 3. View of (a) the meridional (x–y) and (b) equatorial (x–z) planes passing through the centre of the 5 : 1 ellipsoid. The instantaneous z-vorticity and y-vorticity are shown in (a) and (b), respectively (red, positive; blue, negative).

Boundary layer separation is typically estimated to occur at the location where surface shear stress becomes zero. While this criterion provides a reasonable estimate for quasi-two-dimensional (quasi-2-D) flow separation, such as for cross-flow over cylinders with uniform cross-section, it may not always provide a complete description of 3-D separation (Maskell Reference Maskell1955; Wang Reference Wang1972). Separation in three-dimensions may occur through other mechanisms, such as through vortex roll-up observed on prolate spheroids inclined at moderate angles of attack. In such scenarios, the helicity density has been shown to be a useful metric for identifying separation locations (Chesnakas & Simpson Reference Chesnakas and Simpson1997; Plasseraud & Mahesh Reference Plasseraud and Mahesh2025). Notably, for large angles of attack, separation leads to the formation of a large circulation bubble as is typically observed in quasi-2-D cross-flow over circular cylinders (Wang Reference Wang1972). The angle of attack for the ellipsoids in the present work is effectively 90° with respect to the major axis, and a large separation bubble forms in the wake as opposed to the rollup type separation and reattachment observed for moderate angles of attack. Moreover, the boundary layer detaches along a separation line or curve that varies continuously with the surface curvature of the ellipsoids, instead of separating at a nearly constant mean elevation angle as in the case of cylinders.

Wetzel et al. (Reference Wetzel, Simpson and Chesnakas1998) carried out a detailed investigation of 3-D separation using various experimental techniques, including oil flow visualization, surface pressure measurements, skin-friction measurements and laser Doppler velocimetry. They demonstrated that even for complex 3-D separation scenarios, minima in wall shear stress magnitude correlate well with separation locations. Thus, separation is identified in the present work by locating minima in plots of wall (or surface) shear stress magnitude. For this purpose, the instantaneous value of shear stress was calculated as $\tau _{i}=T_{i}-(T_{j}n_{j})n_{i}$ . Here, $T_{i}$ denotes the traction vector calculated using the deviatoric part of the stress tensor ( $T_{i}=\mu ({\partial u_{i}}/{\partial x_{j}}+{\partial u_{j}}/{\partial x_{i}})n_{j}$ ), $n_{j}$ denotes the surface normal unit vector and the wall shear stress magnitude is given as $\tau _{w}=| \tau _{i}|$ . The wall shear stress was time-averaged for a duration of $400U/D$ , and resulting plots for $\tau _{w}$ are shown in figure 4 for both the meridional and equatorial planes for all five ellipsoids. The corresponding numerical values for elevation (α) and azimuthal (θ) separation angles, calculated using the minima in the $\tau _{w}$ curves, are provided in table 1. We observe that the elevation separation angle decreases with decreasing AR ratios, which indicates earlier separation of the boundary layer for increasing body symmetry (or correspondingly, for lower body anisotropy). The sphere undergoes earliest separation at elevation angle α ≈ 85° in the meridional plane, whereas the most elongated 5 : 1 ellipsoid undergoes later separation at α ≈ 89.6°.

Figure 4. The magnitude of surface shear stress ( $\tau _{w}$ ) in (a) the meridional plane (x–y) and (b) the equatorial plane (x–z). Here α and θ represent the elevation and azimuthal angles, respectively, and they vary from 0⁰ to 180⁰.

Table 1. Separation angle from the upstream stagnation point (in degrees) in the meridional and equatorial planes, and time-averaged drag coefficients for various AR values.

To better understand the observed separation behaviour, the pressure coefficient $C_{p}=(p-p_{\infty })/(1/2 \rho U^{2} )$ in the meridional (x–y) plane along the upper surface of the ellipsoids is shown in figure 5(a). For smaller AR ratios (5 : 3, 5 : 4 and 1 : 1), the flow experiences a favourable pressure gradient (FPG) along the surface up until α ≈ 75^°; a gradual decrease in pressure is observed up to this point followed by a gradual increase (adverse pressure gradient (APG)). Due to the FPG, the flow in the boundary layer gets accelerated, leading to an initial increase in surface shear stress in the meridional plane (figure 4 a). The subsequent APG region (figure 5 a) causes the flow to decelerate, and the surface shear stress decreases (figure 4 a). Due to this APG, the boundary layer undergoes separation at some point, and a constant surface pressure is observed downstream of this location (i.e. beyond α ≈ 90° in figure 5 a).

Figure 5. Surface pressure coefficient in (a) the meridional plane (x–y) and (b) the equatorial plane (x–z). Here α and θ correspond to the elevation and azimuthal angles, respectively, as shown in the schematic in figure 4.

For larger AR ratios (i.e. 5 : 1 and 5 : 2) the pressure initially decreases along the surface up until α ≈ 80° (figure 5 a). Beyond this point the flow undergoes an abrupt rise in pressure. The steep decrease and increase in pressure for $\alpha \gt 60^{\circ}$ , when compared with more gradual pressure changes observed for smaller AR ratios (5 : 3, 5 : 4 and 1 : 1), are a result of the sharp curvature near the upper and lower poles of the 5 : 1 and 5 : 2 ellipsoids. This geometric feature first gives rise to a large FPG and subsequently to a large APG. The APG leads to boundary layer separation at some point, beyond which the pressure undergoes minor fluctuations. These fluctuations are not observed for the smaller AR ratios (5 : 3, 5 : 4 and 1 : 1) discussed previously, for which the postseparation surface pressure asymptotes to constant values. The plots also suggest a correlation between wall shear stress magnitude and surface curvature. Regions of high curvature, e.g. $\alpha \gt 60^{\circ}$ for the 5 : 1 and 5 : 2 ellipsoids, experience stronger flow acceleration and hence higher wall shear stress magnitude (figure 4 a), whereas regions of low curvature experience comparatively lower shear stress due to weaker flow acceleration. Overall, these observations provide an initial indication of how the severity of body anisotropy influences boundary layer behaviour, separation angle and postseparation characteristics in the meridional plane.

For considering separation characteristics in the equatorial plane, we note that the equatorial cross-sections for all five ellipsoids are geometrically identical regardless of the AR ratio, since the minor-axis diameter is kept constant. This geometric similarity may explain why the profiles of $\tau _{w}$ plotted against the azimuthal angle $\theta$ do not show significant variation with AR in figure 4(b). However, the azimuthal separation angle $\theta$ , determined using the locations of minima in figure 4(b), exhibits a measurable increase with decreasing AR ratio (table 1). This is contrary to the trend observed in the meridional plane where separation occurs earlier for lower body anisotropy (i.e. for the sphere), but separation in the equatorial plane occurs earlier for higher body anisotropy (i.e. for the most elongated ellipsoid). The 5 : 1 ellipsoid undergoes separation at θ ≈ 76.4° in the equatorial plane, whereas the sphere undergoes separation at θ ≈ 85.2°. This observation can be related to $C_{p}$ reaching its minima earlier for higher AR ratios in figure 5(b). Unlike the abrupt pressure changes observed in the meridional plane for high AR values, the pressure varies gradually in the equatorial plane for all values of AR (figure 5 b). This can be related to the surface curvature remaining constant in the equatorial plane, whereas rapid curvature changes occur in the meridional plane near the poles for the 5 : 1 and 5 : 2 ellipsoids.

In addition to considering separation angles in the meridional and equatorial planes, flow separation over the entire surface of the ellipsoids can be visualized in three-dimensions by locating minima in wall shear stress. Figure 4 indicates that the wall shear stress is close to zero at the minima locations. Thus, contours of zero wall shear stress (excluding the forward stagnation point) were used to define separation curves on the ellipsoids’ surface. To determine whether these contours are consistent with separation locations, limiting streamlines were calculated for the 5 : 1 ellipsoid and the corresponding results are shown in Appendix A. The zero-stress contour aligns closely with the convergence of the limiting streamlines, which demonstrates that the separation curves provide a reasonable indication of boundary layer detachment.

To better understand the influence of body anisotropy on boundary layer detachment, the separation curves for the 5 : 1 and 5 : 4 ellipsoids are compared in figure 6. These curves appear as narrow bands instead of lines since a small range of shear stress values, rather than a single zero value, was used to locate the minima. This approach helps mitigate the effects of noise near the separation location. In addition to the 3-D representations shown in figures 6(a) and 6(b), 2-D projections of the points that generate the separation curves are shown in figures 6(c) and 6(d). These projections were obtained by considering various horizontal cross-sections parallel to the equatorial plane and identifying the azimuthal separation angle in each such cross-section. To depict the resulting information on the polar plot, the distance of each cross-section from the upper pole of the ellipsoid (normalized by D) was used to mark radial distance on the polar plot, and the azimuthal angle was used to mark the polar angle. The resulting plots show a continuous increase in azimuthal separation angle for both ellipsoids as we move from the equator to the poles. However, higher body anisotropy corresponds to comparatively earlier separation for horizontal cross-sections farther from the pole. We observe more noticeable differences between the two ellipsoids closer to the equator than to the poles, which is somewhat unexpected since the equatorial cross-sections are geometrically identical between the 5 : 1 and 5 : 4 ellipsoids. The early separation closer to the equator for the 5 : 1 ellipsoid is an indirect consequence of a larger wake signature in the height dimension (major axis), which in turn leads to a larger overall wake envelope in three dimensions and consequently to a broader wake in the width dimension (minor axis). The broader wake width corresponds to the boundary layer separating earlier as we move closer to the equator.

Figure 6. (a,b) Separation curves and (c,d) polar distribution of the corresponding points in various cross-sections parallel to the equatorial plane, shown for (a,c) the 5 : 1 ellipsoid and (b,d) the 5 : 4 ellipsoid. The radial distance in the polar plots indicates non-dimensional vertical distance from the major axis poles of the respective ellipsoids. The polar angle is measured with respect to the upstream stagnation point on the ellipsoids.

Quantitative differences in wake-height and -width for the two ellipsoids are examined in figure 7. Although the shear layers form a 3-D surface around the periphery of the ellipsoids, 2-D cross-sections can offer valuable insights regarding the influence of body anisotropy on the shape and width of the wake. The wake size in the two orthogonal planes is examined using the plane-normal vorticity components along cross-stream line cuts located 1.5D downstream of the ellipsoid centres. The maxima/minima in vorticity occur near the core of the shear layers, and the distance between the extrema is smaller for the 5 : 4 ellipsoid than for the 5 : 1 ellipsoid. This indicates that the shear layer of the 5 : 4 ellipsoid diverges away from the wake centreline to a smaller degree, resulting in a narrower wake.

Figure 7. (a) Time-averaged non-dimensional z-vorticity in the meridional plane and (b) y-vorticity in the equatorial plane, taken from line cuts located at x = 1.5D downstream of the ellipsoid centres.

3.3. Wake velocity and drag coefficient

As the boundary layer separates and the 3-D shear layer is formed, a recirculation region develops in the near-wake of the ellipsoids. Velocity profiles taken along line cuts in the wake provide an indication of the size of these recirculation regions as well as of other general wake characteristics, albeit only in a single dimension along the line cut. The time-averaged streamwise velocity profiles ( $\langle u_{1}\rangle$ ) shown in figure 8(a) indicate that for all five ellipsoids considered, there is an initial flow reversal corresponding to negative velocity values, followed by a recovery to various levels that are lower than the free stream velocity. The size of the recirculation region can be estimated using the location where the velocity first recovers to zero. It is evident that the size increases monotonically with increasing aspect ratio, and that the most elongated 5 : 1 ellipsoid generates a recirculation bubble that is substantially larger than the rest. The velocity profiles also provide an indication of the intensity of recirculation, with more negative minima corresponding to stronger recirculation bubbles. The trend in this regard suggests an increase in recirculation intensity with increasing aspect ratio.

Figure 8. (a) Time-averaged streamwise velocity and (b) velocity defect along the centreline in the wake (y = 0, z = 0). The horizontal axis indicates non-dimensional distance from the rear surface of the ellipsoids. The dot symbols correspond to reference LES data from Rodriguez et al. (Reference Rodriguez, Lehmkuhl, Soria, Gómez, Domínguez-Pumar and Kowalski2019). The intersection of the dashed line in (a) with the solid lines provides an indication of the size of the recirculation bubble.

The velocity defect along the wake centreline, defined as $u_{d}=(U-u_{1})/U$ , is shown in figure 8(b). Two distinct decay rates of the time-averaged $\langle u_{d}\rangle$ are observed for all cases except the 5 : 1 ellipsoid. For the sphere with AR = 1 : 1, the decay of $\langle u_{d}\rangle$ exhibits a $x^{-2}$ power law in the near-wake region for $1.33D\leq x\leq 3D$ (with x measured from the downstream surface), followed by a slower decay of $x^{-3/4}$ for $x\geq 4D$ . A rapid decay in velocity defect indicates faster recovery towards free stream conditions, and it is typically associated with a smaller wake and lower drag. A slower initial decay rate is indicative of slow recovery, a large wake, and higher drag. The decrease in decay rate beyond $x\geq 4D$ for the sphere is a result of the wake velocity having already recovered substantially within $x\leq 3D$ , which is followed by a very gradual asymptote towards the free stream value (figure 8 a). The 5 : 1 ellipsoid displays a consistent $x^{-1.13}$ decay rate for the majority of the wake, and the remaining asymmetric ellipsoids exhibit decay rates that lie between this case and that of the sphere. This trend in velocity defect suggests that the drag coefficient is higher for increasing body anisotropy, which is confirmed by the calculated drag values provided in table 1. The transition from two distinct decay rates to a single decay rate with increasing aspect ratio can be attributed to differences in the wake structure. As discussed previously, the 5 : 1 ellipsoid gives rise to a larger wake not just in the meridional plane but also in the equatorial plane. This can also be observed in Supplementary movie 1, which depicts the normal vorticity components in these two orthogonal planes. The influence of flow recirculation persists farther downstream for the 5 : 1 ellipsoid, resulting in a consistent but weaker decay of the velocity defect. For the 5 : 4 ellipsoid, which closely resembles the sphere but retains some degree of anisotropy, the influence of the recirculation region is limited in comparison, as can be observed in Supplementary movie 1. Thus, the streamwise velocity recovers rapidly beyond $x\gt 1.5D$ (figure 8 a) giving rise to a faster initial decay rate, followed by a slow asymptote towards the free stream value.

Apart from the wake velocity profiles discussed above, the variation in drag coefficient with body anisotropy is also examined. The time-averaged drag coefficient $C_{D}=F_{D}/(({1}/{2})\rho U^{2}A)$ for all five ellipsoids is provided in table 1. Here $F_{D}$ represents the net drag force and A represents the projected frontal area ( $A=\pi DH/4$ ). The contributions of viscous (or skin-friction) drag and pressure-induced drag to the total drag coefficient have also been listed in table 1. We observe that pressure drag shows an increasing trend for higher AR values. This is a result of the gauge pressure (not shown) in the near-wake being considerably more negative for higher values of AR than for lower values. A larger pressure difference between the upstream and downstream surface leads to higher pressure drag. The wake pressure can also be related indirectly to the velocity defect decay rates discussed previously; a slower decay rate signifies a larger wake, and consequently, an extensive region of low wake pressure. Some variations are observed in the skin-friction drag coefficient, and its contribution is an order of magnitude smaller than that of pressure drag at the given Reynolds numbers. The magnitude of the pressure drag being dominant leads to a monotonic increase in total C _D with increasing values of AR.

3.4. Shear layer characteristics and flow topology

The separated shear layers are examined using instantaneous vorticity contours in the meridional and the equatorial planes in figure 9. The current and subsequent sections focus primarily on the 5 : 1 and 5 : 4 ellipsoids since these two shapes represent the highest and lowest non-zero body anisotropy, respectively. For both the ellipsoids, the separated shear layers were observed to remain coherent for a small distance downstream before undergoing roll-up, with subsequent breakdown giving rise to a broad range of spatial and temporal scales in the wake. The roll-up locations were identified using instantaneous vorticity contours, where a roll-up event was defined as the point where the shear layers in the two orthogonal plane cuts transitioned from a relatively flat profile to forming a discrete rotating material element. The centre coordinates of these discrete material elements were identified prior to pinch-off, and the values were averaged over a duration of 10 non-dimensional time units to obtain the roll-up locations shown in figure 10. We observe that the roll-up process in the meridional plane starts at comparable distances downstream of the ellipsoid centres for both AR values, specifically between x/D ≈ 0.8–1 on average. Moreover, there is very little variation of the roll-up location in the y direction, as is evident from the barely perceptible vertical error bars for the red symbols. In the equatorial plane, roll-up occurs earlier for the 5 : 1 ellipsoid at around x/D ≈ 0.3 on average, but the shear layer stays coherent farther downstream for the 5 : 4 ellipsoid, up to approximately x/D ≈ 0.8. There is significant variation in both the streamwise (x) and cross-stream (z) location of roll-up in the equatorial plane, which matches qualitative observations from Supplementary movie 1.

Figure 9. (a,b) Instantaneous z-vorticity in the meridional plane, and (c,d) instantaneous y-vorticity in the equatorial plane. Panels (a) and (c) correspond to the 5 : 1 ellipsoid, and (b) and (d) correspond to the 5 : 4 ellipsoid. The vorticity has been non-dimensionalized by dividing with $U/D$ , and the axes have been non-dimensionalized using D. An animation of this figure is provided in Supplementary movie 1.

Figure 10. Shear layer roll-up locations for the 5 : 1 ellipsoid (circles) and the 5 : 4 ellipsoid (squares), recorded for multiple roll-up events in the equatorial plane (black) and in the meridional plane (red). The values have been non-dimensionalized using D, and they indicate distance from the ellipsoids’ centres. The cross-stream coordinate locations are plotted along the left-hand vertical axis for the equatorial plane ( $| z|$ ), and along the right-hand vertical axis for the meridional plane ( $| y|$ ). The mean values are plotted using open symbols, and the error bars indicate the standard deviation in both directions.

To examine strong rotational structures that form as the shear layers evolve in space and time, coherent structures in the wakes of the 5 : 1 and 5 : 4 ellipsoids are visualized using the Q-criterion in figure 11 (also shown in Supplementary movie 2). Larger Q values are associated with stronger rotational structures, and the strongest structures appear closer to the equator than to the poles for both the ellipsoids. The near-wakes of the ellipsoids are largely devoid of Q structures, although additional structures would become visible if contours with smaller Q values were used. The absence of Q structures suggests that the strongest tube-like rotational structures are typically absent in the near-wake where the shear layers remain cohesive. However, early breakdown of the shear layer in the equatorial plane for the 5 : 1 ellipsoid (discussed previously) leads to strong near-wake vortex tubes forming close to the equator in figure 11(a). The strongest vortex tubes (i.e. contours with larger Q values) form at a certain standoff distance, approximately 2.5D downstream of the ellipsoid centres for both the cases shown in figure 11. This can be related to enstrophy production reaching its maximum around this downstream distance, which will be discussed later in § 3.5.

Figure 11. Isosurfaces of the Q-criterion for (a) the 5 : 1 ellipsoid and (b) the 5 : 4 ellipsoid. The colours shown correspond to normalized values of $Q/(U/D)^{2}=50$ (grey), 100 (blue) and 200 (red). An animation of this figure is provided in Supplementary movie 2.

The characteristics of the resolved small-scale wake structures are examined in figure 12 using joint probability density functions (PDF) of the second and third invariants (Q and R) of the velocity gradient tensor (Chong, Perry & Cantwell Reference Chong, Perry and Cantwell1990; Soria et al. Reference Soria, Sondergaard, Cantwell, Chong and Perry1994; Ooi et al. Reference Ooi, Martin, Soria and Chong1999). Additional details regarding the relationship of Q and R to the local flow topology are provided in Appendix B. We observe in figure 12 that the Q–R joint PDF for the 5 : 1 and 5 : 4 ellipsoids bear some similarities, and both display a characteristic ‘teardrop’ shape which is typically observed in a variety of turbulent flow configurations. The teardrop shape was observed consistently across all five ellipsoids examined in this study, with minor influence of body anisotropy on the distribution of points. This indicates that despite notable differences in geometry, the resolved scales in the wakes display characteristics that are typical of universal small-scale behaviour in turbulent flows. Some noticeable differences between the two cases shown in figure 12 are a broader spread of the distribution for the 5 : 1 ellipsoid, as well as a higher probability of encountering larger values of Q and R. The spread can be quantified using the standard deviation values ${\unicode[Arial]{x03C3}} _{Q}$ and ${\unicode[Arial]{x03C3}} _{R}$ , which are 1.47 times and 1.53 times larger for the 5 : 1 ellipsoid than for the 5 : 4 ellipsoid. The mean value of Q is positive for both cases, with $\mu _{Q}$ being 1.78 times larger for the 5 : 1 ellipsoid. This suggests that vorticity dominates over strain in the wakes of both ellipsoids, but that the vortex tubes are considerably stronger for the 5 : 1 ellipsoid. Overall, these observations provide a quantitative indication that higher body anisotropy in cross-flow leads to the formation of stronger rotational coherent structures, which matches qualitative observations from figure 11.

Figure 12. Joint PDF of Q and R for (a) the 5 : 1 ellipsoid and (b) the 5 : 4 ellipsoid. Both Q and R have been normalized using three times their respective standard deviations ( ${\unicode[Arial]{x03C3}} _{Q}$ and ${\unicode[Arial]{x03C3}} _{R})$ .

3.5. Enstrophy production

The spatial and temporal evolution of wake structures are linked closely to the production and dissipation of vorticity, and thus, to the production of enstrophy. Enstrophy is defined as $\xi =0.5| \boldsymbol{\nabla }\times \boldsymbol{u}| ^{2}$ and it provides a quantitative measure of vorticity strength. The transport equation for enstrophy is provided in Appendix C, and the main focus of this section is on the enstrophy production term $\omega _{i}S_{\textit{ij}}\omega _{j}$ . Using eigendecomposition of the symmetric matrix $S_{\textit{ij}}$ , this term can be decomposed into contributions from the three eigenvalues s _i of the strain-rate tensor, and the alignment between the corresponding eigenvectors and the vorticity vector (Buxton & Ganapathisubramani Reference Buxton and Ganapathisubramani2010):

(3.2) \begin{equation}\xi _{prod}=\omega _{i}S_{\textit{ij}}\omega _{j}=\omega ^{2}s_{i}\left(\hat{\boldsymbol{\omega }}.\boldsymbol{e}_{i}\right)^{2}.\end{equation}

Here, $\omega$ represents the magnitude of vorticity, $\hat{\boldsymbol{\omega }}$ is the unit vector in the direction of vorticity and $\boldsymbol{e}_{\boldsymbol{i}}$ are the three eigenvectors of the strain-rate tensor. When s _i is arranged in ascending order, the first eigenvalue $s_{1}$ is always negative and the third eigenvalue $s_{3}$ is always positive. The sum of the three eigenvalues must be zero to satisfy the continuity equation. For the enstrophy production term $\omega ^{2}s_{i}(\hat{\boldsymbol{\omega }}.\boldsymbol{e}_{\boldsymbol{i}})^{2}$ in (3.2), we observe that only the term $s_{i}$ can be negative. Hence, the positive eigenvalue $s_{3}$ always contributes to enstrophy generation whereas the negative eigenvalue $s_{1}$ always contributes to enstrophy depletion. The intermediate eigenvalue $s_{2}$ can contribute locally either to enstrophy generation or to depletion, depending on whether its value is positive or negative. A similar eigendecomposition was used by Watanabe et al. (Reference Watanabe, Sakai, Nagata, Ito and Hayase2014) in DNS of a planar turbulent jet to demonstrate that positive enstrophy production becomes less dominant along the upstream-oriented (trailing) edge of the turbulent/non-turbulent interface. This results in a marked reduction in enstrophy generation along the trailing edge, due to an increased prevalence of vortex compression events associated with local negative enstrophy production. Nevertheless, the net enstrophy production along the trailing edge was observed to remain weakly positive.

The expression for enstrophy production in (3.2) can be simplified further as $\xi _{prod}=\omega ^{2}(s_{1}\cos ^{2} \phi _{1}+s_{2}\cos ^{2} \phi _{2}+s_{3}\cos ^{2} \phi _{3})$ , where $\cos \phi _{i}=\hat{\boldsymbol{\omega }}.\boldsymbol{e}_{i}$ are the direction cosines of the vorticity vector in the eigenframe. It is known that the vorticity vector aligns preferentially with the intermediate eigenvector $\boldsymbol{e}_{2}$ in homogeneous isotropic turbulence, and that it is orthogonal to the compressive eigenvector $\boldsymbol{e}_{1}$ and shows no preferential alignment with the extensive eigenvector $\boldsymbol{e}_{3}$ (Ashurst et al. Reference Ashurst, Kerstein, Kerr and Gibson1987). Preferential alignment with $\boldsymbol{e}_{2}$ is observed almost universally across a broad range of turbulent flow configurations. While this leads to the direction cosine $\cos \phi _{2}$ being the largest among the three, the intermediate eigenvalue $s_{2}$ is known to be the smallest in magnitude. Thus, it is difficult to attribute the contribution of any one of the three eigenvalues or eigenvectors as being dominant to enstrophy production. The behaviour of the three direction cosines and the corresponding eigenvalues is analysed later in this section to better understand their individual roles.

Figure 13 shows the normalized time-averaged enstrophy production $P_{\xi }=\lt \omega _{i}S_{\textit{ij}}\omega _{j}\gt /(U/D)^{3}$ in the meridional and equatorial planes for the four asymmetric ellipsoids used in this study. The most prominent differences appear in the meridional planes due to greatest dissimilarity of the corresponding cross-section profiles. For low body anisotropy, i.e. for the 5 : 3 and 5 : 4 ellipsoids, we observe shear layers from opposite poles merging into a single region of high enstrophy production approximately 2.5D downstream of the ellipsoid centres. With increasing body anisotropy, the high production region begins to separate out towards the poles, and two distinct regions are visible for the 5 : 1 ellipsoid in the meridional plane. Another notable difference in the meridional planes is the presence of negative production regions for the highly anisotropic 5 : 1 and 5 : 2 ellipsoids close to the poles. However, $P_{\xi }$ remains positive beyond x/D> 1.5 for these ellipsoids, indicating that the flow primarily undergoes vortex stretching in the intermediate- and far-wake regardless of body anisotropy. It is known that positive enstrophy production dominates overall in turbulent flows due to the prevalence of vortex stretching, which is related to the forward energy cascade. However, enstrophy production can be negative in localized regions, corresponding to strong vortex compression (Buxton & Ganapathisubramani Reference Buxton and Ganapathisubramani2010; Bechlars & Sandberg Reference Bechlars and Sandberg2017). As discussed in Appendix C, negative production indicates a local reduction in enstrophy and a corresponding suppression of strong vortex tubes. This can be confirmed qualitatively from figure 11 and from Supplementary movie 2, where large Q-value contours are absent close to the upper and lower poles of the 5 : 1 ellipsoid.

Figure 13. Time-averaged enstrophy production $P_{\xi }=\lt \omega _{i}S_{\textit{ij}}\omega _{j}\gt /(U/D)^{3}$ in (a–d) the meridional planes, and (e–h) the equatorial planes for the four asymmetric ellipsoids. The values for enstrophy production have been normalized by dividing with $(U/D)^{3}$ and the axes have been normalized using D.

Unlike in the meridional plane, enstrophy production is primarily positive in the equatorial plane for all four ellipsoids, indicating a prevalence of vortex stretching regardless of the degree of body anisotropy. Similar to the behaviour observed in the meridional plane, the shear layers merge into a single high production region for low body anisotropy in the equatorial plane, but they remain separated into two distinct regions for higher body anisotropy. Overall, there are significant differences visible in the equatorial plane among the four ellipsoids, even though all four equatorial cross-sections are geometrically identical. The magnitude of enstrophy production appears to increase noticeably with increasing body anisotropy, which can be correlated qualitatively to the presence of strong rotational structures (i.e. contours of large Q-value) near the equator in figure 11(a). Higher enstrophy production also appears to correlate well qualitatively with the shear layer being more susceptible to instabilities resulting in roll-up and breakdown, as can be seen in figure 9(c).

While the 2-D plane cuts shown in figure 13 provide a preliminary indication of the spatial distribution of enstrophy production, examining 3-D structures can provide a more comprehensive picture of the origin of negative enstrophy production and of large positive production. Figure 14 shows various contours of the positive time-averaged enstrophy production for the 5 : 1 and 5 : 4 ellipsoids. Three distinct contour levels with non-dimensional production values $P_{\xi }=$ 10, 100 and 200 have been visualized in the figure. We observe that for both ellipsoids, regions of highest positive $P_{\xi }$ (red contours) are generally concentrated in the core of the wake, approximately 2.5D downstream from the ellipsoids’ centres. It also appears that the positive $P_{\xi }$ contours exhibit shell-like structures, with high magnitude $P_{\xi }$ contours nested within low magnitude $P_{\xi }$ contours. The strongest $P_{\xi }$ shells do not connect to the ellipsoids’ surface, indicating that positive enstrophy production reaches its maximum in the intermediate wake. As the value of $P_{\xi }$ decreases, the shells extend towards the ellipsoids’ surface, with the lowest $P_{\xi }=10$ shells wrapping around the surface, enveloping the boundary layer and the separated shear layer. The most prominent difference between the highly anisotropic 5 : 1 ellipsoid and the low anisotropy 5 : 4 ellipsoid is the separation of the $P_{\xi }=200$ shell into two distinct lobes closer to the poles for the 5 : 1 ellipsoid, which was discussed previously using 2-D plane cuts in figure 13. Overall, figure 14 indicates that positive enstrophy production remains relatively low in the near-wake region regardless of body anisotropy. This is followed by a gradual increase, reaching a maximum around x = 2.5D for the cases shown here, followed by a gradual decline farther downstream.

Figure 14. (a) Contours showing positive time-averaged enstrophy production in the wake of the 5 : 1 ellipsoid. The contours have been cut away along the meridional plane to show the internal 3-D structure of positive production, and the contours correspond to non-dimensional values of $P_{\xi }= 10$ (white), $100$ (pink) and $200$ (red). (b) Positive enstrophy production for the 5 : 4 ellipsoid, with the same contour levels as shown in (a).

The contours of negative time-averaged enstrophy production for the 5 : 1 and 5 : 4 ellipsoids are shown in figure 15. For both the ellipsoids, $P_{\xi }\leq -10$ is constrained to thin regions around the periphery. For the 5 : 1 ellipsoid, the negative production region takes the form of a hood-like structure, with a large and prominent canopy forming near the poles that diminishes towards the equator. The finger-like striations observed closer to the equator arise from disproportionate amplification of discretization errors when computing enstrophy production, since $\omega _{i}S_{\textit{ij}}\omega _{j}$ is cubic in velocity gradients. However, such oscillations are not present in the hood-like canopy close to the pole of the ellipsoid. Similarly to the positive production shells shown in figure 14, the negative production contours with higher $P_{\xi }$ magnitude are nested within lower magnitude $P_{\xi }$ contours. However, differently from positive production, negative time-averaged production attains its maximum magnitude in the near-wake close to the ellipsoid surface, and not in the intermediate wake. It is also evident from figure 15 that negative production regions for the 5 : 4 ellipsoid are essentially absent. This suggests that high body anisotropy may be more conducive to negative enstrophy production being sustained over extended periods of time in localized spatial regions. Sustained negative production is not a common occurrence in incompressible turbulent flows since vortex stretching is expected to be dominant in general. Thus, the source of negative production for the 5 : 1 ellipsoid is investigated in further detail.

Figure 15. (a) Side view of regions of negative time-averaged enstrophy production for the 5 : 1 ellipsoid, with contour levels corresponding to $P_{\xi } = -10$ (light blue), $-50$ (blue) and $-100$ (violet). (b) Negative production for the 5 : 4 ellipsoid, with the same contour levels as in (a). Panels (c) and (d) show downstream views of the images in (a) and (d).

To determine whether the length scales associated with persistent negative production regions are substantially larger than the LES filter scale, an instantaneous velocity snapshot was filtered a posteriori using a Gaussian kernel of filter width $16{\unicode[Arial]{x0394}}$ . Here, ${\unicode[Arial]{x0394}}$ is the local cell size defined as the cube root of each cell’s volume. Instantaneous negative enstrophy production regions were calculated using the filtered velocity field, and the hood-like canopy near the pole of the 5 : 1 ellipsoid remained clearly visible, indicating that this feature resides on length scales substantially larger the LES filter width. This large-scale nature also suggests that the persistent negative production region is only weakly influenced by the subgrid model, which is consistent with the minimal model contribution observed in the near-wake (Appendix A).

To better understand the cause of sustained negative production near the poles of the 5 : 1 ellipsoid, streamlines were generated near the upper pole using an instantaneous velocity snapshot. A spanwise line of tracer particles was released upstream of the ellipsoid, and the resulting streamlines are shown coloured using the normalized instantaneous enstrophy production $\xi _{prod}/(U/D)^{3}$ in figures 16(a) and 16(b). A relatively large cohesive region of instantaneous negative production is visible in the near-wake region in figure 16(b), and it has been highlighted using a green rectangle. As described earlier, enstrophy production can be decomposed into the three contributing terms, $\omega ^{2}s_{i}\cos ^{2} \phi _{i}$ , and only the compressive and intermediate eigenvectors $s_{1}$ and $s_{2}$ can give rise to negative production since $s_{3}$ is always positive. Figure 16(c) shows PDFs of the magnitude of the three direction cosines, $\cos \phi _{i}$ , within the highlighted rectangular region. It is evident that the eigenvectors $\boldsymbol{e}_{1}$ and $\boldsymbol{e}_{3}$ are both almost perfectly orthogonal to the vorticity vector in this region, whereas the intermediate eigenvector $\boldsymbol{e}_{2}$ shows near perfect alignment with vorticity. This indicates that enstrophy production in the highlighted region is determined almost entirely by the intermediate eigenvector term, $\omega ^{2}s_{2}\cos ^{2} \phi _{2}$ . Thus, the intermediate eigenvalue is examined in figure 16(d) by colouring the streamlines with $s_{2}$ . The positive and negative values of $s_{2}$ show a near perfect correspondence with the sign of enstrophy production (figure 16 b) in the highlighted region. This is expected, since $\omega ^{2}s_{2}\cos ^{2} \phi _{2}$ is the only term with a meaningful contribution to enstrophy production, as discussed above.

Figure 16. (a) Side view of streamlines generated by tracer particles released upstream of the 5 : 1 ellipsoid in a single velocity snapshot. The streamlines have been coloured using the instantaneous normalized enstrophy production $\xi _{prod}/(U/D)^{3}$ – colourbar shown in (b). (b) Top view of the streamlines from (a). A green rectangle highlights the region where a large cohesive region of negative production is observed. (c) The PDFs of the magnitude of the three direction cosines in the highlighted rectangular region, with the line for $| cos \phi _{1}|$ shown in blue, $| cos \phi _{2}|$ shown in black and $| cos \phi _{3}|$ shown in red. (d) Top view of the streamlines coloured with the intermediate eigenvector $s_{2}$ normalized using U/D.

To assess whether the long-term behaviour of the direction cosines, $\cos \phi _{i}$ , matches their instantaneous behaviour, time-averaged values of $|\lt \cos \phi _{i}\gt |$ were examined in the negative enstrophy region defined by $P_{\xi }\leq -10$ in figure 15. Similarly to the instantaneous case, the magnitudes of time-averaged compressive and extensive direction cosines were found to be significantly smaller than that of the intermediate direction cosine. This indicates that over the long term, vorticity remains nearly perfectly aligned with the intermediate eigenvector within the highlighted region of interest, and it remains orthogonal to both the compressive and extensive eigenvectors. Consequently, the long-term enstrophy production is expected to be governed by the same $\omega ^{2}s_{2}\cos ^{2} \phi _{2}$ term that dominates instantaneous production as described above.

The eigenvalue $s_{2}$ being negative in the region of interest (figure 16 d) indicates that the flow experiences compression in two orthogonal directions (since $s_{1}\lt 0$ always) and stretching in one orthogonal direction (since $s_{3}\gt 0$ always). In figure 16, as the horizontal streamlines approach the highly anisotropic upper pole of the 5 : 1 ellipsoid, they experience severe deflection in the vertical (y) direction (figure 16 a) in addition to lateral deflection in the spanwise (z) direction (figure 16 b). Immediately downstream of the ellipsoid, the streamlines are pulled back in sharply from multiple directions into the near-wake, before evolving into convoluted trajectories; the streamlines at the upper pole bend downward in the vertical direction and they contract inward in the spanwise direction. This pronounced cross-stream wake contraction manifests as the flow experiencing two orthogonal compressive eigenvalues, i.e. as $s_{1}$ and $s_{2}$ both being negative in the highlighted region of interest. The severity of streamline bending was quantified by calculating the curvature of streamlines passing close to the upper pole as $\boldsymbol{\kappa }=\boldsymbol{t}\boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{t}$ , where $\boldsymbol{t}=\boldsymbol{u}/| \boldsymbol{u}|$ is the unit velocity vector. The $\kappa _{y}$ component of curvature was observed to be strongly negative within the highlighted rectangular region shown in figure 16, which provides quantitative confirmation that the streamlines bend downward sharply, resulting in wake contraction along the major axis. On average, the magnitude of the $\kappa _{z}$ component was 0.53 times that of $\kappa _{y}$ within the highlighted region, indicating comparatively weaker bending in the lateral direction. The spatial distribution of positive and negative $\kappa _{z}$ values was consistent with a spanwise contraction.

Having established the central role of the intermediate eigenvalue and eigenvector in negative enstrophy production for the 5 : 1 ellipsoid, we now examine the local flow topology using data from an instantaneous snapshot. The relationship between the velocity invariants Q and R and local flow topology is described briefly in Appendix B. Figure 17 shows the prevalence of the four canonical topologies, namely, unstable focus/compressing (UF/C), stable focus/stretching (SF/S), stable node/saddle/saddle (SN/S/S) and unstable node/saddle/saddle (UN/S/S), conditioned on various thresholds of the instantaneous enstrophy production. Probability distributions for negative values of enstrophy production are shown in figure 17(a), and those for positive production are shown in figure 17(b). The results in figure 17(a) indicate that negative enstrophy production is overwhelmingly dominated by the UF/C topology, with its prevalence increasing as production becomes more negative. The remaining three topologies contribute only marginally to negative production. Positive enstrophy production (figure 17 b) is dominated by the SF/S topology, along with a notable presence of the UN/S/S topology. As production becomes more positive, the prevalence of the SF/S topology increases while that of the UN/S/S topology decreases.

Figure 17. Probability distribution of the four incompressible flow topologies for the 5 : 1 ellipsoid, conditioned on (a) negative and (b) positive instantaneous enstrophy production. In (a) the colours correspond to regions where $\xi _{prod}/(U/D)^{3}\leq -1$ (blue), $\leq -1$ 0 (red), $\leq -5$ 0 (green) and $\leq -1$ 00 (black). In (b) the colours correspond to regions where $\xi _{prod}/(U/D)^{3}\geq 1$ (blue), $\geq 1$ 0 (red), $\geq 5$ 0 (green) and $\geq 1$ 00 (black).

The correspondence between negative enstrophy production and local flow topology was examined by Buxton & Ganapathisubramani (Reference Buxton and Ganapathisubramani2010) using stereoscopic particle image velocimetry of an axisymmetric turbulent jet. They found that enstrophy production was predominantly negative within flow regions characterized by the UF/C topology, and that regions associated with other flow topologies displayed a smaller presence of negative production, which is consistent with the observations in figure 17(a). Bechlars & Sandberg (Reference Bechlars and Sandberg2017) also reported similar results from a turbulent boundary layer DNS, where persistent negative production was observed primarily in regions of UF/C topology. Overall, the results presented here indicate that the negative production region is dominated by the unstable focus/compressing topology, where the flow spirals outwards from a focus point and streams inwards along the axis of the spiral. This topological representation is consistent with earlier discussion of the role of the intermediate eigenvalue and eigenvector in negative enstrophy production; the physical interpretation of vorticity being aligned with the intermediate eigenvector, together with the intermediate eigenvalue being negative, is a vortex that is undergoing compression, which is consistent with the UF/C topology. We note that while the observations presented here hold for Re _D = 10 000, persistent negative enstrophy production was also observed near the poles of the 5 : 1 ellipsoid in a separate simulation at a lower Re _D = 5000 (data not shown). Whether or not local negative production persists at significantly higher Reynolds numbers in bluff-body flows remains an open question.