From the beginning of global atmospheric modeling, solution of the governing equations on the sphere has been a challenging problem. Today, most global atmospheric models use either the "spectral" method, based on expansions in terms of triangularly truncated spherical harmonics (Jarraud and Simmons, 1983) or finite differences implemented on regular

General Circulation Model Development

Copyright © 2000 by Academic Press. All rights of reproduction in any form reserved.

grids in spherical (longitude-latitude) coordinates. Each of these methods has strengths and weaknesses. Finite-difference models based on latitude-longitude grids must deal with the convergence of the meridians at the poles, which demands a short time step for computational stability; this is the well known "pole problem." Filtering techniques (Arakawa and Lamb, 1977) are typically employed to allow a moderately long time step. In addition, fine zonal spatial resolution near the poles entails additional computations (e.g., for parameterized physical processes) that contribute little to the overall accuracy of the solution. On the other hand, increased zonal resolution in middle and high latitudes may be useful, particularly when the radius of deformation decreases near the poles.

The spectral method based on spherical harmonics with triangular truncation (e.g., Jarraud and Simmons, 1983) elegantly eliminates the pole problem. On the other hand, it gives poor results for the advection of strongly varying non-negative scalars such as the mixing ratios of water vapor and cloud water (e.g., Williamson and Rasch, 1994). The problem of negative water is sufficiently serious that most spectral models have now been modified to use semi-Lagrangian techniques for advection, rendering the modified models distinctly less spectral in character. Spectral computation of mass convergence and the horizontal pressure gradient force allows fast and easy implementation of semi-implicit time-differencing schemes, which permit a relatively long time step; in contrast, semi-implicit schemes are both more expensive and more complicated in finite-difference models, although they have the potential to eliminate the need for polar filtering.

Several efforts have been made to devise alternative discretizations of the sphere, to permit the use of finite-difference methods without the need for polar filtering. These include grids that skip some longitudinal grid points near the poles (Kurihara, 1965; Halem and Russell, 1973), matched and/or patched polar stereographic grids (Phillips, 1957; Browning et al., 1989), and grids based on polyhedra (Sadourny et al., 1968; Sadourny and Morel, 1969; Williamson, 1968, 1970; Sadourny, 1972; Thacker, 1978; Augenbaum and Peskin, 1985; Baumgardner and Frederickson, 1985; Masuda and Ohnishi, 1986; Nickovic, 1994; Popovic et al., 1996; McGregor, 1996; Rancic et al., 1996; Purser and Rancic, 1998; Thuburn, 1997). As recognized in early work by Arakawa and colleagues, grids based on icosahedra offer an attractive framework for simulation of the global circulation of the atmosphere; their advantages include almost uniform and quasi-isotropic resolution over the sphere. Such grids are termed geodesic because they resemble the geodesic domes designed by Buckmin-ster Fuller. Geodesic grids caught the attention of Arakawa during the 1960s (Sadourny et al., 1969).

Figure 1 (a) An ieosahedron inscribed in a unit sphere, (b) Bisect each edge forming four new faces, (c) Project the new vertices onto the unit sphere. The process is repeated in panels (d)-(f).

A geodesic grid is constructed by starting with an ordinary ieosahedron inscribed inside a unit sphere, as shown in Fig. la. The ieosahedron has 12 vertices. As a first step in the construction of a spherical geodesic grid, each face of the ieosahedron is subdivided into four new faces by bisecting the edges. The result of this process is shown in Fig. lb. Next, the new vertices are "popped out" onto the unit sphere, creating the polyhedron shown in Fig. lc. This recursive process can be continued indefinitely, yielding arbitrarily fine meshes. When a finer geodesic grid is derived by subdividing a coarser one, the finer one "nests" within the coarser one. A geodesic grid of any resolution consists of exactly 12 pentagons, corresponding to the 12 vertices of the original icosahedron, and a (typically much) larger number of hexagons. Such grids are quasi-homogeneous in the sense that the area of the largest cell is only a few percent greater than the area of the smallest cell. For more details, see Heikes and Randall (1995a).

Hexagonal grids are also quasi-isotropic. As is well known, only three regular polygons tile the plane: equilateral triangles, squares, and hexagons. Figure 2 shows planar grids made up of each of these three possible polygonal elements. On the triangular grid and the square grid, some of the neighbors of a given cell lie directly across cell walls, while others lie across cell vertices. As a result, finite-difference operators constructed on these grids tend to use "wall neighbors" and "vertex neighbors" in different ways. For example, the simplest second-order finite-difference approximation to the gradient, on a square grid, uses only wall neighbors; vertex neighbors are ignored. Although it is certainly possible to construct finite-difference operators on square grids (and triangular grids) in which information from all neighboring cells is used (e.g., the Arakawa Jacobian, as discussed by Arakawa, 1966), the essential anisotropics of these grids remain, and are unavoidably manifested in the forms of the finite-difference operators. Hexagonal grids, in contrast, have the property that all neighbors of a given cell lie across cell walls; there are no vertex neighbors. As a result, finite-difference operators constructed on hexagonal grids treat all neighboring cells in the same way; in this sense, the operators are as symmetrical and isotropic as possible. The 12 pentagonal cells also have only wall neighbors; there are no vertex neighbors anywhere on the sphere.


Winninghoff (1968) and Arakawa and Lamb (1977) studied the simulation of geostrophic adjustment on a variety of logically rectangular staggered grids for the solution of the shallow water equations. The variables that they considered were the zonal and meridional components of the wind, and the mass. They discussed five distinct staggered grids, which they labeled A through E. Randall (1994) defined the Z Grid, which is an unstaggered grid for the vorticity, divergence, and mass. With a rectangular mesh, the Z Grid corresponds to the Arakawa-Lamb C Grid for the divergent component of the wind, and to the Arakawa-Lamb D Grid for

12 neighbors, 3 wall neighbors

8 neighbors, 4 wall neighbors

6 neighbors, 6 wall neighbors

Figure 2 Cells neighboring a given cell (shaded) on triangular, square, and hexagonal grids. A "wall neighbor" is a neighbor that lies directly across a cell wall.

the rotational component. The spatial arrangements of the variables on the B, C, and Z Grids are shown in Fig. 3.

The analyses of Winninghoff (1968) and Arakawa and Lamb (1977) showed that geostrophic adjustment is best simulated on the C Grid, provided that the grid spacing is smaller than the radius of deformation. Although a modeler will try to choose a model's grid spacing to be smaller than the radius of deformation, this is not always possible; in atmospheric models, high internal modes can have small radii of deformation, and j+1

h u

h u

,v h u

v h u




Figure 3 The spatial arrangements of the horizontal wind components (u and v) and the mass {h) on B Grid and C Grid. Also shown are the spatial arrangements of the vorticity, divergence, and mas on the Z Grid.

small radii of deformation are even more prevalent in ocean models because the ocean is weakly stratified. These problems are exacerbated near the poles, where the radius of deformation tends to be particularly small because the Coriolis parameter is large. When the grid spacing is larger than the radius of deformation, the C Grid gives poor results (Randall, 1994).

This is illustrated in Fig. 4, which shows the gravity-inertia wave frequency as a function of the horizontal wave numbers for the continuous shallow-water equations, and for the B, C, and Z Grids. The results are shown for a case in which the grid size is smaller than the radius of deformation, and for a case in which it is larger. The C Grid gives realistic results when the grid spacing is finer than the radius of deformation, but behaves very badly when the radius of deformation is smaller than the grid size. As shown by Randall (1994) and as illustrated in Fig. 4, the Z Grid behaves well regardless of the ratio of the grid size to the radius of deformation.

Figure 4 The nondimensional frequency (on the vertical axis), plotted as a function of the normalized horizontal wave numbers (the horizontal axes) for the continuous shallow-water equations, and for the B, C, and Z Grids, for the case in which the grid spacing is twice the radius of deformation (left column) and for the case in which the grid spacing is one-tenth of the radius of deformation (right column).

0 0

Post a comment