## Abstract

The presence of intense recirculations in the abyssal oceans has been revealed by both observations and modeling studies. A suggested mechanism is the interactions between geostrophic eddies and the mean circulation in the presence of variable bottom topography. Here such interactions are studied using an idealized numerical model, consisting of a single active abyssal layer overlying variable bottom topography. An initial ensemble of eddies quickly organize themselves to generate a mean anticyclonic circulation around seamounts. During this rearrangement energy is approximately conserved while potential enstrophy is dissipated, consistent with previous studies of geostrophic turbulence.

A parameterization of geostrophic eddies is formulated in terms of an eddy-induced transport, **U***. Based on results from the eddy-resolving experiments it is hypothesized that **U*** should dissipate potential enstrophy while conserving energy. Using variational methods, a solution for **U*** is found that dissipates potential enstrophy most efficiently, subject to energy being conserved. For a shallow water layer, it is shown that **U*** = *κ*[**∇**(*Q*^{2}/2) + *λ***∇***B*], where *Q* is the potential vorticity, *B* is the Bernoulli potential, *κ*(*x, y*) is an arbitrary function that describes the efficiency with which eddies can rearrange fluid within a layer, and *λ* is a Lagrange multiplier that is determined through the energetic constraint. Results from a coarse resolution version of the model using the parameterization compare favorably with those obtained using the eddy-resolving model.

## 1. Introduction

The traditional view of the abyssal circulation is of a weak, poleward meridional motion driven by upwelling (Stommel and Arons 1960). This view has recently been challenged by new observations and modeling studies. For example, data from the Brazil Basin Experiment reveal an intense basin-scale recirculation in the large-scale structure (Speer and Zenk 1993; de Madron and Weatherly 1994). The presence of mean circulations around isolated seamounts has also been observed by, for example, Brink (1995).

There is considerable evidence from a variety of modeling studies (e.g., Bretherton and Haidvogel 1976; Treguier 1989; Chapman and Haidvogel 1992; Beckmann and Haidvogel 1993, 1997; Haidvogel et al. 1993; Beckmann 1995) that mean anticyclonic circulations are generated around seamounts, often referred to as “cold domes” (Holloway 1997). One suggested mechanism for the generation of these intense recirculations is the interactions between geostrophic eddies and the mean circulation in the presence of variable bottom topography (Holloway 1987, 1992; Spall 1994). Our aim is to understand how the interactions of geostrophic eddies with the mean circulation in the presence of variable bottom topography can lead to the generation of cold domes, and how these effects may be parameterized.

Recent advances in the eddy parameterization problem have stemmed from the recognition that the ocean interior is adiabatic to leading order (Gent and McWilliams 1990). Incorporating adiabatic properties into eddy parameterizations has led to an impressive list of improvements in ocean general circulation models, including improved global temperature distribution, improved poleward and surface heat fluxes, and improved occurrence of deep convection (Danabasoglu et al. 1994). However, the Gent and McWilliams (1990) parameterization spins the ocean down toward a state of rest. This is at odds with the findings from the observations and modeling studies mentioned above.

Holloway (1987, 1992, 1997) argues that eddy parameterizations should instead drive the ocean toward finite, topography following currents. His argument is that the ocean is driven toward a state of maximum entropy—the “Neptune effect.” However, entropy is not easily computed from the variables carried in a standard ocean model. It may be possible to obtain similar results by incorporating additional constraints into the Gent and McWilliams parameterization. For example, Bretherton and Haidvogel (1976) and Sadournay and Basedevant (1985) have suggested that geostrophic eddies dissipate potential enstrophy while conserving total energy.

To address these issues we develop a simple, shallow-water model consisting of a single moving abyssal layer overlying variable bottom topography. While highly idealized, the model contains the essential ingredients to study how geostrophic eddies interact with the mean circulation in the presence of variable bottom topography. Using the results from the numerical model we diagnose those quantities that are conserved and those that are dissipated. These findings are then used to develop a simple parameterization of geostrophic eddies in our model.

The layout of the paper is as follows. In section 2 we review three eddy parameterizations and discuss their conservation properties. In section 3 we describe the numerical experiments. In section 4 we develop the parameterization and illustrate its use in a coarse resolution model. In section 5 we outline how the parameterization might be extended to two or more layers. Finally, in section 6 we discuss the wider implications of our work and suggest some future developments.

## 2. Eddy parameterizations

Here we review three eddy parameterizations—downgradient diffusion, the Gent and McWilliams, and a potential vorticity closure—that are currently in use in numerical ocean models. We reinterpret these parameterizations in terms of the quantities that they dissipate and the quantities that they conserve.

### a. Downgradient diffusion

Traditionally eddies have been parameterized in terms of a downgradient diffusion of tracers and momentum (e.g., Stone 1972). Thus, the eddy flux of density is parameterized

where **u** is the horizontal velocity, *ρ* is the density, and *κ*_{D} is an eddy transfer coefficient; overbars denote time mean, and primes denote time varying components. This closure spins the ocean down to a state of rest in which the isopycnals are flat, but does not conserve the volume of fluid between any two isopycnals. Consequently downgradient diffusion leads to spurious diapycnal mixing and erosion of the thermocline (Veronis 1975; Böning et al. 1995; Roberts and Marshall 1998).

### b. Gent and McWilliams

Gent and McWilliams (1990) proposed an alternative parameterization that is adiabatic; that is, the volume of fluid between any two isopycnals is exactly conserved. To achieve this, their scheme both mixes tracers along isopycnals and rearranges fluid within isopycnal layers by an eddy-induced transport,

where *h* is the thickness of an isopycnal layer (Gent et al. 1995).

To parameterize **U***, Gent and McWilliams assume that eddies act to flatten isopycnals, driving the ocean toward a state of rest, mimicking the effects of baroclinic instability. For example, consider a two-layer ocean overlying variable bottom topography, as shown in Fig. 1, where *h*_{1} is the thickness of the upper layer. Gent and McWilliams parameterize the eddy-induced transport **U*** in terms of the slope of the isopycnals. In two layers this gives

We now show that (4) represents the form of **U***, which dissipates potential energy most efficiently. Here the potential energy *P* is defined to within an arbitrary constant as

where *g*′ = *g*Δ*ρ*/*ρ*_{0} is the reduced gravity and *ρ*_{0} is a reference density.

The change in thickness of the upper layer due to **U*** in a time Δ*t* is

The associated change in potential energy is

Here we have assumed that the normal component of **U*** vanishes at solid boundaries. Additionally we impose a constraint on the magnitude of **U***,

This constrains the rate at which the eddies are able to rearrange fluid parcels within the layers. Without (8) the rearrangement would occur instantaneously.

Now varying **U***, we maximize the dissipation of potential energy subject to the constraint (8), giving

Here *γ* is a Lagrange multiplier determined by the constraint (8). Since the variations *δ***U*** are arbitrary, we obtain

where *κ*_{GM} = *ρ*_{0}*g*′Δ*t*/*γa* is an eddy transfer coefficient, which may vary spatially.

Therefore, in a two-layer model, the Gent and McWilliams parameterization represents the eddy-induced transport that dissipates potential energy most efficiently, regardless of the form of the bottom topography.

### c. Potential vorticity closure

In classical turbulence theory it is only materially conserved properties that are fluxed in a downgradient manner. The potential vorticity,

is materially conserved by fluid parcels, where *f* is the Coriolis parameter, *ζ* is the relative vorticity, and *h* is the isopycnal layer thickness. Green (1970), Marshall (1981), Rhines and Young (1982), and others have therefore argued that geostrophic eddies act to mix potential vorticity along isopycnals. Assuming that the Reynolds stresses are negligible on large scales (Treguier et al. 1997; Greatbatch 1998), we can relate the eddy-induced transport to the eddy flux of potential vorticity:

(Treguier et al. 1997; Killworth 1997). Eddy-resolving numerical calculations with zonal baroclinic jets have generally supported the potential vorticity, rather than the Gent and McWilliams closure (Lee et al. 1997; Treguier 1999; Marshall et al. 1999). The use of a potential vorticity closure to parameterize the formation of cold domes has been advocated by Greatbatch and Li (2000).

The Gent and McWilliams parameterization can be regarded as the closure that dissipates potential energy most efficiently; likewise, the potential vorticity closure can be viewed as the scheme that dissipates potential enstrophy most efficiently. Within an isopycnal layer, the change in potential enstrophy due to **U*** in a time Δ*t* is

where we have integrated by parts as in the previous example, and note Δ**u** = Δ*ζ* = 0.

Now varying **U***, we maximize the dissipation of potential enstrophy subject to the same constraint (8), giving

Since the variations *δ***U*** are arbitrary, we obtain

where *κ*_{Q} = −*Q*^{2}Δ*t*/(*γa**h*) is an eddy transfer coefficient.

Therefore the potential vorticity closure represents the eddy-induced transport that dissipates potential enstrophy, or equivalently mixes potential vorticity, most efficiently within isopycnal layers.

### d. Implications for the abyssal circulation

What are the implications of these closures for an abyssal layer? Consider the thought experiment, shown in Fig. 2, in which we have an initial baroclinic front overlying variable topography. The Gent and McWilliams closure (and also the downgradient diffusive closure) acts to slump the isopycnals and drive the ocean toward a state of rest, as sketched in Fig. 2a. In contrast the potential vorticity closure acts to homogenize the potential vorticity within isopycnal layers. Within an abyssal ocean this causes isopycnals to rise excessively over topographic features, as shown in Fig. 2b. Aside from generating unphysically large flows, this complete homogenization of potential vorticity requires a large input of energy. In reality, cold domes are observed in the deep ocean and also in numerical models, but the potential vorticity is not made uniform in the abyssal layers, as sketched in Fig. 2c. This is more consistent with a hypothesis proposed by Bretherton and Haidvogel (1976) and Sadournay and Basedevant (1985) in which geostrophic eddies dissipate potential enstrophy while conserving energy.

## 3. Numerical experiments

We now present some numerical results to illustrate how geostrophic eddies interact with the mean circulation in the presence of variable bottom topography. The following sections outline the experimental setup and the results.

### a. Model formulation

The numerical experiments are conducted using an “inverted” shallow water model (Fig. 3) on a *β* plane. The motion is confined to a single abyssal layer, overlying a variable bottom topography; the upper ocean is assumed to be at rest. While highly idealized, this model contains the most important ingredients necessary to address the interactions of geostrophic eddies with the mean circulation in the presence of variable bottom topography. At the same time, the model is simple enough that a wide variety of integrations can be performed and interpreted.

The model is configured in a square basin. Finite amplitude topography is totally contained within the depth of the lower layer, in the form of a single seamount, with a maximum at the center of the model domain. An initial ensemble of eddies is introduced by perturbing the height of the density interface between the two layers with an imposed geostrophic velocity condition. From this initial state the model is integrated forward in time.

The shallow-water equations are solved on a C grid:

where overbars indicate an average between adjacent grid points in the stated direction.

Here (*u, υ*) is the fluid velocity, *h* is the layer thickness,

is the Bernoulli potential,

is the relative vorticity,

is the height of the density interface separating the two layers, *h* is the layer thickness, and *H* is the height of the topography above a reference level (Fig. 3). Additionally, *f* = *f*_{0} + *βy* is the Coriolis parameter, *g*′ = *g*(Δ*ρ*)/*ρ*_{0} is the reduced gravity and *ρ*_{0} is a reference density. Here *A* is the coefficient of biharmonic friction and *r* is the coefficient of bottom friction.

The boundary conditions are no normal flow,

and free slip,

We also require a higher-order boundary condition for the biharmonic friction, which is

Here the subscripts denote components perpendicular and parallel to the boundary.

### b. Illustrative experiment

We now present an experiment in which we take a single topographic seamount defined as

where *h*_{max} = 500 m is the maximum height of the seamount. The lateral boundaries are located at *x* = 0, *L* and *y* = 0, *L.* The width of the basin is *L* = 1000 km. We have considered a range of more complex topographies, and generally obtain results similar to those reported here.

The mean height of the density interface separating the two layers is 750 m. The maximum initial perturbation of the interface height field is ±60 m. Other physical parameters are *f*_{0} = 0.7 × 10^{−4} s^{−1}, *β* = 2 × 10^{−11} m^{−1} s^{−1}, *ρ*_{0} = 10^{3} kg m^{−3}, *g*′ = 0.02 m s^{−2}. Our value of *g*′ is somewhat large for the abyssal ocean but has been chosen to increase the size of the Rossby deformation radius, such that we can better resolve the turbulent cascades. For the parameters given above the Rossby deformation radius is *O*(50 km).

In this integration, the horizontal grid spacing is 5 km, allowing a turbulent cascade to be explicitly resolved, and *A* = 2.5 × 10^{8} m^{4} s^{−1}, giving a spindown time of approximately 4 weeks at the grid scale. There is no explicit bottom friction; that is, *r* = 0. The experiments are integrated with a time step of 300s using a leapfrog scheme. A Robert–Asselin time filter is applied to prevent alternate time steps from diverging.

First, we consider the evolution of the interface height during the integration (Fig. 4). Interface height contours are streamlines for the geostrophic flow, and also give an indication of the spatial structure of the energy field. Figure 4a shows the initial field, containing a series of geostrophic eddies. The mean interface height is perturbed by ±*O*(60 m), with maximum geostrophic velocities of *O*(30 cm s^{−1}). Figure 4b, after 4 model months, shows the clear signature of the interface raising over the sea mount at the center of the domain. Examination of the data at times before and after this period show that this is a persistent feature, and not merely a transient eddy passing over the sea mount. Figure 4c, after 5 model years, shows considerable doming of the interface over the sea mount, and associated mean anticyclonic flow. The interface is raised over the sea mount by *O*(120 m) with geostrophic velocities approaching 40 cm s^{−1}. The amplitude of the doming of the interface is proportional to the energy available in the initial eddy field. In addition, there is a tendency for anticyclonic eddies to accumulate at the northern boundary and for cyclonic eddies to accumulate at the southern boundary. This is a standard property of closed-basin experiments (Wang and Vallis 1994). For comparative purposes results from an integration with a flat bottom are shown in Fig. 5. After 5 model years the solution takes the form of two “Fofonoff gyres,” free modes where the flow is along potential vorticity contours (Fofonoff 1954). It should be noted that, if the seamount experiments are performed on an *f* plane, then the solution is characterized by a doming of the interface over the topography and a series of large-scale eddies around the perimeter of the basin (not shown).

Returning to the experiment with the seamount (Fig. 4), even after 40 model years, Fig. 4d shows the doming of the interface and associated anticyclonic circulation remain. There is no suggestion that the circulation is spinning down toward a state of rest. Note that throughout the experiment, the interface height variations, and hence the energy field, remain trapped at large spatial scales.

Now consider the evolution of the potential vorticity field (Fig. 6). Initially, Fig. 6a shows how variations in the potential vorticity are dominated by the topography, with maximum values of potential vorticity found over the seamount. After 4 months, Fig. 6b shows that the eddy field has stirred the potential vorticity, creating a series of filaments. There is a tendency for filaments of high potential vorticity fluid to be transferred off the seamount, and conversely for filaments of low potential vorticity to be transferred onto the seamount. These filaments cascade to ever smaller spatial scales and are ultimately dissipated by the biharmonic momentum dissipation. After 5 model years, Fig. 6c shows that the filamentation is reduced over the seamount, although there is still mixing of potential vorticity in the surrounding fluid. This is reinforced in Fig. 6d after 40 model years. The turbulent cascade mixes the potential vorticity, which is consequently reduced over the seamount. As large-scale variations in the potential vorticity are dominated by variations in layer thickness, then decreasing the potential vorticity over the seamount increases the layer thickness and consequently raises the interface.

These changes in interface height and potential vorticity over the seamount are clearly shown in the scatterplots in Fig. 7. There is a tendency for the two quantities to collapse onto a simple functional form. The maximum value of potential vorticity is reduced, whereas the maximum value of interface height is increased.

### c. Energetics and potential enstrophy

We now consider the evolution of the energy and potential enstrophy for a series of integrations with differing dissipation coefficients.

Time series of the available energy integrated over the basin are shown in Fig. 8 for different values of biharmonic momentum dissipation. A value of zero corresponds to the minimum energy state in which the fluid is at rest and the interface is flat; this is the end state predicted by the Gent and McWilliams parameterization. Figure 8 shows that the amount of energy dissipated is correlated with the magnitude of the biharmonic momentum dissipation coefficient. In all cases the energy of the system never approaches the zero line, showing that the model does not reach a state of rest over the 5 model years. For the smallest value of dissipation, less than 5% of the available energy is dissipated. These results suggest that in the ocean, where the dissipation is truly weak, geostrophic eddies should *conserve,* rather than dissipate, energy.

Time series of potential enstrophy are shown in Fig. 9. Here the upper horizontal line indicates the potential enstrophy for the minimum energy state with a flat interface, the end state predicted by the Gent and McWilliams parameterization. The zero line indicates the potential enstrophy for a state of uniform potential vorticity, the end state predicted by the potential vorticity closure. The potential enstrophy is very quickly dissipated below the level predicted by a state of rest, but remains well above the level predicted for uniform potential vorticity. These results are remarkably insensitive to the level of biharmonic momentum dissipation.

These results contradict both the Gent and McWilliams and the potential vorticity closures, but support the hypothesis of Bretherton and Haidvogel (1976) that geostrophic eddies dissipate potential enstrophy, subject to the constraint of conserving total energy.

### d. Sensitivity to bottom friction

In the integrations discussed above no bottom friction is included. We now show, in Figs. 10 and 11, the impact of including bottom friction on our results. We consider two values, *r* = 1 × 10^{−7} s^{−1} and *r* = 1 × 10^{−8} s^{−1}; in each case the biharmonic momentum dissipation coefficient is *A* = 2.5 × 10^{8} m^{4} s^{−1}. With finite bottom friction there is considerable dissipation of energy throughout the integration. However, the amount of potential enstrophy dissipated over first model year remains relatively insensitive to the bottom friction. The eddies appear to interact with the mean circulation on much faster timescales than does the bottom friction. We therefore feel justified in neglecting bottom friction in our analysis.

## 4. A new parameterization

### a. Formulation

We now discuss how the numerical results described in the previous section can be exploited to develop a parameterization for the effects of geostrophic eddies on the mean circulation in the presence of variable bottom topography. As in Gent et al. (1995), we work with an eddy-induced transport **U***. We imagine that **U*** acts to rearrange fluid parcels within an isopycnal layer, leading to the doming of the interface over the seamount as sketched in Fig. 12. To ensure that mass is conserved within the system, a no-flux condition is imposed on **U*** at the boundaries (i.e., **U**^{*}_{⊥} = 0).

Our numerical results suggest that we should choose **U*** to ensure that potential enstrophy is dissipated, while conserving total energy. To achieve this, we use variational methods to find the **U*** that will *most efficiently* dissipate potential enstrophy subject to energy being conserved. While there is no reason why the dissipation should be the most efficient, this ensures that the parameterization drives the model toward a plausible end state. The change in interface height due to **U*** in one time step is given as

Note that **U*** has no direct effect on *u* and *υ*; that is,

The *dissipation* of potential enstrophy in the same time is therefore given by

The total energy is

The change in total energy in one time step due to **U*** is therefore

where *B* is the Bernoulli potential as defined in (19).

As in section 2 we also impose a constraint on **U***:

which limits the rate at which the eddies are able to rearrange the fluid within the layer. Here *a*(*x, y*) can be an arbitrary function. Note that if the constraint is removed, equivalent to setting *a* = 0, we obtain the solution that **U*** is infinite, corresponding to an instantaneous rearrangement toward the end state.

Using the calculus of variations to obtain the **U*** that gives the maximum decrease in potential enstrophy, we find

where *λ* and *γ* are Lagrange multipliers.

Since the variations *δ***U*** are arbitrary, we therefore obtain

Here *λ* is uniquely determined by the energy constraint, (30); *κ* = −Δ*t*/*γa* is analogous to an eddy transfer coefficient and is determined by the constraint on the magnitude of **U***, (31).

Thus the eddy-induced transport is related to *both* the slope of the gradient of potential vorticity and the gradient of Bernoulli potential. On large scales, the Bernoulli potential is dominated by the layer interface height; our **U*** therefore contains elements of *both* the Gent and McWilliams parameterization and the potential vorticity closure. Indeed, the scheme predicts that **U*** vanishes under two limiting cases.

In the limit that the fluid is at rest and the layer interface is flat, we have

**∇***B*= 0. The energetic constraint also ensures that*κ*= 0, and thus**U*** vanishes. Physically, this occurs because we are in a*minimum energy state.*There is no rearrangement of fluid parcels available without increasing the total energy.In the limit that the potential vorticity is uniform, we have

**∇***Q*= 0. The energetic constraint here ensures that*λ*= 0, and thus**U*** again vanishes. Physically, we are in a*minimum potential enstrophy state.*There is no rearrangement of fluid parcels available without increasing the potential enstrophy.

More generally, both *κ* and *λ* are finite, and we obtain partial, but not total, mixing of potential vorticity.

### b. Example 1

We now illustrate how our parameterization can generate plausible doming of the interface over a seamount in a coarse resolution simulation, where the turbulent cascade is not resolved. We present results from parallel integrations, at 5-km and 40-km resolution. In these experiments we initialize the integrations with a series of geostrophic eddies, as shown in Fig. 13. The coarse resolution version of the model is integrated from this initial state using both the resolved dynamics and the parameterization of the unresolved geostrophic eddies; that is

For comparison we also show results obtained without the parameterization at coarse resolution. The same biharmonic momentum dissipation coefficient, *A* = 5 × 10^{9} m^{4} s^{−1}, is used in each case; there is no bottom friction. In our parameterization we set *κ* = 1.5 × 10^{18} m^{5} s. This value for *κ* has been selected such that the timescale for the doming of the interface is comparable in the high resolution and coarse resolution integrations. By analogy with (12), we can convert this value into an effective eddy transfer coefficient for potential vorticity by multiplying by *f*^{2}/*h*^{3}; over the seamount this gives *O*(500 m^{2} s^{−1}), with somewhat smaller values off the seamount.

The instantaneous interface height and potential vorticity fields are shown in Fig. 14 for each of the three cases after 30 model weeks. In both the high resolution case and the coarse resolution case with the parameterization, there is a clear doming of the interface over the seamount (Figs. 14a,c, respectively). In contrast in the coarse resolution case without the parameterization, there is relatively little doming of the interface (Fig. 14e). This behavior is consistent with the potential vorticity fields. In the first two cases, the potential vorticity is reduced over the seamount (Figs. 14b,d) by the resolved and parameterized mixing respectively. In contrast, the potential vorticity retains a much larger value over the seamount in the coarse resolution case without the parameterization (Fig. 14f), in which there is relatively little mixing.

### c. Example 2

To understand how the high-resolution and coarse-resolution models evolve toward the end state it is helpful to consider a further, highly idealized example. The model in initialized with a single anticyclonic eddy sited away from the center of the seamount, as shown in Fig. 15. Here the biharmonic dissipation coefficients are *A* = 2.5 × 10^{8} m^{4} s^{−1} in the high resolution (5 km) case, and *A* = 5 × 10^{9} m^{4} s^{−1} in the coarse resolution (40 km) cases. We use the same value for *κ* as in the previous example above. In Fig. 16 we show the interface height after 5 model years in the three different cases.

In the high resolution case (Fig. 16a) we see that the interface is raised over the seamount (*η*_{max} = 784 m), giving an anticyclonic circulation, as in section 3, with evidence of Fofonoff gyres forming at the northern and southern boundaries. The doming of the interface occurs through the same turbulent cascade process as seen in the earlier experiments. Even if we initialize the model with a cyclonic eddy, we still obtain an anticyclonic circulation over the seamount (not shown).

In the coarse resolution case with the parameterization included (Fig. 16b) we obtain approximately 50% more doming over the seamount (*η*_{max} = 806 m) than in the high resolution case. In this highly idealized experiment the precise level and timescale of the doming is sensitive to the value of *κ.* Also note that the parameterization appears to transfer the energy predominantly into the anticyclonic eddy over the seamount, rather than into Fofonoff gyres.

Finally, in the coarse resolution case without the parameterization (Fig. 16c) we obtain very little doming over the seamount, although there is slight evidence of a Fofonoff gyre forming at the southern boundary.

While there are qualitative similarities between the end states obtained at high resolution, and at coarse resolution using the parameterization, there are significant differences between the routes that the fluid takes onto the seamount. In the high resolution case the eddy initially circulates around the boundary of the domain before breaking up and slowly being mixed onto the bump from all directions. In the coarse resolution case the eddy takes the shortest route onto the sea mount.

These differences are reinforced by looking at the evolution of the total energy over the two integrations (Fig. 17). Due to the decreased resolution, there is slightly less initial energy in the coarse resolution model. There are differences in the timescale over which the energy is dissipated. This is inevitable given the larger levels of dissipation required to keep the coarse resolution model stable. However, significantly more energy is retained in the coarse resolution model with the parameterization than without the parameterization.

Nevertheless, while the details of how the interface domes over the seamount differ between high and coarse resolution integrations, we are encouraged that we obtain qualitatively similar end states by observing the appropriate integral properties.

## 5. Extending the parameterization to two layers

We now briefly outline how our parameterization might be extended to two, or more, layers. For simplicity, the following analysis is developed for the two-layer model sketched in Fig. 1. The upper-layer thickness is *h*_{1} and the lower-layer thickness is *h*_{2}. Assuming a rigid lid, the eddy-induced transport **U*** must be equal and opposite in each layer, as indicated in the figure.

The change in *h*_{1} and *h*_{2} in a time interval Δ*t* due to **U*** is

We now define a measure of the net potential enstrophy as

The constants *α*_{1} and *α*_{2} allow the potential enstrophy to be weighted differently in each layer. The change in Λ due to **U*** is

The total energy, to within an arbitrary constant, is

and the change in *E* due to **U*** is

As previously we also impose the constraint on **U***,

to limit the rate at which the eddies can rearrange the fluid within the isopycnal layers.

Now varying **U*** to maximize the dissipation of potential enstrophy, −ΔΛ, we find

Again *λ* is uniquely determined by the energy constraint; *κ* = −Δ*t*/*γa* is analogous to an eddy transfer coefficient and is determined through the constraint on **U***.

With two active layers, the eddy-induced transport is therefore dependent on the potential vorticity in *both* layers and the difference in Bernoulli potential between the two layers. On large scales, the latter is proportional to the height of the interface separating the two layers, as in the Gent and McWilliams parameterization.

The analysis can be extended to multiple layers through the same procedure, except that it is necessary to impose an additional constraint that the depth-integrated eddy-induced transport vanishes.

While this parameterization in multiple layers ensures that the *net* potential enstrophy is dissipated, it is possible that the potential enstrophy may increase within any individual layer. This may appear unphysical because it implies unmixing of potential vorticity. However, such a two-layer calculation has been performed by Merryfield and Holloway (1999). They find that, while the eddy fluxes of potential vorticity are downgradient in the lowest layer, the potential vorticity fluxes can be either downgradient or upgradient in the upper layer, dependent on the form of dissipation used in the model. This issue clearly needs further study.

There are also severe problems with (42) when a layer outcrops since the potential enstrophy is then infinite. We have side stepped this issue in our shallow-water model by ensuring that the topography is contained entirely within the lower layer. However, the mixing of “delta sheets” of potential vorticity between the boundaries and the ocean interior is a real physical process (Hallberg and Rhines 1997) and needs to be incorporated into a complete parameterization. It may be possible to work instead with the inverse of potential enstrophy, which is well behaved at outcrops.

## 6. Summary and discussion

We have studied the interactions between geostrophic eddies and the mean circulation in the presence of variable topography using an idealized shallow-water model. Consistent with previous studies we have found that the eddies lead to doming of the layer interface over seamounts and an associated anticyclonic circulation. During this process we find that the total potential enstrophy is dissipated whereas energy is conserved, in line with geostrophic turbulence theory.

Using these results we have developed a parameterization for the geostrophic eddies in terms of an eddy-induced transport. Using variational methods, we have solved for the eddy-induced transport that dissipates potential enstrophy most efficiently, subject to energy being conserved. We find that this eddy-induced transport depends on *both* gradients in the potential vorticity and the Bernoulli potential. Our parameterization therefore contains elements of both the Gent and McWilliams and potential vorticity closures. Consistent with the numerical results, our parameterization drives the ocean toward a state in which potential vorticity is mixed but not made uniform. We do not obtain a state of rest. This contrasts with both the Gent and McWilliams scheme that drives the ocean toward a state of rest and the potential vorticity closure that drives the ocean toward a state of uniform potential vorticity. The latter has recently been advocated by Greatbatch and Li (2000) and Dewar (1998). However, we believe that the potential vorticity closure is unphysical over topography as it requires a large input of energy, violating one of the basic principles of physics.^{1}

Our approach is complementary to the statistical mechanics approach of Holloway (1992). Indeed Brands (1998) has recently shown that there are close similarities between increasing entropy and the “selective decay” of potential enstrophy, subject to an energy constraint. However there are also important differences; for example, Holloway’s “Neptune” parameterization will drive an ocean initially at rest towards a moving end state, whereas our parameterization leaves a resting ocean at rest. We also believe that our parameterization may be easier to implement in practice as it is based on variables that are already carried explicitly in numerical ocean models.

Our parameterization is also closely related to the “anticipated potential vorticity method” developed by Sadournay and Basedevant (1985). This method works by applying an additional force perpendicular to the velocity in the momentum equations. This force is chosen to ensure that potential enstrophy is dissipated; energy is automatically conserved *locally* because this force is perpendicular to the velocity. In contrast in our parameterization energy is conserved *globally.* The relative merits of the two approaches requires further study.

It is interesting to speculate how our parameterization behaves in the slumping of a baroclinic front. The classical paradigm of baroclinic instability is spindown to a state of rest. However, in our parameterization the potential energy released through the slumping of the baroclinic front must be converted to kinetic energy. We believe this is physical, as energy remains trapped at large scales in geostrophic turbulence. Numerical studies (e.g., Rhines 1975) suggest that geostrophic turbulence on a *β* plane generates multiple, zonal jets. We are currently exploring the possibility that our parameterization may generate such jets.

In this paper we have considered freely decaying turbulence in which forcing and dissipation are neglected. However, it should be straightforward to apply our parameterization to a forced scenario, provided that there is no interaction between the eddies and the forcing. We believe this to be a reasonable, zero-order assumption over large-scale topographic variations. However, over small-scale topographic variations, the topography can transfer energy to smaller horizontal scales, in which case energy is likely to be dissipated (Treguier and Hua 1988).

Our ultimate goal is to develop a parameterization for use in ocean general circulation models. In extending our parameterization to multiple layers it will be necessary to incorporate outcropping. There is also the issue of how the eddy transfer coefficient should vary with the background mean flow. Nevertheless, in this paper we have developed a framework for incorporating ideas from geostrophic turbulence theory into eddy parameterization schemes. We have also provided an important first step in developing a parameterization for a shallow water layer.

## Acknowledgments

We wish to thank Alberto Alvarez, Bob Hallberg, Greg Holloway, Brian Hoskins, and Alan O’Neill for useful discussions. We are also grateful for the useful comments from two anonymous reviewers. Funding was provided by the UK Natural Environment Research Council, GT4/96/256/HAS, and GR9/3760.

## REFERENCES

## Footnotes

*Corresponding author address:* Dr. Susan T. Adcock, Department of Meteorology, University of Reading, P.O. Box 243, Reading RG6 6BB, United Kingdom.

Email: susan@met.rdg.ac.uk

^{1}

Killworth (1997) and others have suggested that one can make the eddy transfer coefficient a function of the background flow such that the eddy fluxes vanish when the flow becomes stable. While this may prevent the potential vorticity from being homogenized, this approach does not ensure that energy is conserved.