2 Set-up

We consider an unsteady velocity field $\boldsymbol{v}(\boldsymbol{x},t)$ , defined on a possibly time-dependent spatial domain $U(t)\subset \mathbb{R}^{3}$ over a finite time interval $[t_{0},t_{1}].$ We assume that $U(t)$ is invariant under the fluid flow generated by the velocity field (cf. (2.7) below). Thus, $U(t)$ is either a physical domain with an impermeable boundary or a material domain formed by a set of evolving trajectories of $\boldsymbol{v}(\boldsymbol{x},t)$ .

We write the velocity gradient $\boldsymbol{{\rm\nabla}}\boldsymbol{v}$ as

(2.1) $$\begin{eqnarray}\boldsymbol{{\rm\nabla}}\boldsymbol{v}(\boldsymbol{x},t)=\unicode[STIX]{x1D63F}(\boldsymbol{x},t)+\unicode[STIX]{x1D652}(\boldsymbol{x},t),\end{eqnarray}$$

where $\unicode[STIX]{x1D63F}=(\boldsymbol{{\rm\nabla}}\boldsymbol{v}+[\boldsymbol{{\rm\nabla}}\boldsymbol{v}]^{\text{T}})/2$ is the rate-of-strain tensor and $\unicode[STIX]{x1D652}=(\boldsymbol{{\rm\nabla}}\boldsymbol{v}-[\boldsymbol{{\rm\nabla}}\boldsymbol{v}]^{\text{T}})/2$ is the spin tensor. We recall that the vorticity ${\it\omega}=\boldsymbol{{\rm\nabla}}\times \boldsymbol{v}$ of the fluid satisfies

(2.2) $$\begin{eqnarray}\unicode[STIX]{x1D652}\boldsymbol{e}=-{\textstyle \frac{1}{2}}{\bf\omega}\times \boldsymbol{e},\quad \forall \boldsymbol{e}\in \mathbb{R}^{3}.\end{eqnarray}$$

We will also use the instantaneous spatial mean $\bar{{\bf\omega}}$ of the vorticity over $U(t)$ , defined as

(2.3) $$\begin{eqnarray}\bar{{\bf\omega}}(t)=\frac{\displaystyle \int _{U(t)}{\bf\omega}(\boldsymbol{x},t)\,\text{d}V}{\text{vol}\,(U(t))},\end{eqnarray}$$

where $\text{vol}(\cdot )$ denotes the volume for three-dimensional flows and the area for two-dimensional flows. Accordingly, $\text{d}V$ refers to the volume or area element respectively in $U(t)$ .

Under general observer changes of the form (1.1), the spin tensor and the vorticity in the new $y$ coordinate frame take the forms

(2.4a,b ) $$\begin{eqnarray}\tilde{\unicode[STIX]{x1D652}}(\boldsymbol{y},t)=\unicode[STIX]{x1D64C}^{\text{T}}(t)\unicode[STIX]{x1D652}(\boldsymbol{x},t)\unicode[STIX]{x1D64C}(t)-\unicode[STIX]{x1D64C}^{\text{T}}(t)\dot{\unicode[STIX]{x1D64C}}(t),\quad \tilde{{\bf\omega}}(\boldsymbol{y},t)=\unicode[STIX]{x1D64C}^{\text{T}}(t){\bf\omega}(\boldsymbol{x},t)+\dot{\boldsymbol{q}}(t),\end{eqnarray}$$

with the vector $\dot{\boldsymbol{q}}$ defined uniquely by the relation $\dot{\unicode[STIX]{x1D64C}}\unicode[STIX]{x1D64C}^{\text{T}}\boldsymbol{e}=(\dot{\boldsymbol{q}}\times \boldsymbol{e})/2$ for all $e\in \mathbb{R}^{3}$ (see, e.g., Truesdell & Rajagopal Reference Truesdell and Rajagopal2009). Formula (2.4) shows that the spin tensor and the vorticity are not objective quantities: the eigenvalues and eigenvectors of $\unicode[STIX]{x1D652}$ change in rotating frames, and so do the direction and the magnitude of ${\bf\omega}$ . Thus, neither $\unicode[STIX]{x1D652}$ nor ${\bf\omega}$ is, by itself, suitable for defining distinguished material sets in the flow. This is because material sets are tied to evolving fluid particles without any reference to coordinates, and hence are inherently frame-invariant. More generally, statements about the material response of a moving continuum cannot depend on the observer and hence should be objective (Gurtin Reference Gurtin1982).

Fluid particle trajectories generated by $\boldsymbol{v}(\boldsymbol{x},t)$ are solutions of the differential equation

(2.5) $$\begin{eqnarray}\dot{\boldsymbol{x}}=\boldsymbol{v}(\boldsymbol{x},t),\end{eqnarray}$$

defining the flow map

(2.6) $$\begin{eqnarray}\mathscr{F}_{t_{0}}^{t}:\boldsymbol{x}_{0}\mapsto \boldsymbol{x}(t;\boldsymbol{x}_{0}),\quad t\in [t_{0},t_{1}],\end{eqnarray}$$

as the mapping from initial particle positions $\boldsymbol{x}_{0}\in U(t_{0})$ to their later positions $\boldsymbol{x}(t;\boldsymbol{x}_{0})\in U(t)$ . The assumed invariance of $U(t)$ over the time interval $[t_{0},t_{1}]$ can now be conveniently expressed as

(2.7) $$\begin{eqnarray}\mathscr{F}_{t_{0}}^{t}(U(t_{0}))=U(t),\quad t\in [t_{0},t_{1}].\end{eqnarray}$$

We will also use the notion of a material surface, which is a smooth codimension-one time-dependent surface family $\mathscr{M}(t)$ advected by the flow, i.e. $\mathscr{M}(t)=\mathscr{F}_{t_{0}}^{t}(\mathscr{M}(t_{0}))\in [t_{0},t_{1}]$ .

The deformation gradient

(2.8) $$\begin{eqnarray}\unicode[STIX]{x1D641}_{t_{0}}^{\,t}(\boldsymbol{x}_{0})=\boldsymbol{{\rm\nabla}}\mathscr{F}_{t_{0}}^{t}(\boldsymbol{x}_{0})\end{eqnarray}$$

is a linear map, taking initial infinitesimal perturbations to the fluid trajectory $\boldsymbol{x}(t;\boldsymbol{x}_{0})$ at time $t_{0}$ to their later positions at time $t$ . Although often believed otherwise, $\unicode[STIX]{x1D641}_{t_{0}}^{\,t}$ is not objective: its eigenvectors and eigenvalues depend on the frame of reference (see, e.g., Liu Reference Liu, Man and Fosdick2004). Therefore, the invariants of $\unicode[STIX]{x1D641}_{t_{0}}^{\,t}$ do not provide an objective indication of the rotational component of the deformation.

3 Finite material rotation from the dynamic polar decomposition

We seek to identify coherent Lagrangian vortices as the union of tubular material surfaces in which fluid elements exhibit the same bulk material rotation over a finite time interval of interest. Individual material fibres based at an initial point $\boldsymbol{x}_{0}$ in a deforming continuum, however, all rotate around different axes and by different angles. In recent work (Farazmand & Haller Reference Farazmand and Haller2016), we used the classic polar rotation angle (PRA) from finite strain theory to identify pointwise bulk material rotation in a moving fluid systematically. The use of the PRA, however, also leaves several challenges unaddressed, as we discuss in appendix A. Most notable of these are an inconsistency of the PRA with experimentally observed dynamic rotation angles of spherical tracers in fluids, and its lack of objectivity in three dimensions.

To address these challenges, we use here the recently developed dynamic polar decomposition (DPD) to identify a dynamically consistent and fully frame-invariant rotational component in the finite deformation of fluid elements (Haller Reference Haller2016). This decomposition gives a unique time-evolving factorization of the deformation gradient into the product of two deformation gradients: one for a purely rotational flow with zero rate of strain and one for a purely straining flow with zero vorticity.

Specifically, the unique right DPD of $\unicode[STIX]{x1D641}_{t_{0}}^{\,t}$ at $\boldsymbol{x}_{0}$ can be written as

(3.1) $$\begin{eqnarray}\unicode[STIX]{x1D641}_{t_{0}}^{\,t}=\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}\unicode[STIX]{x1D648}_{t_{0}}^{\,t},\quad t\in [t_{0},t_{1}],\end{eqnarray}$$

where the proper orthogonal dynamic rotation tensor $\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}=\partial _{\boldsymbol{ a}_{0}}\boldsymbol{a}(t)$ is the deformation gradient of the purely rotational flow

(3.2) $$\begin{eqnarray}\dot{\boldsymbol{a}}=\unicode[STIX]{x1D652}(\boldsymbol{x}(t;\boldsymbol{x}_{0}),t)\boldsymbol{a},\end{eqnarray}$$

and the non-degenerate right dynamic stretch tensor $\unicode[STIX]{x1D648}_{t_{0}}^{\,t}=\partial _{\boldsymbol{ b}_{0}}\boldsymbol{b}(t)$ is the deformation gradient of the purely straining flow

(3.3) $$\begin{eqnarray}\dot{\boldsymbol{b}}=\unicode[STIX]{x1D64A}_{t}^{\,t_{0}}\unicode[STIX]{x1D63F}(\boldsymbol{x}(t;\boldsymbol{x}_{0}),t)\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}\boldsymbol{b}.\end{eqnarray}$$

The linear velocity field for $\boldsymbol{a}(t)$ in (3.2) is strainless because its coefficient matrix $\unicode[STIX]{x1D652}$ is skew-symmetric. Similarly, the linear velocity field for $\boldsymbol{b}(t)$ in (3.3) is irrotational because its coefficient matrix $\unicode[STIX]{x1D64A}_{t}^{\,t_{0}}\unicode[STIX]{x1D63F}(\boldsymbol{x}(t;\boldsymbol{x}_{0}),t)\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}$ is symmetric. Unlike the classic polar rotation tensor (cf. appendix A), the dynamic rotation tensor $\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}$ is dynamically consistent, i.e. satisfies the fundamental superposition property of solid-body rotations:

(3.4) $$\begin{eqnarray}\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}=\unicode[STIX]{x1D64A}_{s}^{\,t}\unicode[STIX]{x1D64A}_{t_{0}}^{s},\quad s,t\in [t_{0},t_{1}].\end{eqnarray}$$

This follows because $\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}$ is the fundamental matrix solution of a classical linear differential equation and hence satisfies the process property noted in (3.4) (cf. Dafermos Reference Dafermos1971; Arnold Reference Arnold1978). In contrast, $\unicode[STIX]{x1D648}_{t_{0}}^{\,t}$ is the fundamental matrix solution of a non-classical linear differential equation with memory, i.e. with explicit dependence on the initial time $t_{0}$ . Such fundamental solutions do not obey the process property indicated in (3.4). The reason behind the dynamical inconsistency (A 3) of polar rotations is a similar memory effect to that in (A 2).

The decomposition in (3.1) is a right-type decomposition, i.e. the dynamic stretch tensor precedes the dynamic rotation tensor from the right. Just as for the classic polar decomposition, a left-type version of the DPD is also available (Haller Reference Haller2016).

4 Lagrangian-averaged vorticity deviation

Despite its dynamical consistency, the dynamic rotation tensor $\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}$ is not objective. Its frame dependence is the consequence of the frame dependence of the spin tensor $\unicode[STIX]{x1D652}(\boldsymbol{x},t)$ appearing in the differential equation (3.2). The single remaining challenge out of those listed in appendix A is, therefore, to identify an objective part of the rotation described by $\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}$ which also preserves the dynamical consistency of $\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}$ . Below we recall further results from Haller (Reference Haller2016), and use them to address this challenge.

The dynamic rotation tensor $\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}$ can further be factorized into two deformation gradients: one for a spatially uniformly rotating flow and one for a flow that describes deviations from this uniform rotation. Specifically, we have

(4.1) $$\begin{eqnarray}\unicode[STIX]{x1D64A}_{t_{0}}^{\,t}=\unicode[STIX]{x1D731}_{t_{0}}^{t}{\it\bf\Theta}_{t_{0}}^{t},\end{eqnarray}$$

where the proper orthogonal relative rotation tensor $\unicode[STIX]{x1D731}_{t_{0}}^{t}=\partial _{{\bf\alpha}_{0}}{\bf\alpha}(t)$ is dynamically consistent, serving as the deformation gradient of the relative rotation flow

(4.2) $$\begin{eqnarray}\dot{{\bf\alpha}}=[\unicode[STIX]{x1D652}(\boldsymbol{x}(t;\boldsymbol{x}_{0}),t)-\bar{\unicode[STIX]{x1D652}}(t)]{\bf\alpha}.\end{eqnarray}$$

In contrast, the proper orthogonal mean-rotation tensor ${\it\bf\Theta}_{t_{0}}^{t}=D_{{\bf\beta}_{0}}{\bf\beta}(t)$ is the deformation gradient of the mean-rotation flow

(4.3) $$\begin{eqnarray}\dot{{\bf\beta}}=\unicode[STIX]{x1D731}_{t}^{t_{0}}\bar{\unicode[STIX]{x1D652}}(t)\unicode[STIX]{x1D731}_{t_{0}}^{t}{\bf\beta}.\end{eqnarray}$$

The mean-rotation tensor ${\it\bf\Theta}_{t_{0}}^{t}$ is not dynamically consistent because (4.3) exhibits the same memory effect as discussed for (3.3).

The dynamic consistency of $\unicode[STIX]{x1D731}_{t_{0}}^{t}$ implies that the total angle swept by this tensor around its own axis of rotation is dynamically consistent. This angle ${\it\psi}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ , called the intrinsic rotation angle (see figure 1), therefore satisfies

(4.4) $$\begin{eqnarray}{\it\psi}_{t_{0}}^{t}(\boldsymbol{x}_{0})={\it\psi}_{s}^{t}(\boldsymbol{x}_{0})+{\it\psi}_{t_{0}}^{s}(\boldsymbol{x}_{0}),\quad s,t\in [t_{0},t_{1}].\end{eqnarray}$$

In addition, as shown in Haller (Reference Haller2016), ${\it\psi}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ is objective in both two and three dimensions. In two dimensions, even the tensor $\unicode[STIX]{x1D731}_{t_{0}}^{t}$ itself turns out to be objective, not just its associated scalar field ${\it\psi}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ .

Figure 1.

The geometry of the intrinsic rotation angle ${\it\psi}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ along a material trajectory $\boldsymbol{x}(t)$ in a three-dimensional deforming continuum.

Using the results obtained in Haller (Reference Haller2016), the intrinsic dynamic rotation ${\it\psi}_{t_{0}}^{t}(x_{0})$ can be computed as

(4.5) $$\begin{eqnarray}{\it\psi}_{t_{0}}^{t}(\boldsymbol{x}_{0})={\textstyle \frac{1}{2}}\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0}),\end{eqnarray}$$

with the LAVD defined here as

(4.6) $$\begin{eqnarray}\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0}):=\int _{t_{0}}^{t}|{\bf\omega}(\boldsymbol{x}(s;\boldsymbol{x}_{0}),s)-\bar{{\bf\omega}}(s)|\,\text{d}s.\end{eqnarray}$$

The objectivity of ${\it\psi}_{t_{0}}^{t}$ and LAVD can be confirmed directly from formula (2.4a,b ). Indeed, under a Euclidean observer change $\boldsymbol{x}=\unicode[STIX]{x1D64C}(t)\boldsymbol{y}+\boldsymbol{b}(t)$ , the transformed vorticity $\tilde{{\bf\omega}}(\boldsymbol{y},t)$ satisfies

(4.7) $$\begin{eqnarray}\displaystyle |\tilde{{\bf\omega}}(\boldsymbol{y}(s),s)-\tilde{\bar{{\bf\omega}}}(s)| & = & \displaystyle |\unicode[STIX]{x1D64C}^{\text{T}}(s){\bf\omega}(\boldsymbol{x}(s),s)+\unicode[STIX]{x1D64C}^{\text{T}}(t)\dot{\boldsymbol{q}}(t)-(\unicode[STIX]{x1D64C}^{\text{T}}(s)\bar{{\bf\omega}}(s)+\unicode[STIX]{x1D64C}^{\text{T}}(t)\dot{\boldsymbol{q}}(t))|\nonumber\\ \displaystyle & = & \displaystyle |\unicode[STIX]{x1D64C}^{\text{T}}(s)[{\bf\omega}(\boldsymbol{x}(s),s)-\bar{{\bf\omega}}(s)]|\nonumber\\ \displaystyle & = & \displaystyle |{\bf\omega}(\boldsymbol{x}(s),s)-\bar{{\bf\omega}}(s)|,\end{eqnarray}$$

because the rotation matrix $\unicode[STIX]{x1D64C}^{\text{T}}(s)$ preserves the length of vectors. We summarize the results of this section in a theorem.

Theorem 1. For an infinitesimal fluid volume starting from $\boldsymbol{x}_{0}$ , the $\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ field is a dynamically consistent and objective measure of bulk material rotation relative to the spatial mean rotation of the fluid volume $U(t)$ . Specifically, $\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ is twice the intrinsic dynamic rotation angle generated by the relative rotation tensor $\unicode[STIX]{x1D731}_{t_{0}}^{t}$ . The latter tensor is obtained from the dynamically consistent decomposition

(4.8) $$\begin{eqnarray}\unicode[STIX]{x1D641}_{t_{0}}^{\,t}=\unicode[STIX]{x1D731}_{t_{0}}^{t}{\it\bf\Theta}_{t_{0}}^{t}\unicode[STIX]{x1D648}_{t_{0}}^{\,t},\end{eqnarray}$$

with the deformation gradient ${\it\bf\Theta}_{t_{0}}^{t}$ of a pure rigid-body rotation and with the deformation gradient $\unicode[STIX]{x1D648}_{t_{0}}^{\,t}$ of a unique purely straining flow.

For detailed proofs of all statements in Theorem 1, we refer the reader to Haller (Reference Haller2016). Importantly, this theorem enables the extraction of an objective and dynamically consistent material rotation component from the deformation gradient without carrying out the differentiation with respect to the initial conditions in the definition (2.8) of $\unicode[STIX]{x1D641}_{t_{0}}^{\,t}$ .

5 Rotationally coherent Lagrangian vortices

We now use the LAVD to identify objectively material tubes along which small fluid volumes experience the same bulk rotation over $[t_{0},t_{1}]$ relative to the mean rigid-body rotation of the fluid. By Theorem 1, the time $t_{0}$ positions of such material tubes are tubular level surfaces of the scalar function $\text{LAVD}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ . By a tubular set, we mean here a convex, cylindrical, cup-shaped or toroidal set in three dimensions, and a closed convex curve in two dimensions. We require convexity for tubular surfaces, motivated by the near-circular cross-section generally observed for stable vortices.

If the gradient $\partial _{\boldsymbol{x}_{0}}\text{LAVD}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ is non-zero along a tubular LAVD level surface, then this level surface is surrounded by a continuous nested family of tubular level surfaces (Milnor Reference Milnor1963). The singular centre of such a nested sequence of tubes, with inward-increasing LAVD values, gives a definition of a Lagrangian vortex centre. Similarly, the largest convex member of such a nested tube family defines the boundary of a Lagrangian vortex. We summarize these concepts in the following definition, with its geometry illustrated in figure 2.

We refer to the evolving positions $\mathscr{L}(t)$ of the tubular level sets

(5.1) $$\begin{eqnarray}\mathscr{L}(t_{0})=\left\{\boldsymbol{x}_{0}\in U(t_{0}):\text{LAVD}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})=\text{const.}\right\}\end{eqnarray}$$

as rotational Lagrangian coherent structures (rotational LCSs), as indicated in figure 2. These LCSs give a foliation of the evolving Lagrangian vortex $\mathscr{V}(t)$ into tubes along which material elements complete the same intrinsic dynamic rotation ${\it\psi}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ .

Figure 2.

Initial and current positions of a rotationally coherent Lagrangian vortex $\mathscr{V}(t)$ , a Lagrangian vortex boundary $\mathscr{B}(t)$ and a Lagrangian vortex centre $\mathscr{C}(t)$ . Also shown within $\mathscr{V}(t)$ is a rotational LCS $\mathscr{L}(t)$ , i.e. a material surface along which volume elements exhibit the same intrinsic dynamic rotation over the time interval $[t_{0},t_{1}].$ As a consequence, $\mathscr{V}(t)$ experiences only tangential material filamentation without global breakaway.

For rotational LCSs, as well as Lagrangian vortices, their boundaries and centres are material objects by definition. Therefore, their time $t$ position is uniquely determined by Lagrangian advection of their initial positions:

(5.2a-c ) $$\begin{eqnarray}\mathscr{L}(t)=\mathscr{F}_{t_{0}}^{t}(\mathscr{L}(t_{0})),\quad \mathscr{B}(t)=\mathscr{F}_{t_{0}}^{t}(\mathscr{B}(t_{0})),\quad \mathscr{C}(t)=\mathscr{F}_{t_{0}}^{t}(\mathscr{C}(t_{0})),\quad t\in [t_{0},t_{1}].\end{eqnarray}$$

Recent stretching-based definitions of Lagrangian vortices allow for no filamentation in the boundary of the vortex (Haller & Beron-Vera Reference Haller and Beron-Vera2013; Blazevski & Haller Reference Blazevski and Haller2014). In contrast, the LAVD-based definition of a rotational LCS allows for tangential material filamentation. The filamented part of the material surface, however, still rotates together with the LCS without global breakaway.

The above definitions capture Lagrangian vortices with the simplest (i.e. convex) geometry at time $t_{0}$ . More generally, one may allow for small tangential filamentation to be a priori present in the vortex boundary even at time $t_{0}$ . This involves the relaxation of the convexity of $\mathscr{L}(t_{0})$ and $\mathscr{B}(t_{0})$ to material surfaces with small convexity deficiency, as discussed along with other numerical aspects in § 9.

In geophysical flows over a rotating planet, rotationally coherent Lagrangian vortices can be directly computed from the flow induced in the curvilinear coordinate space instead of the curved surface of the planet (cf. appendix B). By construction, the resulting vortices and their centres are invariant with respect to time-dependent rotations and translations within the space of curvilinear coordinates. (Frame invariance cannot be defined restricted to a curvilinear surface, as Euclidean frame changes take the observer off the surface.) With this approach, one simply computes the classic Euclidean vorticity of the longitudinal and latitudinal coordinate speeds, as opposed to computing the vorticity in curvilinear coordinates.

6 Rotationally coherent Eulerian vortices

Over a short time interval $[t_{0},t+s]$ with $|s|\ll 1$ , we can Taylor expand the LAVD field as

(6.1) $$\begin{eqnarray}\text{LAVD}_{t_{0}}^{t+s}(\boldsymbol{x}_{0})=\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})+\text{IVD}(\boldsymbol{x}(t;\boldsymbol{x}_{0}),t)\cdot s+O(s^{2}),\end{eqnarray}$$

with the instantaneous vorticity deviation (IVD) defined as

(6.2) $$\begin{eqnarray}\text{IVD}(\boldsymbol{x},t):=|{\bf\omega}(\boldsymbol{x},t)-\bar{{\bf\omega}}(t)|.\end{eqnarray}$$

By the calculation (4.7), the IVD field is objective. By (6.1), $\text{IVD}(\boldsymbol{x}(t;\boldsymbol{x}_{0}),t)$ describes the rate of change of the LAVD field at an initial condition $\boldsymbol{x}_{0}$ under increasing integration time.

Using the IVD, we now introduce the instantaneous notion of a rotationally coherent Eulerian vortex by taking the limit $t_{0},t_{1}\rightarrow t$ in Definition 1. At a time $t\in [t_{0},t_{1}],$ such an Eulerian vortex is composed of tubular surfaces along which the intrinsic rotation rates $\dot{{\it\psi}}_{t}^{t}$ of fluid elements are equal. Indeed, by formula (4.5), we have $\dot{{\it\psi}}_{t}^{t}(\boldsymbol{x})=(\text{IVD}(\boldsymbol{x},t))/2=\text{const.}$ along these tubular surfaces. The following definition summarizes the details for this objective Eulerian vortex concept.

A rotational Eulerian coherent structure (rotational ECS) is then just a level surface

(6.3) $$\begin{eqnarray}E(t)=\left\{\boldsymbol{x}\in U(t):\text{IVD}(\boldsymbol{x},t)=C_{0}(t)\right\}\end{eqnarray}$$

along which material elements experience the same dynamic rotation rate $\dot{{\it\psi}}_{t}^{t}(x)=C_{0}(t)/2$ . Unlike rotational LCSs, individual rotational ECSs are instantaneous quantitates without a well-defined evolution, unless $C_{0}(t)$ is specifically selected as constant over time. We illustrate the geometry of Definition 2 in figure 3.

Figure 3.

A rotationally coherent Eulerian vortex $\mathscr{V}(t)$ , with boundary $\mathscr{B}(t)$ and vortex centre $\mathscr{C}(t)$ . Also shown within $\mathscr{V}(t)$ is a rotational ECS $\mathscr{E}(t)$ , i.e. a surface along which volume elements exhibit the same intrinsic dynamic rotation rate instantaneously at time $t$ . No material coherence is guaranteed for the Eulerian vortex $\mathscr{V}(t)$ under passive advection.

Given that

(6.4) $$\begin{eqnarray}\frac{\text{d}}{\text{d}t}\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})=\text{IVD}(\boldsymbol{x}(t;\boldsymbol{x}_{0}),t),\end{eqnarray}$$

rotationally coherent Eulerian vortices are effectively the derivatives of rotationally coherent Lagrangian vortices with respect to the length of the extraction time of the latter. Setting $t=t_{0}$ in (6.1) and observing that $\text{LAVD}_{t_{0}}^{t_{0}}(\boldsymbol{x}_{0})=0,$ we conclude that the $\mathscr{B}(t_{0})$ initial positions of rotationally coherent Lagrangian vortex boundaries evolve precisely from their Eulerian coherent counterparts $B(t_{0})$ as the Lagrangian extraction time $t_{1}-t_{0}$ increases from zero (see figure 4).

Figure 4.

A rotationally coherent Eulerian vortex is the limit of a rotationally coherent Lagrangian vortex under vanishing extraction time (i.e. integration time).

Physically, these Eulerian vortices are built of tubular surfaces showing instantaneous coherence in the rate of their bulk material rotation. This instantaneous coherence in rotation rates generally does not imply sustained coherence in the finite rotation of material trajectories released from these surfaces. Furthermore, the interior of an evolving vortex $V(t)$ is not a material domain.

By the objectivity of the IVD, Definition 2 nevertheless gives an objective definition of an Eulerian vortex, its centre and boundary. To the best of our knowledge, no other objective three-dimensional Eulerian vortex definition has been proposed in the literature. Specifically, none of the Eulerian criteria reviewed or proposed in Jeong & Hussein (Reference Jeong and Hussein1995), Chakraborty et al. (Reference Chakraborty, Balachandar and Adrian2005) and Haller (Reference Haller2005) are invariant under time-dependent rotations and translations of the observer. Since truly unsteady flows have no distinguished frames of reference (Lugt Reference Lugt, Müller, Riesner and Schmidt1979), non-objective Eulerian vortex criteria do not yield well-defined material vortices for unsteady flows.

A relevant discussion can be found in Jeong & Hussein (Reference Jeong and Hussein1995) about the $|{\bf\omega}|$ -criterion, by which $|{\bf\omega}|$ must exceed a preselected threshold within a vortex. This approach is found to be intuitive but inadequate in Jeong & Hussein (Reference Jeong and Hussein1995) for several reasons. We agree with this general assessment, because the $|{\bf\omega}|$ -criterion is threshold-dependent and not objective. In contrast, our rotationally coherent Eulerian vortex definition in Definition 2 is based on a threshold-independent and objective assessment of the level-surface topology of $|{\bf\omega}-\bar{{\bf\omega}}|$ .

Importantly, for flows with zero mean vorticity in their frame of definition, the outermost convex tubular level sets of the vorticity magnitude $|{\bf\omega}|$ or of the enstrophy $|{\bf\omega}|^{2}$ coincide with rotationally coherent Eulerian vortices. Therefore, when properly interrogated, the vorticity and enstrophy distributions of zero-mean-vorticity flows do reveal objective structures that can be viewed as derivatives of rotationally coherent Lagrangian vortices.

7 Rotationally coherent vortices in planar flows

In a flow defined in the $(x_{1},x_{2})$ plane, the LAVD (4.6) takes the simple form

(7.1) $$\begin{eqnarray}\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})=\int _{t_{0}}^{t}|{\it\omega}_{3}(\boldsymbol{x}(s;\boldsymbol{x}_{0}),s)-\bar{{\it\omega}}_{3}(s)|\,\text{d}s,\end{eqnarray}$$

with ${\it\omega}_{3}$ referring here to the $x_{3}$ component of the vorticity vector ${\bf\omega}(\boldsymbol{x},t)=(0,0,{\it\omega}_{3}(\boldsymbol{x},t))$ , and with the mean vorticity $\bar{{\it\omega}}_{3}(s)$ .

By Definition 1, a rotational LCS evolves over the time interval $[t_{0},t_{1}]$ via advection from a closed and convex level curve of the $\text{LAVD}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ computed in (7.1). By the same definition, Lagrangian vortex boundaries are outermost members of such rotational LCS families. Similarly, Lagrangian vortex centres are advected positions of isolated maxima of $\text{LAVD}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ .

By Definition 2, a rotational ECS at time $t$ is a closed and convex level curve of

(7.2) $$\begin{eqnarray}\text{IVD}(\boldsymbol{x},t)=|{\it\omega}_{3}(\boldsymbol{x},t)-\bar{{\it\omega}}_{3}(t)|,\end{eqnarray}$$

around one of its local maxima. Accordingly, rotationally coherent Eulerian vortices are outermost members of such nested curve families with outwards non-increasing instantaneous IVD values.

For two-dimensional flows only, Haller (Reference Haller2016) shows that the relative dynamic rotation angle

(7.3) $$\begin{eqnarray}{\it\phi}_{t_{0}}^{t}(\boldsymbol{x}_{0}):=\frac{1}{2}\int _{t_{0}}^{t}[{\it\omega}_{3}(\boldsymbol{x}(s;\boldsymbol{x}_{0}),s)-\bar{{\it\omega}}_{3}(s)]\,\text{d}s\end{eqnarray}$$

is also an objective and dynamically consistent measure of rotation. It measures the net rotation angle generated by the relative rotation tensor around the $x_{3}$ axis, with sign changes in the rotation accounted for. (In contrast, ${\it\psi}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ measures the rotation about the evolving instantaneous axis of rotation of $\unicode[STIX]{x1D731}_{t_{0}}^{t}$ , which points in the $-x_{3}$ direction when ${\it\omega}_{3}-\bar{{\it\omega}}_{3}$ is negative.)

While ${\it\psi}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ , as the total angle swept by the relative rotation tensor over the time interval $[t_{0},t_{1}]$ , is always positive, the sign of the angle ${\it\phi}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ is unrestricted. In vortical regions preserving the pointwise sign of the Lagrangian vorticity ${\it\omega}_{3}(\boldsymbol{x}(s;\boldsymbol{x}_{0}),s)-\bar{{\it\omega}}_{3}(s)$ , contours of ${\it\phi}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ , contours of $\text{LAVD}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ and contours of the Lagrangian-averaged vorticity (LAV)

(7.4) $$\begin{eqnarray}\text{LAV}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0}):=\int _{t_{0}}^{t_{1}}{\it\omega}_{3}(\boldsymbol{x}(s;\boldsymbol{x}_{0}),s)\,\text{d}s\end{eqnarray}$$

all coincide with each other, albeit they generally correspond to different values of the underlying scalar fields. (It should be noted that $\text{LAV}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ is not an objective quantity, but its level curves are objectively defined.) In such regions, the sign of ${\it\phi}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ will also carry objective information about the direction of the relative rotation.

Differences in the contours of $\text{LAVD}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ , ${\it\phi}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ and $\text{LAV}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ will arise in regions where the sign of ${\it\omega}_{3}(\boldsymbol{x}(s),s)-\bar{{\it\omega}}_{3}(s)$ crosses zero, i.e. near the boundaries of regions with a well-defined sign in their deviation from the mean vorticity of the flow. In such regions, connected level curves of $\text{LAV}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ can group together initial conditions with substantially different global rotational histories, as long as their final net rotations are equal. The use of $\text{LAVD}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ eliminates such coincidental agreement in the net rotation angles.

8 Geostrophic Lagrangian vortex centres are attractors for inertial particles

Consider a small spherical particle of radius $r_{0}$ and density ${\it\rho}_{part}$ in a geostrophic flow of density ${\it\rho}$ and viscosity ${\it\nu}$ . Under the ${\it\beta}$ -plane approximation, let $f$ denote the Coriolis parameter (twice the local vertical component of the angular velocity of the Earth). Applying a slow-manifold reduction to the Maxey–Riley equations (Maxey & Riley Reference Maxey and Riley1983) in the limit of small Rossby numbers, Beron-Vera et al. (Reference Beron-Vera, Olascoaga, Haller, Farazmand, Triñanes and Wang2015) showed that the inertial particle motion satisfies

(8.1) $$\begin{eqnarray}\dot{\boldsymbol{x}}=\boldsymbol{v}(\boldsymbol{x},t)+{\it\tau}({\it\delta}-1)f\unicode[STIX]{x1D645}\boldsymbol{v}(\boldsymbol{x},t)+O({\it\tau}^{2}),\quad \unicode[STIX]{x1D645}=\left(\begin{array}{@{}cc@{}}0 & -1\\ 1 & 0\end{array}\right),\end{eqnarray}$$

where

(8.2a,b ) $$\begin{eqnarray}{\it\delta}=\frac{{\it\rho}}{{\it\rho}_{part}},\quad {\it\tau}:=\frac{2r_{0}^{2}}{9{\it\nu}{\it\delta}}.\end{eqnarray}$$

Provenzale (Reference Provenzale1999) considered the Maxey–Riley equation in the same physical setting, but without a slow-manifold reduction. His second-order differential equation also included additional terms that either vanish along the ${\it\beta}$ -plane or appear at higher order in the reduced first-order equation (8.1) (cf. Beron-Vera et al. (Reference Beron-Vera, Olascoaga, Haller, Farazmand, Triñanes and Wang2015) for more detail).

Remarkably, in the limit of vanishing Rossby numbers, cyclonic attractors for light particles $({\it\delta}>1)$ and anticyclonic attractors for heavy particles $({\it\delta}<1)$ in (8.1) turn out to be precisely the rotationally coherent vortex centres defined in Definition 1. The same Lagrangian vortex centres act as cyclonic repellers for heavy particles and anticyclonic repellers for light particles. We state these results in more detail as follows.

Theorem 2. Assume that $\boldsymbol{x}_{0}^{\ast }\in U(t_{0})$ is the initial position of a rotationally coherent Lagrangian vortex centre whose relative rotation keeps constant sign, i.e.

(8.3) $$\begin{eqnarray}\text{sign}[{\it\omega}_{3}(\boldsymbol{x}_{0}^{\ast },t)-\bar{{\it\omega}}_{3}(t)]={\it\mu}(\boldsymbol{x}_{0}^{\ast }),\quad t\in [t_{0},t_{1}],\end{eqnarray}$$

for an appropriate sign constant ${\it\mu}(\boldsymbol{x}_{0}^{\ast })\in \{-1+1\}$ . Then, for ${\it\tau}>0$ small enough, the following hold:

1. (i) in a cyclonic ( ${\it\mu}(\boldsymbol{x}_{0}^{\ast })f>0$ ) rotationally coherent vortex, there exists a finite-time attractor (repeller) for light (heavy) particles that stays $O({\it\tau})$ close to the vortex centre $\mathscr{C}(t)$ ;

2. (ii) in an anticyclonic ( ${\it\mu}(\boldsymbol{x}_{0}^{\ast })f<0$ ) rotationally coherent vortex, there exists a finite-time attractor (repeller) for heavy (light) particles that stays $O({\it\tau})$ close to the vortex centre $\mathscr{C}(t)$ .

Theorem 2 provides an independent experimentally verifiable justification for defining vortex centres as in Definition 1 for geostrophic flows. Specifically, positively buoyant drifters or floating debris released well inside a cyclonic oceanic eddy will spiral onto the evolving Lagrangian vortex centre identified from Definition 1 (figure 5). We will illustrate this effect using simulated inertial particle motion on satellite-based ocean velocities in § 10.5.

Figure 5.

By Theorem 2, cyclonic Lagrangian eddy centres evolving from LAVD maxima act as observed attractors for small positively buoyant particles, such as drifters and floating debris.

9 Numerical aspects

Computation of a rotation angle from the classic polar decomposition requires the computation of the deformation gradient $\unicode[STIX]{x1D641}_{t_{0}}^{\,t}$ (cf. (A 1) in appendix A). This is achieved either by the numerical differentiation of fluid trajectories with respect to their initial conditions or by solving the equation of variations $(\text{d}/\text{d}t)\unicode[STIX]{x1D641}_{t_{0}}^{\,t}=\boldsymbol{{\rm\nabla}}\boldsymbol{v}(\boldsymbol{x}(t;t_{0},\boldsymbol{x}_{0}),t)\unicode[STIX]{x1D641}_{t_{0}}^{\,t}$ , whose solutions typically grow exponentially. Either way, computation of polar rotation has the same long-time numerical sensitivity that arises in computing the invariants of the Cauchy–Green strain tensor.

In contrast, computation of rotationally coherent vortices based on Definition 1 only requires the integration of the normed vorticity deviation along fluid trajectories. This is simplest to do simultaneously with trajectory integration, solving the extended system of differential equations

(9.1) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\displaystyle \frac{\text{d}x}{\text{d}t}=\boldsymbol{v}(\boldsymbol{x},t),\\ \displaystyle \frac{\text{d}}{\text{d}t}\text{LAVD}_{t_{0}}^{t}=|\boldsymbol{{\rm\nabla}}\times \boldsymbol{v}(\boldsymbol{x},t)-\overline{\boldsymbol{{\rm\nabla}}\times \boldsymbol{v}}(t)|\end{array}\right\} & & \displaystyle\end{eqnarray}$$

over a grid of initial conditions over the time interval $[t_{0},t_{1}]$ . In our experience, however, features of the $\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ turn out to be sharper when the trajectory ordinary differential equation (ODE) $\dot{\boldsymbol{x}}=\boldsymbol{v}(\boldsymbol{x},t)$ is solved first, and the vorticity is subsequently integrated along trajectories. This is because adaptive ODE solvers make different decisions about time steps when the trajectory ODE is amended with the ODE for the LAVD field. This is especially so when both ODEs are solved simultaneously over large grids of initial conditions.

The computational domain for solving (9.1) is just $U(t)\equiv U$ in the case of a closed flow with a fixed boundary. In open flows, the focus may be on vortices on a smaller domain. In that case, the domain should still be chosen large enough so that the averaged vorticity $\overline{\boldsymbol{{\rm\nabla}}\times \boldsymbol{v}}(t)$ is representative of the overall mean rotation of the fluid mass under study. In geophysical flows, this mean rotation is expected to be zero, which is confirmed by our calculations even for domains of the size of a few degrees.

For two-dimensional flows, we first identify local maxima of the LAVD field, then extract nearby closed LAVD level curves. In all of our computations, we use the level-set function of MATLAB for this purpose, and identify the closedness of a level curve by probing the output from this function. Definition 1 then requires the identification of the outermost convex LAVD level curve around an LAVD maximum as vortex boundary. This convexity requirement is conservative, ensuring that the material vortex starts out unfilamented at the initial time $t_{0}$ . At later times, our approach allows for tangential filamentation in the advected LCS due to local strain, but disallows large-scale filamentation arising from differences in the bulk rotation along material filaments. As a consequence, filaments developed by rotational LCSs rotate together with the main body of the underlying material vortex.

In actual computations, one reason to relax strict convexity for closed LAVD level surfaces is that they are numerically represented by discrete polygons. The more vertices such a polygon has, the more likely it is that the polygon is not convex, even if the approximated level curve is. A second reason for relaxing convexity is to remove the conceptual asymmetry of figure 2, allowing for small tangential filamentation even in the initial positions of vortex boundaries. A third reason for relaxing strict convexity in multi-scale data sets is the presence of smaller-scale vortices near the perimeter of a larger-scale vortex. This necessitates the use of an appropriately coarse-grained notion for the boundary of the larger vortex.

In two dimensions, addressing the finite grid, the initial tangential filamentation and the multi-scale challenges can be achieved by allowing for a small convexity deficiency in the LAVD level curves (Gonzalez and Woods Reference Gonzalez and Woods2008; Batchelor & Whelan Reference Batchelor and Whelan2012). Here, we define the convexity deficiency of a closed curve in the plane as the ratio of the area between the curve and its convex hull to the area enclosed by the curve. In figure 6, we show cases of closed LAVD level sets with small convexity deficiency. In three dimensions, we require small convexity deficiency for the one-dimensional intersection of LAVD level surfaces with a family $\{\mathscr{P}_{i}\}$ of planes transverse to the expected vortex centre curve $\mathscr{C}(t_{0})$ .

Figure 6.

Cases of closed LAVD level sets benefiting from the relaxation of convexity to small convexity deficiency. The shaded area indicates the area difference between the level set and its convex hull. (a) A convex LAVD level set with a non-convex approximating polygon arising from discretization. (b) A closed LAVD level set with minor initial filamentation in the tangential direction. (c) A closed LAVD level set that is only convex after a filtering of smaller-scale vortices.

In our computations, we set the convexity deficiency bound to $10^{-3}$ or lower, which we find to produce robust results for the first two cases covered in figure 6. Applications to multi-scale data sets (figure 6 c) will probably require a somewhat higher bound. In general, increasing the convexity deficiency bound produces larger eddies that also tend to filament more.

Small-scale local maxima of the LAVD function also arise due to numerical or observational noise in the velocity data. To eliminate the resulting artificial vortex candidates, we choose to ignore closed LAVD contours whose arclength falls below a minimal threshold. In a given application, this threshold should be chosen below the minimal vortex perimeter that is expected to be reliably resolved by the data set.

We summarize the extraction procedure for rotationally coherent Lagrangian vortices in two- and three-dimensional flows in Algorithms 1 and 2 respectively. The extraction of their rotationally coherent Eulerian counterparts follows the same steps (ii)–(iii), but applied to the function $\text{IVD}(\boldsymbol{x},t)$ instead of $\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})$ .

The procedure outlined in Algorithm 2 is mathematically well defined by the level-surface topology of the LAVD. Step (2) of the algorithm describes one possible numerical extraction scheme for these level surfaces, using their intersections with planes transverse to the anticipated vortex centre $\mathscr{C}(t_{0})$ . (A MATLAB implementation of this algorithm is available under https://github.com/LCSETH.) Complex flow geometries with multiple vortices and a priori unknown vortex orientations will require more involved numerical approaches to level-set extraction.

To implement Algorithms 1 and 2 in the forthcoming examples, we use a variable-order Adams–Bashforth–Moulton solver (ODE113 in MATLAB) for trajectory advection. The absolute and relative tolerances of the ODE solver are chosen as $10^{-6}$ or higher. In §§ 10.3–10.5 and 11.2, we use cubic and bilinear interpolation schemes for computing pointwise vorticity values for the IVD and LAVD functions respectively. The lower-order interpolation suffices for the LAVD calculation, because trajectory advection smoothes out the LAVD contours in our experience.

10 Two-dimensional examples

10.1 Planar Euler flows

On any solution of the two-dimensional Euler equation on the $(x_{1},x_{2})$ coordinate plane, the scalar vorticity value ${\it\omega}_{3}(\boldsymbol{x},t)$ is preserved along trajectories. Therefore, the LAVD defined in (4.6) simplifies to

(10.1) $$\begin{eqnarray}\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})=(t-t_{0})|{\it\omega}_{3}(\boldsymbol{x}_{0},t_{0})-\bar{{\it\omega}}_{3}(t_{0})|=(t-t_{0})\,\text{IVD}(\boldsymbol{x}_{0},t_{0}).\end{eqnarray}$$

Consequently, at time $t_{0}$ , the boundaries and centres of all rotationally coherent Lagrangian vortices coincide with those of rotationally coherent Eulerian vortices. Specifically, ${\mathcal{B}}(t_{0})=B(t_{0})$ are outermost, closed and convex level curves of $|{\it\omega}_{3}(\boldsymbol{x}_{0},t_{0})-\bar{{\it\omega}}_{3}(t_{0})|$ , encircling a local maximum $\mathscr{C}(t_{0})=C(t_{0})$ of $|{\it\omega}_{3}(\boldsymbol{x}_{0},t_{0})-\bar{{\it\omega}}_{3}(t_{0})|$ .

If the planar Euler flow is steady, then

(10.2) $$\begin{eqnarray}\frac{\text{d}{\bf\omega}(\boldsymbol{x}(t))}{\text{d}t}=\boldsymbol{{\rm\nabla}}{\bf\omega}(\boldsymbol{x}(t))\boldsymbol{\cdot }\boldsymbol{v}(\boldsymbol{x}(t))=0.\end{eqnarray}$$

Consequently, as long as the vorticity is not constant on open sets, the vorticity contours coincide with streamlines. In that case, both the Lagrangian and the Eulerian vortices are bounded by outermost closed streamlines, with their centres marked by a centre-type fixed point, as expected.

Despite these formal LAVD calculations for planar Euler flows, one must remember that the eternal conservation of vorticity in these flows creates a high degree of degeneracy for the rotation of fluid parcels. Indeed, despite the temporal and spatial complexity of a planar Euler velocity field $\boldsymbol{v}(\boldsymbol{x},t)$ , a fluid element travelling along a fluid trajectory $\boldsymbol{x}(t)$ will keep its initial angular velocity (i.e. half of its vorticity at $t_{0}$ ) for all times $t$ . This property locks the initial conditions of initially rotationally coherent fluid parcels to the same LAVD level curve for all times, no matter how much these parcels separate or deform in the meantime. Specifically, even if material filaments break away transversely from a vortical region and undergo high stretching and global filamentation, they will eternally keep the exact same pointwise angular velocity they initially acquired near the vortex (see, e.g., the Eulerian vortex interaction simulation of Dritschel & Waugh (Reference Dritschel and Waugh1992) for a striking example).

This rotational degeneracy of the planar Euler equation is lost under the addition of the slightest viscous dissipation, compressibility or three-dimensionality. Under any of these regularizing perturbations, LAVD-based vortex detection can be applied to reveal the Lagrangian signature of vortex interactions in the perturbed flow. We illustrate this in a viscous perturbation of the contour-dynamics simulation of Dritschel & Waugh (Reference Dritschel and Waugh1992) in Section 10.4.

10.2 Irrotational vortices

Although physically unrealizable, swirling flows with regions of zero vorticity are important theoretical models. The classic irrotational vortex flow is given by the two-dimensional circularly symmetric velocity field

(10.3) $$\begin{eqnarray}\boldsymbol{v}(\boldsymbol{x})=\left(\begin{array}{@{}c@{}}\displaystyle \frac{-{\it\alpha}x_{2}}{x_{1}^{2}+x_{2}^{2}}\\ \displaystyle \frac{{\it\alpha}x_{1}}{x_{1}^{2}+x_{2}^{2}}\end{array}\right).\end{eqnarray}$$

The vorticity of this flow is identically zero, implying $\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})\equiv 0$ and $\text{IVD}(\boldsymbol{x})\equiv 0$ for any choice of $t_{0}$ and $t$ . This may seem to be at odds with the fact that all particles move on circular orbits, and hence line elements tangent to the trajectories exhibit one full rotation over the period of the circular orbit.

However, the tangents are special line elements and are not representative of the overall bulk rotation of fluid elements. Material fibres transverse to the trajectory tangents all rotate at different speeds. The average rate of rotation for all line elements emanating from a given point is zero, as was already pointed out by Helmholtz (Reference Helmholtz1858) in a related debate with Bertrand (Reference Bertrand1873) (cf. Truesdell & Rajagopal Reference Truesdell and Rajagopal2009; Haller Reference Haller2016). Therefore, an infinitesimal fluid volume has no experimentally measurable bulk rigid-body rotation component in its deformation, as indeed demonstrated by Shapiro (Reference Shapiro1961).

The irrotational vortex flow (10.3) is composed of a single degenerate LAVD and IVD level set $\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})=\text{IVD}(\boldsymbol{x})\equiv 0$ , as opposed to a nested family of codimension-one tubular level sets with outward-decreasing LAVD and IVD values. Therefore, the flow (10.3) does not satisfy our Definitions 1 and 2, and hence does not qualify as a rotationally coherent Lagrangian or Eulerian vortex. A similar conclusion holds on the irrotational domains of the Rankine and Lamb–Oseen vortices (Majda & Bertozzi Reference Majda and Bertozzi2002), or, more generally, in a recirculation region of any potential flow. We note that irrotational vortices are also either explicitly excluded by most systematic vortex criteria or fail the test for being a vortex according to these criteria.

10.3 Direct numerical simulation of two-dimensional turbulence

We solve numerically the forced Navier–Stokes equation

(10.4) $$\begin{eqnarray}\partial _{t}\boldsymbol{v}+\boldsymbol{v}\boldsymbol{\cdot }\boldsymbol{{\rm\nabla}}\boldsymbol{v}=-\boldsymbol{{\rm\nabla}}p+{\it\nu}{\it\Delta}\boldsymbol{v}+\boldsymbol{f},\quad \boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\boldsymbol{v}=0,\end{eqnarray}$$

for a two-dimensional velocity field $\boldsymbol{v}(\boldsymbol{x},t)$ with $\boldsymbol{x}=(x_{1},x_{2})\in U=[0,2{\rm\pi}]\times [0,2{\rm\pi}]$ . We use a pseudo-spectral code with viscosity ${\it\nu}=10^{-5}$ . We evolve a random-in-phase velocity field in the absence of forcing ( $\boldsymbol{f}=0$ ) until the flow is fully developed, then turn on a random-in-phase forcing (cf. Farazmand & Haller Reference Farazmand and Haller2013). We identify this latter time instance with the initial time $t_{0}=0$ , and run the simulation until the final time $t_{1}=50.$

To construct the LAVD field, we advect trajectories from an initial grid of $1024\times 1024$ points over the time interval $[0,50]$ . We integrate the vorticity deviation norm separately (as opposed to solving the combined ODE (9.1)), using 1200 vorticity values along each trajectory, equally spaced in time. Figure 7 shows the results superimposed on the contours of $\text{LAVD}_{0}^{50}(\boldsymbol{x}_{0})$ . In this computation, we have set the minimum arclength $l_{min}=0.3$ and convexity deficiency bound $d_{max}=10^{-4}$ .

Figure 7 shows the Lagrangian vortex boundaries extracted from Definition 1 using Algorithm 1 at the initial time $t_{0}=0$ . In figure 7, we confirm the Lagrangian rotational coherence of these vortex boundaries by advecting them to the final time $t_{1}=50$ . As expected, the vortex boundaries do not give in to the general trend of exponential stretching and folding observed for generic material lines. Instead, they display only local (tangential) filamentation. The complete advection sequence over the time interval $[0,50]$ is illustrated in the online supplementary movie M1 available at http://dx.doi.org/10.1017/jfm.2016.151.

Figure 7.

(a) Three-dimensional contour plot of $\text{LAVD}_{0}^{50}(\boldsymbol{x}_{0})$ . Lagrangian vortex boundaries (outermost rotational LCSs) extracted from Algorithm 1 are shown in red. (b) Vortex boundaries shown in the flow domain at initial time $t=0$ . The contour plot of $\text{LAVD}_{0}^{50}(\boldsymbol{x}_{0})$ is shown in the background for reference. (c) Advected Lagrangian vortex boundaries at the final time $t_{1}=50$ (see the online supplementary movie M1 for the complete advection sequence of the vortex boundaries).

In figure 8(a), we show the time $t_{0}$ positions of rotationally coherent Eulerian vortices in green, along with their Lagrangian counterparts in red. The passively advected positions of the Eulerian vortices are shown in figure 8(b), displaying substantial material fingering into their surroundings. While this lack of full material coherence over time is expected for the Eulerian vortices, some of them are impressively close to their Lagrangian counterparts at the initial time $t_{0}$ , even though they are extracted from just an instantaneous analysis. At the same time, other Eulerian vortices without nearby Lagrangian counterparts disintegrate completely under advection, as seen in figure 8(c).

Figure 8.

(a) Lagrangian (red) and Eulerian (green) rotationally coherent vortices at time $t_{0}$ in the two-dimensional turbulence example. (b) Rotationally coherent Eulerian vortex boundaries and centres at time $t_{0}$ , with the $\text{IVD}(\boldsymbol{x},t_{0})$ field shown in the background. (c) The same objects advected passively to the final time $t_{1}=50$ , with $\text{IVD}(\boldsymbol{x},t_{1})$ shown in the background.

10.4 Interaction of two unequal vortices

Motivated by the inviscid contour-dynamics simulation of Dritschel & Waugh (Reference Dritschel and Waugh1992), we seek to identify Lagrangian vortex cores in the interaction of two vortices of unequal strength (cf. our discussion in § 10.1). We solve the Navier–Stokes equation (10.4) in two dimensions with viscosity ${\it\nu}=2\times 10^{-6}$ . As initial velocity field, we let

(10.5a,b ) $$\begin{eqnarray}\boldsymbol{v}(\boldsymbol{x},0)=\boldsymbol{v}_{1}(\boldsymbol{x})+\boldsymbol{v}_{2}(\boldsymbol{x}),\quad \boldsymbol{v}_{i}(\boldsymbol{x})={\it\Gamma}_{i}\frac{1-\exp (-r_{i}^{2}/{\it\delta}_{i}^{2})}{2{\rm\pi}r_{i}^{2}}\,\unicode[STIX]{x1D64D}_{{\rm\pi}/2}\,(\boldsymbol{x}-\boldsymbol{x}_{i}),\quad i=1,2,\qquad \quad\end{eqnarray}$$

where $r_{i}^{2}=|\boldsymbol{x}-\boldsymbol{x}_{i}|^{2}$ , and the matrix $\unicode[STIX]{x1D64D}_{{\rm\pi}/2}$ refers to counter-clockwise rotation by ${\rm\pi}/2$ . Each $v_{i}$ defines a vortex centred at a point $\boldsymbol{x}_{i}$ with a Gaussian profile. The parameters ${\it\Gamma}_{i}$ and ${\it\delta}_{i}$ control the strength and width of each vortex respectively. In our simulation, we chose

(10.6a-c ) $$\begin{eqnarray}{\it\Gamma}_{1}=0.01,\quad {\it\delta}_{1}=0.05,\quad \boldsymbol{x}_{1}=({\rm\pi},{\rm\pi}+0.2),\end{eqnarray}$$

(10.7a-c ) $$\begin{eqnarray}{\it\Gamma}_{2}=0.12,\quad {\it\delta}_{2}=0.20,\quad \boldsymbol{x}_{2}=({\rm\pi},{\rm\pi}-0.2).\end{eqnarray}$$

The vortex described initially by $\boldsymbol{v}_{1}$ is, therefore, weaker and smaller, and hence is expected to deform significantly under the strain field created by the second vortex initially described by $\boldsymbol{v}_{2}$ . The computational domain is $\boldsymbol{x}\in [0,2{\rm\pi}]\times [0,2{\rm\pi}]$ , discretized into $512\times 512$ grid points. During the simulation, the vortices stay far enough from the boundary so that the boundary effects on their evolution are negligible. A relevant Reynolds number can be computed as $Re={\it\Gamma}_{2}/{\it\nu}=6\times 10^{4}$ (cf. Kevlahan & Farge Reference Kevlahan and Farge1997). For this Reynolds number, the flow still remains fairly close to its inviscid limit, with the standard deviation of the vorticity from its initial value staying below $0.5\,\%$ along fluid trajectories.

We have run the computation from the non-dimensional time $t_{0}=0$ up to $t=30.$ Over this time scale, we could preserve the smoothness of the velocity field without additional numerical effort. We applied Algorithm 1 with the tight convexity deficiency $d_{max}=10^{-3}$ and the arclength filter $l_{min}=0.1$ . We show the rotationally coherent vortex boundaries obtained in this fashion in figure 9, along with two other closed material lines initiated around the vortex cores. The latter two circles approximate closely the Eulerian vorticity boundaries inferred from the instantaneous vorticity at time $t_{0}=0$ . They, therefore, mimic the role of the initially circular vorticity jump contours advected in the inviscid simulation of Dritschel & Waugh (Reference Dritschel and Waugh1992).

Figure 9.

(a) Rotationally coherent Lagrangian vortex boundaries (red) at time $t_{0}=0$ in the vortex interaction example. Also shown are two circles (black) approximating the Eulerian vortex boundaries inferred at time $t_{0}$ from the vorticity distribution. The $\text{LAVD}_{0}^{30}(\boldsymbol{x}_{0})$ is shown in the background for reference. (b) The advected positions of the Lagrangian vortex boundaries and the two closed material curves at time $t=30.$ The advected LAVD field $\text{LAVD}_{0}^{30}(\boldsymbol{x}(30;0,\boldsymbol{x}_{0}))$ is shown in the background (see the online supplementary movie M2 for the complete advection sequence of the vortex boundaries).

As seen in figure 9, the stronger Eulerian vortex creates a sizable rotationally coherent Lagrangian vortex in the centre. This vortex shows only tangential filamentation, as confirmed by its advected position at $t=30$ . The weaker Eulerian vortex has a substantially smaller coherent Lagrangian footprint that orbits around the larger vortex.

This weaker satellite vortex gradually reaches a maximal elongated perimeter at approximately $t=24$ , then preserves its arclength in an approximate rigid-body rotation for the remaining $20\,\%$ of the simulation time. No transverse filamentation occurs in this case either: the satellite vortex remains rotationally coherent in the sense of our Definition 1.

Therefore, as figure 9 illustrates, two unequal viscous vortices may interact strongly and still preserve their rotationally coherent material cores during a finite time interval of their interaction. Such coherent cores remain hidden in contour-dynamics studies focused on the advection of Eulerian vortex boundaries inferred from large initial vorticity gradients. Indeed, the two black material curves of figure 9, mimicking the role of vortex-bounding inviscid vorticity contours, develop substantial transverse filamentation during the same time interval.

10.5 Two-dimensional Agulhas eddies in satellite altimetry

Here, we illustrate the detection of rotationally coherent eddy boundaries in velocity data derived from satellite-observed sea-surface heights under the geostrophic approximation. In this approximation, the satellite-measured sea-surface height ${\it\eta}({\it\varphi},{\it\theta},t)$ serves as a non-canonical Hamiltonian for surface velocities in the $({\it\varphi},{\it\theta})$ longitude–latitude coordinate system.

The evolution of fluid particles satisfies

(10.8) $$\begin{eqnarray}\displaystyle & \displaystyle \dot{{\it\varphi}}({\it\varphi},{\it\theta},t)=-\frac{g}{R^{2}f({\it\theta})\cos {\it\theta}}~\partial _{{\it\theta}}{\it\eta}({\it\varphi},{\it\theta},t), & \displaystyle\end{eqnarray}$$

(10.9) $$\begin{eqnarray}\displaystyle & \displaystyle \dot{{\it\theta}}({\it\varphi},{\it\theta},t)=\frac{g}{R^{2}f({\it\theta})\cos {\it\theta}}~\partial _{{\it\varphi}}{\it\eta}({\it\varphi},{\it\theta},t), & \displaystyle\end{eqnarray}$$

where $g$ is the constant of gravity, $R$ is the mean radius of the Earth and $f({\it\theta})\equiv 2{\it\Omega}\sin {\it\theta}$ is the Coriolis parameter, with ${\it\Omega}$ denoting the Earth’s mean angular velocity. The publicly available AVISO sea-surface height data base for ${\it\eta}({\it\varphi},{\it\theta},t)$ is given at a spatial resolution of $1/4^{\circ }$ and a temporal resolution of $7$ days.

We select the computational domain in the longitudinal range $[-4^{\circ },9^{\circ }]$ and the latitudinal range $[-35^{\circ },-28^{\circ }]$ , which falls inside the region of the Agulhas leakage in the Southern Ocean. A Lagrangian analysis of coherent mesoscale eddies is particularly important in this context, as the amount of warm and salty water carried from the Indian Ocean to Atlantic Ocean has relevance for global circulation and climate (Beal et al. Reference Beal, De Ruijter, Biastoch and Zahn2011).

Recent two-dimensional Lagrangian studies of the Agulhas leakage have used geodesic LCS theory to locate perfectly non-filamenting (black-hole-type) material eddies (Beron-Vera et al. Reference Beron-Vera, Wang, Olascoaga, Goni and Haller2013; Haller & Beron-Vera Reference Haller and Beron-Vera2013). Here, we use a relaxed notion of coherence, which allows for tangential filamentation, but not for global breakaway of material from the eddy. The computational cost in the present method is substantially lower: the calculation of the deformation gradient and the search for limit cycles bounding the black-hole eddies are absent.

We consider the AVISO data set ranging from the initial time $t_{0}=11\,\text{November}\,2006$ to the final time $t_{1}=t_{0}+90~\text{days}$ . We select an initial grid of particles with step size ${\it\Delta}x_{0}=1/50^{\circ }$ . As an additional filter to eliminate LAVD and IVD maxima due to resolution coarseness, we ignore from the start LAVD and IVD maxima that are closer than the submesoscale distance $0.2^{\circ }$ (approximately 20 km). As a convexity deficiency bound, we select $d_{max}=10^{-3}$ . As an arclength threshold, we fix $l_{min}=2{\rm\pi}r$ , with $r=20$  km selected again as a lower bound on mesoscale structures reliably resolved by altimetry.

Beyond executing Algorithm 1 to extract rotationally coherent eddies in this setting, we also use this example to illustrate the predictions obtained from Theorem 2 for the attractor role of rotationally coherent Lagrangian vortex centres. Figure 10(a) shows the rotationally coherent eddies, while figure 10(b) confirms that their boundaries only develop tangential filamentation under Lagrangian advection, as expected. This second plot also confirms that LAVD-based vortex centres are precisely the observed attractors for light particles released in cyclonic eddies and for heavy particles released in anticyclonic eddies.

Figure 10.

(a) Rotationally coherent Lagrangian vortices at time $t_{0}=11~\text{November}~2006$ , identified from Algorithm 1 using the contours of $\text{LAVD}_{t_{0}}^{t_{0}+T}(\boldsymbol{x}_{0})$ with $T=90~\text{days}$ . The contour plot of $\text{LAVD}_{t_{0}}^{t_{0}+T}(\boldsymbol{x}_{0})$ is shown in the background for reference. (b) Initial (red) and final (black) positions of the Lagrangian vortex boundaries at time $t_{0}+T$ , along with representative inertial particle trajectories. Heavy particles (blue) converge to the centres of anticyclonic (clockwise) eddies. Light particles (green) converge to the centres of cyclonic (clockwise) eddies. Here, ${\it\delta}_{heavy}=0.99$ , ${\it\delta}_{light}=1.01$ , $r_{0}=1$  m. (See the online supplementary movie M3 for the complete advection sequence of the vortex boundaries and the inertial particles.)

We show a comparison of the initial positions of rotationally coherent Lagrangian and Eulerian vortices in figure 11(a). Three of the Eulerian eddies (green) are close to Lagrangian eddies (red), and will accordingly show some coherence under advection by the end of the observational period in figure 11(c). The remaining Eulerian eddies show major filamentation and disintegrate under advection. We note the large number of false positives for coherent eddies based on the instantaneous Eulerian prediction at the initial time, even though this prediction is frame-invariant.

Figure 11.

(a) Lagrangian (red) and Eulerian (green) rotationally coherent vortices at time $t_{0}$ in the satellite altimetry example. (b) Rotationally coherent Eulerian vortex boundaries and centres at time $t_{0}$ , with the $\text{IVD}(\boldsymbol{x},t_{0})$ field shown in the background. (c) The same objects advected passively to the final time $t_{1}=t_{0}+T$ , with $\text{IVD}(\boldsymbol{x},t_{1})$ shown in the background.

In appendix D, we also use this example to illustrate differences between LAVD-based vortex detection and two other objective Lagrangian tools for two-dimensional flows: the geodesic LCS approach of Haller & Beron-Vera (Reference Haller and Beron-Vera2013) and the ellipticity-time diagnostic of Haller (Reference Haller2001).

We find that the geodesic LCS approach generally identifies similar material vortex regions to the LAVD approach. By construction, however, geodesic vortex boundaries enclose perfectly non-filamenting vortex cores, and hence miss the rotationally coherent (but tangentially filamenting) outer annuli of Lagrangian vortices. The identification of the perfectly coherent black-hole-type cores by the geodesic LCS approach also comes with a higher computational cost and does not extend to three dimensions (cf. appendix D for details).

We also find that the ellipticity-time diagnostic highlights the general vicinity of the coherent material vortex regions, but offers no well-defined procedure for defining vortex boundaries. Focused on individual trajectory stability rather than global coherence, the ellipticity time has equally high values in some regions that do not remain coherent as a whole. Conversely, the ellipticity time indicates predominantly hyperbolic (i.e. non-vortical) behaviour in the outer tangentially filamenting parts of rotationally coherent vortices.

In summary, despite their higher computational costs, neither globally uniform-stretching material lines nor elliptic fluid trajectories are able to provide the large coherent material vortex boundaries obtained from the LAVD.

11 Three-dimensional examples

11.1 Strong Beltrami flows

Strong Beltrami flows are steady flows whose vorticity field is a constant scalar multiple of the velocity field (Majda & Bertozzi Reference Majda and Bertozzi2002), i.e.

(11.1) $$\begin{eqnarray}{\bf\omega}(\boldsymbol{x})={\it\lambda}\boldsymbol{v}(\boldsymbol{x}).\end{eqnarray}$$

In this case, formula (4.6) gives

(11.2a,b ) $$\begin{eqnarray}\text{LAVD}_{t_{0}}^{t}(\boldsymbol{x}_{0})=|{\it\lambda}|\int _{t_{0}}^{t}|\boldsymbol{v}(\boldsymbol{x}(s))-\bar{\boldsymbol{v}}|\,\text{d}s,\quad \bar{\boldsymbol{v}}:=\frac{1}{\text{vol}\,(U)}\int _{U}\boldsymbol{v}(\boldsymbol{x})\,\text{d}V.\end{eqnarray}$$

Thus, for any strong Beltrami flow with ${\it\lambda}\neq 0$ , a rotationally coherent Lagrangian vortex boundary is a locally outermost, closed and convex level surface of $\int _{t_{0}}^{t}|\boldsymbol{v}(\boldsymbol{x}(s))-\bar{\boldsymbol{v}}|\,\text{d}s,$ the trajectory-averaged deviation of the velocity field from its spatial mean. For general 3D steady flows, level sets of asymptotically Lagrangian-averaged observables have been noted to approximate ergodic components (Budišić & Mezić Reference Budišić and Mezić2012). The new result here is that tubular level sets of the finite-time average of the normed velocity deviation in strong Beltrami flows define rotationally coherent Lagrangian structures (surfaces of constant net bulk rotation) in an objective fashion.

By formula (11.2), rotationally coherent Eulerian vortices in strong Beltrami flows are composed of convex tubular level surfaces of the normed velocity deviation $\text{IVD}(\boldsymbol{x})=|\boldsymbol{v}(\boldsymbol{x})-\bar{\boldsymbol{v}}|$ .

As an example, we consider the $ABC$ flow whose velocity field is given by

(11.3) $$\begin{eqnarray}\boldsymbol{v}(\boldsymbol{x})=\left(\begin{array}{@{}c@{}}A\sin x_{3}+C\cos x_{2}\\ B\sin x_{1}+A\cos x_{3}\\ C\sin x_{2}+B\cos x_{1}\end{array}\right)\end{eqnarray}$$

over the triply periodic domain $U=[0,2{\rm\pi}]^{3}$ . A direct calculation gives $\bar{v}(\boldsymbol{x})\equiv 0$ and ${\it\lambda}=1$ . We select the parameter configuration $A=1,$ $B=\sqrt{2/3}$ and $C=\sqrt{1/3}$ . For our Lagrangian advection, we select an initial grid of $200^{3}$ evenly distributed particles, which we advect over the interval $[t_{0},t_{1}]=[0,50]$ .

Figure 12.

(a) Three representative rotationally coherent Lagrangian vortex boundaries in the ABC flow, obtained as outermost convex tubular level sets of $\text{LAVD}_{0}^{50}(\boldsymbol{x}_{0})$ . Black curves indicate trajectories launched on these boundaries. (b) Two representative $\text{IVD}(\boldsymbol{x})$ level surfaces (with blue denoting the rotationally coherent Eulerian vortex boundary) in the single Eulerian coherent vortex obtained from the IVD field for the ABC flow for the parameter values considered here.

Figure 12(a) shows three representative Lagrangian vortex boundaries obtained as outermost convex tubular level surfaces of $\text{LAVD}_{0}^{50}(\boldsymbol{x}_{0})$ . Due to the simplicity of this steady flow, relaxation of the convexity was not necessary (the convexity deficiency is zero). Also shown in figure 12 are trajectory segments launched along the Lagrangian vortex boundaries, illustrating their Lagrangian invariance. Figure 12(b) shows the single IVD-based Eulerian vortex obtained for the ABC flow. While this Eulerian structure is near one of the Lagrangian vortices, the shape of the IVD-based vortex is not representative of the true Lagrangian vortex in this example. Furthermore, no IVD-based Eulerian vortices arise near the remaining Lagrangian vortices. This might seem puzzling first, but one must remember the spatial anisotropy of the ABC flow for unequal $A,B$ and $C$ parameter values.

11.2 Three-dimensional Agulhas eddies in a data-assimilating circulation model

Here, we apply LAVD-based vortex extraction to a three-dimensional unsteady velocity field set obtained from the Southern Ocean state estimation (SOSE) model (Mazloff, Heimbach & Wunsch Reference Mazloff, Heimbach and Wunsch2010). The domain of the data set is again in the area of the Agulhas leakage in the Southern Ocean, representing a three-dimensional extension of the Agulhas eddy extraction study of § 10.5.

Our Lagrangian study covers a period of $T=120~\text{days}$ , ranging from $t_{0}=15~\text{ May}~2006$ to $t=12~\text{September}~2006$ . The selected computational domain in the South Atlantic ocean is bounded by longitudes $[11^{\circ }\text{E},16^{\circ }\text{E}]$ , latitudes $[37^{\circ }\text{S},33^{\circ }\text{S}]$ and depths $[7,2000]$  m. We compute the LAVD and IVD fields over a uniform grid of $350\times 350\times 350$ points, and identify rotationally coherent Lagrangian and Eulerian vortices using Algorithm 2. As in our two-dimensional turbulence example, we integrate the vorticity deviation norm separately (as opposed to solving the combined ODE (9.1)), using 1000 vorticity values along each trajectory, equally spaced in time. For the plane family $\mathscr{P}_{i}$ featured in Algorithm 2, we consider horizontal planes along nodes of the initial grid, starting from 28 m below the sea-surface level to eliminate noise due to boundary effects at the surface. In these planes we use the arclength threshold $l_{min}=0.1$ and the maximal convexity deficiency $d_{max}=10^{-4}.$ The spatial mean vorticity $\bar{{\bf\omega}}(t)$ , as the practically observed mean vorticity for a large enough fluid mass in the ocean, is taken to be zero.

Figure 13(a) shows the initial position of a rotationally coherent Lagrangian eddy boundary (yellow) and its centre (red), extracted as level sets of $\text{LAVD}_{t_{0}}^{t_{1}}(\boldsymbol{x}_{0})$ by Algorithm 2. Also shown is a nearby LAVD level surface outside the eddy boundary, illustrating the complexity of the near-surface mixing region enclosing the eddy. Figure 13(b) gives a full view of the Lagrangian eddy, whereas figure 13(c) shows the materially advected position of the eddy at the final time, $120$ days later. As expected, there is mild tangential filamentation in the material eddy boundary, but strictly no breakaway from the rotating water mass. Given the complexity of material mixing in the surrounding waters, this high degree of material coherence illustrates well the accuracy of LAVD-based vortex extraction.

Figure 13.

(a) Representative level surfaces of $\text{LAVD}_{t_{0}}^{t_{0}+T}(\boldsymbol{x}_{0})$ with $t_{0}=15~\text{May}~2006$ and $T=120$ days.  The yellow surface is extracted using Algorithm 2, marking the vortex boundary for a mesoscale rotationally coherent Lagrangian eddy, extending from $28$  m down to $646$  m in depth. The green surface is a nearby level surface of $\text{LAVD}_{t_{0}}^{t_{0}+T}(\boldsymbol{x}_{0})$ outside the Lagrangian vortex region. The red curve marks the coherent vortex centre, as defined in Definition 1. (b) Full view of the Lagrangian eddy boundary and its centre at the initial time $t_{0}$ . (c) The advected eddy boundary and vortex centre $120$ days later, extending from 49 m down to 726 m below the surface (see the online supplementary movie M4 for the complete advection sequence of the eddy boundary).

The corresponding rotationally coherent Eulerian eddy extracted at the same initial location, then materially advected for 120 days, is shown in figure 14. In the Eulerian computation, the noise in the level-surface computation is more moderate than in the Lagrangian case. As a result, the plane family $\mathscr{P}_{i}$ featured in Algorithm 2 can be selected to start from as high as 7 m below the surface.

This Eulerian eddy has approximately the same diameter near the sea surface as its Lagrangian counterpart, but maintains this diameter and reaches substantially larger depths. Its vertical size is further increased under advection, with the bottom forming a sharp tip. Overall, the advected surface shows high ribbing and filamentation. Large-scale material breakaway is absent in this example, but this cannot be guaranteed a priori for an Eulerian eddy (see our two-dimensional computations in § 10.5). A comparison of figures 13 and 14 suggests that the Lagrangian eddy forms a smooth coherent centre region that transports water without observable filamentation inside the vortical Eulerian feature.

Figure 14.

(a) The rotationally coherent Eulerian counterpart of the Lagrangian eddy shown in figure 13(a), ranging from 7 m down to 1007 m below the sea surface. (b) The materially advected position of the Eulerian eddy 120 days later, ranging from 15 to 1163 m below the sea surface (see the online supplementary movie M5 for the complete advection sequence of the eddy boundary).

With a larger number of dedicated control points, the computation of LAVD level surfaces can be extended closer to the ocean surface. Figure 15 shows a higher-resolution computation of the initial and advected positions of the top slices of IVD-based and LAVD-based rotationally coherent vortex boundary surfaces. In this computation, the initial top slice of the LAVD-based vortex boundary is located only 15 m below the ocean surface, as opposed to the 28 m distance used in figure 13(c). Even this Lagrangian boundary slice remains more coherent than the Eulerian one, but the Eulerian slice still performs well under material advection for this eddy and for this advection time.

Figure 15.

Initial top slices (dashed lines, at time $t_{0}=15~\text{May}~2006$ ) and their advected positions (solid lines, $T=120~\text{days}$ later) of the LAVD- and IVD-based rotationally coherent vortex boundaries.