## Abstract

The mean ocean circulation near 1000-m depth is estimated with 100-km resolution from the Argo float displacements collected before 1 January 2010. After a thorough validation, the 400 000 or so displacements found in the 950–1150 dbar layer and with parking times between 4 and 17 days allow the currents to be mapped at intermediate depths with unprecedented details. The Antarctic Circumpolar Current (ACC) is the most prominent feature, but western boundary currents (and their recirculations) and alternating zonal jets in the tropical Atlantic and Pacific are also well defined. Eddy kinetic energy (EKE) gives the mesoscale variability (on the order of 10 cm^{2} s^{−2} in the interior), which is compared to the surface geostrophic altimetric EKE showing *e*-folding depths greater than 700 m in the ACC and northern subpolar regions. Assuming planetary geostrophy, the geopotential height of the 1000-dbar isobar is estimated to obtain an absolute and deep reference level worldwide. This is done by solving numerically the Poisson equation that results from taking the divergence of the geostrophic equations on the sphere, assuming Neumann boundary conditions.

## 1. Introduction

In the ocean, the differences in the geopotential values of given pressure surfaces are calculated from density through (vertical) integration of the hydrostatic equation. Horizontal differentiation of the geopotential differences then provides the geostrophic velocity at one level relative to another (e.g., Gill 1982; Wunsch 2008). From the very beginning of physical oceanography, in the late nineteenth century, the quest for a level of known (absolute) motion has thus been at the heart of large-scale ocean circulation studies (e.g., Sverdrup et al. 1942; Defant 1961).

If we knew the geopotential on a given pressure surface (up to an arbitrary constant), we would obtain the absolute geostrophic velocity on this surface, whence everywhere in the ocean. Obviously direct worldwide measurement of the absolute flow at a given level would also solve the reference level problem.

The World Ocean Circulation Experiment (WOCE) float program was meant to supply such a reference (near 1000-m depth) and thanks to the newly developed Autonomous Lagrangian Circulation Explorer (ALACE) float (Davis et al. 1992), it succeeded (at least at a 300-km resolution) in the central and South Pacific (Davis 1998) and the Antarctic Circumpolar Current (ACC) (Gille 2003; Davis 2005). A better resolution but with a much more limited geographical coverage was also made possible by acoustic (RAFOS type) floats in the equatorial Atlantic and in the Brazil Basin (Ollitrault 1999; Ollitrault et al. 2006, Legeais et al. 2012), the North Atlantic subpolar gyre (Bower et al. 2002), and around South Africa (Boebel et al. 2003).

In the absence of large-scale velocity mapping, the determination of the reference level has been addressed indirectly using a variety of hypotheses: 1) adjustment of bottom velocities so that the flow is qualitatively consistent with the distribution of tracers along isopycnal surfaces [Reid 1989, 1994, 1997, 2003; Worthington (1976), and the review of Schmitz (1996a,b); Schmitz and Mc Cartney 1993]; 2) conservation of mass in inverse models applied to the North Atlantic by Wunsch (1978), Wunsch and Grant (1982), Olbers et al. (1985), and Provost and Salmon (1986), to the North Pacific by Roemmich and McCallister (1989), and more recently to the World Ocean by Ganachaud and Wunsch (2000) using WOCE hydrographic data; 3) inverse models using a combination of hydrographic and velocity data (SOFAR floats) by Mercier et al. (1993); 4) restoring of temperature–salinity in ocean general circulation models (OGCMs) to observed data (the models simply adjust their velocity field to the frozen tracer field; Sarkisyan and Ivanov 1971; Holland and Hirschman 1972; Ezer and Mellor 1994; Huck et al. 2008); 5) conservation of large-scale potential vorticity (if the flow is conservative), which has led to the beta spiral method of Stommel and Schott (1977); and 6) the method of Park and Guernier (2001), which assumes isopycnal bottom flow.

Now that the Argo program collects velocity data on the global scale, it has become possible to envision direct determinations of velocities at a reference level (here on an isobaric surface) worldwide and with finer spatial resolution. One can either use these velocities and the thermal wind relation to infer velocity at any other depth or estimate the geopotential Φ on the reference pressure level and use the relation , where *ρ* is density and *p* is pressure, to obtain the geopotential elsewhere in the water column. It is the latter method that we favor here because the geopotential gives major information on the circulation, which is harder to extract from the noisier velocities. The optimal interpolation (OI) of Bretherton et al. (1976) has been widely used in the mesoscale experiments of the 1970s and 1980s to obtain the geopotential. The geopotential is a streamfunction in quasigeostrophic theory (because the velocity field is horizontally non divergent) and may be obtained from the current measurements using the hypotheses of homogeneity and isotropy of the underlying turbulent eddy field. The method relies on model spatial correlation functions valid over the mesoscale array. This method has been adapted and used for the mean large-scale circulation by Davis (1998, 2005) in the Pacific and Indian Oceans, by Gille (2003) in the ACC, and by Bower et al. (2002), Lavender et al. (2000), and Willis and Fu (2008) in the North Atlantic. The elaborate OI method proposed by Davis (2005), also recently used by Katsumata and Yoshinari (2010), is able to enforce the condition of no normal flow at the continental boundaries, but still rely on ad hoc assumptions (e.g., Gaussian covariance functions) to resolve the World Ocean circulation, which has observed scales that vary strongly between the midlatitude interior, the tropics, the continental boundaries, or the proximity of topographic obstacles. Note that a least squares method has also been used by Niiler et al. (2003) to obtain the mean sea level from near-surface drifter velocities.

We propose to obtain the geopotential on the mean 1000-dbar surface directly through the use of an equation of balance, the divergence of the geostrophic momentum equation, which relates the Laplacian of geopotential and vorticity. The result is freed from the a priori fixing of scales for the field to be reconstructed because the method just finds the scales of the mean circulation. Assuming a constant in situ density at 1000 dbar (true within 3 parts per thousand) one also obtains the mean geostrophic pressure on the (constant) geopotential surface close to 1000 dbar.

A careful reprocessing of the Argo float data collected before 1 January 2010 (mainly after 1 January 2000) has been carried out to produce a high-quality deep displacement database covering the World Ocean (Ollitrault and Rannou 2013). A particular attention was paid to the float drifting pressure with the objective that no float should depart more than 100 dbar from the true value. Almost 400 000 displacements are available around 1000 m, from which mean currents were mapped as spatial averages over *O*(100 km)-size boxes and over the whole 10-yr period (see section 2).

We present in section 3 the mean absolute circulation near 1000 m, with a focus on the zonal current component, including seasonal variations close to the equator. Then, we estimate in section 4 the geopotential of the 1000-dbar surface, using a straightforward numerical method. This is followed by a general discussion of the intermediate ocean circulation in section 5. Eddy kinetic energy per unit mass (EKE) estimated from the ANDRO dataset permits to give error bars for the mean currents, which are used in turn to estimate the uncertainty on the computed geopotential. A comparison with EKE at the surface (estimated from altimetry) sheds light on the baroclinic structure of the oceanic variability in the upper 1000 m.

## 2. The ANDRO dataset

Beginning in the late 1990s, the Argo program accumulated, by the end of 2012, more than one million profiles and associated deep displacements. This unprecedented database is freely available on either of the two Global Data Assembly Center (GDAC) websites (http://www.coriolis.eu.org or http://www.usgodae.org/argo/). Argo data come as four different netCDF files known as .meta, .tech, .prof, and .traj files (containing general information, technical ones, temperature and salinity profiles, and surface positions and dates, respectively). However, this is not very user friendly. It is probably why Yoshinari et al. (2006) produced a welcomed deep displacement (or velocity) dataset as a one ASCII file (named YoMaHa’07), regularly updated (Lebedev et al. 2007), containing—besides the surface positions and dates—the drifting depth, the mean velocities (with error estimates), and other useful information. However it appeared that some 10% of the displacements had erroneous parking depths, mainly because only the nominal depth ascribed to a float for its whole mission is preserved in YoMaHa’07. In fact, a float mission may consist of several hundred cycles and the depths at which a given float drifts may vary within its successive cycles. Furthermore, the nominal drifting depth of a float that is found (only) in its .meta file (from which it is copied in YoMaHa’07) may be erroneously entered by the DAC, the in situ measured pressure may also be erroneously decoded, the float pressure sensor may malfunction, or the float may ground.

To eliminate or correct if necessary these parking pressure problems, a new deep displacement dataset, named ANDRO (after a traditional dance of Brittany meaning a swirl), has been produced in the spirit of YoMaHa’07, but containing representative and validated parking pressures (generally an average of the measured values during float drift at depth, and for every cycle). ANDRO ASCII files and the format description are available online (http://wwz.ifremer.fr/lpo/). A description of the work done to obtain ANDRO is given by Ollitrault and Rannou (2013). Presently the ANDRO dataset contains Argo float displacements from before 1 January 2010, and from all Argo data prior to that date collected from the 10 following DACs: Atlantic Oceanographic and Meteorological Laboratory (AOML), Coriolis, Japan Meteorological Agency (JMA), Commonwealth Scientific and Industrial Research Organisation (CSIRO), British Oceanographic Data Centre (BODC), Marine Environmental Data Service (MEDS), Indian National Centre for Ocean Information Services (INCOIS), Korea Ocean Research and Development Institute (KORDI), Korea Meteorological Administration (KMA), and China Second Institute of Oceanography (CSIO). This amounts to a total of 6212 floats contributing 606 119 displacements (29 November 2012 update).

In this paper, we use only the 394 035 displacements found in the 950–1150 dbar layer and whose parking times may vary between 4 and 17 days (see Figs. 1a,b). Here, the parking time is simply given as the difference between the date of first surface fix and the date of last surface fix of the previous cycle. This is always greater than the time actually spent at parking pressure. If we estimate the cumulative float-day numbers within 1° × 1° boxes, worldwide (between 65°S and 75°N), only half of the ocean surface is covered with more than 90 days of data per box. However, the coverage samples more or less uniformly all the oceans, except for south of 60°S around Antarctica and north of 70°N in the Arctic Ocean. If instead, we use a covering by 150-km-diameter disks (centered on half-integer degrees, thus overlapping) 80% of the disks now contain 90 days of data or more (Fig. 2). Although the oldest displacement in ANDRO dates back to September 1995, only in 2000 did the number of floats drifting at a time in the World Ocean began to rise (from 100), reaching almost 3000 in 2007 and thereafter. Thus, the presently selected displacements can be considered to represent the ~1000-dbar circulation over the 10-yr period from January 2000 to December 2009 but note that only 1/10 of the floats sample the 5-yr period from 2000 to 2004.

## 3. Mean currents

Most Argo floats cycle approximately every 10 days and stay at the surface on the order of 10 h to be positioned and transmit data (this is true for the floats using Argos satellite system, not so for the few recent floats using Iridium or OrbCom). Deep float displacements are estimated from the last surface position before sinking and the first surface position after surfacing. Mean deep velocities are calculated as displacements over parking times (Fig. 1b shows that the mean parking time is around 9.5 days).

We now discuss the instrumental errors on the ~9.5-day velocities. Generally ascent float velocity is on the order of 10 cm s^{−1} implying that it takes 2–3 h to rise from 1000 m, but the downward velocity varies widely during descent and between float models. We can nevertheless assume that the sum of the time taken to descend to and ascend from 1000 m is not greater than 10 h, which may contribute a 1 cm s^{−1} error in a 0.5 m s^{−1} uniform shear. Davis and Zenk (2001) argue that the main shear is found in the upper 400 m, implying only a few millimeters per second error. These shear errors, provided in ANDRO, are calculated with uniform shears estimated from the surface velocities and the deep velocities [ appendix A, with more details in Ollitrault and Rannou (2013)]. Here, 90% of the shear errors on the zonal component are less than or equal to 5 mm s^{−1} and 95% are less than or equal to 1 cm s^{−1}. Similar proportions (95% and 99%) are obtained for the shear errors on the meridional component (Fig. 3).

Most floats, whose surface fixes are obtained through the Argos satellite system, wait for 80 ± 60 min on average, at the surface, before being positioned (statistics obtained with 750 floats representative of various conditions). This corresponds to an average 1.5 ± 1.5 km of surface displacement (Ollitrault and Rannou 2013). Because a similar delay and surface displacement occur before diving, their root square sum induces a second error on the order of a few millimeters per second. Of course, the new Iridium floats do not present this error, because they are positioned within a few minutes as soon as they have surfaced, but their number is presently very small in ANDRO.

Figure 4 shows the general circulation current vectors as the averages of the (mostly 9.5 day) mean deep velocities, over the 150-km-diameter disks. To avoid a broadening of the boundary current width, the average current vectors are positioned at the center of the mass of the individual velocity positions found within the disks (the weights being the parking times). Although this is not totally satisfying, because of the overlapping disks and of the boundary current widths [on the order of 50 km for example near Brazil, see Legeais et al. (2012)], this enables us to map 50% (respectively 80%) of the World Ocean near 1000 m with a minimum of 9 (respectively 4) degrees of freedom. This comes from the fact that the Lagrangian integral time scale is on the order of 10 days (e.g., Ollitrault and Colin de Verdière 2002) below thermocline depths and 50% (respectively 80%) of the disks contain more than 180 (respectively 90) float-day numbers. Knowing the velocity variances within each disk, we can estimate the errors on the average current vectors as the velocity standard deviations divided by the square root of the number of degrees of freedom. This is the error caused by mesoscale variability. Velocity variances (or EKE) at 1000 m are on the order of 10 cm^{2} s^{−2} in the interior, reaching on the order of 100 cm^{2} s^{−2} near boundaries or in the ACC (Fig. 5), implying an average current error of 1–3 cm s^{−1} (with 9 degrees of freedom). For comparison, the conservative 1 cm s^{−1} instrumental error estimated previously leads to a 3 mm s^{−1} instrumental error on the average current vectors. Figures 6a–c show the histograms of instrumental errors, mesoscale errors, and their ratios, respectively, justifying the neglect of instrumental errors in the following. Appendix A details the calculation of these errors.

Figures 7–9 are blow ups of Fig. 4 detailing the general circulation for the North Atlantic subpolar gyre, the South Atlantic downstream of Drake Passage, and the ACC flow to the south of Australia and New Zealand. On these three figures, currents with speeds less than 1, between 1 and 5, and greater than 5 cm s^{−1} are shown as displacements over 120 (blue), 45 (cyan), and 15 days (yellow), respectively. Mesoscale error ellipses (with a 0.39 probability inside), centered on arrows, are given for the yellow displacements (and only if the ellipses do not contain the whole displacement arrows). Other blow ups, covering the whole ocean, are available online (http://wwz.ifremer.fr/lpo/Produits/ANDRO/).

The ACC stands out as the strongest mean current near 1000-m depth. Many boundary currents are also revealed clearly: the Agulhas Current and the East Madagascar Current, the East Australian and the Falkland Currents in the Southern Hemisphere, the Kuroshio, the Alaska Stream flowing westward along the Aleutian Islands arc, the Gulf Stream, and the Labrador and the East and West Greenland Currents. The Oyashio and East Kamchatka Currents do not show up in Fig. 4 because, in much of that region, there are less than 90 float days of data within cells of 150-km diameter. Equatorial mean currents are slightly weaker and mainly zonal. Topographic highs like the Reykjanes Ridge (Fig. 7) or lows like the Bounty Trough (southeast of New Zealand between the Campbell Plateau and the Chatham Rise, Fig. 9) are seen to channel the strong flow of the North Atlantic subpolar gyre and the ACC Subantarctic Front, respectively. In Fig. 8, the Zapiola anticyclonic eddy (centered at 45°S, 45°W) stands out while the ACC is deviated northward as it crosses the Mid-Atlantic Ridge (MAR) near the Bouvet Plateau. A similar northward deflection is also visible as the ACC crosses the East Pacific Rise near 50°S, 116°W, and south of Tasmania where the ACC flows north between 60° and 55°S (600 km) before resuming a southeast course (Fig. 9).

The Intermediate Western Boundary Current (IWBC) along the Brazilian coast, which transports Antarctic Intermediate Water (AAIW) into the equatorial Atlantic, is not observed because of its shallower depth near 800 m (Legeais et al. 2012). The Somali Current is not found because it reverses with the monsoons thus giving an almost zero average (Schott and McCreary 2001). Note the local divergence of the flow around 57°N, 52°W in Fig. 7, to be related with deep convection in the Labrador Sea that may reach 1000 m or deeper there (Väge et al. 2009).

### Mean zonal currents

Figure 10 shows the mean zonal velocity component and makes clear the presence of alternate zonal bands between 20°S and 20°N, in the three oceans (although not so well defined in the Indian Ocean) with the strongest signature in the Pacific Ocean. Figures 11 and 12 detail these alternating zonal-mean flows in the Atlantic and the Pacific. Between 5°S and 5°N, ubiquitous zonal-mean flows are shown to extend across the whole width of the Atlantic and Pacific. Note also the eastward flow at 20°–22°S in the Atlantic. In contrast, there is a large area of almost uniform eastward flow, north of the Subarctic Front (found near 40°N) in the Pacific. This corresponds to the North Pacific Current at the surface. Between 35° and 25°S in the Atlantic, the flow is also quasi uniform and westward. This corresponds to the subtropical gyre return circulation.

Figure 13 compares the zonal near-equatorial currents, averaged over the basin width of both the Pacific and Atlantic: extrema are found every 1.5° and 2° of latitude, respectively (shaded areas cover 95% of the variability). In the Indian Ocean, there is some tendency for zonality but the mean zonal currents are hidden by the variability. Note that at the equator proper, the mean zonal currents are not very different from zero (with an annual average), whereas there is actually a strong seasonal cycle as shown in Fig. 14 (maximum monthly zonal velocities can reach 15 cm s^{−1}). In the Indian Ocean, the main period is semiannual because of the monsoon regime. In the Pacific and Indian Oceans the pattern propagation is westward (on the order of 0.5 m s^{−1}), while it is unclear in the Atlantic (with the data at hand) yet. These results agree with those of Cravatte et al. (2012) in the Pacific, except for the Lower Equatorial Intermediate Current (LEIC) that is stronger and better defined in their Fig. 3. It is possible that we do not have enough data in the present ANDRO (Argo data from Cravatte et al. cover until August 2011, ANDRO data until January 2010), but it may be that some floats are drifting at different depths. Actually, both problems have possibly happened in Ollitrault et al. (2006) where the zonal means at 1000 dbar in the Atlantic are a bit stronger than in this study, but the zonal mean was also found ill defined at the equator. This is substantiated by lowered acoustic Doppler current profiler (LADCP) sections (Gouriou et al. 2006) showing instantaneous currents at the equator, westward or eastward flowing, depending on the time of the year (and agreeing with Fig. 14).

Spatial patterns (not shown) of the mean meridional component have scales on the order of a few degrees and do not show coherent structure (except on western boundaries, of course).

Highest values of EKE (see Fig. 5) are found in the retroflection to the south of Cape Town, South Africa, and the confluence region to the east of Mar del Plata, Argentina (reaching 300 cm^{2} s^{−2}). The ACC path is clearly delineated by EKE values greater than 50 cm^{2} s^{−2}. Generally boundary currents are associated with high EKE, but the Labrador and East Greenland Currents and the current loop in the Bounty Trough are exceptions: they are quasi stationary and strongly constrained by the topography.

Now the signature of the Somali Current variability is visible (EKE reaching 200 cm^{2} s^{−1}) and there is a band of high EKE along the equator (associated with the seasonal variations of the Equatorial Intermediate Current, see Fig. 14). One may also notice the ring of high EKE encircling the Zapiola eddy. This almost perfectly barotropic (Saunders and King 1995) and stationary eddy is linked to the Zapiola Rise culminating a few hundred meters above the 5000–5500-m-deep Argentine Basin. Over the rest of the ocean interior, EKE values are on the order of 10 cm^{2} s^{−2}.

Mean kinetic energy over eddy kinetic energy (MKE/EKE) ratios on the 1° × 1° grid (not shown) are smaller than one in the interior. However, the nearly stationary current features mainly near boundaries or constrained by the bottom topography are significant with MKE/EKE greater than one: the Labrador and Greenland Currents, the Zapiola eddy, the eastward North Pacific Current, and the Subantarctic Front.

It is interesting also to compare the surface EKE as deduced from altimetry (thus comprising only the geostrophic part of the surface currents), with EKE at 1000 m (where the currents are expected to be almost geostrophic except at the equator and near the boundaries). The altimetric surface EKE (from data covering the period from January 1993 to September 2011) was kindly provided by Dibarboure et al. (2011). The comparison is meaningful because of similar sampling rates of ~10 days at the surface and 1000 m. It is striking to disclose very similar patterns in both EKE maps, although there are exceptions: a strong surface variability in the Indian Ocean, the Indonesian Throughflow region, the poleward limbs of the South and North Equatorial Currents in the western half of the Pacific, all of which totally absent at 1000 m. This baroclinicity is quantified by computing the *e*-folding depth of EKE in meters shown in Fig. 15. Negative values correspond to EKE at 1000 m greater than EKE at the surface. Equatorial, tropical, and subtropical regions appear very baroclinic with *O*(400 m) *e*-folding depths. However in the South Atlantic around 20°S and in the North Atlantic near western Europe, the variability at 1000 m is equal to that at the surface. In the ACC and subpolar regions, the *e*-folding depth exceeds 1000 m and can reach the bottom occasionally.

## 4. The geostrophic circulation

### a. The estimation of the geopotential at 1000 dbar

We now proceed to obtain the geopotential (on a given pressure surface) or the pressure (on a given geopotential surface) from the mean velocity observations. We use geostrophy because the Rossby number estimated as the ratio of relative vorticity to Coriolis parameter *f* (for latitudes outside three degrees of the equator) is found to be less than 0.1 everywhere and less than 0.01 for 95% of the grid points. The problem is of a global nature in the sense that the geopotential field at one point depends on the geopotential and velocity at other points on the given pressure (or geopotential) surface, respectively. It is a slightly complex operation but well worth the effort because the end product is the deep level of reference required in oceanography.

If *u* and *υ* are the eastward and northward components of the mean velocity field, the geostrophic relationship on a rotating spherical earth reads

Here, Φ is the geopotential of the relevant pressure surface (here typically 1000 db), *λ* is longitude, *θ* is latitude, *a* is the earth radius assumed constant (*a* = 6370 km), and *f* = 2Ω sin*θ*. In this isobaric formulation, the derivatives on the right are carried out at constant pressure. If Φ(*λ*, *θ*) is known at *p*_{0} = 1000 dbar oceanwide, one also obtains, with the hydrostatic relation ∂Φ/∂*p* = −1/*ρ*, *δp* = −*ρδ*Φ on the constant geopotential Φ_{0} (because the vertical excursions of pressure surfaces do not exceed 1 m, *ρ* is the local density). The Φ_{0} may be defined as the mean of Φ over the World Ocean surface, *δ*Φ = Φ(*λ*, *θ*) − Φ_{0} and *δp* = *p*_{0} − *p*(*λ*, *θ*).

with the divergence in spherical coordinates given by

Expanding the divergence in Eq. (2):

with *β* = . The geostrophic flow on a rotating earth is divergent. Equation (3) is the planetary vorticity equation, one of the cornerstones of large-scale ocean circulation theories (see Pedlosky 1996). Although the momentum Eqs. (1a) and (1b) have Rossby number accuracy, this is not the case of the vorticity Eq. (3). Advection of relative vorticity is *O*(*U*/*βL*^{2}) relative to the *β* term and at 1° resolution it may not be small for swift boundary currents. For example, the horizontal divergence mentioned at the Labrador deep convection site is *O*(30 × 10^{−7} s^{−1}) while *βυ*/*f* = 10^{−7}*υ* with *υ* in meters per second (although part of the divergence may be due to nonuniform sampling).

Consider now how to get the geopotential from the velocity. By differentiating Eq. (1a) × cos*θ* with respect to *λ* and Eq. (1b) with respect to *θ*, one obtains

This equation is nothing but a special form of the divergence of the horizontal momentum equation. In quasigeostrophic models, a similar equation (with *f* constant) is used to obtain the streamfunction “by inverting” relative vorticity. An important difference with quasigeostrophy, however, is that here the geopotential is not a streamfunction because the flow is horizontally divergent [see Eq. (3)]. Equation (4) is the equation of balance for planetary geostrophic dynamics. The regimes of validity of this equation have been discussed by Gent and McWilliams (1983) who show that it conserves the Rossby number accuracy of the momentum equations. We do not use Eq. (3) but only Eq. (4) in what follows to get Φ from **u**. This Poisson equation requires boundary conditions at the edges of the domain under consideration. Near solid boundaries, the geostrophic relation becomes a poor approximation for the velocity component normal to the coast. Because it goes to zero at the coast, there is no Coriolis force left to balance the alongshore pressure gradients and nonlinear terms must be considered. On the other hand, the velocity along the coast can still be assumed to be geostrophic to a good approximation. A simple scaling of the momentum equation in the direction parallel to the coast shows that the advective term is on the order of relative to the Coriolis term, whereas it is in the direction normal to the coast where and are horizontal scales of motion respectively normal and parallel to the coast and is the (large) velocity parallel to the coast. Because one expects to be small, the advection terms of the momentum equation in the direction normal to the coast are smaller by the ratio . Hence, the alongshore velocity may remain geostrophic to a high degree of accuracy. Indeed, the Florida and Gulf Stream transports have been estimated early using cross-shore hydrographic sections (e.g., Wüst 1924). Thus, if **n** is the outward unit normal at the edge of the domain and **s** the tangential vector such that **s** = **k** × **n**, (with **k** the upward unit vector), the appropriate geostrophic boundary condition to be used in conjunction with Eq. (4) reads

with and = (*fu*, *fυ* cos*θ*).

The condition [Eq. (5)] is a Neumann-type boundary condition for the elliptic Poisson problem [Eq. (4)]. Observed along-coast velocities are used to find the normal pressure gradient that is imposed to the interior. Imposing instead a Dirichlet condition, that is a constant geopotential would assume that geostrophy remains valid for velocities normal to the coast, which is generally not true.

The solution is defined up to an arbitrary constant provided a consistency relation exists between the forcing [Eq. (4)] and the boundary condition [Eq. (5)]. This consistency relationship is obtained by integrating Eq. (4) over the domain *D*. In the two-dimensional (*λ*, *θ*) Cartesian space Eq. (4) can be rewritten as

Using the divergence theorem on the left-hand side and Green’s theorem on the right-hand side, the domain integral becomes

Use of the boundary condition [Eq. (5)] therefore guarantees that the compatibility condition [Eq. (7)] is satisfied. In practice, the interior of *D* will be defined by those points where the Cartesian vorticity has been estimated, which are adjacent either to a solid boundary or to a data hole (details in appendix B). Note that the problem has no singularity at the equator and can therefore be solved directly for the World Ocean. Of course as the equator is approached, the Coriolis acceleration vanishes and so does the pressure gradient. The pressure gradient that the equatorial flows demand has to come from the neglected terms namely the advection of momentum but this special case is postponed to a later study.

The implementation of the numerical method to solve Eq. (4) with boundary conditions [Eq. (5)] for the World Ocean is fully described in appendix B.

The geopotential (m^{2} s^{−2}) so estimated on the mean 1000-dbar surface is given in Fig. 16 as *δz* = 100Φ/*g*_{c}, that is, in centimeter height (with *g*_{c} = 9.806 65 m s^{−2}). If preferred, this can also be viewed as the geostrophic pressure in millibars on the (constant) geopotential surface close to 990-m depth [a better approximation is *p* − *p*_{0} (mbar) = 1.013*δz* (cm), using a constant *ρ*_{0} = 1033 kg m^{−3}]. The “projection” of the Argo-derived mean velocities on the geostrophic mode is efficient as the correlation between the Argo velocities and the geostrophic velocities reaches 0.86 in the meridional direction (Fig. 18a) and 0.96 in the zonal direction (Fig. 18b). Of course this is true only outside an equatorial band (here 3°S–3°N). This geostrophic “projection” has a further advantage: the Eulerian-mean velocities reconstructed from the Argo float displacements can have large statistical errors in strong eddy regions, which in turn degrade the expected geostrophic properties of a true mean circulation. The hypothesis that we put forward is that by projecting this noisy mean on the geostrophic mode, the sampling errors decrease significantly.

### b. Results

One surprising result is the flatness of the pressure surface equatorward of a line originating at 20° of latitude on the western side and ending at 30° on the eastern side in either hemisphere. One could conclude that the mean circulation is very weak in these regions, a statement that includes all basins. However, outside the 3°S–3°N band the alternating jets shown in Figs. 11 and 12 are still present but not revealed by the contouring interval used. Of course the surface must be flat close to the equator but the model is invalid there and the zonal jets described previously at low latitudes in the Pacific and Atlantic are therefore excluded from the present discussion. The whole of the north Indian Ocean has a very weak circulation.

Note the large-scale anomaly on the order of 5 cm between the tropical Pacific and Atlantic Oceans. The sign is similar to what is found at the surface (Niiler et al. 2003; Maximenko et al. 2009) but the amplitude is less. It is a remarkable fact that velocity observations capture this geopotential height difference between the two oceans.

Outstanding features are the western boundary current’s southern recirculations of the North Atlantic and North Pacific subtropical gyres. Both have a pressure high of about 10 cm. Their zonal extent is roughly similar but the meridional one is smaller in the North Atlantic (30°–37°N) than in the North Pacific (25°–35°N). The low pressure found to the north of the Gulf Stream defines its northern recirculation (Hogg 1983). Note how the method captures the strong front between the Gulf Stream and the waters of Labrador origin. Farther north, the Atlantic subpolar gyre is two times stronger than its Pacific counterpart (15 versus 7 cm troughs). In the Atlantic (see Fig. 17) one also discerns a small gyre near 40°N, 42°W, the anticyclonic flow around the Azores Plateau (note the absence of the Azores Current) and an extended cyclonic circulation (5 cm) associated with the Mediterranean Outflow.

The South Atlantic and South Pacific Oceans show similar broad, weak subtropical gyres centered about 35° and 40°S, respectively. They are western intensified with a 15-cm pressure high, but there is no boundary current along the Brazilian coast, which would transport water at intermediate depth into the equatorial Atlantic. These gyres are bounded to the south by the northward front of the ACC at 45°S in the Atlantic and 50°S in the Pacific. Actually, the strongest subtropical gyre is found in the western half of the South Indian Ocean bounded to the South by the southern branch of the Agulhas retroflection (near 37°S). The East Madagascar Current feeds directly into the swift Agulhas Current. The strong anticyclonic Zapiola gyre in the South Atlantic, centered at 45°S, 45°W shows a high pressure anomaly that reaches 30 cm.

Of course the distinctive element of Fig. 16 is the strong Antarctic Circumpolar Current. The pressure drop across the current can reach 90 cm at selected longitudes. Its zonal geostrophic transports (per unit depth, at 1000 dbar) from the most southern Argo observations to the tips of South Africa, South Australia, and South America all tend toward a maximum of 8 ± 0.2 × 10^{4} m^{2} s^{−1} (reached at 21°E, 136°E, and 62°W). This constancy is indicative of rather weak communication from the Southern Ocean to the Atlantic, Indian, and Pacific at that depth. However, the meridional geostrophic transports that enter the Atlantic, Indian, and Pacific Oceans from the south are not well defined and cannot be used as a check. Their determinations should improve when all western boundary currents are well sampled by the Argo floats. Starting from 45°S at the position of the Zapiola gyre in the South Atlantic, the ACC runs eastward avoiding northward the Bouvet Plateau (near 8°W) and Crozet Island (50°E). From there, it moves in the southeast direction (passing north of the Kerguelen Plateau). Downstream of the Campbell Plateau (to the southeast of New Zealand), the northern flank of the ACC veers strongly northward. In the Pacific, it continues in the southeastward direction meandering once more northward over the East Pacific Rise to reach 60°S at Drake Passage. Downstream of Drake Passage the northern part of the ACC branches out as the Falkland Current, which turns sharply to the north from 55° to 40°S. Actually, these features of the ACC were already observed in Fig. 4, but the “geostrophic” pressure enables to continuously follow the ACC flow at 1000 m. The southeastward orientation of the ACC found in the Indian and Pacific Oceans and the two abrupt northward excursions, east of the Campbell Plateau and Drake Passage, are major features of the circulation at 1000 m.

The geopotential height confirms the strong topographic coastal influences that are observed in subpolar regions. There the steady flows are channeled in the direction of coastal Kelvin waves (Greenland and Labrador coasts, Aleutian Island arc, Reykjanes Ridge, or Bounty Trough). The interaction of subpolar currents with topographic obstacles in the interior is best observed in the ACC. The mean flow often tries to avoid obstacles such as the Crozet and Kerguelen Plateaus or aims to find gaps in the ridges. When forced across a ridge, the mean ACC shows significant northward excursions toward greater depths (Mid-Atlantic Ridge, Macquarie Ridge, and East Pacific Rise).

## 5. Discussion

### The mean circulation

The early work on the middepth circulation of the World Ocean has been reviewed by Reid (1981). Many flow patterns inferred from density and tracer fields and relative dynamic height are confirmed by the present direct measurements that bring more spatial resolution. The weaker flows, however, are often different.

#### 1) North Atlantic

The circulation found herein can be compared with the many basinwide and regional circulation patterns. Schmitz and McCartney (1993) have synthesized these. The tight recirculation south of the Gulf Stream is a well-known feature described early by Worthington (1976) for the lower thermocline layer (4°–7°C). To rationalize his scheme Worthington (1976) argued that the Sargasso Sea gyre at that depth must be confined to the north of the tongue of the Mediterranean Water (MW). What is also apparent here is the lack of a Gulf Stream signature south of 30°N. Martel and Wunsch (1993) show a 1000-m recirculation similar to ours in geographic extent but with a pressure high half as large. The circulation of Figs. 16 and 17 shows that the northward-flowing branch connects the Gulf Stream to the anticyclonic northern gyre between 40° and 45°N at variance with Worthington (his Fig. 24), which presents a northern gyre disconnected from the Gulf Stream system. The 300-km-wide anticyclone centered at 40°N, 42°W is the deep signature of the Mann (1967) eddy. The Gulf Stream then rejoins the strong eastward-flowing Polar Front Current at 50°N that delineates the southern branch of the subpolar gyre. Clarke et al. (1980), Schmitz and Mc Cartney (1993), and Reid (1994) have shown similar branching of the Gulf Stream with the northern gyre from hydrography and tracers, and Wunsch and Grant (1982) and Martel and Wunsch (1993) from inverse models. The Polar Front’s flow in Reid (1994) is weaker than here. South of the Gulf Stream recirculation, the scheme of Figs. 16 and 17 conforms with Worthington’s figure because no mean circulation appears! Worthington (1976) stressed from the nearly zonal distribution of salinities in the tropical Atlantic (for instance, the 34.7-psu contour at 10°N) the weakness of exchanges with the South Atlantic for his lower thermocline layer. However, in the western part the isohalines show some northward excursions and this is associated here with the existence of a nearshore northward mean flow and western-intensified eddy activities between 0 and 10°N. Such a northward motion present in Reid’s picture is only visible with the original currents (Fig. 4). The associated transport appears at 5 Sv (1 Sv ≡ 10^{6} m^{3} s^{−1}) for the layer 7°–12°C in Schmitz and McCartney (1993).

Reid (1994) shows a large-scale, weak cyclonic cell oriented to the northeast along the axis of the MW tongue. There is also a low pressure observed here near the latitude of Gibraltar. The associated surrounding cyclonic flow first follows the Portuguese slope [deviated there by the Coriolis force, Lacombe and Tchernia (1960)] but it does not extend beyond the longitude of the Azores at the difference of Reid (1994). Such a flow coincident with the MW salinity anomaly suggests that the MW is not a passive tracer. Such a hypothesis has been considered early by Kawaze and Sarmiento (1986) and Arhan (1987) who invoked the existence of double diffusive convection.

Previous float experiments have been carried out in the subpolar gyre by Bower et al. (2002) and Lavender et al. (2000) but none covered the present 1000-m depth. The subpolar gyre is best described on Fig. 16 as a main trough at 54°N, 32°W surrounded by the signature of the strong Polar Front to the south and cyclonic currents intensifying near the coasts. Because of the grid spacing (1°) the flow cuts partly through the Reykjanes Ridge. Using float velocities from several different depths, Lavender et al. (2000) were able to reconstruct velocity information at a given depth using geostrophic shears from climatology. Their geostrophic pressure at 700 m (reconstructed with optimal interpolation) shows similarities with our 1000-m map. The circulation reconstructed by Bower et al. (2002), which comes from isopycnal floats drifting between 600 m and the surface (along the 27.5 kg m^{−3} potential density surface), shows nevertheless the same general flow characteristics of the subpolar gyre. This confirms simply the rather weak depth dependence of the subpolar flow (see Fig. 15).

#### 2) North Pacific

Here, the Kuroshio is well defined at 35°N and its extension and its southern recirculation organize all the circulation between 20° and 40°N and up to 150°W—features that are also present in Reid (1997). He also includes currents not found here such as an eastward current branch between 20° and 30°N, a western cyclonic cell to the south of 30°N, and several elongated zonal gyres to the south of 20°N. Reid (1997) also shows a western connection with the South Pacific through a northward flow, which is not found here. Although the pressure signature of zonal flows weakens as the equator is approached, the alternating zonal jets are still captured by the geostrophic method outside the 3°S–3°N equatorial band.

Both North Pacific and North Atlantic circulations at that depth are dominated by the narrow western-intensified recirculations of subtropical gyres. Narrow recirculations at depth are the major signatures of quasigeostrophic models forced by Ekman pumping at the surface. Rhines and Holland (1979) and Holland and Rhines (1980) show that the eddies generated by the baroclinic instability of the wind-driven circulation in turn drive the deep mean gyres. The similarities of the observations with such results suggest similar causes that have been examined previously by Hogg (1983) and Lozier (1997). With the full velocity field now available in such regions, further comparisons can be expected to become more exhaustive.

#### 3) The tropical and subtropical oceans

We speculate here that the weakness of the mean circulation in all oceanic basins at tropical latitudes is consistent with the idea of stagnant shadow zones of Luyten et al. (1983). The early Rooth et al. (1978) argument is adapted to the present situation as follows: for the zonal flow to vanish at an eastern boundary, the geostrophic pressure must be constant along that boundary. In the oceanic interior at that depth and range of latitudes (less than 20°–30°), the potential vorticity (PV) *f*/*h* (with *h* the thickness between two density surfaces) is very much governed by *f* as made clear by Keffer (1985) for his layer *D*. Now, because the geostrophic pressure of a conservative flow is constant on a PV contour, it must be constant in the whole region swept by PV contours connected to the eastern boundary. Therefore, the absence of circulation in these regions is consistent with a hypothesis of conservative geostrophic flow properties. Tracer observations in the North Atlantic by Jenkins et al. (1988) and Kawaze and Sarmiento (1986) support the idea that the intermediate waters are only ventilated to the north of 30°N. Notice the absence of data around 15° of latitude on the eastern boundary of the Atlantic and Pacific Oceans, which indicates that very few floats drifted there, further supporting the idea of shadow zones.

Note, however, that our observed circulation does not completely vanish at low latitudes as zonal flows are found (if the geopotential is differentiated) with extended zonal coherence in the Atlantic and Pacific. These year-mean alternating jets are also seen in high- resolution GCMs (Richards et al. 2006). A physical explanation was proposed for the strongest equatorial jets through the destabilization of mixed Rossby gravity waves (Hua et al. 2008), but is still lacking for the extra equatorial ones (beyond a few degrees from the equator).

#### 4) Indian Ocean

Both Davis (2005) and Reid (2003) agree on the existence of weak cyclonic gyres on either side of India. We have very little motion in the North Indian Ocean with geopotential signals less than 0.5 cm but data are lacking there. An anticyclonic tongue just to the west of northern Australia centered along 15°S is almost identical to that found in Reid (2003). South of 20°S, all three studies agree on a northwestward flow that is the northern limb of the subtropical gyre. On impinging on Madagascar, the flow does not flow north as in Reid (2003) but passes south to join the Agulhas which closes the gyre. The pressure maximum of the subtropical gyre is two times higher than for the corresponding South Atlantic or the South Pacific gyres in agreement with Davis (2005). South of this gyre lies the strong branch of the ACC coming from the Agulhas retroflection.

#### 5) South Pacific

The main gyre encompasses the whole Pacific Ocean, meeting Australia around 25°S in Fig. 16, and is very similar to the circulation in Davis (2005) from ALACE floats. Both datasets find a thin southward East Australian Current that connects later with the Indian Ocean. The Reid (1997) gyre does not extend as far to the east. North of 20°S is the domain of the weak zonal jets in our study (see Fig. 12).

#### 6) South Atlantic

The main gyre starts near 25°S, a more southern position than in the South Pacific and reaches its maximum at 35°S, (15 cm pressure difference) consistent with Reid (1989). North of that gyre (from 25°S to 5°S), Reid 1989 shows a large-scale cyclonic gyre, which is not observed in Fig. 16. Instead, a blow up of the region (Fig. 17) shows a weak (3 cm) anticyclonic circulation between 10° and 20°S, separated from the subtropical gyre to the south by the eastward flow at 20°S shown in Fig. 11. The East Falkland Current and its retroflection in the confluence region is well captured in both studies but Reid (1989) does not show the strong Zapiola gyre presumably because of its weak density signature (the float observed temperature variance is less than 0.1°C^{2}). Both the weak eastward flow at 20°S and the strong Zapiola gyre were first mapped as stationary features by subsurface floats (Boebel et al. 1999). They are firmly confirmed by this analysis.

#### 7) Antarctic Circumpolar Current

The geostrophic component of the ACC found from ALACE floats in Gille (2003) and Davis (2005) is very similar to what is found here with the difference that our method offers better resolution of the smaller scales. The present observations emphasize the clear southeastward path of the ACC in the Indian and Pacific Oceans interrupted by two violent northward excursions downstream of the Campbell Plateau and Drake Passage. A poleward component for the depth-integrated geostrophic transport was speculated long ago by Stommel (1957) to be linked to the upwelling required by the divergence of Ekman transport to the south of the latitude of maximum westerlies. The problem was addressed by Gill (1968) and it is tempting to compare the present 1000-m flow with these idealized model solutions but many questions remain to validate such a scheme.

When the ACC flows over ridges (Mid-Atlantic Ridge, Macquarie Ridge, and East Pacific Rise), equatorward excursions are systematically observed. However, as isobars cross readily the PV (=*f*/*H*_{ocean}) contours, these excursions fall short of accurately conserving PV. Gille (2003) has shown that her 900-m flow is consistent with an ACC reaching the bottom and she explains the smaller equatorward excursions by the equivalent barotropic structure of the ACC with a 700-m *e*-folding depth that weakens the PV constraint. We have noted previously that the ACC remains sensitive to the bottom topography aiming at avoiding obstacles or finding gaps. It is after all quite natural for a geostrophic flow to minimize bottom vertical velocities, hence divergence.

## 6. Conclusions

We have described a direct method to estimate the mean absolute geostrophic flow on a reference pressure surface, free from statistical assumptions on the scales of the flow. It needs however an original mapping of the currents, here by subsurface floats. Of course Eulerian measurements by current meters would do as well but the shorter Lagrangian time scales give a definite advantage to the former. Strictly speaking, the present solution is not valid within ±3° of the equator, even though the geopotential is smoothly continuous at the equator. A further study is necessary to take into account the nonlinear terms neglected in this “geostrophic” formulation.

Comparisons of our mean circulation with other float experiments are generally good but often differ depending on the data coverage. Comparisons with Reid’s hydrographic studies are good away from a band of about 20° on either side of the equator. This demonstrates the insightful analysis of hydrographic data that Reid has carried out over the years to solve this elliptic problem.

The global use of the Argo float data has allowed us to discuss the vertical structure of the oceanic circulation. The mean circulation penetrates to 1000 m predominantly in subpolar latitudes showing boundary-trapped mean flows circulating in the Kelvin wave sense. The signature of subtropical gyres is very much confined to western recirculations (at least in the North Atlantic and North Pacific). There is in our view qualitative support for the theoretical ideas of (i) deep eddy-driven mean recirculations and (ii) middepth energetic regions restricted to western parts of oceanic basins. The depth dependence of the mesoscale variability found here will be useful to better understand the links between the observed geostrophic turbulence and its generation mechanism, that is, baroclinic instability. The Argo data will continue to provide testing grounds for theories of the ubiquitous zonal jets found in the 20°S–20°N band in the Atlantic and Pacific.

There is thus a need of better resolution of the tropical and equatorial deep jets and of the boundary currents by Argo floats and it is hoped that the inclusion of the remaining 40% of the Argo data (as of 1 January 2013) in ANDRO (covering the period January 2010–December 2012) will be possible soon.

## Acknowledgments

The float data were collected and made freely available by the International Argo Program and the national programs that contribute to it (http://argo.jcommops.org). The Argo program is part of the Global Ocean Observing System. The present global mapping of the mean currents at 1000 dbar was made possible only because of the quality and world coverage of the then available Argo data and we acknowledge warmly the time and effort of all participants in the International Argo Program. We thank further Richard Schopp, Guillaume Roullet, and Olivier Marchal for useful comments, and two anonymous reviewers for their very careful reading and clarifying suggestions.

### APPENDIX A

#### Estimation of the Average Current Errors

Shear errors on individual velocity components are provided in the ANDRO dataset as

Here, PD is the parking depth (m), *W* is the vertical ascent rate (m s^{−1}), and Δ*t* is the parking time (s) for cycle number *n*. The descent rate is assumed equal to the ascent rate (which is acceptable because floats are subject to strong shears mainly in the top 400 m while their descents may stop temporarily before reaching 1000 m but generally near or deeper than 400-m depth). The and are the surface velocities estimated (linear least squares fits) over the last and first 6 h while the float is at the surface during the previous and current cycles, respectively. They are also given in ANDRO. Further details may be found in Ollitrault and Rannou (2013).

The “waiting time” errors at the surface (for floats using the Argos satellite system) are estimated as because statistics done on 750 representative floats give an average of 1500 ± 1500 m surface displacement. Actually, we assume a similar error both before sinking and after surfacing.

If now we have *N* individual deep velocity estimates within a specified area, the instrumental error on the average velocity is given as .

Note that if *ε*_{shear} is not available in ANDRO, it is defaulted to 1 cm s^{−1}. Obviously, the above expressions apply to both zonal and meridional components.

Given the velocity variances (on *u* and *υ*) within a specified area, the mesoscale errors are estimated as and , respectively. Here *T*_{Lu} and *T*_{Lv} are the Lagrangian integral time scales (Tennekes and Lumley 1972) that are taken constant at 10 and 5 days, respectively (based on past determination and provisory results not given in this paper).

### APPENDIX B

#### Implementation of the Method

To reduce the data gaps (see Fig. 2), an increase of the grid size is one option but this would be at the cost of spatial resolution near the boundaries where the energetic mean flows are to be found. Furthermore, we want to work on a regular 1° grid in latitude and longitude, while the positions of the mean velocity vectors are not given on such a grid (due to our center of mass averaging).

We use first optimal estimation to create (on a regular 1° grid) an interpolated velocity field, by using a locally nondivergent field of the vector = (*fu*, *fυ* cos*θ*). From observations of the mean field at positions **x**_{k} in the two dimensional horizontal plane, we want to find a nondivergent field **û** as close as possible in the least squares sense to the . By definition, the estimated **û** is linked to an estimated streamfunction by

In the vicinity of the arbitrary point **x**_{A}, is decomposed over a polynomial basis:

Limiting ourselves to polynomials of degree two and writing *δλ* = *λ* − *λ*_{A} and *δθ* = *θ* − *θ*_{A}, the above becomes

To obtain the five coefficients *C*_{mn}, we use a moving least squares algorithm to minimize the cost function:

where *k* is the index of data points surrounding the point **x**_{A} at which the estimator is required. In the moving least squares procedure, the weights *W*_{k} are used to enforce a local fit. Here, we choose

where and are the Euclidian zonal and meridional distances. The scales *L*_{λ} and *L*_{θ} of the weights allow a trade-off between resolution and statistical robustness. We have tested several choices and found appropriate values close to the scale of the mean velocity grid (*L*_{λ} = 139 km and *L*_{θ} = 111 km) to preserve spatial resolution. To minimize computation, only those data points that satisfy |*δλ*|/*L*_{λ} ≤ 1.75 and |*δθ*|/*L*_{θ} ≤ 1.75 are considered in the cost function [Eq. (B3)] resulting in a maximum of 32 neighboring data points. To obtain the least squares solution, write = 0 for the five *m*, *n* pairs, yielding a system of five equations with five unknowns (the *C*_{mn} coefficients), which is solved at each point of the reconstruction grid. If the number of data points is too small, the matrix is singular. By requiring at least five data points in the neighborhood of a point to be estimated, the singularities were limited to a few tens over more than 30 000 points. In such cases no estimation is done, thus resulting in a data hole. All grid points with depths shallower than 950 m are excluded from the domain *D*. A smoothed version of 2-Minute Gridded Global Relief Data (ETOPO2v2) was used to compute the bottom depth (http://www.ngdc.noaa.gov/mgg/global/etopo2.html).

Once the *C*_{mn} coefficients are known, the optimal estimate **û** is obtained from Eqs. (B1) and (B2) as

This formulation insures the nondivergence of the **û** field in the neighborhood of the point **x**_{A}.

The staggered *D* grid (Mesinger and Arakawa 1976) is the natural choice to implement geostrophy:

Note that this grid is isotropic in the (*λ*, *θ*) space that we choose. Vorticity is estimated at *ϕ* points and we use Eq. (B4) to estimate the vector over the points of the *D* grid that map the domain.

To solve Eq. (4) with boundary conditions [Eq. (5)], the Laplacian and vorticity are discretized using space centered finite differences. With no concern for computer time, the solution of Eq. (4) is found as the asymptotic limit of the following inhomogeneous diffusion equation:

Simple Euler forward differencing in time insures a stable numerical solution (Roache 1985). The context of the diffusion equation allows us to understand better the need for the consistency relation [Eq. (7)]: a steady solution of the diffusion equation can be found only if the net heat source (here net vorticity) over the domain can be expelled through the boundaries by adequate outward heat fluxes (here tangential velocities). The evolution toward steady state of Eq. (B5) was monitored by the ratio size of the term on the left (i.e., the residual) over the last term on the right (the forcing). The time for convergence varies as the square of the spatial scale in a diffusion problem. If spin up is achieved rapidly at the mesoscale, this is not so at planetary scale. We have found that the settling of the geopotential difference between the Pacific and Atlantic requires a considerable time of integration on the order of 10^{6} time steps Δ*t* [with Δ*t* = 0.1(Δ*x*)^{2} as required for stability]. The solution shown in this paper ensures that the above ratio is less than 10^{−4} at any point in the domain (actually it is less than 10^{−7} for 97% of the points). To find out the errors on Φ that originate from the errors on the mean velocity data , a Monte Carlo method was used. We add to each component of the vector at point **x**_{k}, a perturbation equal to the mesoscale error on that component times a random number between −1 and +1 drawn from a uniform probability distribution. The whole calculation is then restarted to generate a “perturbed” geopotential Φ_{p} and this is repeated 25 times (convergence is obtained after 10 trials). It is the rms value of Φ_{p} − Φ at each grid point that gives us the statistical error on Φ that we are looking for. Figure B1 gives the rms error (cm) as a function of position. The errors are usually less than ±1 cm but can increase up to ±2 cm in eddy-rich regions (the largest errors are found in the ACC to the southwest of Africa).

## REFERENCES

*Deep Sea Res.,*

**50,**

*Ocean Circulation and Climate,*G. Siedler, J. Church, and W. J. Gould, Eds., International Geophysics Series, Vol. 77, Academic Press, 123–139.

*J. Phys. Oceanogr.,*

*Philos. Trans. Roy. Soc. London,*

**A**

**12,**527–548.

**85,**109–126,

*GARP Publ. Ser.,*

**17,**1–64.

*International WOCE Newsletter,*Vol. 34, WOCE International Project Office, Southampton, United Kingdom, 7–10.

*Evolution of Physical Oceanography: Scientific Surveys in Honor of H. Stommel,*B. A. Warren and C. Wunsch, Eds., MIT Press, 70–110.

*Prog. Oceanogr.,*

**23,**149–244.

*Prog. Oceanogr.,*

**33,**1–92.

*Prog. Oceanogr.,*

**39,**263–352.

*Prog. Oceanogr.,*

**56,**137–186.

*Prog. Oceanogr.,*

**22,**171–204.

*A First Course in Turbulence*. MIT Press, 300 pp.

*Oceanographic Studies*. The John Hopkins University, 110 pp.

*Geophys. Monogr.,*No. 173, Amer. Geophys. Union, 53–74.

*Prog. Oceanogr.,*

**11,**1–59.