3.1.1. Unsteady Bickley jet

For our first example we consider the unsteady Bickley jet. This is a 2-D geophysical model of a quasi-periodic zonal jet with adjacent migrating eddies whose Lagrangian dynamics have been studied in depth (see, e.g. Del-Castillo-Negrete & Morrison Reference Del-Castillo-Negrete and Morrison1993; Rypina et al. Reference Rypina, Brown, Beron-Vera, Koçak, Olascoaga and Udovydchenkov2007). This flow was also used as a benchmark for previous LCS comparisons (Hadjighasem et al. Reference Hadjighasem, Farazmand, Blazevski, Froyland and Haller2017). Its time-dependent stream function is given by

(3.1)\begin{equation} \psi(x,y,t) ={-}UL\tanh\left(\frac{y}{L}\right) + UL\,{\rm sech}^2\left(\frac{y}{L}\right)Re\left[\sum_{n=1}^3 \epsilon_n\exp(-{\rm i}k_nc_nt)\exp({\rm i}k_nx)\right]. \end{equation}

We use the parameters from Rypina et al. (Reference Rypina, Brown, Beron-Vera, Koçak, Olascoaga and Udovydchenkov2007) found in table 1.

Table 1.

Bickley jet model parameters.

In figure 1 we compare ${\mathrm {r}\alpha }$ and ${\rm r}\overline {\mathrm {TSE}}$ with benchmarks of material rotation and stretching, the LAVD (Haller et al. Reference Haller, Hadjighasem, Farazmand and Huhn2016) and the FTLE (Haller Reference Haller2015) on a dense grid of trajectories ($3600\ell ^{-2}$, $\ell ={\rm \pi}$), for an integration time of 30 days. This corresponds with approximately six full rotations and significant translation for the advected eddies. We find that ${\mathrm {r}\alpha }$ is able to match the prominent features of the LAVD field by effectively identifying eddies boundaries as closed convex contours and separating them from the central jet. The ${\rm r}\overline {\mathrm {TSE}}$ ridges identify the edges of the central jet and edges of vortices as regions of significant stretching, similar to FTLE. Deviations primarily exist in the centres of the vortex cores where $\boldsymbol {v}_d(\boldsymbol {x}(t))$ does not evolve in a material way.

Figure 1.

Evaluation of $\mathrm {r}\alpha$ and $\mathrm {r}\overline {\mathrm {TSE}}$ (a,d) with material rotation (LAVD) and stretching (FTLE) (b,e), and traditional quasi-objective single trajectory diagnostics $\overline {\mathrm {TRA}}$ and $\overline {\mathrm {TSE}}$ (c,f), in the unsteady Bickley jet. Note the distinct central jet with rotating eddies on either side. Here $\overline {\mathrm {TRA}}$ and $\overline {\mathrm {TSE}}$ have difficulty accurately capturing material behaviour (stretching and rotation) due to the large bulk advection in the flow and temporal variability in the velocity field.

Previous work has shown how $\overline {\mathrm {TRA}}$ and $\overline {\mathrm {TSE}}$ can effectively reproduce the same LCS as LAVD and FTLE in slowly varying flows (Haller et al. Reference Haller, Aksamit and Bartos2021). In figure 1, we reveal how necessary this slowly varying assumption truly is. With the parameters in table 1, the Bickley jet is highly unsteady. For meaningful quasi-objective calculations, one requires $|\boldsymbol {v}_t|/|\ddot {\boldsymbol {x}}|<<1$, but only $41\,\%$ of trajectories experience $|\boldsymbol {v}_t|/|\ddot {\boldsymbol {x}}|<1$ and only $0.03\,\%$ of particle paths show $|\boldsymbol {v}_t|/|\ddot {\boldsymbol {x}}|<0.1$. This suggests that the coherent structures in this flow are likely evolving or travelling much faster than a Lagrangian particle is able to trace them out. In fact, $\overline {\mathrm {TRA}}$ designates the eddies in figure 1 as exhibiting relatively low rotation. In addition, $\overline {\mathrm {TSE}}$ also suggests nearly everything outside the central jet is undergoing significant stretching, in contrast to the thin hyperbolic regions seen as FTLE ridges. For this degree of unsteadiness and strong eddy advection, ${\mathrm {r}\alpha }$ and ${\rm r}\overline {\mathrm {TSE}}$ are more suitable than $\overline {\mathrm {TRA}}$ and $\overline {\mathrm {TSE}}$ for quantifying fluid rotation and stretching while requiring much less-structured data than LAVD and FTLE.

We now progressively and randomly downsample the number of trajectories to test the robustness of our objective diagnostics at low trajectory densities. In figure 2, we calculate ${\mathrm {r}\alpha }$ and ${\rm r}\overline {\mathrm {TSE}}$ as well as the particle trajectory length (Mancho et al. Reference Mancho, Wiggins, Curbelo and Mendoza2013) and initial particle speed $|\boldsymbol {v}(t_0)|$. Although trajectory length and particle speed are not objective, they are either easily computed or commonly used diagnostics for experimental LPT studies (see, e.g. Fu, Biwole & Mathis Reference Fu, Biwole and Mathis2015; Tauro, Piscopia & Grimaldi Reference Tauro, Piscopia and Grimaldi2017; Rosi & Rival Reference Rosi and Rival2018). After calculating each diagnostic, we reconstruct the full-resolution diagnostic field using radial basis functions, similar to Encinas-Bartos et al. (Reference Encinas-Bartos, Aksamit and Haller2022). In addition, we include the outermost closed convex LAVD contours from figure 1 as a visual reference to aid in the comparisons.

Figure 2.

Comparison of relative rotation (${\mathrm {r}\alpha }$) and stretching ($\mathrm {r}\overline {\mathrm {TSE}}$) with non-objective trajectory length and initial particle speed diagnostics. Each diagnostic was progressively and randomly downsampled with the full-resolution product created using radial basis functions. For each panel, we include the outermost closed convex LAVD contours to highlight relative agreement of each sparse diagnostic with ground-truth material rotation.

At a trajectory concentration of 360$\ell ^{-2}$, ${\mathrm {r}\alpha }$ looks relatively similar to its full resolution sampling, but with added noise. Five concentrations of high ${\mathrm {r}\alpha }$ are still distinguishable from the surrounding flow. The thin hyperbolic structures identified by ${\rm r}\overline {\mathrm {TSE}}$ in the Bickley jet, however, are already beginning to disappear. This degradation of hyperbolic coherent structures ${\rm r}\overline {\mathrm {TSE}}$ is similar to the findings of Haller et al. (Reference Haller, Aksamit and Bartos2021) for $\overline {\mathrm {TSE}}$ and Mowlavi et al. (Reference Mowlavi, Serra, Maiorino and Mahadevan2022) for $\overline {\mathrm {TSE}}$ and FTLE. The trajectory length reveals fluid particles in the central jet have the longest trajectories over this time window when compared with nearby particles in the adjacent eddies. The initial velocity magnitude $|\boldsymbol {v}|$ also appears largest in the centre of flow domain, and suggests some oscillating feature may be present.

Dropping the trajectory resolution by a factor of 100 to 3.5$\ell ^{-2}$, we again reconstruct the original diagnostic resolution and include the location of trajectory positions with white crosses. At this concentration, ${\mathrm {r}\alpha }$ still suggests five distinct local rotation maxima corresponding to the locations of the five eddies. The interpolated ${\rm r}\overline {\mathrm {TSE}}$ field no longer recreates any of the features seen at full resolution and hence provides no meaningful information. The trajectory length is also contradicting features that were revealed with 100 times the number of trajectories. At this resolution, we are testing both the ability of the metric to provide useful flow descriptors, as well as the interpolation scheme. We believe, however, that there is limited influence from the interpolation scheme because $|\boldsymbol {v}|$ is practically constant across all resolutions. We further test the impact of interpolation schemes in the Appendix, and refer the reader to Encinas-Bartos et al. (Reference Encinas-Bartos, Aksamit and Haller2022) for a deeper investigation of those schemes in a similar context. This trajectory concentration is well below the limits of the previous sparse trajectory studies of Lekien & Ross (Reference Lekien and Ross2010) and Mowlavi et al. (Reference Mowlavi, Serra, Maiorino and Mahadevan2022).

From this initial investigation it appears that $|\boldsymbol {v}|$ may be the best option for identifying structures with low-resolution data. The dominant flow feature suggested by $|\boldsymbol {v}|$ is, in fact, misleading, even at full resolution. Figure 3(a) shows an overlay of $|\boldsymbol {v}|$ contours with the FTLE field in greyscale. The distribution of $|\boldsymbol {v}|$ suggests two distinct peaks with a distribution minimum, $s_\delta$, separating the central and outer regions. The fast moving central core is highlighted as the region inside the red $|\boldsymbol {v}|=s_\delta$ contour in figure 3(a). Examination of figure 3(a) show that the jet and eddy structures suggested by the FTLE field have boundaries that are actually transverse to all of the velocity contours.

Figure 3.

Comparison of particle speed, FTLE and the evolution of fluid particles advected from inside a given velocity level set.

In figure 3(b) we plot the final position of fluid particles from only the velocity core ($|\boldsymbol {v}|\ge s_\delta$), coloured by their $|\boldsymbol {v}(t_0)|$ value, after 30 days of advection. Particles corresponding to the red boundary have been advected as well. It is clear that the velocity-based jet is not actually a coherent structure as the proposed feature has been stretched inside and around the eddies, as well as down the eastward jet, with initial particle speed giving no indication of what feature a particle should be attached to.

For comparison, we show the evolution of the ${\mathrm {r}\alpha }$ structures, those that were similarly resilient to downsampling. In figure 4(a) we plot ${\mathrm {r}\alpha }$, with rotationally coherent structure boundaries identified as the outermost closed convex contours of ${\mathrm {r}\alpha }$, as has been previously used in other elliptic LCS methods (Haller et al. Reference Haller, Hadjighasem, Farazmand and Huhn2016; Encinas-Bartos et al. Reference Encinas-Bartos, Aksamit and Haller2022). Figure 4(c) shows the location of the same ${\mathrm {r}\alpha }$-coloured fluid particles after 30 days of advection. The ${\mathrm {r}\alpha }$ eddy boundaries have minimally deformed after being advected downstream more than 20 times their diameter downstream. This time also corresponds with approximately six full rotations of the eddy, ample time for inaccurate eddy boundaries to experience filamentation. This strong coherence can be attributed to the mathematical definition of ${\mathrm {r}\alpha }$ as a measure of relative fluid rotation and suggest we have effectively and objectively identified the true rotationally coherent structures in our flow field.

Figure 4.

Outermost closed convex ${\mathrm {r}\alpha }$ level-set contours overlaid on the ${\mathrm {r}\alpha }$ field. Advected position of ${\mathrm {r}\alpha }$ contours and surrounding fluid, coloured by initial ${\mathrm {r}\alpha }$ values.

The ability of ${\mathrm {r}\alpha }$ to identify these same eddies at a random lower data concentration ($3.5\ell ^{-2}$) is tested in figure 4(b,d). We again identify eddies as the outermost closed convex contours of ${\mathrm {r}\alpha }$ at $t=0$ (figure 4b) and advect them for 30 days. Comparing figures 4(a) and 4(b), we see the leftmost eddy boundary extending into the central jet and differences arising in the exact eddy boundary location. Testing these boundaries under advection, we find the red boundary particles have remain primarily inside the eddies we sought to identify (figure 4d), but with a much higher degree of filamentation.

The ability to identify structure boundaries is as much of a test of the interpolation scheme used as of ${\mathrm {r}\alpha }$ to maintain meaningful values at lower resolution. When there are fewer than one trajectory per coherent structure, we will inevitably misrepresent some flow features in interpolated fields and ${\mathrm {r}\alpha }$ and ${\rm r}\overline {\mathrm {TSE}}$ are more reliably interpreted as point measurements of local behaviour. As mentioned previously, we investigate interpolation scheme impacts in the Appendix. Here, we investigate the ability of low-concentration data to provide accurate reference frame values $\bar {\boldsymbol {v}}$, $\overline {\boldsymbol {\omega }_e}$ and $\bar {\boldsymbol {x}}$ and meaningful pointwise ${\mathrm {r}\alpha }$, independent of the interpolation scheme, with the following experiment.

Starting with full-resolution ($3600\ell ^{-2}$) Bickley jet trajectories, we first calculate ground-truth values for ${\mathrm {r}\alpha }(\boldsymbol {x}_0)_{{full}}$. We then randomly subsample the full resolution trajectory $10^5$ times at a given concentration. For each selection of trajectories, we recalculate ${\mathrm {r}\alpha }(\boldsymbol {x}_0)_{{low}}$ using the lower-resolution data, and calculate the correlation coefficient ($R^2$) of ${\mathrm {r}\alpha }(\boldsymbol {x}_0)_{{low}}$ and ${\mathrm {r}\alpha }(\boldsymbol {x}_0)_{{full}}$ at the corresponding locations. We then obtain a distribution of $R^2$ values that details how closely we can approximate the full-resolution measurements for a large number of potential trajectory orientations. We then repeat this process at progressively lower concentrations and witness the rate at which the accuracy degrades.

The findings from this experiment are summarised in table 2. The mean value of $R^2$ slowly decreases with increasing sparsity but stays above 0.99, with a standard deviation $\le$0.008 down to $16\ell ^{-2}$. The mean of $R^2$ reduces more dramatically around the extremely sparse sampling threshold of $1\ell ^{-2}$. This is precisely the trajectory concentration at which we suggest using time-averaged $\boldsymbol {v}_{{ref}}$, $\boldsymbol {\omega }_{{ref}}$ and $\boldsymbol {x}_{{ref}}$ values. The impact of this choice is seen in the final column where mean $R^2$ values nearly return to the $4\ell ^{-2}$ value. This suggests a significant improvement in pointwise ${\mathrm {r}\alpha }$ accuracy, and meaningful rotation diagnostics for extremely sparse data.

Table 2.

Bickley jet downsampling effect on ${\mathrm {r}\alpha }$ accuracy. Mean ($\mu$) and standard deviation ($\sigma$) of $R^2$ values correlating downsampled ${\mathrm {r}\alpha }$ with full-resolution calculations. Raw ${\mathrm {r}\alpha }$ values maintain a strong correlation until the extremely sparse sampling threshold (1$\ell ^{-2})$, at which point steady time-averaged reference values (‘ref’, § 2.5) allow for a significant improvement in accuracy again.

3.1.2. AVISO

We now consider the identification of elliptic Lagrangian flow structures from 2-D ocean satellite altimetry data provided by AVISO which has been the focus of several coherent structure studies (see, e.g. Haller et al. Reference Haller, Hadjighasem, Farazmand and Huhn2016, Reference Haller, Aksamit and Bartos2021). The zonal and meridional component $\boldsymbol {v} = (\mathrm {v}_1, \mathrm {v}_2)$ of the ocean currents are derived from the sea-surface height profile

(3.2)\begin{gather} f\mathrm{v}_2(\boldsymbol{x},t) = \dfrac{1}{\rho} \dfrac{\partial}{\partial x}p(\boldsymbol{x},t), \end{gather}

(3.3)\begin{gather}f\mathrm{v}_1(\boldsymbol{x},t) = \dfrac{1}{\rho} \dfrac{\partial}{\partial y}p(\boldsymbol{x},t), \end{gather}

(3.4)\begin{gather}-g = \dfrac{1}{\rho}\dfrac{\partial}{\partial z}p(\boldsymbol{x},t), \end{gather}

where $p(\boldsymbol {x},t)$ is the pressure, $\rho$ is the fluid density, $f$ is the Coriolis parameter and $g$ is the Earth's acceleration. The daily-gridded velocity data are freely available from the Copernicus Marine Environment Monitoring Service. Although this is an observational product, the processed nature of the product prevents many of the complications that we will see with our experimental examples in the following section. Our analysis focuses on the North Atlantic Gulf Stream between longitudes $70\,^{\circ }$W and $55\,^{\circ }$W and latitudes $30\,^{\circ }$N and $45\,^{\circ }$N, spanning September and October 2006. We start with an initial grid of trajectories consisting of $150 \times 150$ points. Using a characteristic length scale of mesoscale eddies ($\ell =100$ km), this density roughly corresponds to $100\ell ^{-2}$.

Elsewhere in the Atlantic, Haller et al. (Reference Haller, Aksamit and Bartos2021) and Encinas-Bartos et al. (Reference Encinas-Bartos, Aksamit and Haller2022) have shown that the AVISO velocity field is slowly varying, with a vanishing spatially averaged vorticity. Furthermore, these authors showed that $\overline {\mathrm {TRA}}$ effectively identifies ocean eddies, and outperforms other sparse trajectory rotation diagnostics. As such, in figure 5 we compare the ${\mathrm {r}\alpha }$ (figure 5a) with $\overline {\mathrm {TRA}}$ (figure 5b), using LAVD as the ground truth of material rotation (figure 5c). The white dashed box highlights the meanders of the Gulf Stream. Three eddy boundaries are identified as outmost closed convex contours of LAVD in figure 5(a–c).

Figure 5.

Overview of rotation diagnostics for the AVISO data in the Gulf Stream at full resolution ($100 \ell ^{-2}$). The white contours are the outermost convex contours of the LAVD field. The dashed white box highlights the meanders of the Gulf Stream.

When compared with LAVD, $\overline {\mathrm {TRA}}$ faithfully approximates the LAVD field, identifying the meanders of the Gulf Stream, high rotation in the adjacent eddies and relatively low rotation away from these central features. At this trajectory density, figure 5(a) reveals that ${\mathrm {r}\alpha }$ can match the ability of $\overline {\mathrm {TRA}}$ to visualise the correct flow structure topology. This stems from the similarity in $\overline {\mathrm {TRA}}$ and ${\mathrm {r}\alpha }$ definitions for small bulk advection and rotation, as highlighted in § 2.5.

We compare the sparse data performance of ${\mathrm {r}\alpha }$ against $\overline {\mathrm {TRA}}$ in figure 6 by progressively and randomly subsampling using the same methodology as in § 3.1.1. We again interpolate back to full resolution using a radial basis function for each successive subsampling. We also include the major eddy boundaries from figure 5.

Figure 6.

Comparison for the AVISO data between objective rotational sparse data diagnostic ${\mathrm {r}\alpha }$, $\mathrm {r}\alpha _{{ref}}$ and its quasi-objective single-trajectory counterpart $\overline {\mathrm {TRA}}$ under varying trajectory densities ($10 \ell ^{-2}, 1 \ell ^{-2}, 0.1 \ell ^{-2}$). The white markers denote the trajectory endpoints.

At a trajectory concentration of $10\ell ^{-2}$, $\overline {\mathrm {TRA}}$ (figure 6a) and ${\mathrm {r}\alpha }$ (figure 6d) appear qualitatively similar, with both suggesting relatively high rotation in the Gulf Stream meanders, and three distinct elliptic rotational maxima aligned with the eddy boundaries. For a trajectory density of $1 \ell ^{-2}$, ${\mathrm {r}\alpha }$ (figure 6e) is qualitatively indistinguishable from $\overline {\mathrm {TRA}}$ (figure 6b) with two of the three mesoscale LAVD eddies still visible as local maxima in both diagnostic fields. For this region of the ocean, we often have ocean drifter observations [e.g. Global Drifter Program (GDP)] at concentrations between $1 \ell ^{-2}$ and $0.1 \ell ^{-2}$, as discussed in § 3.2.1.

As we pass the extremely sparse sampling threshold to $0.1 \ell ^{-2}$, the leftmost eddy is still identifiable as a $\overline {\mathrm {TRA}}$ and ${\mathrm {r}\alpha }$ maximum (figures 6c and 6f, respectively) due to the presence of the single trajectory originating inside its boundary (white cross). The remaining eddies are no longer visible for either metric due to a lack of data, but $\overline {\mathrm {TRA}}$ still highlights a relatively higher rotation in the Gulf Stream. The local maximum $\overline {\mathrm {TRA}}$ within the dashed white box may be classified incorrectly as an eddy, but it actually corresponds to a region of high material rotation induced by shear in the Gulf Stream. In contrast, ${\mathrm {r}\alpha }$ has created multiple local maxima near the coast of Nova Scotia not present at higher concentrations.

At this extremely sparse resolution, we again replace $\boldsymbol {\bar {\omega }}(t)$, $\boldsymbol {\bar {v}}(t)$, and $\boldsymbol {\bar {x}}(t)$ with their time-averaged steady approximations, $\boldsymbol {v}_{{ref}}, \boldsymbol {\omega }_{{ref}}, \boldsymbol {x}_{{ref}}$ from (2.17)–(2.19). As with the Bickley jet statistics, this improves the accuracy of ${\mathrm {r}\alpha }$ with ${\mathrm {r}\alpha }_{{ref}}$ displayed in figure 6(i). With the steady reference frame, ${\mathrm {r}\alpha }_{{ref}}$ and $\overline {\mathrm {TRA}}$ again show nearly indistinguishable flow topologies. The coastal local maxima disappear, and the local maximum appears in the shear dominated meanders of the Gulf Stream. Furthermore, if we apply this steady reference frame change to higher concentrations (figure 6g and h), we obtain flow structure visualisations comparable to ${\mathrm {r}\alpha }$. This suggests that obtaining meaningful results using the extremely sparse steady reference-frame approach is resilient over a range of concentrations.

This effect stems from the effect of undersampling on bulk rotation and translation estimates. In figure 7, we display the $\boldsymbol {\bar {\omega }}(t)$ and $\boldsymbol {\bar {v}}(t)$ estimates for each trajectory concentration. At $100 \ell ^{-2}$ the mean flow properties are nearly constant over 15 days. This is consistent with the expectation from earlier studies that Eulerian features do not change significantly over such relatively small time scales in the ocean. As we progressively reduce the trajectory density, however, $\boldsymbol {\bar {\omega }}(t)$ and $\boldsymbol {\bar {v}}(t)$ start displaying strong oscillations in time. At extremely low trajectory densities ($<1 \ell ^{-2}$), the spatial averages are no longer representative of the underlying mean flow properties (dashed red curves in figure 7). These oscillations contribute unphysical fluctuations in the direction of $\boldsymbol {v}_{d}$ and hinder our ability to measure rotation with ${\mathrm {r}\alpha }$.

Figure 7.

Data-based estimates for bulk translation velocity $|\boldsymbol {\bar {v}}(t)|$ and rotation $\boldsymbol {\bar {\omega }}(t)$ under varying trajectory density ($100 \ell ^{-2}, 10 \ell ^{-2}, 1 \ell ^{-2}, 0.1 \ell ^{-2}$).

To further quantify the impact of sparse sampling on ${\mathrm {r}\alpha }$, we repeat the $R^2$ distribution analysis that was performed for the unsteady Bickley jet. Using AVISO surface current data to generate an initial baseline at $100\ell ^{-2}$, we run $10^4$ random subsamplings at a range of concentration levels. The summary statistics of this experiment can be found in table 3. A consistent decrease in the $R^2$ mean, and increase in $R^2$ standard deviation is found as the data concentration levels decrease, though the level of ${\mathrm {r}\alpha }$ accuracy remains quite high. Even at a low concentration level of $1\ell ^{-2}$, we obtain a mean $R^2$ value of 0.92.

Table 3.

Downsampling effect on ${\mathrm {r}\alpha }$ accuracy for AVISO Gulf Stream surface currents. Mean and standard deviation of $R^2$ values relating downsampled ${\mathrm {r}\alpha }$ with full-resolution calculations for $10^4$ random subsamplings at each resolution. Raw ${\mathrm {r}\alpha }$ values maintain a strong correlation until the extremely sparse threshold $(1 \ell ^{-2})$, at which point time-averaged reference values (ref) allow for a significant improvement in accuracy again.

The most precipitous drop in ${\mathrm {r}\alpha }$ accuracy is found from $1\ell ^{-2}$ to $0.1\ell ^{-2}$. One such example of this poor performance was highlighted in figure 6. We further validate the improvement in ${\mathrm {r}\alpha }$ accuracy at this extreme sparsity level in the last two columns of table 3. After transitioning to time-averaged reference frame values, there is a significant increase in accuracy, with the mean of $R^2$ jumping from 0.92 to 0.95 for $1\ell ^{-2}$, and from 0.45 to 0.89 for $0.1\ell ^{-2}$. With this simple reference-frame methodology we find $r\alpha$ is a promising Lagrangian flow rotation diagnostic for extremely sparse data.

3.1.3. Stationary concentrated vortex

Our third example considers the identification of elliptic coherent structures inside the stationary concentrated vortex model (SCVM) of Onishchenko et al. (Reference Onishchenko, Fedun, Horton, Pokhotelov, Astafieva, Skirvin and Verth2021). The SCVM represents finite-size cyclones in the Earth's atmosphere as steady axially symmetric solutions to the Euler equations. Both the vertical and radial extent of the vortex is controlled by model parameters. The SCVM fundamentally differs from the previous flows as it is a highly rotational flow with no bulk translation. Our investigation of the SCVM helps confirm that ${\mathrm {r}\alpha }$ and ${\rm r}\overline {\mathrm {TSE}}$ are robust diagnostics for approximating the shape and location of LCS and are not tailored to a specific kind of flow experiment (e.g. both high and low swirl number; Gupta, Lilley & Syred Reference Gupta, Lilley and Syred1984).

The SCVM consists of inner and outer regions, an internal upwards motion and external downwards motion. A distinct central torus was also identified by Aksamit (Reference Aksamit2023) as a LCS. The cyclone is further divided vertically into a centripetal and centrifugal flow with upward-moving fluid recirculated from the top of the vortex back down.

The SCVM was originally derived as an axisymmetric flow from a stream function $\varPsi (r,v)$ and can be written in cylindrical coordinates as

(3.5)\begin{gather} \mathrm{v}_r ={-}\mathrm{v}_0 \frac{r}{L} \left( 1-\frac{z}{L} \right)\exp\left({-\frac{z}{L}-\frac{r^2}{r_0^2}}\right), \end{gather}

(3.6)\begin{gather}\mathrm{v}_\theta ={\pm} \mathrm{v}_{\theta_0}\frac{z}{L}\exp\left({-\frac{z}{L}-\frac{r^2}{r_0^2}}\right), \end{gather}

(3.7)\begin{gather}\mathrm{v}_z = 2\mathrm{v}_0 \frac{z}{L} \left( 1-\frac{r^2}{r_0^2} \right)\exp\left({-\frac{z}{L}-\frac{r^2}{r_0^2}}\right), \end{gather}

where $r_0$ is the radial extent separating inner and outer vortex structures, $L$ is the height of maximum vertical velocity at $r=0$ and $\mathrm {v}_0$ and $\mathrm {v}_{\theta _0}$ are characteristic velocities. For our example, we use $\mathrm {v}_0 = 1$, $L = 10$, $r_0 = 1$ and $\mathrm {v}_{\theta _0}=9$. We consider our characteristic length of the flow to be equal to the separation of inner and outer motions, $\ell =r_0=1$, and the approximate radius of the main central torus.

In figure 8 we compare the objective rotational diagnostics ${\mathrm {r}\alpha }$ and $\mathrm {LAVD}$ advected for $\Delta t = 500$ with the vorticity magnitude, $|\boldsymbol {\omega }|$ and the widely used lambda-2 criterion of Jeong & Hussain (Reference Jeong and Hussain1995). If $\boldsymbol {J}$ is the velocity gradient tensor, and $\boldsymbol {S}$ and $\boldsymbol {\varOmega }$ are the symmetric and antisymmetric parts, the lambda-2 criterion of Jeong & Hussain (Reference Jeong and Hussain1995) suggests that regions with negative intermediate eigenvalues of $\boldsymbol {S}^2 + \boldsymbol {\varOmega }^2$ (labelled $\lambda _2$) are part of a vortex core. Stationarity of the flow suggests $\overline {\mathrm {TRA}}$ would also be a suitable comparison, but in light of the strong performance of $r\alpha$ against $\overline {\mathrm {TRA}}$ in the previous two examples, we have left it out for compactness of presentation.

Figure 8.

Comparison of rotation diagnostics ${\mathrm {r}\alpha }$, $\mathrm {LAVD}$, $\lambda _2$ and $|\boldsymbol {\omega }|$ for the SCVM in the $y=0$ plane. White contour lines are contours of the SCVM stream function $\varPsi$, indicating circulation patterns in the model. Note the agreement between ${\mathrm {r}\alpha }$, $\mathrm {LAVD}$, ${\rm r}\overline {\mathrm {TSE}}$ and $\mathrm {FTLE}$ features with the stream function structures, and the dissimilarity of $|\boldsymbol {\omega }|$ and $\lambda _2$ features.

From initial high spatial resolution ($2000\ell ^{-3}$) calculations, a 2-D slice along the $y=0$ plane is shown for each field in figure 8. We overlay these plots with the level-set contours of the SCVM stream function $\varPsi (r,v)$. As this is a stationary axisymmetric flow, $\varPsi$-contours provide us a ground truth of transport barrier structure. A qualitative contour comparison of the sparse metrics (${\mathrm {r}\alpha }$ and ${\rm r}\overline {\mathrm {TSE}}$) with $\varPsi$ show a comparable level of structure agreement as one can find from LAVD and FTLE structures. We find that ${\mathrm {r}\alpha }$ and $\mathrm {LAVD}$ contours suggest very similar rotational behaviour, with both indicating the same location of a central torus vortex surrounded by recirculating flow. Both r$\overline {\mathrm {TSE}}$ and FTLE show strong spiraling into the centre of the vortices identified by ${\mathrm {r}\alpha }$ and $\mathrm {LAVD}$. One notable difference between sparse and spatial-gradient-reliant diagnostics exists in the centre of the toroidal vortex core in ${\rm r}\overline {\mathrm {TSE}}$ which again have unnaturally low values when compared with FTLE, as we also saw in the unsteady Bickley jet.

We see a stark difference, however, when comparing $\varPsi$ with $|\boldsymbol {\omega }|$ and $\lambda _2$ fields. Contours of $|\boldsymbol {\omega }|$ and $\lambda _2$ are often orthogonal to $\varPsi$-contours, suggesting a disagreement in structure identification. We further investigate these qualitative findings by testing the ability of our rotation diagnostics to identify the dominant elliptic coherent structures in this vortex flow. Potential vortex boundaries are identified as the set of outermost closed convex contours for ${\mathrm {r}\alpha }$, LAVD and $|\boldsymbol {\omega }|$. Following the derivation of the $\lambda _2$ criterion, we use $\lambda _2<0$ as the designation of strong rotational motion (Jeong & Hussain Reference Jeong and Hussain1995). Fluid parcels from inside the boundaries of these features are then advected for twice the integration time used for the Lagrangian diagnostics, $\Delta t = 1000$. This approach allows us to test both the Lagrangian and Eulerian diagnostics to identify features that maintain their coherence beyond the time horizon of data used for the initial boundary calculation.

These fluid particle paths are shown in figure 9 with trajectories coloured by the initial diagnostic field value. We include the corresponding figure 8 diagnostic field on the $y=0$ plane as a reference. The vortex core identified by ${\mathrm {r}\alpha }$ and $\mathrm {LAVD}$ maintain coherence over time, with seemingly no deviation of fluid from the domain identified as the central torus vortex (figure 9a,b). Given that ${\mathrm {r}\alpha }$ does not rely on the same 3-D velocity gradient data as $\mathrm {LAVD}$, ${\mathrm {r}\alpha }$ shows exceptional performance in identifying this LCS. Fluid initially inside $\lambda _2$ and $|\boldsymbol {\omega }|$ structures have instead maintained no coherence, and are sufficiently dispersed in the flow. The high $|\boldsymbol {\omega }|$ (figure 9d) fluid parcels are now mixed throughout the entire inner SCVM ‘beehive’ domain, inhabiting regions with both high and low $|\boldsymbol {\omega }|$ values. The $\lambda _2\!<\!0$ fluid is recirculating in both the inner and outer regions of the SCVM (figure 9c), with the $\lambda _2$-estimated vortex core parcels inhabiting the domain of both negative and positive $\lambda _2$ values. This is distinct from LAVD and ${\mathrm {r}\alpha }$ vortex trajectories that stay confined to the regions of high LAVD and ${\mathrm {r}\alpha }$, respectively.

Figure 9.

Trajectories corresponding to coherent structures identified by the four diagnostics, calculated for the same integration times.

The ability of ${\mathrm {r}\alpha }$ to identify these same LCS at increasing sparsity is shown in figure 10. As our only 3-D numerical example, we utilise a random set of initial particle trajectories in three dimensions and extract isosurfaces from a reconstructed ${\mathrm {r}\alpha }$ field using natural neighbour interpolation. We utilise the same ${\mathrm {r}\alpha }$ value for 3-D level-set extraction as the outermost closed convex contour used to test vortex boundaries (figure 9). Plots of reconstructed ${\mathrm {r}\alpha }$ fields on the $y=0$ plane (figure 10a–c) reveal that small deviations from this range would not significantly change the topology of the surfaces in figure 10(d–f). Note the difference of $z$-value ranges between figures 10(a–c) and 10(d–f).

Figure 10.

(a–c) Values of ${\mathrm {r}\alpha }_0^{500}$ on the $y=0$ plane reconstructed from randomly selected trajectories at decreasing particle density. The inner vortex is extracted at each concentration as the closed contours corresponding to the level-set value from figure 8. (d–f) Comparison of the inner vortex surface in ${\mathrm {r}\alpha }_0^{500}$ fields reconstructed from randomly selected trajectories at decreasing particle density. The inner vortex is extracted in each scenario as closed contours corresponding to the outermost closed convex contour in figure 9.

In figure 10, ${\mathrm {r}\alpha }$ at $10\ell ^{-3}$ looks only moderately different from $2000\ell ^{-3}$. Two distinct intersections of the central torus with the $y=0$ plane are still evident in white, and the annular shape is clearly visible in the associated isosurface. At $1\ell ^{-3}$, the ‘doughnut hole’ above $x=0$ becomes less obvious, but the horizontal and vertical extent of the central torus is still accurate. For $0.1\ell ^{-3}$, we use ((2.17)–(2.18)) to implement the same bulk time-averaged reference values for $\boldsymbol {v}_{{ref}}, \boldsymbol {\omega }_{{ref}}, \boldsymbol {x}_{{ref}}$. At this extremely sparse sampling, the correct extent of strong central rotation is still evident, as is a suggestion of the inner beehive feature, but the vortex surface is now sufficiently deformed that it no long resembles a torus.

3.2.1. GDP database

For our first experimental data set we consider a set of drifters from the GDP (Lumpkin & Centurioni Reference Lumpkin and Centurioni2019) located in the North Atlantic Ocean. The GDP data set contains more than 40 000 drifters spanning the last four decades. Roughly 1200 drifters are currently active worldwide and report their position every 6 hours (Lumpkin & Pazos Reference Lumpkin and Pazos2007).

Specifically, we focus here on a subset of drifters active over 15 days from 20 September to 4 October 2006 in the Gulf Stream (figure 11) during the same time period of analysis as in § 3.1.2. In figure 11, we reconstruct diagnostic fields from the true GDP tracks with a drifter density of $0.63 \ell ^{-2}$, for $\ell = 100\,\mathrm {km}$. We are thus in an extremely sparse data setting ($< 1 \ell ^{-2}$) and we again compute ${\mathrm {r}\alpha }$ by replacing the time-dependent mean flow properties with time-averaged $\boldsymbol {\omega }_{{ref}}, \boldsymbol {v}_{{ref}}$ and $\boldsymbol {x}_{{ref}}$ from the drifter data.

Figure 11.

Drifter-based $\overline {\mathrm {TRA}}$ (a), ${\mathrm {r}\alpha }$ (b), $\overline {\mathrm {KE}}$ (c) and LAVD (d) in the North Atlantic Gulf Stream. The black crosses indicate drifter positions.

Following the sparse eddy identification study by Encinas-Bartos et al. (Reference Encinas-Bartos, Aksamit and Haller2022), we compare two sparse drifter-based rotation diagnostics, $\mathrm {r}\alpha _{{ref}}$ in figure 11(a) and $\overline {\mathrm {TRA}}$ in figure 11(b), with commonly used trajectory-averaged kinetic energy (KE; figure 11c) (e.g. Lumpkin Reference Lumpkin2016), and the full-resolution LAVD field using AVISO data from § 3.1.2. Though the AVISO current data may not be representative of the full spatial scale of flow features affecting the drifters, we include LAVD as a mesoscale comparison for general flow behaviour.

We do not compare with a direct ‘looping’ count as a user-defined threshold is inevitably required and differentiating between looping and non-looping trajectories is user-dependent and non-trivial (Encinas-Bartos et al. Reference Encinas-Bartos, Aksamit and Haller2022). Instead, we extract looping segments for the reader to see from the drifter trajectories as proposed by the frame-dependent algorithm from Lumpkin (Reference Lumpkin2016). The white trajectory segments indicate cyclonic loopers which frequently arise in connection with ocean eddies (see, e.g. Griffa, Lumpkin & Veneziani Reference Griffa, Lumpkin and Veneziani2008; Dong et al. Reference Dong, Liu, Lumpkin, Lankhorst, Chen, McWilliams and Guan2011). The full black curves in figure 11(a,b,d) denote the approximate vortex boundaries respectively obtained from drifter-based $\overline {\mathrm {TRA}}$ and ${\mathrm {r}\alpha }_{ref}$, and AVISO-based LAVD, respectively.

Both $\mathrm {r}\alpha _{{ref}}$ and $\overline {\mathrm {TRA}}$ are able to detect rotational features just to the south of the Gulf Stream from the extremely sparse GDP data. Cyclonic loopers are also observed in the regions surrounding the principal local maxima of $\mathrm {r}\alpha _{{ref}}$ and $\overline {\mathrm {TRA}}$. The leftmost eddy contains two looping drifters that display high rotational values (figure 11a,b). In contrast, $\overline {\mathrm {KE}}$ fails to clearly identify the leftmost eddy since it displays sharp gradients between the two drifters inside the eddy (figure 11c). Both $\mathrm {r}\alpha _{{ref}}$ and $\overline {\mathrm {TRA}}$ display a weaker rotational feature at roughly $32\,^{\circ }$N latitude and $64^{\circ }$W longitude. No loopers are identified in this region, but the yellow trajectories in the high ${\mathrm {r}\alpha }$ and $\overline {\mathrm {TRA}}$ zones show greater circulating behaviour with respect to the surrounding green trajectories (figure 11, bottom insets), suggesting stronger rotational motion in the surface currents. There is no feature that resembles any vortical structure in the $\overline {\mathrm {KE}}$ field in this region.

In this ocean buoy experiment, we were able to identify eddies adjacent to the Gulf Stream in both $\overline {\mathrm {TRA}}$ and ${\mathrm {r}\alpha }$ fields that were not identified using existing sparse data diagnostics. In contrast with the $\overline {\mathrm {TRA}}$, which requires additional assumptions on the frame of reference, the ${\mathrm {r}\alpha }$ returns structures that are valid for all observers. This suggests that ${\mathrm {r}\alpha }$ may be beneficial for systematic eddy analysis of the GDP data set that can complement remote-sensing based tools that lack submesoscale motion resolutions. In addition, ${\mathrm {r}\alpha }$ shows great promise as a buoy diagnostic for focused buoy deployments (e.g. LASER and GLAD) where buoy concentrations higher than $1\ell ^{-2}$ will result in more refined descriptions of the flow.

3.2.2. Large-scale LPT experiment

Our last example involves a novel, large-volume 3-D LPT velocimetry data set from an industrial-scale wind tunnel (Hou et al. Reference Hou, Kaiser, Sciacchitano and Rival2021). Utilising glare-point spacing on large soap bubbles, Hou et al. (Reference Hou, Kaiser, Sciacchitano and Rival2021) were able to investigate the vortical wake behind a tractor–trailer model at a $9^{\circ }$ yaw angle. For details on the methodology, see Kaiser & Rival (Reference Kaiser and Rival2023). We compare two wind-tunnel experiments, one with the tractor–trailer model obstructing the flow, and a uniform flow without the tractor–trailer.

Each experiment provides around 300 s of trajectory data, sampled at 150 Hz with an average windspeed of 10 m s$^{-1}$. The bubble generator and camera setup produce time-varying particle concentrations at the very lower end of our computation spectrum. For measurements of the tractor–trailer wake, trajectory concentrations oscillate between 0.2 and 1.5$\ell ^{-3}$ with a mean concentration of 0.6$\ell ^{-3}$. The uniform flow measurements have a slightly higher resolution, ranging from 0.5 to 3$\ell ^{-3}$ and a mean concentration of 1.3$\ell ^{-3}$. To avoid significant oscillations in relative Lagrangian velocities, we again utilise the time-averaged reference frame determined by the bubble trajectories. Furthermore, as we are using trajectory data that quickly transits through the observation volume, we investigate 10-frame increments and require all trajectories to be present for every frame. At these extremely low concentrations, a single particle that disappears and reappears during the integration window may have a substantial effect on the deformation velocity of other fluid particles. A 10-frame integration time was chosen as a balance between the length of integration time and the number of tracks available to calculate meaningful mean values.

At the low trajectory concentration available in this turbulence experiment, meaningful time-resolved spatial derivatives are unattainable. Instead of relying on a spatial interpolant to visualise strong regions of rotation, we now demonstrate a complementary statistical approach using this sparse data set. In figure 12 we show the discrete probability histograms for common LPT flow diagnostics in addition to ${\mathrm {r}\alpha }$ analysis. The distribution of the particle speed over the entire 300 s (figure 12b) shows a slight speed up of the flow for the trailer wake, but the bubble trajectories travel for a shorter distance on average (figure 12c). The fluctuating streamwise components of bubble velocities were also nearly identical for the two flows.

Figure 12.

Comparison of discrete probability histograms of frame-dependent metrics for the truck wake vortex and object-free uniform flow experiments at $\sim 1\ell ^{-3}$.

Sensor development for large-scale LPT is still an active field of research. Given the uncertainty in signal noise from nascent methodologies, such as the glare-point method used to generate for these data, flow diagnostics that maintain their utility for unknown noise levels can help to expedite the development of large-scale measurement tools. With the available information from the industrial-scale wind tunnel, none of the diagnostics in figure 12(b–d) present a clear picture of the wake dynamics. In fact, $\boldsymbol {v}_x'$ cannot distinguish between the two flows.

To complement these approaches, we compute 10-frame ${\mathrm {r}\alpha }$ values beginning at every frame available. This gives us a distribution of ${\mathrm {r}\alpha }_t^{t+10}(\boldsymbol {x}_0)$ for a wide range of $t$ and $\boldsymbol {x}_0$. From a set of trajectories at each time $t$, we then calculate the median of ${\mathrm {r}\alpha }$ ($\widetilde {{\mathrm {r}\alpha }}$) to quantify the general behaviour of the flow in each time window. The discrete probability histogram of $\widetilde {{\mathrm {r}\alpha }}$ for the two flows is shown in figure 12(a). The empty tunnel data are concentrated around a low peak, suggesting a relatively uniform flow domain and an accumulation of small turbulent oscillations, but no major rotational features. In contrast, the trailer wake data are distinguished by a larger range of values, with a heavy tail towards considerable rotation and high peak $\widetilde {{\mathrm {r}\alpha }}$. The physical meaning of this data suggests a more complex flow, with a higher degree of rotation, in this case because of the vortical wake. Even at this low concentration, and with this experimental methodology, we are able to identify distinguished rotating features in the flow.

We can begin to distinguish the spatial features of the flow by analysing the long-time trends in ${\mathrm {r}\alpha }$. In figure 13 we project all initial particle locations to the $y$–$z$ plane, and bin them in a $100\times 100$ grid. We again calculate $\widetilde {{\mathrm {r}\alpha }}$ for each bin and plot this in figure 13(a,b). The trailer wake in figure 13(a) has high rotation near its centre, and relatively higher ${\mathrm {r}\alpha }$ distributed throughout the flow domain when compared with the quiescent centre and lower values in figure 13(b). We quantify the spatial dependence of rotation on flow geometry by binning ${\mathrm {r}\alpha }$ in 10 bins along the radial distance from the centre of the flow. We show this radial axis as a dot and line in figure 13(a,b). Again, the median ${\mathrm {r}\alpha }$ values from these bins shows the trailer wake generates a much greater degree of rotation on the fluid, starting from a minimum where little data are available, with a local maximum in the vortex core, then gradually increasing as we move away in the flow. At nearly every location, the empty tunnel generates limited fluid rotation with lower ${\mathrm {r}\alpha }$ values.

Figure 13.

(a,b) Binned median value of ${\mathrm {r}\alpha }_0^t(\boldsymbol {x}_0)$ after projecting in the streamwise direction onto the $y$–$z$ plane. The trailer wake has high rotation in its core, and relatively higher $\widetilde {{\mathrm {r}\alpha }}_0^t$ values distributed across the flow. The empty tunnel, in contrast has a low rotational core, with higher $\widetilde {{\mathrm {r}\alpha }}_0^t$ appear only on the edges of the measurement domain. This is confirmed by binning the ${\mathrm {r}\alpha }_0^t$ values by radial distance from the flow centre, and calculating median values ($\widetilde {{\mathrm {r}\alpha }}_0^t(r)$, panel c). Note the radial centres and axes are marked by blue dots and lines in panels (a,b).

The results here provide a good baseline for experimental LPT data that were not designed with ${\mathrm {r}\alpha }$ analysis in mind. As future experiments may be designed that consider ${\mathrm {r}\alpha }$ requirements, further improvement on the identification and analysis of coherent rotational flow features is expected. Furthermore, given the connections between ${\mathrm {r}\alpha }$, curvature (e.g. § 2.4) and rotation times, concerted future effort may help connect classical notions of geometry, eddy turnover times, spectral analysis and sparse trajectory data in an objective manner, in addition to designating strongly and weakly rotating regions of a flow.