## Abstract

Recent work on the global overturning circulation and its energetics assumes that processes caused by nonlinearities of the equation of state of seawater are negligible. Nonlinear processes such as cabbeling and thermobaricity cause diapycnal motion as a consequence of isopycnal mixing. The nonlinear equation of state also causes the helical nature of neutral trajectories; as a consequence of this helical nature, it is not possible to define a continuous “density” surface that aligns with neutral tangent planes. The result is an additional diapycnal advection, which needs to be accounted for in water-mass analysis. In this paper, the authors take advantage of new techniques in constructing very accurate continuous density surfaces to more precisely estimate isopycnal and diapycnal processes caused by the nonlinear equation of state. They then quantify the diapycnal advection due to each of these nonlinear processes and show that they lead in total to a significant downward diapycnal advection, particularly in the Southern Ocean. The nonlinear processes are therefore another source of dense water formation in addition to high-latitude convection. To maintain the abyssal stratification in the global ocean, these dense water masses have to be brought back toward surface layers, and this can occur by either diabatic or adiabatic processes. Including these nonlinear processes into the advection–diffusion balance, the authors show that observed diapycnal diffusivities are unlikely to be able to support the amount of dense water produced in the global ocean, thus placing more importance on the adiabatic way of bringing the deep waters back to the surface.

## 1. Introduction

In most existing work on ocean circulation, the causes for water-mass transformation in the ocean are assumed to be either surface forcing or small-scale turbulent mixing. It has been known for decades that nonlinearities in the equation of state of seawater can lead to water-mass transformation processes (e.g., cabbeling and thermobaricity), but these processes have usually been assumed to be negligible in much recent work. There have been a few estimates of diapycnal velocities caused by these nonlinear processes, but none of these has been made in the more general context of the global overturning circulation.

Oceanographers tend to use the “density” conservation equation to describe water-mass transformation processes in the ocean. This density conservation equation can easily be derived from the conservation equations of salt and temperature, with the only “nuisance” being that this leads to additional terms describing nonlinear processes such as cabbeling and thermobaricity. In most previous work, including work on using this density conservation equation to understand ocean energetics, these terms have been dropped.

An additional issue arises because of extra terms that need to be taken into account when using this density conservation equation on a continuous density surface. This is due to the nonexistence of continuous density surfaces that globally connect neutral tangent planes (i.e., the local direction of mixing in the absence of diapycnal mixing processes) everywhere on this surface (McDougall and Jackett 1988). Continuous density surfaces can be associated with a density variable, such as potential density or neutral density, which is constant on the surface. These surfaces are mathematically well defined (i.e., continuous), but water parcels moving along these continuous density surfaces encounter buoyant restoring forces.

This ambiguity in defining continuous density surfaces is due to neutral helicity being nonzero in the ocean. Neutral helicity is caused by the pressure dependence of *α*^{Θ}/*β*^{Θ}, where *α*^{Θ} is the thermal expansion coefficient and *β*^{Θ} the saline contraction coefficient. This thermobaric nonlinearity can be written as the thermobaric coefficient *T _{b}*

^{Θ}=

*β*

^{Θ}(

*α*

^{Θ}/

*β*

^{Θ})

*. We can write neutral helicity as*

_{p}where **∇*** _{n}* is the two-dimensional gradient along a neutral tangent plane and Θ is the conservative temperature (McDougall 2003).

Using continuous density surfaces therefore leads to errors associated with the density surfaces having a different slope than the local neutral tangent plane. These slope errors cause a fictitious diapycnal diffusivity (Klocker et al. 2009a,b) and an additional diapycnal advection (Klocker and McDougall 2010), which have to be taken into account when using the density conservation equation on a continuous density surface.

Recently, an algorithm has been introduced that significantly minimizes the slope errors between neutral tangent planes and density surfaces compared to continuous density surfaces previously used, producing what is called an *ω* surface (Klocker et al. 2009a,b). Because these *ω* surfaces describe the isopycnal direction more accurately than other continuous density surfaces, the isopycnal gradients of temperature and pressure, which are necessary for estimating diapycnal velocities due to cabbeling and thermobaricity, are also evaluated more accurately. It is also shown that these *ω* surfaces minimize the additional diapycnal advection due to neutral helicity (Klocker and McDougall 2010). The need for this additional diapycnal velocity is due to the lateral velocity being incorrectly assumed to be along continuous density surfaces instead of along neutral tangent planes.

In this paper, we will quantify the diapycnal advection caused by nonlinearities in the equation of state of seawater (i.e., cabbeling, thermobaricity and the diapycnal motion due to neutral helicity) and show that they cause significant water-mass transformation, particularly in the Southern Ocean. Mostly these processes produce denser water masses and are therefore sometimes described as the “densification” of mixing. In comparison to the formation of Antarctic Bottom Water (AABW) or North Atlantic Deep Water, which are due to surface fluxes, cabbeling and thermobaricity are due to isopycnal mixing leading to a diapycnal velocity due to nonlinearities in the equation of state of seawater. The diapycnal advection due to neutral helicity is a consequence of using continuous density surfaces as a reference frame for water-mass analysis; that is, this diapycnal advection needs to be taken into account when using a surface that has a different slope than that of neutral tangent planes (which is the case with any continuous surface in an ocean with nonzero neutral helicity), but it does not exist when analyzing water-mass transformation locally with respect to neutral tangent planes. All these nonlinear effects are independent of any surface forcing and can therefore be seen as a different path of dense water production in the global overturning circulation.

All the dense water masses produced in the global ocean (i.e., dense water masses formed by high-latitude convection and by nonlinearities in the equation of state of seawater) have to be brought back toward surface layers; otherwise, the ocean would simply fill up with dense fluid (Munk and Wunsch 1998). This can be done diabatically by small-scale turbulent mixing [the classic hypothesis as used in Munk (1966) and Munk and Wunsch (1998)] or adiabatically by isopycnal advection along density surfaces (Toggweiler and Samuels 1995; Webb and Suginohara 2001). The rates of diapycnal advection and diapycnal diffusion are of central importance to the dynamics of the deep ocean (Stommel and Arons 1960a,b; Gargett 1984; Davis 1994).

Here, we will revisit the work by Munk and Wunsch (1998) to show the changes in diapycnal diffusivity profiles and therefore the estimates of the energy necessary to support these diffusivity profiles when taking into account the nonlinearities of the equation of state of seawater. We use these results to argue that the dense waters are unlikely to be brought back to surface layers by diabatic processes alone.

## 2. The density conservation equation

The aim of this section is to use the conservation equations for salinity and conservative temperature (which is proportional to potential enthalpy and represents the “heat content” per unit mass of seawater; see McDougall 2003) to derive the density conservation equation, which can also be written as an equation for the diapycnal advection. The conservation equations for salinity *S* and conservative temperature Θ are (McDougall and Jackett 2009)

These equations have been written in the advective form and with respect to neutral tangent planes so that *e* is the vertical velocity through neutral tangent planes and **V** is the thickness-weighted horizontal velocity obtained by temporally averaging the horizontal velocity between closely spaced approximately neutral tangent planes. Similarly, the salinity and conservative temperature are the thickness-weighted values obtained by averaging between closely spaced pairs of neutral tangent planes (McDougall and McIntosh 2001). In these equations, *h*(*x*, *y*) is the mean thickness between two closely spaced approximately neutral tangent planes. The mixing processes that appear on the right-hand sides are simply isopycnal turbulent mixing of passive tracers (with isopycnal diffusivity *K*) along density surfaces and diapycnal turbulent mixing (with diapycnal diffusivity *D*) across density surfaces; that is, we have not considered double-diffusive convection or double-diffusive interleaving.

By cross-multiplying Eqs. (2) and (3) by *β*^{Θ} and *α*^{Θ}, respectively, and substracting, we obtain what can loosely be called the density conservation equation with respect to neutral tangent planes or the equation for the dianeutral velocity *e*: namely,

where the cabbeling and the thermobaric coefficients (McDougall 1984, 1987b) are

and represent the nonlinear nature of the equation of state of seawater. The first two terms on the right-hand side of Eq. (4) are due to small-scale turbulent mixing and the other terms are due to cabbeling and thermobaricity. The diapycnal advection due to cabbeling can be written as

and the diapycnal advection due to thermobaricity can be written as

This density conservation equation [Eq. (4)] has the advantage compared to the conservation equation of salinity and conservative temperature that the lateral advection does not affect density so that this form of the conservation equation can be used to infer diapycnal advection (which cannot be directly measured in the ocean) from large-scale property distributions (Davis 1994). Note that *e*^{cab} and *e*^{therm} describe diapycnal velocities (which are caused by isopycnal mixing) and can therefore not be seen in dissipation measurements that are used to infer the diapycnal diffusivity *D* [as in Eq. (24)]. That is, the diapycnal velocities *e*^{cab} and *e*^{therm} are independent of the dissipation of kinetic energy in the ocean.

From McDougall and Jackett (1988), we know that it is impossible to connect neutral tangent planes globally to form a well-defined surface in three-dimensional space. Because of this helical nature of neutral trajectories, which is a consequence of neutral helicity being nonzero in the ocean, we have to deal with the additional diapycnal advection *e*^{hel}. The diapycnal advection *e*^{hel} is a consequence of the lateral velocity **V** being assumed to be along a well-defined surface rather than a neutral tangent plane, and it can be written as (McDougall and Jackett 1988; Klocker et al. 2009a,b; Klocker and McDougall 2010)

where **V** is the lateral velocity and **s** is the slope error between the neutral tangent plane and the continuous density surface used,

where **∇*** _{n}* is the gradient along a neutral tangent plane and

**∇**

*is the gradient along a continuous density surface.*

_{a}The term *e*^{hel} would only be zero everywhere if the equation of state of seawater were linear or if neutral helicity in the ocean were zero. We know that the equation of state is nonlinear. Neutral helicity can be close to zero in some regions of the global ocean but can be large in others. Therefore, depending on which regions of the global ocean we work in, *e*^{hel} can be large.

Knowing that it is impossible to construct a continuous density surface that describes neutral tangent planes everywhere, the best solution is to minimize the error involved in using a continuous density surface. This is done by the algorithm described in Klocker et al. (2009a; see also Klocker et al. 2009b), who solve a least squares problem minimizing the slope errors on a continuous density surface, leading to what is called an *ω* surface. Using these surfaces minimizes both the fictitious diapycnal diffusivity *D ^{f}* and

*e*

^{hel}in such a way that

*D*and

^{f}*e*

^{hel}are only due to neutral helicity, whereas in other surfaces

*D*and

^{f}*e*

^{hel}are due to neutral helicity and also errors arising from the way these surfaces are defined. These

*ω*surfaces allow us to calculate the accurate isopycnal gradients of pressure and temperature that are necessary to calculate diapycnal velocities caused by cabbeling and thermobaricity. The diapycnal advection

*e*through an approximately neutral surface is given by (Klocker and McDougall 2010)

^{a}The diapycnal velocity *e* is the one through neutral tangent planes; *e* appears in Eq. (4). In more detail, Eq. (11) can be written as

where the parts involving the diapycnal diffusivity *D* are due to diapycnal diffusion. The last term in Eq. (12), the vertical advection due to temporal changes in the ocean’s hydrography *e*^{tmp}, can be ignored because this velocity is approximately an order of magnitude smaller than the other terms (Klocker and McDougall 2010).

Note that the previous equations describing the diapycnal advection change when the lateral mixing is done along potential density surfaces, as is done in most layered models of the ocean (which usually use a reference pressure of *p _{r}* = 2000 dbar). For these layered models, the equation for the diapycnal advection through potential density surfaces [which is the equivalent of Eq. (4) for the neutral tangent planes] can be written as (IOC et al. 2010)

where *e ^{σ}* is the diapycnal velocity through a potential density surface and

*K*is the isopycnal diffusivity along the potential density surfaces. The

^{σ}*α*

^{Θ}(

*p*),

_{r}*β*

^{Θ}(

*p*), and

_{r}*C*

_{b}^{Θ}(

*p*) are the thermal expansion coefficient, the saline contraction coefficient, and the cabbeling coefficient calculated at the reference pressure

_{r}*p*. The most serious difference between Eqs. (4) and (13) is the absence of thermobaricity in Eq. (13). As we will show later, thermobaricity can cause a substantial amount of water-mass transformation in the global ocean, and the error in ignoring this contribution in layered ocean models using potential density surfaces as their vertical coordinate might be significant. The absence of thermobaricity as a contributor to the diapycnal velocity in layered ocean models has previously been noted by Iudicone et al. (2008). Also missing from layered ocean models is the mean diapycnal advection

_{r}*e*

^{hel}due to the helical nature of neutral trajectories in the ocean. The absence of both thermobaricity

*e*

^{therm}and the diapycnal advection due to neutral helicity

*e*

^{hel}from isopycnal coordinate ocean models is due to the direction of lateral mixing in these models not coinciding with neutral tangent planes.

The density conservation equation is different again when we use potential density surfaces to analyze hydrographic data from an ocean that mixes along neutral tangent planes. According to McDougall (1991), the diapycnal velocity *e ^{σ}* can be written as

where *h ^{σ}*(

*x*,

*y*) is the thickness between two closely spaced potential density surfaces and

*μ*is defined as

The stability ratio of the water column *R _{ρ}* and

*r*are defined as

When hydrographic data are analyzed in potential density coordinates, it is common to see just the terms on the first line of Eq. (14) being used. The last two lines of Eq. (14) arise from (incorrectly) writing the lateral diffusion in the first line of Eq. (14) as occurring along potential density surfaces rather than along neutral tangent planes. Even at the reference pressure when *p* = *p _{r}* and

*r*=

*μ*= 1, the last two lines of Eq. (14) do not reduce to zero but rather to

*T*

_{b}^{Θ}

*K*

**∇**

_{n}Θ ·

**∇**

*, showing that the thermobaric effect remains. Because of this difference between conservation equations for density [Eqs. (4), (13), and (14)], care has to be taken that the correct version is used for one’s application.*

_{n}p## 3. Cabbeling, thermobaricity, and the dianeutral advection due to neutral helicity

The lateral diffusive fluxes of Θ and *S* along neutral tangent planes, −*K***∇*** _{n}*Θ and −

*K*

**∇**

*, have equal but opposite effects on the locally referenced potential density; that is, −*

_{n}S*Kα*

^{Θ}

**∇**

_{n}Θ = −

*Kβ*

^{Θ}

**∇**

*. To illustrate the cause of the diapycnal advection of cabbeling and thermobaricity, consider a situation where the lateral diffusive flux of salt is spatially constant, implying that*

_{n}S**∇**

*· (*

_{n}*K*

**∇**

*) = 0. However, the neutral relationship*

_{n}S*α*

^{Θ}

**∇**

_{n}Θ =

*β*

^{Θ}

**∇**

*implies that the Θ field must obey*

_{n}S**∇**

*· [(*

_{n}*α*

^{Θ}/

*β*

^{Θ})

*K*

**∇**

_{n}Θ] = 0 so that the nonlinearity in the equation of state (i.e., the spatial variation of

*α*

^{Θ}and

*β*

^{Θ}) implies that in this situation there is a nonzero divergence of the lateral heat flux: that is,

**∇**

*· (*

_{n}*K*

**∇**

_{n}Θ) ≠ 0. It is then not surprising to learn that there is a net diapycnal flow caused by these nonlinearities of the equation of state.

### a. Quantifying cabbeling and thermobaricity and the diapycnal advection due to neutral helicity in the global ocean

Here, we will quantify the diapycnal advection caused by nonlinearities in the equation of state of seawater [the diapycnal advection processes on the right-hand side of Eq. (12); i.e., *e*^{cab}, *e*^{therm}, and *e*^{hel}]. For accurate estimates of these processes, it is necessary to compute accurate isopycnal gradients of pressure and temperature. This is achieved by using *ω* surfaces to calculate these gradients. The terms *C _{b}*

^{Θ},

*T*

_{b}^{Θ}, and

*N*

^{2}can easily be calculated on any density surface.

This leaves the isopycnal diffusivity *K*. There has been much discussion about isopycnal diffusivities in the ocean, with values ranging from *O*(10^{2} m^{2} s^{−1}) to *O*(10^{4} m^{2} s^{−1}). These estimates come from tracer release experiments (e.g., Ledwell et al. 1998), Lagrangian statistics of surface drifters (e.g., Sallée et al. 2008), geostrophic flow estimates that were used to advect a numerical tracer to calculate eddy diffusivities (e.g., Marshall and Radko 2006), and a study investigating the relationship of *K* and *D* to the strength of the Southern Ocean meridional overturning circulation (Zika et al. 2009). Estimates such as those by Sallée et al. (2008) and Marshall and Radko (2006) are only valid close to the surface. Tracer release experiments that are used to estimate eddy diffusivities in the deeper ocean only exist for small regional areas, and there is no knowledge about the vertical structure of *K*. Particularly in the Southern Ocean, in which isopycnal eddy diffusivities are clearly dynamically important, there is a lack of reliable estimates of *K* below the mixed layer. The Southern Ocean is also the region with large isopycnal gradients of temperature and pressure due to the outcropping of density surfaces, leading to large diapycnal velocities caused by cabbeling and thermobaricity. Abernathey et al. (2010) and A. C. Naveira Garabato et al. (2009, unpublished manuscript) have also suggested that thermohaline fronts arise where there is little isoneutral mixing, which would lead to an overestimate of the dianeutral velocities due to cabbeling–thermobaricity along fronts. However, according to Abernathey et al. (2010) and A. C. Naveira Garabato et al. (2009, unpublished manuscript), the eddy diffusivity close to fronts is approximately 10^{3} m^{2} s^{−1}, with values several times that farther away from these fronts. In a global sense, our use of *K* = 10^{3} m^{2} s^{−1} would then lead to an underestimates of the dianeutral velocity due to cabbeling and thermobaricity.

Because of this great variance in estimates for *K*, we use a value of *K* = 10^{3} m^{2} s^{−1} (this is also the value that the ocean model used here uses). This choice of *K* is our largest uncertainty in our calculations of *e*^{cab} and *e*^{therm}. However, until further observations are made that improve our understanding of isopycnal diffusivities, we will not be able to decrease these errors in our calculations.

For all following calculations, we will be using model output from the Modular Ocean Model version 4 (MOM4; Griffies et al. 2004) using a standard configuration apart from the inclusion of conservative temperature (McDougall 2003) instead of potential temperature as its conserved temperature variable. This change is not central to the results in this paper. It is necessary to use model data to calculate *e*^{hel} because of the need for lateral velocities. We will also compare these estimates of cabbeling and thermobaricity in MOM4 to estimates from climatological data to confirm that the model physics produce similar results to observational data.

One of the issues with using model output is that, because of a lack of knowledge about realistic values for the isopycnal diffusivity, ocean models are run with an isopycnal diffusivity set in a way to improve the model output. If the choice of the isopycnal diffusivity *K* used in MOM4 is too large, we can expect the isopycnal gradients of properties to be too small and vice versa. If we then assume that the lateral flux of Θ, −*K***∇**_{n}Θ, is independent of the model’s value of *K*, then in fact the strength of the thermobaric diapycnal advection would be independent of *K*, whereas that of cabbeling (being proportional to *K***∇**_{n}Θ · **∇**_{n}Θ) would actually increase as *K* is decreased. Note also that the use of an isopycnal diffusivity *K* is only a valid approach in coarse data/model output in which eddies are not resolved (in these cases, *K* describes the mixing an eddy would cause if it would be resolved). In high-resolution data/model output, the isopycnal diffusive property fluxes should be estimated as −**V**′Φ′, where Φ is any tracer and **V**′ and Φ′ are the perturbations from the mean values caused by eddies.

The diapycnal advection caused by cabbeling *e*^{cab} on an *ω* surface with an average pressure of 1400 dbar is shown in Fig. 1a. This figure shows that cabbeling causes diapycnal velocities on the order of −1 × 10^{−7} m s^{−1} in large areas in the Southern Ocean. The strongest diapycnal velocities due to cabbeling are located in confluence zones of the Antarctic Circumpolar Current (ACC). These regions, such as the Brazil Basin in the Atlantic, the Agulhas retroflection in the Indian Ocean, and the confluence zone off New Zealand, are regions where subtropical water masses and water masses from the Southern Ocean have frontal boundaries, creating large isopycnal gradients of salinity, temperature, and pressure. These large isopycnal gradients then lead to large diapycnal advection caused by cabbeling and thermobaricity. Small regions of enhanced *e*^{cab} are also evident in the northern North Atlantic along the outcropping density surfaces and in the Arctic Ocean. The zonal mean of the diapycnal velocity due to cabbeling between the 27.4 kg m^{−3 }*ω* surface and the 28.1 kg m^{−3 }*ω* surface (averaged along *ω* surfaces), as seen in Fig. 2a, shows large values south of 40°S. The diapycnal velocities due to thermobaricity (on the same surface as used for cabbeling earlier) are shown in Fig. 1b, and the zonal mean of this diapycnal velocity is shown in Fig. 2b. These figures show downward diapycnal velocities of *O*(10^{−7} m s^{−1}) in a continuous band around the Antarctic Circumpolar Current, south of 40°S. Similarly to cabbeling, enhanced values of *e*^{therm} can be seen in the northern North Atlantic. Compared to cabbeling, which always produces negative (downward) diapycnal velocities, thermobaricity also produces positive values (an upward diapycnal velocity), mainly between 20° and 40°S. The diapycnal advection due to neutral helicity *e*^{hel} is shown in Fig. 1c for the *ω* surface with the average pressure of 1400 dbar and in Fig. 2c for the zonal mean. As for cabbeling and thermobaricity, the largest values are south of 40°S, with the main difference to cabbeling and thermobaricity being that, instead of the diapycnal velocity being mainly downward, *e*^{hel} produces patches of banded positive and negative values, both occupying similar areas. The only exception to this is the Pacific sector, where *e*^{hel} is upward.

The sum of the diapycnal advection due to cabbeling, thermobaricity, and the diapycnal advection due to neutral helicity on this same *ω* surface is shown in Fig. 1d. The sum of the zonal mean of *e*^{cab}, *e*^{therm}, and *e*^{hel} is shown in Fig. 2d. In both of these figures, the largest values for the diapycnal velocities caused by these nonlinear processes are surrounding the Antarctic Circumpolar Current and to a lesser extent along the outcropping density surfaces in the northern North Atlantic. Integrating these diapycnal velocities globally, we derive diapycnal transports. A vertical profile of these transports is shown in Fig. 3a for cabbeling (blue), thermobaricity (red), the diapycnal advection due to neutral helicity (green), and their sum (black). Cabbeling produces a downward diapycnal transport of approximately 2 Sv (1 Sv ≡ 10^{6} m^{3} s^{−1}), which decreases with depth, whereas thermobaricity increases from approximately 2 Sv on the shallower surfaces to more than 6 Sv on the denser surfaces. The transport due to *e*^{hel} is smaller than 1 Sv through most of the water column. Note that the vertical structure of the *e*^{therm} and *e*^{cab} diapycnal transports is only due to the changing **∇*** _{n}p* and

**∇**

_{n}Θ with depth (because the isopycnal diffusivity is chosen to be constant). By contrast, the vertical structure of

*e*

^{hel}is caused by the changing lateral velocities with depth and the spatial structure of the neutral helicity field.

The global integral of these diapycnal transports in the ocean model caused by cabbeling, thermobaricity, and the diapycnal advection due to neutral helicity results in a downward diapycnal transport of approximately 6 Sv through all the density surfaces shown here. This diapycnal transport is mainly located in the Southern Ocean. Compared to estimates of AABW production of 8.1–9.4 Sv (using a chlorofluorocarbon budget) or 12.3 Sv (using a mass budget; Orsi et al. 1999), this is a substantial amount of dense water production that should not be ignored. Note that one of the main differences between the production of AABW and the dense water production due to the nonlinear equation of state is that the AABW production is due to surface forcing, whereas the nonlinear dense water production is due to isopycnal mixing and therefore independent of surface forcing. However, even though these dense water masses are formed differently, both the Antarctic Bottom Water and the nonlinear diapycnal processes add to the estimates of the northward transport of Antarctic Bottom Water and Lower Circumpolar Deep Water north of the Antarctic Circumpolar Current. It is important to remember that this estimate of 6 Sv of dense water formed by these nonlinear processes in the ocean model is linearly dependent on the isopycnal diffusivity *K*, which is uncertain and likely to be spatially and temporally varying (for our estimates, *K* is taken constant).

To validate the results from the model output that we used for our calculations, we also calculated the global integral of the transports caused by the nonlinearities in the equation of state for the World Ocean Circulation Experiment (WOCE) climatology (Gouretski and Koltermann 2004), again using *K* = 10^{3} m^{2} s^{−1}. These transports are plotted in Fig. 3b, showing similar results to the model output. The main differences are that thermobaricity seems to be less important on most density layers relative to cabbeling, opposite to the results from the model output. This could be due to the model getting the isopycnal gradients of temperature wrong (for an estimate of these errors, see Sloyan and Kamenkovich 2007) or due to the averaging done in the climatology. The sum of the nonlinear diapycnal transports is similar for both cases, starting with a downward diapycnal transport of about 5 Sv on the lighter density layers and increasing to about 8 Sv on the denser layers. Note that, because of the averaging done to construct the global climatology from sparse measurements, the tracer fields are likely to be smoother than in reality, which would then lead to an underestimate of these nonlinear diapycnal velocities (assuming we are using a correct value for the isopycnal diffusivity).

Here, it is important to note that all these calculations have to be done on accurate density surfaces, where accurate means that the density surface has to be close to describing the direction of neutral tangent planes. Transports due to cabbeling, thermobaricity, and *e*^{hel} are usually larger on density surfaces other than *ω* surfaces. Examples for the change of *e*^{hel} when using different potential density surfaces are shown in Klocker and McDougall (2010). Using the *γ ^{n}* surfaces of Jackett and McDougall (1997), which are the most common density surfaces used recently, care has to be taken when labeling data that have a significantly different water-mass structure to the reference hydrography used in the algorithm by Jackett and McDougall (1997). Because of this,

*γ*performs well for observational data but less well for model data, which sometimes have quite different water masses to today’s ocean. The need for accurate continuous density surfaces for the calculation of diapycnal velocities caused by thermobaricity and cabbeling has been mentioned first by Iudicone et al. (2008), even though in this work the authors use the

^{n}*γ*variable to estimate these processes in the Southern Ocean of an ocean model, leading to diapycnal transports by these processes that are substantially larger than the estimates of this work.

^{n}### b. A simple scale analysis of cabbeling and thermobaricity

Earlier, we quantified diapycnal advection caused by nonlinearities in the equation of state in a three-dimensional hydrography, showing that these processes cause a substantial amount of water-mass transformation. The main uncertainty in these calculations is the isopycnal diffusivity *K*. Now, we will use a simple scale analysis to show the importance of thermobaricity and cabbeling in the Southern Ocean relative to isopycnal diffusion. For this simple scale analysis, the absolute value of the isopycnal diffusivity *K* is not important.

Using the advective conservation statements for *S* and Θ [see Eqs. (2) and (3)] and the relations for the variation of *S* and Θ on the neutral tangent plane (McDougall 1987a),

one finds

Here, *DgN*^{−2}Θ_{z}^{3}*β*^{Θ}(*d*^{2}*S*/*d*Θ^{2}) is the total effect of diapycnal mixing on water-mass transformation: that is, the change of Θ on a neutral tangent plane (McDougall 1987b). This water-mass transformation is therefore dependent on the diapycnal diffusivity *D* multiplied with the curvature of the Θ−*S* diagram.

At a thermoclinic front, the magnitude of the isopycnal diffusion scales as the isopycnal gradient divided by the half-width of the front *L*: that is, |∇_{n}^{2}Θ| ≈ |**∇**_{n}Θ|/*L*. Setting half the isopycnal Θ and *p* contrasts across the front equal to ΔΘ and Δ*p*, respectively (ΔΘ = *L*|**∇**_{n}Θ| and Δ*p* = *L*|**∇*** _{n}p*|), we find from Eq. (19) that

In the Antarctic Circumpolar Current, where (*R _{ρ}* − 1) is of order 1, ΔΘ is about 1 K, and Δ

*p*is about 500 dbar (and in the same sense as ΔΘ), Eq. (20) is about 0.5, implying that the contribution of thermobaricity and cabbeling to water-mass transformation is 50% of the contribution of the isopycnal diffusion term. However, thermobaricity and cabbeling are much more important at causing water-mass transformation in the ACC region than this comparison suggests, because the isopycnal diffusion term changes sign across the front and therefore averages close to zero in the frontal region, whereas cabbeling and thermobaricity contribute a term of the same sign across the whole front. This can be seen in Fig. 4, in which Fig. 4a shows the temperature Θ across a front, Fig. 4b shows its gradient Θ

*, and Fig. 4c shows its second derivative Θ*

_{y}_{yy}. From this figure, it is obvious that the temperature gradient is always one sign, whereas its second derivative has a positive extreme and a negative extreme that cancel when integrated over the frontal width. This discussion indicates that, averaged over the region of the Antarctic Circumpolar Current, cabbeling and thermobaricity are much more effective at causing water-mass transformation (e.g., changes of Θ on density surfaces) than is Laplacian lateral diffusion.

It is also important to note that the amount of cabbeling across a front is not dependent on the sharpness of the front (i.e., **∇*** _{n}*Θ) but on the temperature flux

*F*

^{Θ}across it. As an example to show this, we use the Antarctic Circumpolar Current. In the ACC, we assume the temperature flux along the zonal direction to be zero. The temperature flux in the meridional direction can be written as

and we will take *F*^{Θ} to be directed southward and independent of latitude. Therefore, the amount of cabbeling across this front is proportional to (using *e*^{cab} ∝ *K***∇*** _{n}*Θ ·

**∇**

*Θ)*

_{n}showing that ∫*e*^{cab }*dy* is only dependent on the temperature change ΔΘ and the temperature flux *F*^{Θ} across the front but is not dependent on the isopycnal gradient of temperature **∇*** _{n}*Θ. Therefore, if the temperature change across the front ΔΘ and the heat flux

*F*

^{Θ}are known, the amount of cabbeling should be estimated correctly. The same is true for thermobaricity, based on the known pressure change across a front.

When cabbeling and thermobaricity are analyzed by considering the mixing of two fluid parcels, one finds that the density change is proportional to the square of the property (Θ and/or *p*) contrasts between the two fluid parcels. This leads to the thought that, if an ocean front is split up into a series of many smaller fronts, such as occurs when eddies create filaments along fronts as described by Smith and Ferrari (2009), then perhaps the effects of cabbeling and thermobaricity will be reduced by the square of the number of such fronts. The previous argument shows that this is not the case. Rather, from Eq. (22) we see that the total dianeutral transport due to cabbeling and thermobaricity across the whole frontal region is proportional to the lateral heat flux through the region times either ΔΘ or Δ*p*, the total isopycnal change in temperature or pressure across the frontal region.

## 4. Maintaining the abyssal stratification: The complicated relation between upwelling and dissipation

In the previous section, we saw that in the ocean model approximately 6 Sv; in the climatology, approximately 8–10 Sv of dense water is produced by diapycnal advection as a consequence of nonlinearities in the equation of state of seawater. This nonlinear diapycnal advection is described by several terms of the density conservation equation on a continuous density surface [Eq. (12)]: that is, *e*^{cab}, *e*^{therm} and *e*^{hel}. We now turn our attention to the terms involving the diapycnal diffusivity *D* on the right-hand side of Eq. (12).

There has been substantial effort to quantify the diapycnal diffusivity *D* in the global ocean. These efforts include tracer release experiments (e.g., Ledwell et al. 1998, 2000), microstructure measurements (e.g., St. Laurent and Simmons 2006, and references therein), and CTD strain–ADCP shear measurements (Sloyan 2005; Kunze et al. 2006). Most of these methods give a wide range of results. These measurements are also confined to small regions, certain density layers, or a few transects, making it difficult to extrapolate these diapycnal diffusivities to ocean basins.

The estimates of dense water masses produced and diapycnal diffusivities are very uncertain; however, to maintain the abyssal stratification as it exists in today’s ocean, both sides of the density conservation equation on a continuous density surface [Eq. (12)] have to balance.

If high-latitude convection would exist without interior diapycnal mixing, the ocean would turn into a stagnant pool of cold salty water (Munk and Wunsch 1998). Because the pathways of how fluid returns back to the surface are still mainly unknown, it is easiest to assume a one-dimensional advection–diffusion balance to understand the amount of mixing necessary to upwell a certain amount of dense water and the associated energy needed to sustain this mixing (note that the amount of energy needed to sustain a certain amount of mixing is highly dependent on the stratification of the region where this mixing happens; Sloyan 2006). This one-dimensional advection–diffusion has been used in the seminal papers by Wyrtki (1962), Munk (1966), and Munk and Wunsch (1998). Munk (1966) has fit this balance to density and radiocarbon data to arrive at a diapycnal diffusivity of *D* = 10^{−4} m^{2} s^{−1}, which is necessary to upwell 25 Sv of bottom water. Munk and Wunsch (1998) expanded this work to calculate the power required to sustain the mixing necessary to upwell 30 Sv of bottom or deep water, which they estimated as 2.1 TW(1 TW = 10^{12} W) (using an imprecisely defined set of deep and abyssal ocean layers and not to the upper part of the permanent pycnocline). This balance of both sides of the density conservation equation on a continuous density surface [Eq. (12)] does not include the role of adiabatic motion (i.e., advection along a continuous density surface) and is therefore only one extreme case. Because of our main interest here being the correct description of nonlinear processes leading to diapycnal advection, we will look at the purely diabatic case first and discuss its importance relative to the adiabatic case later.

To discuss the diabatic case, which we also call the classic hypothesis for historical reasons, we will briefly revisit the work by Munk and Wunsch (1998) and show how their results change when calculating budgets along accurate isopycnal surfaces and taking into account the full nonlinear equation of state [using Eq. (12)]. Until now, most work on the one-dimensional advection–diffusion balance has ignored these nonlinear terms, some stating that the equation of state is not far from linear (e.g., Wunsch and Ferrari 2004). Accurate estimates of the power/energy supplied to mixing are necessary, because it is clear that the meridional overturning circulation (MOC) and its associated transports are determined not only by the high-latitude buoyancy forcing but also by the energy available for mixing (Munk and Wunsch 1998; Huang 1999; Gnanadesikan et al. 2005).

### a. The classic hypothesis

Here, we will follow the assumptions of Munk and Wunsch (1998) that 30 Sv of dense water has to be brought back to the surface layers. In their “classic” hypothesis, all this dense water is brought back to the surface layers by diabatic processes (i.e., small-scale turbulent mixing). As in Munk and Wunsch (1998), our region of interest is the water column between 1000 and 4000 dbar. Above 1000 dbar, other processes (e.g., Ekman pumping) clearly become important, causing a more complicated structure of density surfaces, and 4000 dbar is just above the depth of the densest water masses. A schematic of this hypothesis can be seen in Fig. 5, in which *Q*_{BW} is the amount of deep and bottom water that has to be brought back toward surface layers and *Q* is the diapycnal transport through approximately neutral surfaces (*γ ^{a}* surfaces). The downward dianeutral transport due to the nonlinear equation of state implies a corresponding larger

*Q*than in the absence of these nonlinear processes. All the following calculations will be on density surfaces, in contrast to Munk and Wunsch (1998), who did their calculations on isobaric surfaces. We use

*ω*surfaces from

*ω*= 27.4 kg m

^{−3}(which are at a depth of about 1000 dbar in the midlatitudes) to

*ω*= 28.1 kg m

^{−3}(which are at a depth of about 4000 dbar in the midlatitudes) for our calculations.

From Eq. (4), ignoring for a moment thermobaricity, cabbeling, and the diapycnal advection due to neutral helicity, we know that the diapycnal diffusivity *D* causes diapycnal advection *e* at the rate

The last term in this equation is the term used by Munk and Wunsch (1998), which assumes the equation of state to be linear. The difference between these two expressions should not be ignored. Using the relation between the diapycnal diffusivity *D* and the turbulent kinetic energy dissipation rate *ε* (Osborn 1980),

where Γ is the ratio of the turbulent buoyancy flux to the turbulent kinetic energy dissipation rate (which we call the mixing efficiency), we find that Eq. (23) can be written as (McDougall 1988)

In our calculations we will use a mixing efficiency of Γ = 0.2, which is an upper bound for Γ (Osborn 1980). It has been suggested that a constant value for Γ is an oversimplification (Smyth et al. 2001) and that a value of Γ ≈ 0.11 is a better estimate for patchy turbulence and should be used for mixing in oceanic patches (Arneborg 2002); however, because of such a range of estimates of Γ and to simplify comparison to earlier work, we will use Γ = 0.2.

Rearranging Eq. (25) and including the effects of hypsometry (i.e., the inclusion of the area of the continuous density surface into the density conservation equation), we get

or, after rearrangement,

In Eqs. (26) and (27), *A* is the area of the density surfaces. The decreasing area *A* with depth in the ocean implies that we need higher diapycnal diffusivities to bring a certain amount of dense water through a small area on a density surface than to bring the same amount of dense water through a larger area of the density surface, as can be seen from the term (*Aε*)* _{z}* in Eqs. (26) and (27). In the Brazil Basin, the dissipation rate is observed to increase toward the bottom; when examined locally on a cast-by-cast basis, the observed increase of

*ε*toward the seafloor suggests that the diapycnal velocity is downward [see Eq. (25)]. However, we know that the mean diapycnal velocity is positive (i.e., upward) in this area, and we suggest that this seemingly inconsistent result may be explained by the rate of increase in the area of isopycnals with height

*A*in the Brazil Basin [i.e., it is explained by the use of Eq. (27) rather than Eq. (25)]. For the derivation of Eq. (27), see the appendix.

_{z}From Eq. (27), we can see that given the hydrography and the hypsometry, the solution of Eq. (27) for (*Aε*) is subject to an unknown constant of integration; that is, to calculate a vertical profile of turbulent kinetic energy dissipation *ε* or diapycnal diffusivity *D*, which is necessary to bring 30 Sv back to surface layers, we have to set an initial value somewhere in the water column. Estimates of *ε* and *D* have large error bars everywhere in the water column, and we will therefore show different cases of vertical profiles of *ε* and *D* for different initial conditions at *ω* = 28.1 kg m^{−3} (which is at a pressure of about 4000 dbar in midlatitudes).

Results for three of these cases, all upwelling 30 Sv from the dense layers toward surface layers, can be seen in Fig. 6a for the diapycnal diffusivity *D* and Fig. 6b for the turbulent kinetic energy dissipation *ε*. In case 1 (the red curve), *D* is chosen to be zero at the bottom (which is not very realistic but gives a minimum boundary for *D* through the water column that is necessary to bring 30 Sv back to surface layers). In case 2 (the blue curve), *D* is approximately constant throughout the water column. This is the case Munk and Wunsch (1998) chose for their estimates. In case 3 (the black curve), we chose a profile that shows an intensification of *D* at the denser levels.

From these vertical profiles of *D* and *ε*, we calculate the energy *E*, which is necessary to sustain this diapycnal mixing, using

For the three cases described earlier, we calculate *E* ≈ 0.6, 0.7, and 0.8 TW for the red, blue, and black profiles, respectively. These amounts of energy would therefore be necessary to upwell 30 Sv from dense layers to the bottom of the thermocline purely by diabatic processes: that is, by small-scale turbulent mixing. The profile needing 0.7 TW of mechanical energy to be maintained is the case that has to be compared to the estimates by Munk and Wunsch (1998) of 2.1 TW, meaning that in our case only a third of the energy is necessary to bring 30 Sv of dense water back to the surface. Note that our estimate of 0.7 TW does not include the effects of dianeutral advection due to cabbeling, thermobaricity, and neutral helicity. It can be argued that the 6 Sv of downwelling caused by the nonlinear nature of the equation of state acts to affect our best estimate of the 30 Sv of total downwelling in the ocean. On one hand, it can be argued that the 6 Sv of nonlinear downwelling means that the diapycnal turbulent mixing must by itself account for a total of 36 Sv of upwelling. On the other hand, if one takes 30 Sv as a best estimate of the observed flow of bottom water close to their sources, then the 6 Sv of nonlinear downwelling can perhaps be interpreted as implying that only 24 Sv of bottom water is due to air–sea interaction at the poles. The first inference seems the more obvious to the authors. In the case of 36 Sv of dense water, an additional 0.2 TW (i.e., 0.9 TW instead of 0.7 TW) would be necessary to upwell this amount of dense water. However, because the 30 Sv used in Munk and Wunsch (1998) is only a crude estimate and because here we are only interested in showing an accurate methodology for estimating the relation between upwelling and dissipation, we will ignore the dianeutral advection due to the nonlinear nature of the equation of state of seawater for now.

Here, we explain the reasons for the difference of a factor of 3 between the estimates of Munk and Wunsch (1998) and ours, in an order of decreasing importance. The reasons for this difference are as follows:

- The depth-independent diapycnal diffusivity
*D*required to sustain a certain profile*ε*(*z*) of turbulent kinetic energy is given by with Δ*ρ*being the top (1000 dbar) to bottom (4000 dbar) density difference, In contrast to Eq. (30), Munk and Wunsch (1998) used differences of potential density referenced to 4000 dbar between 4000 and 1000 dbar in their calculations instead of Eq. (30), producing a slight underestimate of*E*compared to the equivalent calculation using “adiabatic leveling.” Even if we assume the errors due to adiabatic leveling to be negligible, from Fig. 2 of Munk and Wunsch (1998) this Δ*ρ*would seem to be Δ*ρ*≈ 0.5 kg m^{−3}, whereas in their calculations they use (for no reason we can determine) Δ*ρ*≈ 1 kg m^{−3}, giving a factor of 2 difference in the results between Munk and Wunsch (1998) and ours. Munk and Wunsch (1998) exclude regions poleward of 50°N and 40°S for their calculations. In our calculations, we only exclude the Arctic Ocean and the marginal seas. Because of our choice of

*ω*surfaces, we exclude the regions of dense water production (these dense water masses have*ω*> 28.1 kg m^{−3}). The Southern Ocean is the region in which the highest values of diapycnal diffusivity have been measured (Naveira Garabato et al. 2007), and it is therefore likely that dense water masses are brought back to surface layers in those regions. If all the 30 Sv of dense water would be brought back to surface layers south of 40°S, approximately an order of magnitude less energy would be needed than when bringing these water masses back to surface layers north of 40°S. If we exclude regions as in Munk and Wunsch (1998),*E*is overestimated by approximately 0.2 TW compared to our estimates.We have carried out our analysis on density surfaces, whereas Munk and Wunsch (1998) have done their analysis on isobaric surfaces. When excluding regions poleward of 48°N and 40°S, these differences are small because of density surfaces being relatively flat, but when extending this analysis to high latitudes the error in such an analysis would increase rapidly because of the steep slopes of isopycnals in these regions.

Ignoring the nonlinear terms in the equation of state of seawater [using the last term in Eq. (23) instead of the second term, or equivalently ignoring the second term on the right-hand side of Eq. (25); see McDougall 1988] leads to an underestimate of dissipation necessary to sustain a certain amount of upwelling below 1500 dbar (under typical conditions in the ocean, the second term on the right-hand side of Eq. (25) equals zero at a depth of about 1500 dbar). The difference in diapycnal diffusivity profiles between the linear and nonlinear case can be seen in Fig. 7. Here, the nonlinear case is shown as a black line and two linear cases are shown as red lines. The linear cases are chosen such that the boundary conditions at the lightest–densest density surface align with the nonlinear case. From this figure, one can see that the inclusion of these nonlinearities can make a substantial contribution to the profile of diapycnal diffusivities, depending on the boundary conditions chosen.

Hypsometry does not cause a large change in this analysis, because the density surfaces on the global scale have very similar surface areas.

Some of these points might be more/less important when using this one-dimensional advection–diffusion balance for regional estimates. For example, in some deep-ocean studies involving the upwelling of AABW, hypsometry might play a significant role because of complicated bathymetry (see the appendix for a discussion of the Brazil Basin). Additionally, the nonlinear terms might play a larger role in some regions in the ocean. These modifications to the calculations by Munk and Wunsch (1998) are simple to include in any work using the one-dimensional advection–diffusion balance, leading to more accurate results, and should therefore not be neglected from such.

### b. Possible upwelling scenarios from an observational point of view

In the last section, the one-dimensional advection–diffusion balance has been discussed, assuming that 30 Sv has to be brought back diapycnally from the deep ocean to the bottom of the seasonal thermocline. This is what we call the classic hypothesis because of the seminal work by Munk (1966) and Munk and Wunsch (1998). This classic hypothesis assumes that these dense water masses are entirely brought back to the surface layers diabatically by small-scale turbulent mixing, ignoring adiabatic advection, which is thought to be able to bring these dense water masses back to surface layers along outcropping isopycnals (Toggweiler and Samuels 1995; Webb and Suginohara 2001; Sloyan and Rintoul 2001). When allowing for adiabatic advection along density surfaces, small-scale turbulent mixing processes would only have to bring these dense water masses to the depth of the isopycnals that have a surface expression. Here, we will use estimates of global mean diffusivity profiles to show how much dense water can be upwelled with these diffusivity profiles. If the amount of dense water we can upwell with these diffusivity profiles is less than the amount of dense water produced, we can infer that processes other than small-scale turbulent mixing must play a role.

The global mean diffusivity profiles we use here are as follows:

Estimates inferred from lowered ADCP shear and CTD strain profiles (Kunze et al. 2006): These estimates are likely to be on the lower end of possible diffusivity values because of their technique only parameterizing turbulence due to internal waves and therefore not accounting for the dissipation in overflows, canyons, sills, and other boundary processes.

Estimates from an inverse study of the global ocean (Lumpkin and Speer 2007).

Using these diffusivity profiles and Eq. (27), we estimate how much dense water can be upwelled toward surface layers by diabatic processes only. The best-fit diffusivity profiles for these estimates, calculated from the WOCE climatology (Gouretski and Koltermann 2004), are shown in Fig. 8. The red line shows the best fit for the Kunze et al. (2006) estimates and the black line shows the best fit for the Lumpkin and Speer (2007) estimate. To find how much diapycnal advection is consistent with these estimates of the diapycnal diffusivity, we first assume simplicity and use profiles of constant diapycnal transports through the water column and later we use profiles for diapycnal transports that are given a certain value at the densest layer and linearly decreases to zero transport through the lightest density layer (as used to form the two curves in Fig. 8). The best-fit solution for the Kunze et al. (2006) profile was a constant profile of 4-Sv diapycnal transport or a profile with 8 Sv through the densest layer decreasing to 0 Sv through the lightest layer, whereas the best-fit solution for the Lumpkin and Speer (2007) profile was a constant profile of 15-Sv diapycnal transport or a profile with 30 Sv through the densest layer decreasing to 0 Sv through the lightest layer. Both the constant profile and the profile with a linear decrease toward the surface give very similar results. The energy necessary to maintain these profiles is 0.13 and 0.8 TW, respectively, compared to the 2.1 TW estimated by Munk and Wunsch (1998).

These results therefore show the following:

When fitting diapycnal diffusivity profiles from observations and an inverse model to the hydrography of the WOCE climatology, it is unlikely that there is enough small-scale mixing in the global ocean to upwell all dense water masses back toward surface waters (under the assumption that approximately 30 Sv of dense water have to be upwelled). This shows that it is very likely that other processes, such as adiabatic advection along density surface (Toggweiler and Samuels 1995; Webb and Suginohara 2001) or entrainment into sinking regions (Hughes and Griffiths 2006), play a significant role in the global overturning circulation. This idea is also supported by radiocarbon constraints on upwelling pathways (Gnanadesikan et al. 2007), highlighting the importance of the Southern Ocean upwelling pathway: that is, the adiabatic advection along density surfaces outcropping in the Southern Ocean, which is particularly important for the lower branch of the meridional overturning circulation in the Southern Ocean as shown in an inverse model by Sloyan and Rintoul (2001).

The energy necessary to upwell these dense waters back toward surface waters is likely to be much less than estimated by Munk and Wunsch (1998). This difference is due to (i) the assumption of Munk and Wunsch (1998) that all dense water is upwelled toward surface waters by small-scale mixing (which is only one extreme scenario, giving an upper estimate of energy needed to upwell dense water masses); (ii) the unexplainable choice of top to bottom density difference in Munk and Wunsch (1998); and (iii) several simplifications made by Munk and Wunsch (1998), such as ignoring effects of the nonlinear equation of state on the advection–diffusion equation used in their estimates, using isobaric instead of isopycnal surfaces, etc.

Note that, even though we believe that the results derived in this section give a more accurate picture on the complicated relation between upwelling and dissipation in the global oceans, the one-dimensional advection–diffusion balance used here is an oversimplification of the actual problem because of the complicated and largely unknown paths of water masses, with many built-in assumptions that could significantly change the results, such as the assumption of the mixing efficiency, which could easily change results by a factor of 2. Even though the results for the global advection–diffusion problem discussed here are uncertain, the advection–diffusion equations derived here, taking into account the nonlinear equation of state of seawater, making use of accurate density surfaces (i.e., *ω* surfaces), hypsometry, etc. are definitely more accurate than previously used version of these equations (and not hard to implement) and should therefore be used for oceanographic applications in which an advection–diffusion balance is used to close budgets.

## 5. Summary and conclusions

We have derived a “density” conservation equation from the conservation equations of salinity and conservative temperature. This derivation leads to diapycnal advection caused by the thermobaric and cabbeling nonlinearities in the equation of state. Most oceanographic analysis is done in density space and therefore this extra diapycnal advection has to taken into account when using density as a reference frame.

Using the reference frame of continuous density surfaces to describe water-mass transformation, an additional diapycnal advection process has to be taken into account. This additional diapycnal advection process occurs even in the absence of the dissipation of mechanical energy and is due to the ill-defined nature of neutral trajectories, causing flow through any continuous density surface.

These three nonlinear processes cause a downward diapycnal transport on the order of 6 Sv in the global (model) ocean, with the largest fraction of these transports surrounding the Antarctic Circumpolar Current. This extra amount of dense water produced in the global ocean is not insignificant compared to the amounts of deep and bottom waters produced by high-latitude convection. Although deep and bottom water production is thought of as being due to surface forcing, cabbeling and thermobaricity are due to isopycnal mixing causing diapycnal motion and the diapycnal advection due to the ill-defined nature of neutral surfaces is a consequence of using a continuous density surface as a coordinate in a nonlinear ocean.

Having derived this density conservation equation for a nonlinear ocean, we revisit work done by Munk and Wunsch (1998), who use the one-dimensional advection–diffusion balance to estimate that 2.1 TW is necessary to maintain diapycnal mixing in the ocean to upwell 30 Sv of dense water. Our estimate, using similar assumptions but the full nonlinear equation of state of seawater, is 0.7 TW, a third of what Munk and Wunsch (1998) estimated [to facilitate the comparison to Munk and Wunsch (1998), these estimates ignore the diapycnal advection due to thermobaricity, cabbeling, and the diapycnal motion due to neutral helicity; i.e., we assume that 30 Sv have to be upwelled, as in Munk and Wunsch (1998)]. Most of this difference is due to an error in the calculations of Munk and Wunsch (1998), which has been corrected here, and a smaller part of this error is due to the different assumptions in Munk and Wunsch (1998) and this work. However, even though our estimate includes the full nonlinear equation of state of seawater and is done on density surfaces, the uncertainties are still large. These uncertainties include, among others, the boundary conditions necessary for the calculation of diapycnal diffusivity profiles and the mixing efficiency. A simple one-dimensional model as described here is also an oversimplification of the complicated pathways of dense water back toward surface layers. However, even though the uncertainties in these global estimates of diapycnal mixing are large, the equations for the advection–diffusion balance described here are more accurate than similar equations used before because of the inclusion of nonlinear effects, hypsometry, etc. and might lead to a large difference in results when used for regional oceanographic applications where the pathways of different water masses are better constrained than in the global case.

In addition to estimating the energy necessary to upwell 30 Sv from the deep ocean back toward surface layers, we use diapycnal diffusivity profiles from observations and from an inverse model to see if there is actually enough mixing to upwell the dense water produced in the ocean. From these data, we infer that it is unlikely that all the dense water in the global ocean is upwelled via diapycnal mixing, suggesting that adiabatic processes, as shown by Toggweiler and Samuels (1995) and Webb and Suginohara (2001), or the entrainment into sinking regions, as shown by Hughes and Griffiths (2006), are likely to play a significant role in the global overturning circulation.

Even a decade after the seminal paper by Munk and Wunsch (1998), little more is known about ocean energetics and the pathways of dense water produced in the ocean. Similarly, our knowledge of cabbeling and thermobaricity is limited by the uncertain strength of isopycnal diffusion. We have shown that the nonlinear diapycnal advection processes, cabbeling, thermobaricity, and the diapycnal advection due to the ill-defined nature of neutral surfaces play a significant role in the global overturning circulation and should not be ignored in water-mass analysis. More accurate estimates of the magnitude and spatial variation of eddy diffusivities, which are needed to improve the estimates of dense water produced by these nonlinear processes, will hopefully become available with mixing experiments planned in the near future.

## Acknowledgments

We thank Louise Bell and Lea Crosswell for help with the figures in this paper and Russ Fiedler for his help with the model runs. We also thank Bernadette Sloyan and Jean-Baptiste Sallée for valuable comments on the manuscript. The first author was supported by a joint CSIRO-UTAS Ph.D. scholarship in quantitative marine science (QMS) and a top-up CSIRO Ph.D. stipend (funded from Wealth from Oceans National Research Flagship). This research was supported by the Antarctic Climate and Ecosystems Cooperative Research Centre (ACE CRC). This work contributes to the CSIRO Climate Change Research Program and has been partially supported by the CSIRO Wealth from Oceans National Research Flagship.

## REFERENCES

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

### APPENDIX

#### Hypsometry and Vertical Profiles of Dissipation

The usual relationship between the diapycnal velocity *e* and mixing processes [e.g., Eq. (4)] is written as a local balance. There are some situations where the quantity of interest is the total diapycnal transport across an isopycnal, integrated laterally to the boundaries (side walls) of an ocean basin. For example, it is relevant to understanding the relationship between the upwelling across isopycnals and the dissipation of mechanical energy: for example, as has been done for the deep Brazil Basin by St. Laurent et al. (2001). This appendix extends the usual density conservation, Eqs. (4) and (12), to find expressions for the area-integrated diapycnal advection and the area-integrated dissipation of mechanical energy. We will find that it is possible to have the mean diapycnal advection to be upward while, on many of the casts in the interior, the diapycnal velocity is downward. The key extra ingredient is the proportional rate of change with height of the area of an isopycnal.

The three basic (Boussinesq) conservation equations (continuity, conservation of salinity, and conservation of conservative temperature) can be written as the volume integral of the divergence form over the control volume (see Fig. A1),

Here, the superscripts *u* and *l* refer to the upper and lower interfaces bounding each layer and *h* is the vertical distance (layer thickness) between these bounding density interfaces. The overbar is an average over the horizontal area *A* of the control volume or over the upper or lower approximately neutral surfaces bounding the control volume. The **n** is the outward pointing unit vector in two dimensions, pointing upstream from the “choke point,” where the current meters are located. The diapycnal diffusion terms are actually the volume integral of the general diffusive terms **∇** · (*D***∇***S*) and **∇** · (*D***∇**Θ) for the isotropic diffusivity *D*. The *A ^{u}* and

*A*are the areas (measured in the projected horizontal

^{l}*x*–

*y*plane) of the upper and lower approximately neutral surfaces. In this, we have ignored double-diffusive convection in that the only diapycnal mixing is done with the same diffusivity

*D*for salinity as for conservative temperature. The lateral advection could represent (see Fig. A1) the situation where water is exchanged between a sinking plume and the interior or where fluid flows over a sill.

We multiply Eq. (A1) by the constant values *S*_{1} and Θ_{1} (which we will take to be the average of the salinity and conservative temperatures on the upper and lower interfaces) and subtract these equations from Eqs. (A2) and (A3), respectively, obtaining exactly the same form of equations but now for the anomaly variables *S*′ ≡ *S* − *S*_{1} and Θ′ ≡ Θ − Θ_{1}: namely (the primes are not needed on the diffusion terms, because spatial gradients of the constants are zero),

Now we multiply Eq. (A5) by the thermal expansion coefficient *α*^{Θ} and Eq. (A4) by the saline contraction coefficient *β*^{Θ} and subtract these equations to eliminate the unsteady terms: *α*^{Θ} and *β*^{Θ} are the values averaged over both the upper and lower approximately neutral surfaces. The resulting equation is

The left-hand side of Eq. (A6) is the difference between the upper and lower values of the area integral of the product of the dianeutral velocity *e* and the “perturbation buoyancy” (*α*^{Θ}Θ′ − *β*^{Θ}*S*′), and along each of these surfaces the perturbation buoyancy is very nearly constant, so any correlation between spatial variations of *e* and (*α*^{Θ}Θ′ − *β*^{Θ}*S*′) can be ignored, so that the left-hand side of Eq. (A6) can be approximated as (recall the overbar is the area averaging operator)

Our choice of the constant offset salinity and conservative temperature *S*_{1} and Θ_{1} means that the second line of Eq. (A7) is zero. Because the thickness of a layer is reduced toward zero, we write the average of the dianeutral transport across the upper and lower interfaces 0.5(*A ^{u}*

*e*+

^{u}*A*

^{l}*e*) in Eq. (A7) as simply

^{l}*A*

*e*so that Eq. (A6) becomes

The lateral diffusion term in Eq. (A8) can be shown to be a good approximation to the area-integrated diapycnal transport due to cabbeling and thermobaricity, but the last term in Eq. (A8) is difficult to discount as being unimportant, because it also scales as the cabbeling and thermobaric terms in Eq. (A8). The key to resolving this is to realize that, because we have used constant values of the thermal expansion and saline contraction coefficients, the velocity *e* on the left-hand side of Eq. (A8) is not actually the diapycnal velocity through an approximately neutral surface (this diapycnal velocity involves multiplying by the local values of these coefficients everywhere).

Hence, we need to make the simplification that we treat the equation of state as linear so that the lateral diffusive and lateral advection terms in Eq. (A8) become zero. Equation (A8) is now divided by the thickness *h* and the limit as *h* tends to zero is taken obtaining

On dividing by the area-averaged *N* ^{2}/*g*, we arrive at the form of the equation for the area-averaged diapycnal transport *e* that corresponds to Eq. (4) for the dianeutral velocity *e* on an individual cast: namely,

Apart for the area averaging operator, the only difference between this equation and the diapycnal diffusive terms in Eq. (4) is the presence of the area *A* appearing inside the vertical gradient operator of the diapycnal mixing terms. Perhaps this is more obvious in the expression for the area-average diapycnal velocity *e* obtained by simply dividing (A10) by *A* (ignoring possible correlations between the diapycnal diffusivity and the vertical property gradients),

which is identical to the diapycnal diffusion terms in Eq. (4), except for the additional term *DA _{z}*/

*A*.

Continuing to ignore the influence of the nonlinear equation of state and using Osborn’s (1980) expression that relates the diapycnal diffusivity *D* to the dissipation rate of mechanical energy *ε*, (A11) implies that the area-averaged diapycnal velocity follows from

where we have allowed the mixing efficiency Γ to vary in the vertical.

This simple relation may hold the key to explaining a conundrum that arose in the deep measurements of dissipation in the Brazil Basin where on individual casts the dissipation rate of mechanical energy was observed to decrease with height (*ε _{z}*/

*ε*≤ 0) and yet it was clear that water must be moving upward through isopycnals. Perhaps the proportional rate of the vertical increase of the area of the isopycnals

*A*/

_{z}*A*is sufficiently positive there to ensure that

*e*is positive in Eq. (A12).

## Footnotes

* Current affiliation: Massachusetts Institute of Technology, Cambridge, Massachusetts

*Corresponding author address:* Andreas Klocker, Massachusetts Institute of Technology, 77 Massachusetts Ave., Cambridge, MA 02139. Email: aklocker@mit.edu