The motions of the atmosphere are governed by physical laws which are expressed in the form of the equations of hydrodynamics and thermodynamics. These equations are, in principle, well known, but some of the physical processes in the atmosphere, such as small-scale turbulence and convection, cloud formation, and precipitation mechanisms, are still difficult to express quantitatively.
Let us write the basic system of atmospheric equations as follows:
In these equations, x, y, and z are the three Cartesian space coordinates, which are directed eastward, northward, and upward, respectively. The time coordinate is denoted by t. The velocity components in the x, y, and z coordinates, respectively, are denoted by u, v, and w. The horizontal velocity IV denotes ui + vj, where i and j are unit vectors in the x and y coordinates, respectively. ρ represents density, p pressure, s specific entropy, and T temperature. f = 2Ω sin φ and f ' = 2Ω cos φ, where Ω denotes the angular speed of the earth's rotation, and φ denotes latitude. Fx, Fy and Fz, are the three components of the frictional force per unit mass. Q is the rate of heating per unit mass. Cp and Cv, are the specific heat at constant pressure and constant volume, and R is the gas constant of dry air equal to Cp-Cv. g denotes the acceleration due to gravity. For simplicity the equations for predicting the water vapor field and describing the condesation process have been omitted.
In the system of equations (1) through (7), x, y, z, and t are called the coordinate variables. u, v, w, ρ, and s are called prognostic variables because there are equations for the time rate-of-change of these variables. p and T are called diagnostic variables because they are computed from Eq. (6) and Eq. (7), which do not include time derivatives.
If Fx, Fy, Fxz and Q were expressed in terms of u, v, w, p, ρ, s, and T, these equations, together with the proper boundary conditions, would be complete in the sense that there are the same number of equations as of unknowns.
The usual method of solving such a system of non-linear partial differential equations is to replace the partial derivatives by finite differences in both space and time so that the unknowns, u, v, w, ρ, p, s, and T, are defined only at the points of a four-dimensional lattice in space and time. If the variables u, v, w, ρ, and and s are known at t = t0 at all grid points in space, the right-hand sides of Eqs. (1) through (5) may be computed by appropriate space differences. The resulting point values for the local time derivatives of u, v, w, ρ, and s may then be multiplied by a suitable time increment, Δt, and added to the values of u, v, w, ρ, and s at t= t0 in order to define them at the time t0+Δt. Repetition of this process will, in principle, yield a forecast for any later time. [1] [2] [3] [4]
It turns out, however, that a constraint on the value of Δt exists when these equations are solved by explicit difference schemes. This constraint is known as a computational stability condition, and was first discussed by Courant, Friedrichs, and Lewy [5] in 1928. Roughly speaking, the time increment, Δt, must satisfy the condition
where Δs represents one of the space increments, Δx, Δy, or Δz, and Cm denotes the maximum magnitude of the characteristic velocity of the system. (If implicit difference schemes are used, such a constraint may not appear. However, implicit schemes are usually not economical to use.)
The maximum characteristic speed in this system is the speed of sound, 300 m/sec. If we choose 300 km as a representative horizontal grid distance, we find that Δt is less than 1000 sec, or 15 min. However, sound waves in this system propagate not only in the horizontal direction but also in the vertical. If we choose 3 km as a representative vertical grid distance, we find that Δt must be less than 10 sec! To use such a small time step for numerical integration of a system designed to predict weather patterns for a period of ten days or so is impractical, even with modern high-speed computers. We must, therefore, consider modifications to the basic system of the atmospheric equations in order to eliminate sound waves from the system.
The hydrostatic approximation is defined by
If we substitute Eq. (9) into (3), we find that sound waves are filtered out, but that gravity waves undergo little change (Monin and Obukhov [6]). The maximum speed of gravity waves can be as large as that of sound waves, but the fact that gravity waves propagate only horizontally in the hydrostatic system relaxes the severe computational stability condition mentioned earlier.
The use of the hydrostatic approximation based on scale analysis can also be justified for large- to medium-scale motions in the atmosphere (as defined later). It should be borne in mind that w is no longer a prognostic variable since the equation for the time rate of change of w does not exist.
Because of the hydrostatic equation (9), one can no longer calculate ∂ρ/∂t and ∂s/∂t (or ∂p/∂t) independently from Eqs. (4) and (5) [or from Eq. (5) together with (6)]. This constraint leads to a diagnostic equation for w. In other words, w is so calculated that the computations of ∂ρ/∂t from (4) and ∂s/∂t from (5) satisfy the hydrostatic condition (9) all the time. This important constraint was first noted in 1922 by Lewis F. Richardson [7] in his pioneering work, Weather Prediction by Numerical Process. Therefore, this diagnostic equation for w is now called the Richardson equation. The Eulerian hydrodynamic equations modified by the assumption of hydrostatic equilibrium are called, in meteorology, the primitive equations.
Richardson actually tried numerical integrations of the primitive equations using observed meteorological data available to him at that time. It is no surprise that his forecast was unsuccessful. In retrospect, explaining why Richardson's forecast failed at that time seems equivalent to describing the history of numerical weather prediction. Several major technological and theoretical breakthroughs were necessary before today's successful weather forecasts could be produced. The real value of Richardson's work is that he laid down a sound scientific scheme for weather prediction by assembling some of the diverse knowledge about various physical processes which was available in his day.
For a long time after what Richardson called his glaring error, no one attempted to follow his objective approach to weather forecasting. I would like to quote a famous sentence from Richardson's book:
"Perhaps some day in the dim future it will be possible to advance the computations faster than the weather advances and at a cost less than the savings to mankind due to the information gained. But that is a dream."
(The reader may refer to an interesting review of Richardson's book by Platzman.[8])
In the mid-1940s, John von Neumann, of the Institute of Advanced Studies at Princeton, began to build an electronic computer, primarily for the purpose of weather prediction. At the Institute, a meteorological research group was established to attack the problem of numerical weather prediction through a step-by-step investigation of models designed to approximate more and more closely the real properties of the atmosphere. In 1950, results of the first successful numerical prediction were published by Charney, Fjortoft, and von Neumann [9] Their atmospheric model was not based upon the primitive equations, but was rather a quasi-geostrophic model, in which even gravity waves are filtered out. I shall explain the quasi-geostrophic model later.
So far, I have not been specific about the scale of atmospheric motions with which we are concerned. As a matter of fact, in the days of Richardson, meteorologists did not have a clear idea about the scale of weather systems and the principal physical processes involved in a particular scale of motions. Perhaps this is why Richardson was reluctant to make even minor approximations and why he considered practically all the physical processes involved, including the effects of condensation nuclei, evaporation from soil and transpiration from plants, even for a one-day forecast. Recent successes in numerical weather prediction have depended upon recognition of the different scales of motion which govern weather systems, and upon simplification of the Navier-Stokes hydrodynamic equations by mathematical approximations in order to describe motions of a particular scale. Table I shows the scales of motion associated with different atmospheric phenomena.
| Classification | Scale(km) | Typival Phenomena |
|---|---|---|
| large scale | 10,000 | very long waves |
| large scale | 5,000 | cyclones and anticyclones |
| medium scale | 1,000 | frontal cyclones |
| medium scale | 500 | tropical cyclones |
| meso scale | 100 | local severe storms, squall lines |
| meso scale | 50 | hailstorms, thunderstorms |
| small scale | 10 | cumulonimbi |
| small scale | 5 | cumuli |
| small scale | 1 | tornados, waterspouts |
| micro scale | 0.5 | dust devils, thermals |
| micro scale | 0.1 | dust devils, thermals |
| micro scale | 0.05 | transport of the heat, momentum, and water vapor |
| micro scale | 0.01 | transport of the heat, momentum, and water vapor |
Motions of a scale greater than a few thousand kilometers are known as large-scale motions. Figure 1 shows a pressure map of the Northern Hemisphere for a particular winter day. The lines are isobars drawn at 5mb intervals at a height of 6 km. The symbols H and L denotes the high and low pressure areas, called cyclones (or troughs) and anticyclones (or ridges), depending upon the configuration of the patterns. For this scale of motion, the relation between the pressure field and the field of motion is given, to a first approximation, by the geostrophic wind law. This law states that the air moves parallel to the isobars, clockwise around a high pressure area and counterclockwise around a low pressure area in the Northern Hemisphere, and in the opposite direction in the Southern Hemisphere. The pressure force is, of course, directed at right angles to the isobars, from higher to lower pressure. In the mid-tropo sphere, friction may be neglected, and the pressure force is ordinarily balanced by the Coriolis force (the deflecting force which results from the earth's rotation and which acts at rght angles to the motion).
Since the pressure generally decreases from the equator to the North Pole, the wind generally blows from the west to the east. It is called westerly when the flow is from west to east and easterly when the flow is from east to west. Figure 1 shows the prevailing mid-latitude westerlies meandering around the earth. The portion of the zonal flow where the intensity is maximum is called the jet stream.
Large-scale meandering motions, of the type seen in Figure 1, are called long waves. They were discovered from the upper air analyses during World War II. The work of J. Bjerknes and C. Rossby led to the important conclusions that (1) long waves are closely related to weather systems, and (2) long waves are characterized by a state of quasi-balance between the horizontal wind and pressure fields. The latter observation gave Charney [10] and Eliassen [11] a basis for deriving the concept of a quasi-geostrophic approximmation to formulate a prediction model which was later used successfully by the Princeton meteorological group. This was the beginning of a new era, not only for weather forecasting, but also for meteorology in general.
In 1955, the National Meteorological Center, U.S. Weather Bureau, began to issue numerical forecasts on an operational basis by means of an electronic computer. Until recently, their forecasting model has been of the quasi-geostrophic type mentioned earlier . [12][13]
There were two major problems which dynamic meteorologists began to worry about: (1) Although numerical forecasts of the mid-tropospheric flow patterns were superior to the subjective forecasts by skilled forecasters, this was not true of the numerical forecasts of surface flow patterns. (2) Very long waves tended to move too rapidly westward in the numerical forecasts unless the magnitude of a certain parameter in the fore-models was increased beyond an acceptable range. [14]
Figure 2 shows a sea-level pressure map of the Northern Hemisphere for the same date as Figure 1. The lines are isobars, drawn at 5 mb intervals. Compared with the flow pattern at 6 km of Figure 1, the sea-level pressure map shows the presence of smaller scale motions. As shown in Table I, motions of a scale on the order of 1000 km are called medium-scale motions. Typical phenomena are tropical cyclones in the tropics (Figure 3) and frontal cyclones in middle to high latitudes.
Figure 4 shows a schematic model of frontal cyclones (after J. Bjerknes and H. Solberg [15]). Frontal cyclones form on fronts separating a cold air mass to the north from a warm air mass to the south. Figure 5 shows diagrammatically the relation between fronts and cyclones at sea level and the long waves aloft (after Palmen [16]).
For motions of medium scale and smaller, the deviation of wind from a geostrophic balance gradually becomes important. One cause for this deviation is surface friction. Since these motions are dominant in the lower part of the troposphere, they are influenced by conditions at the earth's surface. When a current of air moves over an uneven surface, the motion usually becomes irregular or turbulent. Because of surface friction, the kinetic energy of the flow is dissipated by eddy actions. Not only does the surface of the earth act as a sink of kinetic energy, but it also plays a role as an energy source. When cold air flows over a warm ocean, heat is transported from the ocean to the atmosphere by the turbulent actions of conduction and convection.
Another important heat source for motions of medium scale and smaller is the latent heat released due to the condensation of water vapor. This heat source is particularly important for the development of tropical cyclones, hurricanes, and typhoons. Local severe storms, such as hailstorms and thunderstorms, are called mesoscale motions; and the latent heat released in gigantic cumulonimbi is a major energy source for such systems. It is clear that we must incorporate various heat sources and energy sinks in a forecasting model of flow patterns near the surface of the earth.
Coming back to the other difficulty, namely why the propagation of very long waves is overpredicted by quasi-geostrophic models, it was later found that there are two types of geostrophic motion in the earth's atmosphere.[17] As noted earlier, motions on the order of 5000 km can ordinarily be in a state of geostrophic balance. This is also true in the case of very long waves; but because of the large scale, variation of the Coriolis parameter with latitude becomes very important. The upshot was that quasi-geostrophic models were inaccurate in describing the motions of very long waves and that the spherical geometry of the earth must be taken into account. This led to the conclusion that the primitive equations on spherical coordinates must be used in predicting the global circulation system days in advance. Thus, we have returned to the position of L. F.Richardson a half century ago.
With the advent of high-speed electronic computers and the development of numerical weather prediction techniques, dynamic meteorologists have begun looking into the fundamental problem of the general circulation of the earth's atmosphere. Phillips [18] made the first successful attempt to simulate large-scale motions through a step-by-step integration of quasi-geostrophic equations. Despite many simplifications employed in setting up the model, the experiment demonstrated that certain gross properties of atmospheric motion, including the energy-transformation processes (which will be discussed later), could be numerically simulated by the model. Later Smagorinsky [19] performed a similar experiment using the primitive equations.
Following these earlier attempts, more elaborate types of general circulation experiments using the primitive equations were designed and carried out, [20][21][22][23] In 1964, when Kasahara and Washington decided to work on a numerical experiment of the general circulation of the atmosphere at the National Center for Atmospheric Research (NCAR), we made a survey of the characteristics of various general circulation models thus far reported, and we found that height as the vertical coordinate (rather than pressure, as customarily done) would be convenient in formulating a prediction model. similar to Richardson's [7]. Besides the different vertical coordinate, the NCAR model includes various physical processes in the atmosphere which have been formulated somewhat differently from those incorporated in other general circulation models.[24]
Let us now consider the basic physical processes involved in the mechanism of the general circulation. All energy of atmospheric motion is ultimately derived from incoming solar radiation. If directly reflected radiation is ignored, the atmosphere receives from the sun 2.0 cal cm-2 min-1, which is known as the solar constant. This gives a mean flux of solar energy perpendicular to the earth's surface of about 0.5 cal cm-2 min-1. (The factor 4 is the ratio of surface area to cross-section for a sphere.) Of this, approximately 35 percent is reflected from the atmosphere (including clouds) and portions of the earth's surface. This percentage of energy reflected back to space is called the albedo of the earth as a whole. Approximately 20 percent is absorbed in the atmosphere, and 45 percent is absorbed at the surface. Since the albedo of the earth is 0.35, an average of 0.65×0.5 or about 0.3 cal cm-2 min-1 is available for direct or indirect heating of the earth and the atmosphere [25].
Since the mean temperature of the atmosphere and the earth do not change appreciably from one year to the next, the amount of energy received by the earth must be sent back to outer space. This ultimate return of energy to space is in the form of low-temperature terrestrial radiation. Figure 6 shows the two curves. One represents the annual mean of solar radiation absorbed by the earth and the atmosphere, and the other the outgoing terrestrial radiation from the atmosphere, both as functions of latitude.[26][27]
At low latitudes the earth-atmosphere system gains more heat energy per unit area by the absorption of short wave radiation from the sun than it loses to space by the emission of long-wave radiation; the reverse happens at high latitudes. The net radiation gained by the earth-atmosphere system at the equator is, for example, about 100 ly/day, while the net loss at the North Pole is about 260 ly/day. At a latitude of 38°, the incoming and outgoing radiation are balanced.
ly is an abbreviation for langley, a unit equal to one gram calorie per square centimeter.
It is clear that the surplus heat energy in the tropics must somehow be carried to the poles. Otherwise, tropical regions would become steadily hotter and polar regions steadily colder. About ten percent of the surplus heat is transported poleward by ocean currents. The remainder is transported by the atmospheric circulation, including both the large-scale eddy motions (which are cyclones and anticyclones) and the mean meridional circulations.
Since atmospheric motions are generated by non-uniform heating of the air, a physical prescription is needed for the heating rate, Q, in the thermodynamic equation (5). In general, the following three processes (expressed as rates) are important in the atmosphere (Figure 9): heating from radiational sources, denoted by Qa; the release of latent heat by condensation, denoted by Qc; and heating/cooling due to eddy diffusion, denoted by Qd.
The rate of heating due to radiational sources, Qa, may be divided into two parts: one is the rate of heating due to absorption of solar insolation in the atmosphere, Qas, and the other is the rate of heating/cooling due to infrared radiation, Qal These terms can be evaluated by integrating the transfer equations, taking into account the distributions of water vapor, carbon dioxide, and ozone, the major absorption gases in the atmosphere.
However, the amounts of Qas, and Qal are considerably altered by the presence of clouds in the atmosphere. The evaluation of Qc the rate of released latent heat of condensation, is, of course, directly dependent on the prediction of cloud formation. Thus, the prediction of clouds is one of the important aspects of a general circulation model.
The moisture field of the atmosphere is continually replenished by evaporation from the earth's surface. Since heat is required to evaporate water, and heat is released by condensation of water vapor, the evaporation process transports latent heat from the earth's surface to the atmosphere. Similarly, sensible heat is transported from the surface to the atmosphere, and vice versa, by the actions of small-scale turbulence and convection.
The direction of energy flow, in the form of latent heat and sensible heat, is determined primarily by the difference between the temperature of the earth's surface and the air temperature directly above it. Thus, it is necessary to calculate the temperature of the earth's surface.
The temperature of the earth's surface is determined by the balance of energy flux at the earth's surface, as shown in Figure 8. The surface of the earth absorbs solar radiation, Sg, and downward infrared radiation, D, from the atmosphere. A portion of the absorbed energy is transported down into the ground by heat conduction in the soils or is carried away by the motion of water in the oceans. This quantity is denoted by M.
Energy is also carried away in the form of sensible heat, Es, and of latent heat, El. The magnitudes of Es and El are dependent on both the thermal and dynamical conditions of the atmosphere and the properties of soil and water.
The remainder of the absorbed energy is re-radiated to space as infrared radiation from the earth's surface. This amount is expressed by ρTg4, where ρ is the Stefan-Boltzmann constant and Tg, is the temperature of the earth's surface which we wish to determine.
The evaluation of M requires knowledge of the mechanism of heat and moisture transports in the lithosphere - and also of the circulation of water in the hydrosphere. Oceanic circulation is important in determining M as well as Es. This is why there is a growing awareness of the role of the ocean's circulation in meteorology.
Perhaps equally important to meteorology is the role of soil. Knowledge of soil physics, particularly the transport processes of heat and moisture, will be essential when we try to explain, for example, the formation of deserts in connection with the long-term effects of the atmospheric circulation upon the earth's surface.
It is clear that what happens in the first hundred meters above the earth's surface is very important in accounting for the sources and sinks of energy for the atmosphere. Micrometeorology is a branch of meteorology dealing with various phenomena in the earth's boundary layer. Figure 9 shows a schematic structure of the lower atmosphere (after Sutton [28]).
In the surface boundary layer, which extends up to 50-100m above the surface, shearing stress of wind and heat flux are both constant with height, and the strueture of the wind is primarily determined by the characteristics of the surface and the vertical gradient of temperature. This region is referred to as the Prandtl layer.
In the planetary boundary layer, which extends approximately 500-1000 m above the ground, the structure of the wind is influenced not only by the pressure gradient, density stratification, and surface friction, but also by the earth's rotation. This region is referred to as the Ekman layer.
After the heating rate, Q, and the frictional terms, Fx, and Fy, are suitably expressed in terms of dependent variables, the primitive equations can be integrated from prescribed initial conditions under well-posed boundary conditions.
Now I would like to discuss some computational aspects of our problem, specifically the six-layer version of the NCAR general circulation model. The model is global. Horizontal grid points are placed at 5 degree longitude and latitude intersections. However, the grid system is staggered in time, as shown in Figure 10. The cross and dot points are at different time levels of one time increment, Δt.
Figure 11 shows a vertical cross-section of the model.
The height increment, Δz, is chosen to be 3 km. The top is located at z = zT (= 6 Δz), which corresponds approximately to the 75 mb pressure level. It is assumed that vertical motion vanishes at the top in order to conserve the mass of the model atmosphere. The lower boundary of the model is a curved surface which takes into account the earth's orography (indicated by shading).
Surface pressure, ps, is evaluated along the mountain surface. Pressure, p, and vertical velocity, w, are placed at z = 3,6,9, 12,and 18 km; and wind components, u and v, temperature, T, density, p, and the mixing ratio of water vapor, q, are placed at the intermediate levels.
The effects of both the surface boundary layer and the planetary boundary layer are included in the model, as shown in the inset of Figure 11. The variables us, vs , Ts, qs, and ps, which correspond to variables in the surface boundary, are evaluated from the lower boundary conditions.
In the model, the earth's surface is divided into three parts: oceanic regions, continents, and any region where the surface temperature is below the freezing point (0°C), be it land or ocean. The third region is, therefore, considered to be covered by either ice or snow.
The ocean surface temperature is prescribed from the distribution of observed climatological mean temperatures for the appropriate season of the year. The surface temperature of the land masses is computed by the method previously described. Over ice covered regions, an upper limit of 0°C is assumed, with the implication that any excess heat flux goes into the latent form [21].
One unique feature of the atmospheric equations is that the time rate-of-change of the variables is about an order of magnitude. less than the rest of the terms in the equations. For example, the acceleration of air in large-scale motions is about one-tenth the magnitude of the Coriolis or the horizontal pressure gradient force. This fact is in striking contrast to other hydrodynamic problems, particularly shock wave problems in which large accelerations are involved. Thus, in predicting large-scale motions, extra care must be taken to compute individual terms invloving space derivatives, since the time tendency is evaluated as the sum of individual terms in the equations.
Another special requirement for the finite difference prediction equations is that the model be run long enough to study long-term fluctuations of weather systems. To produce a 100 day forecast with a time step of 6 min requires over 20,000 iterations. This is why strongly dissipative difference schemes are not preferable to use for numerical weather prediction. Of course, any practical finite difference scheme must be both stable and accurate. For long-term numerical integration of the atmospheric equations, it is, therefore, necessary to design finite difference schemes which conserve certain properties of the original equations, such as momentum and energy.
Even though a finite difference equation is stable in the sense of the von Neumann condition [29], there is no guarantee of stability if such a scheme is applied to non-linear equations. As an example, consider the very simple hyperbolic nonlinear equation.
When the centered scheme is applied in order to approximate both the time and space derivatives, the finite difference analog of the above equation may become unstable, even under the usual stability condition:
Since this instability typically appears in nonlinear equations, it is called nonlinear instability [29][30]. However, a similar kind of instability may occur even in linear equations if the coefficients are not constant. Phillips [30] showed that instability arises because the grid system cannot resolve wavelengths shorter than two grid intervals. When such wavelengths are formed by the nonlinear interaction of longer waves, the grid system incorrectly interprets them as long waves.
This type of instability is inherent in all meteorological equations, and many means to eliminate this difficulty have been proposed. Our general circulation model is based on a scheme which was originally proposed by Lax and Wendroff [32] and later modified by Richtmyer [31] who called it the two-step Lax-Wendroff scheme. The scheme consists of a combination of the first step, called diffusing, and the second step, called leapfrog. The diffusing step is strongly dissipative, but the leapfrog step is neutral in stable conditions. Because of the dissipative nature of the first step, the combined steps provide effective smoothing for shorter waves and suppression of nonlinear instability. We found, however, that the two-step method damps the solution too severely. In order to eliminate this difficulty, we apply the diffusing step less often, so that short waves are not over-damped. By experimenting, we have found it sufficient to apply the diffusing step only after every 45 steps of leapfrog.
In order to give an idea of the magnitude of computation involved in our general circulation experiments, let us consider the amount of data which must be stored in the computer. In each horizontal plane, there are 72×36 = 2592 data points for a 5 degree longitude and latitude mesh. There are seven levels, including the surface, so that the total number of data points in three dimensional space is 2592×7 = 18144. Since the model is staggered in horizontal space, only half of this number (9072) is needed at one time. There are four prognostic variables, u, v, p, and q, which are needed at two consecuive time levels. Here q denotes the mixing ratio of water vapor in the atmosphere. There are three diagnostic variables, w, ρ, and T. We need 15 t0 30 extra fields of heating, friction, and similar terms for diagnostic purposes, in addition to such physical parameters as cloudiness and precipitation, depending on the complexity of the model. Thus, a total of more than two million pieces of information plus about 12,000 arithmetic and logic computer instructions must be stored in the computer memory. Of course, there is no need to store all of these data in the high-speed core memory all the time during a computation. Since computation is done in a sequential manner, only a portion of the data is needed at one time for the calculations.
At NCAR, we use a Control Data 6600 computer which has only 65,000 words of core memory. The grid point data are stored in six magnetic drums, each of which contains about one-half million words. Data are constantly read from the drums to the core and transferred from the core to the drums simultaneously with the arithmetic computations. About 150 floating operations are required at each mesh point, and it takes about 25 sec for one time step. Since a time step of 6 min is used, it takes about 90 min of computing time to produce a one-day forecast.
Because a large amount of data are involved, there is a serious problem in processing the output. The situation is in contrast to, for example, an eigenvalue calculation in which the output may be a single number after hours and hours of computation. Although it is possible to print out a large amount of numbers in a short time, it is practically impossible to digest the results and evaluate the experiment without resort to graphic analyses of numbers to display the patterns. At NCAR, the output data processing is aided by a cathode ray tube (CRT) plotter called the Data Display 80 (dd80). Using this equipment, grid point data are analyzed and graphic contour lines are then superimposed on a suitable map projection.
Even though the output is analyzed in the form of weather maps, it is still difficult to comprehend the evolution of the flow patterns in a short time. Here the idea of motion pictures is very useful. The making of computer-generated movies is relatively new. Welch [33], for example, describes a computer movie of numerical simulation of water waves from a breaking. dam. Washington, O'Lear, Takamine, and Robertson [34] discuss the techniques of making black-and-white and color computer movies of the output of the NCAR general circulation experiments. We have found the use of color very effective in illustrating simultaneously the patterns of different variables for comparsion. Color film is processed from a black-and-white original film by a multi-exposure method using filters of different colors. Eventually, we hope to have an output device to demonstrate three-dimensional patterns. We will then be able to increase greatly the information content of a single graphical display.
Now I would like to talk about the result of one experiment of our two-layer general circulation calculations. This experiment simulated flow patterns of January (Northern Hemisphere winter). The season is determined by the declination of the sun. In addition, we prescribed a climatological mean distribution of sea surface temperature for January. The calculation started from an isothermal atmosphere of 250°K. This particular value is chosen as the equilibrium temperature of the earth-atmosphere system, which is evaluated from the balance requirement of incoming solar radiation with outgoing infrared radiation. This value is lower than the average for the earth's surface, but is close to the average temperature of the atmosphere itself.
Let us now look at the color computer movie which illustrates the evolution of flow patterns. The upper figure shows the surface pressure distribution, and the lower figure shows the temperature distribution at 3 km. The blue lines delineate continents, the yellow contours are isobars at 5 mb intervals, and the red lines are isotherms at 5°C intervals.
As I explained earlier, at low latitudes the earthatmosphere system gains more energy by the absorption of short-wave radiation from the sun than it loses to space by the emission of long-wave radiation; the reverse happens at high latitudes. Because of this radiational imbalance, the north-south temperature gradients develop. The circulation pattern in the early phase of the development is of large-scale monsoons. Large highs develop over the North American continent and over Siberia, and lows form over the warm North Atlantic and to some extent over the North Pacific Figure 12a, b). Meanwhile, the warming in the low latitudes slows down, but the cooling in the higher latitudes, particularly near the North Pole, continues. Since the map is a cylindrical equal-spaced projection, the patterns near the poles are greatly exaggerated.
These large monsoonal circulations are too ineficient to transport heat northward to make up the radiational imbalance. As a result, the north-south temperature gradients kept increasing. After approximately 15 days, the north-south temperature gradients became so large that the large-scale monsoonal circulation broke down into eddy motions. These eddies are nothing but cyclones and anticyclones, the familiar systems that appear on weather maps (Figure 13a, b).
The breakdown of the monsoonal circulations into eddy motions is known in meteorology as baroclinic instability, and these eddies are called baroclinic waves. (The nature of baroclinic instability of long waves was first theoretically analyzed by Charney [35] and Eady [36].)
Having noted some of the principal features of the evolution of global circulation, let us see to what extent those features can be explained. Since there is a net warming in low latitudes and a net cooling in higher latitudes in the troposphere, the pressure aloft is, in general, higher in low latitudes and lower in high latitudes. Thus, the presence of westerly flows can be anticipated in the middle latitudes, and the intensity of westerlies must increase with respect to height in the troposphere.
In addition to transporting heat from low latitudes to high latitudes, the meridional circulation and the superimposed eddy motions also transport momentum and moisture poleward. Thus, we can infer that there must be a source of momentum in the lower latitudes and a sink somewhere in the middle or high latitudes. It turns out that the frictional drag of the surface of the earth gives rise to a source of zonal mean westerly momentum in those latitudes with easterly winds, and a sink in those latitudes with westerly winds. Therefore, the presence of surface easterlies in the low latitudes can be explained from the point of view of momentum balance. In middle latitudes, surface westerlies play the role of momentum sink.
In reality, the easterlies in the low latitudes extend throughout the troposphere, though the latitudinal width is very much reduced in the upper troposphere. The presence of easterlies in the upper troposphere is difficult to explain at present. In our general circulation experiment simulating January climatology, we were able to simulate surface easterlies in the low latitudes, but upper easterlies were not present. This failure may be due to the fact that our two-layer model lacks an efficient mechanism for transporting easterly momentum upward from the earth's surface. There are many indications that the convective activites of cumulonimbi in lower latitudes provide this mechanism.
Here is an interesting point in connection with the poleward transport of momentum. In the mid-troposphere, there is a belt of easterlies in low latitudes and one of westerlies in middle latitudes. Yet westerly momentum is transported poleward. This means that momentum is transported against the north-south gradient of momentum. It is obvious, then, that the mechanism of momentum transport is not one of eddy diffusion, for we would have to assume a negative coefficient of turbulent viscosity in order to permit the counter-gradient flow.
To summarize the mechanisms of the general circulation of the atmosphere, let us look at Figure 14, which shows schematically the energy transformation. Each of three boxes represents a type of energy. The notation I indicates the zonal (east-west) average of internal energy of the atmosphere. Since the kinetic energy is a squared quantity, zonal averaging of the kinetic energy can be divided into zonal kinetic energy, K, and eddy kinetic energy, K'. The zonal average of the potential energy, P, is added into the K box to form the zonal mechanical energy.
In the paper I , K and P have a bar over the letter.
The net energy source, Q, contributes heat to the air columns. As a result of differential heating, cyclones and anticyclones are maintained. Thus, energy flows from the I box to the K' box. The mechanism of energy flow is due to baroclinic instability. A portion of the eddy kinetic energy is lost through the action of dissipation. The eddy motions then feed momentum into the zonal motion indicated by K + P. This process of energy flow is accomplished by the counter-gradient momentum transport mentioned earlier, and is due to barotropic stability. Here the terms instability and stablity are defined according to the responses of the eddy motions to the zonally averaged motions. If eddy motions grow at the expense of the zonally averaged motion, the resulting phenomenon is called instability. If eddy motions feed their energy into the zonal motion, the resulting phenomenon is called stability. The importance of the role of barotropic stability in connection with the mechanism of the general circulation was first pointed out by Kuo [37].
Again there is an energy sink for K + P due to dissipative actions. There is also energy flow between K + P and I. If we consider a water-filled dishpan having uniformly heated sides, we expect radial flow to be generated between the center of the dishpan and its side. If the dishpan is made to rotate, the radial motion of the flow may cause the generation of tangential motion as the result of the deflecting force of the system's rotation. Hadley [38], using this idea, was the first to explain the trade winds system [39]. If zonal kinetic energy is generated as a result of meridional motion, such a meridional circulation is called direct circulation, or Hadley-type circulation. On the other hand, if zonal kinetic energy is consumed in order to maintain meridional circulation, such a circulation is called indirect circulation.
In middle and high latitudes, eddy meridional transport of kinetic energy seems strong enough so that the energy is even driven back to the I box, as first suggested by Starr [40]. However, in tropical latitudes, the conversion of energy from K' to K + P is rather small, and, therefore, direct circulation converts zonal internal energy into zonal kinetic energy [41]. Many investigations have been carried out to determine the magnitude of energy conversion terms. However, scarcity of upper air data over the tropics and the Southern Hemisphere still prevents us from performing satisfactory diagnostic analyses on the heat, momentum, and energy conversions in the global atmosphere.
We have seen that numerical experiments help us to understand many questions about the mechanisms of the general circulation. The degree of understanding may be measured by our ability to simulate, with a model, actual atmospheric motions. A successful model, therefore, must be capable of predicting actual weather changes. Or, put in another way, we judge how good a model is by verifying predicted patterns. For this reason, we have also been working on the use of real data in order to check the forecasting capability of our general circulation model. So far, the results of our forecasts are very encouraging. Incidentally, it is stil difficult to prepare global upper air maps due to the paucity of observational data over the oceans and over the Southern Hemisphere in general. This means that it will be a long time until satisfactory short-range weather predictions are achieved because of the lack of accurate initial conditions for forecasting models.
Every year in the Untied States alone, weather hazards, such as hurricanes, severe storms, floods and drought, kill approximately 1300 people and cause property damage in excess of $11 billion [42]. When this problem is considered on a global basis, it is obvious that improved weather prediction can lead to immense social and economic benefits. In the last few years there has been an encouraging development toward eliminating deficiencies in the present international system of global weather observation. The program is called the World Weather Watch (WWW), and it will be carried out during 1972 under the sponsorship of the World Meteorological Organization (WMO). Another important and related plan is called the Global Atmospheric Research Program (GARP), which is the research-oriented partner of WWW. GARP is sponsored jointly by the International Council of Scientific Unions and the International Union of Geodesy and Geophysics. The goal of GARP is to promote a vastly improved understanding of the general circulation of the atmosphere [43].
One year's operational forecasts at the National Meteorological Center, Weather Bureau, ESSA, have demonstrated that their new six-layer primitive equaition model significantly improves the accuracy of short-range forecasts. Their experience with this model, which is as complicated as a general circulation model [44], clearly indicates that even a short-range forecasting model must take into account detailed physical processes and that only by doing so can an improvement in the accuracy of weather forecasts be expected. It should be borne in mind that such an attempt has become possible only recently through the use of largecapacity and high-speed computers, such as the Control Data 6600 and the IBM 360.
Very recent studies by NCAR staff members and others have demonstrated that the accuracy of numerical weather simulation will increase significantly if a horizontal grid mesh finer than the current 5 degree mesh is used. A 2½ degree mesh will increase computing time by a factor of eight. It is obvious that to run a general circulation model which includes ocean circulation and detailed atmosphere-ground boundary layer calculations with a finer horizontal resolution than that presently used will require a supercomputer 500 times faster than the Control Data 6600 [45]. A few supercomputers are being developed, the best known of which is the University of Illinois' parallel-processing ILLIAC IV. ILLIAC IV will be available by 1972.
Using upper air data gathered during 1972 by the World Weather Watch, together with a supercomputer, we will soon be able to attack the fundamental problem of an objective long-range weather forecast for two weeks and beyond. However, in order to achieve meaningful long-range forecasts, we will first have to learn more about the mechanism of long-term weather changes, an important, yet very poorly understood, problem. Results of our general circulation experiments have already indicated the importance of the condition of the earth's surface, including dynamic and thermal effects, upon atmospheric motions. We have begun to see the important effects of clouds interacting with atmospheric radiation and the equally important effects of small-scale turbulence and convection in the atmospheric boundary layer. These problems must be investigated in even greater detail before the various effects can be properly included in general circulation models.
[1] N A PHILLIPS, Numerical weather prediction, In Advances in Computers ed F L Alt, M Rubinoff, Academic Press New York 1960 Vol 1 43-90
[2] P D THOMPSON, Numerical weather analysis and prediction, The Macmillian Co, New York 1961 170 pp
[3] I A KIBEL, An introduction to the hydrodynamic methods of short period weather forecasting, Translated from the original volume published in 1957 Moscow, The Macmillian Co New York 1963 383 pp
[4] G I MARCHUK, Numerical methods in weather forecasting,In Russian Gidromet Izdat Leningrad 1967 356 pp
[5] R COURANT, K O FRIEDRICHS, H LEWY, Uber die partellen differenzenqleichungen der physik, Math Annalen 100 32-74 1928
[6] A S MONIN, A M OBUKHOV, A note on general classification of motions in a baroclinic atmosphere, Tellus 11 159-162 1959
[7] L F RICHARDSON, Weather prediction by numerical process, Cambridge University Press 1922 236 pp. Reprinted by Dover Publications New York with an introduction by S Chapman 1965
[8] G W PLATZMAN, A retrospective view of Richardson's book on weather prediction, Bull Am Meteorol Soc 48 514 550 1967
[9] J G CHARNEY, R FJORTOFT, J VON NEUMANN, Numerical integration of the barotropic vorticity equation, Tellus 2 237-254 1950
[10] J G CHARNEY, On the scale of atmospheric motions, Geofys Publikasjoner 17 2 17 pp 1948
[11] A ELIASSEN, The quasi-static equations of motion with pressure as independent variable, Geofys Publikasjoner 17 3 44 pp 1949
[12] P D THOMPSON, W L GATES, A test of numerical prediction methods based on the barotropic and two-parameter baroclinic models, J Meterol 13 127-141 1956
[13] G P CRESSMAN, A three level model suitable for daily numerical forecasting, Tech Memo No 22, National Meteorological Center Weather Bureau ESSA 1963 22 pp
[14] G P CRESSMAN, Barotrophic divergence and very long atmospheric waves, Monthly Weather Review 85 293-297 1958
[15] J BJERKNES, H SOLBERG, Life cycle of cyclones and the polar front theory of atmospheric circulation, Geofys Publikasjoner 3 1 18 pp 1922
[16] E PALMEN, The aerology of eztratropical disturbances, In Compendium of Meteorology ed T F Malone, American Meteorological Society 1951 599-620
[17] A P BURGER, Scale consideration oj planetary motions of the atmosphere, Tellus 10 195 205 1958
[18] N A PHILLIPS, The general circulation of the atmosphere: A numerical experiment, Quart J Roy Meteorol Soc 82 123-164 1956
[19] J SMAGORINSKY, General circulation experiments with the primitive equations: I The basic experiment, Monthly Weather Review 91 99-164 1963
[20] C E LEITH, Numerical simulation of the earth's atmosphere, In Methods in Computational Physics Vol 4 ed B Alder et al, Academic Press New York 1965 1-28
[21] Y MINTZ, Very long term global integration of the primitive equations of atmospheric motion, In Proc WMO/IUGG Symposium on the Research and Development Aspects of Long Range Forecasting Boulder Colo 1964 WMO Tech Note No 66 1965 141-167
[22] S MANABE, J SMAGORINSKY, R F STRICKLER, Simulated climatology of a general circulation model with a hydrological cycle, Monthly Weather Review 93 769--798 1965
[23] J SMAGORINSKY, S MANABE, J L HOLLOWAY JR, Numerical results from a nine level circulation model of the atmosphere, Monthly Weather Review 93 727-768 1965
[24] A KASAHARA, W M WASHINGTON, NCAR global general circulation model of the atmosphere, Monthly Weather Review 95 389-402 1967
[25] R M GOODY, Atmospheric Radiation Vol 1 Theoretical Basis, Clarendon Press Oxford, 1964 436 pp
[26] H G HOUGHTON, On the annual heat balance of the northern hemisphere, J Meteorol 11 1-9 1954
[27] J LONDON, A study of the atmospheric heat balance Contract AF 19 122-165 ASTIA No 117227, Department of Meterology and Oceanography New York University 1957 99 pp Final Report
[28] O G SUTTON, Micrometeorology, McGraw-Hill Book Co New York 1953 333 pp
[29] R D RICHTMYER, K W MORTON, Difference methods for initial-value problems, Interscience Publishers New York 1967 405 pp
[30] N A PHILLIPS, An example of non-linear computational instability, In The Atmosphere and the Sea in Motion Rossby Memorial volume ed B Bolin Rockefeller Inst Press New York 1959 501-504
[31] R D RIGHTMYER, A survey of difference methods for'non-steady fluid dynamics, Tech Note 63-2 National Center for Atmospheric Research, Boulder Colo 1963 25 pp
[32] P D LAX, B WENDROFF, Systems of conservation laws, Comm Pure Appl Math 13 217-237 1960
[33] J E WELCH, Computer simulation of water waves, Datamation 12 41-47 1966
[34] W M WASHINGTON, B T O'LEAR, J TAKAMINE, D ROBERTSON, The application of CRT contour analysis to general circulation experiments, To be published in BulAm Meteorol Soc 1968
[35] J G CHARNEY, The dynamics of long waves in a baroclinic westerly current, J Meteorol 4 135-162 1947
[36] E T EADY, Long waves and cyclone waves, Tellus 1 33-52 1949
[37] H L KUO, Dynamic instability of two-dimensional non-divergent flow in a barotropic atmosphere, J Meteorol 6 105-122 1949
[38] G HADLEY, Concerning the cause of the general trade wind, Phil Trans Roy Soc 39 58 62 1735 Reprinted by C Abbe Smithsonian Misc Coll 51 5-71910
[39] E N LORENZ, The circulation of the atmosphere, Am Scientist 54 402-420 1966
[40] V P STARR, An essay on the general circulation of the earth's atmosphere, J Meteorol 5 3943 1948
[41] E PALMEN, General circulation of the tropics, In Proc Symposium on Tropical Meteorology, Rotoru, New Zeland Nov 1963 3--30 1964
[42] R M WHITE, The world weather program, Bull Am Meteorol Soc 48 81-84 1967
[43] W O ROBERTS, The global atmospheric research program, Bull Am Meteorol Soc 48 85-88 1967
[44] F G SHUMAN, J B HOVERMALE, An operational six-layer primitive equation model, J Appl Meteorol 7 525-547 1968
[45] H G KOLSKY, Computer requirements in meteorology, Tech Rept No 38002 IBM Scientific Center Palo Alto Calif 1966158pp