Simulation of secular trends in the middle atmosphere, 1950–2003

2.1. Domain and Resolution

[8] WACCM3 is a global model with 66 vertical levels from the ground to 4.5 × 10⁻⁶ mbar (approximately 145 km geometric altitude). As in CAM3, the vertical coordinate is purely isobaric above 100 mbar, but is hybrid below that level. The vertical resolution is variable: 3.5 km above about 65 km, 1.75 km around the stratopause (50 km), 1.1–1.4 km in the lower stratosphere (below 30 km), and 1.1 km in the troposphere (except near the ground where much higher vertical resolution is used in the planetary boundary layer).

[9] WACCM3 currently supports two standard horizontal resolutions: 1.9° × 2.5° and 4° × 5° (latitude × longitude). The simulations presented in this paper, which encompass the 54-year period 1950–2003 and place very large demands on computational resources, have been carried out at 4° × 5° resolution. At all resolutions, the time step is 1800 s for the physical parameterizations. Within the finite volume dynamical core only, this time step is subdivided as necessary for computational stability.

2.2. Gravity Wave Parameterization

[10] WACCM3 incorporates a parameterization for a spectrum of vertically propagating internal gravity waves based on the work of Lindzen [1981], Holton [1982], Garcia and Solomon [1985], and Sassi et al. [2002]. Orographically generated gravity waves follow the parameterization of McFarlane [1987]. Both the orographic and spectral components of the parameterization take into account the rapid increase with altitude of molecular diffusion, which leads to diffusive separation and becomes the principal dissipation mechanism for upward propagating waves. Details of the implementation of these parameterizations in WACCM3 are given in Appendix A.

2.3. Molecular Diffusion

[11] Molecular diffusion is included in WACCM3 using the formulation of Banks and Kockarts [1973]. Enhanced molecular diffusivity suppresses the breaking of parameterized gravity waves above about 100 km, where wave dissipation occurs mainly via this process. Molecular diffusion also leads to diffusive separation at altitudes where the mean free path becomes large. Since WACCM3 extends only into the lower thermosphere, we avoid the full complexity of the diffusive separation problem by representing the diffusive separation velocity for each constituent with respect to the usual dry air mixture used in the lower atmosphere (mean molecular weight of 28.97 g mol⁻¹).

2.4. Chemistry

[12] The WACCM3 chemistry module is derived from the three-dimensional (3-D) chemical transport Model for Ozone and Related chemical Tracers (MOZART) [Brasseur et al., 1998; Hauglustaine et al., 1998; Horowitz et al., 2003; http://gctm.acd.ucar.edu/mozart]. It solves for 51 neutral species, including all members of the O_X, NO_X, HO_X, ClO_X, and BrO_X chemical families, along with tropospheric “source species” such as the N₂O, H₂O, CH₄, chlorofluorocarbons (CFCs) and other halogenated compounds, etc. Nonmethane hydrocarbons are excluded from this middle atmosphere mechanism, but several ion species important in the MLT (N₂⁺, O₂⁺, N⁺, NO⁺ and O⁺, plus electrons) are taken into account. Heterogeneous processes on sulfate aerosols and polar stratospheric clouds (liquid binary sulfate, supercooled ternary solutions, nitric acid trihydrate, and water ice), as well as aerosol sedimentation, are represented following the approach of Considine et al. [2000]. In almost all cases the chemical rate constants are taken from JPL02-25 [Sander et al., 2003]. A complete listing of species and reactions is given by D. E. Kinnison et al., Sensitivity of chemical tracers to meteorological parameters in the MOZART3 chemical transport model, submitted to Journal of Geophysical Research, 2006, hereinafter referred to as Kinnison et al., submitted manuscript, 2006.

[13] The calculation of photolysis rates in WACCM3 is divided into two regions: 120–200 nm (34 wavelength intervals) and 200–750 nm (67 wavelength intervals). The photolysis rate for each absorbing species is calculated during model execution as a function of the exoatmospheric flux, the atmospheric transmission function, the molecular absorption cross section, and the quantum yield. Details are given by Kinnison et al. (submitted manuscript, 2006). The exoatmospheric flux over the model wavelength intervals is parameterized in terms of the solar 10.7 cm radio flux (f10.7) following Solomon and Qian [2005] for wavelengths shortward of Lyman α, and Woods and Rottman [2002] for wavelengths between Lyman α and 350 nm. Beyond 350 nm, the flux is parameterized by regressing the difference between the total solar irradiance data of Froelich [2000] and the integrated flux up to 350 nm onto the 10.7 cm radio flux.

2.5. Longwave and Shortwave Heating

[14] WACCM3 retains the longwave (LW) formulation used in CAM3 [Kiehl and Briegleb, 1991]. However, modeling of the mesosphere and lower thermosphere requires a suite of LW parameterizations that deal with NLTE of the 15 μm band of CO₂ [Fomichev et al., 1998] and cooling due to NO at 5.3 μm [Kockarts, 1980]. The LW heating/cooling rates produced by these parameterizations are merged smoothly at 65 km with those produced by the standard CAM3 LW code, as recently revised by Collins et al. [2002].

[15] Shortwave (SW) heating in the CAM3 formulation employs the δ-Eddington approximation longward of 200 nm [Briegleb, 1992]. At altitudes higher than ∼70 km, radiation of shorter wavelength must also be included in WACCM3. Heating shortward of 200 nm is obtained from the same wavelength-dependent photolysis module used in the chemistry solver. The bond dissociation energy is subtracted for each O₂ and O₃ photolytic pathway, leaving only localized thermal heating. The additional energy is stored as chemical potential energy and realized later through 24 exothermic reactions, or lost as airglow through the 762 nm O₂(¹Σ) and 1.27 μm O₂(¹Δ) emission lines [Mlynczak and Solomon, 1993].

[16] Solar energy deposition in the EUV (shortward of Lyman α) and X-ray region is handled in a manner similar to longer wavelength ultraviolet radiation, with the spectrum divided into moderate-resolution bands and ionization, dissociation, and heating rates calculated in each band as a function of altitude [Solomon and Qian, 2005]. At EUV wavelengths, energy partitioning is complicated by photoionization, which generates energetic photoelectrons that, in turn, cause additional ionization, dissociation and heating, and become particularly important in the lower ionosphere. WACCM3 uses a high-resolution parameterization based upon the 1-D photoelectron model of Solomon and Qian [2005] to calculate heating rates due to photoelectrons.

[17] The SW heating rates calculated as described above are merged with those obtained with the CAM3 scheme at approximately 65 km. As in the case of photolysis, all heating rates are scaled by the wavelength-dependent exoatmospheric flux.

2.6. Auroral Processes, Ion Drag, and Joule Heating

[18] An auroral parameterization based on existing code from NCAR's Thermosphere-Ionosphere-Mesosphere Electrodynamics General Circulation Model (TIME-GCM) [Roble and Ridley, 1987] has been developed for rapid calculation of the total auroral ionization rate, particle precipitation in the polar cusp, and general polar cap precipitation (“polar drizzle”). The parameterization takes as input the hemispheric power (HP) of precipitating auroral electrons, and outputs total ionization rates and neutral heating. HP itself is parameterized as a function of the K_p geomagnetic index [Maeda et al., 1989], which is allowed to vary based upon observations. Once ionization rates are determined, the production rates for the E region ions N₂⁺, O₂⁺, N⁺, NO⁺ and O⁺ are calculated. Auroral production of NO can then be determined from the reaction of molecular oxygen and N(²D), the latter produced through dissociative recombination and charge exchange.

[19] The effects of momentum forcing by ion drag and of Joule heating associated with electric fields, which are particularly important above 110 km at high geomagnetic latitudes, are implemented in WACCM3 following Dickinson et al. [1981] and Roble et al. [1982], respectively. These models require knowledge of the Earth's electric field, which is parameterized according to the model of Weimer [1995] for high latitudes and that of Richmond et al. [1980] at low and middle latitudes. The Weimer model uses the interplanetary magnetic field as an input; this is estimated in WACCM3 from K_p, which, as in the case of the aurora, is allowed to vary according to observations.

2.7. Boundary Conditions

[20] The upper boundary conditions for momentum and for most constituents are the usual zero flux conditions used in CAM. However, in the energy budget of the thermosphere, much of the SW radiation at wavelengths <120 nm is absorbed above 145 km (the upper boundary of the model), where LW radiation is very inefficient. This energy is transported downward by molecular diffusion to below 120 km, where it can be dissipated more efficiently by LW emission. Imposing a zero flux upper boundary condition on heat omits a major term in the heat budget and causes the lower thermosphere to be much too cold. Instead, we use the Mass Spectrometer-Incoherent Scatter (MSIS) model [Hedin, 1987, 1991] to specify the temperature at the top boundary as a function of season and phase of the solar cycle. The particular version of the MSIS model used in WACCM3 is NRLMSISE-00 (see http://uap-www.nrl.navy.mil/models_web/msis/msis_home.htm).

[21] For chemical constituents, surface mixing ratios of CH₄, N₂O, CO₂, H₂, CFC-11, CFC-12, CFC-113, HCFC-22, H-1211, H-1301, CCl₄, CH₃CCH₃, CH₃Cl, and CH₃Br are specified from observations. The model accounts for surface emissions of NO_X and CO based on the emission inventories described by Horowitz et al. [2003]. The NO_X source from lightning is distributed according to the location of convective clouds based on Price et al. [1997a, 1997b] with a vertical profile following Pickering et al. [1998]. Aircraft emissions of NO_X and CO are included in the model and based on Friedl [1997].

[22] At the upper boundary, a zero-flux upper boundary condition is used for most species whose mixing ratio is negligible in the lower thermosphere, while mixing ratios of other species are specified from a variety of sources. The MSIS model is used to specify the mixing ratios of O, O₂, H, and N; as in the case of temperature, the MSIS model returns values of these constituents as functions of season and phase of the solar cycle. CO and CO₂ are specified at the upper boundary using output from the TIME-GCM [Roble and Ridley, 1994]. NO is specified using data from the Student Nitric Oxide Explorer (SNOE) satellite [Barth et al., 2003], which has been parameterized as a function of latitude, season, and phase of the solar cycle in Marsh et al.'s [2004] Nitric Oxide Empirical Model (NOEM). Finally, a global mean value (typical of the sunlit lower thermosphere) is specified for species such as H₂O, whose abundance near the top of the model is very small under sunlit conditions, but which can be rapidly transported upward by diffusive separation in polar night (since they are lighter than the background atmosphere). In these cases, a zero flux boundary condition leads to unrealistically large mixing ratios at the model top in polar night.

2.8. Specification of Boundary and Initial Conditions for 1950–2003

[23] The boundary conditions in this study are based upon, but not identical to, the specifications for the first reference case (REF1) used in the model intercomparison exercise of Eyring et al. [2005, 2006]. These specifications include surface mixing ratios for GHGs defined by scenario A1B of IPCC [2001]; surface mixing ratios for halogen compounds taken from Table 4B-2 of WMO [2003]; monthly mean sea surface temperatures (SSTs) from the UK Met's Hadley Centre data set; chemical and radiative effects of volcanic aerosols; and 11-year solar cycle irradiance variability parameterized in terms of observed f10.7 radio flux.

[24] In our simulations, the surface area density (SAD) of sulfate aerosols is derived from satellite observations by the Stratospheric Aerosol and Gas Experiment (SAGE, SAGE II) and the Stratospheric and Mesospheric Sounder (SAMS), as described by Thomason et al. [1997] and updated by D. B. Considine [WMO, 2003]. Daily observations of f10.7 (and also of the K_p geomagnetic index) were obtained from the Space Environment Center of the U.S. National Oceanographic and Atmospheric Administration (NOAA) (http://www.sec.noaa.gov). Additional details on the REF1 reference case are given by Eyring et al. [2005, 2006].

[25] The WACCM3 calculations differ from the REF1 specification in several respects: SST are prescribed from the global HadISST data set prior to 1981 and from the Smith/Reynolds data set after 1981 [Hurrell et al., 2006]; a QBO is neither generated spontaneously by the model nor specified externally; heating from volcanic aerosols is not included (although the chemical effects thereof are taken into account, as noted above); chemical kinetics follow JPL02-25 [Sander et al., 2003], as noted in 5; and, in addition to solar cycle variations in photolysis and heating, WACCM3 also calculates changes in ion and NO production in the aurora, and changes in ion drag and Joule heating, as explained in 7.

[26] Note that the treatment of the effect of volcanic aerosols in these model calculations is incomplete in that heating due to absorption of solar radiation by the aerosols is neglected. We did not include aerosol heating because we lacked a suitable parameterization thereof at the time the model runs were begun. We were particularly concerned about the effects of heating at the tropical cold point, which, if not accurately modeled, can lead to unrealistically large water vapor mixing ratios in the air entering the stratosphere. Once in the stratosphere, this excess water can persist for years, and can affect ozone chemistry through catalysis by the HO_X family. On balance, we decided it was preferable not to include aerosol heating than to include a heating distribution that might cause the aforementioned problems. Thus any effects of volcanic eruptions on temperature, tropical circulation, and water vapor in the lower tropical stratosphere due to aerosol heating of the lower stratosphere are not included in these runs.

[27] Three realizations of the period 1950–2003 were carried out using the boundary conditions described above. The realizations start from an equilibrated initial state for 1950, which was obtained by integrating the model for at least 10 years with fixed boundary conditions and solar inputs appropriate for 1950. Independent realizations are obtained by introducing small perturbations in the equilibrated initial state.

4.1. Temperature

[41] Figure 5 shows zonal mean temperature trends for the stratosphere (K/decade) calculated from monthly mean WACCM3 output (Figure 5). The temperature trend is smallest near 70–100 mbar (∼16–17 km) and increases with altitude up to about 1 mbar (∼45 km), where it reaches values of −1.25 to −1.50 K/decade, with minor variations in latitude, except at high latitudes of the Southern Hemisphere. Over the southern polar cap the model computes a strong negative trend of more than −2.5 K/decade centered around 18–20 km, in the region of the ozone hole. In the Northern Hemisphere lower stratosphere, WACCM3 does not produce a significant temperature trend poleward of 70°. These results are broadly consistent with the CCM calculations of Austin and Butchart [2003] for the period 1980–1999, and with several of the models discussed by Shine et al. [2003]. However, they differ from observations [e.g., Pawson and Naujokat, 1999] and certain recent modeling results [e.g., Lahoz, 2000; Braesicke and Pyle, 2004; Dameris et al., 2005] in that Arctic winters in the 1990s are not especially cold in the lower stratosphere, even though the model is driven by observed SSTs. The modeling studies cited suggest that specification of observed SSTs produces Arctic stratospheric temperatures that are cold in the 1990s, in agreement with observations. This behavior is not present in the WACCM3 simulations [cf. Eyring et al., 2006, Figure 4], and contributes to the lack of a significant temperature trend in the Arctic lower stratosphere; it also has consequences for the calculated ozone trends in the Arctic, as discussed below.

[42] Figure 6 compares WACCM3 temperature trends for 1979–2003, averaged between 60°S and 60°N, with 1979–2004 trends computed from satellite data, and from radiosonde observations. The observed trends are obtained from linear regression upon time, omitting the two years after the volcanic eruptions of El Chichón and Mount Pinatubo. The model trends make no allowance for volcanic eruptions because heating by volcanic aerosols was not included in the simulations, as noted in 9. The satellite observations are from the stratospheric sounding unit (SSU) for altitudes between about 20 km and the stratopause, and from the microwave sounding unit (MSU) channel 4 in the lowermost stratosphere. Radiosonde results are from a subset of stations between 60°S and 60°N, described by Lanzante et al. [2003] and updated as described by Randel and Wu [2006]; this subset is chosen to omit stations with large artificial cooling biases, in particular those for which differences between MSU channel 4 and radiosonde trends are greater than 0.3 K/decade. WACCM3 trends are shown individually for each of the model realizations to illustrate the internal variability of the model. The 2σ error bars for the three model realizations overlap at all altitudes; differences with respect to SSU/MSU data are significant around 35 km (∼6 mbar) and near the stratopause (∼1 mbar). The reason for these discrepancies is not known. Trends calculated with other recent CCMs tend to bracket our results. For example, Austin and Butchart [2003] obtained global trends of about −1.6 K/decade at 1 mbar and −0.8 K/decade at 6 mbar for the period 1980–1999, which are similar to the results from WACCM3, whereas Langematz et al. [2003] calculated a trend for 1980–2000 of nearly −2.5 K/decade at 1 mbar, which is actually larger than the SSU/MSU trend. At mesospheric altitudes there are few observations to compare with WACCM3, but results from CCMs that extend beyond the stratosphere are generally consistent with those shown in Figure 6 [see, e.g., Shine et al., 2003, Figure 4].

[43] Attribution of temperature trends to different factors (ozone decrease, increases in GHGs), as done for certain models by Shine et al. [2003], cannot be carried out with WACCM3 because the model has been run with interactive chemistry and radiation, which makes it impossible to separate the influences of each. However, calculations carried out with an earlier, noninteractive version of the model (not shown) yielded conclusions in line with those discussed by Shine et al. [2003], namely, that in the upper and lowermost stratosphere the cooling trend is dominated by the effect of ozone loss (see 13), while in the middle stratosphere (∼10 mbar) the effect of GHGs is most important.

[44] Figure 7 shows the temperature trend for the entire atmosphere up to 135 km calculated from the three model realizations for the entire period of simulation, 1950–2003. The morphology of the trend in the stratosphere is very similar to that of the trend shown in Figure 5, except that the magnitude is smaller. In the lower thermosphere (above 100 km) large trends are computed, peaking at −2.5 K/decade near 120 km. Interestingly, above that altitude the trends become smaller, which appears to be the result of the increasing dominance above about 125 km of IR cooling by the 5.3 μm emission of NO, a gas whose abundance remains essentially constant through the period 1950–2003. It bears repeating here that all model results, including those shown in Figure 7, are displayed in (isobaric) log pressure altitude. This should be kept in mind when comparing these results against observations made at geometric altitudes, especially in the thermosphere (above ∼100 km [see, e.g., Akmaev and Fomichev, 2000]).

[45] Near the mesopause, at 80–90 km, the calculated temperature trend is either insignificant or very small. The lack of a temperature trend in a range of altitude where CO₂ is the main infrared emitter is puzzling, but appears to be consistent with available observations. For example, Beig et al. [2003] have compiled estimates of mesospheric temperature trends obtained from a variety of observations (ground-based, rocketsonde, satellite, etc.); they point out that the majority of the data sets examined, including the most reliable ones, show no significant temperature trend in this range of altitude. Note that the decrease in the geometric altitude of isobaric surfaces due to cooling of the atmosphere over the period 1950–2003 is less than 800 m near the mesopause, so temperature trends at 80–90 km would be essentially the same as shown in Figure 7 had they been calculated at constant geometric altitude.

[46] The reason for the lack of a temperature trend near the mesopause in the WACCM3 simulations is the subject of current investigation. However, Schmidt et al. [2006] have recently used the HAMMONIA CCM in a 2 × CO₂ experiment to show that the temperature change is small, and even statistically insignificant, at many locations near the mesopause. They attribute this behavior to compensating changes in dynamical heating by the mean meridional circulation. In contrast, Manzini et al. [2003] used the MAECHAM/CHEM model to perform time slice simulations for 1960 and 2000 conditions, and derived a temperature decrease between 1969 and 2000 of −4 K in the upper mesosphere, which is much larger than obtained by us or by Schmidt et al. This is perhaps due to the fact that the top boundary in MAECHAM/CHEM is located at 0.01 mbar, i.e., near the mean altitude of the mesopause, which may preclude compensating effects by the mean meridional circulation.

4.2. Ozone

[47] Figure 8 compares calculated ozone trends with results obtained from satellite measurements by the Stratospheric Aerosol and Gas Experiment (SAGE I and SAGE II), supplemented with ozonesonde observations at high latitudes [Randel and Wu, 2007]. The observations are regressed upon ESSC, QBO, and solar cycle indices, whereas model results omit regression on the QBO, which is absent in WACCM3. The SAGE data cover the latitude range ±55°, from 20 to 50 km for SAGE I and from the tropopause to 50 km for SAGE II; ozonesonde observations are available in the polar regions at Syowa (69°S) and Resolute (75°N) from the tropopause to 30 km. Note therefore that the values shown in Figure 8b above 30 km in the polar regions are extrapolated from lower latitudes. Note also that the SAGE/ozonesonde results are expressed as the net change over the period 1979–2005, whereas WACCM3 trends with respect to EESC are computed for 1979–2003, since the simulations end in 2003. Finally, because WACCM3 results are expressed in percent change per unit of EESC, the values in Figure 8a should be multiplied times 1.5 (the change in EESC between 1979 and 2003) in order to compare them with the SAGE/ozonesonde net changes shown in Figure 8b.

[48] Model results the observations are consistent in most regions of the stratosphere. Thus the largest ozone trends in the upper stratosphere are found around 40 km, and are about −8% per unit of EESC (or −12% change over the period 1979–2003), which agrees well with the change derived from the SAGE/ozonesonde data set. In both WACCM3 and the data there are regions of slight, albeit statistically insignificant, ozone increase in the tropics, at ∼25 km and immediately above the tropopause, and a region of small trends between 25 and 30 km in extratropical latitudes whose significance is marginal. In the data there is region of strong negative trends in the tropics, centered near 18 km, which is not present in the WACCM3 results; however, as noted by Randel and Wu [2007], SAGE trends in this region are of questionable validity.

[49] In the region of the ozone hole, WACCM3 calculates a maximum trend of −28% per unit of EESC at ∼17 km (or −42% between 1979 and 2003, smaller than the −60% obtained from the Syowa ozonesonde data at ∼15 km). On the other hand, negative trends over Antarctica extend throughout the stratosphere in WACCM, whereas those computed from the data are actually positive (although statistically insignificant) between 25 and 30 km. In the Arctic lower stratosphere, the discrepancy is larger: WACCM3 does not produce a significant trend, whereas the data indicate a net change of about −8%. This is consistent with the temperature results shown in Figure 5, where WACCM3 trends at high latitudes are insignificant poleward of 70°N. It is possible that observed trends in both temperature and ozone are influenced by the series of very cold Northern Hemisphere winters that occurred in the mid-1990s [cf. Eyring et al., 2006, Figures 4 and 15], behavior that is not present in the WACCM3 results, as noted above. On average, WACCM3 zonal mean temperatures in the Northern Hemisphere polar cap in winter are warmer than observed by about 5 K, which can affect the efficiency of chlorine activation and, consequently, both ozone depletion and the trend thereof.

[50] Figure 9 compares results for the total ozone column as a function of season, with monthly mean resolution. The data are from the Total Ozone Mapping Spectrometer (TOMS) and the Solar Backscattered Ultraviolet (SBUV) satellite instruments for 1979–2005, as recently analyzed by Randel and Wu [2007]. The ozone data are regressed upon ESSC, QBO, and solar cycle indices, and the results are expressed as net change in Dobson units (DU) over the period in question; WACCM3 results omit regression on the QBO, are shown in terms of DU per unit of EESC, and are computed for 1979–2003. As in Figure 8, the WACCM3 numbers should be multiplied by 1.5, the change in EESC between 1979 and 2003, in order to compare with the data.

[51] WACCM3 column trends are generally consistent with observations: They are small in the tropics and increase toward high latitudes, maximizing during the Northern and Southern Hemisphere spring seasons, when ozone depletion is largest. The model trends over the southern polar cap in October are −65 DU/EESC unit, or −97 DU over the period 1979–2003, very similar to what is obtained from the data for 1979–2005. Note, however, the persistence of large model ozone trends into January (−50%/EESC unit, or a change of −75% over 1979–2003), whereas the percentage change in the data drops rapidly after November, to about −30% in January. This is a manifestation of the cold bias that develops in the model during southern spring in the polar lower stratosphere, where temperatures remain cold and westerlies persist into January (see 11).

[52] In the Northern Hemisphere, WACCM3 trends are considerably smaller than those observed; they are largest over the Arctic in February (−15 DU/EESC unit, or −22.5 DU for the period 1979–2003), but only −5 to −10 DU/EESC unit (−7.5 to −15 DU change from 1979 to 2003) poleward of 60°N in March, compared to as much as −20 DU at 60°N seen in the data for that month. Further, aside from the months of February and March, high-latitude trends in the Northern Hemisphere are very small (and statistically insignificant) in WACCM3, whereas the observed change remains larger than −8 DU throughout the entire period (April–July) when the satellite instruments are able to observe the northern polar cap. This is a reflection of the discrepancy between the ozone loss in the Arctic lower stratosphere calculated with WACCM3 and the larger observed loss, as noted in connection with Figure 8.

[53] The calculation of trends or net changes in ozone, as in Figures 8 and 9, provides a compact description of the evolution of ozone in the last few decades of the 20th century. However, because ozone changes in this period did not occur at a constant rate, it is useful to examine the secular evolution of the ozone column to see whether WACCM3 can capture the salient features of the observational record. Figure 10 shows the behavior of the ozone column anomaly since 1964 in WACCM3 and in observations, averaged between ±60° (Figure 10, top) and globally (Figure 10, bottom), as reported by WMO [2003]. In all cases the anomalies are calculated with respect to the mean ozone column for the period 1964–1980, as was done in the WMO [2003] Ozone Assessment [see also Fioletov et al., 2002]. The data come from a variety of sources (TOMS, SBUV and ground-based instruments), and exhibit a high degree of consistency among data sets. For WACCM3, the results are shown for each of the three realizations and for their average. The dates of the eruptions of El Chichón and Mount Pinatubo are indicated in the plots.

[54] The evolution of the ozone column anomalies derived from WACCM3 agrees in most respects with the observations. In particular, the magnitude of the decrease is very similar for the globally averaged column (−5% in both cases), although slightly smaller in WACCM3 (−4%) than in the data (−4.5 to −5%) when averaged over ±60°. A sharp dip is seen in the model and the data after the eruption of Mount Pinatubo (June 1991), with column ozone reaching minimum values toward the end of 1992 whether averaged globally or over ±60°. This behavior has been attributed to the effect of the aerosol load introduced into the stratosphere by the eruption [Brasseur and Granier, 1992; Hofmann et al., 1994]. Note, on the other hand, that there is little indication of a large column decrease after El Chichón (spring of 1982) in either the model or the observations. Note also that in both model and observations, ozone column anomalies show no consistent decrease in the period after about 1994. This is consistent with the fact that the value of EESC is nearly constant after 1994–1995, so no ozone change would be predicted by our ozone-EESC regression (Figure 8a). Whether this represents the beginning of ozone recovery cannot be answered definitively by the present calculations. However, WACCM3 calculations of the period 2000–2050, to be presented elsewhere, indicate that recovery, in the sense of a sustained increase in ozone column, does not begin until approximately 2005. Thus the behavior calculated (and observed) since the mid-1990s is perhaps best characterized as a “slowdown” of ozone loss [Newchurch et al., 2003].

[55] The behavior of the global ozone column in one of the WACCM3 realizations, shown in gold in Figure 10, is remarkable in that the anomaly averaged over ±90° is almost zero in model year 1991 (the two other realizations have anomalies of about −2 to −3% with respect to the 1964–1980 average, which is typical of the early 1990s). This behavior is not seen in the column anomalies averaged over ±60°, which indicates it must be due to the contribution from the polar regions. Further examination reveals that the behavior can be isolated to the southern polar cap, where 1991 is a highly anomalous model year, with column ozone amounts in Antarctic spring (not shown) reaching values similar to those calculated for the 1960s and 1970s, before the development of the Antarctic ozone hole.

[56] The disappearance of the ozone hole in model year 1991 is reminiscent of the behavior observed in Antarctica in 2002, when a major stratospheric sudden warming virtually eliminated the ozone hole in September [see, e.g., Charlton et al., 2005; Hio and Yoden, 2005; Krüger et al., 2005; Richter et al., 2005; Roscoe et al., 2005; Stolarski et al., 2005] However, the behavior in WACCM3 is different in several important respects: Disturbances of the Antarctic vortex due to large-amplitude planetary wave events occurred early in the winter season (as early as July) and several times thereafter; this prevented the normal development of cold temperatures, the activation of catalytic chlorine species, and the formation of the ozone hole. In addition, the zonal mean zonal wind, although exhibiting very large negative anomalies with respect to climatology (as much as −55 m s⁻¹), did not meet at any time during winter or spring the usual criterion for a major sudden warming (easterlies present at 10 mbar and 60°). These aspects of model year 1991 are similar to the Antarctic winter of 1988, as reported by Kanzawa and Kaguchi [1990] and Hirota et al. [1990], although that winter appears to have been somewhat less disturbed (large planetary wave events occurred later in the season, starting in August 1988 rather than in July, and they reduced but did not eliminate the ozone hole). A description of the unusual model year of 1991 will be presented elsewhere.

[57] Figure 11 shows the whole atmosphere ozone trend over the entire period of simulation, 1950–2003. Note that because the period shown begins well before the time of largest anthropogenic emissions of chlorine and bromine compounds, and because trends outside the stratosphere are due to factors others than chlorine/bromine catalysis of ozone, the results in Figure 11 are presented as actual temporal trends (in percent change per decade) rather than as changes per unit of EESC, as was done in Figures 8 and 9. In the troposphere, stratosphere and lower mesosphere the morphology of the trends is very similar to that obtained for 1979–2003 (Figure 8), although changes are smaller, reflecting the fact that the main factors that influence ozone, temperature and halogen abundance, have changed most rapidly in the last 2–3 decades of the 20th century. At higher altitudes, there is a negative trend up to ∼100 km as a result of increasing water vapor (and hence HO_X, which dominates ozone destruction in the mesosphere); the water vapor increases are attributable, in the model, to increases in methane (as shown in 14). In the lower thermosphere ozone actually increases (by as much as 6–8% per decade near 120 km), as a result of the large, negative temperature trend in this region (compare Figure 7). Finally, in the troposphere there is a net increase in ozone of about 2–3% per decade that can be attributed to increases in methane, whose oxidation leads to the production of ozone [Crutzen, 1973].

[58] As noted at the beginning of this section, WACCM3 trends for ozone are obtained from multiple regressions that include the 10.7 cm solar flux as a predictor, so it is appropriate to make note of the sensitivity displayed by WACCM3 to 11-year solar variability. Figure 12 shows the latitude dependence of the regression coefficient of column ozone on f10.7. At most latitudes the regression coefficient is between 2.5 and 3 DU per 100 units of f10.7. This is consistent with values derived from ground-based instruments and from BUV/SBUV satellite observations [WMO, 2003]; however, the values are substantially smaller than those obtained recently by Stolarski et al. [2006] using the Goddard Space Flight Center's chemistry transport model, and those derived by the same authors from combined TOMS/SBUV observations, which are about 4–6 DU per 100 units of f10.7 at most latitudes. The reason for this discrepancy is not clear, although it should be pointed out that the estimates shown in Figure 12 are based upon the entire 1950–2003 simulation period, whereas Stolarski et al.'s results are for 1979–2003; in addition Stolarski et al. also included volcanic aerosol loading as a predictor in their multiple regression, something that was not done in the WACCM3 analysis. It is also noteworthy that when regressions from WACCM3 output are calculated for the period 1979–2003, the regression coefficient for f10.7 (not shown) increases to about 4 DU per 100 units of f10.7.

4.3. Water Vapor

[59] While calculated temperature and ozone trends are generally consistent with observations, as shown in 12 and 13, this is not true of trends in water vapor, a constituent whose evolution has been carefully documented from hygrosonde observations made in Boulder, Colorado, over the last two decades, and in the satellite record provided by the HALOE instrument onboard UARS since 1992 [Randel et al., 2004b]. Figure 13 shows a comparison of the trends calculated with WACCM3 (Figure 13a) and from HALOE data (Figure 13b) for the period 1992–2002. The latter are obtained from regression upon time and a QBO index, as explained by Randel et al. While HALOE shows declining water in the lower stratosphere together with strong increases (as much as 6% per decade) in the middle and upper stratosphere, WACCM3 trends are of the opposite sign in the lower stratosphere (2–4% per decade increases in the tropics), and positive but substantially smaller than HALOE trends in the middle and upper stratosphere.

[60] Agreement between model and observations is no better for the Boulder hygrosonde data, available since 1980 [Oltmans et al., 2000], which are also computed from regression upon time and QBO index and are shown in Figure 14. These observations have been cited in recent works as evidence that stratospheric water vapor is undergoing a rapid increase (as much as 10% per decade), which may imply significant changes in tropical tropopause temperatures, or even dynamics [see, e.g., Randel et al., 2004b]. Comparison with WACCM3 trends over 1992–2002 shows that the latter are about one half, or even less, of those derived from the Boulder hygrosondes. Further, model trends are rather variable among the three realizations, which are shown separately in Figure 14. Also shown in Figure 14 is the 1992–2002 trend at the latitude of Boulder derived from HALOE data; this trend does not agree with either the hygrosondes or the model and, indeed, it is of the opposite sign below about 23 km.

[61] While it is possible that water vapor has indeed undergone large secular changes over the last 10–20 years, an alternative interpretation of the results shown in Figures 13 and 14 is that trends calculated over decadal timescales are unstable; that is, they are not indicative of long-term change but instead reflect the presence of low-frequency variability inherent in the behavior of water vapor. Such variability can resemble a trend, even a statistically significant one, but may disappear when a sufficiently long time series is analyzed. In support of this interpretation Figure 15 shows WACCM3 trends computed for two arbitrary 10-year periods, 1975–1985 and 1980–1990. In the first period, trends are positive throughout most of the stratosphere and exceed 6% per decade in certain regions; in the second, trends are substantially smaller in the upper stratosphere, and negative in the lower stratosphere. Only over Antarctica, where water vapor abundance is strongly influenced by dehydration associated with the ozone hole, are the trends similar between the two periods shown in Figure 15. It is clear that the trends in either simulation cannot be representative of long-term behavior but instead must reflect the influence of low-frequency variability.

[62] For the WACCM3 simulations, it is necessary to compute trends over three or more decades in order to obtain stable results. When WACCM3 trends are calculated for the entire period of simulation 1950–2003, as shown in Figure 16, a stable pattern emerges in the stratosphere and mesosphere, where the long-term trend can be attributed to the increase in methane between 1950 and 2003 (about 0.6 ppmv, which implies an increase of 1.2 ppmv in water vapor, sufficient to account for the maximum 4% per decade trend in water vapor seen in the upper stratosphere and mesosphere). In addition, a negative trend is calculated over Antarctica even in this long record, indicative of the growth of the ozone hole (and thus of colder temperatures) in the last 20 years; and a positive trend is calculated in the troposphere, as a result of the increase in relative humidity allowed by the small but significant tropospheric temperature trend due to the greenhouse effect (compare Figure 7).

[63] These results raise the question what are the sources of low-frequency variability in the behavior of water vapor in the middle atmosphere. Because water vapor in the middle atmosphere is very strongly influenced by the temperature of the entry region at the tropical cold point tropopause, it is to be expected that any processes that affect the cold point temperature will exert a strong influence on water vapor abundance. The most obvious such processes that operate on timescales of several years are the QBO (through the adiabatic effect of downwelling and upwelling associated with the QBO secondary circulation) [Plumb and Bell, 1982; Giorgetta and Bengtsson, 1999; Randel et al., 2004b]; volcanic eruptions (which heat the tropopause region when solar radiation is absorbed by volcanic aerosols) [Stenchikov et al., 2002; Joshi and Shine, 2003]; and El Niño–Southern Oscillation (ENSO, which modifies upper tropospheric and lower stratospheric temperatures) [Geller et al., 2002; Calvo Fernández et al., 2004; Fueglistaler and Haynes, 2005].

[64] The effect of ENSO on the simulation of water vapor in WACCM3 is illustrated in Figure 17, which shows the behavior of water vapor anomalies in the tropics from 1950 to 2003 (Figure 17a), and their correlation with the Niño-3.4 (N3.4) index as a function of altitude and time lag (Figure 17b). N3.4 is a conventional indicator of ENSO intensity based upon sea surface temperature anomalies in the region 170°W–120°W, 5°S–5°N. It is clear from Figure 17a that enhanced water vapor in the troposphere often accompanies the occurrence of ENSO events, and this is followed by an increase in the lowermost stratosphere that propagates upward at the speed of the tropical “tape recorder.” Figure 17b shows that this effect is reflected in the lag correlation between water vapor anomalies and N3.4, which is as large as 0.7 in the troposphere a few months after the maximum of N3.4; in the stratosphere the correlation maximizes at lags that increase with time (as expected from transport), and reach values between ±0.2 and ±0.4. These results are consistent with those of Scaife et al. [2003], who documented the impact of ENSO events on stratospheric water vapor using the UK Met's Unified Model.

[65] Because WACCM3 does not generate a QBO, and because heating by volcanic aerosols was not included in the present simulations, it is not surprising that short-term water vapor trends, such as those shown in Figures 13 and 14, do not agree with observations. It also appears that successful simulation of observed water vapor trends will require at a minimum the inclusion of a realistic QBO and aerosol heating rates in addition to ENSO effects. Even then, it may be necessary to examine trends over perhaps as long as three decades before an unambiguous attribution can be made, a conclusion consistent with the findings of Fueglistaler and Haynes [2005]. Finally, it is also apparent that the trends derived from HALOE and shown in Figure 14, which represent a zonal mean at the latitude of Boulder, cannot be reconciled with the local hygrosonde data. This suggests that the behavior of water in the lower stratosphere over Boulder may be influenced by local processes not captured in the HALOE observations, or that there are unknown errors in one or both sets of observations.

4.4. Stratospheric Age of Air

[66] The stratospheric age of air (AOA) is a measure of the strength of the stratospheric circulation obtained by comparing the mixing ratio of a steadily increasing “age of air tracer” at some reference point to the mixing ratio elsewhere in the stratosphere. The reference point is usually chosen to be in the upper tropical troposphere, in the vicinity of the entry region of air into the global stratosphere [see Hall and Plumb, 1994; Hall et al., 1999]. Specifically, the AOA may be defined as

where t_χ(y, z) is the time at which a certain mixing ratio, χ, of the AOA tracer is reached at some location (x, y) in the meridional plane, and t_χ(y₀, z₀) is the (earlier) time when the same mixing ratio occurs at the reference point (x₀, y₀). AOA is determined from observations of long-lived, steadily increasing trace gases with sinks present only at very high altitudes, like CO₂ or SF₆. In WACCM3 we use an ad hoc, conservative AOA tracer whose concentration increases linearly in time with a constant surface flux.

[67] Figure 18 shows the evolution of τ_A(1 mbar), averaged over ±22 ° for each of the three WACCM3 realizations. Note that AOA is only shown starting in 1963 because the tracer is initialized to zero everywhere in the model domain at the start of each simulation and, as a consequence, it takes about a dozen years before its mixing ratio throughout the meridional plane equilibrates to values representative of the AOA. The results are remarkably consistent among the realizations and indicate that in the 40 years displayed in Figure 18, τ_A(1 mbar) decreases by about 4 months, from a little under 4 years to a bit more than 3.5 years, or about 8.25% with respect to its initial value.

[68] The strengthening of the stratospheric (Brewer-Dobson) circulation in response to increases in GHGs has been documented previously by Hansen et al. [2005], who compared an ensemble of General Circulation models. Similar results for AOA have been obtained by Austin and Li [2006] with the AMTRAC model of the Geophysical Fluid Dynamics Laboratory, although in that model the decrease of AOA from 1960 to 2005 at 1 mbar (averaged between ±20°) is about 8 months, about twice the value than we obtain. This would seem to imply a more sensitive response of the global circulation to climate change in the last half of the 20th century in AMTRAC compared to WACCM3, the reasons for which remain moot at this time. In any case, we show below that the AOA change in WACCM3 is consistent with the trend in tropical upwelling.

[69] Figure 19 shows the trend in the tropical average (±22°) of the zonal mean vertical velocity as a function of altitude for the combined ensemble of WACCM3 simulations, along with individual time series at selected levels. We use the conventional Eulerian zonal mean, , which is not expected to be significantly different from the Transformed Eulerian mean in the tropics [Andrews et al., 1987]. The model results have been smoothed with a 12-month running mean to remove short-term variability and emphasize the secular behavior. The WACCM3 ensemble of shows significant positive trends at all altitudes between 15 and 50 km; typical values between 20 and 35 km are a bit less than 5 × 10⁻⁶ m s⁻¹ per decade, with somewhat larger values near the tropopause and a much larger increase between 40 and 45 km. (In all cases, these trends are <10% of the time mean at the corresponding altitude.) A very simple test of consistency between these results and the AOA changes shown in Figure 18 can be carried out as follows: Ignoring lateral mixing in the “tropical pipe” region [Plumb, 1996] near the equator, we assume that the AOA at some distance, ΔZ, above the tropopause is given approximately by the traveltime,

so that

At 1 mbar (∼45 km), the distance above the tropopause is ΔZ ∼ 30 km, while δ ∼ 2 × 10⁻⁵ m s⁻¹, or about 2 m d⁻¹. The last number is arrived at by taking an estimate of 0.5 × 10⁻⁶ m s⁻¹ per decade as representative of the long-term trend in throughout much of the stratosphere and multiplying times 4 decades (the period over which the AOA shown in Figure 18 is computed). With these numbers, and a value of τ(1 mbar) ∼4 years, or about 1400 days, we have from (4)

which is consistent with the 4 month decrease in AOA seen in Figure 19. This result, which is also consistent with the findings of Austin and Li [2006], suggests that the decrease in AOA can indeed be interpreted as a strengthening of the global circulation insofar as the tropical average of is a measure of the global upwelling in the stratosphere.

[70] One might ask whether such changes in the global stratospheric circulation play a role in the secular trends calculated by WACCM3. A simple way to approach the question is to ask how changes in the flux of trace species due to an enhancement of the circulation compare to changes induced by other processes, in particular by trends in the tropospheric mixing ratio of the tracer in question. If the tropical average of the advective flux of χ entering the stratosphere is given by

then fractional changes in the flux may be expressed as

In WACCM3, the fractional change attributable to changes in tropical upwelling, δ/, is less than 0.1, while the fractional change due to changes in mixing ratio, δ/, can be quite large over the period 1950–2003 for CFCs and other halogens; even for methane and CO₂, δ/ has values of 0.5 and 0.2, respectively. Thus it appears that the strengthening of the stratospheric circulation documented in Figure 19 plays a minor role in altering the stratospheric abundance of trace gases whose tropospheric sources have undergone large changes in the period of simulation.

[71] There remains the question whether changes of 5–10% in the strength of tropical upwelling might influence temperatures in the tropical cold point region and therefore the water vapor mixing ratio of air entering the stratosphere. According to Figure 19, the net change in at the tropical tropopause over the period 1950–2003 is ∼5 × 10⁻⁵ m s⁻¹. From the thermodynamic equation one can estimate the change in temperature implied by a change in vertical velocity as

Using values α = 0.01 d⁻¹ [Randel et al., 2001] for the radiative relaxation rate and Γ = 2 K km⁻¹ for the lapse rate in the upper tropical troposphere in equation (7) gives δT = −1 K. However, the temperature trends over 1950–2003, shown in Figure 7, do not indicate a temperature change at the model's cold point (85 mbar, or ∼17 km) nearly as large as this. The model cold point temperature change between 1950 and 2003 averaged over ±22° is smaller than −0.25 K, and its statistical significance is marginal. Similarly, water vapor trends in the lower stratosphere (Figure 16) do not indicate a change commensurate with a 1 K decline in temperature (the average water vapor change over the period 1950–2003, averaged over ±22 is nearly zero, and statistically insignificant in any case). This suggests that other processes (e.g., changes in the radiative budget due to changes in GHGs, or changes in tropical convection) overwhelm the impact of changes in tropical upwelling on the cold point temperature (and hence on the water vapor mixing ratio of air entering the stratosphere) in the WACCM3 calculations.