# Local Mixing Events in the Upper Troposphere and Lower Stratosphere. Part I: Detection with the Lyapunov Diffusivity

Francesco d'Ovidio, Emily Shuckburgh, Bernard Legras

► To cite this version:

Francesco d'Ovidio, Emily Shuckburgh, Bernard Legras. Local Mixing Events in the Upper Troposphere and Lower Stratosphere. Part I: Detection with the Lyapunov Diffusivity. Journal of the Atmospheric Sciences, 2009, 66, pp.3678. 10.1175/2009JAS2982.1. hal-00759990

HAL Id: hal-00759990

https://hal.science/hal-00759990

Submitted on 22 Oct 2021

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# Local Mixing Events in the Upper Troposphere and Lower Stratosphere.
## Part I: Detection with the Lyapunov Diffusivity

FRANCESCO D'OVIDIO

Laboratoire de Météorologie Dynamique, Ecole Normale Supérieure, and ISC-PIF, Paris, France

EMILY SHUCKBURGH

British Antarctic Survey, Cambridge, United Kingdom

BERNARD LEGRAS

Laboratoire de Météorologie Dynamique, Ecole Normale Supérieure, and CNRS, Paris, France

(Manuscript received 16 October 2008, in final form 5 June 2009)

ABSTRACT

A new diagnostic (the “Lyapunov diffusivity”) is presented that has the ability to quantify isentropic mixing in diffusion units and detects local mixing events by describing latitude–longitude variability. It is a hybrid diagnostic, combining the tracer-based effective diffusivity with the particle-based Lyapunov exponent calculation. Isentropic mixing on the 350-K surface shows that there is significant longitudinal variation to the strength of mixing at the northern subtropical jet, with a strong mixing barrier over Asia and the western Pacific, a weaker mixing barrier over the western Atlantic, and active mixing regions at the jet exits over the eastern Pacific and Atlantic.

# 1. Introduction

It is now well established that the distribution of tracers in the upper troposphere and the lower stratosphere (UTLS) strongly depends on the transport and mixing properties of the flow. It is also well established that the dominant isentropic motion induces a chaotic type of tracer advection, giving rise to strongly inhomogeneous stirring (and thus, in the presence of diffusion, inhomogeneous mixing).<fnref n="1" /> This segregates tracers into distinct well-mixed reservoirs separated by regions of relatively weaker isentropic mixing. For simplicity, in the following we drop the term “isentropic,” and we will use the term “mixing barrier” to denote a coherent region

across which the isentropic mixing strength presents a local minimum. It should be noted, however, that some irreversible mixing does occur across such regions, which is evident in the mixing lines seen in observed tracer–tracer correlations (e.g., Hoor et al. 2002).

Motivated by the observation that the segregation into distinct reservoirs is mainly zonal, several diagnostics based on zonal averages have been developed and applied to large wind datasets provided by analysis and reanalysis from operational weather centers. Calculation of the stretching rates of potential vorticity (PV) contours have been used to show that the edge of the stratospheric polar vortex corresponds to a minimum of this quantity (e.g., Pierce and Fairlie 1993), and this technique has been applied to the UTLS (Bithell and Gray 1997; Scott et al. 2003), providing some evidence of a minimum at the subtropical jet (at least in winter). Effective diffusivity, a diagnostic introduced by Nakamura (1996) and applied to the UTLS by Haynes and Shuckburgh (2000b) and Scott et al. (2003), goes one step further by considering the evolution of a passive tracer in the wind field. This diagnostic reduces the advection and diffusion of a tracer on an isentropic surface to a pure (enhanced) diffusion in

<footnote n="1">¹ To avoid ambiguity, we refer to the reversible process of tracer filament formation by advection as “stirring,” and we call “mixing” the enhanced effect of background diffusion by stirring.</footnote>

<footnote>*Corresponding author address:* Francesco d'Ovidio, Institut des Systèmes Complexes Paris Ile-de-France, 57-59 rue Lhomond F-75005, Paris, France.
E-mail: francesco.dovidio@locean-ipsl.upmc.fr</footnote>
a new coordinate system that follows the contours of the tracer. A key advantage of this latter approach is that it provides a parameterization of isentropic mixing in terms of a contour-averaged diffusion coefficient that measures quantitatively the enhancement of small-scale diffusion by the stirring process.

For the large-scale atmospheric circulation, which is characterized by a dominant zonal component, the contours of advected tracers are close to latitude circles. Therefore, diagnostics based on either PV or tracer contour averages only vary with latitude and do not provide any information about the longitudinal variability of mixing. However, in the UTLS, mixing does not occur uniformly over latitudinal circles but can be induced by localized events like synoptic-scale systems (Polvani and Esler 2007). As noted by Bithell and Gray (1997), in this case there are frequently strong PV gradients near the edges of the synoptic systems, but there are not generally strong gradients close to a single PV contour value at all longitudes (unlike in the winter stratosphere). The investigation of such longitudinal variability is clearly not amenable to methods where a single value represents all mixing events along an entire PV or tracer contour. This problem is especially relevant in the interpretation of mixing barriers for which the contour-based diagnostics are not able to distinguish whether the barrier is continuous over the latitude circle; is modulated in longitude; or even breaks at some locations, giving way to latitudinal exchanges.

An appropriate way to represent stretching in fluid mechanics is the local Lyapunov exponent (Pierre-humbert 1991), which measures the dispersion of fluid particles along a Lagrangian trajectory. It pinpoints transport properties accurately in space and time and is not limited to contour averages. Therefore, unlike the contour-based diagnostics already mentioned, it is able to detect mixing events in both latitude and longitude. Extrema of Lyapunov exponents have been used to map coherent structures (Haller and Yuan 2000; Haller 2002; Lapeyre 2002) or identify mixing regions (d'Ovidio et al. 2004), and it can be shown that lines of extrema often act as kinematic transport barriers (Shadden et al. 2005) over a bounded range of times. However, the use of the local Lyapunov exponent alone as a quantitative measure of mixing presents a problem. As we shall see, its latitudinal variations are significantly different from—and in same cases may be even anticorrelated with—the Nakamura effective diffusivity. Understanding the discrepancies between the variability of effective diffusivity and Lyapunov exponent and the relation between stretching and mixing will be the starting point of our work and will lead us to the definition of a hybrid diagnostic that is able to quantify mixing in diffusion units and to localize mixing

events geographically. This Lyapunov-based diffusivity diagnostic combines information concerning the stretching intensity (from the Lyapunov exponents) and geometry (from the Lyapunov vectors) with the Nakamura effective diffusivity, and it can be used to map the local properties of atmospheric mixing.

We focus here on a climatological study of two-dimensional isentropic mixing at 350 K. This provides, for the first time, a quantification of the variability of mixing in latitude, longitude, and time. A more comprehensive analysis, including the relation to the major climatic modes—El Niño–Southern Oscillation (ENSO) and the North Atlantic Oscillation (NAO)—will be presented in Shuckburgh et al. (2009a, hereafter Part II).

This paper is organized in the following manner: in section 2, we describe the data [a subset of the 40-yr European Centre for Medium-Range Weather Forecasts (ECMWF) Re-Analysis (ERA-40) winds] and recall the definitions of effective diffusivity and finite-size Lyapunov exponents (FSLEs) with some details. Readers familiar with these concepts can skip section 2 and jump directly to section 3, where we compare Lyapunov calculations and effective diffusivity by looking at their correlation in space and time and analyzing their advantages and limitations as diagnostics of transport and mixing. Theoretical arguments guide the construction of a hybrid diagnostic, the Lyapunov diffusivity. This novel diagnostic is shown to be able to pinpoint and quantify mixing events in both latitude and longitude. It is applied to the UTLS region, and a climatology of 10 years is presented for the 350-K isentropic level. Finally, section 4 briefly discusses possible theoretical developments and further use of the Lyapunov diffusivity to more general geophysical cases. Technical details of the theoretical analysis are left for the appendix.

## 2. Data and method

The calculations of the effective diffusivity and of the FSLEs are performed on the 350-K isentropic surface using the ERA-40 winds provided by the European Centre for Medium-Range Weather Forecasts (Uppala et al. 2005) over the period 1981–2001. The choice of this period is dictated by the relative homogeneity of the assimilated data, notably the satellite data provided by Television and Infrared Observation Satellite (TIROS) Operational Vertical Sounder (TOVS) instruments, cloud-wind products, and surface data from buoys. The two diagnostics are computed every two days for the entire period and aggregated in monthly means.

### a. Effective diffusivity

Effective diffusivity has been introduced by Nakamura (1996) and proceeds from a transformation of the
advection-diffusion equation for a tracer with concentration $c$,

$$ \dot{c} = \frac{\partial c}{\partial t} + \mathbf{u} \cdot \nabla c = \nabla \cdot (\kappa \nabla c), \quad (1) $$

to tracer-based coordinates where the contour $\Gamma(C)$ with $c(x, y) = C$ is used as a coordinate line. Embedded contours are labeled by their area $A_C$. When the tracer has a mean pole to equator gradient, the area on the sphere is related to the equivalent latitude $\phi_e$ by $A_C = 2\pi R^2(1 - \sin\phi_e)$, where $R$ is the radius of the earth. Hence, $C$ can be considered as a function of $\phi_e$ and $t$. The minimum length $L_e(\phi_e)$ of a tracer contour with area $A_C(\phi_e)$ is $L_e = 2\pi R \cos\phi_e$, corresponding to the case in which the contour is aligned with a latitude circle.

It has been shown by Nakamura (1996) and Haynes and Shuckburgh (2000b) that, when the velocity $\mathbf{u}$ is nondivergent, (1) can be transformed to

$$ \frac{\partial C(\phi_e, t)}{\partial t} = \frac{1}{R^2 \cos\phi_e} \frac{\partial}{\partial \phi_e} \left[ \kappa \frac{L_{\text{eff}}^2(\phi_e, t)}{L_e^2(\phi_e)} \cos\phi_e \frac{\partial C(\phi_e, t)}{\partial \phi_e} \right], \quad (2) $$

where

$$ L_{\text{eff}}^2 = \oint_{\Gamma} |\nabla c| \, dl \oint_{\Gamma} \frac{dl}{|\nabla c|}. \quad (3) $$

It is assumed that $\nabla c$ does not vanish over the contour (this is generic if there are no saddle points). The value of $L_{\text{eff}}$ is larger than the length $L$ of the contour, but in practice the two quantities are fairly close (Shuckburgh and Haynes 2003; Scott et al. 2003). Both are generally much larger than $L_e$, because the contour is folded many times by chaotic advection, generating tracer filaments. The actual length results from a balance between continuous generation of small-scale filamentary structures and smoothing by dissipation. Here, dissipation $\kappa$ is a simple uniform diffusion that is meant to represent subgrid-scale three-dimensional turbulence (Legras et al. 2005). The effective diffusivity takes into account the enhancement above the minimum length of the tracer contour length by the effects of stirring and mixing, and it is given from Eq. (2) by

$$ \kappa_{\text{eff}} = \kappa \frac{L_{\text{eff}}^2}{L_e^2}. \quad (4) $$

The effective diffusivity is obtained daily by analyzing the contours of a chosen “probe” tracer evolving according to an advection-diffusion equation where advection is by the ERA-40 winds (see Haynes and

Shuckburgh 2000a,b). The tracer evolution calculation is performed by a pseudospectral method using a T159 spherical harmonics representation and associated collocation grid. The tracer is initialized with a concentration proportional to the sine of the latitude. The parameter $\kappa$ is a constant numerical diffusivity (taken as $1.5 \times 10^5 \text{ m}^2 \text{ s}^{-1}$) used to parameterize the irreversible effect of small-scale turbulence. It has been shown that the effective diffusivity is largely independent$^2$ of the value of $\kappa$. A repeat calculation with a numerical diffusivity of $5 \times 10^5 \text{ m}^2 \text{ s}^{-1}$ indicated little difference in the results. The equivalent contour length $L_{\text{eff}}$ is calculated from the tracer each day, after a spin-up time of one month,$^3$ as a function of the equivalent latitude $\phi_e$ with spacing of approximately $1^\circ$. From this, $\kappa_{\text{eff}}$ is then obtained.

### b. Finite-size Lyapunov exponents

In the absence of diffusion and sources, the evolution of an infinitesimal line element $\delta\mathbf{x}$ and of a passive tracer gradient $\nabla c$ follow very similar equations, because the variation $\delta c = \delta\mathbf{x} \cdot \nabla c$ is preserved by the motion. We have

$$ \frac{D}{Dt} \delta \mathbf{x}_i = \frac{\partial u_i}{\partial x_j} \delta \mathbf{x}_j \quad \text{and} \quad \frac{D}{Dt} \frac{\partial c}{\partial x_i} = - \frac{\partial u_j}{\partial x_i} \frac{\partial c}{\partial x_j}, \quad (5) $$

where $D/Dt$ means the time derivation along a given trajectory $\mathbf{x}(t)$, along which the motion is linearized. Over a time interval $[t_0, t_0 + \tau]$, a linear application $\mathbf{M}(\mathbf{x}_0|t_0, t_0 + \tau)$ is defined that maps a line element or a tracer from its initial value at $\mathbf{x}_0 = \mathbf{x}(t_0)$ to its final value at $\mathbf{x}(t_0 + \tau)$:

$$ \delta \mathbf{x}(t_0 + \tau) = \mathbf{M}(\mathbf{x}_0|t_0, t_0 + \tau) \delta \mathbf{x}(t_0) \quad \text{and} \quad (6) $$

$$ \nabla c(t_0 + \tau) = - \mathbf{M}^T(\mathbf{x}_0|t_0, t_0 + \tau) \nabla c(t_0). \quad (7) $$

When computing $\delta\mathbf{x}(t_0 + \tau)$ for finite time, one has the choice of prescribing the integration time $\tau$ or the ratio $\delta\mathbf{x}(t_0 + \tau)/\delta\mathbf{x}(t_0)$ between the initial and final separation. The first choice leads to the finite-time Lyapunov exponent calculation, whereas the second choice defines the finite-size calculation.

More precisely, FSLE is defined as (Ott 1993; Aurell et al. 1997)

<footnote>$^2$See Shuckburgh and Haynes (2003) and Marshall et al. (2006) for a sensitivity study of the effect of $\kappa$ on the effective diffusivity.</footnote>

<footnote>$^3$Haynes and Shuckburgh (2000b) showed that this spin-up time gives results independent of the initial conditions.</footnote>
$$ 
\begin{aligned}
\lambda[\mathbf{x}(t_0), \delta(t_0), \delta] &= \max \frac{1}{\tau} \ln \frac{|\mathbf{M}\delta\mathbf{x}|}{|\delta\mathbf{x}|} \\
&= \max \frac{1}{\tau} \ln \frac{|\mathbf{M}^T\mathbf{V}_c|}{|\mathbf{V}_c|} = \ln \frac{\delta}{\delta_0} \max \frac{1}{\tau}, \quad (8)
\end{aligned}
$$

where $\tau$ is the time at which the initial perturbation of size $\delta_0$ has grown to $\delta$, and the maximum is taken over the initial orientation of the line element or the tracer gradient. For small $\delta_0$ and $\delta$ but large $\delta/\delta_0$ (or large $\tau$, for the case of finite-time calculations) and if the flow is ergodic, $\lambda[\mathbf{x}(t_0), \delta(t_0), \delta]$ tends to a unique $\bar{\lambda}$ as a result of the Oseledec theorem (Oseledec 1968). In the other cases, $\lambda$ exhibits large spatial and temporal variations, which represent the variations of stretching properties of the flow. In particular, the length of a contour $\Gamma$ should grow exponentially with time, on average, for a length of time $t$ (i.e., $L = e^{\lambda t}$, where $\bar{\lambda}$ is the most probable value of $\lambda$ in the domain spanned by the contour).

The calculation of Lyapunov exponents is well-defined forward and backward in time. As we shall see, both calculations contain valuable and complementary information on mixing.

The finite-size Lyapunov exponents are computed on a quasi-regular latitude-longitude grid, with a 0.5° spacing in latitude. Winds from the ERA-40 spectral T159 representation (as used for the effective diffusivity calculation) are interpolated in time using cubic splines and bilinearly in space from an intermediate projection to a 480 × 242 Gaussian grid. Trajectories of Lagrangian parcels are constructed both forward and backward in time. At each grid point, the trajectories starting at the four cardinal corners of a small square with a diagonal $\delta_0 = 0.1°$ are followed until a prescribed relative separation $\delta = 20\delta_0$ is reached for one of the diagonals. Other ratios of $\delta/\delta_0$ in the same order of magnitude (5, 10, and 50) have been also tested and provide very similar results. The finite-size Lyapunov exponents and vectors (which, in the following, will be referred to as Lyapunov exponents and Lyapunov vectors for simplicity; see appendix A for more details) are obtained from the singular values and directions of the matrix that transforms the initial square into a parallelepiped (Ott 1993; Koh and Legras 2002). Other details are in Mariotti et al. (1997) and Joseph and Legras (2002).

### 3. Lyapunov diffusivity

Our aim is to develop a new diagnostic, based on Lyapunov calculations and able to parameterize the effect of mixing as a diffusion process, to provide a quantitative map of mixing events in both latitude and longitude. Here, we will attempt to deduce such a diagnostic by dynamical considerations, some analytical relations obtained for the case of a shear-dominated flow, and the direct comparison of independent calculations of effective diffusivity and Lyapunov-derived quantities.

The construction of a Lyapunov diffusivity will proceed in three steps. Key to the success of the effective diffusivity diagnostic as a parameterization of mixing is its ability to characterize the complexity of tracer filaments. Thus, in the first step, we incorporate information concerning the structure of tracer filaments from the Lyapunov vectors to form a modified Lyapunov exponent. Then, in the second step, we compare the modified Lyapunov exponents and the effective diffusivity, with the aim of generating a mixing diagnostic in units of diffusion by fitting free parameters of a functional form. Because the effective diffusivity only depends on equivalent latitude, this comparison is done with modified Lyapunov exponents averaged along contours of equivalent latitude. In the third step, the fitted parameters are used with the functional form to provide a local estimation of mixing that represents the variability in both latitude and longitude.

#### a. Relation between effective diffusivity and Lyapunov-derived quantities

In a turbulent system, a passive tracer is stirred in filamentary patches. In the early stage of their life cycle, filaments intensify local gradients by creating fingered intrusions between regions of different tracer concentrations. Such a gradient intensification mechanism is a precondition to mixing, because the redistribution of tracers in progressively thinner filaments eventually enhances the irreversible dispersion effect of small-scale turbulence.

The effective diffusivity and the Lyapunov exponent calculation provide two different approaches for detecting filamentation events. The effective diffusivity probes the velocity field with a passive tracer and a background diffusion. The enhanced dispersion effect resulting from the stirring is measured indirectly by averaging the filamentation events occurring along tracer contours. A contour grows during the initial stage of filamentation and retracts when its (thin) content is dispersed by the background diffusion. If the turbulent field has reached (or is close to) steady-state conditions, then growth and retraction events over an entire contour are at statistical equilibrium and the contour length evolves toward an asymptotic value. Longer contours correspond to a stronger filamentation process and to a more effective enhancement of the background diffusion. Therefore, by quantifying the contour length, an estimate of mixing is obtained in terms of an effective diffusion. However, the contour averaging nature of the effective diffusivity is the main limitation of this technique, because mixing variability along an individual contour cannot be resolved.
<figure><img src="image_2.png" /><figcaption>FIG. 1. (a) In a shear flow, where stable and unstable directions of deformation are parallel, a passive tracer initially orientates its gradient parallel to the shear direction, but thereafter the gradient is unaffected by the flow. (b) Conversely, when stable and unstable directions of deformation intersect, the gradient of a passive tracer is exponentially intensified, enhancing the diffusion effect of small-scale turbulence at any time.</figcaption></figure>

Lyapunov exponents and Lyapunov vectors detect local stretching properties of a velocity field with both longitudinal and latitudinal variability. In the case of filamentation events, forward Lyapunov calculations provide intensity and future direction of tracer gradient intensification, whereas backward calculations estimate the present orientation of a tracer gradient that has been passively advected. Strong filamentation is typically associated with large Lyapunov exponents because of the presence of strong chaotic regions. On the contrary, weakly mixing structures, such as the cores of slowly evolving vortices, are often associated with stagnation, low particle separation, and therefore low values of Lyapunov exponents (Joseph and Legras 2002; see also Garny et al. 2007).

Hence, one might expect the effective diffusivity to be large when the Lyapunov exponent is large. In fact, a relationship can be derived from the following heuristic model for blob dilution (Shuckburgh and Haynes 2003). Consider a finite blob of fluid in a domain of size $L_0 \times L_0$. As time $\tau$ goes on, under the repeated action of stretching and folding, the fate of the blob is to invade the whole domain. It does so by elongating as $e^{\lambda\tau}$, and its lateral size is bounded by $L_d$. When the blob fills the domain, its length $L$ is such that $LL_d = L_0^2$ (i.e., $L^2 = \lambda L_0^4/\kappa$); hence, the effective diffusivity is $\kappa_{eff} \approx \lambda L_0^2$.

The actual calculation of effective diffusivity involves replenishing of the tracer by a large-scale gradient but follows a similar scaling<fnref n="4" /> if the contour fills a significant portion of the accessible domain.

b. Particle dispersion and mixing in a shear-dominated flow

There are, however, cases in which particle separation does not always act together with folding, notably in shear-dominated regions within jet cores. In this case, the use of the Lyapunov exponent as a diagnostic of mixing can be misleading.

In a chaotic region, the stirring of a passive tracer is always active and eventually leads to a uniform tracer concentration. In contrast, the asymptotic effect of a pure shear is not the suppression of the tracer gradient but its orientation parallel to the shear direction (Fig. 1). Calling $\phi$ the angle between the gradient direction and the shear direction, calling $\Lambda$ the shear intensity, and solving the advection equation for the case of a linear tracer gradient, one finds

$$ \phi(t) = \arctan \left\{ \frac{1}{\frac{1}{\tan[\phi(0)]} + \Lambda t} \right\}. \quad (9) $$

Once the angle between the tracer gradient and the shear direction is zero, the tracer is in a stationary condition. No mixing occurs at this stage, and this mechanism maintains asymptotically a large-scale tracer gradient. In this sense, shear is more consistent with a transport barrier than a mixing region, even if dissipation of the tracer gradient can occur during the transient or in response to perturbations (e.g., induced by shear instabilities).

On the other side, particles that do not strictly align on a streakline are separated by the shear. There is no mixing associated with this separation if the particles carry the same value of tracer and are thus indistinguishable. In a pure shear flow, the separation is linear and leads asymptotically to a zero Lyapunov exponent. However, over finite time this separation can be large. Moreover, when small random strain with amplitude $\gamma$ is added to the shear, it is shown in appendix B that the Lyapunov exponent no longer vanishes but scales as $\Lambda^{2/3}\gamma^{1/3}$, which is dominated by the large value of the shear. For this reason, Lyapunov exponents may exhibit large values, and a maximum in Lyapunov exponents may appear where the effective diffusivity has a minimum. Examples will be shown later, and other cases of failed relation between effective diffusivity and Lyapunov exponent of this type can be found in simple analytic flows (Shuckburgh and Haynes 2003).

<footnote n="4">⁴ The fact that effective diffusivity is independent of $\kappa$ is fairly well satisfied in practice (Shuckburgh and Haynes 2003).</footnote>
<figure><img src="image_3.png" /><figcaption>FIG. 2. Scatterplot of monthly means of effective diffusivity and (a) Lyapunov exponents $\lambda$ and (b) modified Lyapunov exponents $\lambda \sin^2\alpha$ averaged along equivalent latitude contours for the period 1981–2001 at 350 K. Points are colored in blue and red for tropical (±15°) and extratropical (15°–80°) latitudes, respectively. Density contours enclosing 50% and 25% of the points are indicated.</figcaption></figure>

### c. Particle dispersion from a dynamical systems perspective

The ambiguous interpretation of large values of Lyapunov exponents and mixing is also apparent by looking at filaments from a dynamical systems perspective. By identifying the physical space with a phase space, tracer filament formation has been shown to correspond to the folding of the stable and unstable manifolds of hyperbolic structures that are embedded in the velocity fields. This folding process is typically associated with the presence of homoclinic and heteroclinic intersections of the manifolds of hyperbolic points (see, e.g., Beige et al. 1991; Wiggins 2005; Haller 2000; Lapeyre et al. 2001; Mancho et al. 2004; Shadden et al. 2005). These manifolds are often used as templates of transport properties. In fact, under general conditions, it can be shown that a passive tracer tends to align its gradient orthogonally to the unstable manifolds. When the manifolds fold in thin lobes, as in the case of a perturbed homoclinic and heteroclinic intersection, passive tracers are redistributed in filament patches, enhancing the background diffusion and eventually leading to mixing. Manifolds of hyperbolic points have typically larger stretching rates than the ones found in nearby regions. On the basis of this heuristic assumption, these manifolds are often unveiled as local maxima (ridges) in a map of Lyapunov exponents.<fnref n="5" /> Although the Lyapunov calculation also provides estimates of manifold orientation through the Lyapunov vectors, the Lyapunov

exponent depends entirely on the stretching, not on the geometry of the folding. This complicates the use of the exponent alone for mixing, because the formation of longer and thinner filaments is not necessarily associated with larger Lyapunov exponents. The problem can be exemplified by considering the two cases of unperturbed

<figure><img src="image_4.png" /><figcaption>FIG. 3. Comparison of time-averaged effective diffusivity ($\kappa_{eff}$; solid black line) with the Lyapunov exponents ($\lambda$; dashed line) and modified Lyapunov exponents ($\lambda'$; solid gray line) averaged in time and along equivalent latitude contours for the period 1981–2001 at 350 K. The modified Lyapunov exponents have been multiplied by a factor 2.5 to fit on the same scale. Note that the anticorrelation between the unmodified Lyapunov exponent and the effective diffusivity at the latitudes of the subtropical jets, where transport is dominated by shear.</figcaption></figure>

<footnote n="5">⁵ See, for instance, Koh and Legras (2002) for the manifolds encircling the stratospheric polar vortex.</footnote>
<figure><img src="image_5.png" /><figcaption>FIG. 4. Comparison of time-averaged effective diffusivity ($\kappa_{eff}$; solid black line) with the Lyapunov exponents ($\lambda$; dashed line) and modified Lyapunov exponents ($\bar{\lambda}$; solid gray line) averaged in time and along equivalent latitude contours for each season over the period 1981–2001 at 350 K. The modified Lyapunov exponents have been multiplied by a factor 2.5 to fit on the same scale.</figcaption></figure>

and perturbed homoclinic structures (see, e.g., Beige et al. 1991, their Fig. 1a; Joseph and Legras 2002, their Fig. 1). A tracer released in the vicinity of these structures tends to orientate its gradients orthogonally to the manifolds (with a relaxation time given by the Lyapunov exponent). These dynamics result in two very different tracer redistributions: transport by lobes and formation of longer interface surfaces (with enhanced mixing in presence of background diffusion) for the perturbed case and segregation for the unperturbed case. However, the values of Lyapunov exponents in the two cases do not have to be different. Similarly, the values of the Lyapunov exponents along a manifold are not necessarily larger where the manifold folds and mixing events occur than when the manifolds intersect at low angles, and a tracer gradient can be maintained for longer times. These observations suggest that the Lyapunov exponent, in order to correctly detect mixing events, needs to be complemented by some additional information (the

Lyapunov vectors) for describing the intersection of stable and unstable manifolds. In the following, we will see how to obtain an analytical relation that combines Lyapunov exponents and vectors. This will be done by comparing climatologies of Lyapunov exponents to the effective diffusivity and by studying in more detail the case of shear flows that corresponds to the case in which stable and unstable manifolds intersect at low angles.

### d. Modified Lyapunov exponents

If we assume that the effective diffusivity provides a good parameterization of mixing, then the fact that stretching alone is not always a correct proxy of mixing is apparent in Fig. 2a, which shows no compact relation between monthly means of effective diffusivity and monthly means of Lyapunov exponents averaged around equivalent latitude contours. From Figs. 3 and 4, one can
<figure><img src="image_6.png" /><figcaption>FIG. 5. (left) Effective diffusivity (black) and Lyapunov exponent averaged around equivalent latitude contours (gray). (right) (top) Gridded Lyapunov exponent $\lambda$ and (bottom) modified Lyapunov exponent $\tilde{\lambda} = \lambda \sin^2\alpha$. The time mean is for December 2000–February 2001, at 350 K, the streamlines averaged over the same time period are overplotted.</figcaption></figure>

see that the Lyapunov exponents (dashed line) and effective diffusivity (black solid line) are strongly anti-correlated in the vicinity of the subtropical jets, where transport is dominated by shear (with $\kappa_{eff}$ indicating a minimum, but $\lambda$ indicating a maximum). An examination of gridded maps of Lyapunov exponents (Fig. 5) also indicates that a stretching maxima appears to be located in correspondence with the strongest zonal segment of the subtropical jet. This confirms that the Lyapunov exponent alone can be a misleading diagnostic of mixing in shear-dominated regions.

However, the Lyapunov calculation also contains Lyapunov vectors that provide information on the ge-

ometry of the stretching that a passive tracer undergoes.<fnref n="6" /> As noted earlier, the eigenvector associated with the principal Lyapunov exponent of the backward calculation is known to indicate the orientation of the gradient of a passively advected tracer. Conversely, the eigenvector obtained for the Lyapunov calculation forward in time is the stable direction relative to future gradient amplification. A simple condition for the presence of

<footnote n="6">⁶ Note that the Lyapunov vectors are always obtained from the Lyapunov calculation as a by-product when the exponent is computed as described in section 2, even if this information is often neglected.</footnote>
<figure><img src="image_7.png" /><figcaption>FIG. 6. Time evolution of the two quantities used to fit Eq. (11): (top) the quantity based on the Lyapunov calculations, $L_{eddy}^2 \lambda \sin^2 \alpha$, averaged around contours of equivalent latitude and (bottom) the effective diffusivity $\kappa_{eff}$.</figcaption></figure>

shear is therefore given by the angle $\alpha$ between these two Lyapunov vectors. This is because in shear-dominated regions and large times these two directions tend to coincide, whereas they differ in the case of mixing regions (see appendix A for more details). As a first step toward a Lyapunov diffusivity, we propose a quantity that depends on the particle dispersion in the future (the forward Lyapunov exponent $\lambda$) and also on the angle $\alpha$ between forward and backward Lyapunov vectors, so that the values of strong Lyapunov exponent in shear regions are suppressed:

$$ \tilde{\lambda} = \lambda \sin^2 \alpha. \qquad (10) $$

The factor $\sin^2 \alpha$ modulates the growth rates of the tracer gradient [see Eq. (A4)]. It can be justified analytically considering the dynamics of a fluid parcel in a shear-

dominated flow with random perturbations, as discussed in appendix B. The angle $\alpha$ increases when the thinning of filaments develops. The multiplication by the factor $\sin^2 \alpha$ thus emphasizes the Lyapunov values where the filaments are thinner and suppresses the values where shear dominates; hence, it is expected to have a better correlation with the effective diffusivity.

Returning to Figs. 3 and 4, the modified Lyapunov exponent $\tilde{\lambda}$ (gray line) correctly represents a mixing minimum corresponding with the wintertime subtropical jet (although there are still differences with respect to the effective diffusivity; in particular, the weaker minimum in correspondence with the summertime subtropical jet is not well resolved). An examination of gridded maps of this quantity (Fig. 5) indicates that the inversion of the subtropical extrema is collocated with the most zonally active part of the subtropical jet (cf. the minimum of
<figure><img src="image_8.png" /><figcaption>FIG. 7. Linear coefficients obtained by fitting to Eq. (11): (a) additive coefficient $κ_a$, (b) multiplicative coefficient $b$, and (c) correlation coefficient $r$. Outside the tropics, the correlation is always greater than $r = 0.4$.</figcaption></figure>

$\tilde{\lambda}$ and the maximum of $\lambda$ over the west Pacific). Figure 2b shows that there is a more compact relation in the extratropics (red dots) between the monthly mean modified Lyapunov exponent and effective diffusivity compared to the unmodified version in Fig. 2a. In the extratropics, the correlation coefficient is $r = 0.73$. In the tropics (blue dots), the correlation is still weak (see later comments).

### e. The Lyapunov diffusivity

Supported by the improved correlation between $\tilde{\lambda}$ and the effective diffusivity as well as by theoretical arguments (see appendix B), we define a Lyapunov diffusivity $D_{\lambda}(\theta, \phi)$ as a quantity with the units of diffusion and a linear dependence on $\tilde{\lambda}$:

$$D_{\lambda} \equiv \kappa_a + b L_{\text{eddy}}^2 \tilde{\lambda}, \quad (11)$$

where $κ_a$ and $b$ are proportionality coefficients and $L_{\text{eddy}}$ is a length scale representing the most energetic perturbations. As a first approximation, outside the tropics, one might expect $L_{\text{eddy}}$ to scale as the Rossby deformation radius. However, there are times (e.g., when the flow is dominated by the Asian monsoon anticyclone) when this may be a very poor approximation. In light of this, we set $L_{\text{eddy}}$ to a generic value of 1000 km, reducing at high latitudes (above 70°) as sine of latitude (i.e., with the deformation radius). We propose to find the coefficients $κ_a$ and $b$ by taking the average around equivalent latitude contours of Eq. (11) and fitting it to the effective diffusivity $κ_{\text{eff}}$. Figure 6 shows the temporal evolution of $L_{\text{eddy}}^2 \lambda \sin^2 \alpha$ averaged around equivalent

latitude contours, along with $κ_{\text{eff}}$. A similar temporal evolution can be seen in both quantities, reflecting the good correlation noted earlier.

When evaluating the coefficients in Eq. (11), we compensate for the fact that we are not properly representing the latitudinal variations in $L_{\text{eddy}}$ by determining the values of $κ_a$ and $b$ separately for each equivalent latitude or in other words by allowing these values to represent the unaccounted variation in latitude. Specifically, we take the monthly means of $L_{\text{eddy}}^2 \lambda \sin^2 \alpha$ averaged around equivalent latitude contours and the monthly means of $κ_{\text{eff}}$ for the entire time period; at each equivalent latitude, we use a linear fit to obtain $κ_a(\phi_e)$ and $b(\phi_e)$. It should be noted that in the tropics, where divergent motion forced by (in the analyzed wind field, parameterized) convection is important, the length scale likely gets very small, possibly to that of the grid size, and the Lyapunov exponent may not be very meaningful. In addition, the weak temporal variability of both $L_{\text{eddy}}^2 \lambda \sin^2 \alpha$ and $κ_{\text{eff}}$ in the tropics (see Fig. 6) makes the fitting of Eq. (11) unreliable.

Figure 7 presents the results of the linear fit procedure. At each equivalent latitude outside the tropics ($\pm 15^\circ$), there is a good fit between $L_{\text{eddy}}^2 \lambda \sin^2 \alpha$ and $κ_{\text{eff}}$, with a correlation coefficient always greater than 0.4. For these extratropical equivalent latitudes, the multiplicative coefficient $b$ is typically ~0.5 and the additive coefficient $κ_a$ is typically small. The exception is in the Northern Hemisphere subtropics where, as previously noted, the scaling we have chosen for $L_{\text{eddy}}$ may be erroneous because of the dominance of the
<figure><img src="image_9.png" /><figcaption>FIG. 8. The effective diffusivity $\kappa_{\text{eff}}$ as a function of equivalent latitude (black) and the Lyapunov diffusivity $D_{\lambda}$ averaged around contours of equivalent latitude (gray) for each season averaged over the period 1980–2001.</figcaption></figure>

Asian monsoon anticyclone. The correlation between $L_{\text{eddy}}^2 \lambda \sin^2 \alpha$ and $\kappa_{\text{eff}}$ is weak in the tropics; as previously discussed, the Lyapunov calculation may have less relevance there. Correspondingly, the $D_{\lambda}$ values are dominated by the additive coefficient $\kappa_{\alpha}$ in this region, and the contribution from $L_{\text{eddy}}^2 \lambda \sin^2 \alpha$ is small.

Figure 8 presents the results of the fit (gray) to the effective diffusivity data (black) for each season averaged over the period 1980–2001. It can be seen that, at all times, the linear fit reproduces the effective diffusivity well.

### f. A climatology of local mixing at 350 K

The effective diffusivity for the UTLS region has been analyzed previously by Haynes and Shuckburgh (2000b) and Scott et al. (2003). Here, we focus on the ability of the Lyapunov diffusivity to detect longitudinal variability of isentropic mixing. Figure 9 presents the Lyapunov diffusivity as a gridded map that is time averaged for December 2000–February 2001 and as a climatology for the entire period 1981–2001.

The picture for December 2000–February 2001 demonstrates that there is considerable longitudinal variability in the mixing that arises from the longitudinal variability in the Lyapunov exponents (Fig. 5). At all latitudes, local mixing regions may be characterized by intensities several times stronger than their zonal mean.

The climatological picture highlights some interesting systematic features of the longitudinal variability. The effective diffusivity detects reduced mixing, indicating a barrier to eddy transport at the latitudes of the subtropical jets. The Lyapunov diffusivity uncovers a more complex picture, with sections of mixing minima separated by synoptic-scale regions of vigorous mixing along the subtropical jets.

In the Northern Hemisphere, the jet velocities are strongest in a climatological sense in the western Pacific and North Africa–Middle East and the barrier strength is also strongest in these regions. The strongest mixing in the region of the subtropical jet ($2 \times 10^6 \text{ m}^2 \text{ s}^{-1}$ in the climatological mean) is at the jet exits in the central and
<figure><img src="image_10.png" /><figcaption>FIG. 9. (left) Effective diffusivity (black) and fitted Lyapunov diffusivity averaged around equivalent latitude contours (gray). (right) Gridded Lyapunov diffusivity $D_\lambda$. The time means for (top) December 2000–February 2001 and (bottom) the period 1981–2001 at 350 K. The streamlines averaged over the same time periods are overplotted.</figcaption></figure>

eastern sides of the ocean basins, extending into North America and northern Europe, respectively. These regions coincide with regions of enhanced Rossby wave breaking (Postel and Hitchman 1999). Poleward of the subtropical jet is a region of rather homogeneous mixing at all longitudes. In the Southern Hemisphere, mixing appears more uniform and relatively weaker than in the Northern Hemisphere in the climatological mean. A strong subtropical barrier is observed at all longitudes (strongest over Australia). There is also evidence of a barrier associated with the polar jet at 60°S over the Indian and Atlantic Oceans. Regions of enhanced mixing in the Pacific subtropics and at all longitudes along the coast of Antarctica roughly correspond to regions of

enhanced Rossby wave breaking (Postel and Hitchman 1999; Berrisford et al. 2007). All these longitudinal variations in mixing will be examined in more detail in Part II.

4. Discussion

In this paper, we have introduced a new diagnostic of transport and mixing, the Lyapunov diffusivity. This diagnostic uses the effective diffusivity for quantifying mixing as a function of equivalent latitude and both Lyapunov exponents and vectors for measuring the intensity and geometry of local stretching. The combination of the two diagnostics overcomes the main limitation of
<figure><img src="image_11.png" /><figcaption>FIG. 10. The evolution of a blob of tracer in forward and backward time under the linear application **M**. The directions **G**<sup>+</sup>, **G**<sup>-</sup>, **H**<sup>+</sup>, and **H**<sup>-</sup> are indicated at different times.</figcaption></figure>

the effective diffusivity; it allows mixing variability to be resolved along equivalent latitude contours (which, for the atmosphere, corresponds broadly to the longitudinal direction).

The geometric information coming from Lyapunov vectors is especially relevant in shear-dominated regions, where stretching is orthogonal to the gradient developed by a passively advected tracer. In this case, without a correction by a geometrical factor, the Lyapunov exponents are not correlated in space and/or time with the effective diffusivity. By studying the case of a shear-dominated flow perturbed by a random deformation, we have found an analytical relation for obtaining a local diffusivity from Lyapunov calculations. Fitting this relation to values of effective diffusivity allows a “Lyapunov diffusivity” to be defined. This quantity has the complementary advantages of Lyapunov calculations and effective diffusivity: the ability to resolve local mixing events in both latitude and longitude, along with the parameterization of mixing in diffusion units.

The theoretical analysis we have performed is based on the study of a shear-dominated flow perturbed by a random and weak deformation. This case is especially well suited for the subtropics, but the link between Lyapunov calculations and effective diffusivity that we find can be interpreted in the general case in terms of filament dynamics and is supported by the good correlation with effective diffusivity outside the tropics. The linear fit we perform has the effect of redistributing the effective diffusivity over a longitudinal contour, which is proportional to the value of $L_{eddY}\lambda \sin^2\alpha$. A future goal is to develop the theoretical framework to provide a parameterization of the proportionality coefficients.

By applying the Lyapunov diffusivity to the 350-K isentropic surface, we have been able to characterize local transport and mixing structures in the UTLS region. In Part II, we conduct a detailed examination of the seasonal and interannual variability of this local transport and mixing structure. A future study will consider other regions of the atmosphere: in particular, isentropic exchanges across the tropical pipe in the stratosphere (Plumb 2002; Shuckburgh et al. 2001).

Both the effective diffusivity (Marshall et al. 2006; Shuckburgh et al. 2009b) and the Lyapunov exponents and vectors (Abraham and Bowen 2002; d'Ovidio et al. 2004; Waugh et al. 2006; Lehahn et al. 2007) have been shown to provide robust diagnostics of ocean circulations. This offers the possibility of using the Lyapunov diffusivity to quantify local eddy mixing in the ocean. This would be particularly interesting, because studies using the effective diffusivity based on restricted patches of the ocean (Shuckburgh et al. 2009c) have indicated that there is very significant local variability in the eddy diffusivity.

*Acknowledgments.* F. d'O. has been supported during this research by the Marie-Curie Fellowship (024717-DEMETRA). EFS is supported by a NERC post-doctoral fellowship. We also acknowledge the support of the EU integrated project SCOUT (Contract 505390-GOCECT-2004).

## APPENDIX A

### Lyapunov Exponents and Lyapunov Vectors

We consider here the evolution of an infinitesimal blob of tracer that is assumed to be circular at time $t_0$ (see Fig. 10). Equation (6) shows that this evolution is entirely defined by the linear application **M**, which satisfies, from (5),

$$ \frac{D}{D\tau} \mathbf{M}(x_0|t_0, t_0 + \tau) = \nabla \mathbf{u}(x_0, t_0 + \tau) \mathbf{M}(x_0|t_0, t_0 + \tau), \quad (A1) $$

with $\mathbf{M}(x_0|t_0, t_0) \equiv \mathbf{I}$.

From its definition in (8), we see that the Lyapunov exponent is $\lambda = 1/2 \tau \ln \sigma$, where $\sigma$ is the largest eigenvalue of $\mathbf{M}^T \mathbf{M}$. For an incompressible flow, the other eigenvalue is $1/\sigma$, defining a second Lyapunov exponent as $\lambda_2 = -\lambda$. These two eigenvalues are associated with two normalized orthogonal eigenvectors $\mathbf{G}^+(t_0)$ and $\mathbf{H}^+(t_0)$, satisfying for the former
$$ \mathbf{M}^{\mathrm{T}}(t_0, t_0 + \tau)\mathbf{M}(t_0, t_0 + \tau)\mathbf{G}^{\mathrm{+}}(t_0) = e^{-2\lambda(t_0, t_0+\tau)\tau}\mathbf{G}^{\mathrm{+}}(t_0). \quad (A2) $$

The eigenvector $\mathbf{G}^{\mathrm{+}}(t_0)$ is called the stable direction at time $t_0$. A line element along $\mathbf{G}^{\mathrm{+}}(t_0)$ is shrunk by the flow between $t_0$ and $t_0 + \tau$ and gets aligned with the transverse direction of the blob at time $t_0 + \tau$ given by

$$ \mathbf{G}^{\mathrm{+}}(t_0 + \tau) = e^{\lambda(t_0, t_0+\tau)\tau}\mathbf{M}(t_0, t_0 + \tau)\mathbf{G}^{\mathrm{+}}(t_0). $$

The orthogonal direction $\mathbf{H}^{\mathrm{+}}(t_0 + \tau)$ is the elongation axis of the blob at time $t_0 + \tau$ and the image of $\mathbf{H}^{\mathrm{+}}(t_0)$ that we refrain, however, to call the unstable direction for the following reason; if a line element is initially projected onto $\mathbf{G}^{\mathrm{+}}(t_0)$ and $\mathbf{H}^{\mathrm{+}}(t_0)$ as $\delta\mathbf{x}(t_0) = a\mathbf{G}^{\mathrm{+}}(t_0) + b\mathbf{H}^{\mathrm{+}}(t_0)$, then its counterpart at time $t_0 + \tau$ will be

$$ \delta\mathbf{x}(t_0 + \tau) = ae^{-\lambda\tau}\mathbf{G}^{\mathrm{+}}(t_0 + \tau) + be^{\lambda\tau}\mathbf{H}^{\mathrm{+}}(t_0 + \tau). \quad (A3) $$

Hence, all line elements but those aligned initially with $\mathbf{G}^{\mathrm{+}}(t_0)$ tend to become aligned with $\mathbf{H}^{\mathrm{+}}(t_0 + \tau)$ as $\tau$ increases. In other words, the only distinctive direction at time $t_0$ for future evolution is the stable direction $\mathbf{G}^{\mathrm{+}}(t_0)$, because all other directions are unstable.

Applying the same reasoning for past evolution, which maps the blob at $t_0$ to its previous shape at $t_0 - \tau'$ by $\mathbf{M}^{-1}(t_0 - \tau', t_0)$, a new direction $\mathbf{G}^{-}(t_0)$ can be defined for which a line element is shrunk for backward evolution in time and becomes aligned with $\mathbf{G}^{-}(t_0 - \tau')$ at $t_0 - \tau'$. Correspondingly, $\mathbf{H}^{-}(t_0)$ and $\mathbf{H}^{-}(t_0 - \tau')$ are defined as the axis of expansion for past evolution. We call $\mathbf{G}^{-}(t_0)$ the unstable direction at time $t_0$, because it is the direction along which all tracer elements at time $t_0 - \tau'$, but those aligned with $\mathbf{H}^{-}(t_0 - \tau')$ tend to align when $\tau'$ gets large. All the other directions at time $t_0$ are images of direction close to $\mathbf{H}^{-}(t_0 - \tau')$ at time $t_0 - \tau'$. Hence, $\mathbf{G}^{-}(t_0)$ is the only distinctive relation for past evolution at time $t_0$.

Notice that $\mathbf{G}^{\mathrm{+}}(t_0)$ and $\mathbf{G}^{-}(t_0)$ depend on $\tau$ and $\tau'$, but they converge exponentially in time (Goldhirsch et al. 1987; Legras and Vautard 1996), unlike the Lyapunov exponent, which converges as $1/\sqrt{\tau}$. Notice also that the necessary distinction between past and future evolution is often obscured in the mathematical literature by the fact that Lyapunov exponents are generally introduced for a time-periodic dynamical system, a restriction that is not necessary for the Oseledec theorem (Oseledec 1968) and does not apply to geophysical flows. An extended version of this discussion for N-dimensional flows and examples can be found in Legras and Vautard (1996).

Under suitable conditions, the extensions of $\mathbf{G}^{\mathrm{+}}(t_0)$ and $\mathbf{G}^{-}(t_0)$ at finite distance generate the stable and unstable manifolds of the Lagrangian trajectories containing all the points converging to it in the future or past evolution, respectively. Haller (2000, 2002) discusses how this notion applies to finite durations.

The roles of $\mathbf{G}$ and $\mathbf{H}$ in this discussion are swapped if one considers gradients instead of line elements. As a result, the gradient of a passive tracer tends to align with $\mathbf{H}^{-}(t_0)$ at time $t_0$, and the stable direction for future evolution is given by $\mathbf{H}^{\mathrm{+}}(t_0)$. This alignment is also preconditioning the evolution for $t > t_0$. Assuming alignment with $\mathbf{H}^{-}(t_0)$ at time $t_0$ and using (A3) modified for the gradient by swapping $\mathbf{G}$ and $\mathbf{H}$, the growth rate is then

$$ \frac{d}{dt} |\nabla_C|^2 = \lambda \sin^2 \alpha e^{2\lambda\tau}, \quad (A4) $$

where $\alpha$ is the angle between the stable and the unstable direction (see Fig. 10).

# APPENDIX B

## Stretching and Effective Diffusivity in a Shear-Dominated Flow

A pure shear $U(y) = \Lambda y$ is characterized by a linear growth of the distance between two parcels; hence, the Lyapunov exponent vanishes at large time. However, it is clear that, for a flow with large shear region, the separation between two parcels can grow quickly within the shear region, as long as the parcels are not aligned with the wind. The shear tends to produce such alignment with the wind (or an alignment of the gradient perpendicular to it), but this alignment can be destroyed by an added perturbation. The basic mechanism is that, under pure shear, the alignment is a stable equilibrium for a rotation of the segment formed by two parcels against the shear, whereas it is unstable for a rotation in the same direction as the shear. When the segment is approaching the equilibrium from the stable side, a small perturbation can send it on the unstable side, where the shear tumbles it by a $\pi$ angle toward equilibrium again. It is shown here that the repeated action of this process generates a mean growth of the distance between the two parcels and that it provides a relation between the Lyapunov exponent and the shear, which is the main ingredient for relating the effective diffusivity to Lyapunov calculations.

We consider the case of an incompressible two-dimensional flow dominated by a fixed shear $U = \Lambda y$ with added fluctuating strain and rotation. The instantaneous velocity gradient is
$$F = \nabla u = \Omega \begin{pmatrix} 0 & -1 \\ 1 & 0 \end{pmatrix} + \gamma \begin{pmatrix} \sin 2\phi & -\cos\phi \\ -\cos 2\phi & -\sin 2\phi \end{pmatrix} + \Lambda \begin{pmatrix} 0 & 1 \\ 0 & 0 \end{pmatrix},$$

where $\Omega$, $\gamma$, and $\phi$ are the rotation, strain, and orientation of the strain axis, respectively.

A small blob of tracer can be described by its inertial matrix

$$I(t) = \frac{1}{2N} \int \mathbf{x}\mathbf{x}^\mathrm{T} C(\mathbf{x}, t) d^2\mathbf{x},$$

where $N = \int C(\mathbf{x}, t) d^2\mathbf{x}$ is the total concentration. Under advection by the flow plus diffusion, the inertial matrix evolves according to (Balkovsky and Fouxon 1999)

$$\frac{dI}{dt} = \kappa + FI + IF^\mathrm{T}.$$

but the basic behavior can be readily seen directly from (B1) and (B2); the angle $\Theta$ tends to 0 or $\pi$, and $\rho$ increases to infinity but less than linearly in time. Consequently, the Lyapunov exponent tends asymptotically to zero, whereas it can still be large at finite time. However, the two equilibria are only stable on one side. As the angle approaches 0 or $\pi$ from upper values, a small perturbation induced by additional strain or rotation can induce a jump to the other side of the equilibrium, which is followed by tumbling and stretching of the patch by the shear. Thus, for a perturbed shear flow, the Lyapunov exponent no longer vanishes, and the stretching is mainly performed by the shear, whereas the perturbation triggers the jumps across $\Theta = 0$ and $\pi$. Adding a small diffusion does not change the picture very much, because the equilibrium is then displaced to a small $O(\kappa)$ angle and a finite $\rho$, but with the same stability properties.

To investigate further the behavior of the perturbed flow, we consider the large $\rho$ limit and a random short-lived perturbation modeled as a Wiener process. In this case, the angular equation is

$$d\Theta = 2\tilde{\gamma}^{1/2} dw - \Lambda \sin^2\Theta dt, \quad (B3)$$

with $\langle dw(t) dw(t') \rangle = \delta(t - t') dt^2$, where $\langle \cdot \rangle$ means a statistical average and we have assumed that $\phi$ is random and that the strain and the rotation have equal

A convenient representation of $I$ is $I = R(\Theta)S(\rho)R(-\Theta)$, where $R(\Theta)$ is the rotation of angle $\Theta$ and $S(\rho)$ is a diagonal matrix with diagonal $(e^{2\rho}, e^{-2\rho})$. Here, $\Theta$ is the angle of the main inertial direction of the blob with the x axis, and the equations for $\Theta$ and $\rho$ are

$$\frac{d\rho}{dt} = -\gamma \sin 2(\Theta - \phi) + \frac{1}{2} \Lambda \sin 2\Theta + \frac{\kappa}{2} e^{-2\rho} \quad \text{and} \quad (B1)$$

$$\frac{d\Theta}{dt} = \Omega - \gamma \cos 2(\Theta - \phi) \coth 2\rho + \frac{\Lambda}{2}(\cos 2\Theta \coth 2\rho - 1). \quad (B2)$$

The equation for a line element is obtained for $\kappa = 0$ and in the limit of large $\rho$ [i.e., letting $\coth 2\rho = 1$ in (B2)]. In this case, the angular equation decouples from the equation of the elongation. In the following, we will neglect $\kappa$, but it can be easily reintroduced in the final result. In the simplest case, $\Omega = \gamma = 0$, the full solution is simply given by

$$S[\rho(t)] = R[-\Theta(t)] \begin{pmatrix} 0 & \Lambda t \\ 0 & 0 \end{pmatrix} R(\Theta_0) S(\Theta_0) R(-\Theta_0) \begin{pmatrix} 0 & 0 \\ \Lambda t & 0 \end{pmatrix} R[\Theta(t)],$$

variance. The Fokker-Planck equation for the probability distribution $P(\Theta, t)$ is then

$$\frac{\partial P}{\partial t} = \Lambda \frac{\partial P \sin^2\Theta}{\partial \Theta} + \tilde{\gamma} \frac{\partial^2 P}{\partial \Theta^2}. \quad (B4)$$

The stationary solution for this equation is (Turitsyn 2007)

$$P(\Theta) = K f(\Theta) \int_{\Theta-\pi}^{\Theta} \frac{f(\Theta)}{f(s)} ds, \quad \text{where} \quad (B5)$$

$$f(\Theta) = \exp \left[ -\frac{\chi}{2} \left( \Theta - \frac{1}{2} \sin 2\Theta \right) \right]$$

and the shear to perturbation ratio is $\chi = \Lambda/\tilde{\gamma}$. The periodicity condition $P(\Theta + \pi) = P(\Theta)$ is satisfied by (B5) and the constant $K$ is obtained by the normalization condition

$$2 \int_{0}^{\pi} P(\Theta) d\Theta = 1.$$

For large $\chi$, we can perform explicit calculations using the steepest descent method, and we have $K = c\chi^{2/3}$, where $c = 3^{1/3}2^{-1}\Gamma^{-2}(1/3) = 0.1005...$. In this limit, we write

$$P(\Theta) = c\chi^{2/3} \int_{\Theta-\pi}^{\Theta} \exp\left[-\frac{\chi}{3}(\Theta^3 - s^3)\right] ds. \quad (B6)$$
It is visible from (B6) that the distribution $P(\Theta)$ peaks at a small angle with a mean value scaling as $\chi^{-1/3}$ and the Lyapunov exponent $\bar{\lambda}$ is obtained from

$$ \bar{\lambda} = \left\langle \frac{d\rho}{dt} \right\rangle = \Lambda \int_{0}^{\pi} \sin(2\Theta) P(\Theta) d\Theta = c_3 \Lambda \chi^{-1/3}, \quad (B7) $$

where the constant $c_3$ is numerically obtained as $c_3 \approx 0.37$. In other words,

$$ \bar{\lambda} = c_3 \Lambda^{2/3} \tilde{\gamma}^{1/3}. \quad (B8) $$

Equation (B8) can be used to get a relation between the effective diffusivity and the Lyapunov exponents. In a shear flow, the stable and unstable directions defined in appendix A tend to align with the flow, and the tracer gradients tend to be orthogonal to the flow. In a pure shear flow, when the two directions are aligned and the stable and unstable manifold are the same, a perfect mixing barrier is generated. This barrier breaks under perturbation, and the width of the stochastic zone or the angle between the stable and the unstable manifolds depends on the amplitude of the perturbation, according to Melnikov theory (Wiggins 2003).

In the simple model considered here, the unstable direction aligns with the main inertial axis, whereas the stable direction is obtained by reverting the temporal integration from some time in the future. The stable direction then lies at an angle scaling as $\chi^{-1/3}$ on the opposite side of the shear axis with respect to the unstable direction. Consequently, the angle between the stable and unstable directions is small and scales as $\chi^{-1/3}$.

The tracer gradient is basically orthogonal to the unstable direction, and this is favored in the atmosphere by the existence of a mean pole-to-equator difference for most species and a mean flow along the parallels. To establish the link between stretching and effective diffusivity, we also need to assume a typical spatial scale $L_f$ for the random deformation. In this case, one can argue that the effective diffusivity scales as $\tilde{\gamma}L_f^2$. Because the Lyapunov exponent $\lambda$ scales as $\Lambda^{2/3}\tilde{\gamma}^{1/3}$, we can write the relation as $\kappa_{\text{eff}} \sim \lambda(\tilde{\gamma}/\Lambda)^{2/3}L_f^2$. Taking into account that the angle $\alpha$ between the stable and unstable directions is $\sim(\tilde{\gamma}/\Lambda)^{1/3}$, we obtain

$$ \kappa_{\text{eff}} \sim \lambda \alpha^2 L_f^2. $$

This result is obtained for a small ratio, $\tilde{\gamma}/\Lambda$. To generalize it to the whole flow and all circumstances, we propose that the relation between Lyapunov exponent and effective diffusivity should be

$$ \kappa_{\text{eff}} = \kappa_a + b L_f^2 \bar{\lambda} \sin^2 \alpha, \quad (B9) $$

where $\kappa_a$ and $b$ are experimentally determined quantities.

## REFERENCES

Abraham, E. R., and M. M. Bowen, 2002: Chaotic stirring by a mesoscale surface-ocean flow. *Chaos*, **12**, 373–381.

Aurell, E., G. Boffetta, A. Crisanti, G. Paladin, and A. Vulpiani, 1997: Predictability in the large: An extension of the concept of Lyapunov exponent. *J. Phys.*, **30A**, 1–26.

Balkovsky, E., and A. Fouxon, 1999: Universal long-time properties of Lagrangian statistics in the Batchelor regime and their application to the passive scalar problem. *Phys. Rev.*, **60E**, 4164–4174.

Beige, D., A. Leonard, and S. Wiggins, 1991: Chaotic transport in the homoclinic and heteroclinix regions of quasiperiodically forced two-dimensional dynamical systems. *Nonlinearity*, **4**, 775–819.

Berrisford, P., B. J. Hoskins, and E. Tyrlis, 2007: Blocking and Rossby wave breaking on the dynamical tropopause in the Southern Hemisphere. *J. Atmos. Sci.*, **64**, 2881–2898.

Bithell, M., and L. J. Gray, 1997: Contour lengthening rates near the tropopause. *Geophys. Res. Lett.*, **24**, 2721–2724.

d'Ovidio, F., V. Fernández, E. Hernández-García, and C. López, 2004: Mixing structures in the Mediterranean sea from finite-size Lyapunov exponents. *Geophys. Res. Lett.*, **31**, L17203, doi:10.1029/2004GL020328.

Garny, H., G. E. Bodeker, and M. Dameris, 2007: Trends and variability in stratospheric mixing: 1979–2005. *Atmos. Chem. Phys.*, **7**, 5611–5624.

Goldhirsch, I., P.-L. Sulem, and S. A. Orszag, 1987: Stability and Lyapunov stability of dynamical systems: A differential approach and a numerical method. *Physica D*, **27**, 311–337.

Haller, G., 2000: Finding finite-time invariant manifolds in two-dimensional velocity fields. *Chaos*, **10**, 99–108.

—, 2002: Lagrangian coherent structures from approximate velocity data. *Phys. Fluids*, **14**, 1851–1861.

—, and G. Yuan, 2000: Lagrangian coherent structures and mixing in two-dimensional turbulence. *Physica D*, **147**, 352–370.

Haynes, P., and E. Shuckburgh, 2000a: Effective diffusivity as a diagnostic of atmospheric transport. 1. Stratosphere. *J. Geophys. Res.*, **105** (D18), 22 777–22 794.

—, and —, 2000b: Effective diffusivity as a diagnostic of atmospheric transport. 2. Troposphere and lower stratosphere. *J. Geophys. Res.*, **105** (D18), 22 795–22 810.

Hoor, P., H. Fischer, L. Lange, J. Lelieveld, and D. Brunner, 2002: Seasonal variations of a mixing layer in the lowermost stratosphere as identified by the CO-O$_3$ correlation from in situ measurements. *J. Geophys. Res.*, **107**, 4044, doi:10.1029/2000JD000289.

Joseph, B., and B. Legras, 2002: Relation between kinematic boundaries, stirring, and barriers for the Antarctic polar vortex. *J. Atmos. Sci.*, **59**, 1198–1212.

Koh, T., and B. Legras, 2002: Hyperbolic lines and the stratospheric polar vortex. *Chaos*, **12**, 382–394.

Lapeyre, G., 2002: Characterization of finite-time Lyapunov exponents and vectors in two-dimensional turbulence. *Chaos*, **12**, 688–698.

—, B. Hua, and P. Klein, 2001: Dynamics of the orientation of active and passive scalars in two-dimensional turbulence. *Phys. Fluids*, **13**, 251–264.

Legras, B., and V. Vautard, 1996: A guide to Liapunov vectors. *Proc. Predictability Seminar*, Reading, United Kingdom, ECMWF, 143–156.
—, I. Pisso, G. Berthet, and F. Lefèvre, 2005: Variability of the Lagrangian turbulent diffusion in the lower stratosphere. *Atmos. Chem. Phys.*, **5**, 1605–1622.

Lehahn, Y., F. d'Ovidio, M. Lévy, and E. Heifetz, 2007: Stirring of the northeast Atlantic spring bloom: A Lagrangian analysis based on multisatellite data. *J. Geophys. Res.*, **112**, C08005, doi:10.1029/2006JC003927.

Mancho, A. M., D. Small, and S. Wiggins, 2004: Computation of hyperbolic trajectories and their stable and unstable manifolds for oceanic flows represented as data sets. *Nonlinear Processes Geophys.*, **11**, 17–33.

Mariotti, A., M. Moustaoui, B. Legras, and H. Teitelbaum, 1997: Comparison between vertical ozone soundings and reconstructed potential vorticity maps by contour advection with surgery. *J. Geophys. Res.*, **102** (D5), 6131–6142.

Marshall, J., E. Shuckburgh, H. Jones, and C. Hill, 2006: Estimates and implications of surface eddy diffusivity in the Southern Ocean derived from tracer transport. *J. Phys. Oceanogr.*, **36**, 1806–1821.

Nakamura, N., 1996: Two-dimensional mixing, edge formation, and permeability diagnosed in an area coordinate. *J. Atmos. Sci.*, **53**, 1524–1537.

Oseledec, V., 1968: A multiplicative ergodic theorem. Lyapunov characteristic numbers for dynamical systems. *Trans. Moscow Math. Soc.*, **19**, 197.

Ott, E., 1993: *Chaos in Dynamical Systems*. Cambridge University Press, 478 pp.

Pierce, R. B., and T. D. Fairlie, 1993: Chaotic advection in the stratosphere: Implications for the dispersal of chemically perturbed air from the polar vortex. *J. Geophys. Res.*, **98** (D10), 18 589–18 595.

Pierrehumbert, R., 1991: Large-scale horizontal mixing in planetary atmospheres. *Phys. Fluids*, **3A**, 1250–1260.

Plumb, R. A., 2002: Stratospheric transport. *J. Meteor. Soc. Japan*, **80**, 793–809.

Polvani, L. M., and J. G. Esler, 2007: Transport and mixing of chemical air masses in idealized baroclinic life cycles. *J. Geo-phys. Res.*, **112**, D23102, doi:10.1029/2007JD008555.

Postel, G. A., and M. H. Hitchman, 1999: A climatology of Rossby wave breaking along the subtropical tropopause. *J. Atmos. Sci.*, **56**, 359–373.

Scott, R. K., E. Shuckburgh, J.-P. Cammas, and B. Legras, 2003: Effective diffusivity as a diagnostic of atmospheric transport. *J. Geophys. Res.*, **108**, 4394, doi:10.1029/2002JD002988.

Shadden, S. C., F. Lekien, and J. E. Marsden, 2005: Definition and properties of Lagrangian structures from finite-time Lyapunov exponents in two-dimensional aperiodic flows. *Physica D*, **212**, 271–304, doi:10.1016/j.physd.2005.10.007.

Shuckburgh, E., and P. H. Haynes, 2003: Diagnosing transport and mixing using a tracer-based coordinate system. *Phys. Fluids*, **15**, 3342–3357.

—, W. Norton, A. Iwi, and P. Haynes, 2001: Influence of the quasi-biennial oscillation on isentropic transport and mixing in the tropics and subtropics. *J. Geophys. Res.*, **106**, 327–337.

—, F. d'Ovidio, and B. Legras, 2009a: Local mixing events in the upper troposphere and lower stratosphere. Part II: Seasonal and interannual variability. *J. Atmos. Sci.*, **66**, 3695–3706.

—, H. Jones, J. Marshall, and C. Hill, 2009b: Robustness of effective diffusivity diagnostic in oceanic flows. *J. Phys. Oceanogr.*, **39**, 1993–2009.

—, —, —, and —, 2009c: Understanding the regional variability of eddy diffusivity in the Pacific sector of the Southern Ocean. *J. Phys. Oceanogr.*, **39**, 2011–2023.

Turitsyn, K. S., 2007: Polymer dynamics in chaotic flows with a strong shear component. *J. Exp. Theor. Phys.*, **105**, 655–664, doi:10.1134/S1063776107090245.

Uppala, S. M., and Coauthors, 2005: The ERA-40 Re-Analysis. *Quart. J. Roy. Meteor. Soc.*, **131**, 2961–3012.

Waugh, D. W., E. R. Abraham, and M. M. Bowen, 2006: Spatial variations of stirring in the surface ocean: A case study of the Tasman Sea. *J. Phys. Oceanogr.*, **36**, 526–542.

Wiggins, S., 2003: *Introduction to Applied Nonlinear Dynamical Systems and Chaos*. Springer-Verlag, 843 pp.

—, 2005: The dynamical systems approach to Lagrangian transport in oceanic flows. *Annu. Rev. Fluid Mech.*, **37**, 295–298.