Abstract. The numerical implementation of an ocean model based on the incompressible Navier Stokes equations which is designed for studies of the ocean circulation on horizontal scales less than the depth of the ocean right up to global scale is described. A "pressure correction" method is used which is solved as a Poisson equation for the pressure field with Neumann boundary conditions in a geometry as complicated as that of the ocean basins. A major objective of the study is to make this inversion, and hence nonhydrostatic ocean modeling, efficient on parallel computers. The pressure field is separated into surface, hydrostatic, and nonhydrostatic components. First, as in hydrostatic models, a two-dimensional problem is inverted for the surface pressure which is then made use of in the three-dimensional inversion for the nonhydrostatic pressure. Preconditioned conjugate-gradient iteration is used to invert symmetric elliptic operators in both two and three dimensions. Physically motivated preconditioners are designed which are efficient at reducing computation and minimizing communication between processors. Our method exploits the fact that as the horizontal scale of the motion becomes very much larger than the vertical scale, the motion becomes more and more hydrostatic and the threedimensional Poisson operator becomes increasingly anisotropic and dominated by the vertical axis. Accordingly, a preconditioner is used which, in the hydrostatic limit, is an exact integral of the Poisson operator and so leads to a single algorithm that seamlessly moves from nonhydrostatic to hydrostatic limits. Thus in the hydrostatic limit the model is "fast," competitive with the fastest ocean climate models in use today based on the hydrostatic primitive equations. But as the resolution is increased, the model dynamics asymptote smoothly to the Navier Stokes equations and so can be used to address smallscale processes. A "finite-volume" approach is employed to discretize the model in space in which property fluxes are defined normal to faces that delineate the volumes. The method makes possible a novel treatment of the boundary in which cells abutting the bottom or coast may take on irregular shapes and be "shaved" to fit the boundary. The algorithm can conveniently exploit massively parallel computers and suggests a domain decomposition which allocates vertical columns of ocean to each processing unit. The resulting model, which can handle arbitrarily complex geometry, is efficient and scalable and has been mapped on to massively parallel multiprocessors such as the Connection Machine (CM5) using data-parallel FORTRAN and the Massachusetts Institute of Technology data-flow machine MONSOON using the implicitly parallel language Id. IntroductionDetails of the numerical implementation of a model which has been designed for the study of dynamical processes in the ocean from the convective, through the geostrophic eddy, up to global scale are set out. The "kernel" algorithm solves the incompressible Navier Stokes equations on the...

Abstract. Ocean models based on consistent hydrostatic, quasi-hydrostatic, and nonhydrostatic equation sets are formulated and discussed. The quasi-hydrostatic and nonhydrostatic sets are more accurate than the widely used hydrostatic primitive equations. Quasi-hydrostatic models relax the precise balance between gravity and pressure gradient forces by including in a consistent manner cosine-of-latitude Coriolis terms which are neglected in primitive equation models. Nonhydrostatic models employ the full incompressible Navier Stokes equations; they are required in the study of smallscale phenomena in the ocean which are not in hydrostatic balance. We outline a solution strategy for the Navier Stokes model on the sphere that performs efficiently across the whole range of scales in the ocean, from the convective scale to the global scale, and so leads to a model of great versatility. In the hydrostatic limit the Navier Stokes model involves no more computational effort than those models which assume strict hydrostatic balance on all scales. The strategy is illustrated in simulations of laboratory experiments in rotating convection on scales of a few centimeters, simulations of convective and baroclinic instability of the mixed layer on the 1-to 10-km scale, and simulations of the global circulation of the ocean.

Height coordinate ocean models commonly represent topography as a ''staircase'' of discontinuous steps that are fitted to the model grid. Here the ramifications of an alternative approach are studied in which ''shaved cells'' are used to represent irregular topography. The problem is formulated using the finite-volume method and care is taken to ensure that the discrete forms have appropriate conservation properties. Two representations of topography, ''partial step'' and ''piecewise linear,'' are considered and compared with the staircase approach in some standard problems such as the topographic ␤ effect and flow over a Gaussian bump. It is shown that shaved cells are clearly more accurate than the conventional staircase representation. The use of partial steps, although not as accurate as the piecewise linear approach, is seen to be superior to the staircase approach. Moreover, partial steps can be readily implemented in existing height coordinate models.

This paper describes the sea ice component of the Massachusetts Institute of Technology general circulation model (MITgcm); it presents example Arctic and Antarctic results from a realistic, eddy-admitting, global ocean and sea ice configuration; and it compares B-grid and C-grid dynamic solvers and other numerical details of the parameterized dynamics and thermodynamics in a regional Arctic configuration. Ice mechanics follow a viscous-plastic rheology and the ice momentum equations are solved numerically using either line-successive-over-relaxation (LSOR) or elastic-viscous-plastic (EVP) dynamic models. Ice thermodynamics are represented using either a zero-heat-capacity formulation or a two-layer formulation that conserves enthalpy. The model includes prognostic variables for snow thickness and for sea ice salinity. The above sea ice model components were borrowed from currentgeneration climate models but they were reformulated on an Arakawa C grid in order to match the MITgcm oceanic grid and they were modified in many ways to permit efficient and accurate automatic differentiation. Both stress tensor divergence and advective terms are discretized with the finite-volume method. The choice of the dynamic solver has a considerable effect on the solution; this effect can be larger than, for example, the choice of lateral boundary conditions, of ice rheology, and of ice-ocean stress coupling. The solutions obtained with different dynamic solvers typically differ by a few cm s −1 in ice drift speeds, 50 cm in ice thickness, and order 200 km 3 yr −1 in freshwater (ice and snow) export out of the Arctic.

Abstract. We first describe the principles and practical considerations behind the computer generation of the adjoint to the Massachusetts Institute of Technology ocean general circulation model (GCM) using R. Giering's software tool Tangent-Linear and Adjoint Model Compiler (TAMC). The TAMC's recipe for (FORTRAN-) line-by-line generation of adjoint code is explained by interpreting an adjoint model strictly as the operator that gives the sensitivity of the output of a model to its input. Then, the sensitivity of 1993 annual mean heat transport across 29øN in the Atlantic, to the hydrography on January 1, 1993, is calculated from a global solution of the GCM. The "kinematic sensitivity" to initial temperature variations is isolated, showing how the latter would influence heat transport if they did not affect the density and hence the flow. Over 1 year the heat transport at 29øN is influenced kinematically from regions up to 20 ø upstream in the western boundary current and up to 5 ø upstream in the interior. In contrast, the dynamical influences of initial temperature (and salinity) perturbations spread from as far as the rim of the Labrador Sea to the 29øN section along the western boundary. The sensitivities calculated with the adjoint compare excellently to those from a perturbation calculation with the dynamical model. Perturbations in initial interior salinity influence meridional overturning and heat transport when they have propagated to the western boundary and can thus influence the integrated east-west density difference. Our results support the notion that boundary monitoring of meridional mass and heat transports is feasible. IntroductionThe impending need to synthesize, basin-wide and globally, ocean data such as the entire World Ocean Circulation Experiment (WOCE) data set including altimetry and the surfaceforcing data obtained from weather centers and scatterometers makes imperative the use of sophisticated ocean general circulation models (GCMs) to (1) interpolate in space and time between the observations and (2) diagnose unobservable but important quantities such as vorticity and heat transports. Conversely (and, indeed, prior to all these interpretations), the data stream must be used to test and improve the GCMs (the stated goal 1 of International WOCE). A very powerful and general approach to synthesis is the use of optimization methods; we will concentrate on the particular flavor that has become known as the "adjoint approach" in meteorology and oceanography [e.g., Talagrand