## Abstract

Frontal collisions of mesoscale baroclinic dipoles are numerically investigated using a three-dimensional, Boussinesq, and *f*-plane numerical model that explicitly conserves potential vorticity on isopycnals. The initial conditions, obtained using the potential vorticity initialization approach, consist of twin baroclinic dipoles, balanced (void of waves) and static and inertially stable, moving in opposite directions. The dipoles may collide in a close-to-axial way (cyclone–anticyclone collisions) or nonaxially (cyclone–cyclone or anticyclone–anticyclone collisions). The results show that the interacting vortices may bounce back and interchange partners, may merge reaching a tripole state, or may squeeze between the outer vortices. The formation of a stable tripole from two colliding dipoles is possible but is dependent on diffusion effects. It is found that the nonaxial dipole collisions can be characterized by the interchange between the domain-averaged potential and kinetic energy. Dipole collisions in two-dimensional flow display also a variety of vortex interactions, qualitatively similar to the three-dimensional cases.

## 1. Introduction

Mesoscale vortical structures are frequent phenomena in the oceans and atmosphere, and the vortex dipole, also called vortex pair, vortex couple, double-vortex, or mushroomlike vortex, is one of the simplest among these mesoscale structures (Carton 2001). In nonrotating fluids, for example, dipoles are spontaneously generated in Von Kármán wakes in two-dimensional turbulence (Couder and Basdevant 1986). A dipolar solution in nonrotating two-dimensional fluids, different from the Lamb dipole, was recently found by Juul Rasmussen et al. (1996).

Vortex dipoles have been observed in many places in the ocean (e.g., Fedorov and Ginzburg 1986; Munk et al. 1987), including the Alaska Coastal Current (Ahlnäs et al. 1987), along the Vancouver Island coast (Ikeda et al. 1984), along the California coast (Sheres and Kenyon 1989; Simpson and Lynn 1990), south of Madagascar (de Ruijter et al. 2004), in Tartar Strait (Ginzburg and Fedorov 1984), in the Norwegian Coastal Current (Johannessen et al. 1989), and in the Oyashio Front (east of Japan; Vastano and Bernstein 1984).

Laboratory experiments in rotating tanks show that barotropic dipoles can be generated from an impulsive jet (Kloosterziel et al.1993), become important transporters of fluid (Eames and Flór 1998), and finally decay—for example, because of bottom friction effects (Sansón et al. 2001). Theoretically, Stern (1975) derived an exact dipolar solution, for nondivergent barotropic flow on the *β* plane, called the modon (Flierl et al. 1983; McWilliams 1983; Berson and Kizner 2002, and references therein). Numerically, the generation of oceanic dipoles was investigated using a two-layer, *f* plane, shallow-water model by Mied et al. (1991).

Since the dipole is a coherent vortical structure with a propagation speed, it may approach a coast or obstacle and rebound (Carnevale et al. 1997; Danaila 2004), interact with a sloping boundary (Kloosterziel et al. 1993), or interact with a jet (Vandermeirsh et al. 2002). The dipoles may also generate or interact with inertia–gravity waves (Afanasyev 2003; Godoy-Diana et al. 2006). Furthermore, interactions between dipoles, or collisions, are also possible. These interactions may occur in a number of ways, namely, frontal or head-on axial collisions, head-back (dipole overtaking or merging) axial collisions, oblique collisions, or frontal nonaxial collisions (Voropayev and Afanasyev 1994, chapter 4). Laboratory experiments, in nonrotating platforms, have reproduced axial head-on dipole collisions, and subsequent vortex partner interchange (van Heijst and Flór 1989), as well as oblique dipole collisions (Voropayev and Afanasyev 1992). Two-dimensional numerical simulations, in nonrotating fluids, of oblique dipole collisions were carried out by Couder and Basdevant (1986) in the context of two-dimensional turbulence. Head-on and overtaking dipole collisions of barotropic equivalent modons on the *β* plane were numerically investigated by McWilliams and Zabusky (1982). Laboratory experiments in rotating fluids on a *β* plane, confirmed by point-vortex numerical simulations, have shown that nonaxial dipole collisions lead to a large mass exchange between the dipoles or the ambient fluid (Velasco Fuentes and van Heijst 1995).

Here, we take a step forward in our understanding of dipole interactions by numerically investigating nonaxial frontal collisions of mesoscale baroclinic dipoles. Since the vortices are easily characterized by their conserved potential vorticity (PV), we use, as a primary model, a three-dimensional numerical algorithm that explicitly conserves PV on isopycnals (described in section 2). The initial conditions, obtained using the PV initialization approach, consist of twin baroclinic dipoles moving in opposite directions (section 3a). The translating dipoles are balanced (void of waves) and remain always static and inertially stable. We describe six different classes of possible collisions, which depend on the initial distance between the dipoles axes of motion. These cases include the axial cyclone–anticyclone collision (section 3b) and the nonaxial anticyclone–anticyclone (section 3c) and cyclone–cyclone (section 3d) collisions. We show that, besides the familiar close-to-axial collision, there is a rich variety of vortex interactions during the dipole collisions, including vortex merging and splitting, vortex bouncing, vortex squeezing, and tripole formation. These processes involve interchange between the domain averaged kinetic and potential energy. Dipole collisions in two-dimensional flow are also numerically investigated (section 3e), resulting also in a variety of vortex interactions, qualitatively similar to the three-dimensional cases.

## 2. Numerical model and parameters

The time evolution of the three-dimensional (baroclinic) dipoles is simulated using a triply periodic, volume-preserving, nonhydrostatic numerical model with the Boussinesq and *f*-plane approximations (Dritschel and Viúdez 2003) initialized using the PV initialization approach (Viúdez and Dritschel 2003). The vertical displacement 𝒟 of isopycnals with respect to the reference density configuration is 𝒟(**x**, *t*) ≡ *z* − *d*(**x**, *t*), where *d*(**x**, *t*) ≡ (*ρ*(**x**, *t*) − *ρ*_{0})/ϱ* _{z}* is the depth, or vertical location, that an isopycnal located at

**x**at time

*t*has in the reference density configuration defined by

*ρ*

_{0}+ ϱ

*, where*

_{z}z*ρ*is the density, and

*ρ*

_{0}> 0 and ϱ

*< 0 are given constants. Static instability occurs when 𝒟*

_{z}*≡ ∂𝒟/∂*

_{z}*z*> 1. The Rossby number ℛ ≡

*ζ*/

*f*, and the Froude number ℱ ≡

*ω*/𝒩, where

_{h}*ω*and

_{h}*ζ*are the horizontal and vertical components of the relative vorticity

**, respectively,**

*ω**f*is the constant Coriolis frequency, and 𝒩is the total Brunt–Väisälä frequency. We simulate here baroclinic dipoles in the regime of static (𝒟

*< 1) and inertial stability (|ℛ| < 1) so that the flow remains largely in hydrostatic balance.*

_{z}The dimensionless PV anomaly ϖ ≡ Π − 1, where Π ≡ (** ω**/

*f*+

**k**) ·

**∇**

*d*is the dimensionless total PV, is represented by contours lying on isopycnals, and PV material conservation (

*d*ϖ/

*dt*= 0) is made explicit by PV contour advection. The state variables are the components of the vector potential

**= (**

*ϕ**ϕ*,

*ψ*,

*ϕ*), which provide the velocity

**u**= (

*u*,

*υ*,

*w*) = −

*f*

**∇**×

**and the vertical displacement of isopycnals 𝒟 = −ɛ**

*ϕ*^{2}

**∇**·

**, where ɛ ≡**

*ϕ**c*

^{−1}≡

*f*/

*N*is the ratio of

*f*to the mean Brunt–Väisälä frequency

*N*. Since

**is triply periodic,**

*ϕ***u**and 𝒟 are triply periodic as well, while obviously

*d*(or

*ρ*) is not. This numerical model (referred to as the Aℬϖ model) integrates the horizontal components of the dimensionless ageostrophic vorticity

**A**

*≡ (𝒜, ℬ) ≡*

_{h}

*ω**/*

_{h}*f*−

*c*

^{2}

**∇**

*𝒟 (appendix A). The horizontal potentials*

_{h}

*ϕ**= (*

_{h}*ϕ*,

*ψ*) are obtained every time step by inverting ∇

^{2}

*ϕ**=*

_{h}**A**

*, while the vertical potential*

_{h}*ϕ*is obtained by inverting the definition of ϖ(

*ϕ*,

*ψ*,

*ϕ*).

The total energy *E _{T}* ≡

*E*+

_{K}*E*, where

_{P}*E*≡

_{K}*u*

^{2}+

*υ*

^{2}+

*w*

^{2}and

*E*≡

_{P}*N*

^{2}𝒟

^{2}are the kinetic and potential energy, respectively. The domain-averaged total energy 〈

*E*〉 is conserved (appendix B), where

_{T}is the domain average of any field ℰ, and *n _{x}*,

*n*, and

_{y}*n*are the dimensions of the numerical grid.

_{z}We use an *n _{x}* ×

*n*×

_{y}*n*= 128

_{z}^{3}grid, with

*n*= 128 isopycnal layers, in a domain of vertical extent

_{l}*L*= 2

_{z}*π*(which defines the unit of length) and horizontal extents

*L*=

_{x}*L*=

_{y}*cL*, with

_{z}*c*= 100. We take the (mean) buoyancy period (

*T*

_{bp}≡ 2

*π*/

*N*) as the unit of time by setting

*N*≡ 2

*π*. One inertial period (

*T*

_{ip}≡ 2

*π*/

*f*) equals 100

*T*

_{bp}. The time step

*δt*= 0.05, and initialization time

*t*= 5

_{I}*T*

_{ip}. The initialization time is the minimum time required for the fluid to reach its initial perturbed state with minimal generation of inertia–gravity waves.

To avoid the generation of grid-size noise, a biharmonic hyperdiffusion term − *μ*(∇^{4}_{q}𝒜, ∇^{4}_{q}ℬ), where **∇*** _{q}χ* ≡ (∂

*χ*/∂

*x*, ∂

*χ*/∂

*y*, ɛ∂

*χ*/∂

*z*) is the gradient operator in the vertically (quasigeostrophic) stretched space, is added to the equations for the rate of change of

**A**

*. The hyperviscosity coefficient*

_{h}*μ*is chosen by specifying the damping rate (

*e*-folding,

*e*) of the largest wavenumber in spectral space per inertial period, which was set constant to 10.

_{f}As a secondary numerical model we use the full pseudospectral version of the hybrid 𝒜ℬϖ algorithm (appendix A). This model, referred to as the 𝒜ℬ𝒞 model, uses the same grid-based procedures as the 𝒜ℬϖ model except that it does not use any of the procedures involving the PV contours. The prognostic variables are **A** ≡ (𝒜, ℬ, 𝒞) ≡ ** ω**/

*f*−

*c*

^{2}

**∇**𝒟. Hence, there is no inversion of PV, but only a Poisson equation for all three components of the vector potential, ∇

^{2}

**=**

*ϕ***A**, which is inverted spectrally every time step. All the parameter settings are the same except for the hyperviscous damping rate, which was significantly increased to maintain numerical stability.

## 3. Numerical results

### a. Initial dipole configuration

The initial dipole consists of two baroclinic vortices, each one defined by a three-dimensional ellipsoidal distribution of PV anomaly ϖ, constant on ellipsoidal surfaces, which varies linearly with the ellipsoidal volume, with ϖ = 0 on the outermost surface, and ϖ_{min} = −0.75 and ϖ_{max} = 0.75 at the center of the cyclone and anticyclone, respectively (Fig. 1a). The PV distribution is discretized by placing a number of PV contours within each isopycnal surface crossing through the vortex. The middle isopycnal surface (*i _{l}* = 65) has the maximum number of contours (

*n*= 10). The initial PV contours are ellipses with a ratio of major (

_{c}*a*) to minor (

_{x}*a*) semiaxes lengths

_{y}*a*/

_{x}*a*= 1.5. The largest ellipse, located on the middle isopycnal surface, has

_{y}*a*= 0.6

_{x}*c*and

*a*= 0.4

_{y}*c*. The vortex ellipsoidal shape is initially prescribed only to facilitate the transition toward an equilibrium shape, largely independent on the initial conditions, reached at the end of the initialization time (explained below).

The vertical semiaxes are *a*^{+}_{z} = 0.4 and *a*^{−}_{z} = 0.27 for the cyclone (ϖ > 0) and the anticyclone (ϖ < 0), respectively. These values were chosen so that the dipole described a straight trajectory (Figs. 1b,c). This asymmetry in the vertical extent of the vortices is due to the PV anomaly being prescribed in the reference configuration (flat isopycnals) at the beginning of the initialization time (*t*′ = 0). During the initialization period (from *t*′ = 0 to 5*T*_{ip}) the isopycnals stretch (in the anticyclone) and shrink (in the cyclone) to reach a balanced state so that the final adjusted state of the dipole is not exactly antisymmetric. For the initial conditions (*t* = 0) of the collision simulations we use the state of the dipole at *t*′ = 19*T*_{ip} (Fig. 1c, hereinafter just referred to as the initial dipole configuration). At this stage the PV vortices have long time ago deformed from their initial ellipsoidal configuration. Thus, this large period of time (19*T*_{ip}) assures that the initial PV configuration of the dipoles has been adjusted to a steady PV distribution.

In the initial dipole configuration the horizontal velocity **u*** _{h}* in the anticyclone is somewhat larger than in the cyclone (Fig. 2a), with the maximum horizontal velocity, reaching

*u*

_{max}≡ max{|

**u**

*|} = 0.77, located at the center of the dipole. The vertical velocity*

_{h}*w*(Fig. 2b) is 10

^{4}times smaller than |

**u**

*|, and presents the quadrupolar pattern typical of mesoscale quasigeostrophic balanced dipoles (Pallàs-Sanz and Viúdez 2006). The isopycnal displacement 𝒟 (𝒟 < 0 in the anticyclone; Fig. 2c) and the vertical vorticity*

_{h}*ζ*(

*ζ*< 0 in the anticyclone; Fig. 2d) have similar vertical extents in both vortices. The relative vorticity changes sign on the northern and southern sides of the dipole with both |

*ζ*| and |𝒟| in the anticyclone larger than in the cyclone. The dipole is static and inertially stable with 𝒟

_{z}_{max}≡ max{𝒟

*(*

_{z}**x**)} = 0.36, ℛ

_{min}≡ min{ℛ(

**x**)} = −0.62, ℛ

_{max}≡ max{ℛ(

**x**)} = 0.41, and ℱ

_{max}≡ max{ℱ(

**x**)} = 0.32 at

*t*= 0.

Last, the initial PV configuration of the two-dipole frontal collision simulations is obtained by adding, to the dipole described above (left dipole, *L* in Fig. 3), a *y*-inverted copy of dipole *L* (right dipole, *R* in Fig. 3) separated from *L* by Δ*X* ≃ *πc* (fixed) and Δ*Y* (variable). The initially prescribed *y*-offset Δ*Y* makes the initial location of the dipoles change along lines *L* and *R* in Fig. 3. Specification of Δ*Y* defines, therefore, the different numerical simulations described below. Since the two-dipole PV configuration is now different from the previous single dipole configuration, the numerical simulations are initialized from *t* = 0 to 5*T*_{ip}, using again the PV initialization approach. The actual *y*-offset between the dipoles, that is, at the end of the initialization period (*t* = 5*T*_{ip}), is defined as 𝒴 ≡ 𝒴* _{R}* − 𝒴

*, where the*

_{L}*y*coordinate of the

*L*and

*R*dipole are the geometric centers of the pair of vortices,

where the *y* coordinates 𝒴^{±}_{R} and 𝒴^{±}_{L} of every vortex are defined as

where the integration volumes 𝒱^{±}_{L} and 𝒱^{±}_{R} comprise the grid locations with |ϖ| > 0.1. We define the dimensionless *y*-offset *δ* as 𝒴per unit of dipole *y* extent,

After the initialization time, the dipoles explicitly conserve PV and approach each other almost steadily along the axes *y* = 𝒴* _{L}* (dipole L) and

*y*= 𝒴

*(dipole R).*

_{R}### b. P–N collision

We describe first the familiar frontal, close-to-axial dipole collision. As an example we show a case with *δ* = −0.27 (Fig. 4). The dipoles collide in such way that the interaction occurs between vortices of opposite PV sign; that is, we have *L*^{−}–*R*^{+} and *L*^{+}–*R*^{−} interactions. We refer to this case as a P–N collision (positive–negative PV collision). The dipole interaction occurs between *t* = 14*T*_{ip} and 20*T*_{ip}, when the dipoles interchange vortex partners cleanly, with little PV mixing between vortices. After collision, the new dipoles leave the impinging region and propagate with a straight trajectory, which can be approximately obtained following the relative extrema of ϖ during the numerical simulation (Fig. 5a). Owing to the negative *δ*, the new axis of motion is however not perpendicular to the original one (along constant *y*). In this case, the angle between old and new dipoles trajectories (*α*, defined from the trajectory of vortex *R*^{−}) is *α* < *π*/2. An initial *δ* = 0 (perfectly axial collision) results in new dipole trajectories perpendicular to the original ones (*α* = *π*/2), while *δ* > 0 results in new dipole trajectories with *α* > *π*/2 (not shown). We show later that the P–N collision is typical for a range of values *δ*_{min} < *δ* < *δ*_{max}, where −0.52 < *δ*_{min} and *δ*_{max} < 0.20.

In this P–N case, the flow before, during, and after the collision remains always static and inertially stable with |𝒟* _{z}*| < 1 and |ℛ| < 1. Extreme Rossby and Froude numbers remain quite constant during the dipoles interaction (Table 1). Both

*u*

_{max}and maximum vertical velocity

*w*

_{max}≡ max{

*w*} decrease and 𝒟

_{max}increases slightly before the dipoles collision, and the opposite changes happen after collision (Figs. 6 and 7a). Thus, some deceleration (acceleration) of the flow, at least in the extreme values, occurs before (after) collision. Numerically, 〈

*E*〉 is well conserved (Fig. 8a), and there is little interchange between 〈

_{T}*E*〉 and 〈

_{K}*E*〉 in the flow (Fig. 8b).

_{P}### c. N–N collision

We describe next two different classes of dipoles collisions in which the dipole interaction occurs mainly between the anticyclonic vortices (N–N collisions).

#### 1) Case N–N(1)

In this case *δ* = −0.52 is smaller than in the previous case so that only the anticyclones collide (collision *L*^{−}–*R*^{−}; Fig. 9). The interaction between the anticyclones resembles, in some way, an elastic collision between solid bodies. There is some mass transfer between the interacting anticyclones (between *t* ≃ 16–24*T*_{ip}), although their core remain largely isolated.

The anticyclones do not merge due to the influence of the cyclones. The gradient pressure force exerted by the outer cyclones induces a linear momentum, in the new anticyclone partners, opposite to the initial dipole momentum, preventing the merging of the two anticyclones and causing the change in direction of the new dipoles. After collision, the anticyclones bounce back and join the cyclone of the companion dipole, which experiences no trajectory curvature change (Fig. 5b). Thus, this case may be considered as a limit case of the previously described class P–N since here, qualitatively speaking, the percentage of positive–negative PV collision may be taken as zero, the percentage of negative–negative PV collision as 100%, and *α* = 0.

As in case P–N, *u*_{max} reaches a minimum (Fig. 6a), and 𝒟 a maximum (Fig. 7a), during the interaction period, with 𝒟_{max} larger than in case P–N. Here, however, *w* reaches a maximum during this period (Fig. 6b); 𝒟_{z}_{max} reaches a significant maximum when the cyclones collide (Fig. 7b). This maximum is due to the fact that the colliding anticyclones (which partially merge) have larger 𝒟* _{z}* (Fig. 2c) than the cyclones (which do not merge). In this case ℛ

_{max}and ℱ

_{max}are very similar to those in the case P–N (Table 1). The ratio 〈

*E*〉/〈

_{P}*E*〉 reaches a maximum during the collision, and therefore 〈

_{T}*E*〉/〈

_{K}*E*〉 reaches a minimum (Fig. 8b). Thus, there is an interchange between 〈

_{T}*E*〉 and 〈

_{K}*E*〉. This interchange seems to be due to the deceleration of the flow and deepening of isopycnals during the collision.

_{P}The 〈*E _{T}*〉 experiences, however, an increment Δ〈

*E*〉 ≃ 0.02/2.07 ≃ 0.01 = 1%, relative to its initial value (Fig. 8a). This nonconservative change is related to the intrinsic numerical diffusion associated to the discretization of the PV field (which implies contour merging during the anticyclones partial merging). We address this numerical process in the next case where the vortex merging is larger and, therefore, its effects on 〈

_{T}*E*〉 are more important.

_{T}#### 2) Case N–N(2)

If *δ* is decreased by a small amount relative to the previous case N–N(1), setting *δ* = −0.56 (i.e., a 7% decrease), the evolution of the dipoles collision is very different (Fig. 10). Initially (*t* ≃ 16–20*T*_{ip}) the anticyclones collide in a way similar to the case N–N(1) (Fig. 9). However, here the anticyclones almost fully merge so that the vortical structure transforms into a tripole. The tripole rotates (*t* = 22–34*T*_{ip}) with a negative (anticyclonic) phase speed, completing at least a rotation of 90° during about 12 inertial periods (Fig. 5c). The tripole is unstable and eventually the anticyclone splits (*t* = 34–36*T*_{ip}), and the two-dipole system is recovered with an axis of motion rotated about 90° relative to the initial axis. The rotation angle of the axis of motion depends on *δ*, first increasing with decreasing *δ*, reaching a maximum value of about 150°, and decreasing afterwards. A stable tripole was not found.

In this case, *u*_{max} decreases when the dipoles approach and collide [as happens in the two previous cases, P–N and N–N(1)], then increases (reaching a maximum at *t* = 24*T*_{ip}) and decreases during the tripole rotation (Fig. 6a). Then *u*_{max} increases again after the tripole splits. The time evolution of *w*_{max} is similar to the previous case N–N(1) except after the vortex merging (*t* > 22*T*_{ip}) since a new relative maximum (*t* = 35*T*_{ip}) is reached when the tripole splits (Fig. 6b). Here 𝒟_{max} experiences the largest maximum, compared to all cases, increasing during the vortex merging and generally decreasing during the tripole and vortex splitting episodes (Fig. 7a). Furthermore, 𝒟_{z}_{max} reaches an overall maximum at the end of the collision period, oscillates during the tripole episode, and finally decreases when the vortex splits (Fig. 7b). Though still relatively small, this case experienced the largest change in 𝒟_{z}_{max} (Table 1). These changes in the extremum values are quite consistent with the changes in 〈*E _{K}*〉 and 〈

*E*〉. The ratio 〈

_{P}*E*〉/〈

_{P}*E*〉 increases during the vortex merging, oscillates during the tripole rotation, and decreases after the splitting (Fig. 8b). Thus, there is an interchange between 〈

_{T}*E*〉 and 〈

_{K}*E*〉, with kinetic energy being converted into potential energy during the vortex merging and converted back into kinetic energy during the vortex splitting.

_{P}As in the previous and following cases, 〈*E _{T}*〉 increases with time during the vortex merging (Fig. 8a). In this respect, the behavior of the discrete numerical model differs from the continuous theoretical model (which conserves 〈

*E*〉). The change in 〈

_{T}*E*〉 is related to the discretization of the PV field, which necessarily implies PV contour merging and splitting during strong vortex interactions. In the case P–N, described above, the vortices do not merge and 〈

_{T}*E*〉 is very well conserved (Fig. 8a).

_{T}To assess the importance of the above nonconservative effects we carried out a series of numerical simulations using the 𝒜ℬ𝒞model, which does not make any use of PV contours nor PV advection and inversion. The numerical diffusion, acting now on the three-dimensional biharmonic operator ∇^{4}_{q} (in the quasigeostrophic stretched or computational domain) of **A** = (𝒜, ℬ, 𝒞), had to be significantly increased, to maintain numerical stability, to *e _{f}* = 200 (

*e*= 180 developed grid-size noise). The 𝒜ℬ𝒞 simulations were initialized with the potential

_{f}**provided by the 𝒜Bϖ simulation [case N–N(2)] at the end of the initialization time (**

*ϕ**t*= 5

*T*

_{ip}). The other numerical parameters remained unchanged.

The results are similar to the 𝒜ℬϖ simulation during the vortex merging except that, due to the larger diffusion, PV is not conserved and the dipoles displacement speed is smaller (Fig. 11). The anticyclones merge and a tripole is formed. In this case the tripole remains stable and only decays by diffusion effects. Consistent with the large diffusion dumping rate, 〈*E _{T}*〉 decreases with time (Fig. 12). However, besides this diffusion effect, the 〈

*E*〉 and 〈

_{K}*E*〉 behave in a way similar to the 𝒜ℬϖ simulation during the vortex merging, where 〈

_{P}*E*〉/〈

_{P}*E*〉 increases. Further 〈

_{T}*E*〉–〈

_{K}*E*〉 interchange is related to oscillations in the cyclones trajectories around the anticyclone. Note that the change of 〈

_{P}*E*〉 in the 𝒜ℬϖ simulation amounts for a change (per

_{T}*T*

_{ip}) of

*χ*≃ 0.1% of its initial value, while

*χ*≃ −0.3% in the 𝒜ℬ𝒞 simulation. Thus, the process of reaching a stable tripole from a two-dipole collision is related to diffusion effects typical of an irreversible dynamics since the reverse process (namely a stable tripole splitting into two dipoles) is a contradiction. The formation of an unstable tripole, followed by tripole splitting, and recovery of the two-dipole system is more consistent with a reversible dynamics.

### d. P–P collision

We describe next three classes of frontal dipoles collisions in which *δ* > 0 so that the dipoles interact through the positive PV vortices.

#### 1) Case P–P(1)

In this case *δ* = 0.20, the cyclones collide, partially merge, but finally split and, in a way similar to case N–N(1), bounce back interchanging the anticyclonic partner (Fig. 13). Though there is considerable vortex merging, the PV core of the cyclones remains isolated (Fig. 14a). The PV filamentation is, however, larger than in the case N–N(1), and the cyclones after the interaction become weaker than the anticyclones, so that the resulting dipoles have negative trajectory curvature.

Similar to the previous cases, *u*_{max} decreases during the vortex merging (*t* = 15–19*T*_{ip}) and increases during the vortex splitting (*t* = 23–31*T*_{ip}) (Fig. 6a). Contrary to the N–N cases, 𝒟_{max} decreases during the vortex merging and increases during the vortex splitting (Fig. 7a). These changes of 𝒟 are probably related to the large deformation of the cyclones, in comparison to the deformation of the anticyclones in the cases N–N, during the vortex interaction (cf. Fig. 13 with Figs. 9 and 10). Consequently, 〈*E _{P}*〉/〈

*E*〉 decreases during the cyclone merging and increases during the cyclone splitting (Fig. 8b).

_{T}#### 2) Case P–P(2)

In this and the following case, the interchange of the vortex partners is not so clear as in the previous case. Here *δ* = 0.24; that is, 8% larger than in P–P(1), a percentage similar to that between cases N–N(1) and N–N(2). The cyclones collide and fully merge forming a transitory tripole (Fig. 15). This process is, thus, similar to N–N(2) except that here the vortex merging occurs between the cyclones and the tripole splits faster into two dipoles with little change in the anticyclones trajectories (Fig. 14b). These differences are probably due to the amount of PV (i.e., the volume integration of ϖ) in the cyclone being smaller than in the anticyclones (Figs. 2c,d).

As is characteristic of the P–P collision, the deformation of the cyclone is larger than the deformation of the anticyclones of the N–N collisions. This implies that 𝒟_{max} and 〈*E _{P}*〉/〈

*E*〉 decrease during the cyclone merging (

_{T}*t*= 15–21

*T*

_{ip}) and increase during the cyclone splitting (Figs. 7a and 8b).

#### 3) Case P–P(3)

For the last case, we slightly increase *δ* = 0.36, that is, a 26% increase relative to P–P(1). In this case the cyclones collide but do not fully merge (Fig. 16). Thus, there is no interchange of vortex partners as happens in case P–P(1) (Fig. 13). This case is very similar to the previous P–P(2) except that here the cyclonic cores remain mostly isolated, so a tripole episode is not completely reached (cf. Figs. 14b,c). This case is the only one in which *u*_{max} increases during the cyclones impingement (Fig. 6a). This increase of *u*_{max} is probably related to the fact that the interacting cyclones have to squeeze between the anticyclones, which, having larger amount of PV anomaly, remain more stable and static.

### e. Two-dimensional dipole collisions

Numerical simulations with the ϖ2*D* model, the two-dimensional (2D) version of the 𝒜ℬϖ model (appendix A), have been also carried out for comparison with the 3D results. The process to obtain the initial conditions for the 2D dipole was similar to the one in the 3D simulations, using a 2D dipole with the same initial PV value (ϖ_{max} = −ϖ_{min} = 0.75) and geometry, and in the 2D case selecting the adjusted state at *t*′ = 15*T*_{ip}. Though a case-to-case comparison is not strictly possible since a 2D vortex is an entity different from a baroclinic 3D vortex, and the PV criterion is not the unique one to relate both cases, we found that, depending on the *y* offset, the dipoles may interchange partners or squeeze between the outer vortices as in the 3D case. A transient 2D tripole state, as long-lasting as in the 3D case, was not found (with ϖ_{max} = 0.75). This is probably due to the larger propagation speed of the dipoles relative to the 3D ones with the same ϖ. However, it was observed as a new temporary state, a five-pole structure (Fig. 17). During the dipole collision the inner vortices (the cyclones in this case) merge and a tripole is temporarily formed. The anticyclones are however strong enough to tear apart the cyclone, which splits in three small cyclones, forming two dipoles with a central cyclone. The dipoles are asymmetric and have a large trajectory curvature.

Longer-lasting temporary tripole states are also possible with 2D dipoles of vortices with smaller ϖ. The trajectories of the outer vortices, which are easier to follow since they do not merge, are displayed in Fig. 18 for dipoles with ϖ_{max} = 0.25, 0.50, and 0.75, respectively. In the three cases the dipoles behave qualitatively in a similar way. For small *δ* the vortices exchange partners and change 90° from their trajectory. As *δ* increases, so does the scattering angle *α*, up to a maximum value *α*_{max} after which it diminishes since the dipoles no longer interact. Note that *α* increases faster with *δ* for small ϖ_{max} (there is a smaller number of trajectories pointing northward and southward in Fig. 18a than in Fig. 18c), which is due to the fact that the dipoles with small ϖ_{max} have slower speed and interact, therefore, during longer times than the dipoles with large ϖ_{max}.

The scattering angle *α* as a function of *δ* is shown in Fig. 19 for the three PV cases. The scattering angle is computed by fitting the first and last five points of every trajectory in Fig. 18 to a straight segment. The results confirm that *α* ≃ 90° for *δ* = 0, that *α* increases with *δ* faster for smaller ϖ_{max}, and that a maximum *α*_{max} > 180° is reached after which *α* decreases with increasing *δ* toward *α* → 180°. The *α*_{max} is larger for small ϖ_{max}. The largest vortex interaction and merging occur for *α* ∈ [180°, *α*_{max}] (tripole and five-pole transient states are indicated in Fig. 19). Note that in these extreme cases the outcome vortices may have large trajectory curvatures and therefore the scattering angle *α* is time dependent.

## 4. Concluding remarks

The frontal collision of baroclinic dipoles is a complicated process that includes different vortex interactions, depending on the horizontal dimensionless *y*-offset *δ* between the colliding dipoles. The cyclone–anticyclone collisions (as in case P–N here) occur for a range of values of *δ* close to zero. The axial (*δ* = 0) cyclone–anticyclone causes a new axis of motion perpendicular to the original one. For larger *δ*, the dipoles may interact in cyclone–cyclone or anticyclone–anticyclone nonaxial collisions, and partially or fully merge. In these cases the interacting vortices may bounce back and interchange partners [cases N–N(1) and P–P(1)], may merge reaching a tripole state [cases N–N(2) and P–P(2)], or just squeeze between the outer vortices and continue without interchanging partners [case P–P(3)]. The nonaxial dipole collisions may be characterized by the interchange between the domain averaged potential and kinetic energy. No significant energy interchange occurs, at least within the dipoles parameters used here, in the close-to-axial type of dipole collision.

The formation of a tripole from two colliding dipoles is possible but is dependent on diffusion effects. For very small diffusivity the tripole state is only transient [cases N–N(2) and P–P(2) using the 𝒜ℬϖ model], consistent with reversible dynamics. For larger diffusivity the impinging dipoles may form a stable, though decaying, tripole, which is an irreversible process [case N–N(2) using the 𝒜ℬ𝒞 model].

This study has shown that new interesting phenomena on dipoles collision are possible, though it has obviously not exhausted the complete parameter space, which is very large. We have only explored the effect of changing the initial horizontal offset between the colliding dipoles. Other possible variables are the PV anomaly, the size of the vortices, the dipoles trajectory curvature, and the angle of impingement (oblique collisions). Also, it is possible that other phenomena, like the spontaneous generation of small-scale inertia–gravity waves, may occur during the collision of balanced (void of waves) mesoscale dipoles with larger PV anomalies. This topic is left for future research.

## Acknowledgments

We are very thankful to two anonymous reviewers for their valuable comments. We acknowledge partial support from the Spanish Ministerio de Educación y Ciencia (Grant CGL2005-01450/CLI).

## REFERENCES

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

### APPENDIX A

#### The Prognostic Equations

Let *χ̃* ≡ *χ*/*f* for any quantity *χ*. From the vorticity equation

and mass conservation equation

we obtain the equation for **A** ≡ ** ω̃** −

*c*

^{2}

**∇**𝒟, that is, the prognostic equation of the 𝒜ℬ𝒞 model,

The 𝒜ℬϖ model integrates the horizontal components of (A3), that is, the equation for the rate of change of the dimensionless horizontal ageostrophic vorticity **A**_{h} ≡ *ω̃*_{h} − *c*^{2}**∇**_{h}𝒟,

The third prognostic equation is the explicit conservation of potential vorticity *d*ϖ/*dt* = 0.

Last, the ϖ2D model is the two-dimensional (horizontal) version of the 𝒜ℬϖ model. In this case *ω _{h}* = 𝒟 = 𝒜 = ℬ =

*ϕ*=

*ψ*= 0, and the 𝒜ℬϖ model degenerates in the PV conservation equation for the PV anomaly ϖ, which becomes identical to the conservation of the dimensionless vertical vorticity ζ̃ ≡

*ζ*/

*f*,

### APPENDIX B

#### Derivation of ∂〈ET〉/∂t = 0

Multiplying by **u** the three-dimensional Boussinesq momentum equation

where

and the density anomaly *ρ*′(**x**, *t*) ≡ *ρ*(**x**, *t*) − ϱ* _{z}z* −

*ρ*

_{0}, we obtain

where we have used the mass conservation equation

Thus, in a triply periodic domain

## Footnotes

* Additional affiliation: École Nationale Supérieure de Techniques Avancées, Paris, France

*Corresponding author address:* Sara Dubosq, Clayrac, 81630 Salvagnac, France. Email: sara.dubosq@free.fr