Eunok Yim

Paul Billant

On the mechanism of the Gent-McWilliams instability of a columnar vortex in stratified rotating fluids

. This bending instability, called "Gent-McWilliams instability" herein, is distinct from the shear, centrifugal, or radiative instabilities. The goals of the paper are to better understand the origin and properties of this instability and to explain why it operates only in stratified rotating fluids. Both numerical and asymptotic stability analyses of several velocity profiles have been performed for wide ranges of Froude number F h = Ω 0 /N and Rossby number Ro = 2Ω 0 /f , where Ω 0 is the angular velocity on the vortex axis, N the Brunt-Väisälä frequency and f the Coriolis parameter. Numerical analyses restricted to the centrifugally stable range show that the maximum growth rate of the Gent-McWilliams instability increases with Ro and is independent of F h for F h 1. In contrast, when F h > 1, the maximum growth rate decreases dramatically with F h . Long axial wavelength asymptotic analyses for isolated vortices prove that the Gent-McWilliams instability is due to the destabilization of the long-wavelength bending mode by a critical layer at the radius r c where the angular velocity Ω is equal to the frequency ω: Ω(r c ) = ω. A necessary and sufficient instability condition valid for long wavelengths, finite Rossby number and F h 1 is that the derivative of the vertical vorticity of the basic vortex is positive at r c : ζ (r c ) > 0. Such critical layer r c exists for finite Rossby and Froude numbers because the real part of the frequency of the long-wavelength bending mode is positive instead of being negative as in homogeneous non-rotating fluid (Ro = F h = ∞). When F h > 1, the instability condition ζ (r c ) > 0 is necessary but not sufficient because the destabilizing effect of the critical layer r c is strongly reduced by a second stabilizing critical layer r c2 existing at the radius where the angular velocity is equal to the Brunt-Väisälä frequency. For non-isolated vortices, numerical results show that only finite axial wavenumbers are unstable to the Gent-McWilliams instability.

Introduction

Vortices can be subjected to many instabilities in stably stratified and rotating fluids because of their own structure or because of interactions with other vortices. If we exclude the latter possibility by restricting the problem to the idealized configuration of a single vertical columnar axisymmetric vortex, there are still several instability mechanisms that can disrupt a vortex. The shear and centrifugal instabilities are not specific to stratified-rotating fluid and are well-known since [START_REF] Rayleigh | On the stability, or instability, of certain fluid motions[END_REF][START_REF] Rayleigh | On the dynamics of revolving fluids[END_REF]. The former is a two-dimensional instability that can destabilize non-axisymmetric perturbations with azimuthal wavenumbers m 2 if the vertical vorticity derivative ζ (r) vanishes for some radius r [START_REF] Rayleigh | On the stability, or instability, of certain fluid motions[END_REF][START_REF] Howard | On the hydrodynamic and hydromagnetic stability of swirling flows[END_REF][START_REF] Carton | Barotropic and baroclinic instabilities of axisymmetric vortices in a quasigeostrophic model[END_REF]Carnevale Eunok Yim and Paul Billant & Kloosterziel 1994). This necessary condition is equivalent to the inflection point condition for parallel shear flows. On the other hand, the centrifugal instability is intrinsically three-dimensional and most unstable for short axial wavelength and axisymmetric perturbations m = 0 in inviscid fluids. The classical Rayleigh criterion [START_REF] Rayleigh | On the dynamics of revolving fluids[END_REF] for the centrifugal instability can be extended to stratified-rotating inviscid fluid [START_REF] Kloosterziel | An experimental study of unstable barotropic vortices in a rotating fluid[END_REF][START_REF] Billant | Generalized rayleigh criterion for non-axisymmetric centrifugal instabilities[END_REF] by replacing the circulation by the absolute circulation, i.e. a necessary and sufficient instability condition is that the square of the absolute circulation decreases for increasing radius, regardless of the stratification. Even if the axisymmetric mode is dominant in inviscid fluid, it should be mentioned that non-axisymmetric modes can be the most unstable in viscous fluids owing to the combined effects of stratification and viscosity [START_REF] Billant | Instabilities of a vortex pair in a stratified and rotating fluid[END_REF] or confinement in shallow layers [START_REF] Lahaye | Centrifugal, barotropic and baroclinic instabilities of isolated ageostrophic anticyclones in the two-layer rotating shallow water model and their nonlinear saturation[END_REF] which damp short axial wavelengths.

Recently, a different type of instability has been shown to occur specifically on vortices in stratified rotating fluid: the radiative instability [START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF]; [START_REF] Schecter | Damping and pumping of a vortex Rossby wave in a monotonic cyclone: Critical layer stirring versus inertia-buoyancy wave emission[END_REF]; Billant & Le Dizès (2009); Le Dizès & Billant (2009); [START_REF] Riedinger | Viscous stability properties of a lamb-oseen vortex in a stratified fluid[END_REF][START_REF] Riedinger | Radiative instability of the flow around a rotating cylinder in a stratified fluid[END_REF]; [START_REF] Park | Radiative instability of an anticyclonic vortex in a stratified rotating fluid[END_REF][START_REF] Park | Instabilities and waves on a columnar vortex in a strongly stratified and rotating fluid[END_REF]). This instability is due to a coupling between the waves sustained by the vortex and inertial-gravity waves in the surrounding fluid. The mechanism of the instability can be interpreted as a wave over-reflection at the critical radius where the azimuthal phase velocity matches the angular velocity of the base flow, a radius which exists only for non-axisymmetric waves. The radiative instability is most unstable when the fluid is strongly stratified and non-rotating. As the background rotation increases, its growth rate decreases and vanishes in the quasigeostrophic limit (strongly stratified and rapidly rotating fluid).This regime pertains to large scale atmospheric and oceanic vortices (See [START_REF] Schecter | Conditions that inhibit the spontaneous radiation of spiral inertiagravity waves from an intense mesoscale cyclone[END_REF] for a summary of the typical characteristics of various geophysical vortices).

A fourth type of instability has been evidenced by [START_REF] Gent | The instability of barotropic circular vortices[END_REF] on vortices with zero-circulation in quasi-geostrophic fluids. In contrast to the shear or centrifugal instabilities, it is most unstable for a finite axial wavenumber and occurs only for the azimuthal wavenumber m = 1. This instability bends and fragments the vortex into lenticular vortices [START_REF] Gent | The instability of barotropic circular vortices[END_REF][START_REF] Hua | The internal barotropic instability of surface-intensified eddies. part i: Generalized theory for isolated eddies[END_REF][START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF]. [START_REF] Gent | The instability of barotropic circular vortices[END_REF] called this instability "internal instability",i.e. threedimensional as opposed to the "external instability", i.e. the two-dimensional shear instability. Here, we shall call the internal instability "Gent-McWilliams instability" in order to distinguish it from the other types of three-dimensional instability. The Gent-McWilliams (GMW) instability is the most dangerous instability in quasi-geostrophic fluids for sufficiently steep vortex profiles such as for example the [START_REF] Carton | Barotropic and baroclinic instabilities of axisymmetric vortices in a quasigeostrophic model[END_REF] vortex (Gaussian angular velocity profile). Even if the instability is three-dimensional, [START_REF] Gent | The instability of barotropic circular vortices[END_REF] have also shown that the necessary conditions for the shear instability (Rayleigh and Fjortoft's criteria) also apply to the GMW instability in quasi-geostrophic fluid. [START_REF] Flierl | On the instability of geostrophic vortices[END_REF] has shown that the GMW instability also occurs on vortices with piecewise profile with two levels of non-zero uniform vorticity in quasi-geostrophic fluids. A necessary condition for instability is that (ζ i -ζ o )ζ o < 0, where ζ i and ζ o are the inner and outer vorticities. This is equivalent to the condition that the vorticity gradient changes sign between the two vorticity jumps. When the vortices are isolated (i.e. zero circulation), the GMW instability is of long-wavelength nature, i.e. its growth rate is positive as soon as the axial wavenumber is non-zero. In contrast, only finite axial wavenumbers (i.e. wavenumbers above a non-zero wavenumber cut-off) are unstable when the circulation is non-zero. [START_REF] Flierl | On the instability of geostrophic vortices[END_REF] has also provided a simple explanation of the instability mechanism. Two different processes are involved: self-induced motion when the vortex is [START_REF] Flierl | On the instability of geostrophic vortices[END_REF] for piecewise vortex profiles in quasi-geostrophic fluid. As sketched in (a), we consider a columnar vortex with two levels of constant vorticity which is perturbed by a long-wavelength bending perturbation. The solid and dotted lines show the perturbed and unperturbed vorticity contours, respectively. The vorticity is assumed positive in the core region and negative in the surrounding annulus. The vorticity contours in the horizontal cross-section A are represented in (b). The initial displacements of the inner and outer contours are chosen to be in the positive x direction and equal (∆x1 = ∆x2 > 0). Because of the self-induction, the bent vortex will then rotate in the positive sense about its unperturbed position. As indicated in (c), this leads to displacements in the positive y direction (∆y1 > 0, ∆y2 > 0) in the cross-section A. However, the two contours are not displaced by the same amount: the inner contour turns out to be more displaced than the outer one (∆y1 > ∆y2 > 0). In turn, since the vorticity is not of the same sign in the inner and outer regions, this differential displacement creates positive (⊕) and negative ( ) vorticity anomalies. As sketched in (d), this dipole of vorticity anomalies tends to translate further the vortex in the positive x direction, i.e. in the same direction as the initial displacement. The same reasoning can be carried out for the motions in the cross-section B. There is therefore a positive feedback which will make the initial bending perturbation to grow.

slightly bent along the vertical and two-dimensional advection by the vorticity anomalies generated when the two vorticity contours are shifted relative to each other. This mechanism is illustrated in figure 1 and explained in details in the caption. [START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF] have studied the stability of the Carton & McWilliams vortex beyond the quasi-geostrophic regime. Most interestingly, they have shown that the GMW instability exists only in presence of both stratification and rotation. However, this observation remains still unexplained. The Rayleigh and Fjortoft's stability criteria [START_REF] Gent | The instability of barotropic circular vortices[END_REF] indeed apply only to quasi-geostrophic fluids.

Hence, the main questions that we would like to address in this paper are: why does the GMW instability occur only in stratified-rotating fluids? and is it possible to derive an instability condition valid for arbitrary rotation and stratification? These questions are intimately linked to the preliminary question: what is the general mechanism of the GMW instability? Even if a mechanism has been highlighted by [START_REF] Flierl | On the instability of geostrophic vortices[END_REF] for piecewise vorticity profiles in quasi-geostrophic fluids, it is indeed not straightforward to extend it to arbitrary rotation and stratification and/or to continuous vorticity profiles.

Interestingly, [START_REF] Reasor | Three-dimensional alignment and corotation of weak, tc-like vortices via linear vortex rossby waves[END_REF] and [START_REF] Reasor | A new look at the problem of tropical cyclones in vertical shear flow: Vortex resiliency[END_REF] have shown that certain vortices in quasi-geostrophic fluid tend to align when they are tilted, i.e. a behaviour opposite to the GMW instability. This alignment process has been attributed to a critical layer damping of the bending mode when ζ (r c ) < 0, where r c is the critical radius where the phase velocity of the wave ω/m is equal to the angular velocity of the vortex Ω(r c ) [START_REF] Briggs | Role of landau damping in cross-field electron beams and inviscid shear flow[END_REF][START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF][START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF]. In viscous fluids, this damping comes from viscous effects that smooth the singularity while in inviscid fluids, it can be understood from the conservations of angular momentum [START_REF] Schecter | Inviscid damping of asymmetries on a two-dimensional vortex[END_REF] or wave pseudo-momentum [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF].

Even if they focus on the alignment process and vortex profiles for which ζ (r c ) < 0, [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF] and [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF] have also predicted theoreti-Eunok Yim and Paul Billant cally for profiles with a weak outer vorticity gradient that the critical layer should destabilize the bending mode when ζ (r c ) > 0 in strongly stratified and rotating fluids. Such critical layer instabilities have been first reported in the case of the two-dimensional stability of vortices [START_REF] Briggs | Role of landau damping in cross-field electron beams and inviscid shear flow[END_REF][START_REF] Dizès | Non-axisymmetric vortices in two-dimensional flows[END_REF][START_REF] Schecter | Inviscid damping of asymmetries on a two-dimensional vortex[END_REF] or parallel shear flows in shallow water [START_REF] Balmforth | Shear instability in shallow water[END_REF][START_REF] Riedinger | Critical layer and radiative instabilities in shallow-water shear flows[END_REF]. [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF] and [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF] did not make a connection with the GMW instability but subsequently [START_REF] Reasor | A new look at the problem of tropical cyclones in vertical shear flow: Vortex resiliency[END_REF], while further studying vortex alignment, have also made a linear numerical simulation for a different vortex profile which leads to an exponential growth of the vortex tilt. Interestingly, [START_REF] Reasor | A new look at the problem of tropical cyclones in vertical shear flow: Vortex resiliency[END_REF] pointed out that the critical radius is located in the region of positive vorticity gradient and they speculated on the possible link between the instability observed by [START_REF] Gent | The instability of barotropic circular vortices[END_REF] on isolated vortices and a destabilization of the bending mode by a critical layer when ζ (r c ) > 0 predicted by [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF] and [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF] for profiles with a weak outer vorticity gradient. Here, we shall prove this conjecture by means of long-wavelength stability analyses of the bending mode of columnar vortices in stratified-rotating fluids.

Such long-wavelength approach differs from the analyses of [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF] and [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF] where a formal expression for the growth rate is derived under the assumption of a weak vorticity gradient at the critical radius and under the quasi-geostrophic or asymmetric balance approximations. These derivations are based on the conservation of angular momentum or wave pseudo-momentum and are equivalent to the perturbative approach of [START_REF] Briggs | Role of landau damping in cross-field electron beams and inviscid shear flow[END_REF]. The resulting formula shows that the growth rate is of the same sign as the vorticity gradient at the critical radius if the vorticity gradient is negative throughout the vortex core. However, this formula is implicit since it requires prior knowledge of the eigenfunction and frequency of the mode in the vortex core. These can be computed analytically only for the Rankine vortex. Furthermore, this formula cannot predict formally the sign of the growth rate for the Carton & McWilliams vortex since the vorticity gradient changes sign in the vortex core. It may be further noted that small vorticity gradients are encountered for large radius for many vortex profiles, meaning that the assumption of a weak vorticity gradient at the critical radius implies that the critical radius should be large. This situation occurs for small frequency, i.e. long-wavelength, in the case of the bending mode. Here, we shall use the long-wavelength assumption from the outset and this will enable us to derive explicit analytical formula for the growth rate of the bending mode for several vortex profiles including the Carton & McWilliams vortex for any finite Froude and Rossby numbers, i.e. beyond the quasi-geostrophic or asymmetric balance regimes.

These asymptotic analyses will also allow us to understand why the bending mode can be unstable only in stratified-rotating fluid. The explanation can be anticipated from what is known on the behaviours of the bending mode for a vortex with a non-zero circulation. For such profiles, long-wavelength asymptotic analyses have shown that the bending mode is neutral with a negative frequency in homogeneous non-rotating fluid [START_REF] Widnall | Theoretical and experimental study of the stability of a vortex pair[END_REF][START_REF] Saffman | Vortex Dynamics[END_REF] whereas it is neutral with a positive frequency in stratified rotating fluids when the Froude number is lower than unity: F h = Ω 0 /N < 1 [START_REF] Billant | Zigzag instability of vortex pairs in stratified and rotating fluids. part 1. general stability equations[END_REF], where Ω 0 is the angular velocity on the vortex axis and N the Brunt-Väisälä frequency. In the latter case, a critical radius r c where ω = Ω(r c ), thus exists in contrast to homogeneous non-rotating fluids. For a vortex with non-zero circulation, such critical radius has no effect in the long-wavelength limit in stratified rotating fluids because the derivative of the vorticity gradient ζ (r c ) is negligible. In this paper, we will show that the critical radius r c destabilizes the long-wavelength bending mode when ζ (r c ) is not negligible and positive. Of course, the prerequisite will be that there exists a critical radius r c , i.e. the real part of the frequency of the long-wavelength bending mode has to be positive.

Before carrying out these asymptotic analyses, we will conduct a numerical stability analysis of the m = 1 azimuthal wavenumber for the Carton & McWilliams vortex. This analysis will allow us to extend the previous results of [START_REF] Gent | The instability of barotropic circular vortices[END_REF] and [START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF] on the GMW instability by showing the separate effects of the Rossby and Froude numbers. In particular, we will show that the growth rate of the GMW instability decreases as soon as the Froude number is larger than unity whatever the Rossby number. These numerical results will also serve as a basis to check the asymptotic results. At the end of the paper, the stability of other profiles that are not unstable in the long-wavelength limit will be also studied numerically.

The paper is organized as follows: the general stability problem is described in §2. In §3, the stability of the m = 1 azimuthal wavenumber is computed numerically for the Carton & McWilliams vortex in stratified-rotating fluids. In section 4, a long-wavelength asymptotic stability analysis is first carried for the Carton & McWilliams vortex. It is shown that the frequency and growth rate of the bending mode are always positive for long-wavelength in stratified-rotating fluid while the bending mode is neutral with a negative frequency in homogeneous non-rotating fluid. To investigate the role of the sign of ζ (r c ), we next consider an angular velocity profile with the non-dimensional form Ω a n /r 2n for large radius r, where a n and n are positive constants. It will be proved that the critical radius r c can destabilize or stabilize the long-wavelength bending mode when n > 1 and n < 1, respectively. In section 5, the latter instability condition in the long-wavelength limit is generalized to ζ (r c ) > 0 for arbitrary angular velocity profiles with a weak vorticity gradient ζ . Using the numerical results of section 3 for the Carton & McWilliams vortex, we further show that this instability condition turns out to be valid not only for small vertical wavenumber but also for any vertical wavenumber. In section 6, the stability of the m = 1 azimuthal mode is computed numerically for a family of profiles possessing regions where ζ > 0 together with a non-zero circulation. In this case, we shall see that the GMW instability occurs only for finite axial wavenumbers.

Problem formulation

We consider an axisymmetric vortex with velocity u(r, θ, ẑ) = [0, Ûθ (r), 0] in cylindrical coordinates (r, θ, ẑ) where ẑ is the vertical coordinate. The fluid is assumed incompressible, stably stratified with constant Brunt-Väisälä frequency N and rotating about the vertical axis at rate f /2 where f is the Coriolis parameter. In the following, we nondimensionalize length by the radius of the vortex R and time by 1/Ω 0 , where Ω 0 is the angular velocity on the axis. Quantities without a hat will denote non-dimensional variables.

The vortex is assumed to be perturbed by infinitesimal perturbations (denoted with tilde) of velocity [ũ r , ũθ , ũz ], pressure p, and density ρ. Since the basic flow is axisymmetric and uniform along the vertical, they are written as normal modes, [ũ r , ũθ , ũz , p, ρ] = [u r (r), u θ (r), u z (r), p(r), ρ(r)]e i(kz+mθ-ωt) + c.c.,

(2.1)

where k is axial wavenumber, m the azimuthal wavenumber and ω the frequency. We consider that k and m are positive since negative wavenumbers can be retrieved by the symmetry: ω(-k, m) = -ω * (-k, -m) = ω(k, m). Under the Boussinesq approximation and in the inviscid limit, the linearized governing equations can be reduced to a single Eunok Yim and Paul Billant

equation for ϕ = ru r d 2 ϕ dr 2 - 1 r + Q Q dϕ dr - m 2 r 2 + k 2 F 2 h φ -s 2 1 -s 2 F 2 h + m rs ζ - 2 Ro + ζ 2 r + Q Q ϕ = 0, (2.2) where Q = m 2 /r 2 -k 2 F 2 h s 2 /(1 -s 2 F 2 h ), s = mΩ -ω is the Doppler shifted frequency, φ = (2Ω + 2/Ro)(ζ + 2/Ro) is the Rayleigh discriminant, ζ = (1/r)∂(rU θ )/
∂r is the vertical vorticity and Ω = U θ /r is the non-dimensional angular velocity of the basic vortex. The Froude number F h and Rossby number Ro

F h = Ω 0 N , Ro = 2Ω 0 f , (2.3)
measure the effect of the stratification and rotation, respectively. The boundary conditions impose that the perturbation is non-singular at r = 0 and decays or corresponds to a wave propagating outward for r → ∞. This implies that ϕ ∼ r m as r → 0 for m 1.

To impose the outer boundary condition, we use the two following asymptotic solutions of (2.2) for r 1

ϕ = K m (βkr) + ωRo m(2 + ωRo) βkrK m-1 (βkr), (2.4) ϕ = H (n) m (|β|kr) - ωRo m(2 + ωRo) |β|krH (n) m-1 (|β|kr), (2.5)
with

β 2 = (4/Ro 2 -ω 2 )/(1/F 2 h -ω 2 )
, where K m is the modified Bessel function of order m of the second kind and H (n) m is the Hankel function of order m of kind n. The solution (2.4) ensures that the perturbation decays as r → ∞ except when ω i = 0 and β 2 < 0. In the latter case, we use the solution (2.5) which describes an outward propagating wave if the kind is set to n = 1 when ω r (Ro 2 /4 -F 2 h ) > 0 and n = 2 when ω r (Ro 2 /4 -F 2 h ) < 0, where ω r = Re(ω). When ω i > 0, the asymptotic solution (2.4) corresponds also to the solution which propagates outward.

Equation (2.2) has been generally solved by a shooting method. The numerical integration is started from a small radius r 1 and a large radius r 2 toward a radius r f using the asymptotic solutions of (2.2) and an initial guess for ω. This leads to two couples of values at

r f : [ϕ 1 (r f ), ϕ 1 (r f )] and [ϕ 2 (r f ), ϕ 2 (r f )]. The value of ω for which the Wron- skian ϕ 2 (r f )ϕ 1 (r f ) -ϕ 1 (r f )ϕ 2 (r f )
vanishes is then searched by an iterative scheme. The path of integration is deformed in the complex plane in order to avoid the singular radius r c where s = mΩ(r c ) -ω = 0 or s = ±1/F h . Since the vortex profiles that will be studied have Ω (r c ) < 0, the path is deformed in the upper complex plane in order that the inviscid solution be the proper limit of the viscous solution [START_REF] Lin | The theory of hydrodynamics stability[END_REF][START_REF] Dizès | Viscous critical-layer analysis of vortex normal modes[END_REF].

The stability problem in presence of viscous and diffusive effects has been also solved by a Chebyshev pseudo-spectral collocation method [START_REF] Antkowiak | Transient energy growth for the lamb-oseen vortex[END_REF]. This code has allowed us to check the results of the shooting code and to compute the structure of the eigenmodes when critical layers are present since they are smoothed by viscosity and diffusion. The Reynolds number is defined as Re = Ω 0 R 2 /ν and the Schmidt number as Sc = ν/D, where ν is the kinematic viscosity and D the molecular diffusivity of the stratifying agent.

Numerical results for the Carton & McWilliams vortex

We first study numerically the stability of the azimuthal wavenumber m = 1 for the [START_REF] Carton | Barotropic and baroclinic instabilities of axisymmetric vortices in a quasigeostrophic model[END_REF] vortex whose non-dimensional angular velocity is

Ω = e -r 2 . (3.1)
The stability of this profile has been already analysed by [START_REF] Gent | The instability of barotropic circular vortices[END_REF] in the quasi-geostrophic regime and by [START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF] for finite Ro and F h but only for some particular ratios: F h /Ro = ∞, 0.5, 0.005, 0. Here, we shall investigate separately the effects of F h and Ro.

Figure 2 shows the growth rate and frequency of the azimuthal wavenumber m = 1 for a fixed Froude number F h = 1 and several positive Rossby numbers from Ro = 1 to Ro = 10. The growth rate curves for Ro 5 exhibit a bell shape as reported by [START_REF] Gent | The instability of barotropic circular vortices[END_REF] and [START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF]. When Ro increases from Ro = 1 to Ro = 5, the most amplified wavenumber increases and the maximum growth rate slightly increases. The unstable wavenumber band also widens when Ro increases but, remarkably, it corresponds to a fixed frequency range: 0 ω r 0.135 (figure 2b). Note that the critical value 0.135 is very close to exp(-2). The reason for this will be seen in §5.

The vertical vorticity of the most unstable perturbation for Ro = 1 is shown in figure

Eunok Yim and Paul Billant 3a. It exhibits an inner dipolar structure surrounded by outer spiraling perturbations. The inner and outer perturbations have a phase shift and are localized in the regions where the basic vorticity ζ is negative and positive, respectively. Thus, the perturbations will tend to translate the positive and negative vorticity regions of the basic vortex in different directions. This is a key feature of the mechanism of the GMW instability as discussed by [START_REF] Flierl | On the instability of geostrophic vortices[END_REF] and illustrated in figure 1.

For higher Ro, as exemplified in figure 2 by Ro = 10, the growth rate is no longer maximum for a finite axial wavenumber k but increases monotonically with k (figure 2a). The vertical vorticity of the eigenmode for large k is then strongly localized at a particular radius (figure 3b). As shown by [START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF], this corresponds to the centrifugal instability which is most unstable for k → ∞. Indeed, the minimum of the Rayleigh discriminant φ = (2Ω + 2/Ro)(ζ + 2/Ro) is negative when Ro > exp(2) = 7.39.

For negative Rossby number, the centrifugal instability is also present when Ro < -1 since min(φ) is negative. This can be seen in figure 4 where the growth rate and frequency are displayed for several negative Rossby numbers still for F h = 1. The growth rate increases monotonically with k for large k when Ro < -1, whereas for Ro -1, the growth rate is maximum for a finite wavenumber k. Interestingly, there are two distinct and independent branches for Ro = -3 and Ro = -5 (figure 4): the GMW instability branch at low k and the centrifugal instability branch at large k. As seen in figure 4b, the frequencies of these two branches largely differ. In contrast, for Ro = -10 and Ro = -15, the two branches are merged into a single continuous one with similar frequencies. Remarkably, the maximum frequency of the GMW instability branch is still ω r = 0.135 independently of Ro for Ro -5 (figure 4b).

The effect of the Froude number for a fixed Rossby number Ro = 1 is now displayed in figure 5. The maximum growth rate is almost independent of the Froude number when F h 1 (figure 5a). The unstable frequency range remains also the same 0 ω r 0.135 even if the upper vertical wavenumber cut-off varies (figure 5b). However, when F h is increased above unity, the maximum growth rate decreases abruptly and the unstable frequency range shrinks. As illustrated in figure 3c, the vertical vorticity of the most unstable perturbation for F h = 2 still exhibits a dipolar structure like for F h 1 but there is an additional circular region of rapid variations at a particular radius. This corresponds to a singularity of (2.2) which occurs at the radius where s(r c2 ) = 1/F h , i.e. r c2 = 0.8 for F h = 2. Since such critical layer singularity is regularized in presence of viscous and diffusive effects, the eigenmode in figure 3c and Sc = ν/D = 1. In the next section, we will show that the growth rate is strongly reduced when F h > 1 because of the damping due to this critical layer.

The separate effects of the Rossby number for F h = 1 and the Froude number for Ro = 1 on the maximum growth rate are summarized in figure 6. In figure 6a, we see that the maximum growth rate continues to decrease when the Froude number is large but it remains positive even for the largest Froude number investigated: F h = 20. The maximum growth rate increases monotonically with Ro (figure 6b): it is approximately doubled from Ro = -1 to Ro = 7.39. Only the Rossby number range: -1 < Ro < 7.39 is shown in figure 6b since outside this range, the vortex is centrifugally unstable. The centrifugal instability is then dominant over the GMW instability (figures 2a, 4a) except close to the thresholds Ro = -1 and Ro = 7.39.

Finally, the combined effects of Ro and F h on the maximum growth rate and the most amplified wavenumber k of the GMW instability are depicted in figure 7. The effect of the Froude number is the same independently of Ro: the maximum growth rate is independent of F h for F h 1 and decreases abruptly with F h when F h > 1. The monotonic increase of max(ω i ) with Ro is also observed whatever F h . The most amplified wavenumber has been scaled by F h /|Ro| in figure 7b. With this scaling, we see that it varies little and is around kF h /|Ro| ≈ 0.25 ∼ 0.45 regardless of F h and Ro.

Eunok Yim and Paul Billant

Long-wavelength asymptotic analyses

As seen in §3, the GMW instability for the Carton & McWilliams profile (3.1) starts at k = 0, i.e. it is of long-wavelength nature. For k = 0, an exact solution of (2.2) for m = 1 which derives from the translational invariance is ϕ = r(Ω -ω) [START_REF] Michalke | On inviscid instability of certain 2-dimensional vortex-type flows[END_REF]. Using this solution at leading order, the frequency and eigensolution of the azimuthal wavenumber m = 1 can be computed for small wavenumber k by means of an asymptotic expansion. A similar asymptotic analysis has been performed in stratified rotating fluids for vortices with non-zero circulation [START_REF] Billant | Zigzag instability of vortex pairs in stratified and rotating fluids. part 1. general stability equations[END_REF]. Using the present non-dimensionalisation and definition of F h and Ro (see section 2), the frequency has been found to be

ω = 2a kβ 0 2 2 -ln k|β 0 | 2 + δ(F h , Ro) a 2 -γ e , (4.1)
up to order k 2 , where a = Γ/(2πR 2 Ω 0 ) with Γ the dimensional circulation, β 0 = 2F h /Ro and γ e = 0.5772 is the Euler constant. The constant δ is given by

δ(F h , Ro) = A(F h )Ro 2 4 + B(F h )Ro + D(F h ), (4.2) with A(F h ) = ∞ 0 ξ 3 Ω 4 1 -F 2 h Ω 2 dξ, (4.3) B(F h ) = ∞ 0 ξ 3 Ω 3 1 -F 2 h Ω 2 dξ, (4.4) D(F h ) = lim η0→∞ η0 0 ξ 3 Ω 2 1 -F 2 h Ω 2 dξ -F 3 (η 0 ), (4.5)
where

F 3 (η 0 ) = η0 ξ 3 Ω2 1 -F 2 h Ω2 dξ, (4.6)
where Ω denotes the asymptotic form of Ω for large radius. For vortices with non-zero circulation, we have Ω = a/r 2 , giving F 3 (η 0 ) a 2 ln η 0 for η 0 √ F h . The factor a appears in (4.1) because time is non-dimensionalized by 1/Ω 0 instead of 2πR 2 / Γ in [START_REF] Billant | Zigzag instability of vortex pairs in stratified and rotating fluids. part 1. general stability equations[END_REF].

The frequency (4.1) is purely real when F h 1 whatever Ro. When F h > 1, the integrands in (4.3-4.5) are singular at the radius where Ω(r c2 ) = 1/F h . This singularity is regularized in presence of viscous and diffusive effects and can be taken into account in the inviscid limit by bypassing the critical radius in the upper complex plane. The constant δ is then complex with a negative imaginary part implying that the mode m = 1 is damped. Therefore, the long-wavelength bending mode is at best neutral but never unstable when the circulation is non-zero. However, the formula (4.1) breaks down when a = 0, i.e. when the angular velocity does not behave as Ω a/r 2 for r 1. For this reason, specific asymptotic analyses have been carried out for the Carton & McWilliams profile (3.1) and for a profile of the form Ω = a n /r 2n for r 1, where a n are n constants but n is different from unity. For both profiles, two regions need to be considered in order to solve (2.2) asymptotically: an inner region where r 1/k and an outer region such that r 1. However, the analysis for each profile is different and so they are presented separately.

Carton & McWilliams vortex

A specificity of the profile (3.1) concerns the location of the critical radius r c where Ω(r c ) = ω. When ω 1 and Ω a/r 2 , the critical radius is located at r c a/ω. Hence, since (4.1) shows that ω = O(k 2 ), we have r c ∼ O(1/k) meaning that r c is located in the outer region. In contrast, the critical radius for the profile (3.1) is located at r c = √ -lnω. If we anticipate that ω O(k) for the long-wavelength bending mode of the profile (3.1), we can deduce that the critical radius is located in the inner region:

r c = O( √ -lnk)
. This asymptotic problem is solved in details in appendix A and only briefly summarized here. The inner and outer solutions are expanded in the form

ϕ = ϕ 0 + F 2 h k 2 ϕ 2 + • • • , (4.7)
with ϕ 0 = r(Ω -ω) in the inner region. Even if ω will be eventually small, it is simpler to consider it arbitrary to solve the inner problem. This is the reason why ω appears in the leading order inner solution ϕ 0 . The expansion (4.7) is written in power of kF h for convenience since k always appears multiplied by F h in (2.2). The second order inner solution ϕ 2 and the outer solution are computed in appendices A.1 and A.2. The behaviour of the inner solution for small ω and for r 1 is then determined in appendix A.3. The matching between the inner and outer solutions for 1 r 1/k and ω 1 is performed in appendix A.4. It leads to the dispersion relation:

ω 2 = δk 2 β 2 0 2 - 1 r 2 c + 1 r 4 c (1 + iπγ) + 1 r 6 c 2π 2 3 -2iπγ - 2π 2 r 8 c + O 1 r 10 c
, k 2 β 2 0 ln(kβ 0 ) , (4.8) where δ is the constant defined in (4.2) and γ = sgn(ω i ). The terms inside the square brackets in (4.8) correspond to an expansion in inverse power of r 2 c = -ln ω. Several orders have been computed since the series converges slowly for ω 1 because of the logarithm. The complex terms iπγ/r 4 c and -2iπγ/r 6 c come from the presence of the critical radius r c .

The dispersion relation (4.8) is implicit since the critical radius depends on ω. To determine whether or not there exist some solutions, it is first convenient to consider only the first two leading orders in 1/r 2 c :

ω 2 = δk 2 β 2 0 2 - 1 r 2 c + 1 r 4 c (1 + iπγ) . (4.9)
By writing ω = ρe iθ and δ = µe iα , the real and imaginary parts of (4.9) read

ρ 2 cos(2θ -α) = µk 2 β 2 0 2 (-ε ρ + ε 2 ρ ), (4.10) ρ 2 sin(2θ -α) = µk 2 β 2 0 2 ε 2 ρ (πγ -θ), (4.11)
where ε ρ = -1/ ln ρ 1. Combining (4.10) and (4.11) gives

tan(2θ -α) = -ε ρ (πγ -θ) + O(ε 2 ρ ). (4.12)
This relation is satisfied if

θ = α 2 + qπ 2 - ε ρ 4 (2πγ -α -qπ) + O(ε 2 ρ ), (4.13)
where q is an integer. Inserting (4.13) into (4.10) then gives an implicit relation for ρ = |ω|

ρ 2 (-1) q = µk 2 β 2 0 2 (-ε ρ + ε 2 ρ + O(ε 3 ρ )), (4.14)
which can be satisfied at leading order only if q is odd. Hence, the implicit dispersion relation (4.9) can be solved explicitly at leading orders by first choosing a small value for ρ. Then, the corresponding value of θ and k can be obtained directly from (4.13) and (4.14), respectively. In order to exhibit the behaviour of the solutions, we first focus on the simplest case F h 1 for which there is no critical layer r c2 where Ω(r c2 ) = 1/F h so that δ is purely real and positive, i.e. α = 0. Then, (4.13) yields two solutions: θ = π/2 -ε ρ π/4 for q = 1 and θ = -π/2 + ε ρ π/4 for q = -1. Other values of q are not relevant since -π θ π. The dispersion relation (4.9) has therefore two solutions complex conjugate of each other, one unstable and the other stable, which can be written in the form:

ω = ρe iθ = ρ ε ρ 4 π + i , (4.15) ω = ρ ε ρ 4 π -i . (4.16)
where ρ (and ε ρ = -1/ ln ρ) is related to k through (4.14). The frequency ω r is positive and O(ε ρ ) smaller than the growth rate.

We now consider the case F h > 1 for which δ is complex with a negative imaginary part, i.e. -π < α 0 because of the critical layer r c2 where Ω(r c2 ) = 1/F h . In this case, the two solutions (4.15) and (4.16) become at leading order:

ω = ρ -sin α 2 + ε ρ 4 (π -α) cos α 2 + i cos α 2 + ε ρ 4 (π -α) sin α 2 , (4.17) ω = ρ sin α 2 + ε ρ 4 (π + α) cos α 2 -i cos α 2 + ε ρ 4 (π + α) sin α 2 , (4.18)
where ρ can be still determined as a function of k by means of (4.14). The two solutions are no longer complex conjugate of each other. If we consider a fixed value of ρ, the relation (4.17), which is illustrated in figure 8, shows that the growth rate of the unstable solution decreases when α decreases from zero meaning that the instability is damped by the critical layer at r c2 . From (4.13), we can deduce that the growth rate goes to zero (i.e. θ = 0) when α -π + πε ρ (dotted line in figure 8). When α decreases further below this critical value, (4.13) has no longer a solution for q = 1, i.e. the solution (4.17) does not exist anymore because the integration path has been assumed to be on the real axis. However, the solution (4.17) could be continued by deforming the integration path in the upper complex plane above the critical radius r c . The mode (4.17) would be then non-regular and damped since the critical radius would be located between the integration path and the real axis. In contrast, the stable solution (4.18) (thin lines in figure 8) exists in the whole range -π < α < 0. The damping rate of this mode is reduced when α decreases and vanishes for α = -π. Its frequency becomes negative as soon as α < -ε ρ π/2. Figure 9 shows α as a function of F h and Ro. We see that α decreases and tends to α = -π for large Froude number for any Rossby number. This decay occurs slower when Ro is negative and moderate i.e. -5 Ro 0. Hence, for a given value of ε ρ = -1/ ln |ω|, the growth rate of the unstable mode (4.17) should vanish for large Froude number when However, it should be mentioned that the dispersion relation (4.8) (and (4.9)) is no longer valid when F h = ∞ when Ro is finite. The behaviour of the inner solution for large radius is indeed different in this case so that the matching performed to derive (4.8) breaks down. Nevertheless, it remains valid when Ro = ∞ i.e. for a stratified non-rotating fluid and this is considered in appendix B. The evolution of the two solutions of (4.8) is studied as a function of F h . In particular, it is shown that in the limit F h = ∞, there remains only one solution with a negative and purely real frequency ω as predicted by (4.17) and (4.18).

We now compare the prediction of the full long-wavelength dispersion relation (4.8) for the unstable mode to the numerical results. Figure 10a shows the asymptotic and numerical growth rates and frequencies for small Froude and Rossby numbers (F h = 0.1, Ro = 0.1) approaching the quasi-geostrophic limit. We see that they agree but only in a limited range of small wavenumbers: kF h /Ro 0.05. This relatively rapid divergence is due to the fact that the asymptotic dispersion relation (4.8) has been partly obtained as an expansion in power of the parameter 1/r 2 c = -1/ ln ω. Because of the logarithm, this parameter can be actually not very small even if ω is small so that the convergence of the asymptotics is slow. Other comparisons for finite Froude and Rossby numbers are displayed in figure 10b,c,d. The slopes of the frequency and growth rate at k = 0 predicted by the asymptotics agree with the numerical results for all the values of F h and Ro investigated. As predicted by (4.17), the growth rate remains positive for sufficiently small k even for the largest Froude number investigated F h = 5 (figure 10b). We can notice that the slope of the growth rate is smaller than the one of the frequency in figure 10b,d in contrast to figure 10a,c. As shown by (4.17) and figure 8,this occurs for α -π/2, i.e. when F h is sufficiently larger than unity.

In summary, we have found that the azimuthal wavenumber m = 1 for the Carton & McWilliams vortex (3.1) is always unstable in the long-wavelength limit for finite Ro when F h 1. Thus, even if the GMW instability is dominant only in the centrifugally stable range -1 < Ro < 7.39, it exists for long-wavelength outside this range of Rossby

F h = 0.1; (b) Ro = -1, F h = 5; (c) Ro = 0.2, F h = 1; (d) Ro = 5, F h = 2.
numbers. When F h is increased above unity, the GMW instability continues to exist in the long-wavelength limit but is strongly damped by the critical layer where Ω(r c2 ) = 1/F h . In the limit Ro = F h = ∞, the long-wavelength bending mode is neutral with a negative frequency.

Vortices with algebraic decay of the angular velocity

We now investigate an angular velocity profile which behaves like Ω a n /r 2n for large radius, where a n and n are constants. For example, this is the case of the profile

Ω = 1/(1 + r 2 ) n , (4.19)
considered by [START_REF] Gent | The instability of barotropic circular vortices[END_REF]. The interest for this class of profiles lies in the fact that the vertical vorticity derivative ζ 4a n n(n -1)/r 2n+1 for large radius can be positive or negative depending if n > 1 or n < 1, respectively. This will allow us to highlight the role of the sign of ζ at the critical radius r c . However, this asymptotic analysis can be carried out analytically only when n = 1 + ε with |ε| 1, i.e. when ζ (r c ) is small. As mentioned before, [START_REF] Briggs | Role of landau damping in cross-field electron beams and inviscid shear flow[END_REF], [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF] and [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF] also use this assumption to derive their expression of the growth/decay rate due to the critical layer. In contrast to section 4.1, the critical radius r c is located in the outer region for small wavenumber since r c (a n /ω) 1/2n O(1/k 1/n ), where it has been anticipated that ω O(k 2 ) at leading order. This asymptotic analysis is carried out in details in appendix C. It gives the frequency up to orders k 2 and ε in the form:

ω = ω (0) + εω (1) + O(ε 2 ), (4.20)
where

ω (0) = 2a n kβ 0 2 2 -ln k|β 0 | 2 + δ a 2 n -γ e , (4.21
)

ω (1) = 2a n kβ 0 2 2 1 -E -γ e - δ a 2 n + ln 2 + ln (k|β 0 |) 1 -2γ e - 2δ a 2 n + 2 ln 2 , (4.22)
and E is a constant defined by

E = lim η0→0 η 0 k ∞ 4K 2 1 (β 0 kr) r(1 -ω (0) r 2 /a n ) dr + 2 β 2 0 η 2 0 -2 ln (η 0 β 0 ) -1 + 2γ e + 2ω (0) k 2 β 2 0 a n + ln η 0 β 0 4 . (4.23)
The leading term ω (0) in (4.20) corresponds to (4.1). When ε is non-zero and F h 1, all the terms in (4.20) are real except the constant E. As seen in (4.23), this constant can be indeed complex since the integrand is singular at the critical radius where Ω(r c ) -ω (0) ≡ a n /r 2 c -ω (0) = 0. Assuming that the integration path in (4.23) is deformed in the upper complex plane since Ω (r c ) < 0, the imaginary part of E is then

E i = imag(E) = πa -1 , (4.24)
with a -1 the residue at r c ,

a -1 = -2K 2 1 (kβ 0 r c ). (4.25)
Hence, we see that E i is negative implying that the growth rate ω i = -εa n k 2 β 2 0 E i /2 is positive or negative, when ε > 0 or ε < 0, respectively. Note that there exists also a decaying mode when ε > 0 which is the complex conjugate of the unstable mode. It can be obtained by integration in the lower complex plane. This stable mode is regular since the integration path can be deformed continuously from the real axis without encountering the critical radius r c unlike the stable mode when ε < 0.

When F h > 1, the additional critical radius Ω(r c2 ) = 1/F h is present like for the Carton & McWilliams vortex ( §4.1). Its effect appears again in (4.20) through δ i the imaginary part of the constant δ. The destabilizing effect of the critical radius Ω(r c ) = ω (0) when ε > 0 is then in competition with the stabilizing effect due to r c2 . Assuming that δ i is small, the growth rate is then at leading order

ω i = a n k 2 β 2 0 2 δ i a 2 n -εE i . (4.26)
Therefore, there will be an instability in the long-wavelength limit only when

δ i a 2 n -επa -1 > 0. (4.27)
From (4.3-4.5), we can obtain

δ i = πr 3 c2 2Ω (r c2 ) Ro 2F h + 1 2 1 F 3 h . (4.28)
This shows that δ i depends only on (F h , Ro) for a given velocity profile. In contrast, the residue a -1 depends also on kr c and therefore varies with (ω (0) , k). Thus, the condition (4.27) depends on the four parameters (F h , Ro, ω (0) , k). However, if we consider a very small wavenumber k, we have from (4.21): ω (0) -a n β 2 0 k 2 ln(β 0 k/2)/2 at leading orders for finite Rossby number and Froude number. This implies that kβ 0 r c = kβ 0 a n / √ ω (0) [-ln(β 0 k/2)/2] -1/2 showing that kβ 0 r c tends to zero as k tends to zero. Hence, we can estimate the residue as

a -1 = -2K 2 1 (kβ 0 r c ) ln β 0 k 2 . (4.29)
This shows that |a -1 | tends to infinity when k vanishes. In contrast, (4.28) is independent of the wavenumber. Therefore, if ε > 0, we can deduce that for any finite Rossby and Froude numbers, the condition (4.27) will be satisfied for sufficiently small wavenumber.

In other words, there should always exist an instability in the long-wavelength limit for finite Froude and Rossby numbers like for the Carton & McWilliams vortex ( §4.1). However, when either Ro = ∞ or F h = ∞, the previous condition does not apply and it would be necessary to carry out specific analyses for these limits. In the particular limit Ro = F h = ∞, it is known that the long-wavelength bending mode when n = 1 is neutral with a negative frequency [START_REF] Widnall | Theoretical and experimental study of the stability of a vortex pair[END_REF][START_REF] Saffman | Vortex Dynamics[END_REF]) like for the Carton & McWilliams vortex (appendix B). Hence, no critical radius exists.

Figure 11 shows some comparisons between the predictions of (4.20) and numerical results for the angular velocity profile (4.19) for n = 0.9 (a,c) and n = 1.1(b, d) for two Froude numbers F h = 0.5 (a,b) and F h = 1.3 (c,d) for a fixed Rossby number Ro = 1. In each case, a good agreement is found for small wavenumber. When F h = 0.5, there is an instability in the long-wavelength limit when n = 1.1 (figure 11b) while the case n = 0.9 (figure 11a) is stable since ε is positive and negative, respectively. When F h = 1.3, the value n = 0.9 (figure 11c) is still stable. The value n = 1.1 for F h = 1.3 (figure 11d) seems also stable but if we zoom on the small wavenumber region, we see that the growth rate is actually positive as predicted by (4.27). The maximum growth rate is however very small. When F h is increased further, the maximum growth rate continues to decrease but, strictly speaking, it should remain positive for sufficiently small k when F h is finite.

General instability condition

The long-wavelength asymptotic analysis performed in section 4.2 can be generalized to angular velocity profiles differing from the form Ω = a n /r 2(1+ε) for large radius but with still a weak gradient of vorticity ζ for large radius i.e. Ω = a n /r 2 + εΩ 1 (r) so that ζ = O(ε). As shown in the appendix C.6, the growth rate ω i is then at leading order

ω i = a n k 2 β 2 0 2 δ i a 2 n - πζ (r c )K 2 1 (kβ 0 r c ) Ω (r c )
.

(5.1)

Since δ i 0, a general necessary condition for instability of the bending mode in the longwavelength limit is that there exists a critical radius r c where Ω(r c ) = ω with ζ (r c ) > 0 since Ω (r c ) -2a n /r 3 c < 0. This instability condition is sufficient when

F h 1 since δ i = 0. When F h > 1, it can be also sufficient if |δ i |/a 2
n is smaller than the destabilizing term (last term of (5.1)) like for the profiles considered in section 4.2. However, this depends on each particular vortex profile and therefore the instability condition will not be always sufficient for F h > 1.

When F h < 1, i.e. δ i = 0, the expression (5.1) closely resembles to the growth rate formula derived for weak vorticity gradient at the critical radius by [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF] under the quasigeostrophic approximation and [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF] under the asymmetric and hydrostatic balance approximations. However, the formula of [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF] and [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF] are implicit since they require prior numerical computation of the eigenfunction and frequency for zero outer vorticity gra-Eunok Yim and Paul Billant dient. These formula can become explicit only for the Rankine vortex profile with skirt (weak outer vorticity field) for which an analytic solution can be found. In this case, these formula in the limit k 1 become identical to (5.1). In contrast, (5.1) is always explicit for any vortex profiles with small vorticity gradient at large radius. The critical radius is indeed given by r c = a n /ω (0) where the frequency at leading order ω (0) is given explicitly by (4.21). However, the formula (5.1) is restricted to small wavenumbers k but, actually the formula of [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF] and [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF] are also implicitly restricted to small wavenumbers since the assumption of weak vorticity gradient at the critical radius applies when the critical radius is large, i.e. when the frequency and wavenumber are small for the bending mode. In addition, we emphasize that (5.1) is valid for any finite Froude and Rossby numbers and, therefore, is not restricted to the regime of validity of the asymmetric and hydrostatic balance approximations in contrast to the formula of [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF].

It is also worth pointing out that the instability condition ζ (r c ) > 0 is equivalent to the Rayleigh and Fjortoft necessary criteria for instability in the case of a quasigeostrophic fluid and a vortex with a monotonically decreasing angular velocity [START_REF] Gent | The instability of barotropic circular vortices[END_REF]. Indeed, the Rayleigh criterion requires that there exists an inflection point r I where ζ (r I ) = 0 while the Fjortoft criterion demands that ζ > 0 for r > r I and ζ < 0 for r < r I when Ω (r) < 0 for all r. [START_REF] Montgomery | Generalized charney-stern and fjortoft theorems for rapidly rotating vortices[END_REF] have extended these criteria to flows under the asymmetric and hydrostatic balance approximations. In contrast,the present asymptotic analysis and resulting instability condition are valid for any finite Froude and Rossby numbers However, their validity are restricted to small k and small vorticity gradient ζ (r c ) in contrast to the Rayleigh and Fjortoft criteria. ) at the location of the critical radius on the real axis rc.

In practice, we can show from the numerical results obtained for the profile (3.1)(section 3) that ζ (r c ) > 0 is a necessary condition for instability of the bending mode for any vertical wavenumber and magnitude of ζ (r c ) for finite Froude and Rossby numbers. Figure 12a shows again the frequency and growth rate for Ro = F h = 1 for the profile (3.1) (i.e. corresponding to the bold line in figure 2). The derivative of the vorticity at the critical radius ζ (r cr ) is also plotted. Note that only the real part of r c is considered since the imaginary part of r c is small because ω i is small. As seen in figure 12a, the growth rate ω i and ζ (r cr ) are positive in the same range of vertical wavenumbers. Hence, the cut-off wavenumber k c F h /|Ro| 0.8 corresponds to the limit where ζ (r c ) = 0, i.e. when the critical radius r c becomes equal to the inflection point r I = √ 2. The maximum value of the frequency corresponds therefore to ω r = Ω(r I ) = 0.135 as observed in section 3. Figure 13a, which shows the growth rate as a function of ω r /Ω(r I ), confirms that the GMW instability occurs if and only if 0 ω r /Ω(r I ) 1 for F h = 1 for the different Rossby numbers investigated: -0.5 Ro 5 for which only the GMW instability exists.

Figure 12b shows now the growth rate and frequency for a Froude number larger than unity: F h = 3. In this case, the unstable wavenumber band corresponds only to a portion of the range where ζ (r cr ) > 0. This means that the instability condition ζ (r c ) > 0 is necessary but not sufficient. As shown in figure 13b, the range of rescaled frequency ω r /Ω(r I ) for which the growth rate is positive is 0 ω r /Ω(r I ) 1 only when F h 1. When F h > 1, the upper frequency cut-off is lower than Ω(r I ) because of the presence of the other critical radius Ω(r c2 ) = 1/F h which is always stabilizing. This frequency cut-off decreases when F h increases but the asymptotic dispersion relation (4.8) has shown that the instability is never suppressed for sufficiently small frequency ω for finite F h and Ro for the profile (3.1).

Non-isolated vortices

In the previous sections, we have considered vortex profiles that turn out to be unstable to the GMW instability only when their circulation is zero, i.e. when they are isolated. However, this is not generic and it is possible to have profiles with a non-zero circulation satisfying the GMW instability condition, i.e. a positive vorticity gradient for some radius. The purpose of this section is to investigate numerically the stability of the m = 1 azimuthal wavenumber for these vortices. To this end, we consider a fam-Eunok Yim and Paul Billant ily of vortices with a non-dimensional angular velocity combining a Lamb-Oseen and a Carton-McWilliams profiles:

Ω(r) = a 1 -e -r 2 r 2 + (1 -a)e -r 2 , (6.1)
where a is proportional to the non-dimensional circulation Γ = 2πa. When a = 1, the vortex is a pure Lamb-Oseen vortex while a = 0 corresponds to a pure Carton-McWilliams vortex. The angular velocity on the vortex axis is independent of a owing to the non-dimensionalization. The corresponding vorticity ζ is plotted in figure 14a for different values of a in the range 0 a 1. As a increases, the negative minimum of vorticity is reduced and disappears for a = 1. Accordingly, the domain of existence of the centrifugal instability (i.e. φ = (2Ω + 2/Ro)(ζ + 2/Ro) < 0) varies with a: the upper critical Rossby number increases from Ro c = exp(2) for a = 0 to Ro c = ∞ for a = 1 (figure 14b). In contrast, the lower critical Rossby number remains Ro c = -1 regardless of a. As in section 3, we shall study the GMW instability only for the Rossby numbers that are centrifugally stable.

Effects of the circulation parameter a and Froude number F h

Figure 15 shows the growth rate and frequency of the azimuthal wavenumber m = 1 for F h = 1, Ro = 1 for various values of the circulation parameter a. The most striking feature is that the instability no longer starts at k = 0 when a is non-zero. The maximum growth rate and the upper wavelength cut-off also decreases when a increases.

kF h /|Ro| ωr (b) 0 0.1 0.2 0.3 kF h /|Ro| (c) 0 0.1 0.2 0.3 kF h /|Ro|
The stabilization of small wavenumbers k comes from the fact that the angular velocity now behaves like Ω a/r 2 for large r. Indeed, the critical radius r c for a given small value of ω is much larger, r c = a/ω, than for a = 0 for which r c = √ -lnω. This implies that the vorticity gradient at r c is O(e -a/ω ω -3/2 ), i.e. exponentially small instead of being O(ω(lnω) 3/2 ) for a = 0. Since the growth rate for small wavenumber k is directly proportional to ζ (r c ) (see (5.1)), it remains exponentially small. The frequency and growth rate are compared to the asymptotic frequency (4.1) and growth rate (5.1) in figure 16 for various values of a. The asymptotic frequency (4.1) and the exact frequency agree for small k and then differ as soon as the growth rate ω i begins to rise. Remarkably, the asymptotic growth rate (5.1) also predicts well the lower wavenumber cut-off k c at which the exact growth rate rises. However, it quickly diverges afterward. As seen in figure 16 and 15a, the wavenumber cut-off k c increases monotonically with a. This behaviour can be understood from the fact that the growth rate becomes non-negligible when ζ (r c ) is no longer exponentially small, i.e. approximately when r c 3 whatever a. Thus, the corresponding frequency ω c a/r 2 c increases linearly with a. Since the formula (4.1) turns out to vary weakly when a varies (see figure 16), the corresponding wavenumber k c , i.e. ω(k c ) = a/r 2 c , increases with a. Figure 17a further shows the growth rate as a function of the rescaled frequency ω r /Ω(r I ), where r I is the inflection point (ζ (r I ) = 0) determined for each value of a. This figure demonstrates that the growth rate vanishes when ω r = Ω(r I ) regardless of Eunok Yim and Paul Billant a. Therefore, the instability occurs only when ζ (r c ) > 0 as already observed for a = 0 in section 5. The growth rate is however negligible for small rescaled frequency ω r /Ω(r I ) because ζ (r c ) is then exponentially small as explained previously.

Figure 18 shows now the growth rate and frequency for a Froude number above unity: F h = 1.2, still for Ro = 1. The trends when a increases are the same as for F h = 1 (figure 15). The new feature is that the growth rate is negative for small k and then suddenly increases like for F h = 1 but the maximum growth rate becomes positive only for a 0.7. Figure 17b shows that the upper wavenumber cut-offs for a 0.7 correspond to the same value of the rescaled frequency ω r /Ω(r I ). This critical value is very close to the threshold ω r /Ω(r I ) = 1 in contrast to figure 13b where the critical rescaled frequencies ω r /Ω(r I ) were clearly smaller than unity for a = 0. This difference comes from the fact that the Froude number was larger in figure 13b, F h 1.5, than in figure 17b. Another difference is that the growth rate is negative for small wavenumber (i.e. small ω r /Ω(r I )) in figure 17b while for the Carton & McWilliams profile in figure 13b, the growth rate remains positive for small k even when F h > 1. This stabilization increases with a but also with the Froude number. For example, figure 19 displays the growth rate and frequency for a fixed value of the circulation parameter (a = 0.7), but different Froude numbers. The maximum growth rate decreases when F h is above unity and becomes negative for F h = 1.3. Thus, the instability is completely suppressed when the Froude number F h is above a critical value that depends on the circulation parameter a. ; F h = 1 ; F h = 1.2 ; F h = 1.3 .

Effect of the Rossby number

The effect of the Rossby number for a fixed F h is shown in figure 20. The Rossby number is varied between 0.8 and 50, i.e. within the centrifugally stable range for a = 0.7 (-1 < Ro < 254, see figure 14b). The maximum growth rate increases with Ro like for the Carton & McWilliams vortex (figure 2a) but the most amplified wavenumber kF h /|Ro| decreases with Ro in contrast to figure 2a.

Figure 21a shows the growth rate as a function of ω r /Ω(r I ). All the curves for Ro 10 go to zero very close to the critical value ω r /Ω(r I ) = 1 as observed before. However, the growth rate for Ro = 20 and Ro = 50 remain clearly positive beyond this cut-off. Two examples of eigenmode for Ro = 5 and Ro = 50 are depicted in figure 21b. The eigenmode for Ro = 50 exhibits oscillations outside the vortex core in contrast to the one for Ro = 5. Indeed, when the Rossby number is large, the eigensolution for large radius can be wavelike (see (2.5)) since the parameter β 2 = (4/Ro 2 -ω 2 )/(1/F h 2 -ω 2 ) is negative when |ω r | > 2/Ro (assuming ω i = 0 and F h < |Ro|/2). Such outward radiation of inertia gravity waves lead to the radiative instability as shown by [START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF], [START_REF] Schecter | Damping and pumping of a vortex Rossby wave in a monotonic cyclone: Critical layer stirring versus inertia-buoyancy wave emission[END_REF], Billant & Le Dizès (2009), Le Dizès & Billant (2009), [START_REF] Riedinger | Viscous stability properties of a lamb-oseen vortex in a stratified fluid[END_REF], [START_REF] Riedinger | Radiative instability of the flow around a rotating cylinder in a stratified fluid[END_REF][START_REF] Park | Radiative instability of an anticyclonic vortex in a stratified rotating fluid[END_REF][START_REF] Park | Instabilities and waves on a columnar vortex in a strongly stratified and rotating fluid[END_REF].

Since the maximum frequency of the GMW instability is max(ω r ) = Ω(r I ) = 0.163 for a = 0.7 independently of the Rossby and Froude numbers, the critical Rossby number at which the eigenmode become radiative for some wavenumbers can be estimated to be Ro c = 2/max(ω r ) 12. Above this Rossby number, the radiation of inertia-gravity waves provides a second source of instability in addition to the critical layer at r c . This explains why the mode m = 1 remains unstable beyond the cut-off frequency ω r = Ω(r I ) for Ro = 20 and Ro = 50 in figure 21a.

Conclusion

In this paper, we have investigated the stability of the azimuthal wavenumber m = 1 of an axisymmetric columnar vortex in a stratified and rotating fluid. As shown by [START_REF] Gent | The instability of barotropic circular vortices[END_REF] and [START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF], this azimuthal wavenumber can be unstable to an instability which bends the vortex and leads to the formation of lenticular vortices when the fluid is both stratified and rotating. This three-dimensional instability, called "Gent-McWilliams (GMW) instability" herein, can occur in the cen-Eunok Yim and Paul Billant ) and (Ro = 5, kF h /|Ro| = 0.15:

).

trifugally stable regime and in the quasi-geostrophic limit. These properties make this instability clearly different from all other known instabilities of columnar axisymmetric vortices: shear instability, centrifugal instability and radiative instability. We have carried out both numerical and asymptotical stability analyses of several vortex profiles to better understand the origin of this instability, to explain why it occurs only in stratifiedrotating fluids and to derive an instability condition valid outside the quasi-geostrophic regime.

We have first recovered and extended the numerical results of [START_REF] Gent | The instability of barotropic circular vortices[END_REF] and [START_REF] Smyth | Instability of an axisymmetric vortex in a stably stratified, rotating environment[END_REF] for the stability of the azimuthal wavenumber m = 1 for the Carton & McWilliams vortex. In the centrifugally stable regime -1 Ro 7.39, we have shown that the maximum growth rate of the GMW instability increases with Ro and is independent of F h for F h 1. In contrast, when F h > 1, the maximum growth rate abruptly decreases with F h and the eigenmode exhibits a critical layer singularity at the radius where Ω(r c2 ) -ω = 1/F h . Thereby, the maximum growth rate becomes very small for large Froude number.

These numerical results for the Carton & McWilliams vortex have been completed by an asymptotic analysis of the long-wavelength bending mode. Such mode is a generic mode which derives from the translational invariances in the two-dimensional limit.

Hence, its frequency ω vanishes for k = 0 and can be computed by an asymptotic expansion for small wavenumber k. We have found that there exist two long-wavelength bending modes with complex conjugate frequency ω. The real part of the frequency ω r and growth rate ω i of the unstable mode are always positive for finite Rossby and Froude numbers for sufficiently small wavenumber and in good agreement with those of the GMW instability computed numerically. This shows that the GMW instability exists for long-wavelength not only in the centrifugally stable regime: -1 Ro 7.39, but also outside. As the Froude number is increased above unity, the GMW instability continues to exist but the slope of the growth rate with k rapidly decreases and becomes negligible because of the damping by the critical layer r c2 where Ω(r c2 ) = 1/F h . The asymptotic results have also allowed us to study the transition from a stratified to an homogeneous fluid in the case of a non-rotating fluid (Ro = ∞). We have found that the range of wavenumbers of the unstable long-wavelength bending mode shrinks to zero as the Froude number becomes infinite because of the damping by the critical layer r c2 . In contrast, the stable long-wavelength bending mode continues to exist as F h increases. It becomes neutral with a negative frequency as F h → ∞ like for the vortices with a non-zero circulation [START_REF] Widnall | Theoretical and experimental study of the stability of a vortex pair[END_REF].

This asymptotic analysis for the Carton & McWilliams vortex shows that the GMW instability is intimately linked to the existence of a critical radius r c where Ω(r c ) = ω as conjectured by [START_REF] Reasor | A new look at the problem of tropical cyclones in vertical shear flow: Vortex resiliency[END_REF]. In order to prove the role played by the sign of the vorticity gradient ζ (r c ) at the critical radius r c , a second long-wavelength asymptotic analysis has been carried out for vortices exhibiting an algebraic decay of the angular velocity Ω ∼ a n /r 2n for large radius, with n close to unity. The sign of ζ is positive if n > 1 and conversely if n < 1. In the particular case n = 1 corresponding to a vortex with a constant circulation, the long-wavelength bending mode for finite Rossby and Froude numbers is known to have a positive frequency ω r with zero growth rate if F h 1 and a negative growth rate if F h > 1 because of the critical layer r c2 where Ω(r c2 ) = 1/F h [START_REF] Billant | Zigzag instability of vortex pairs in stratified and rotating fluids. part 1. general stability equations[END_REF]. Here, we have shown that the critical layer r c where Ω(r c ) = ω has a positive or negative contribution to the growth rate if n > 1 or n < 1, respectively. For finite Rossby and Froude numbers, this positive contribution overcomes the negative contribution of the critical layer r c2 for sufficiently small wavenumber. Thereby, the bending mode is always unstable for sufficiently long wavelength if n > 1.

For general vorticity profiles, the latter instability condition becomes ζ (r c ) > 0 in agreement with the condition derived by [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF] and [START_REF] Schecter | On the symmetrization rate of an intense geophysical vortex[END_REF] under the quasi-geostrophic and asymmetric balance approximations, respectively. This condition is also equivalent to the Rayleigh inflection point and Fjortoft necessary criteria derived by [START_REF] Gent | The instability of barotropic circular vortices[END_REF] in quasi-geostrophic fluids and by [START_REF] Montgomery | Generalized charney-stern and fjortoft theorems for rapidly rotating vortices[END_REF] under the asymmetric balance approximation. However, the necessary instability condition derived herein is not restricted to the quasigeostrophic or asymmetric balance regimes but is valid for any finite Froude and Rossby numbers. Furthermore, it is necessary and sufficient if F h 1. Nevertheless, we stress that it is restricted to small vorticity gradient ζ and has been derived in the long-wavelength limit. The numerical results for the [START_REF] Carton | Barotropic and baroclinic instabilities of axisymmetric vortices in a quasigeostrophic model[END_REF] vortex show that this necessary instability condition is in practice valid for any axial wavenumber k and magnitude of ζ for finite Froude and Rossby numbers.

Finally, we have investigated the stability of a family of vortices combining a Lamb-Oseen and Carton & McWilliams profiles. These vortices exhibit both positive vorticity gradient ζ (r) > 0 for some radius and a non-zero circulation in contrast to the two profiles studied before. We have shown that these non-isolated vortices are also unstable to the GMW instability in stratified rotating fluid but only for finite axial wavenumbers Eunok Yim and Paul Billant where δ is defined in (4.2). In order to further integrate (A 15), it is first convenient to write ω = ρe iθ and to use the change of variable x = e -η 2 /ρ. This yields

J r 1 = - 2δε 2 ρ ω 2 Ro 2 e -r 2 /ρ 1 x(1 -bx) 2 (1 -ε ρ ln x) 2 dx, ( A 16)
where b = e -iθ and ε ρ = -1/ ln ρ.

Since ω is assumed to be small, the integral (A 16) can be computed asymptotically in power of the small parameter ε ρ when ε

ρ | ln x| 1, i.e. when ρ x 1 ρ . ( A 17)
The integration range in (A 16) is therefore splitted into three intervals:

e -r 2 /ρ = x1 + x2 x1 + e -r 2 /ρ x2 , ( A 18)
where the bounds (x 1 , x 2 ) are chosen such that ρ x 2 1 and 1

x 1 1/ρ. The three corresponding integrals are denoted J 1 , J 2 and J 3 .

In the first interval, we have x 1 so that we can write

J 1 = - 2δε 2 ρ ω 2 Ro 2 x1 1 x(1 -bx) 2 (1 -ε ρ ln x) 2 dx = - 2δε 2 ρ ω 2 Ro 2 x1 1 + O(1/x) b 2 x 3 (1 -ε ρ ln x) 2 dx, (A 19) This yields J 1 = 2δ ω 2 b 2 Ro 2 2e -2/ερ Ei(2(ε -1 ρ -ln x 1 )) - ε ρ x 2 1 (1 -ε ρ ln x 1 ) + O 1 x 3 1 , ( A 20)
where Ei is the exponential integral [START_REF] Abramowitz | Handbook of Mathematical Functions[END_REF]. Similarly, we have x 1 throughout the third interval so that

J 3 = - 2δε 2 ρ ω 2 Ro 2 e -r 2 /ρ x2 1 x(1 -bx) 2 (1 -ε ρ ln x) 2 dx = - 2δε 2 ρ ω 2 Ro 2 e -r 2 /ρ x2 1 + O(x) x(1 -ε ρ ln x) 2 dx, (A 21) giving J 3 = - 2δε 2 ρ ω 2 Ro 2 1 ε 2 ρ r 2 - 1 + O(x 2 ) ε ρ (1 -ε ρ ln x 2 ) . (A 22)
Finally, in the intermediate interval, we have ε ρ | ln x| 1 leading to

J 2 = - 2δε 2 ρ ω 2 Ro 2 x2 x1 1 x(1 -bx) 2 (1 -ε ρ ln x) 2 dx = - 2δε 2 ρ ω 2 Ro 2 x2 x1 1 + 2ε ρ ln x + 3ε 2 ρ (ln x) 2 + O(ε 3 ρ ) x(1 -bx) 2 dx. ( A 23)
This gives

J 2 = - 2δε 2 ρ ω 2 Ro 2 1 1 -bx -ln 1 -bx x + ε ρ ln x ln x -2 + 2 1 -bx -2 ln(1 -bx) + 2 ln(1 -bx) -2Li 2 (bx) + ε 2 ρ ln x 3bx ln x 1 -bx + (ln x) 2 + 6 ln(1 -bx) -3 ln x ln(1 -bx) -6ε 2 ρ ((-1 + ln x)Li 2 (bx) -Li 3 (bx)) + O(ε 3 ρ ) x2 x1 , ( A 24)
where Li n is the polylogarithm function of order n [START_REF] Abramowitz | Handbook of Mathematical Functions[END_REF]. We can also expand J 1 and J 3 with ε ρ since ε ρ ln x 1 1 and ε ρ ln x 2 1. This gives

J 1 = δε 2 ρ ω 2 b 2 x 2 1 Ro 2 1 + ε ρ (1 + 2 ln x 1 ) + 3 2 ε 2 ρ 2 ln x 1 + 2(ln x 1 ) 2 + O(ε 3 ρ ) , ( A 25)
J 3 = - 2δε 2 ρ ω 2 Ro 2 1 ε 2 ρ r 2 - 1 ε ρ -ln x 2 -ε ρ (ln x 2 ) 2 -ε 2 ρ (ln x 2 ) 3 + O(ε 3 ρ ) . (A 26)
The integral J 2 can be also simplified using the fact that x 2 1 and x 1 1:

J 2 = - 2δε 2 ρ ω 2 Ro 2 1 + ln(-b) + ln x 2 + 1 2b 2 x 2 1 -ε ρ π 2 3 + 2 ln(-b) + (ln(-b)) 2 -(ln x 2 ) 2 - 1 + 2 ln x 1 2b 2 x 2 1 + ε 2 ρ (ln x 2 ) 3 + π 2 + π 2 ln(-b) + 3 ln(-b) 2 + ln(-b) 3 + 3(1 + 2 ln x 1 + 2(ln x 1 ) 2 ) 4b 2 x 2 1 (A 27)
Collecting the three integrals, we obtain

J r 1 = J 1 + J 2 + J 3 = - 2δ ω 2 Ro 2 1 r 2 -ε ρ + ε 2 ρ (1 + ln(-b)) -ε 3 ρ π 2 3 + 2 ln(-b) + (ln(-b)) 2 + ε 4 ρ π 2 + π 2 ln(-b) + 3 ln(-b) 2 + ln(-b) 3 + O(ε 5 ρ ) . (A 28)
Since the small parameter depends on the logarithm of |ω|: ε ρ = -1/ ln |ω|, the convergence of (A 28) is quite slow for small ω. This is the reason why we have computed J r 1 up to order ε 4 ρ . It should be noted that (A 28) can be simplified by using as small Eunok Yim and Paul Billant As before, the inner and outer solutions are expanded in the form

ϕ = ϕ 0 + F 2 h k 2 ϕ 2 + • • • . (C 2)
The inner solution will be determined for arbitrary n while the outer solution can be computed analytically only for n = 1 + ε with |ε| 1.

C.1. Inner region

In the inner region, the solutions ϕ 0 and ϕ 2 have the same form as in (A 2) and (A 5) except that, since ω = ω 2 k 2 F 2 h + • • • , s can be approximated by Ω and the term -rωC in the leading order solution (A 2) appears in the second order solution ϕ 2 :

ϕ 0 = CrΩ, ( C 3)
ϕ 2 = ϕ 0 r 0 η Ω(Ω + 2/Ro) 1 -F 2 h Ω 2 dη + r 0 1 η 3 Ω 2 Ĩ1 (η) + 4 Ĩ2 (η) Ro + 4 Ĩ3 (η) Ro 2 dη -Cω 2 r, ( C 4)
where

Ĩp (η) = η 0 ξ 3 Ω 5-p 1 -F 2 h Ω 2 dξ, ( C 5)
with p = 1, 2 or 3.

C.2. Outer region

Like in appendix A.2, the outer solution can be determined by defining a rescaled radius U = kr with U = O(1). However, the angular velocity Ω and vorticity ζ are now not negligible since Ω ∼ a n /r 2n for r 1. Indeed, the rescaling of (2.2) gives for such profile:

d 2 ϕ dU 2 + 1 U [1 + M ] dϕ dU -ϕ 1 U 2 + β 0 2 + O(F 2 h ω 2 , F 2 h k 4n , F 2 h k 2n ω) - 4n(1 -n) U 2 (1 -ωU 2n ank 2n ) + M U 2n-1 (k 2n -ωU 2n /a n ) 2 Ro + 2(1 -n)k 2n U 2n = 0, ( C 6)
where M = O(F 2 h ω 2 , F 2 h k 4n , F 2 h k 2n ω) and β 0 = 2F h /Ro. Hence, since ω = ω 2 k 2 F 2 h , (C 6) reduces at leading order to:

d 2 ϕ dU 2 + 1 U dϕ dU -ϕ 1 U 2 + β 2 0 - 4n(1 -n) U 2 [1 -k 2(1-n) F 2 h ω 2 U 2n /a n ] = 0. (C 7)
This equation can be solved asymptotically when n = 1 + ε with |ε| 1. To do so, we expand the outer solution and ω 2 with ε ϕ out = ϕ (1) 0

dU 2 + 1 U dϕ (1) 0 dU -ϕ (1) 0 1 U 2 + β 2 0 = 4ϕ (0) 0 U 2 (1 -F 2 h ω (0) 2 U 2 /a n ) . ( C 11)
The solution is found by reduction of order ϕ

(1)

0 = K 1 (β 0 U ) U ∞ 1 ηK 2 1 (β 0 η) η ∞ 4K 2 1 (β 0 ξ) ξ(1 -F 2 h ω (0) 2 ξ 2 /a n ) dξdη, ( C 12)
where the integration constants have been chosen such that ϕ

(1) 0 tends to zero for U → ∞.

C.3. Behaviour of the inner solution for large r

In preparation for the matching, we determine here the behaviour of the inner solution for r 1. We start by determining the behaviours of the integrals Ĩp (η) defined in (C 5). Because Ω ∼ a n /r 2n as r → ∞, we have for η (a n F h ) 1/2n Ĩ1 (η) ∼ A(F h ) + O(1/η 8n-4 ), (C 13) Ĩ2 (η) ∼ B(F h ) + O(1/η 6n-4 ), (C 14) Ĩ3 (η) ∼ D(F h ) + F 3 (η), (C 15) where A, B, D are defined in (4.3), (4.4), (4.5) and F 3 in (4.6). Using the asymptotic form of Ω for large radius Ω = a n /r 2n gives the function:

F 3 (η) = η Ω2 ξ 3 1 -F 2 h Ω2 dξ = a 2 n η 4-4n -1 4 -4n . (C 16)
The behaviour of the second order inner solution (C 4) for r max(F 1/2n h , 1) is therefore

ϕ 2 ∼ C 4 a n Ro 2 δr 2n-1 4n -2 + a n 2Ro 2 (1 -n) r 3-2n - r 2n-1 2n -1 -ω 2 r , ( C 17)
where δ is defined in (4.2). The behaviour of the inner solution for large r is then In order to match the inner and outer solutions, we need also to determine the behaviour of (C 12) for U 1, i.e. for r 1/k. To do so, we first rewrite (C 12) in the form ϕ

ϕ in ∼C a n r 2n-1 + F 2 h k 2 4 a n Ro 2 δr 2n-1 4n -2 + 1 2Ro 2 (1 -n) r 3-2n - r 2n-1 2n -1 -ω 2 r + O(F 4 h k 4 ). ( C
(1) 0 = M (U )K 1 (β 0 U ) where

M (U ) = U ∞ χ(η) ηK 2 1 (β 0 η) dη, ( C 20)
and where E is a constant defined as

χ(η) = η ∞ 4K 2 1 (β 0 U ) U (1 -F 2 h ω (0) 2 U 2 a -1 n ) dU. ( C
E = η0 ∞ 4K 2 1 (β 0 U ) U (1 -F 2 h ω (0) 2 U 2 a -1 n ) dU + 2 β 2 0 η 2 0 -2 ln(η 0 β 0 ) -1 + 2γ e + 2ω (0) 2 F 2 h β 2 0 a n + ln η 0 β 0 4 , ( C 23)
with 0 < η 0 1. Using (C 22), we can now determine the behaviour of M for small U . We obtain

M =G -2 ln(U β 0 ) + β 2 0 U 2 E 2 - F 2 h ω (0) 2 β 2 0 a n + ln(U β 0 ) -1 + 2γ e + 2F 2 h ω (0) 2 β 2 0 a n -2 ln 2 + β 2 0 U 2 (ln U β 0 ) 2 + O(U 4 ), ( C 24)
where G is a constant. This constant can be set to zero without loss of generality by a simple rescaling of the leading order solution. This implies that ϕ

(1) form: C.6. Generalisation to arbitrary profiles with a weak gradient of vorticity for large radius

ϕ (1) 0 = M K 1 (β 0 U ) = - 2 ln(U β 0 ) β 0 U + β 0 U 2 E - 2F 2 h ω (
In this section, we consider an angular velocity profile different from the form Ω a n /r 2(1+ε) but with still a weak gradient of vorticity for large radius. In other words, its angular and vorticity profiles for large radius can be written in the general form Ω = Ω 0 (r) + εΩ 1 (r), ζ = εζ 1 (r), (C 31) where Ω 0 (r) = a n /r 2 , ε 1 and (Ω 1 , ζ 1 ) are of order unity. In this case, the outer solution has the same form as found in §C.2 except that (C 12) is replaced by

ϕ (1) 0 = K 1 (β 0 U ) U ∞ 1 ηK 2 1 (β 0 η) η ∞ ζ 1 (ξ/k)K 2 1 (β 0 ξ) k(Ω 0 (ξ/k) -ω (0) )
dξdη, (C 32)

Figure 1 .

 1 Figure1. Mechanism of the GMW instability as explained by[START_REF] Flierl | On the instability of geostrophic vortices[END_REF] for piecewise vortex profiles in quasi-geostrophic fluid. As sketched in (a), we consider a columnar vortex with two levels of constant vorticity which is perturbed by a long-wavelength bending perturbation. The solid and dotted lines show the perturbed and unperturbed vorticity contours, respectively. The vorticity is assumed positive in the core region and negative in the surrounding annulus. The vorticity contours in the horizontal cross-section A are represented in (b). The initial displacements of the inner and outer contours are chosen to be in the positive x direction and equal (∆x1 = ∆x2 > 0). Because of the self-induction, the bent vortex will then rotate in the positive sense about its unperturbed position. As indicated in (c), this leads to displacements in the positive y direction (∆y1 > 0, ∆y2 > 0) in the cross-section A. However, the two contours are not displaced by the same amount: the inner contour turns out to be more displaced than the outer one (∆y1 > ∆y2 > 0). In turn, since the vorticity is not of the same sign in the inner and outer regions, this differential displacement creates positive (⊕) and negative ( ) vorticity anomalies. As sketched in (d), this dipole of vorticity anomalies tends to translate further the vortex in the positive x direction, i.e. in the same direction as the initial displacement. The same reasoning can be carried out for the motions in the cross-section B. There is therefore a positive feedback which will make the initial bending perturbation to grow.

Figure 2 .Figure 3 .

 23 Figure 2. (a) Growth rate (ωi) and (b) frequency (ωr) as function of the vertical wavenumber k for the profile (3.1) for m = 1 and F h = 1 for different Rossby numbers: Ro = 1 ; Ro = 2 ; Ro = 3 ; Ro = 5 ; Ro = 10 . The horizontal dotted line in (b) shows the cut-off frequency ωr = 0.135.

Figure 4 .

 4 Figure 4. (a) Growth rate (ωi) and (b) frequency (ωr) as function of the vertical wavenumber for the profile (3.1) for m = 1 and F h = 1 for different Rossby numbers: Ro = -0.5 ; Ro = -1 ; Ro = -3 ; Ro = -5 ; Ro = -10 ; Ro = -15 . The horizontal dotted line in (b) shows the cut-off frequency ωr = 0.135.

Figure 5 .Figure 6 .

 56 Figure 5. (a) Growth rate (ωi) and (b) frequency (ωr) as function of the vertical wavenumber k for the profile (3.1) for m = 1 and Ro = 1 for different Froude numbers: F h = 0.25 ; F h = 0.5;F h = 1 ; F h = 1.5 ; F h = 2. The horizontal dotted line in (b) shows the cut-off frequency ωr = 0.135.

Figure 7 .

 7 Figure 7. (a) Maximum growth rate and (b) scaled most amplified axial wavenumber kF h /|Ro| as a function of F h and Ro for the profile (3.1) for m = 1. The contour interval is 0.002 in (a) and 0.007 in (b)

Figure 8 .

 8 Figure8. Frequency ωr (dashed lines) and growth rate ωi (solid lines) predicted by (4.17) (thick lines) and (4.18) (thin lines) as a function of α =arg(δ) for a given value of ρ = 0.01. The dotted line shows the critical value for which the solution (4.17) disappears.

Figure 9 .

 9 Figure 9. Contours of α = arg(δ) as a function of F h and Ro for the profile (3.1). The contour interval is 0.3

Figure 10 .

 10 Figure10. Comparison (thick lines: ωi; thin lines: ωr) between numerical results (solid lines) and asymptotic results (dashed lines) for the profile (3.1) for m = 1 for several combinations of F h and Ro: (a)Ro = 0.1, F h = 0.1; (b) Ro = -1, F h = 5; (c) Ro = 0.2, F h = 1; (d) Ro = 5, F h = 2.

Figure 11 .

 11 Figure 11. Frequency ωr (thin lines) and growth rate ωi (thick lines) for m = 1 for the profile (4.19), for n = 0.9 (a,c) and n = 1.1 (b,d) for F h = 0.5 (a,b) and F h = 1.3 (c,d) for Ro = 1. The solid and dashed lines show the numerical and asymptotic results, respectively. The insert in (d) displays a close-up view for small kF h /|Ro|.

Figure 12 .

 12 Figure 12. Growth rate ( : ωi) and frequency, ( : ωr) for the profile (3.1) for m = 1 as a function of kF h /|Ro| for (a) Ro = F h = 1 and (b) Ro = 1, F h = 3. The associated vorticity gradient ζ (rc) () at the location of the critical radius on the real axis rc.

Figure 13 .Figure 14 .

 1314 Figure 13. Growth rate as a function of the rescaled frequency ωr/Ω(rI ) for m = 1 for the profile (3.1) for (a) F h = 1 fixed and Ro varying: Ro = -0.5 ; Ro = 1 ; Ro = 2 ; Ro = 3 ; Ro = 5 and (b) fixed Ro = 1 fixed and F h varying: F h = 0.5 ; F h = 1 ; F h = 1.5; F h = 2 ; F h = 3 .

Figure 15 .

 15 Figure 15. (a) Growth rate and (b) frequency as a function of kF h /|Ro| for the profile (6.1)for m = 1 for F h = 1 and Ro = 1 for increasing value of a by step of 0.1.

Figure 16 .

 16 Figure 16. Comparison (thin lines: ωr; thick lines: ωi) between the numerical results (solid lines) and the asymptotic frequency (4.1) (thin dashed lines) and growth rate (5.1) (thick dashed lines) for the profile (6.1) for m = 1 for different values of a for Ro = F h = 1: (a) a = 0.3; (b) a = 0.5; (c) a = 0.8.

Figure 17 .Figure 18 .

 1718 Figure17. The rescaled growth rate ωi as a function of the rescaled frequency ωr/Ω(rI ) for the profile (6.1) for m = 1 by step of 0.1 and for Ro = 1 for increasing values of a, (a) F h = 1 and (b) F h = 1.2.

Figure 19 .

 19 Figure 19. (a) Growth rate and (b) frequency as a function of kF h /|Ro| for the profile (6.1) for m = 1 for a = 0.7, Ro = 1 and varying F h : F h = 0.8;F h = 1 ; F h = 1.2 ; F h = 1.3 .

Figure 20 .Figure 21 .

 2021 Figure 20. Growth rate and frequency as a function of kF h /|Ro| for the profile (6.1) for m = 1, a = 0.7, F h = 1.2 and varying Ro: Ro = 50 ; Ro = 20 ; Ro = 10 ; Ro = 5 ; Ro = 2 ; Ro = 1 ; and Ro = -0.5 (a) growth rate, (b) frequency. The horizontal dotted line in (b) shows the cut-off frequency ωr = 0.163.

Figure 22 .

 22 Figure 22. Parameters A( ), B( ) and D( ) as a function of F h for the Carton & McWilliams vortex (3.1): (a) real part and (b) imaginary part.

Figure 23 .

 23 Figure 23. Frequency ωr (thin lines) and growth rate ωi (thick lines) as a function of k for the azimuthal wavenumber m = 1 for the Carton & McWilliams vortex (3.1) for (a) F h = 3, (b) F h = 4, (c) F h = 4.5 and (d) F h = ∞ for Ro = ∞. The solid lines show the numerical results while the dashed line represents the asymptotic frequency (B 1).

  Behaviour of the outer solution for small U

									18)
	If we now further assume that n = 1 + ε with |ε|	1, we obtain
	ϕ in = Ca n	1 r	-2ε	ln r r	+ k 2 F 2 h r	2δ a 2 n Ro 2 +	2 ln r -1 Ro 2 -	(0) 2 ω a n
		+ ε -	ω 2 (1) a n	4δ a 2 n Ro 2 (ln r -1) +	1 Ro 2 (2 -2 ln r)	+ O(ε 2 , k 4 F 4 h ). (C 19)
			C.4.

  The expressions of the inner and outer solutions (C 19) and (C 26) match ifThe first relation gives the constant C while the two others give the dispersion relation

													0)
													2
													β 2 0 a n
									+ ln(U β 0 ) -1 + 2γ e +	4F 2 h ω β 2 0 a n (0) 2	-2 ln 2	+ O(U 3 ).	(C 25)
	We can then deduce that the outer solution (C 9) becomes for small U = kr
	ϕ out =	1 β 0 kr	+	β 0 kr 2	ln(β 0 kr) 2	+ γ e -	1 2	+ ε -	2 ln(β 0 kr) β 0 kr
		+	β 0 kr 2	E -	2F 2 h ω β 2 0 a n (0) 2	+ ln(β 0 kr) -1 + 2γ e +	4F 2 h ω 2 (0) β 2 0 a n	-2 ln 2	+ • • • .
													(C 26)
													C.5. Matching
	Ca n =	1 β 0 k	(1 -2ε ln(kβ 0 )),		(C 27)
	2Ca n k 2 F 2 h Ro 2 1 + ε		2δ a 2 n	-1		=	β 0 k 2	+	εβ 0 k 2	-1 + 2γ e -2 ln 2 +	4F 2 h ω 2 (0) β 2 0 a n	,
													(C 28)
	Cak 2 F 2 h	2δ a 2 n Ro 2 -	1 Ro 2 -	ω 2 (0) a n	+ ε -	4δ a 2 n Ro 2 +	2 Ro 2 -	ω 2 (1) a n	=	β 0 k 2	ln	β 0 k 2
	+γ e -	1 2			+	εβ 0 k 2	E -	2F 2 h ω 2 (0) β 2 0 a n	+ ln(kβ 0 ) -1 + 2γ e +	4F 2 h ω 2 (0) β 2 0 a n	-2 ln 2	. (C 29)
	ω =F 2 h k 2 (ω	(0) 2 + εω 2 ) = (1)	a n β 2 0 k 2 2	-ln	β 0 k 2	-γ e +	δ a 2 n	+ εa n	β 2 0 k 2 2	1 -E -γ e -	δ a 2 n
	+ ln 2 + ln(β 0 k) 1 -2γ e -		2δ n a 2	+ 2 ln 2 + O(ε 2 , k 4 ln(k)).	(C 30)

Eunok Yim and Paul Billant as previously observed by [START_REF] Flierl | On the instability of geostrophic vortices[END_REF] for piecewise vorticity profiles in quasi-geostrophic fluids. The long-wavelength limit k 1 is stable (or marginally unstable) because the vorticity gradient at the critical radius ζ (r c ) is exponentially small. This result comes from the fact that the critical radius r c is very large when the circulation (Γ = 2πa) is non-zero and the frequency is small: r c (a/ω) 1/2 . Another particularity of the GMW instability for these non-isolated vortices is that the damping effect of the critical radius r c2 where Ω(r c2 ) = 1/F h can completely overcome the destabilizing effect of the critical radius r c when the Froude number is sufficiently large.

In summary, the origin of the GMW instability can be traced to a destabilization of the long-wavelength bending mode of a columnar vortex by a critical layer r c where Ω(r c ) = ω/m when ζ (r c ) > 0. This mode is specific to the azimuthal wavenumber m = 1 since it derives from the translational invariances. However, more generally, threedimensional for other azimuthal wavenumbers m 2 could be also destabilized by a critical layer if their frequency range lies in the range 0 < ω/m < Ω 0 . Some vortex waves for these azimuthal wavenumbers may indeed have a frequency in this range, both in homogeneous fluid [START_REF] Kelvin | Vibrations of a columnar vortex[END_REF][START_REF] Fabre | The Kelvin waves and the singular modes of the Lamb-Oseen vortex[END_REF]) and stratified-rotating fluid [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF][START_REF] Park | Instabilities and waves on a columnar vortex in a strongly stratified and rotating fluid[END_REF], and it would be interesting to investigate their stability for vortex profiles satisfying the condition ζ (r c ) > 0.

We thank M.T. Montgomery and W.D. Smyth for their comments to improve this manuscript.

Appendix A. Long-wavelength analysis for the Carton & McWilliams vortex in stratified rotating fluids

In this appendix, the equation (2.2) is solved asymptotically for small kF h and m = 1 for the Carton & McWilliams vortex (3.1) using the expansion (4.7). The frequency ω is first considered arbitrary for simplicity. Two regions are considered: an inner region where r 1/k and an outer region where r 1. The matching of the solutions in these two regions will lead to the dispersion relation.

A.1. Inner region

At leading order in kF h , (2.2) reduces to

The solution which is non-singular at r = 0 is

where C is a constant and we recall that s = Ω -ω. At second order, we have

where

The solution can be found by reduction of order in the form

where

and

with p = 1, 2, or 3. The inner solution is therefore

. Outer region

In the outer region, we define a rescaled radius such that U = kr with U = O(1). Since Ω and ζ are exponentially small for r = O(1/k) 1, (2.2) becomes at leading order

where

and

Assuming that β 2 > 0 and ω i = 0 or ω i = 0, the solution which satisfies the boundary condition at infinity is

where E is a constant.

A.3. Behaviour of J for large r

Before matching the inner and outer solutions, we need to determine the behaviour of the function J defined in (A 6) for r 1. To this end, we assume now that ω is small. We first determine the behaviour of the integrals I p (η) defined in (A 7) for η 1. We have

provided that 1 η 1/ √ ω and ωF h 1. The constants (A, B, D) are defined in (4.3), (4.4) and (4.5) and the function F 3 in (4.6). For the Carton & McWilliams profile (3.1), the angular velocity decays exponentially with r so that the asymptotic form of Ω for large r can be taken as Ω = 0, leading to F 3 (η 0 ) = F 3 (η) = 0. Accordingly, we have

1. The constants (A, B, D) are plotted as a function of the Froude number in figure 22. We see that they are real and positive for F h < 1 and complex with a negative imaginary part when F h > 1. As mentioned in section 4.1, this is because the integrands are singular at the critical radius r c2 where Ω(r c2 ) = 1/F h and the integration path is deformed in the upper complex plane. From (A 12 -A 14), we deduce that the behaviour of J for large r is at leading order

Eunok Yim and Paul Billant parameter ερ = -1/ ln ω, i.e. equivalently ερ = 1/r 2 c . It yields

(A 29) where γ = 1 when 0 < θ < π and γ = -1 when -π < θ < 0, if the integration is performed on the real axis. All the complex terms in (A 29) are included in the parameter ερ except two. These two terms correspond to the residue at the singular radius r c = √ -ln ω.

A.4. Matching

The inner and outer solutions should match in the overlap region 1 r 1/k where they are both valid. To this end, we consider the behaviours of the inner solution (A 8) for r 1 and the outer solution (A 11) for U = kr 1. Using the behaviour of the function J determined in the previous section for r 1 and ω 1, the inner solution for r 1 is

In turn, the outer solution (A 11) becomes at leading order for kr 1 and ω 1

where β 0 = 2F h /Ro. The expressions (A 30) and (A 31) match at leading order if

(A 33) The latter relation yields a dispersion relation valid for small k:

We can check a posteriori that it is legitimate to neglect the right-hand side of (A 33).

In fact, (A 34) simply imposes that the inner solution (A 30) decays with r and this condition could have been obtained from (A 30) without solving the outer problem. An alternative way of writing (A 34) is When Ro = ∞, the matched asymptotic analysis conducted in appendix A and section 4.1 for the Carton & McWilliams vortex (3.1) remains valid whatever the Froude number F h . In this case, we have therefore the opportunity to study the transition from a stratified fluid to an homogeneous fluid F h = ∞.

In the limit Ro = ∞, the dispersion relation (4.8) can be rewritten in the convenient form (see (A 35))

where

where A is defined in (4.3).

In the limit F h → ∞, Λ tends to a finite value which is real and negative: Λ = -∞ 0 ξ 3 Ω 2 dξ. Hence, it is straightforward to see that (B 1) has a solution where ω is purely real and negative (i.e. Θ = 0)

In fact, this is the only solution of (B 1). Like for vortices with a non-zero circulation [START_REF] Widnall | Theoretical and experimental study of the stability of a vortex pair[END_REF], the long-wavelength bending mode of the Carton & McWilliams vortex in homogeneous non-rotating fluid is therefore neutral with a negative frequency. This implies that there is no critical radius Ω(r c ) = ω.

Figure 23 shows the two solutions of (B 1) for several Froude numbers from F h = 3 to F h = ∞ and compare them to numerical results. We see that the branch with positive growth rate predicted by (B 1) exists only for small wavenumbers when F h = 4 (figure 23b). It disappears when the growth rate vanishes because the integration path is on the real axis. For F h = 4.5, this branch shrinks further and when F h = ∞, it has completely disappeared. In contrast, the other branch continues to exist as F h increases. Its frequency remains negative but its damping rate diminishes as F h increases. It becomes neutral in the limit F h = ∞. For all the Froude numbers investigated, (B 1) is in good agreement with the numerical results for small wavenumber.

Appendix C. Long-wavelength analysis for vortices with angular velocity decaying algebraically in stratified rotating fluid

In this appendix, a long-wavelength asymptotic analysis similar to the asymptotic analysis of appendix A is carried out for vortices with an angular velocity behaving like Ω a n /r 2n for large radius. In contrast to appendix A, it is more convenient to assume that ω is small from the outset and expanded with k

Eunok Yim and Paul Billant

This implies that the frequency will be of the same form as in (4.20) with ω (0) given by (4.21) and ω (1) by

where

dr + g 2 (k, η 0 , ω (0) , F h , Ro), (C 34) where the functions g 1 and g 2 will be different for each particular profile considered and would have to be determined by matching the inner and outer solutions. However, without determining the explicit expression of these functions, it is easy to see that they should be real if ω (0) is real. Hence, if δ i is small, we can deduce that the growth rate will be at leading order

since the residue of the integral in (C 34) at r c is a -1 = ζ 1 (r c )K 2 1 (kβ 0 r c )/Ω 0 (r c ). Equation (C 35) can be also rewritten in terms of the full profiles (C 31) 36)