# A 3-D minimum-enstrophy vortex in stratified quasi-geostrophic flows

Yan Barabinot, Jean N Reinaud, Xavier J Carton, Charly de Marez, Thomas Meunier

► To cite this version:

Yan Barabinot, Jean N Reinaud, Xavier J Carton, Charly de Marez, Thomas Meunier. A 3-D minimum-enstrophy vortex in stratified quasi-geostrophic flows. Journal of Fluid Mechanics, 2024, 986 (4), pp.R1. 10.1017/jfm.2024.336. hal-04841016

HAL Id: hal-04841016

https://univ-paris8.hal.science/hal-04841016v1

Submitted on 16 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# A 3-D minimum-enstrophy vortex in stratified quasi-geostrophic flows

Yan Barabinot<fnref n="1" />,†, Jean N. Reinaud<fnref n="2" />, Xavier J. Carton<fnref n="3" />, Charly de Marez<fnref n="4" /> and Thomas Meunier<fnref n="3" />

(Received 6 February 2024; revised 2 April 2024; accepted 3 April 2024)

Applying a variational analysis, a minimum-enstrophy vortex in three-dimensional (3-D) fluids with continuous stratification is found, under the quasi-geostrophic hypothesis. The buoyancy frequency is held constant. This vortex is an ideal limiting state in a flow with an enstrophy decay while energy and generalized angular momentum remain fixed. The variational method used to obtain two-dimensional (2-D) minimum-enstrophy vortices is applied here to 3-D integral quantities. The solution from the first-order variation is expanded on a basis of orthogonal spherical Bessel functions. By computing second-order variations, the solution is found to be a true minimum in enstrophy. This solution is weakly unstable when inserted in a numerical code of the quasi-geostrophic equations. After a stage of linear instability, nonlinear wave interaction leads to the reorganization of this vortex into a tripolar vortex. Further work will relate our solution with maximal entropy 3-D vortices.

**Key words:** quasi-geostrophic flows, vortex dynamics, variational methods

## 1. Introduction

### 1.1. *Brief review of minimum-enstrophy vortex and maximum-entropy theories in 2-D fluids*

In homogeneous incompressible fluids, free-decay two-dimensional (2-D) turbulence is characterized by the spontaneous appearance of vortices which grow by merging with

† Email address for correspondence: yan.barabinot@lmd.ipsl.fr

<footnote n="1">¹LMD/ENS, 24 rue Lhomond, 75231 Paris Cedex 05, France</footnote>

<footnote n="2">²School of Mathematics and Statistics, University of St Andrews, North Haugh, St Andrews KY16 9SS, UK</footnote>

<footnote n="3">³LOPS/IUEM, Rue Dumont D'Urville, 29280 Plouzane Cedex, France</footnote>

<footnote n="4">⁴Institute of Earth Sciences, University of Iceland, Sturlugata 7, 102 Reykjavik, Iceland</footnote>
their homostrophic partners, and shed filaments to conserve angular momentum during the merging events (McWilliams 1984, 1990b). This increase in vortex sizes is the physical manifestation of the upscale cascade of energy in wavenumber space, while the production of thin filaments is associated with the downscale cascade of enstrophy (variance of vertical vorticity). While kinetic energy is only weakly (or not) dissipated at large scales, enstrophy is erased by viscosity at the dissipation wavenumber that often corresponds to the numerical grid size in simulations. In the absence of forcing, dissipation, bottom topography and beta effect, the final state of such a turbulent flow is a pair of contra-rotating vortices. The search for the spatial structure of final vortex states in 2-D turbulence has motivated the original work on minimum-enstrophy vortices (hereafter MEVs) (Leith 1984a,b), and more recently, a search for exact coherent structures in 2-D turbulent flows (Zhigunov & Grigoriev 2023). In his seminal paper, Leith used variational principles to derive the velocity profile of an axisymmetric vortex in 2-D incompressible flows, with minimal enstrophy for a given energy, angular momentum or circulation. He found two different MEVs, one with finite angular momentum (MEV M), one with finite circulation (MEV C), both described by Bessel functions of the radius.

In conclusion, Leith mentioned that final states in numerical simulations of 2-D turbulence may look more like MEV C than MEV M. Furthermore, his MEVs were barotropically unstable to small-scale perturbations, leading to non-axisymmetric (tripolar) vortices with even lower enstrophy.

Further work was done by Kozlov (1994), studying the compensated MEV, that is a MEV with different vorticity distributions, *q*, in the core and in the periphery; he studied a non-axisymmetric vortex. But this MEV had a piecewise-constant vorticity distribution and not a continuously distributed one. Later, efforts to calculate tripolar MEVs in 2-D flows, with continuous vorticity distributions, and having a linear vorticity-streamfunction relation, were not successful. Furthermore, numerical experiments suggested that steady tripolar states have nonlinear $q - \psi$ relations, with $\psi$ the streamfunction. Theoretical work done with such nonlinear $q - \psi$ relations also had to deal with critical layers, leading to very complex algebra.

Approaching the final state of free-decay 2-D turbulence from a statistical mechanics point of view, Smith (1991) determined that the squared-Lorentzian vortex $\nabla^2\psi(r) = \zeta_0/[1+(r/r_0)^2]$, with $\zeta_0$ the amplitude, $r$ the radial coordinate and $r_0$ the characteristic radius, was an exact solution to the mean field theory.

Huang & Driscoll (1994) reported an experiment showing that mixing is not very well realized in the core of the vortex. Their experiment was based on an electron plasma confined by a magnetic field. For high values of the magnetic field, the time evolution of the electron density becomes equivalent to the Euler dynamics of the vorticity and a quasi-inviscid flow is produced. Huang and Driscoll studied the relaxation of an initial hollow column of electrons towards an axisymmetric meta-equilibrium state and compared their result with existing theories of 2-D almost inviscid fluid dynamics. They claimed that their equilibrium density vorticity profile exhibited close agreement with a MEV restricted to a subdomain to avoid spurious negative vorticity values but differed substantially from maximum-entropy predictions.

Further work was done using the statistical theory. Chavanis & Sommeria (1996, 1998) solved the 2-D Euler equations in the limit of the statistical theory where maximum-entropy states are equivalent to minimum-enstrophy states and proposed a convenient classification of the maximum-entropy states for bounded and unbounded domain.
A 3-D minimum-enstrophy vortex in stratified QG flows

Motivated by the experiments of Huang & Driscoll (1994), Brands et al. (1999) compared vortex structures having maximal entropy with those having minimal enstrophy. They described the minimum-enstrophy principle as an approximation of the statistical theory. In the case investigated by Huang and Driscoll, they claimed that this approximation even yielded a better agreement with experiments than the full statistical theory because the relaxation of the system toward the maximum-entropy state was not complete. However, they emphasized the problem of spurious negative oscillations and explained that this agreement cannot hold in more general cases. Nevertheless, they concluded that when relaxing the angular momentum constraint, enstrophy minimization was consistent with the physics of the system except near the maximum accessible energy. Therefore, the MEV remained a useful solution for the statistical theory of the Euler equation, except at high energy.

1.2. Motivation for the work

Following the review in § 1.1 the MEV's ability to represent axisymmetric coherent vortices in 2-D turbulence has been debated. However, two aspects have been left untouched.

First, non-axisymmetric 2-D MEVs with continuous vorticity distributions have never been studied; this is important since vortices are rarely axisymmetric. This point will deserve a future study.

Second, a generalization to the continuously stratified three-dimensional (3-D) geostrophic turbulence including vortex stretching remains to be done. In 3-D geostrophic turbulence, according to Charney's work (Charney 1971), the energy is equipartitioned horizontally and vertically (between kinetic and potential components). This can lead to isotropic vortex structures (McWilliams, Weiss & Yavneh 1991).

Following this statement, we look for isotropic vortices, minimizing potential enstrophy with given total energy and generalized momentum. This is in accordance with the cascades of geostrophic turbulence, associated with horizontal and vertical growth of the vortices and with the formation of small-scale filaments, containing a substantial amount of potential enstrophy, and which are dissipated by viscosity (Rhines 1979; McWilliams 1989, 1990a). Providing the 3-D structure of such MEVs and assessing their robustness are the aim of the present article.

The article is organized as follows. The first part is dedicated to the methodology to construct the functional to minimize, including the calculation of enstrophy and global invariants such as energy, momentum and Casimir's. In the second part, the calculus of variations and the solution are obtained and are discussed. Finally, we assess numerically the stability and evolution of such 3-D MEVs.

## 2. Methodology

Considering rotating stratified flows, we work in the framework of the continuously stratified quasi-geostrophic (QG) equations (Charney 1948; Charney & Phillips 1953). The equations require two conditions to be satisfied by the flow: a small Rossby number ($Ro = V/(f_0 R)$) and an order-one Burger number ($Bu = (NH/(f_0 R))^2 \sim 1$), where $N(Z)$ and $f_0$ are the stratification frequency and the Coriolis parameter, $H$ and $R$ are the vertical and horizontal scales of the vortex ($Z$ is the vertical coordinate) and $V$ the maximum velocity of the flow. In the absence of forcing and of dissipation, the QG model is governed
by the equation of potential vorticity conservation:

$$
\left.
\begin{aligned}
D_t q &= 0, \\
q &= (\partial_x^2 + \partial_y^2)\psi + \partial_Z \left( \frac{f_0^2}{N^2(Z)} \partial_Z \psi \right), \\
D_t &= \partial_t + u_g \partial_x + v_g \partial_y,
\end{aligned}
\right\} \qquad (2.1)
$$

with $(u_g = -\partial_y\psi, v_g = \partial_x\psi)$ the horizontal geostrophic velocities. Usually, we apply the boundary conditions $\partial_Z\psi(Z=0, -H) = 0$. In the special case of constant stratification (that we consider here) $N(Z) = N_0$, and using the constraint on $Bu$, we can rescale the vertical coordinate $Z$ via $z = N_0Z/f_0$ and the relation between the potential vorticity and the streamfunction becomes $q = (\partial_x^2 + \partial_y^2 + \partial_z^2)\psi = \nabla^2\psi$. In these coordinates, the potential vorticity is simply the Laplacian of the streamfunction.

We define a 3-D (local) radius $r^2 = x^2 + y^2 + z^2$ and finally, we scale again it by the vortex radius $s = r/R$. We also define a 2-D radius via $\rho^2 = x^2 + y^2$ scaled as $\sigma = \rho/R$. In this spherical frame of reference, $\theta$ and $\phi$ respectively refer to the latitude (positive towards the North) and the longitude (positive towards the West).

Our goal is to find a spherically symmetric vortex in the QG model, satisfying a minimal
enstrophy constraint while energy (and possibly other integral quantities) is prescribed.
Our vortex should have a finite extent. Therefore, we assume that the potential vorticity
is zero beyond a radius $R$. Thus, for $r > R$, that is for $s > 1$, we have $\nabla_s^2 \psi = 0 =$
$\frac{1}{s^2} \partial_s (s^2 \partial_s \psi)$ where $\nabla_s^2 = (1/s^2) \partial_s (s^2 \partial_s)$ is the non-dimensional Laplacian, leading to
a solution $\psi(s > 1) = K/s$, with $K$ a constant to be determined from continuity of the
velocity field at $s = 1$.

We next determine the internal structure of the vortex (for *s* < 1). We first list the
integral invariants of the unforced, non-dissipative QG equations, in order to build the
functional of streamfunction *F*(*ψ*) the extremum of which we wish to determine. Let us
define *Z* the enstrophy, *E* the energy and *L* = ∫ *ρ*²*q* dV the angular momentum with the
horizontal radius *ρ* = *Rs* cos(*θ*) as well as the associated non-dimensional quantities Ẑ, Ĉ
and Ĉ. We also define the non-dimensional streamfunction ψ̂ such that ψ = ψ̂ψ̂, with ψ
the scale of the streamfunction. We write

$$
\left.
\begin{aligned}
\hat{Z} &= \frac{R}{\Psi} Z = 2\pi \int_0^1 s^2 (\nabla_s^2 \hat{\psi})^2 \, ds, \\
\hat{E} &= \frac{E}{R\Psi^2} = 2\pi \int_0^1 s^2 |\nabla_s \hat{\psi}|^2 \, ds = 2\pi \int_0^1 s^2 (d\hat{\psi}/ds)^2 \, ds, \\
\hat{L} &= \frac{L}{R^3 \Psi} = \frac{8\pi}{3} \int_0^1 \nabla_s^2 \hat{\psi} s^4 \, ds.
\end{aligned}
\right\} (2.2)
$$

Using the radial Laplacian only is justified since $\hat{\psi}$ must not depend on the spherical angles for a centrally symmetric vortex in 3-D. Finally, there is an infinity of Casimir invariants for the unforced, non-dissipative QG equations, among which are the potential vorticity integral $C$ (and its non-dimensional form $\hat{C}$):

$$
\hat{C} = \frac{C}{R\Psi} = 4\pi \int_{0}^{1} s^{2} (\nabla_{s}^{2} \hat{\psi}) \, ds. \qquad (2.3)
$$
A 3-D minimum-enstrophy vortex in stratified QG flows

Nevertheless, it should be noted that, if the vortex is isolated (that is, if it does not advect the whole fluid in rotation), this last integral must vanish. In the following, we drop the hat for non-dimensional quantities.

Following Leith (1984b), we could have performed variations of the functional $F(\psi)$ with respect to the vortex radius $R$, but this would have led to a relation between the Lagrange multipliers, and the enstrophy that we minimize (which is therefore unknown). We refer the reader to Appendix A for more details.

Therefore, we choose the functional to be varied with respect to $\psi$:

$$F(\psi) = Z(\psi) - \lambda E(\psi) - \mu L(\psi), \quad (2.4)$$

with $C(\psi) = 0$, $\mu$ and $\lambda$ two Lagrange multipliers. Physically, and mathematically as we will see below, this leads to $K = 0$ for the external solution.

### 3. Calculus of variation and the solution

To calculate a variation of $F$ with respect to $\psi$ we need to have $\delta\psi$ to appear explicitly in the expressions of $\delta Z$, $\delta E$, $\delta L$, where $\delta$ denotes the variation of the quantity. To do so, we will perform integrals by parts. Again, remember that we work in spherically symmetric (3-D) configuration, using the scaled spherical radius $s$.

#### 3.1. First-order variation

Before we calculate these part integrals, we make explicit the condition $C = 0$. Integrating using the expression of the 3-D radial Laplacian in $s$ leads simply to $d\psi/ds(s = 1) = 0$. This condition will be used later in the calculation.

Now, we calculate $dF/d\psi\delta\psi$ by integration by parts for the enstrophy, energy and momentum. For these integrations by parts, we have the usual conditions $\delta\psi = 0$ in $s = 0, 1$. We have successively the enstrophy, energy and momentum variations:

$$\left.
\begin{aligned}
\delta Z &= 4\pi \int_0^1 \nabla_s^2 q \delta\psi s^2 \, ds, \\
-\lambda \delta E &= 4\pi\lambda \int_0^1 q \delta\psi s^2 \, ds, \\
-\mu \delta L &= -16\pi\mu \int_0^1 \delta\psi s^2 \, ds.
\end{aligned}
\right\}
\quad (3.1)$$

Therefore, cancelling the whole variation of $F$ at first order in $\delta\psi$ (to obtain an extremum of $F$) leads to

$$\nabla_s^2 q + \lambda q - 4\mu = 0. \quad (3.2)$$

This equation has a special solution satisfying both $q = 4\mu/\lambda$ and $\nabla_s^2 q = 0$ simultaneously. This special solution is written as

$$\psi_s = \psi_2 + \frac{2\mu}{3\lambda}s^2 = \psi_2 + \psi_1 s^2, \quad (3.3)$$

where $\psi_2$ is a constant, and $\psi_1 = 2\mu/(3\lambda)$. Then, we have to solve the homogeneous rest of the first-order equation $\nabla_s^2 q + \lambda q = 0$. The general solution of which is

$$q_g(s) = q_0 j_0(\sqrt{\lambda} s), \quad (3.4)$$

where $j_0(x) = \sin(x)/x = \sinc(x)$ is the zeroth-order spherical Bessel function. For a 2-D configuration, Leith (1984) found a general solution as a function of the first-order standard
Bessel function but he calculated the functional variation according to the velocity instead of the streamfunction. Considering our geometrical hypothesis, our expression is consistent with his. Finally, we have the complete solution:

$$ \psi(s) = \psi_0 j_0(\sqrt{\lambda}s) + \psi_1 s^2 + \psi_2, \qquad (3.5) $$

where we have three unknowns $\psi_0$, $\mu$, $\lambda$, since $\psi_1$ is expressed as a function of $\mu$ and $\lambda$. As the streamfunction is defined to a constant, $\psi_2$ does not count as an unknown. However, we have to make a choice as we actually have four constraints:

(i) the boundary condition: $q(1) = \nabla_s^2 \psi(s=1) = 0$,

(ii) the potential vorticity: $d\psi/ds(s=1) = 0$,

(iii) the prescribed energy: $E = E_0$,

(iv) the prescribed momentum: $L = L_0$.

Leaving for now the first condition, accounting only for the three last conditions and after some long algebra, we obtain three equations (written in the same order as the constraints).

$$ \left. \begin{aligned} \psi_1 &= \frac{2\mu}{3\lambda} = \frac{\sqrt{\lambda}}{2} j_1(\sqrt{\lambda}) \psi_0, \\ L_0 &= \frac{16\pi}{15} \sqrt{\lambda} j_1(\sqrt{\lambda}) \left( \frac{5}{\sqrt{\lambda}} \frac{j_2(\sqrt{\lambda})}{j_1(\sqrt{\lambda})} - 1 \right) \psi_0, \\ E_0 &= \frac{7\pi}{5} \lambda j_1^2(\sqrt{\lambda}) \left( 1 - \frac{25}{7\sqrt{\lambda}} \frac{j_2(\sqrt{\lambda})}{j_1(\sqrt{\lambda})} + \frac{5}{7} \frac{j_2^2(\sqrt{\lambda})}{j_1^2(\sqrt{\lambda})} \right) \psi_0^2, \end{aligned} \right\} \qquad (3.6) $$

where $j_1$ and $j_2$ are spherical Bessel functions of order one and two, respectively. Calculations have mostly been done using the results of Bloomfield, Face & Moss (2017) which performed some analytical derivations on spherical Bessel functions.

Now the issue is to isolate the unknowns so as to determine them, which is not that easy as the system is nonlinear. In 2-D, the variation of the functional according to the vortex radius gave a simple relation between $E_0$ and $L_0$ in order to determine $\lambda$. In all cases, the energy must be positive and by similarity to the 2-D case (where $\lambda$ was the solution of $J_2(\sqrt{\lambda}) = 0$, where $J_2$ was the second-order standard Bessel function), we can set $\lambda$ as a solution of $j_2(\sqrt{\lambda}) = 0$. Therefore, we have $\sqrt{\lambda} = \gamma_n$ for some $n = 1, 2, \dots$. It is worth noting that setting this constraint makes both conditions $q(1) = 0$ and $d\psi/ds(s=1) = 0$ equivalent. Therefore, the continuity with the external problem is also verified.

In that case, the system is reduced to

$$ \left. \begin{aligned} \psi_1 &= \frac{2\mu}{3\lambda} = \frac{\sqrt{\lambda}}{2} j_1(\sqrt{\lambda}) \psi_0, \\ L_0 &= -\frac{16\pi}{15} \sqrt{\lambda} j_1(\sqrt{\lambda}) \psi_0, \\ E_0 &= \frac{7\pi}{5} \lambda j_1^2(\sqrt{\lambda}) \psi_0^2. \end{aligned} \right\} \qquad (3.7) $$

As in the 2-D case of Leith (1984), $L_0$ and $E_0$ are thus linked, here, through the relation $L_0^2/E_0 = \frac{256}{315}\pi$; $\psi_0$ and $\mu$ are thus given by $\psi_0 = -(16/21\sqrt{\lambda}j_1(\sqrt{\lambda}))(E_0/L_0)$ and $\mu = -\frac{4}{7}\lambda(E_0/L_0)$.
A 3-D minimum-enstrophy vortex in stratified QG flows

With these considerations, using $r_h = R\sigma$, this streamfunction leads to an azimuthal velocity:

$$V(s) = \frac{\partial \psi}{\partial \sigma} = \frac{16}{21} \left(\frac{E_0}{L_0}\right) \frac{1}{j_1(\sqrt{\lambda})} (j_1(\sqrt{\lambda}s) - j_1(\sqrt{\lambda})s) \cos \theta, \quad (3.8)$$

using the fact that $ds/d\sigma = \sigma/s = \cos(\theta)$. This expression is similar to equation (34) by Leith (1984). Therefore, the velocity field is maximum for $\theta = 0$, the median plane of the vortex, and minimum for $\theta = \pm\pi/2$, at the top and bottom of the vortex.

## 3.2. Second-order variation

Following Leith (1984), the second variation analysis needs to assure that one of the stationary solutions already found is a true minimum. This approach is more easily carried out using a spectral representation. In this part, calculus will be performed using $d\psi/ds$ instead of $\psi$ as it leads to simpler expressions. Initially, the aim was to show that $dF/d\psi = 0$ and $d^2F/d\psi^2 > 0$ with $F$ previously defined. Then, we slightly change our problem to $dF/du = 0$ and $d^2F/du^2 > 0$ with $u = d\psi/ds$. The important thing is to show both conditions with the same variable. Replacing $\psi$ by $u$ in the second-order derivative is not simple but doable. However, replacing in the first-order derivative remains easy. Indeed, if $dF/d\psi = 0$ then $dF/du = (dF/d\psi)(d\psi/du) = 0$. Using the properties of $\gamma_n$, $d\psi/ds$ can be represented in terms of coefficients $A_n$ in the series

$$\frac{d\psi}{ds} = \sum_{n=0}^{\infty} A_n g_n(s), \quad (3.9)$$

where the basis functions are defined by $g_0(s) = \sqrt{5}s$ and $g_n(s) = \sqrt{2}(j_1(\gamma_n s)/j_1(\gamma_n))$ for $n = 1, 2, \dots$ so that the following orthonormality relation is satisfied:

$$\int_0^1 s^2 g_n(s) g_m(s) ds = \delta_{nm} \quad \text{for } 0 \le n, m < \infty, \quad (3.10)$$

when $j_2(\gamma_n) = 0$ for $n = 1, 2, \dots$. Then, we can derive the Laplacian of the streamfunction such that

$$\nabla_s^2 \psi = \sum_{n=0}^{\infty} A_n h_n(s), \quad (3.11)$$

where we have $h_0(s) = 3\sqrt{5}$ and $h_n(s) = 3\sqrt{2}(j_0(\gamma_n s)/j_0(\gamma_n))$ for $n = 1, 2, \dots$. Note that these functions are not orthogonal. The condition on the potential vorticity of an isolated vortex imposes a linear constraint on the coefficients $A_n$ such that $\sqrt{5}A_0 + \sqrt{2}\sum_{n=1}^{\infty} A_n = 0$. Finally, using this condition and (3.11), the angular momentum, the energy and the enstrophy are

$$\left. \begin{array}{l} L = \frac{8\pi}{3}(\sqrt{2} - \sqrt{5})A_0, \\ \\ E = 2\pi \sum_{n=1}^{\infty} A_n^2, \\ \\ Z = 2\pi \sum_{n,m=0}^{\infty} H_{mn} A_n A_m, \end{array} \right\} \quad (3.12)$$
where $H_{mn} = \int_0^1 s^2 h_m(s) h_n(s) ds$. Then, using again the results from Bloomfield et al. (2017) on integrals of spherical Bessel functions, one can find the explicit expression for $H_{mn}$ as follows: $H_{00} = 15$, $H_{n0} = H_{0n} = 3\sqrt{10}$, $H_{nn} = 6 + \gamma_n^2$ for $n = 1, 2, \dots$ and $H_{mn} = 6$ for $n, m = 1, 2, \dots$.

By substituting the previous constraint on $A_n$ coefficients into the entrophy expression, we end up with

$$Z = 2\pi \sum_{n=1}^{\infty} \gamma_n^2 A_n^2. \qquad (3.13)$$

Therefore, using (3.11) we can write $A_0$ as a function of $L$ and express the energy and entrophy. That leads to $A_0 = 3L/8\pi(\sqrt{2} - \sqrt{5})$, $E = (9L^2/32\pi(\sqrt{2} - \sqrt{5})^2)(1 + \sum_{n=1}^{\infty} B_n^2)$ and $Z = (9L^2/32\pi(\sqrt{2} - \sqrt{5})^2) \sum_{n=1}^{\infty} \gamma_n^2 B_n^2$ with $B_n = A_n/A_0$. Therefore, in spectral representation, the first-order variations of $E$ and $Z$ are straightforward and the constrained first variation equation becomes

$$\delta Z - \lambda \delta E = \frac{9L^2}{16\pi(\sqrt{2} - \sqrt{5})^2} \sum_{n=1}^{\infty} (\gamma_n^2 - \lambda)(\delta B_n)B_n = 0. \qquad (3.14)$$

This equation may be satisfied for all $\delta B_n$ which leads to $\gamma_N^2 = \lambda$, for $N = 1, 2, \dots$. In particular, the constraints on $A_n$ coefficients leads to $B_N = -\sqrt{5}/2$ and we have $B_n = 0$ for $0 < n \neq N$. We thus have $E = 63L^2/64\pi(\sqrt{2} - \sqrt{5})^2$ and $Z = (45L^2/64\pi(\sqrt{2} - \sqrt{5})^2)\gamma_N^2$.

Of these various stationary solutions, the most important one with $N = 1$ has the smallest value of $Z$ since $\gamma_1^2 < \gamma_2^2 < \gamma_3^2 < \dots$. But to assure ourselves that it is a true minimum in $Z$, we must examine the second variation. The combined second variation about the stationary solution of interest becomes for $N = 1$

$$\delta^2 Z - \lambda \delta^2 E = \frac{9L^2}{16\pi(\sqrt{2} - \sqrt{5})^2} \sum_{n=2}^{\infty} (\gamma_n^2 - \gamma_1^2)(\delta B_n)^2. \qquad (3.15)$$

Since $\gamma_1^2 < \gamma_2^2 < \gamma_3^2 < \dots$, this second variation is positive for any $\delta B_n$. The solution is thus a true minimum.

### 3.3. Solution

Finally, for $N = 1$ that minimizes the entrophy value, the azimuthal velocity and the potential vorticity field are

$$\left. \begin{aligned} V(s) &= \frac{16}{21} \left(\frac{E_0}{L_0}\right) \frac{1}{j_1(\gamma_1)} (j_1(\gamma_1 s) - j_1(\gamma_1) s) \cos \theta, \\ q(s) &= \frac{16}{7} \left(\frac{E_0}{L_0}\right) \frac{1}{j_0(\gamma_1)} (j_0(\gamma_1 s) - j_0(\gamma_1)), \end{aligned} \right\} \qquad (3.16)$$

with $\gamma_1 \approx 5.76$. Fixing arbitrarily $16E_0/(7L_0) = -1$, the relative amplitude of both quantities are shown in figure 1.
<figure><img src="image_3.png" /><figcaption>Figure 1. Potential vorticity (a), azimuthal velocity (b) for $\theta = 0$.</figcaption></figure>

4. Vortex stability – numerical approach

Here we numerically assess the stability of our MEV. We use two different numerical codes of the continuously stratified QG equations. The first one uses the contour advective semi-Lagrangian (CASL) algorithm (Dritschel, Scott & Reinaud 2005; Özugurlu, Reinaud & Dritschel 2008). The second one is a standard purely Eulerian pseudo-spectral method.

The CASL simulation uses a fine grid to represent the gridded potential vorticity (PV) of resolution 1024³ with a coarser grid of 256³ for the velocity field. The continuous PV field is discretized using $n_l = 52$ uniformly distributed PV levels associated with the PV jumps $\Delta q = (PV_{max} - PV_{min})/n_l$. At $t = 0$ there is a total of 7750 PV jumps discretizing the 3-D MEV. Time integration is performed using a fourth-order Runge-Kutta scheme with a time step set by the maximum PV. The method is essentially inviscid, with contour surgery being the only source of pseudo-diffusion. The initial conditions are unperturbed and we let the small background numerical noise be the only source of perturbation. Results are shown in figure 2. They show that the vortex is slightly unstable with a azimuthal mode $m = 2$, resulting in the formation of a tripolar vortex. The evolution on unstable 3-D QG vortices has already been put in evidence by several authors (Miyazaki, Fujiwara & Yamamoto 2003; Reinaud 2017). The resulting tripolar vortex is, however, robust and does not split. We conduct a similar numerical experiment using a purely Eulerian pseudo-spectral method to check that the instability is not the result of the piecewise-uniform nature of the discretized PV in CASL. Equations are marched in time with a semi-implicit Leapfrog scheme, with a CFL = 0.1. An Euler’s time step is performed every 50 iterations to avoid the divergence of the odd and even time steps. The grid resolution is 512³. A small Newtonian dissipation, $\nu$, is introduced with $\nu = 10^{-8}$. Since the CASL simulation indicates that the vortex initially deforms following the azimuthal mode $m = 2$, we introduce an initial perturbation whereby the vortex is slightly spheroidal at $t = 0$, with a horizontal aspect ratio of 1.01. Figure 3 shows a cross-section of the PV field in the horizontal mid-plane of the vortex at $t = 0, 10$ and 19. The results confirm that the vortex is sensitive to a mode $m = 2$ and reorganizes itself into a robust tripolar vortex.
<figure><img src="image_4.png" /><figcaption>Figure 2. (a,b) Horizontal mid-plane cross-section of the PV jumps discretizing the PV field at $t = 9.55$ and 15.92 for the CASL simulation. Red contours correspond to positive PV jumps, blue ones to negative PV jumps. (c) Orthographic view on the PV jumps at the edge of the vortex at $t = 15.92$ (corresponding to the outermost PV jumps of figure 2). The vortex is viewed at an angle of 60° from the vertical direction. The colour shading indicates the height of the PV jump (the top of the vortex is in light blue, the bottom of the vortex is in dark blue).</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Figure 3. Horizontal, mid-plane cross-section of the PV field for the pseudo-spectral simulation of the perturbed 3-D MEV at (a) $t = 0$, (b) $t = 10$ and (c) $t = 19$.</figcaption></figure>

## 5. Conclusion

The mathematical calculation (with a variational method) determined that 3-D, minimal enstrophy, spherical vortices have opposite-signed QG potential vorticity radially. The internal potential vorticity distribution is continuous with that of the external fluid, which vanishes. Similarly, the rotational (azimuthal) velocity of the vortex increases quasi-linearly for a fifth of its total radius (or half of the maximum velocity radius), followed by a smoothly decaying tail. This is similar to the MEV M in 2-D incompressible fluids (Leith 1984b). This is also in qualitative agreement with detailed measurements of deep vortices at sea (Paillet et al. 1999, 2002; Ayouche et al. 2021). Recent results of very-high-resolution numerical simulations of the Atlantic Ocean (Chouksey 2023), indicate that deep/bottom vortices mostly have a shield in potential vorticity. Note that our solutions (in spherical Bessel functions) bear similarity with those recently derived by Zoeller & Viudez (2023), but from a different point of view (that of 3-D isotropic isolated vortices without the minimal-enstrophy constraint). Finally, it will be of interest to compare this MEV to 3-D maximum-entropy vortices, to be determined in the spirit of
A 3-D minimum-enstrophy vortex in stratified QG flows

the work by Brands et al. (1999), Fine et al. (1995). Further work will be needed to

(i) assess analytically if non-axisymmetric 3-D MEVs can be determined analytically (and then compared with the end states of our numerical solutions),

(ii) extend the present work to non-geostrophic dynamics (rotating shallow-water model) and if possible, to 3-D continuously stratified flows.

**Declaration of interests.** The authors report no conflict of interest.

Author ORCIDs.

● Yan Barabinot https://orcid.org/0009-0003-0628-5334;

● Jean N. Reinaud https://orcid.org/0000-0001-5449-6628;

● Xavier J. Carton https://orcid.org/0000-0002-7849-6611;

● Charly de Marez https://orcid.org/0000-0002-1712-729X;

● Thomas Meunier https://orcid.org/0000-0002-9465-0113.

Appendix A. Variation according to the vortex radius R

Following Leith (1984b), we could have performed variations of the functional F(ψ) with
respect to the vortex radius R, but this would have led to a relation between the Lagrange
multipliers, and the enstrophy that we minimize (which is therefore unknown). Indeed,
with dimensions, the variations of enstrophy δ_RZ, energy δ_RE and angular momentum δ_RL
are

$$
\left.
\begin{aligned}
\delta_R Z &= -\frac{Z}{R} \delta R, \\
\delta_R E &= \frac{E}{R} \delta R, \\
\delta_R L &= \frac{4L}{R} \delta R.
\end{aligned}
\right\}
\qquad \text{(A1)}
$$

Minimizing variations leads to

$$
Z = -\lambda E - 4\mu L. \tag{A2}
$$

Therefore, the enstrophy would have been fixed by the energy and the generalized
momentum provided, which would have led to a contradiction. In the 2-D case, the
variation according to *R* does not involve the enstrophy and helps to find a relation
between the energy and the angular momentum. Contrary to the 2-D case, the enstrophy
is computed using the whole potential vorticity *q* and thus the spherical Laplacian instead
of the vertical vorticity (*∂*<sub>*x*</sub><sup>2</sup> + *∂*<sub>*y*</sub><sup>2</sup>)ψ which gives a cylindrical Laplacian. Using the 2-D
introduced radius *ρ* and a cylindrical frame of reference would have led to a similar
development as that of Leith (1984). The 3-D geometry makes this parameter inconsistent.

REFERENCES

AYOUCHE, A., DE MAREZ, C., MORVAN, M., L'HEGARET, P., CARTON, X., LE VU, B. & STEGNER, A. 2021 Structure and dynamics of the Ras al Hadd oceanic dipole in the Arabian Sea. *Oceans* **2** (1), 105–125.

BLOOMFIELD, J.K., FACE, S.H.P. & MOSS, Z. 2017 Indefinite integrals of spherical bessel functions. arXiv:1703.06428

BRANDS, H., CHAVANIS, P.H., PASMANTER, R. & SOMMERIA, J. 1999 Maximum entropy versus minimum enstrophy vortices. *Phys. Fluids* **11** (11), 3465–3477.
CHARNEY, J.G. 1948 On the scale of atmospheric motions. In *The Atmosphere – A Challenge: The Science of Jule Gregory Charney*, pp. 251–265. Springer.

CHARNEY, J.G. 1971 Geostrophic turbulence. *J. Atmos. Sci.* **28**, 1087–1095.

CHARNEY, J.G. & PHILLIPS, N.A. 1953 Numerical integration of the quasi-geostrophic equations for barotropic and simple baroclinic flows. *J. Atmos. Sci.* **10** (2), 71–99.

CHAVANIS, P.-H. & SOMMERIA, J. 1996 Classification of self-organized vortices in two-dimensional turbulence: the case of a bounded domain. *J. Fluid Mech.* **314**, 267–297.

CHAVANIS, P.H. & SOMMERIA, J. 1998 Classification of robust isolated vortices in two-dimensional hydrodynamics. *J. Fluid Mech.* **356**, 259–296.

CHOUKSEY, A. 2023 Long-lived deep coherent vortices in the Atlantic Ocean. PhD thesis, EDSML, Université de Bretagne Occidentale.

DRITSCHEL, D.G., SCOTT, R.K. & REINAUD, J.N. 2005 The stability of quasi-geostrophic ellipsoidal vortices. *J. Fluid Mech.* **536**, 401–421.

FINE, K.S., CASS, A.C., FLYNN, W.G. & DRISCOLL, C.F. 1995 Relaxation of 2D turbulence to vortex crystals. *Phys. Rev. Lett.* **75** (18), 3277.

HUANG, X.-P. & DRISCOLL, C.F. 1994 Relaxation of 2D turbulence to a metaequilibrium near the minimum enstrophy state. *Phys. Rev. Lett.* **72** (14), 2187.

KOZLOV, V.F. 1994 On compensated eddies with minimal enstrophy. *Oceanol. Russian Acad. Sci.* **34** (1), 10–13.

LEITH, C.E. 1984a Minimum enstrophy vortex. In *AIP Conference Proceedings*, vol. 106, pp. 159–168. American Institute of Physics.

LEITH, C.E. 1984b Minimum enstrophy vortices. *Phys. fluids* **27** (6), 1388–1395.

MCWILLIAMS, J.C. 1984 The emergence of isolated coherent vortices in turbulent flow. *J. Fluid Mech.* **146**, 21–43.

MCWILLIAMS, J.C. 1989 Statistical properties of decaying geostrophic turbulence. *J. Fluid Mech.* **198**, 199–230.

MCWILLIAMS, J.C. 1990a The vortices of geostrophic turbulence. *J. Fluid Mech.* **219**, 387–404.

MCWILLIAMS, J.C. 1990b The vortices of two-dimensional turbulence. *J. Fluid Mech.* **219**, 361–385.

MCWILLIAMS, J.C., WEISS, J.B. & YAVNEH, I. 1991 The vortices of homogeneous geostrophic turbulence. *J. Fluid Mech.* **401**, 1–26.

MIYAZAKI, T., FUJIWARA, T. & YAMAMOTO, M. 2003 Quasigeostrophic confocal spheroidal vortices. *J. Phys. Soc. Japan* **72**, 2786–2803.

ÖZUGURLU, E., REINAUD, J.N. & DRITSCHEL, D.G. 2008 Interaction between two quasi-geostrophic vortices of unequal potential vorticity. *J. Fluid Mech.* **597**, 395–414.

PAILLET, J., CANN, B.L., CARTON, X., MOREL, Y. & SERPETTE, A. 2002 Dynamics and evolution of a northern meddy. *J. Phys. Oceanogr.* **32**, 55–79.

PAILLET, J., LE CANN, B., SERPETTE, A., MOREL, Y. & CARTON, X. 1999 Real-time tracking of a Galician Meddy. *Geophys. Res. Lett.* **26**, 1877–1880.

REINAUD, J.N. 2017 Piecewise uniform potential vorticity pancake shielded vortices. *Geophys. Astrophys. Fluid Dyn.* **111**, 32–64.

RHINES, P.B. 1979 Geostrophic turbulence. *Annu. Rev. Fluid Mech.* **11**, 401–441.

SMITH, R.A. 1991 Maximization of vortex entropy as an organizing principle in intermittent, decaying, two-dimensional turbulence. *Phys. Rev. A* **43**, 1126–1129.

ZHIGUNOV, D.C. & GRIGORIEV, R.O. 2023 Exact coherent structures in fully developed two-dimensional turbulence. *J. Fluid Mech.* **970**, 1–31.

ZOELLER, V.C. & VIUDEZ, A. 2023 Three-dimensional quasi-geostrophic vortices with vanishing amount of potential vorticity anomaly. *Phys. Fluids* **35**, 126607.