A streamlined family photochemistry module reproduces major nonlinearities in the global tropospheric ozone system

Pergamon

00!37-8485(%)00054-2

Compurers Chem. Vol. 20, No. 2, pp. 235-259, 19 8.i Copyright Q 1996 Elsevier Science 1 !d Printed in Great Britain. All rights mened 0097-8485/96 $15.00 + 0.00

A STREAMLINED FAMILY PHOTOCHEMISTRY MODULE REPRODUCES MAJOR NONLINEARITIES IN THE GLOBAL TROPOSPHERIC OZONE SYSTEM SCOTT ELLIOTT,‘* ME1 SHEN,Z C. Y. J. KAO,’ R. P. TURCO* and MARK Z. JACOBSON’ ‘Earth and Environmental Sciences Division, Geoanalysis Group, Los Alamos National Laboratory, Los Alamos, NM 87545, U.S.A. ZAtmospheric Sciences Department, University of California, Los Angeles, CA 90024, U.S.A. lEnvironmenta1 Fluid Mechanics and Hydrology, Civil Engineering Department, Stanford University, Stanford, CA 94305, U.S.A. (Received 28 November 1994; in revised form I7 April 1995) Abstract-Tropospheric photochemistry is central to terrestrial climate change and pollution effects and so will be modelled on global 3-D grids. The chemistry is also complex, numerically stiff and kinetically nonlinear. A packaged family based integrator has been developed specifically to combat the difficulties associated with computational modelling of atmospheric chemistry on the global scale. The present work describes its ability to reproduce the major nonlinear features of tropospheric kinetics-those relating the nitrogen oxides (NO,) and oxidizing organics to ozone. It is shown that the family modules can duplicate typical changes in ozone production as a function of NO, level while consuming a minimum number of mathematical operations. The tests are first performed in the box model mode for a variety of pristine and pollutant scenarios. Zero-dimensional runs are patterned largely after the nonlinearity investigations of Liu and coworkers. The testing is then extended to column representations for vertical mixing of ozone precursors in convective storms. Here the calculations follow the climatology of ozone production enhancements assembled by Pickering and colleagues. Benchmarking is reported for a mechanism containing full inorganic kinetics as well as decomposition sequences for six nonmethane hydrocarbons. Chemical species in the simulation number 92. The operations count is roughly 10,000 per cell step for time increments of 1h or more. The coding should thus enable decadal scale runs on massively parallel processors. Scaling experiments indicate full vectorization has been achieved. The chemistry packages are optimized not only for speed but also for convenience. Modularity and routines automating setup of solutions to the kinetic continuity equations are outlined as incidentals.

INTRODUCTION

The chemistry of the troposphere controls several widespread pollution phenomena such as oxidant dispersion (McKeen et al., 199la, b; Thompson, 1992; Jacob et al., 1993a, b), and is closely tied to earth’s climate through the trace species ozone and methane, which are both radiatively and kinetically active (Ramanathan et al., 1987; Ramanathan, 1988; Hauglustaine et al., 1994). Because of the extensive scales involved, oxidant and greenhouse gas distributions can be effectively modelled on global 3-D grids (Levy et al., 1985; Penner et al., 1991; Crutzen & Zimmerman, 1991; Kao et al., 1992a; Kasibhatla et al., 1993). Indeed, nonlinearities in kinetic relationships between fleeting trace species may require detailed global meshes for quantification (Houghton et al., 1990, 1992; Thompson, 1992; Kanakidou & Crutzen, 1993). At the same time, however, tropospheric chemistry simulations tend to incorporate many dozens of constituents interacting as a stiff * To whom all correspondence should be addressed.

numerical system (Trainer et al., 1987; Donahue & Prinn, 1990; Atherton & Penner, 1990; Kanakidou et al., 1991). The conjunction of full chemistry schemes with global gridding offers a strong challenge to modern supercomputers and programmers (Rood et al., 1987, 1990; Kaye & Rood, 1989; Jacob et al., 1989; Kao et al., 1990; Elliott et al., 1993a). Several distinct approaches to the modelling of atmospheric chemistry have been developed for global meshes, with varying degrees of parameterization (e.g. Jacob et al., 1989; Spivakovsky et al., 1990a, b; Crutzen & Zimmerman, 1991). Our group has experimented with accelerated Gear routines (Turco & Whitten, 1974; Jacobson & Turco, 1994) and also with family based integrators (Turco & Whitten, 1974; Turco, 1985; Kao et al., 1990). The family codes are automated for ease in alteration of reaction mechanism and modularized for attachment to arbitrary transport frameworks (Jacobson & Turco, 1994; Elliott et al., 1993a-c 1994, 1995; Shen et al., 1994). Extensive individual species accuracy tests have been reported for the inorganic chemistry occurring above the tropopause (Elliott et al., 1993a,

236

Scott Elliott er al.

1995) and preliminary global 3-D runs have been conducted within the stratosphere of a general circulation model (GCM; Kao et al., 1990; Zhao et al., 1993, 1994; Elliott et al., 1995). Several of the coauthors on this paper (MS and MZJ; Shen et al., 1994) have recently constructed a complete tropospheric version of the packaged family chemistry, including six nonmethane hydrocarbons which can serve as surrogates for the suite of global organics (Chatfield & Delany, 1990; Pickering et al., 1990, 1991; Elliott et aI., 1993b, 1995). Species by species accuracy verification has been undertaken and will be documented in an upcoming Ph.D. Thesis (MS). The current work reports tests of this tropospheric model at a somewhat higher level within the photochemical hierarchy. We demonstrate that the code can reproduce production rates for the key oxidant and greenhouse gas ozone over a range of tropospheric chemical and meteorological conditions. The nonlinearities in the tropospheric kinetics system which have proven most crucial to understanding ozone distributions regard catalysis of production by the nitrogen oxide radicals (NO, = NO + NOz) during hydrocarbon oxidation (Crutzen, 1973; Chameides & Walker, 1973; Fishman & Crutzen, 1978; Liu et al., 1987; Crutzen, 1988; Lin et al., 1988). Ozone generation is most efficient at intermediate levels of NO and NOz observed between pristine and polluted regimes (Logan, 1983, 1985; Penner et al., 1991). We recreate several prominent box modelling efforts using our specialized integrator to show that it can capture the ozone/NO, relationship in air masses from remote oceanic, remote continental and urban areas (Levy et al., 1985; Trainer et al., 1987; Liu et al., 1987). Mixing of air parcels with differing nitrogen oxide signatures can greatly augment overall ozone production (Liu et al., 1987; Lin et al., 1988; Pickering et al., 1989). Since major anthropogenic nitrogen oxide sources are located at the surface (Logan, 1983; Dignon, 1991; Hameed & Dignon, 199 1; Dignon, 1992) and the species are short lived (Penner et al., 1991), significant gradients often occur in the vertical and are disrupted by convection processes. We have also duplicated 1-D models of the effects of convective storms on tropospheric ozone (Pickering et al., 1990; Chatfield & Delany, 1990; Pickering et al., 1993). Our text begins with some background material regarding the ongoing effort to specialize chemistry integrators for the global scale. The need for fast modules is established by placing hypothetical chemistry schemes upon a set of global grids mimicking GCM codes (Elliott et al., 1993b, 1995). Crude operations counts for traditional integrators of the LMBDF type point to options for streamlining, and further, toward massively parallel processing (Elliott & Jones, 1994; Dabdub & Seinfeld, 1994). Parameterizations of the LMBDF strategy are sketched including matrix free and offline calculations (Byrne & Hindmarsh, 1987; Jacob et al., 1989). The fundamen-

tals of the family approach are then discussed as an alternative (Chapman, 1930; Crutzen, 1970, 1971; Turco & Whitten, 1974). Coupling of species into groups leads to a rise in the lowest time constants. The operations counts thus are more favorable. Runs with lengths comparable to the time constants expected for climate change become conceivable (Houghton et al., 1990, 1992). Parallelization of the family chemistry coding may actually render coupled chemistry-climate calculations convenient (Elliott & Jones, 1994). Next, we describe our own response to the challenges posed by the computational intensities of atmospheric kinetics-our suite of family modules (Kao et al., 1990; Zhao ef al., 1993, 1994, Elliott et al., 1993a-c, 1994, 1995; Shen et al., 1994). Key features of the latest chemistry packages we have constructed are: (1) stabilizing family internal implicit integration for partitioning groups of species; and (2) automation of setup for the solutions, which is necessitated by the continual evolution of reaction schemes. A model format section then provides the reader with details of the outer configuration adopted for tests in the nonlinearity runs. Pure chemical issues are dealt with initially. Some of the requirements for a tropospheric kinetics scheme are determined through examination of key organic oxidation sequences, first for the single carbon species methane and carbon monoxide (Liu et al., 1980; Crutzen, 1988), then for ethane (Aiken et al., 1982, 1983; Kasting & Singh, 1986; Thompson & Cicerone, 1986a, b; Elliott et al., 1994) and generic nonmethane hydrocarbons (Liu et al., 1987). As a byproduct of this preparation for delineation of a reaction sequence, we are able to elucidate the role of NO, as a catalyst in generating ozone. In particular, we note that kinetic titration points distinguish ozone removal and production scenarios (Crutzen, 1973, 1988). The analysis also facilitates generalization of ozone production yields from different types of hydrocarbon, after the fashion of Crutzen (1982) or Jacob et al. (1989). The figures obtained are useful as a check of complete model calculations. Furthermore, our choices for surrogate and representative nonmethane hydrocarbons can be justified through a blend of the yield arguments with global scale reactivities (Singh & Zimmerman, 1992). Sources of the reaction schemes and rate constants we select are ultimately listed. Design of a radiation calculation for the wide variety of tropospheric conditions covered was sufficiently difficult that it is treated as a separate section in and of itself. In essence, the incoming direct solar beam is attenuated through pure Beer’s Law absorption by ozone and molecular oxygen, then adjustments are attached for multiple scattering (Luther & Gelinas, 1976; Isaksen et al., 1977; Madronich, 1987). The core of the text is divided into two segments reporting our replication of nonlinearities in the ozone to nitrogen oxide relationship. We start with the O-D runs in the principal tropospheric trace

231

Streamlined family photochemistry module chemical regimes of remote ocean, free troposphere and continental clean and polluted boundary layer. Although a number of groups have explored the nonlinearities in some depth (Chameides & Walker, 1973; Fishman & Crutzen, 1978), particular emphasis is given to recreating the ozone production figures of S. Liu and coworkers (Levy et al., 1985; Liu et al., 1987; Trainer et al., 1987; Liu & Trainer, 1988). Their studies form a nicely unified set across the critical regions of the troposphere. The role of convective vertical mixing in redistributing minor constituents has been investigated from the single tracer viewpoint by Gidel (1983) Isaac and coworkers (Isaac, 1983; Isaac et al., 1983) Chatfield & Crutzen (1984), Liu et al. (1984), Ching and coworkers (Ching & Alkezweeny, 1986; Ching et al., 1988) and Jacob et al. (1987), among others. Studies coupling photochemistry with convection transport models are less common (Costen et al., 1988; Chatfield & Delany, 1990) but an especially thorough collection is that of Pickering et al. (1989, 1990, 1991, 1992a%, 1993). Mirroring the Pickering climatology, generic simulation is conducted of the ozone production enhancements which follow convective storms. Both the local and column chemistry portions of the analysis consist first of subsections giving alterations required to the baseline model to bring it into alignment with Liu et al. or Pickering et al. scenarios, then a brief interpretation of results. Following the main body of results, a benchmarking section confirms the crude operations estimates for the global mesh, and scaling experiments are described which demonstrate complete vectorization. The family package in fact functions at near the maximum possible speed for an online atmospheric integrator. A final discussion section elaborates on the earlier interpretation of our computations. Reasons for the classic nonlinearities are reviewed as they are manifested in our runs. Areas for improvement of the tropospheric models are itemized. For example, lightning and stratospheric nitrogen oxide sources have been deleted here in the interests of simplicity (Chatfield & Delany, 1990), and photolysis rates are not strictly consistent with the cloud structures implied by convective overturn (Jacob & Wofsy, 1988, 1990). We note that as powerful as the convective climatology may be, its limitations still underscore the need for rapid 3-D chemistry calculations; a coupled column cannot simulate the effect of shear, separation and dilution of parcels into the free troposphere. Our fundamental conclusion is reiterated to end. The critical nonlinearity in the tropospheric kinetic system, that connecting NO, as a catalyst for ozone production with organic oxidations, is successfully and efficiently represented.

Zimmerman, 1991; Kanakidou & Crutzen, 1993). The dearth of simulations can be readily understood through some simple estimates of computational intensities. The operations counting performed here is a much condensed version of treatments by Elliott et al. (1993b, c, 1995) and Elliott & Jones (1994). Even a minimalist and purely inorganic stratospheric chemistry model will include dozens of species and on the order of 100 reactions (Crutzen, 1970, 1971; Crutzen et al., 1978; Brasseur & Solomon, 1984). Within the troposphere, hydrocarbon oxidation chains must be added to the fundamental H-N-O-methane system, and the increase in the chemical data set is manifold (Falls & Seinfeld, 1978; Whitten et al., 1980; Liu et al., 1980; Thompson & Cicerone, 1986a; Finlayson-Pitts & Pitts, 1986; Isaksen & Hov, 1987; Kanakidou et al., 1991). Lower dimensionality models of biologically active areas near the earth’s surface may include 1000 reactions (Atkinson & Lloyd, 1984; Atkinson et al., 1989; Donahue & Prinn, 1990; Atherton & Penner, 1990). The progression of chemical complexity moving downward through the atmosphere is illustrated in Fig. 1. The histogram gives rough species counts for typical low dimensionality models in the major altitude ranges, along with estimates of the number of chemical constituents which may in fact exist (e.g. Solomon, 1988; Rodriguez et al., 1989; Singh & Zimmerman, 1992; Elliott et al., 1995). Emphasis has been placed on our own calculations in constructing the graphic (Turco Kc Whitten, 1977; Turco, 1985; Elliott et al., 1993a-c 1995; Jacobson & Turco, 1994), but reference is also made to the scope of models cited. For purposes of operations counting, let us assume that a tropospheric reaction set interrelating about 100 constituents is adopted. The global grids onto which photochemistry will be placed are likely to derive largely from general circulation model coding. Wind fields and advection numerics are often taken directly from GCMs into chemical tracer models (Prather et al., 1987; Jacob et al., 1987; Kao et al., 1992a, b). One of the first steps

Stratosphere

I

$ 2 $ ,

Free Troposphere

_

I

Boundary Layer I

I

100

1000

I 10,000

I 100,000

Species

OPERATIONS COUNTS FOR GLOBAL CHEMISTRY INTEGRATORS

Fig. 1. Estimates of the total number of species in typical

Global 3-D photochemical models are currently a rather rare commodity (Kao et al., 1990; Crutzen &

low dimensionality models of atmospheric photochemistry at different altitudes (short bars), and of the total number of species which may actually be present (long bars).

1

I

238

Scott Elliott et al.

toward earth system modelling will be coupling of atmospheric kinetics models into GCMs to create what might be called chemical climate codes (Dannevik et al., 1993; Elliott et al., 1995). Our group has already been involved in the construction of two models in which full photochemistry interacts with climate through ozone in the stratosphere (Kao et al., 1990; Zhao et al., 1993, 1994; Elliott et al., 1995). GCM and CTM (chemical tracer model) resolutions have been increasing steadily, from on the order of 10” in the horizontal to on the order of 1” (Arakawa, 1966, 1972; Andrews et al., 1983; Pitcher et al., 1983; Arakawa & Suarez, 1983; Malone et al., 1986; Miyahara et al., 1986; Golombek & Prinn, 1986; Prather et al., 1987; Mahlman & Umscheid, 1987; Hack et al., 1992, 1993; Dannevik et al., 1993). At an intermediate value of about 4”, there may be 100,000 individual cells to deal with. Our own preferred version of the photochemical continuity equations has been developed in Elliott et al. (1993a, b, 1995), inclusive of any discussions of alternate integrators, general family numerics and our specialized family routines. The treatment descends mainly from the works of Gear (1971a, b), Chang et al. (1973), Turco & Whitten (1974) Turco & Whitten (1977), Hesstvedt et al. (1978), Wofsy (1978), McRae et al. (1982a, b), Brasseur & Solomon (1984), Kasting et al. (1984), Kasting & Singh (1986) Rood (1987), Byrne & Hindmarsh (1987), Toon et al. (1988), Kaye & Rood (1989) and Brasseur et al. (1990). It is closely related to theory in other publications from this group (Jacobson & Turco, 1994; Jacobson, 1994). The several dozen expressions derived will not be repeated explicitly here. In their simplest form, the photochemical equations can be time split from dynamics and consist mainly of nonlinear production and loss terms (Yanenko, 1971; McRae et al., 1982a, b; Rood, 1987). The method of lines can also be employed in order to preserve the competition between kinetics and transport in finite differencing form (Chang et al., 1973; Wuebbles & Kinnison, 1989). We will restrict our discussion to the operator split mode because it is likely to be a favored choice in global three dimensions (McRae et al., 1982a, b; Chang et al., 1987; Rood, 1987; Kaye & Rood, 1989). Our own strategy has been to follow splitting procedures in Toon et al. (1988). The chemical system is extraordinarily stiff (Lapidus & Seinfeld, 1971; Lapidus et al., 1974; McRae et al., 1982a, b; Brasseur & Solomon, 1984; Finlayson-Pitts & Pitts, 1986). Reactive atoms may survive only nanoseconds, but still play a crucial role in the decay of species with lifetimes of hundreds of years. The linear multistep backward differentiation formulae (LMBDF) are an established general set of techniques for integrating the stiff chemistry (LMBDF; Gear, 1971a, b; Byrne & Hindmarsh, 1987). The reverse Euler constitutes one instructive example (Wofsy, 1978; Logan et al., 1981; Kasting & Singh, 1986; Byrne & Hindmarsh, 1987; Jacob et al., 1989). In their original form, the

LMBDF involved application of the Taylor series approximation in the chemical vector space, and repetitive construction of the Jacobian with Gauss Jordan solution of the resulting linear system. On the order of n3 operations were required per iteration, where n is the number of chemical species (Oran & Boris, 1981; Boris & Winsor, 1982; Boris, 1989). Figure 2 displays the cubic dependence of the algebra as a function of the constituent load. The line with slope three is the lower limit applying to a first step in convergence on solution to the nonlinear system. Actual counts will fall in the shaded area. The figure is similar to, but not identical with, Fig. 15 from Elliott ef al. (1993a). In the time split mode the calculations are replicated in every grid cell. For the sample configuration of 100 constituents and 10’ cells, a single photochemical time step will consume i x n3 x (cell No.) or i x 10” operations, where i is an iteration counter. Step sizes were usually made adjustable in LMBDF solvers in order to preserve certain predetermined error limits, and were much smaller than the value of roughly 1 h often appearing in the global transport operator (Toon et al., 1988; Kao et al., 1990). Normalized to the dynamical time increment, then, the number of operations per step is i x Y x lo”, where r is the ratio of transport to chemical delta t. The factors i. r and their product will all be larger than unity, and possibly by a wide margin. Even on a recent generation Cray vector supercomputer with roughly Gigaflop throughput (Almasi & Gottlieb, 1989; Levesque & Williamson, 1989), global chemistry calculations using a classical LMBDF solver in

E .g

n’

2 8.

0

//

I RREl(

Species

Families

103

’

t

I

I

I

loo

10’

102

I 100

103

I 10’

102

’ I

Fig. 2. Logarithmic operations relationships for several chemical integration strategies. The several scales along the abscissa assume that the reaction scheme contains an average of three kinetic processes per species and that each .^ .^ .

detmed tamely contains three species

Streamlined

family photochemistry

the time split mode would be difficult to advance beyond a few steps. The climate problem, on the other hand, has a time constant of many decades (Houghton et al., 1990, 1992; Dannevik et al., 1993). It seems clear that successful global scale chemistry solvers should strive to reduce the cubic dependency of the algebra, avoid iterations, increase step sizes or adopt some combination of the strategies. Several viable options for meeting these goals are being pursued within the atmospheric chemistry community. Descendants of the early LMBDF integrators are now sparse or matrix free, but preserve the capability to maintain accuracy control (Brown & Hindmarsh, 1986a, b; Byrne & Hindmarsh, 1987; Brown et al., 1989; Brown & Saad, 1990; Jacobson & Turco, 1994). A group based at Harvard University has developed statistics from a large set of reverse Euler calculations (e.g. Wofsy, 1978; Logan et al., 1981) into a functionalized, offline kinetics database relating production and loss for a few major constituents to local concentrations (Jacob et al., 1989; Spivakovsky et al., 1990a, b; Jacob et al., 1993a, b). The technique derives originally from the concepts of Dunker (1986) and Marsden et al. (1987). Although the parameterizations permit integration over 4 h time steps for a small species set, several limitations can be listed. Calculations are restricted to the chemical space defined by the detailed baseline modelling. A corollary is that the calibration problem will be exhaustive over the concentration regimes spanned by the global troposphere. Major changes in mechanism would require recalibration and alteration of reaction schemes has historically been quite common in atmospheric chemistry. In the stratosphere, half a dozen or so new heterogeneous processes discovered since 1985 have been instrumental in explaining the Antarctic ozone hole (Solomon, 1988; Rodriguez et a/., 1989; Hanson & Ravishankara, 199 la, b; Crutzen et al., 1992; Prather, 1992; Abbatt & Mohna, 1992). They are often represented as gas phase bimolecular equivalents in computations (Solomon et al., 1986; McElroy et al., 1986; Jones ef al., 1989; Prather, 1992). In the troposphere, heterogeneous chemistry is only beginning to appear in global models (Lelieveld & Crutzen, 1990, 1991) but will probably be a source of large uncertainties (Pickering et al., 1990, 1991). Rate constants along organic oxidation sequences are often completely unknown (Trainer et al., 1987; Atkinson ef al., 1989; Demore et al., 1990) so that reactions must sometimes be deleted for lack of laboratory kinetics measurements (Thompson & Cicerone, 1986a, b; Cicerone et al., 1991; Elliott et al., 1994). The functionalized chemistry also sacrifices detail because species must be lumped to reduce the total number of statistical variables (Jacob et al., 1989). Once grouped, the constituents cannot be separated again. Individual hydrocarbons, for example, could not be singled out for tracer studies. The proportion of species within an organic classification is constant, and intermediates along hydro-

module

239

carbon oxidation sequences must be set at steady state, which will not always be a realistic approximation. The family approach has proven useful in facilitating large-scale kinetic models as well (Chapman, 1930; Crutzen, 1970, 1971; Crutzen et al., 1978; Garcia & Solomon, 1983; Brasseur & Solomon, 1984) and has been the most popular global integrator in three dimensions (Rood, 1987; Kaye & Rood, 1989; Crutzen & Zimmerman, 1991; Granier & Brasseur, 1991; Kanakidou & Crutzen, 1993). Certain species in the atmospheric chemical system interchange rapidly enough that they can be considered longer lived tracers collectively. Unstable entities are subsumed, and integration can proceed explicitly, with matrices of order n and iterations both being avoided entirely. Time steps are increased by orders of magnitude over early LMBDF. Experience in both low dimensionality (Crutzen et al., 1978; Garcia & Solomon, 1983; Brasseur et al., 1990) and 3-D (Kaye & Rood, 1989; Kao et al., 1990; Crutzen & Zimmerman, 1991; Kanakidou & Crutzen, 1993; Elliott et al., 1993a; Zhao et al., 1993 1994; Elliott et al., 1995) dictates that family coding can be run at increments near or above those in a global transport operator, so that the factor r approaches unity. Some expert judgement is demanded in family definition (Turco & Whitten, 1974; Farrow & Edelson, 1974; Farrow & Graedel. 1977; Graedel, 1977) and verification through established integrators of known accuracy is advisable (Douglass et al., 1989) so that family solvers can be portrayed as complementary to modern LMBDF. However, speed will approach that of parameterized offline integrators (Elliott et al., 1993a) while preserving full kinetic detail. Calibration exercises for our family models against traditional Gear codes are reported in several earlier works (Turco & Whitten, 1974; Elliott et al., 1993a; Jacobson & Turco 1994) and will be presented for the tropospheric chemistry in thesis form (MS). The operations counting procedure can now be repeated for a hypothetical but efficient family routine. Again we are condensing arguments from Elliott et al. (1993b, c, 1995) and Elliott & Jones (1994). Grouping constituents might reduce the number of real tracers to on the order of 30. We presume that there are 300 reactions in the tropospheric mechanism. Their rates can be constructed for about three multiplicative operations per, or a total of 1000. Production and loss may require summing of 30 or so rates per family, for another 1000 operations, this time additions. First power scaling with number of reactions and chemical entities is indicated on Fig. 2 for both processes. Partitioning in the traditional manner at photochemical steady state and with family external terms excluded tends to reduce to linear algebra (Brasseur & Solomon, 1984; Shimazaki, 1985; Turco, 1985; Elliott er al., 1993a, b, 1995). At three species per assemblage, distribution will consume perhaps an additional 30 operations per. or 1000

240

Scott Elliott et al.

total. The system has been linearized at all levels so that iteration is not imperative, although it can be applied to improve accuracy (Elliott et al., 1993a; Jacobson & Turco, 1994). It would be expected that a properly constructed family integrator could then advance global photochemistry calculations at about lo4 x (cell No.) or lo9 operations per transport increment for the mesh defined above. Year long runs would still be difficult on current generation low CPU number vector supercomputers. Massively parallel machines may soon reach the Teraflop barrier, however, rendering decadal scale family kinetics calculations viable (Almasi & Gottlieb, 1989; Dannevik et al., 1993; Elliott & Jones, 1994). AUTOMATED

AND MODULAR

FAMILY

CODES

Our global photochemistry packages achieve the fullest potential of family based numerics through a judicious combination of automation, vectorization and group internal integrations. The individual techniques employed are generally not recent innovations, but the specialized combinations are somewhat unique. The latest coding versions are described in detail in references such as Elliott et al. (1993a-c, 1995), Jacobson & Turco (1994) and Jacobson (1994). We review the major points here. Setup of solutions to the photochemical continuity equations is automated in our fast kinetics modules so that the user need only edit data input lists of species, family definitions and reactions in order to adjust the kinetic mechanism. The automation programming descends from the 1-D stratospheric models of Turco (Turco & Whitten, 1977; Turco, 1985). Full photochemical rate expressions and kinetic partial derivatives are all constructed transparently. Similar data list reading procedures are undoubtedly utilized by other groups, but are usually not mentioned in the atmospheric chemical literature. For a few interesting exceptions to this rule, however, see Ridler (1977) and Ridler et al. (1977). Vectorization is attained by running all cells in the model along a single geographic index, then placing loops over that index innermost whenever possible. The kineticist’s first instinct for programming in one dimension or greater is often to enter a grid location conceptually and perform photochemical calculations there before moving on (Elliott et al., 1993a, b, 1995). Loops over species or reaction may well be longer than a typical machine vector length (Boris & Winsor, 1982; Levesque & Williamson, 1989), but they are sometimes recursive. Species or family specific rates, for example, constitute reduction functions (Elliott & Jones, 1994). The number of cells in multidimensional atmospheric models will also exceed the vector length, and in the time split mode the calculations are completely grid cell internal. The chemistry may be regarded as an infinite medium calculation. The independence of the time split kinetics will render it exceptionally easy to

parallelize (Boris, 1989; Gear, 1988, 1989 1993; Gear & Xuhai, 1993; Elliott&Jones, 1994). Little interprocessor communication will be required during a chemical time step. On recent generation massively parallel computers incorporating virtual processing, distribution across the nodes is handled transparently. Message passing issues should only arise with respect to load balancing (Weber, 1994; Weber et al., 1994), and so will not pertain to a noniterative family algorithm. Parallelization of our family coding is thus a high priority. Definition of families under the standard photochemical equilibrium approximation (Brasseur & Solomon, 1984) is sometimes limited by species which defy grouping because they have multiple sources. The trouble constituents may also be quite short lived, however. In our experience, for example, the hydroxyl and methyl peroxy radicals can restrict stratospheric family integration step sizes to about the minute range (Kao et al., 1990; Elliott et al., 1993a). We have come to favor a family internal integration approach for circumventing such grouping barriers to raising the photochemical time step. The method was outlined as early as 1974 (Turco and Whitten) and has been used in the models of Crutzen et al. (1978) Gidel et al. (1983) and Chang et al. (1987) among others. In our baseline coding we generalize it across all families. Full photochemical continuity equations are established for each species in an assemblage, written in implicit finite differencing form. Family internal nonlinear terms representing production are linearized into a CrankNicholson-like average (Crank & Nicholson, 1947; Turco et al., 1979a, b; Elliott et al., 1995). Removal terms are linearized in the positive definite form (Elliott et al., 1995). The implicit linear system which results is solved over the time step and yields individual species concentrations which are sufficiently stabilizing that the actual integration for the families can be conducted in the forward mode. The particular tropospheric package tested in the present work will be detailed in an upcoming Ph.D. Thesis (MS), and in a paper dealing with isoprene chemistry (Shen et al., 1994). The family partitioning routines are discussed in Jacobson (1994). We present here a brief sketch of the major program segments. Computations begin with calls to subroutines which initialize the chemistry and radiation calculations and which establish the physical grid. Input constituent and reaction lists are read and arrays constructed for the automated rate and partial derivative computations. Pressure, temperatures and coordinates are given. Within the actual chemistry loop, a radiative transfer program provides the actinic radiation flux at the center of any one grid cell based on time, date, location and the columnar fields of major absorbers and scatterers. Photolysis and kinetic rate constants are then set, photorates calculated, and the family internal implicit integrations performed in order to set the stabilizing individual species concentrations.

Streamlined family photochemistry module Geometric averaging of the stabilizing species concentrations (Sokal & Rohlf, 1981) with a set saved from the previous time step may be inserted as an accuracy enhancing option. After total family concentrations are advanced the species concentrations are normalized to them based on the original ratios. MODEL

DESCRIPTION;

CHEMISTRY

A global tropospheric kinetics scheme must contain a core of inorganic data resembling the chemistry of the stratosphere (Crutzen, 1973; Liu, 1977; Liu et al., 1980; Logan et al., 1981; Demore et al., 1990; Crutzen & Zimmerman, 1991). As a matter of conceptual convenience here we will place the crucial single carbon species methane and carbon monoxide, along with their oxidation products, into the inorganic category. Oxidation of the single carbon molecules produces the pollutant and greenhouse gas ozone, as well as the HO., radical oxidants (OH and HO,; Levy, 1971; Crutzen, 1982; Thompson, 1992) if the NO, concentration lies above a local titration point determined by a competition between NO and other sinks for peroxy radicals (Chameides & Walker, 1973; Liu et ai., 1987; Crutzen, 1988). The threshold for the nitrogen oxides is usually in the vicinity of 10 pptv; at lower concentrations removal of ozone and the hydrogen oxides dominates (Liu et al., 1987; Crutzen, 1988). The competition for peroxy species is one of the main causes of nonlinearities in ozone production as a function of nitrogen oxide levels (Liu & Trainer, 1988). Oxidation of multiple carbon molecules may also lead to ozone generation and must be accounted for in a tropospheric reaction listing (Kasting & Singh, 1986; Finlayson-Pitts & Pitts, 1986; Liu et al., 1987; Crutzen, 1988; Kanakidou et al., 1991; Singh & Zimmerman, 1992; Kanakidou & Crutzen, 1993). The multiple carbon organics will be referred to in the conventional jargon as nonmethane hydrocarbons, or NMHC. Simple mass budgeting arguments for their global import in determining oxidant distributions can be found in Liu et al. (1987), Crutzen (1988), and Jacob et al. (1993a, b). As anthropogenic releases of the nitrogen oxides begin to overlap strong tropical inputs of the natural hydrocarbons, the NMHC could rival single carbon species in integrated ozone sourcing (Crutzen, 1988). Major channels are presented in Table 1 for oxidation of the single carbon species both below and above the nitrogen oxide titration point. The mechanism is a composite from Liu (1977), Liu et al. (1980) and Crutzen (1988). Note that the role of the NO,? species in ozone production is catalytic. Reaction of NO with a peroxy radical transforms it into NO,, which photolyzes to atomic oxygen, which in turn combines with molecular oxygen. The photolysis also recreates NO. A maximum molar ozone yield for destruction of 1 mole of methane is five. In practice significantly less ozone results because formaldehyde

241

sometimes photolyzes to yield molecular hydrogen, H,. Ozone production efficiencies for the NMHC rapidly become difficult to bracket as the number of carbon atoms increases. In a seminal exposition of the global effects of organics, Liu et al. (1987) offered the generalized reaction set in Table 2, in which CARB signifies a byproduct carbonyl. The sequence implies a safe minimum of two ozone molecules per organic. However, it is readily demonstrated that the effective yield should scale with the number of carbon atoms (Crutzen, 1982; Jacob et al., 1989). Our analysis will be based on Table 3, an ethane destruction sequence for the high nitrogen oxide regime, but we will attempt to extend to higher alkanes as well as unsaturated species. Estimates of midday hydroxyl radical (OH) concentrations and photolysis rates are those of Trainer et al. (1987) and Jacob et al. (1989). All hydrogen atoms in an ethane molecule are terminal. Through acetaldehyde (C2H40) the total ozone production is two molecules, and resembles the sequence in Table 1 from methane to formaldehyde. This is the general minimum for CARB from Table 2. Photolysis of acetaldehyde gives no ozone directly, but is a minor channel relative to hydroxyl attack on the aldehydic hydrogen, which tends to be rapid (Atkinson & Lloyd, 1984; Thompson & Cicerone, 1986a, b). Hydrogen abstraction leads to a peroxy radical of the form RCOO,, where R signifies an arbitrary organic moiety. In the ethane case the fragment is a methyl radical. The peroxy reacts with NO to give NO, and so also ozone, but a carbon atom is lost to carbon dioxide in the process. The total ozone yield for the first carbon atom oxidized is thus three molecules. The product methyl enters into the single carbon decay of Table 1. The data employed here indicate that hydroxyl attack will constitute the major loss process for formaldehyde. The average ozone yield per carbon atom during the ethane decomposition is then three molecules. For larger alkanes the initial hydrogen abstraction may not be terminal (Atkinson & Lloyd, 1984; Liu et al., 1987; Trainer et al., 1987). The usual first reaction steps give two ozone molecules on the way from an interior hydrogen to a ketone molecule. Ketones photolyze efficiently to the products ROz and RCOOz (Trainer et al., 1987). The simpler of the two peroxy radicals, R02, behaves precisely as if its genesis had been from a terminal alkane. The carbony1 peroxy results in one final ozone molecule from the original site of the hydroxyl attack, so that the total is again three ozone molecules per carbon atom. Our development suggests that alkanes hold the potential to generate three ozone molecules per carbon during oxidation in the presence of the nitrogen oxides. Jacob et al. (1989) indicate that unsaturated species will behave in a similar fashion. These authors report actual yields obtained from an LMBDF box model of decomposition for a number of different organics. The values cluster around the figure of 1.5 ozone per carbon, but reasons for any discrepancy

242

Scott Elliott er al. Table I. Sink

carbon molecule oxidation

channels

at hinh and low nitroaeo

oxide levels

03+hv+O’D+0, O’D+HIO+OH+OH O,+H,O+20H+O, NO High CO+OH-CO,+H H+07+HOI HO,+NO+OH+NO, NO,+hv~NO+O o+o,-+o,

NO Low CO+OH-+C02+H H+O,-+HO, H02+0,+OH+20,

CH, + OH - CH, + H,O CH, + 0, --t CH,O, CH,O>+NO-+CH,O+NO, CH,O+O,+CH,O+ HO, NO,+hv-+NO+O o+o,-+o, HO,+NO-+OH+NO, NO,+hv-+NO+O oi~o,+o,

CH, + OH - CH, + H,O CH, + 0, -+ CH,O, CH,O, + HO, + CH,O,H CH,O,H+hv-+CH,O+OH CH,0+02+CH20+H02

CH, + 40, -

CO+0,+CO,+O*

+ 0,

Competition CH,O,+HOz+CH,OzH+O, CH,02H+OH+CH,02+Hz0 HO,+OH+H,O+O,

CHI + 0, --+ CH,O + H,O

CH,O + H,O + 20,

CH,O + /iv - CHO + H H+O,+HO, CHO+O,dCO+HO, HO,+NO-+OH+NO, HO,+NO-+OH+NO, NO, + hv - NO + 0 NO,+hv+NO+O o+o,-+o, o+o,-+o,

CHIO+hv-+CHO+H H+O,+HO, CHO+Oz+CO+HO, HOz+0,+OH+202 HO,+O,+OH+20, CH,0+203-+CO+20H+20,

CH,0+402+CO+20H+20, CH,O + OH - CHO + H,O CHO+O,dCO+HO, HO,+NO-+OH+NO, NO,+hv+NO+O 0+0,-o,

CH,O + OH --* CHO + H,O CHO+O,+CO+HO, HO,+O,+OH+20, CH,0+03+CO+01+H20

CH,O + 20, + CO + H,O + 0,

relative to our own calculations can readily be identified. Long lived intermediates such as ketones and alkyl nitrates were removed from the box as they appeared, lowering the total ozone production per molecule in the direction of the Table 2 limit. Jacob et al. (1989) also exclude carbon monoxide from the computations, which effectively deletes one ozone from oxidation of any final alkane methyl radicals. Essentially, the crude yields constructed here agree with calculations in Jacob et al. (1989) and Crutzen (1982). Our strategy for implementing a global tropospheric data set has been to take core inorganic Table 2. A generalized nonmethane hydrocarbon oxidation. Only the ozone producing, high NO, situation is represented NMHC + OH +0*-r RO, RO, + NO + 0, - NO, + HO, + CARB HO,+NO+NO,+OH

Net

NMHC

+ 40, + hv -

20, + CARB

reactions from a stratospheric model (Turco, 1985; Demore et al., 1990; Cicerone et al., 1991; Drdla et al., 1993; Elliott et al., 1994) and then to superimpose lower atmospheric inorganics and key NMHC species. Since per carbon ozone yields will be similar from organic to organic, weightings for different species in the global system can be deduced by estimating their initial reactivities. The calculations in Table 4 are based on data provided in Singh and Zimmerman (1992). The first column gives the concentration for selected ambient organics which would yield the methane reactivity relative to the hydroxyl radical. The second column lists order of magnitude estimates of typical surface concentrations near sources. The third column is simply the quotient of the second and first, and is a fractional CH, reactivity equivalent. In the final column the overall fractional equivalents are multiplied by the number of carbon atoms in each molecule for a reflection of relative ozone production capabilities. We have elected to include the two simplest members of the alkane and

Streamlined family photochemistry module Table 3. A scheme for

estimating ozone yields

from alkane oxidation,

HO, -+ 01

CzH,+OH+C,H,+H,O C,H,+ O,-+C,H,O, C,H,O, + NO+ C,H,O CIH,0+02--rC2H,0+

NO, -+ 0,

C2H,0 + OH + CH,CO + H,O CH,CO f O,-, CH,COO, CH,COO, + NOCH, + NO, + CO,

NO, -+ 0, HO, -+ 0,

CH,+OI-rCH,O, CH,O, + NO + CH,O + NO, CH,O + 02--+ CH,O + HO,

(HO,+O,) HO* -+ 0,

CH,O + OH + CHO + H,O CHO+Oz--CO+HO,

HO, -+ 0,

CO+OH+CO,+H H+02+HOz

NO* * 0,

ethane sequence

I750 40 9.6 3.5 10.4 I.0 0.3 7.0 1.3 0.08

+ hv -+ CH, + CO + H (minor)

C,H,O

CH,COO,

+ NO, -

PAN (semi-minor)

(CH,O-+

2H0, + CO)

CH,O-+

H, + CO

from the following sources: ethane (Thompson & Cicerone, 1986a, b; Elliott et al., 1994) propane (Liu et al., 1987; Trainer et al., 1987; Elliott et al., 1994). ethylene and propylene (Trainer et al., 1987; Pickering et al., 1989, 1990, 1991; Chatfield & Delany, 1990). Understanding of branching and rate constants along isoprene oxidation channels continues to evolve rapidly (Zimmerman et al., 1978; Hov et al., 1983; Isaksen & Hov, 1987; Jacob & Wofsy, 1988; Paulson & Seinfeld, 1992). We follow the data compiled by Paulson & Seinfeld (1992) and applied in Shen et al. (1994). Family definitions can be altered easily within the photochemistry package because of the automated setup routines. The family internal integrations which lead to partitioning additionally confer a great deal of flexibility in the groupings selected. For example, species which are tightly bound to one another in photochemical equilibrium during the daytime can be considered a family assemblage but can effectively separate at night (Elliott et al., 1993a, 1995). Definitions used here are similar to those in Elliott et al. (1993a, 1994) for the inorganic species, and for ethane and propane and their oxidation intermediates. Groupings for the other nonmethane hydrocarbons are handled in an analogous fashion. Generally speaking, alkyl radicals are considered with parent hydrocarbons, alkyl peroxy and alkoxy radicals are taken together, hydroperoxide derivatives are treated as species in and of themselves, aldehydes and derivative acetyl radicals are grouped. Peroxyacetyl radicals, PAN (peroxyacetylnitrate) type substances and

of ozone production ppbv for CH, reactivity

Methane Ethane Propane n-Butane Acetylene Ethylene Propylene Benzene Toluene Isoprene

based on a simplified

+ NO, HO2

alkene homologous series here (ethane, propane, ethylene and propylene), along with acetylene and isoprene. Table 4 indicates that higher alkanes and perhaps toluene might also be appropriate choices. Our six NMHC happen to be the species most commonly used as surrogates for other organics in models of the global ozone system (Liu et al., 1987; Isaksen & Hov, 1987; Chatfield & Delany, 1990; Pickering et al., 1989, 1990, 1991). The foundation stratospheric reaction set is taken from ozone hole publications (Cicerone et al., 1991; Drdla et al., 1993; Elliott et al., 1994) and so includes not only standard chlorine oxide chemistry, but also more exotic species such as OClO, the chloroperoxy radical and the Cl0 dimer, along with bimolecular equivalent representations of the major sulfate particle and polar stratospheric cloud heterogeneous reactions. Full bromine oxide kinetics are enfolded as well. The chlorine and bromine schemes may eventually come into play in our simulations; recent studies suggest that the significance of halogen oxides in lower tropospheric photochemistry has been underestimated (Fan & Jacob, 1992; Keene et al., 1993; Pszenny et al., 1993; Finlayson-Pitts, 1993). Rate constants for most of the inorganic kinetics processes in the model can be found in Demore et al. (1990). Our single carbon species reaction list is quite similar to those in classic studies such as Liu et al. (1987) and many of the crucial reactions are itemized in Table 1. Rate constants have been taken from Atkinson & Lloyd (1984) and Atkinson et al. (1989). Nonmethane hydrocarbon oxidation sequences are adopted Table 4. An estimate

243

capability

GlObal actual I750 2 1.0 0.1 0.1 0.1 0.1 0.1 0.1 0.1

for various

hydrocarbons

Fraction of CH, equivalent

x#C

I

I

0.05 0.10 0.03 0.01 0.10 0.33 0.01 0.10 1.0

0.10 0.30 0.12 0.02 0.2 I .o 0.06 0.7 5.0

244

Scott Elliott et al.

organic nitrate molecules are individual species. Overall the kinetic reaction scheme contains 92 species collected into 40 families interrelated through 275 kinetic and photolytic processes. Detailed reaction listings will appear in an upcoming Ph.D. Thesis (MS); results of the nonlinearity experiments will be emphasized here. The computations in essence constitute accuracy tests, and the satisfactory outcome should provide at least some indication that the kinetics have been configured appropriately. Two peripheral points will be mentioned in this section before moving on to a description of radiative transfer in the coding. We have found it expedient for our purposes to calculate ozone yields from the perspective of hydrocarbon oxidation. There is an alternate viewpoint with some advantages. The yields can also be expressed in terms of nitrogen oxide emissions. Derivations for the technique can be found in Liu et al. (1987) and Lin et al. (1988). It has been applied to continental scale ozone budgeting several times (McKeen et al., 199 1a, b; Jacob et al., 1993a, b) and to the global scale by the originators (Liu et al., 1987). An exposition of the differing methods for estimation of ozone production is offered in the discussion below. We also note that reaction schemes such as the ones in Tables 1-3 can be simplified through judicious redefinitions of the photochemical family into which ozone is placed. In stratospheric models, an odd oxygen group is often composed simply of ozone and the several states of atomic oxygen (Chapman, 1930; Brasseur & Solomon, 1984; Turco, 1985). We have sometimes adopted this definition in our tropospheric modelling. Conceptually, however, it can be beneficial to include other oxygenated substances. Inclusion of the NO, species NO?, for example, eliminates the need for inclusion of its photolysis reaction in the Tables. Attack of NO upon the peroxy radicals becomes the rate determining step in odd oxygen production. If NO, is defined as an odd oxygen member, other nitrogen molecules must be included as well to prevent unwanted losses. Liu and coworkers (Liu, 1977; Liu et al., 1980) add NO,, N,O, and nitric acid to the odd oxygen family. The key rate determining steps in ozone production may then be generalized as X0* + NO, and the key removals are the reaction of the excited OlD atom with water vapor, and of ozone itself with OH and HO, (Liu, 1977; Liu et al., 1980; Pickering et al., 1989). The odd oxygen family has on occasion been extended even further, to embrace HO, species and hydrogen oxide reservoirs (Levy et al., 1985). For example, if the hydroxyl radical OH is defined as half a unit of odd oxygen, then OlD + H,O + 20H is a net null process. Reservoirs such as hydrogen peroxide are included. This may be conceptually convenient under some circumstances, but the grouping becomes somewhat cumbersome and artificial. The odd oxygen photochemical lifetime is increased because removal requires more kinetic intermediates. However, the hydrogen oxides and their reservoirs do

not necessarily share radiative and oxidant properties with ozone, so that the justification for collectivizing becomes obscure. It could be argued that smaller families are more realistic. MODEL DESCRIPTION;

RADIATION

The actinic or spherically integrated fluxes of radiation which must be used in photodissociation rate calculations can be estimated at several levels of sophistication. Early studies of atmospheric radiative transfer were often plane parallel, with isotropic or Eddington diffuse fields (Chandrasekahr, 1960; McElroy & Hunten, 1966; Prather, 1974). The first atmospheric chemical models were operationally intensive relative to the available computers even in low dimensionality (Chang et al., 1973; Luther, 1980) so that only the direct solar beam was treated and it was attenuated only through absorption by oxygen species (Luther & Gelinas, 1976; Isaksen et al., 1977). Parameterizations of diffuse effects appeared in the mid-seventies and later, and have been applied continually since (Luther & Wuebbles, 1976; Luther & Gelinas, 1976; Isaksen et al., 1977; Luther, 1980; Nicolet et al., 1982; Garcia & Solomon, 1983; Solomon & Garcia, 1983; Bruhl & Crutzen, 1989; Elliott et al., 1995). Adjustment factors were often tabulated relative to absorption only values, for multiple scattering molecular atmospheres as a function of surface albedo. The effects of clouds and aerosols were sometimes ignored for simplicity (Luther & Gelinas, 1976), and it has been argued that particles will influence the actinic radiation field most significantly in polluted areas below a few kilometers altitude (Demerjian et al., 1980; Logan et al., 1981; Thompson, 1984; Madronich, 1987). Fully spherical radiative transfer calculations have appeared in order to treat twilight problems (Anderson, 1983) and become crucial with respect to the Antarctic ozone issue because zenith angles are large throughout the springtime photolytic period at high latitudes (Anderson & Lloyd, 1990; Elliott et al., 1994). Other important works we have reviewed but which are not included in the above discussion are listed in Logan et al. (1981), and in the introduction to Nicolet et al. (1982). A particularly useful recent reference has been Madronich (1987), which distinguishes carefully between the actinic fluxes calculated under the assumptions of isotropic and collimated diffuse fields. The potential exists for errors of the order of a factor of two if collimation is taken for granted (McElroy & Hunten, 1966; Isaksen er al., 1977; Luther, 1980; Thompson, 1984; Madronich, 1987). Many of our calculations will focus on remote areas of the troposphere, and we have emphasized computationally intensive 3-D problems as a rule (Kao et al., 1990; Elliott et al., 1993a), so that it is to our advantage to adopt a fast, elementary approach to the calculation of photolysis rates initially. We begin by focussing on cloud free conditions in a clean

Streamlined family photochemistry module where aerosols are not a factor influencing actinic fluxes. This is clearly an idealized situation and we intend to perturb it in the sensitivity test fashion. Ultimately we will incorporate directly calculated radiative transfer. We compute attenuation of the direct beam in the pure ozone absorption mode, then apply the wavelength dependent and zenith angle dependent adjustments from Luther & Gelinas (1976) (see also Luther & Wuebbles, 1976; Isaksen et al., 1977; Luther 1980; Madronich, 1987) to reconstruct a diffuse, multiple molecular scattering field. Tabulations in Nicolet et al. (1982) show similar results. The effects on the ratios of the direct beam, absorption only flux to the diffuse flux are described with respect to simplified collimated two stream approximations in Isaksen et al. (1977). Interpretation in other works is closely related (Luther & Wuebbles, 1976; Luther & Gelinas, 1976; Nicolet et al., 1982). The accuracy of the multiple scattering adjustments can be verified under suitable conditions through ozone or NO2 actinometry (Dickerson, 1980; Dickerson et al., 1980; Parrish et al., 1983). Photoelectric detectors have been used for flux measurements as well (Frederick & Lubin, 1991; Arellano et al., 1994). The tabled adjustments naturally apply strictly only to selected typical ozone fields and not to the specific profiles used in our models. This offline approach was common in 1-D models of atmospheric chemistry (Luther & Wuebbles, 1976; Luther, 1980). Online radiative transfer calculations would be expected to improve our simulations but are not currently feasible. We commonly perform sensitivity tests varying the adjustments on the order of 10% with little effect on overall results. The influence of clouds on the photolytic radiation field is particularly difficult to model. Madronich (1987) has summarized the problems involved. Clouds vary widely in albedo, vertical and horizontal extent, and overlap. Furthermore, their particles usually fall in the regime for which droplet radius exceeds wavelength, so that scattering has a complex forward lobe. Early treatments of clouds in photochemical calculations were highly idealized. Chemists studying in-cloud or in-drop processes tended to seal? clear sky actinic fluxes (Chameides & Davis, 1982; Graedel & Goldberg, 1983; Chameides, 1984; Graedel et al., 1986). Gas phase kinetics modellers used global and temporally averaged cloud properties in two stream models (Logan et al., 1981; Briihl & Crutzen, 1989). Clouds were smeared to reproduce surface or satellite observations of coverages and overall albedo (London, 1952), and were stacked on top of one another. Sometimes they were even combined with the surface (Luther, 1980). Thompson (1984) tested the smearing technique against variable and realistic fields in 1-D photochemical models and suggested that effects on tropospheric ozone were not significant, but she does not give details. Results for individual cloud types may be more extreme (Madronich, l987), but averaged computations from location

245

two stream radiative transfer codes usually indicate above cloud increases and below cloud decreases in actinic fluxes both on the order of tens of percent (Logan et al., 1981; Thompson, 1984, 1992 Madronich, 1987). Our approach here is to assume clear sky conditions whenever possible, thus circumventing the cloud questions entirely in many instances. The special case of chemistry coupled to convective storms is discussed later in its own section. In large-scale models we plan to treat clouds as parameters varying in height and in albedo. We may multiply (above cloud) or divide (below cloud) actinic fluxes by some arbitrary factor ranging from one to two. The method is crude but justifiable as a start up expedient in global 3-D work. Eventually, interactive photochemical climate models must conduct cloud kinetics computations online. LOCAL NONLINEAR OZONE PRODUCTION

In this section we test the ability of our organic containing kinetics module to reproduce the nonlinearities which have been observed in ozone production as a function of NO, level, and also quantitative production rates, across the global troposphere. We begin with remote free tropospheric runs. Comparisons will be made against the model results of Levy et al. (1985) for generalized air parcels at the 500 mbar level. The Levy et al. calculations are descendants of earlier works from the Liu group (Liu, 1977; Liu et al., 1980). The Liu photochemistry has been verified against measurements during open ocean field campaigns (Liu et al., 1983, 1992). Methane and carbon monoxide oxidations in Table 1 are the chief ozone generating processes (Liu et al., 1980; Crutzen, 1988). Concentrations for the major species carbon monoxide, water vapor and ozone are set as in Levy et al. (1985) Table 2. Ozone column integrals and local temperature are derived from the same source. Concentrations of the nitrogen oxides are varied over the range in Levy et al. Levels for nitrogen reservoirs such as nitric acid and PAN are taken from the available measurement and modelling literatures (e.g. Logan, 1983; Penner et al., 1991; Kasibhatla er al., 1993). Since the reservoirs are thus not strictly internally consistent with the computations which follow, ratios of species within the nitrogen oxide groupings sometimes changed during the runs by as much as several tens of percent. Integration was performed for 1 day at 1 h time steps, beginning at midnight. The overall ozone change is calculated for the day and then plotted at half log units in Fig. 3. The pattern of nonlinearity is similar to Levy et al. (1985) for all four seasons. We have simulated the solstices and equinoxes, but display only summer and winter results. The effects are somewhat milder in colder seasons. Levy et al. excluded methane oxidation from their equations in the interests of simplicity. Our goal is to develop a detailed photochemical modelling

246

**l#l!d Scott Elliott et al.

5.0

/

/

,’

I

% %

8

I

,I/

I’ / /

/

/

/ / /

u.0

I

100

-_

I

/

Summer

I

//

200

300

Winter

100

400

200

(x 10)

300

400

NOx (PP~V)

Fig. 3. Changes in ozone as a function of NO, over a 1 day integration period, with the fast chemistry package placed at 500mbar within the free troposphere.

capability, and consequently methane and its decomposition products are present. Through the sort

of yield arguments presented earlier and by examination of Table 1, it would be expected that methane oxidation would give several times as many ozone molecules as the single one obtained from carbon monoxide. Following this point to its logical conclusion would lead through reactivity calculations (Crutzen, 1982; Singh & Zimmerman, 1992) and it can be estimated that methane oxidation channels should raise the Levy et al. ozone production curves about 50%. This would bring our results into close agreement. We have run with methane zeroed and come within about 10% of the Levy ef al. values. The comparison of our own remote tropospheric runs with Levy et al. constitutes an effective accuracy test in addition to demonstrating our ability to pick up the key nonlinearities. The coding has been given only initial concentrations for the long lived species. It then constructs the correct levels for radical families such as HO, which control oxidation rates. Explicit species by species accuracy verification has also been conducted, against classical Gear code type integrators. It will be reported in an upcoming Ph.D. Thesis (MS). We elect to plot our results with methane oxidation included as opposed to making the more direct comparison with Levy et al. because of our emphasis on detail. We have also performed the nonlinearity tests for the continental troposphere. Liu et al. (1987) and Trainer et al. (I 987) combined several data sets from the Niwot Ridge site in the Colorado Rockies (Greenberg & Zimmerman, 1984; Mihelcic et al., 1985; Parrish et al., 1986a. b) to validate their model for the continental

regime, which is sometimes

rich in organ-

its. Niwot Ridge is an ideal location for photochemistry investigations because it receives both clean air from the mountains and polluted parcels from the Denver urban area depending on the direction of the prevailing wind. We have repeated some of the Liu

et al. (1987) and Trainer et al. (1987) experiments with our optimized package. Following the two baseline works, we position the model at 40”N latitude and 700 mbar (about 3 km) The temperature is set at 295 K. Runs begin at midnight on 21 July, and 313 Dobson units of ozone are placed overhead. All simulations are for clear skies. Initial concentration sets are mainly identical to Liu et al. (1987) and are reiterated in Trainer et al. (1987). Conditions are moist, with total water vapor at 1%. Ozone is initialized at 40 ppbv. Carbon monoxide is set proportional to the nitrogen oxides and ranges from 200 to 1OOOppbv. Nitric acid is initialized at one third NO,?. Anthropogenic NMHC ratios are taken from Greenberg & Zimmerman (1984), adjusted upward somewhat to account for organics with four carbon atoms or greater which were not reflected in the measurements. The standard values are ethane 2.5 ppbv, ethylene 0.5 ppbv, propane 1.5 ppbv and propylene 0.2 ppbv. Our model does not contain butane, which Liu et al. place at 1.0 ppbv in the standard set. We have used propane as a surrogate and added an extra amount corresponding to the 1.Oppbv butane reactivity, based on ratios of initial hydroxyl radical attack rate constants (Atkinson & Lloyd, 1984; Atkinson et aI., 1989). All anthropogenic organics are scaled to NO, with the above concentrations applying at 0.8 ppbv; polluted air masses should contain both groups of molecules. Liu et al. acknowledge that the ratios of hydrocarbons will change with time and distance from a source, but this subtlety influences only agreement between model and measurement, not intermodel cross checking. Greenberg and Zimmerman observed about 0.65 ppbv of isoprene and 0.35 ppbv of other terpenoids at the 1 m level at Niwot. These values comprised the bulk of the natural hydrocarbon loading. Based on ratios of surface isoprene to boundary layer averages in models such as Hov et al. (1983) and Trainer et al. (1987) which resolve the lowest few hundred meters, Liu et al. (1987) estimate a boundary layer average of about 0.1 ppbv for isoprene. The terpenoids are represented using isoprene as a surrogate and its total concentration is held in the model at 0.15 ppbv by adjusting an input flux beginning at sunrise. Integration was conducted beginning at midnight and again employed 1 h steps. Liu et al. integrated for 5 days and included boundary layer dilution processes. Ozone changes over the day long runs are shown in Fig. 4, along with the corresponding values from Liu et al. (1987). The ozone increases reported are calculated from morning concentrations, which differed little from initial levels, and from an average over the late afternoon. Our averaging period is that of Trainer et al. (1987). The lowest curve in the figure is in fact an extension of the 500mbar runs from the Levy et al. (1985) comparisons, with water vapor at one part per thousand. The first run conducted under the Niwot conditions excludes all NMHC, so that ozone production

Streamlined family photochemistry module

L

I

I

1.0

I

2.0

NOx(PPW

Fig. 4. Changes in ozone as a function of NO, over a

daytime integration period, with the fast chemistry package placed at the Niwot Ridge experimental site in the Colorado Rocky Mountains.

results solely from carbon monoxide and methane oxidation. Values are a factor of three higher than the Fig. 3 results reprised in the dotted curve because of increased irradiation due to lower ozone column abundances, increased OH due to moist conditions, and other factors. Satisfactory agreement is again obtained with more traditional models. The second set of Niwot experiments adds the isoprene fluxes. Agreement with Liu et al. (1987) is good but may be somewhat fortuitous. Liu et al. do not report ozone production for isoprene alone. We have constructed the increases by adding to the Liu et al. CO, CH, results enhancements due to isoprene from terpene only runs in Trainer et al. (1987). The latter cannot be used directly for our comparison because they are reported in terms of absolute concentration, not differences and morning values cannot be read from the data. Nevertheless the proximity of the curves is encouraging. It seems likely that our treatment of isoprene oxidation is adequate. When runs are conducted with all the anthropogenic organics as well as isoprene, ozone production during the integration period exceeds 45 ppbv but is about 20% below the Liu et al. values. There are several possible explanations for this discrepancy. The longer integration time in Liu et al. may have led to a buildup of oxidation intermediates from day to day, ultimately augmenting ozone production. Accumulation would be minimal over a single day. In fact the short run length is probably more realistic. Liu et al. state that parcel travel time from Denver is around 1 day. Polluted air masses are thus expected to proceed beyond Niwot for multiple day processing. Our own curve in fact falls closer to the Niwot delta ozone measurements (Liu et al., 1987; Trainer et al., 1987). Differences in mechanism could conceivably account for some of the divergence as well. Remember, for

241

example, that we have adopted propane as a butane surrogate, so that the chemical schemes are not identical. Agreement in the isoprene case and for methane/CO suggests that numerical problems are not responsible. We have compared concentrations for the key chain initiators hydroxyl and hydroperoxy radicals wherever possible with the earlier computations. Data presented in Trainer et al. were especially convenient for this purpose. In all cases our noon OH concentrations fell within 30% of their older counterparts, and usually within 20%. Considering the potential for differences in radiative transfer routines, this would seem to be acceptable. Patterns of hydrogen oxide radical concentrations as functions of the nitrogen oxides are similar to Liu et al. (1992) for NO, under IOOpptv and to Trainer et al. (1987) otherwise. CONVECTION AND NONLINEAR PRODUCTION

Rapid convective mixing is an especially critical component of the tropospheric chemical transport system. It is subgrid in scale in regional, continental or global transport models (Byers, 1959; Haltiner & Williams, 1980; Ching & Alkezweeny, 1986; Prather et al., 1987; Jacob et aI., 1993a, b) and so is difficult to represent numerically. However, it is responsible for lifting pollutants away from the surface, thus isolating them from any deposition losses (Jacob et al., 1993a, b; Feichter & Crutzen, 1990) and also for placing them into strong upper level wind systems where they can be moved long distances (Liu et al., 1984; Chatfield & Crutzen, 1984; Pickering et al., 1989, 1992a, c; Feichter & Crutzen, 1990). Convective uplift is a process which often entails considerable vertical mixing in the free troposphere and necessarily involves changes in the photochemical environment of pollutants (Gidel, 1983; Chatfield & Crutzen, 1984). It has been identified as a source of major alterations in ozone production potential of air masses stemming from nonlinearity phenomena similar to those in the Liu et al. works (Pickering et al., 1989, 1990, 1991, 1992a-c; Chatfield & Delany, 1990). As another test of our tropospheric photochemistry code, we will attempt to reproduce ozone production rates which have been modelled following convective events. A brief overview of convection models developed in the fields of global climate and dynamics is in order to provide the reader with proper background material. In GCM codes or 3-D tracer transport calculations derived from them, the simulation of convective vertical transport may be quite crude. Some GCM programs employ a simple energy conserving adjustment to restore stability to air columns with appropriate potential temperature profiles (Manabe et al., 1965; Mahlman & Moxim, 1978; Malone et al., 1986; for chemical tracer transport applications see Levy er al., 1985; Kasibhatla et al.,

248

Scott

Elliott et

al.

Table 5. A selection of convective tracer redistribution events which have been reported in the atmospheric chemical literature. The storms are arranged in order of descending maximun cloud height. Author codes are P = Pickering e/ al., CD = Chatfield and Delany, FC = Feichter and Crutzen, N = Niewiadomski, C = Ching et

ai. Location Amazon Basin Brazilian Cerrado Oklahoma City, U.S.A. Indonesia Amazon Basin Ohio, U.S.A. Montana. U.S.A. Virginia, U.S.A.

Surface

Cloud tops (km)

Land Land Land Sea Land Land Land Land

I5 I4 13 I3 I2 8 7 3

1993). Feichter & Crutzen (1990) note that this manipulation does not allow for the movement of any quantities other than moist static energy (Haltiner & Williams, 1980). In a tracer code based on the GISS general circulation model, a semiarbitrary fraction of one-half of the lowest layer involved in an event is lifted to the highest layer (Hansen et al., 1983; Prather et al., 1987; Jacob et al., 1987) with corresponding subsidence through the column. The method can clearly transport trace chemical species, but it ignores variability in cloud base fluxes and in the locality of detrainment. The GISS model does capture the major vertical features of large storms. Both Walton et al. (1988) and Penner et al. (1991) have commented that convection represented by vertical fluxes can only simulate communication between adjacent layers. In studies with an air chemistry orientation, the effects of individual cloud systems on trace species have sometimes been analyzed in detail. Both measurements and simulations of constituent fields have been performed. Gidel (1983), Chatfield & Crutzen (1984) and Liu et al. (1984) review the earliest works on convective tracer redistribution and apply them in a preliminary sense to issues of global photochemistry. Isaac (1983) [see also Isaac et al. (1983)], Ching and coworkers (Ching & Alkezweeny, 1986; Ching et al., 1988) and Niewiadomski (1986) discuss convective pumping of pollutants in a regional context. The most extensive collection of tracer redistribution data has been compiled and interpreted by Pickering and coworkers (Pickering ef af., 1989, 1990, 1991, 1992a-c, 1993). We have selected from the above works some representative events and listed their main properties in Table 5. The viewpoint is quite chemocentric. We attempt to reduce the complex three-dimensional mixing in storm systems to a shuffling of air masses within a column. This is obviously physically unrealistic, but is useful at the very least for pedagogical purposes. The method may also be useful in practice because it bears resemblance to parameterizations in some major global tracer transport codes (Hansen ef al., 1983; Prather et al., 1987). A similar approach has been outlined by Isaac and coworkers (Isaac 1983; Isaac et al., 1983). The fraction F1 in the table is the portion of boundary layer mixed into the free troposphere during each event, and F2 is the fraction of the air

r Fz

Reference

0.5.0.0 0.33, I .o I .o, 0.0 1.0,o.o 0.25, I.0 0.5,o.o 0.5.0.0 0.1, 1.0

P (1992) CD (1990) P (1992) P (1993) P (1991) FC (1990) N (19X61 \~---r C (1988)

6

originally found in the boundary layer which remains in a well-defined layer at the storm top. Our climatology is constructed crudely on the basis of a set of visual inspections of results which contain considerably more detail, but certain gross conclusions can be drawn nonetheless. It appears that roughly half of the boundary layer is typically transported upward during the observed storms and that detrainment is variable but often spread from cloud base to top. Mixing during convection is intuitively expected; storms generate a considerable amount of turbulence (Pickering et al., 1989, 1990, 1991). The selection is weighted toward land surfaces because of the relative ease of observation and toward deep events because they lead to greater vertical movement. Feichter & Crutzen (1990) have listed global average values for storm parameters and for convective overturn based on the NCAR cloud atlases (Hahn et al., 1982; Warren et al., 1986) and precipitation/humidity data. The results indicate that despite the limited scope of Table 5, the values are somewhat representative. It is interesting to note that the GISS general circulation model convection statistics implied in the Prather et al. (1987) and Jacob et al. (1987) chemical tracer model experiments approach the Feichter and Crutzen averages. Furthermore, it is surely not coincidental that the GISS value of 50% loss from the lowest level matches the climatology in the table. We have included only storms originating from within the boundary layer. Pickering et al. (1990) have argued that while convection is common at higher altitudes it yields little tracer redistribution because the major gradients occur between the boundary layer and the free troposphere. Pickering et al. have not only studied a large number of events but have also run photochemical models of undisturbed and convectively perturbed air columns to determine the effects on ozone production. They stress the conjunction of Liu et al. type nonlinearity issues with storm redistribution of trace constituents. Our strategy for testing our kinetics code in convection situations will be to reproduce the Pickering results in an overall sense. We adopt altitude as the vertical coordinate in a column model, as in Pickering et al. (1990, 1991). Ten one kilometer layers are simulated with the top of the boundary layer set at a crude average of 2 km taken from the Pickering climatology. Initial pressure, temperature

Streamlined Table 6. Initial concentrations

family photochemistry

module

249

of key trace constituents in the tropospheric ozone system for several combinations of boundary laver and free troposoheric orotiles

Species

Units

Low BL

Low free

0, co NO, -4 C,H, C,H, GH,

ppb ppb PPt PPm ppb ppb ppb

25 15 10 1.6 2.0

25 15 10 1.6 2.0

and water vapor profiles are averages from the Pickering mid-latitude works. The water mole fractions derived are typical for low to mid-latitude summer (Mastenbrook, 1966, 1968; Logan et al., 1981). Initial concentrations for key trace species in the tropospheric ozone system are given for several generalized vertical concentration scenarios in Table 6. The word “free” signifies the free troposphere, and the acronym BL stands for boundary layer. Ozone and carbon monoxide levels are fairly consistent throughout the Pickering et al. climatologies at near the values in the table. The nitrogen oxides are relatively short lived and have a patchy input distribution, so that observed and modelled concentrations are geographically quite variable. The figure 10 pptv is typical of very clean conditions such as those found over the Amazon basin during the wet season (Pickering et al., 1991, 1992~; Jacob & Wofsy, 1988, 1990) and also over the remote ocean (Liu et al., 1983, 1992). Urban plumes often contain parts per billion of NO,. Pickering cites observations on this order made over Manaus (Brazil) and Oklahoma City on the U.S. Great Plains (Pickering et al., 1992~). High concentrations for the free troposphere are set at about one-third of the ppbv plume level, but not entirely arbitrarily. The ratio of three is noted by Pickering et al. in both Oklahoma events they have studied (1990, 1992~). Ratios of NOY to NO, and partitioning of NO, into the major reservoirs including nitric acid and peroxyacetyl nitrate are taken from measurements reported in the Pickering studies and are consistent with recent observations and global models (Kondo et al., 1987; Drummond et al., 1988; Penner et al., 1991; Kasibhatla et al., 1991, 1993). Our results are largely insensitive to refinements in NO, partitioning because longer lived species are essentially decoupled from the single day runs. Absolute NO, levels, however, are of course crucial. Pickering et al. have tended to adopt the organic surrogate technique for representation of their hydrocarbons. All alkanes have generally been lumped on a reactivity basis as propane and all alkenes as propylene. Ethane is treated as a separate tracer because of its long lifetime and correspondingly even concentration pattern (Blake & Rowland, 1986; Singh & Zimmerman, 1992). Isoprene is not included in the convection simulations here, which simply means we do not choose to model events occurring over forested areas (Chatfield & Delany, 1990;

High BL

Low free

50 150 1000 1.6 2.0 10.0 0.5

25 15 10 1.6 2.0 0.2 0.1

High BL 50 150 1000 1.6 2.0

10.0 0.5

High free 50 150 300 1.6 2.0 2.0 0.1

Pickering et al., 1992a). The alkane and alkene profiles are derived from gas chromatographic measurements made over the Great Plains (Greenberg & Zimmerman, 1984; Pickering et al., 1990; Singh & Zimmerman, 1992). The concentrations listed in the far right-hand columns of Table 6 are taken directly from Pickering et al. (1990). The lower alkane values for the free troposphere in the central data columns are generated by reducing the Pickering et al. levels by an order of magnitude. This yields concentrations above the boundary layer which roughly mimic global average measurements (Rudolph et al., 1982, 1984a, b; Tille et al., 1985; Kasting & Singh, 1986; Rudolph, 1988; Singh et al., 1988; Singh & Zimmerman, 1992). In the scenario in the table with low NO, concentrations for both the boundary layer and free troposphere, we leave out all organics other than carbon monoxide, methane and ethane. The associated calculations are thus more representative of the remote ocean (Pickering et al., 1991, 1993) than any clean continental site (Pickering et al., 1991). From Table 5 we define two extreme types of mixing within the convective events. They are illustrated here in Fig. 5. It is assumed in a first set of runs that air passing from the boundary layer to the free troposphere is mixed from the top of the new boundary layer height to the maximum cloud top height. Changes in the depth of the boundary layer can result

I

Mixing

TOWar I

---I I I

Fig. 5. A schematic of empirical mixing parameterizations for the convection included in the columnar photochemical model.

Scott Elliott er al.

250

Table 7. A matrix of runs performed to deduce alterations in ozone production following convective events. We adopt a well mixed post cloud free troposphere as baseline. In the case for which concentration gradients are largest we run a low mixing sensitivity simu-

lation Processing

Low-low

High-low

High-high

Mixed Hot tower

x

x x

x

(Pickering et al., 1991) but are ignored here in the interests of simplicity. As a sensitivity test, one run is repeated assuming hot tower type conditions (Riehl 8t Malkus, 1958; Riehl & Simpson, 1979) in which the lifted air mass rises directly to the cloud tops with little turbulent mixing. The matrix of computations is displayed in Table 7. Following the storm redistributions water vapor concentrations often surpassed saturation in the upper few kilometers. Where supersaturation was encountered, water levels were reduced to the vapor pressure. Net precipitation efficiencies generated in this rough fashion match large scale averages from several sources (Braham, 1952; Austin & Houze, 1973; Cheng & Houze, 1979; Feichter & Crutzen, 1990; Chatfield & Delany, 1990). Feichter & Crutzen (1990) discuss uncertainties in the efficiency estimates at some length. Our final water vapor profiles agree with those offered by Pickering et al. in their collection of convective events. Full mass conservation checks were performed post-storm on distributions of all species. All runs begin at midnight on the day before the storm is to take place. The date is set at the summer solstice and the latitude at 45”N. Some of the Pickering et al. work has focussed on tropical systems, but our goal is to model global ozone, so that we select an intermediate latitude and the season of maximum photochemical activity. The step size is 1 h in most of the convection simulations. The time of day during which storms form varies with location and season, and is somewhat earlier over the ocean (Hahn ct al., 1982; Warren et al., 1986; Feichter & Crutzen, 1990). We assume that clouds appear in the 6 h preceding local noon and that they dissipate quickly thereafter. The timing is consistent with much of the Pickering et al. climatology. Photolysis rates are lowered during the morning by a factor of two to simulate cloud effects on photochemical intensities. The 50% reduction is rough but consistent with generalizations of cloud influences presented in Logan et al. (1981), Thompson (1984), Madronich (1987) and Thompson (1992). No wet chemistry is included during the cloud period. Pickering et al. (1991) simulated heterogeneous processes during the morning of their runs by inserting first-order loss rates for solubles. Sensitivity tests indicated that the effects were small in their case study, although Lelieveld & Crutzen (1990, 1991) have recently shown that dissolution of formaldehyde may reduce tropospheric HO,Ylevels if considered on a global average basis. The convective redistribution is imposed on the column instantaneously at noon. In real storms pumping lasts for several hours during the from storms

cloud period (Feichter & Crutzen, 1990; Pickering et al., 1989, 1990, 1991) but the Pickering et al. models have demonstrated reduced kinetic activity during the overturn so that instantaneous transport is adequate for our purposes. Full photolysis rates are restored for the remainder of the run. Changes in ozone concentrations are monitored for 24 h after redistribution, then are totalled and reported as daily values. We now analyze results from our interactive photochemical and convective column models. Ozone production values given the low boundary layer and low free tropospheric concentration data are plotted in Fig. 6. The solid curve indicates the baseline, nonstorm run. The vertical bars and darker line are corresponding profiles at the closest available NO, levels from the works of Pickering et al. The sequence at 5-10 km is from the STEP project south of New Guinea (Pickering et al., 1993). The total prestorm NO, concentration was about 50 pptv over the range. The vertical bars are from Pickering et al. (1991), a study of the early dry season in Amazonia. Nitrogen oxide levels were low because biomass burning had yet to commence after the wet season. The difference between our undisturbed-run ozone generation and Pickering et al. can be explained by our lower NO, levels. Note the rapid change of ozone production rates over the 10-50 pptv regime in the Levy et ai. comparisons (Fig. 3). The integrated daily ozone gain of 0.5-2.0 ppbv is also consistent with the Levy et al. computations. The Pickering et al. profiles in the figure are partial either because NO, was dissimilar at other altitudes or because data simply were not presented. Since our runs contain prescribed NO, in the form of HNO, and PAN, their concentrations are not necessarily at photochemical equilibrium. Thermal dissociation in the relatively warm boundary layer produced some extra NO, there and elevated ozone production. Humidity is also a factor in the

I

-1.0

I

0

1

I

1.o

2.0

A% (PPW

Fig. 6. Ozone production integrated over the day following imposition of the convective event, for low boundary layer and low free tropospheric constituent concentration profiles.

251

Streamlined family photochemistry module high boundary layer values. The storm processed profile is shown as the dashed curve in the figure. Ozone production remains in force over the entire column. This is to be expected even though the analog papers (Pickering et al., 1991, 1993) reported some losses. The removals occurred where NO dropped well below 10 pptv. Loss was actually observed in our model early in the morning under the hypothetical storm cloud, because reduced photolysis rates increased the NO* to NO ratio. However, once photochemical activity is fully restored at noon, NO returns, and in fact, thermal dissociation in the boundary layer prior to the convection gives enough NO, to raise free tropospheric levels. The pattern of the dashed curve is thus one which is prevalent in the Pickering et al. climatology, with a reduction in ozone production in the lower troposphere due to dilution by nitrogen oxide poor mid-tropospheric air, and enhancements in the upper free troposphere ascribed to inward mixing of surface material (Pickering et al., 1990). The overall magnitude of post storm ozone production is in line with Pickering et al. low NO, studies either in the column integral or locally. The reader will notice that we have not included lightning or cross tropopause sources of the nitrogen oxides in these tests of our numerics. In the Pickering et af. New Guinea analysis (1993) enrichment of NO, in the upper troposphere caused the results to be opposite in sign relative to those shown here. Overall our calculations display a net enhancement of ozone generation in the free troposphere. Again this is in genera1 agreement with the Pickering et ul. investigations. Results from the runs with high boundary layer concentrations below a pristine free troposphere are plotted in Fig. 7. One day ozone productions for the unprocessed column are indicated by the darker solid curve. The boundary layer values are consistent with our earlier reproductions of the Liu et al. Niwot

10.0

)

9.0

-

6.0

-

I -

--

Thiswork

a

m

Pickeringetal.

7.0 5

6.0 -

E -0 5.0 P j g

4.0 -

a

3.0 2.0 1.0 01

I 0.0

I 10.0

I

20.0

I

30.0

1

A@i(PPW

Fig. 7. Ozone production integrated over the day following imposition of the convective event, for high boundary layer and low free tropospheric constituent concentration profiles.

10.0 9.0

I

\ L II

1,

-

--

-

- -

Thiswork Pickeringetal.

6.0 7.0 e $

6.0 -

5 B

5.0 -

$ 3 I

4.0 -

a

3.0 2.0 1.0 0

I

I

I

I

0.0

10.0

20..0

30.0

Fig. 8. Ozone production integrated over the day following imposition of the convective event, for high boundary layer and high free tropospheric constituent concentration profiles.

Ridge experiments. The storm processed ozone rates are indicated by the dashed curve. Surface NO, remains near the ppbv level after convection while in the free troposphere the concentrations rise to on the order of 100 pptv. Two of the Pickering et al. computations are quite comparable in terms of their initial organic and nitrogen oxide profiles as well as their final ozone productions. They are summarized as the curves bounding shaded regions. Lighter shading indicates the nonstorm calculations. Rural air from the Oklahoma City area with somewhat larger initial free tropospheric NO, levels than in the scenario here gave the larger ozone productions at 2-4 km. Hypothetical convective mixing of a Brazilian Savannah biomass bum led to the lower curves just above the boundary layer (Pickering et al., 1992a). Shading has intentionally been restricted to the portion of the column lying above overlap between the Oklahoma and Brazil runs in order to emphasize similarities with our own work. Pronounced enhancements in free tropospheric ozone production are predicted in all the Fig. 7 cases. Pickering et al. obtained larger augmentations for two tropical forest atmospheres containing isoprene sources (Pickering et al., 1992a, c). The largest values were from a convected urban plume extending from the Brazilian city of Manaus, the economic heart of Amazonia. High isoprene and NO, concentrations coexisted in the Manaus situation, yielding free tropospheric ozone increases on the order of 25-SOppbv during 24 h periods. The Oklahoma City event gave an average enhancement factor of about 3. Our values in the figure fall within the range defined by the Pickering et al. comparisons. Computations with all concentrations set at maxima from Table 6 are outlined in Fig. 8. Free tropospheric ozone generation is greater than in the other simulations, as expected based on the elevated nitrogen oxide levels. An Oklahoma City urban plume forms the best comparable from the Pickering

252

Scott Elliott et al. are restricted to the top of the column, and vertically integrated effects are reduced in the free troposphere.

BENCHMARKING AND SCALING

8:

,/

;:111, : ‘1)

;; :I

2.0

:I-__

1.0 0

0.0

‘\

--_

'\

---___'.

10.0

,--_

20.0

-_ 30.0

Fig. 9. A comparison of ozone production in the high boundary layer low free tropospheric concentration scenario, following full free tropospheric convective mixing and a hot tower process.

et al. climatology (1992~). The Pickering values are shown as the lighter curves in the figure. Local

enhancements for the Oklahoma event were as high as an order of magnitude and an average for the region above the boundary layer was about 4. Overall our calculated enhancements are quite similar. The large local increases in production from Pickering et al. near 67 km are attributable to a layer of clean air with nitrogen oxide concentrations below 100 ppbv. They cannot be simulated with the uniform mixing ratio profiles input here. A single run conducted with release of the uplifted boundary layer material at the top of its storm column is plotted in Fig. 9. This hot tower type redistribution was performed on the central columns of Table 6 because they maximize NO, concentration gradients. Corresponding results assuming a well mixed free troposphere are included from Fig. 7. The computations are meant to be interpreted solely as a sensitivity test. The release of latent heat by water condensation is not permitted to feed back into the temperature distribution. The saturated water vapor concentrations in the upper few layers are thus lower limits and this will reduce photochemical activity somewhat. Unlike the majority of the Pickering et al. convective events, the one hot tower investigated (Pickering et al., 1991) was not accompanied by both before and after dew point profiles. As a next step in the development of our simple models of kinetics proceeding under convective conditions, we intend to couple latent heat into our transport processes. We will then be capable of assessing the sensitivity of tropospheric chemistry to more extreme mixing algorithms. For the moment we de-emphasize the hot tower redistributions by treating only the concentration data likely to give large photochemical changes. Ozone production potential enhancements

We argued conceptually in the preliminary counting section that an appropriately constructed family module could integrate the tropospheric chemical system in an individual grid cell for about lo4 operations per transport time step. The manipulations presumed that step sizes are near the upper limit set by the photolytic diurnal cycle and that iteration and matrix inversion are avoided. The values derived were thus near theoretical minima for detailed photochemical numerics (Elliott et al., 1993a, 1995). Here we verify the counting estimates by reporting performance tracing. The code employed in the nonlinearity tests was expanded from 10 real to 1000 hypothetical atmospheric layers for the performance runs. It was then advanced at the usual 1 h time steps for 10 atmospheric hours to obtain the statistics in Table 8. The operations totals were first summed for subroutines dealing with partitioning of families into constituents, and also for those involved in assembling photochemical rates. They were then divided by the number of cells and steps taken. Species and family numbers in the actual program are quite close to the round figures given during the crude counting procedures. The results indicate that our code can in fact be run for about lo9 operations per transport step on the hypothetical 100,000 point mesh. The chemistry algorithms described here are nearly positive definite because of the stabilizing nature of the family internal implicit integrations (Elliott et al., 1993a). Occasional occurrence of negative family concentrations required a simple filling strategy to restore positives while maintaining mass conservation. The total computational cost for the filling amounted only to a few tens of percent of the partitioning routines in our runs. If the mass balancing were to become prohibitively expensive, it could be implemented on a sensitivity basis and set aside during production. Data independent loops over a layer index were placed innermost in an attempt to achieve a high degree of vectorization within the chemistry (Boris & Winsor, 1982; Jacobson, 1994; Elliott et al., 1993b, 1995). A scaling experiment was performed to demonstrate the effectiveness of the measure. Identical kinetics calculations corresponding to the summer 500mbar runs with nitrogen oxide levels set to 100 pptv were conducted in total cell numbers varying from 1 to 1000 in half decadal increments.

Table 8. Operations counts and rates for key portions of the photochemistry package

Partitioning Rates Automation

Operations

Mtlops

IO4 IO’ -

100 150 0

_---r Streamlined family photochemistry module

I’

/

,

/

/

/

/

/

/

/

/

c-

Initialize

/’

Set Up

________-____-----

I

0.0

I 1.o

I

2.0

I 3.0

_I

Log,, (cell#)

Fig. 10. Scaling of relative floating point operations rates against variation in the number of cells in a 1-D model, for several groupings of subroutines.

Floating point operation rates were normalized to the maximum and then plotted in Fig. 10. Two samples from among the several subroutines called only once within the program are denoted by the dashed lines. The curve labelled “setup” represents a section of coding which reads in and organizes rate data. Vectorization is not attained because numerous logical statements act as inhibitors and because the dominant loops are over chemical dimensions such as reaction, species or family identification number. Since such routines are never repeated, speedups are not essential. Coding within the actual chemical time steps generally behaved as did the partitioning routine represented by the solid line. Some of the initialization segments vectorized naturally. The upper dashed curve displays the behavior of a subprogram which distributed constant background concentrations over the grid.

DISCUSSION We have shown through a series of box model simulations that our tropospheric chemistry module can successfully capture nonlinearities in the ozone production system in a variety of air types. Conditions investigated include the remote free troposphere, in which major oxidizing species are methane and CO (Levy ef al., 1985) and the continental, including natural and anthropogenic hydrocarbons (Liu et al., 1987; Trainer et al., 1987). A general result which we reproduce is that ozone production is most efficient per unit concentration of the nitrogen oxides at intermediate NO, levels (Liu et nl., 1987). Several groups have recently reviewed the reasons for this behavior (Liu et al., 1987; Crutzen, 1988; Lin et al., 1988; Jacob ef al., 1989; Chatfield & Delany, 1990; Jacob et al., 1993a, b). As the titration points are reached for NO to dominate HO* or ozone in reaction with peroxy radicals, ozone production efficiencies of oxidation sequences maximize. The titration

253

points vary with species, however, causing a smearing in the Liu et al. Niwot results and our simulations thereof. Increasing amounts of organic material enhance overall oxidation rates on both a mass basis (Jacob et al., 1989) and because of raised hydrogen oxide and radical input rates (Crutzen, 1988), and this also causes a stretching of the profiles. Eventually NO, begins to counteract these trends by converting the hydroxyl radical to nitric acid. This effect is directly apparent in the control curve for remote ozone production plotted along with the Niwot runs, because high NO, concentrations are attained for constant methane and carbon monoxide levels. Large quantities of organics will also suppress the hydroxyl and reactive, unsaturated species such as isoprene will themselves consume ozone (Jacob & Wofsy, 1988, 1990). The maximum in production relative to nitrogen oxides has several critical ramifications (Liu et al., 1987; Chatfield & Delany, 1990; McKeen et al., 1991a, b). In highly polluted areas, reductions in NO, may in fact increase ozone levels. This result is quite counterintuitive, and has been recognized for the urban context by smog researchers (Finlayson-Pitts & Pitts, 1986; Hough, 1988; Dodge, 1989). As pollutants spread to rural areas, however, the opposite may be true. The balance between the two regimes can only adequately be assessed in 3-D chemical tracer transport models. Several have been built for the regional and continental scales (Jacob et al., 1989, 1993a, b; McKeen et al., 199la, b; Sillman et al., 1990a, b), but global budgeting will require even broader calculations. High resolutions will definitely be advantageous as hemispheric and global chemistry models appear. Sources of the hydrocarbons and nitrogen oxides are often subgrid in scale. Instantaneous mixing through the volume of a regional or global scale model cell may well lead to artificial dilution and mixing, and so boost ozone production rates (Chatfield & Delany, 1990; Jacob et al., 1993a, b). Liu and coworkers (Liu et al., 1987; Lin et al., 1988) acknowledge the necessity for development of global 3-D chemical models, but have constructed an interesting technique for estimating ozone production rates based on NO, fluxes. It is essentially a partial quantification of the above arguments. Assume that NMHC and NO, levels are proportional to one another in a set of air parcels. As NO, concentrations rise, ozone production P fails to track them in a linear fashion. Production can be normalized to NO,Yas a matter of conceptual convenience, and we will denote the new variable as P ‘. Local nitrogen oxide concentrations can be expressed roughly as the product of emissions E and the lifetime T. The triple product P ‘TE can then be viewed either as production of ozone per unit NO, (P ‘T multiplied by emissions), or as P ’ multiplied by a concentration. The perspectives are distinct, but both are valid. The former permits the estimation of ozone generation over large scales

254

Scott Elliott et al

given the nitrogen oxide input pattern. A related bit of nomenclature which has arisen in the regional ozone literature is worth defining (e.g. McKeen et al., 199la, b; Jacob et al., 1993a, b). Chemistry in a fresh pollutant plume has been described as NMHC limited. A justification for the description is that sufficient NO, will be present to ensure that all organic oxidations lead to ozone. At low nitrogen oxide concentrations this will not necessarily be the case and the chemistry is sometimes termed NO, limited. Ozone production can also be estimated from the perspective of organic oxidations (Crutzen, 1982, 1988; Jacob et al., 1989). This is the angle we have favored in designing our reaction schemes. Both Liu ef al. (1987) and Lin et al. (1988) outline methods for estimating ozone generation rates based on organic oxidation, but the details are not developed. The nonlinearity issues extend to vertical mixing during convective events. Pollution sources of the nitrogen oxides and hydrocarbons are often located at the surface (Liu et al., 1987; Chatfield & Delany, 1990; Singh & Zimmerman, 1992). High pollutant concentrations in the boundary layer suppress ozone production. Mixing into the free troposphere can lower overall nitrogen oxide concentrations and so enhance local ozone generation per unit mass of organic (Liu et al., 1987; Lin et al., 1988; Pickering et al., 1989; Chatfield & Delany, 1990). Uplift also isolates ozone and its precursors from surface loss processes, extending their lifetimes while injecting them into faster upper level wind systems for large scale distribution (Liu et al., 1987; Chatfield & Delany, 1990; Jacob et al., 1993a, b). Pickering et al. (1989, 1990, 1991, 1992a+ 1993) have modelled enhancements in ozone production following an entire climatology of convective storms. Their studies coupled a 2-D cloud physics model with columnar photochemistry. We have used idealized vertical mixing ratio profiles and simple 1-D tracer overturn algorithms to demonstrate that our fast integrators can duplicate essential features of the Pickering et al. calculations. Mixing of nitrogen oxide rich boundary layer air upward can increase local and vertically integrated ozone production rates. Our storm chemistry model could be improved in the future in several ways. For example, we have not included an upper level lightning source of nitrogen oxides. The deletion of lightning is not uncommon in models of this type (Chatfield & Delany, 1990), but the Pickering collection does show that thunderstorms often provide their own, significant sourcing of NO, (Pickering et al., 1989, 1990, 1991). Pickering et al. also note that upper tropospheric NO, can be elevated due to stratospheric inputs (Pickering et al., 1993). We have intentionally omitted the stratosphere here in the interests of simplicity, but it could readily be accounted for through more elaborate initialization profiles. Most of our runs were performed assuming a well-mixed post storm free troposphere. Comparison with earlier results indicates that this is

often an adequate representation. A single run with hot tower type dynamics (Riehl & Malkus, 1958; Riehl & Simpson, 1979) suggests that ozone chemistry will be quite sensitive to the treatment of detrainment in convection parameterizations. The Pickering et al. collection clearly implies that this will be the case. Gidel (1983) and Feichter & Crutzen (1990), among others, investigated the influence of detrainment profiles on individual tracers and arrived at similar conclusions. Our hot tower scenarios must be viewed as preliminary because after storm temperature profiles were not available and so we have not included energy release from condensation. Tracer redistribution algorithms with self-consistent energetics will be constructed. They will enable us to explore the need for flexibility in detrainment computations. The consistency of our photolysis rate calculations is another area for potential improvement of our simulations. Cloud parameterizations and calculations of radiative transfer are sometimes decoupled by photochemists (Jacob & Wofsy, 1988, 1990). This deficiency will ultimately be eliminated in coming generations of global models. Studies of the photochemistry which follows convective overturn of the boundary layer have amply demonstrated that the process influences free tropospheric ozone levels on the regional and even global scales. However, the bulk of the detailed kinetics modelling has thus far been 1-D in character (Pickering et al., 1989; Chatfield & Delany, 1990). The computations are carried on for only a few days following the storm event in question, and then a common theme emerges; in the real atmosphere shear will separate the parcels of a column and mixing will dilute them, so that column representations are in fact quite limited. To determine the ultimate fate of ozone precursors in convected columns will require large scale 3-D simulations. We conceive of our fast tropospheric modules as a step in the direction of detailed large scale kinetics investigations. They reflect the major features of existing models of the chemistry of convective storms, but can also bc ported into massive 3-D grids. SUMMARY There are a number of valid and even pressing justifications for the study of tropospheric chemistry in global and 3-D grids. An intimate interaction of kinetics and dynamics is responsible for the distribution of major oxidative pollutants (Thompson, 1992; Jacob et al., 1993a, b). Several trace species in the lower atmosphere are both radiatively and kinetically active, creating a close connection between tropospheric chemistry and climate (Ramanathan et al., 1987; Hauglustaine et al., 1994). Ozone is a prime example in both categories. The photochemical kinetics of the troposphere present a difficult numerical challenge, however, because of a combination of complexity, nonlinearity and stiffness (Crutzen &

255

Streamlined family photochemistry module Zimmerman, 1991; Jacobson & Turco, 1994; Elliott et al., 1995). Organic chemistry contributes greatly to the complexity (Singh & Zimmerman, 1992), and a coupling of organic oxidations with nitrogen oxide catalysts is a major manifestation of the nonlinearity (Liu et al., 1987). Organic oxidations tend to be most ozone productive at intermediate levels of NO, falling between those of the polluted and pristine regimes (Lin et al., 1988; Pickering et al., 1989). Integrators for simulation of global tropospheric chemical kinetics must accurately represent the ozone to NO, relationship. Several strategies for circumventing the numerical barriers to modelling of stiff tropospheric kinetics have been conceived. Some involve optimization of traditional LMBDF or Gear code solvers through sparse matrix techniques (Byrne & Hindmarsh, 1987; Rotman et al., 1993; Jacobson & Turco, 1994). Order of magnitude speedups have been achieved but still may not be adequate. Others are essentially off-line parameterizations which entail considerable sacrifices in chemical resolution (Marsden et al., 1987; Jacob et al., 1989). We have developed a suite of highly automated and modular chemical solvers based on the family approach (Turco & Whitten, 1974; Kao et al., 1990; Elliott et al., 1993a, 1995; Jacobson, 1994; Shen et al., 1994). They are designed specifically for the global 3-D context. Inorganic versions are currently running interactively with climate and dynamics in GCM coding (Kao et al., 1990; Elliott et al., 1995). Tropospheric extensions of the modules have recently been completed by several of us (MS, MZJ). The new programs encompass key nonmethane hydrocarbons for use as surrogates for the entire tropospheric organic system (Elliott et al., 1993b). The present work describes tests of their ability to reproduce the intricacies of the ozone production/nitrogen oxide nonlinearities. The testing exercise begins in the O-D or box modelling mode. We run simulations paralleling the 500mbar, remote free tropospheric work of Levy et al. (1985) and obtain results which are quite comparable. With methane and carbon monoxide as the chief ozone sources, production rates rise to a maximum at NO, levels of several hundred pptv, then fall steadily to several ppbv. The module is then moved into the continental regime. Parameters are set to measurement constraints for the Niwot Ridge experimental station, which receives air masses from both the pristine and urban surface regimes. Our code duplicates simulations by Liu et al. (1987) of daily integrated ozone production for several combinations of anthropogenic and natural hydrocarbons, at nitrogen oxide levels ranging up to 3 ppbv. Sources of the nitrogen oxides and organics as pollutants tend to occur at the surface boundary in the troposphere, so that large gradients in the vertical are the rule (Logan et al., 1981; Penner et al., 1991; Singh & Zimmerman, 1992). Overturn of the boundary layer and free troposphere in convective storms is

a critical mechanism through which mixing of high and low NO, air masses takes place (Gidel, 1983; Chatfield & Crutzen, 1984). We have also shown that our fast, packaged integrator can reproduce enhancements in ozone production which follow the vertical mixing induced by storms. Generalized scenarios are constructed which closely resemble key events from the Pickering et al. (1989, 1990, 1991, 1992a-c, 1993) convection climatology. Ozone production is integrated as a function of altitude in a 1-D model for the day following imposition of a storm parameterization. Results agree well with the more sophisticated but numerically intensive Pickering et al. work. Crude operations counts are performed to estimate the potential computational efficiency of our family strategy (Elliott et al., 1993a, 1995). We find that a well-designed family code should be able to integrate a full 300 reaction tropospheric scheme for as few as 10,000 operations per cell step, at time increments of 1 h or more. Benchmarking on the same coding used in nonlinearity tests shows that this ideal is closely approached. Vectorization is verified through scaling experiments. The number of cells in the model is steadily increased from one to 1000, and routines within the chemical cycle maximize in terms of their flops rate. The vectorization is achieved by moving data independent physical loops inside of chemical loops over reaction, species or family. The physical cells are numerically independent because chemistry is treated as an infinite medium calculation in the time split mode (Boris & Winsor, 1982; Elliott et al., 1993a, b). They will translate readily to massively parallel processors as the number of dimensions is increased (Elliott & Jones, 1994). At the time steps and operations totals we have obtained, the tropospheric coding should enable decadal scale chemistry runs even in the highest resolution atmospheric GCM gridding (Mahlman & Umscheid, 1987; Dannevik et al., 1993). -Current versions of the fast chemical integrators outlined here can be obtained by contacting the R. P. Turco group, Department of Atmospheric Sciences,

Program availability

University of California, Los Angeles.

REFERENCES Abbatt J. P. D. & Molina M. J. (1992) Geophys. Res. Left. 19, 461464. Aiken A. C., Herman J. R., Maier E. J. & McQuillan C. J. (1982) J. Geophys. Res. 87, 3105-3118. Aiken A. C., Herman J. R., Maier J. R. & McQuillan C. J. (1983) Planet. Space Sci. 31, 1075-1082. Almasi G. S. & Gottlieb A. (1989) Highly-parallel Computing. Benjamin Cummings, Redwood City, CA. Anderson D. E. (1983) Pianef Soace Sci. 31. 1517-1523. Anderson D. E. & Lloid S. A. (i990) J. Geophys. Res. 95, 7429-7434. Andrews D. G., Mahlman J. D. & Sinclair R. W. (1983) J. Armos. Sci. 2768-2777.

Arakawa A. (1966) J. Comput. Phys. 1, 119-143. Arakawa A. (1972) Tech. Rept. 7, Department of Meteorology, University of California at Los Angeles.

256

Scott Elliott et al.

Arakawa A. & Suarez M. J. (1983) Mon. Wea. Rev. 111, 34-45.

Crutzen

P. J. (1973) Pure

Appl.

Geophys.

106-108,

1385-1399.

Arellano J., Duynkerke P. G. & van Weele M. (1994) J. Geophys. Rex 99, 3699-3705.

Atherton C. S. & Penner J. E. (1990) J. Geophys. Rex 95, 14,027-14,038.

Atkinson R. & Lloyd A. C. (1984) J. Phys. Chem. Ref Data. 13, 3 15-398.

Atkinson R., Baulch D. L., Cox R. A., Hampson R. F., Kerr J. A. Jr Troe J. (1989) J. Phys. Chem. Ref. Data 18, 881-1097. Austin P. M. & Houze R. A. (1973) J. Ammos. Sci. 30, 1100-1111. Blake D. S. & Rowland F. S. (1986) Nature, 321, 231-233. Boris J. P. (1989) Ann. Rev. Fluid Mech. 21, 345-385. Boris J. P. & Winsor N. K. (1982) Vectorized computation of reactive flow. In Parallel Computations (Edited by G. Rodrigue), pp. 173-214. Academic Press, New York. Braham R. R. (1952) J. Meteorol. 9, 227-242. Brasseur G. & Solomon S. (1984) Aeronomy of the Middle Atmosphere. Reidel, Dordrecht. Brass& G., Hitchman M. H., Walters S., Dymek M., Falise E. & Pirre M. (1990) J. Geoohvs. Res. 95. . . . 5639-5655. Brown P. N. & Hindmarsh A. C. (1986a) SIAM J. Numer. Anal. 23, 610438.

Brown P. N. & Hindmarsh A. C. (1986b) Appl. Math. Compul. 31, 40-91.

Brown P. N. & Saad Y. (1990) SIAM J. Sci. Sfaf. Compur. 11, 45048 1.

Brown P. N., Byrne G. D. & Hindmarsh A. C. (1989) SIAM J. Sci. Stat. Comput. 10, 1038-1051. Briihl C. & Crutzen P. J. (1989) Geophys. Res. Lett. 16, 703-706.

Byers H. R. (1959) General Meteorology. McGraw-Hill, New York. Byrne G. D. & Hindmarsh C. (1987) J. Comput. Phys. 70, 162. Chameides W. L (1984) J. Geophys. Res. 89, 4739-4755. Chameides W. L. & Walker J. C. G. (1973) J. Geophys. Res. 78, 8751-8760. Chameides W. L. & Davis D. D. (1982) J. Geophys. Res. 87, 48634877. Chandrasekhar S. (1960) Radiative Transfer. Dover, New York. Chang J. S., Hindmarsh A. C. & Madsen N. K. (1973) Simulation of chemical kinetics transport in the stratosphere. In Stiff Differential Systems (Edited by R. S. Willoughby). Plenum Press, New York. Chang J. S., Brost R. A., Isaksen I. S. A., Madronich S., Middleton P., Stockwell W. R. & Walcek C. J. (1987) J. Geophys. Res. 92, 14,681-14,700. Chapman S. (1930) Mem. R. Mereorol. Sot. 3, 103-117. Chatfield R. B. & Crutzen P. J. (1984) J. Geophys. Res. 89, 7111-7132. Chatfield R. B. & Delany A. C. (1990) J. Geophys. _ . Res. 95, 18,473-18,488. Chene C. P. & Houze R. A. (1979) Mon. Weu. Rev. 107. ~ 13’%1381. Ching J. K. S. & Alkezweeny A. J. (1986) J. Climate Appl. I

Meteorol. 25, 1702.

Ching J. K. S., Shipley S. T. & Browell E. V. (1988) Atmos. Environ. 22, 225-242.

Cicerone R. J., Elliott S. & Turco R. P. (1991) Science 254, 1191-l 194. Costen R. C., Tennille G. N. & Levine J. S. (1988) J. Geophys. Res. 93, 15,941-15,954. Crank J. & Nicholson P. (1947) Proc. Camb. Phil. Sot. 43, 5067.

Crutzen P. J. (1970) Q. J. R. Met. Sot. %, 320-333. Crutzen P. J. (1971) J. Geophys. Res. 76, 7311-7327.

Crutzen P. J. (1982) The global distribution of hydroxyl. In Atmospheric Chemistry (Edited by E. D. Goldberg), pp. 313-328. Crutzen P. J. (1988) Tropospheric ozone: an overview. In Tropospheric Ozone (Edited by I. S. A. Isaksen), pp. 3-32. Reidel, Dordrecht. Crutzen P. J. & Zimmermann P. H. (1991) Tellus. 43AB, 136151. Crutzen P. J., Isaksen I. S. A. & McAfee J. R. (1978) J. Geophys. Res. 83, 345-362. Crutzen P. J., Muller R., Bruhl C. & Peter T. (1992) Geophys. Res. Lett. 19, 1113-l 116. Dabdub D. & Seinfeld J. H. (1994) Submitted to Atmos. Environ.

Dannevik W. P., MacCracken M. C. & Mirin A. A. (1993) Toward a high-performance climate systems model. In Computing at the Leading Edge (Edited by A. A. Mirin). National Energy Research Supercomputer Center, Lawrence Livermore National Laboratory, Livermore, CA. Demerjian K. L., Schere K. L. & Peterson J. T. (1980) Advances in Environ. Sci. and Tech. 10 (Edited by J. N. Pitts Jr and R. L. Metcalf), pp. 369459. Wiley, New York. Demore W. B., Sander S. P., Golden D. M., Molina M. J., Hampson R. F., Kurylo M. J., Howard C. J. & Ravishankara A. R. (1990) Jet Propulsion Laboratory Publication 9&l. Dickerson R. R. (1980) Ph.D. Thesis, University of Michigan, Ann Arbor. Dickerson R. R., Stedman D. H., Chameides W. L., Crutzen P. J. & Fishman J. (1979) Geoohvs. Res. Lerr. 6, . 833-835. Dignon J. (1991) UCRL-ID-108511 Lawrence Livermore National Laboratory, October 7. Dignon J. (1992) Atmos. Environ. 26A, 1157-1163. Dodge M. C. (1989) J. Geophys. Res. 94, 5121-5136. Donahue N. M. & Prinn R. G. (1990) J. Geophys. Res. 95, 18,387-18,411. Douglass A. R., Jackman C. H. & Stolarski R. S. (1989) J. Geophys. Res. 94, 9862-9872. Drdla K., Turco R. P. & Elliott S. (1993) J. Geophys. Res. 98, 8965-898 1. Drummond J. W., Ehhalt D. H. & Volz A. (1988) J. Geophys. Res. 93, 15,831-15,849. Dunker A. M. (1986) Atmos. Environ. 20, 479486. Elliott S. & Jones P. (1994) Los Alamos National Laboratory Technical Report. LA-12814-MS, pp. l-35. Elliott S., Turco R. P. & Jacobson M. Z. (1993a) Computers Chem. 17, 91-102. Elliott S., Kao C. Y. J., Turco R. P. & Zhao X. P. (1993b) Los Alamos National Laboratory Technical Report LA-12539-M& pp. l-72. Elliott S., Turco R. P., Zhao X. P., Shen M. & Kao C. Y. J. (1993~) Efficient and modular kinetics packages for global scale photochemical modelling. In The Role of I

._

Meteorology in Managing the Environment of the Nineties (Edited by Amiram Roffman), pp. 73-84.

AWMA, Pittsburgh. Elliott S., Cicerone R. J., Turco R. P., Drdla K. & Tabazadeh A. (1994) .I. Geophys. Res. 99, 3497-3508. Elliott S., Turco R. P., Zhao X. P., Kao C. Y. J. & Shen M. (1995) J. Appl. Meteor. 34, 694-718. Falls A. H. & Seinfeld J. H. (1978) Environ. Sci. Technol. 12, 1398-1406.

Fan S. M. & Jacob D. J. (1992) Nature 359, 522-524. Farrow L. A. & Edelson D. (1974) Inr. J. Chem. Kin. 6, 787-800.

Farrow L. A. & Graedel T. E. (1977) J. Phys. Chem. 81, 248&2485 t Feichter J. & Crutzen P. J. (1990) Tells. 42B, 100-117.

Streamlined family photochemistry module Finlayson-Pitts 235-249. Finlayson-Pitts

B. J. (1993) Res. Chem. Intermed.

19,

B. J. & Pitts J. N. (1986) Atmospheric

Chemistry: Fundamentals and Experimental

Techniques.

Wiley, New York. Fishman J. & Crutzen P. J. (1978) Nature 274, 855-858. Frederick J. E. & Lubin D. (1988) J. Geophys. Res. 93, 3825-3832. Garcia R. R. & Solomon S. (1983) J. Geophys. Res. 88, 1379-1400. Gear C. W. (1971a) Numerical Initial Value Problems in Ordinary Dt&mtial Equations. Prentice Hall, Englewood Cliffs, New Jersey. Gear C. W. (1971b) Inform. Process. 68, 187-193. Gear C. W. (1988) Rpt No. UIUCDCS-R-88-1442. University of Illinois at Urbana-Champaign, Department of Computer Sciences. Gear C. W. (1989) Parallel solution of ODES. In RealTime Integration Metho& for Mechanical System Simulation (Edited by E. J. Hang and R. C. Deyo),

~~$33-248.

NATO ASI Series F 69, Springer Verlag,

Gear C. W. (1993) Appl. Numer. Math. 11, 27-43. Gear C. W. & Xuhai X. (1993) Appl. Numer. Math. 11, 4568. Gidel L. T. (1983) J. Geophys. Res. 88, 65876599. Gidel L. T., Crutzen P. J. & Fishman J. (1983) J. Geophys. Res. 88, 66226640. Golombek A. & Prinn R. G. (1986) J. Geophys. Res. 91, 39854001. Graedel T. E. (1977) J. Phys. Chem. 81, 2372-2374. Graedel T. E. & Goldberg K. I. (1983) J. Geophys. Res. 88, 10,865510,882. Graedel T. E., Mandich M. L. & Weschler C. J. (1986) J. Geophys. Res. 91, 5205-5221. Granier C. & Brasseur G. (1991) J. Geophys. Res. %, 2995-301 I. Greenberg J. P. &Zimmerman P. R. (1984) J. Geoohys. _ _ Res. 89, 47674778. Hack J. J., Boville B. A., Briegleb B. P., Kiehl J. T., Rasch P. J. & Williamson D. L. (1992) NCAR Technical Note TN-382- + STR. National Center for Atmospheric Research, Boulder, CO. Hack J. J., Boville B. A., Kiehl J. T., Rasch P. J. & Williamson D. L. (1993) Submitted to J. Geophys. Res. Hahn C. J., Warren S. G., London J., Chervin R. M. & Jenne R. (1982) NCAR/TN-201 + STR. National Center for Atmospheric Research, Boulder, CO. Haltiner G. J. & Williams R. T. (1980) Numerical Prediction and Dynamic Meteorology. Wiley, New York. Hameed S. & Dignon J. (1991) J. Air & Waste Mgmt Assoc. 42.

Hansen J., Russell G., Rind D., Stone P., Lacis A., Lebedeff S., Ruedy R. & Travis L. (1983) Mon. Wea. Rev. 111, 609662. Hanson D. R. & Ravishankara A. R. (1991a) J. Geophys. Res. %, 5081-5090. Hanson D. R. & Ravishankara A. R. (1991b). J. Geophys. Res. %, 17,307-17,314. Hauglustaine D. A., Granier C., Brasseur G. P. & Megie G. (1994) J. Geophys. Res. 99, 1173-1186. Hesstvedt E., Hov 0. & Isaksen I. S. A. (1978) Int. J. Chem. Kin. 10, 971-994. Hough A. M. (1988) J. Geophys. Rex 93, 37893812. Houghton J. T., Jenkins G. J. & Ephraums J. J. (1990) Climate Change: The IPCC Scientific Assessment. Cam-

bridge University Press, Cambridge. Houghton J. T., Callander B. A. SCVarney S. K. (1992) Climate Change 1992: The Supplementary Report to the IPCC Scientific Assessmenf. Cambridge University Press,

Cambridge. Hov 0.. Schjoldager J. & Wathne B. (1983) J. Geophys. Res. 88, 10.679-10,688.

257

Isaac, G. A. (1983) The Vertical transport and redistribution of pollutants by clouds. In Transaction Specialty Conference on the Meteorology of Acid Deposition, pp. 4965 12. Air Pollution Control Association, Hartford. Isaac G. A., Strapp J. W., Wiebe H. A., Leaitch W. R., Kerr J. B., Anlauf K. G., Summers P. W. & Macpherson J. I. (1983) The role of cloud dynamics in redistributing pollutants and the implications for scavenging studies. In Proc. Fourth Int. Conf on Precipitation Scavenging, Santa Monica, CA. lsaksen I. S. A. & Hov 0. (1987) Tellus. 39B, 271-285. Isaksen. A., Midtbs K. N., Sunde J. & Crutzen P. J. (1977) Geophys. Norv. 31, 11-26. Jacob D. J. & Wofsy S. C. (1988) J. Geophys. Rex 93, 1477-1486. Jacob D. J. & Wofsy S. C. (1990) J. Geophys. Res. 95, 16,737-l 6,754. Jacob D. J., Prather M. J., Wofsy S. C. & McElroy M. B. (1987) J. Geophys. Res. 92, 66146626. Jacob D. J., Sillman S., Logan J. A. & Wofsy S. C. (1989) J. Geovhvs. Res. 94, 8497-8509. Jacob D: J:, Logan J.’ A., Yevich R. M., Gardner G. M., Spivakovsky C. M., Wofsy S. C., Munger J. W., Sillman S., Prather Zimmerman

M. J., Rodgers M. O., Westberg H. & P. R. (1993a) J. Geovhvs. _ _ Res. D8.

14,797-14,816. Jacob D. J., Logan J. A., Gardner G. M., Yevich R. M., Spivakovskv C. M. & Wofsv S. C. (1993b) J. Geoohvs. Res. 98, 14,817-14,826. ~ ’ ’ ’ Jacobson M. Z. (1994) Ph.D. Thesis, developing, coupling, and applying a gas aerosol, transport, and radiation model to study urban and regional air pollution. University of California, Los Angeles. Jacobson M. Z. & Turco R. P. (1994) Atmos. Environ. 28, 273-284.

Jones R. J., Austin J., McKenna D. S., Anderson J. G., Fahey D. W., Farmer C. B., Heidt L. E., Kelly K. K., Murphy D. M., Proffitt M. H., Tuck A. F. & Vedder J. F. (1989) J. Geophys. Res. 94, 11,529-11,558. Kanakidou 787-801.

M. & Crutzen

P. J. (1993) Chemosphere. 26,

Kanakidou M., Singh H. B., Valentin K. M. & Crutzen P. J. (1991) J. Geophys. Res. 96, 15,39515,413. Kao C. Y. J., Glatzmaier G. A., Malone R. C. & Turco R. P. (1990) J. Geophys. Res. 95, 22,495-22,512. Kao C. Y. J., Tie X. & Mroz E. (1992a) Simulations of greenhouse trace gases using the Los Alamos chemical tracer model. In American Institute of Physics Conference Proc. Series on Global Climate Change, pp. 134-139. Kao C. Y. J., Tie X., Mroz E., Cunnold D. & Alyea F. (1992b) J. Geophys. Res. 97, l&827-15,838. Kasibhatla P. S., Levy II H., Moxim W. J. & Chameides W. L. (1991) J. Geophys. Res. 96, 18,631-18,646. Kasibhatla P. S., Levy H. & Moxim W. J. (1993) J. Geophys. Res. 98, 7165-7180. Kasting J. F. & Singh H. B. (1986) J. Geophys. Res. 91, 13,239-13,256. Kasting J. F., Pollack J. B. & Crisp D. (1984) J. Atmos. Chem. 1, 403408. Kaye J. A. & Rood R. B. (1989) J. Geophys. Res. 94, 1057-1083. Keene W. C., Jacob D. T., Pszenny A. A. P., Duce R. A., Schulitz-Tikos J. J. 8c Galloway J. N (1993) J. Geophys. Res. 98, 9047-9049. Kondo U., Matthew W. A., Iwata A., Morita Y. & Takagi M. (1987) J. Atmos. Chem. 5, 37-58. Lapidus L. & Seinfeld J. H. (1971) Numerical Solutions of Ordinary Differential Equations. Academic Press, New York. Lapidus L., Aiken R. C. & Liu Y. A. (1974) The occurrence and numerical solution of physical and chemical systems having widely varying time constants. In StiffDifferential

258

Scott Elliott et al.

Systems (Edited by R. S. Willoughby). Plenum Press,

New York. Lelieveld J. & Crutzen P. J. (1990) Nature 343, 227-233. Lelieveld J. & Crutzen P. J. (1991) J. Atmos. Chem. 12, 229-267.

Levesque J. M. & Williamson J. W. (1989) A Guidebook to FORTRAN on Supercomputers. Academic Press, New York. Levy H. (1971) Science 173, 141-143. Levy H., Mahlman J. D., Moxim W. J. & Liu S. C. (1985) 1. Geophys. Rex 90, 3753-3772. Lin X., Trainer M. & Liu S. C. (1988) J. Geophys. Res. 93, ^ 15,879-I 5,888. Liu S. C. (1977) Geoohvs. Res. Lett. 4. 325-328. Liu S. C. ‘& Trainer‘M. (1988) J. Atmos. Chem. 6, 221233.

Liu S. C., Kley D., McFarland M., Mahlman J. D. & Levy H. (1980) J. Geophys. Res. 85, 75467552. Liu S. C., McFarland M., Kley D., Zafiriou 0. & Huebert B. (1983) J. Geophys. Res. 88, 136&1368. Liu S. C., McAfee J. R. & Cicerone R. J. (1984) J. Geophys. Res. 89, 7291-7297. Liu S. C., Trainer M., Fehsenfeld F. C., Parrish D. D., Williams E. J., Fahey D. W., Hubler G. & Murphy P. C. (1987) J. Geophys. Res. 92, 41914207. Liu S. C., Trainer M, Carroll M. A., Hubler G., Montzke D. D., Norton R. B., Ridley B. A., Walega J. G., Atlas E. L.. Heikes B. G.. Huebert B. J. & Warren W. (1992) \ , J. Geophys Res. 97,‘ 10,463-10,471. Logan J. A. (1983) J. Geophys. Res. 88, 10,785-10,807. Logan J. A. (1985) J. Geophys. Res. 90, 10,463-10,482. Logan J. A., Prather M. J., Wofsy S. C. & McElroy M. B. (1981) J. Geophys. Res. 86, 7210-7254. London J. (1952) 1. Meteorol. 9, 145-151. Luther F. M. (1980) Publ. UCRL-50042-80. Lawrence Livermore Laboratory, Livermore, CA. Luther F. M & Gelinas R. J. (1976) J. Geophys. Res. 81, 1125-l 132. Luther F. M. & Wuebbles D. J. (1976) UCRL-78911 Lawrence Livermore Laboratory, Livermore, CA. Madronich S. (1987) J. Geophys. Res. 92, 9740-9752. Mahlman J. D. & Moxim W. J. (1978) J. Atmos. Sci. 35, 134&1374. Mahlman J. D. & Umscheid L. J. (1987) Comprehensive modeling of the middle atmosphere: the influence of horizontal resolution. In Transport Processes in the Middle Atmosphere (Edited by G. Visconti and R. Garcia), pp. 251-256. Reidel, Dordrecht. Malone R. C., Auer L. H., Glatzmaier G. A, Wood M. C. & Toon 0. B. (1986) J. Geophys. Res. 91, 1039-1053. Manabe S., Smagorinsky J. & Stickler R. F. (1965) Mon. Wea. Reo. 93, 769-798. Marsden A. R., Frenklach M. & Reible D. D. (1987) J. Air Pollut.

Control Assoc. 37, 370-376.

Mastenbrook H. J. (1966) Rep. 6447, Naval Res. Lab., Washington, DC. Mastenbrook H. J. (1968) J. Afmos. Sci. 25, 299-311. McElroy M. B. & Hunten D. M. (1966) J. Geophys. Res. 71, 363553638. McElroy M. B., Salawitch R. J., Wofsy S. C. & Logan J. A. (1986) Nature 321, 159-162. McKeen S. A., Hsie E.-Y. & Liu S. C. (1991a) J. Geophys. Res. 96, 15,377-l 5,394. McKeen S. A., Hsie E.-Y., Trainer M., Tallamraiu R. & Liu S. C. (1991b) J. Geophys. Res. %, IO.8099iO,845. McRae G. J.. Goodin W. R. & Seinfeld J. H. (1982a) Almos. Environ. 16, 619-696.

McRae G. J., Goodin W. R. & Seinfeld J. H. (1982b) J. Comput. Phys. 45, 142.

Mihelcic D., Musgen P. & Ehhalt D. H. (1985) J. Afmos. Chem. 3, 341-361. Miyahara S., Hayashi Y. & Mahlman J. D. (1986) J. Atmos. Sci. 43. 1844.

Nicolet M., Meier R. R. & Anderson D. E. (1982) Planet Space Sci. 30, 935-983. Niewiadomski M. (1986) Armos. Environ. Zo, 139-145. Oran E. S. & Boris J. P. (1981) Prog. Energy Cornbust. Sci. 7, I-72.

Parrish D. D., Murphy P. C., Albritton D. L. & Fehsenfeld F. C. (1983) Atmos. Environ. 17, 1365-1379. Parrish D. D., Fahey D. W., Williams E. J., Liu S. C., Trainer M., Murphy P. C., Albritton D. L. & Fehsenfeld F. C. (1986a) J. Afmos. Chem. 4, 63380. Parrish D. D., Huebert B., Norton R. B., Bollinger M. J., Liu S. C., Murphy P. C., Albritton D. L. & Fehsenfeld F. C. (1986b) J. Geophys. Res. 91, 5379-5393. Paulson S. E. & Seinfeld J. H. (1992) J. Geophys. Res. 97, 20,703.-20,715. Penner J. E., Atherton C. S., Dignon J., Ghan S. J. & Walton J. J. (1991) J. Geophys. Res. %, 959-990. Pickering K. E., Thompson A. M. &Dickerson R. R. (1989) J. Geophys Res. 94, 14,879-14,892. Pickering K. E., Thompson A. M., Dickerson R. R., Luke W. T.. McNamara D. P.. Greenberg J. P. & Zimmerman P. R. ‘(1990) J. Geophys.’ Res. 95, k,O49--14,062. Pickering K. E., Thompson A. M., Scala J. R., Tao W. K., Dickerson R. R., Simpson J. & Garstang M. (1991) J. Geophys. Res. 96, 309993114. Pickering K. E., Thompson A. M., Scala J. R., Tao W. K. & Simpson J. (1992a) J. Atmos. Chem. 14, 297-313. Pickering K. E., Scala J. R., Thompson A. M., Tao W. K. & Simpson J. (1992b) Geophys. Res. Letr. 19, 289-292. Pickering K. E., Thompson A. M., Scala J. R., Tao W. K., Dickerson R. R. & Simpson J. (1992~) J. Geophys. Res. 97, 17,985-18,000. Pickering K. E., Thompson A. M., Tao W. K. & Kucsera T. L. (1993) J. Geoohvs. Res. 98. 8737-8749. Pitcher E. J., Malone k.‘C., Ramanathan V., Blackman M. L., Puri K. & Bourke W. (1983) J. Atmos. Sci. 40, 631642.

Prather M. J. (1974) Astrophys. J. 192, 787-792. Prather M. J. (1992) Nature 355, 534537. Prather M. J., McElroy M. B., Wofsy S. C., Russell G. & Rind D. (1987) J. Geophys. Res. 92, 6579-6613. Pszenny A. A. P., Keene W. C., Jacob D. J., Fan S., Maben J. R., Aetuo M. P., Springeryoung M. & Galloway J. N. (1993) Geophys. Res. Lett. 20, 699-702. Ramanathan V. (1988) The radiative and climatic consequences of the changing atmospheric composition of trace gases. In The Changing Atmosphere (Edited by F. S. Rowland and I. S. A. Isaksen)., pp. 1599186. Wiley, New York. Ramanathan V., Callis L., Cess R., Hansen J., Isaksen 1. S. A., Kuhn W., Lacis A., Luther F., Mahlman J., Reck R. & Schlesinger M. (1987) Rea. Geophys. 25, 144-1482. Ridler G. W., Ridler P. F. & Sheppard J. G. (1977) J. Phys. Chem. 81, 243552437. Ridler P. F. (1977) J. Phys. Chem. 81, 2419-2423. Riehl H. & Malkus J. S. (1958) Geophysics. 6, 503-538. Riehl H. & Simpson J. (1979) Contrib. Armos. Phys. 52, 287-305.

Rodriguez J. M., Lo M. K. W., Aze N. D., Pierce S. D., Anderson J. G.. Fahey D. W., Jakoubek R., Kelly K., Mount G. H., Proffitt M. H., Ravishankara A: R., Schmeltekopf A. L., Tuck A. F., Wahner A.. Farmer C. B., Toon G. C., Coffey M. T., Heidt L. E., Mankin W. G., Chan K. R., Starr W. L., Vedder J. F., Browell E. V., & McCormick M. P. (1989) J. Geophvs. . _ Res. 94. 16,683-16,703. Rood R. B. (1987) Rev. Geophys. 25, 71-100. Rood R. B., Kaye J. A., Nielsen E., Schoeberl M. R. & Geller M. A. (1987) Physica Scrip& 36, 337-354. Rood R. B., Kaye J. A., Douglass A. R., Allen D. J., Steenrod S. & Larson E. M. (1990) 1. Ammos. Sci. 47, 26962709.

Streamlined

family ph lotochemistry

Rudolph J. (1988) J. Geophys. Res. 93, 8367-8377. Rudolph J., Ehhalt D. H., Schmidt U. & Khedim A. (1982) Vertical distribution of some C& hydrocarbons in the nonurban troposphere. In Proc. Second Symp. on Composition sf the Nonurban Troposphere. Williamsburg, VA. Rudolph J., Ehhalt D. H. & Khedim A. (1984a) J. Afmos. Chem. 2, 117.-124. Rudolph J., Jobsen C., Khedim A. & Johnen F. (1984b) Phys. Chem. Behacior of Atmos. Poll., Third European Symp., Varese, Italy. Shen M.. Turco R. P., Paulson S. E. & Jacobson M. Z. (1994) J. Geophys. Res. Shimazaki T. (1985) Minor Constiruenrs in rhe Middle Atmosphere. Reidel, Dordrecht. Sillman S., Logan J. A. & Wofsy S. C. (1990a) J. Geophys. Res. 95, 1837-1851. Sillman S., Logan J. A. & Wofsy S. C. (1990b) J. Geophys. Res. 95, 573 l-5748. Singh H. B. & Zimmerman P. (1992) Atmospheric distribution and sources of non methane hydrocarbons. In Gaseous Pollutants: Characterization and Cycling (Edited by J. 0. Nriagu), pp. 177-235. Wiley, New York. Singh H. B., Viezee W. & Salas L. J. (1988) J. Geophys. Res. 93, 15,861m 15,878. Sokal R. R. & Rohlf F. J. (1981) Biometry. Freeman, New York. Solomon S. (1988) Rezl. Geophys. 26, 131-148. Solomon S. & Garcia R. R. (1983) Eur. Space Agency [Spec. Publ.] ESA SP-183, pp. 117-122. Solomon S., Garcia R. R., Rowland F. S. & Wuebbles D. J. (1986) Na/ure 321, 755--757. Spivakovsky C. M., Yevich R., Logan J. A., Wofsy S. C., McElroy M. B. & Prather M. J. (1990a) J. Geophys. Res. 95, 18.441- 18.472. Spivakovsky C. M., Wofsy S. C. & Prather M. J. (1990b) J. Geophys. Res. 95, 18.433318,440. Thompson A. M. (1984) J. Geophys. Res. 89, 1341-1349. Thompson A. M. (1992) Science 256, 1157-I 165. Thompson A. M. & Cicerone R. J. (1986a) J. Geophys. Res. 91, 10,852~ 10,864. Thompson A. M. & Cicerone R. J. (1986b) Nature 321, 1478-1480. Tille K. J. W., Saelsberg M. & Bachmann K. (1985) Atmos. Emiron. 19. I75 I 1760.

module

259

Toon 0. B., Turco R. P., Westphal D., Malone R. & Liu M. S. (1988) J. Afmos. Sci. 45, 2123-2143. Trainer M, Hsie E. Y., McKeen S. A., Tallamraju R., Parrish D. D., Fehsenfeld F. C. & Liu S. C. (1987) J. Geophys. Res. 92, II,879911,894. Turco R. P. (1985) The photochemistry of the stratosphere. In The Photochemistry of Atmospheres (Edited by J. S. Levine). Academic Press, New York. Turco R. P. & Whitten R. C. (1974) J. Geophys. Res. 79, 3179-3185. Turco R. P. & Whitten R. C. (1977) NASA Tech. Publ. TP-1002. Turco R. P., Hamill P., Toon 0. B., Whitten R. C. & Kiang C. S. (1979a) J. Armos. Sri. 36, 699.-717. Turco R. P., Hamill P., Toon 0. B., Whitten R. C. & Kiang C. S. (1979b) NASA Tech. Publ. TP-1362. Walton J. J., MacCracken M. C. & Ghan S. J. (1988) J. Geoohvs. Res. 93. 8339-8354. Warren ‘s..G., Hahn’ C. J., London J., Chervin R. M. & Jenne R. L. (1986) NCAR/TN-273 + STR, National Center for Atmospheric Research, Boulder, co. Weber J. W. (1994) Univ. of Maryland UM-AERO 94-01. Weber J. W., Anderson J. D., Oran E. S.. Patnaik G. & Whaley R. (1994) High Performance Computing Conference 94; Grand Challenges in Computer Simulation, La Jolla, CA. Whitten G. Z., Hogo H. & Killus J. P. (1980) Enoiron. Sci. Tech. 14, 69&700. Wofsy S. C. (1978) J. Geophys. Res. 83, 364-378. Wuebbles D. J. & Kinnison D. E. (1989) A two-dimensional model study of past trends in global ozone. In Ozone in fhe Atmosphere (Edited by R. D. Bojkov and P. Fabian). Deepak, Hampton, VA. Yanenko N. A. (1971) The Method of Fractional Steps. Springer-Verlag, Germany. 1 Zhao X. P.. Turco R. P.. Kao C. Y. J.. Elliott S. & Newman M. (1993) EOS Trans. AGU, Fall ‘Meeting. Zhao X. P.. Turco R. P. & Lin R. (1994) Submitted to J. Geophys. Res. Zimmerman P. R., Chatfield R. B., Fishman J., Crutzen P. J. & Hanst P. L. (1978) Geophys. Res. Lew 5, 679689.