Model Description and User Guide

for WACCM1b

 

F. Sassi, B.A. Boville, B. Eaton,

R.R. Garcia, and R.G. Roble

 

National Center for Atmospheric Research

Boulder, CO, USA

 

9 September 2003

1. Introduction

The Whole Atmosphere Community Climate Model (WACCM) has been developed at the National Center for Atmospheric Research (NCAR) in Boulder, Colorado, as a collaboration among three NCAR Divisions: the Atmospheric Chemistry Division (ACD), the Climate and Global Dynamics division (CGD), and the High Altitude Observatory (HAO). The first goal of the project was to build a General Circulation Model (GCM) from the ground to the lower thermosphere based on the physics package of the NCAR Community Climate Model, version 3 (CCM3) but without interactive chemistry. The user of this note is referred to Kiehl et al [1996] for a description of CCM3 and its parameterizations. Although WACCM is based on CCM3 physics, it has been revised to use the software framework of the Community Atmosphere Model, version 2 (CAM2), which takes advantage of parallel scalar architectures. CAM2 is the atmospheric component of the Community Climate System Model; a User’s Guide for CAM2 can be found at http://www.ccsm.ucar.edu/models/atm-cam/UsersGuide/.

The current version of WACCM (based on CCM3 physics but using the software framework of CAM2) is known as WACCM1, version b (WACCM1b). The purpose of this note is to document the changes and additions that were made to the CCM3 physics package in order to produce WACCM1b, and to explain the use of the model. WACCM1b is “frozen” as of the time of its release. Because the physics package of CCM3 and is no longer maintained, no further maintenance or upgrades of WACCM1b will take place. The model described in this note is guaranteed to run on the NCAR platforms available at the time of its release. No guarantee is made that the model will run on other platforms.

2. Model Domain and Resolution

WACCM1b has 66 vertical levels from the ground to 5.1 x 10-6 hPa. As in CAM2, the vertical coordinate is purely isobaric above 100 hPa, but it is hybrid below that level. At any model grid point, the local pressure p is determined by

                                                Eq. (2.1)                          (2.1)

where a and b are functions of model level k only, and ps is the surface pressure, which is a function of model longitude and latitude (indexed by i and j).

Eq. (2.1) notwithstanding, in the remainder of this note we refer to the vertical coordinate in terms of log-pressure altitude:

                                                Eq. (2.2)                                                   (2.2)

where p0 is a reference pressure (=1000 hPa), p is the atmospheric pressure, and H=7 km. The value adopted for the scale height, H, is representative of the real atmosphere up to ~100 km; above that altitude temperature increases very rapidly and the typical scale height becomes correspondingly larger. It is important to distinguish z* from the geopotential height, which is obtained from integration of the hydrostatic equation.

In terms of log-pressure altitude, the model top level is found at z*=140 km. It should be noted that the solution in the top 15-20 km of the model is not expected to be realistic because this range of altitude is treated as a “sponge layer” to absorb upward propagating waves. The standard vertical resolution is variable; it is 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). The standard horizontal resolution is T63, with 128 longitude by 64 latitude points in a quasilinear grid [Williamson, 1997]. The dynamical equations are solved using a 2 time-level semi-Lagrangian scheme based on Williamson and Olson [1994], with a time step of 1800 s.

3. Upper Atmosphere Parameterizations

The following, new parameterizations are added to the CCM3 physics package in order to represent correctly the middle and upper atmosphere when the vertical domain is extended to z*=140 km:

  1. In addition to the parameterization of orographic (stationary) gravity waves, a parameterization of a spectrum of traveling gravity waves is implemented. Both parameterizations are based on Lindzen’s [1981] formulation. Aspects of this parameterization that are different from the description in Kiehl et al [1996] are illustrated in Section 3.1.
  2. The CCM3 longwave radiation code does not include non-LTE effects, which become important in mesosphere. These are taken into account by adopting the parameterization of Fomichev et al [1998]. This parameterization calculates longwave heating, including non-LTE effects, in the 15-µm band of CO2 and the 9.6-µm band of O3. Details of this procedure can be found in Section 3.2.
  3. In the upper mesosphere and lower thermosphere, absorption of solar radiation at wavelengths <200 nm by the Schumann-Runge bands and continuum of O2 is the primary shortwave heating mechanism. Since CCM3 was not designed to study the upper atmosphere, it calculates solar heating longward of 200 nm only. Therefore, in WACCM1b the CCM3 shortwave heating rates are merged with heating rates provided by the NCAR TIME-GCM [Roble and Ridley, 1994]. The procedure is similar to that used for longwave heating and is described in greater detail in Section 3.3.
  4. Molecular viscosity is an important process that influences the dynamical and thermal structure of the atmosphere above about z*=100 km. In WACCM1b, a parameterization based on the treatment of Banks and Kockarts [1973] is used. Details of the implementation of this parameterization in WACCM1b are described in Section 3.4.
  1. In the thermosphere, ionized species become abundant and interact with the geomagnetic field resulting in a force (“ion drag”) acting on the neutral wind that becomes significant above z*~120 km. Details of the implementation of ion drag in WACCM1b are discussed in Section 3.5.

The remainder of this Section describes each of the physical parameterizations used in WACCM1b in more detail.

3.1 Gravity waves

WACCM1b includes parameterizations for orographic gravity waves and for a spectrum of “non-orographic” waves. The non-orographic waves are parameterized in a manner similar to the orographic component, which is described in Kiehl et al [1996] (Section 4.1.e); the interested reader is referred to that publication for details of the mathematical derivations.

A spectrum of non-orographic gravity waves is launched at 100 hPa; it includes components with phase velocities in the range [–40, +40] m s-1, at intervals of 10 m s-1. The phase velocities are oriented in the direction of the wind at the source level. The stress at the source level has a Gaussian dependence on phase speed,

                                                Eq. (3.1)                                       (3.1)

where c is the phase speed, and =3.2×10-3 Pa is the stress at the source corresponding to a stationary wave (c = 0). The momentum source (3.1) is independent of time and longitude, but has a meridional structure defined by:

                        Eq. (3.2)                        (3.2)

where is the model latitude in degrees, =40o, and =10o. The constraint >0.2 ensures that there is a non-zero momentum flux in the tropics.

Above the source level, the saturation condition is calculated following Lindzen [1981], and the stress profile is damped using a total effective diffusivity:

                                                Eq. (3.3)                              (3.3)

where u is the local wind component in the same direction as the wind at the source level, c is the wave phase velocity, N is the buoyancy frequency, Fc is the critical Froude number (=0.5), k is the zonal wavenumber (=100 km), H is the scale height, is a rate of thermal dissipation (=10-6 s-1), and Dm is the molecular diffusivity of momentum.

Orographic gravity waves are not permitted to deposit momentum above 30 hPa. In early tests of the gravity wave parameterization, with the non-stationary spectrum and the orographic component both operating throughout the entire model vertical domain, the orographic component produced unrealistically large decelerations of the zonal winds in the upper stratosphere and mesosphere. This suggests that either the orographic source function or the saturation hypotheses are inadequate. Both will be replaced in the next version of WACCM.

WACCM1b does not include transport of heat and constituents (see Section 3.6) by breaking waves.

3.2 non-LTE longwave cooling

The details of the non-LTE parameterization can be found in Fomichev et al [1998]. Only those aspects pertaining to its implementation in WACCM1b are discussed here.

Fomichev’s parameterization can be applied above z*=10 km. However, the CCM3 longwave parameterization is accurate below z*=65 km, and includes the effects of additional gases (see Sect. 3.6). Therefore:

                                                Eq. (3.4)                                 (3.4)

where QLTE is the CCM3 longwave cooling, QNLTE is the cooling from the Fomichev parameterization, Qmerged is the resulting cooling, and is

                                                Eq. (3.5)                                              (3.5)

where z is the model height in km, z0 = 60 km, and dz=5 km.

This procedure is similar to that used by Fomichev and Blanchet [1995], except that they merged the heating rates in the middle stratosphere.

Fomichev’s non-LTE cooling scheme does not include nitric oxide cooling, which is an important term in the heat budget above z*=120 km. However, as noted above, the range of altitudes above 120 km are treated as a sponge layer to absorb upward propagating waves; therefore, disregard of nitric oxide cooling is not a major omission in the present model.

Fomichev’s parameterization requires the distribution of several constituents (O, O2, O3, N2, and CO2). The distributions of these species are obtained from HAO’s Thermosphere-Ionosphere-Mesosphere Electrodynamic GCM (TIME–GCM), as instantaneous global fields at the middle of each month. These distributions are spatially and temporally interpolated to the WACCM1b grid and time of the year, and are cycled every year. A discussion and additional information on the TIME–GCM can be found at http://www.hao.ucar.edu/public/research/tiso/tgcm/tgcm.html.

3.3 Upper atmosphere shortwave heating

WACCM1b does not calculate the distribution of chemical constituents necessary to compute the heating at wavelengths shorter than 200 nm. Therefore, we use heating fields from the TIME-GCM above z*~65 km, where shortwave heating at wavelengths <200 nm becomes important and cannot be calculated by the CCM3 parameterization.

The procedure is similar to what is described in Sect. 3.2 for the longwave component:

The TIME-GCM auroral heating component, which is prominent at high latitudes above about 100 km, has not been included in WACCM1b.

3.4 Molecular viscosity

The formulation of molecular diffusion in WACCM1b follows Banks and Kockarts [1973]. Momentum and temperature are diffused assuming a mean atmosphere where nitrogen and oxygen remain the major constituents with constant mixing ratio throughout the model's range of altitudes. This assumption is not completely justified near the top boundary of the model, but the top 15-20 km serve primarily to control upward wave propagation.

Diffusion of momentum obeys the law,

                                                                      (3.6)

where u is the zonal wind, is the atmospheric density, g is the acceleration of gravity, p is the pressure, and is the kinematic viscosity, which is obtained from an empirical formula:

                                                                                                    (3.7)

Note that we have used hydrostatic equilibrium, , to obtain the second equality in (3.6).

Diffusion of temperature is obtained in a similar fashion:

                                                                    (3.8)

where T is temperature, and is the thermal diffusivity which is obtained from the kinematic viscosity:

                                                                                                                    (3.9)

Pr is the Prandtl number, its value determined by Eucken’s formula,

                                                                                                                  (3.10)

and is the ratio of specific heats,

                                                                                                                         (3.11)

Using values appropriate for a well-mixed atmosphere, =1.4 and Pr =3/4.

Although WACCM1b calculates the distribution of several greenhouse gases (see section 3.6), it does not implement molecular diffusion of these constituents.

3.5 Ion drag

A parameterization of the so-called ion drag, based on the work of Dickinson et al [1975], is employed in WACCM1b. It is based on drag coefficients that are globally uniform, and does not account for the displaced position of the magnetic pole with respect of the geographic pole.

The ion drag tensor is defined as

                                                                                                        (3.12)

where and are globally uniform values,

                                                                                                        (3.13)

and d is the dip angle, which is obtained directly from the model latitude:

                                                                                                   (3.14)

Then, the horizontal wind tendency is obtained from product of the ion drag tensor to the horizontal wind components, that is

                                                                                                                (3.15)

where V=(u,v) is the horizontal wind vector and is the vector of horizontal wind tendencies. The matrix product in Eq. (3.16) yields the following set of horizontal tendencies:

                                                                                         (3.16)

3.6 Water vapor and other greenhouse gases

As in CCM3, WACCM1b uses the greenhouse gases option: H2O, CH4, N2O, CFC11, and CFC12 are advected constituents, with loss frequencies specified from the Garcia-Solomon (GS) two-dimensional model [Garcia et al, 1992; Garcia and Solomon 1994], as functions of latitude, height and month. The molecular hydrogen released by methane loss is assumed to be converted instantaneously to water vapor.

In addition to the processes that are available from CCM3, WACCM1b also includes Lyman-alpha photolysis, which is a major sink of water vapor above the stratopause. Photolysis rates are taken from the GS model. The effect OH chemistry, which leads to production of water vapor in the mesosphere and lower thermosphere, is included crudely by specifying a spatially and seasonally-varying source of water vapor from the results of the GS model. The water vapor production from this source is bounded by requiring that total hydrogen (H2O + CH4) does not exceed the value entering through the tropical tropopause. It is important to bear in mind that this procedure does not ensure conservation of total molecular hydrogen, and cannot be considered a substitute for a complete chemistry module. Although the resulting distribution of water vapor is fairly realistic and useful for radiative calculations, the water vapor budget above the stratopause should not be expected to be correct.

4. Where to get the code and how to run WACCM1b at NCAR

The source code, initial and boundary datasets, as well as scripts to run the model are available from the WACCM website.

Download the files:

·      waccm1b_src_code.tar.gz, which contains the distribution source code.

·      waccm1b_inpd.tar.gz, which contains boundary and initial datasets necessary to run the model.

If running on one of the NCAR IBM platforms, we recommend the following procedure:

  1. Move the file to a directory where you wish to unpack the source code:

·      Un-zip the source code tar-file using the UNIX command

gunzip waccm1b_src_code.tar.gz

This will create the file waccm1b_src_code.tar.

·      Un-tar the file with the command

tar xvf waccm1b_src_code.tar

This will create a local subdirectory named cam1, and place within it all the subdirectories containing the WACCM1b code. It is strongly recommended that files under the cam1 directory not be changed, so that the original source code is preserved intact.

  1. Move to the directory where you wish to unpack and store the datasets:

·      Un-zip the datasets tar-file,

gunzip waccm1b_inpd.tar.gz

This will create the file waccm1b_inpd.tar.

·      Un-tar the file with the command

tar xvf waccm1b_inpd.tar

This will create the subdirectory inputdata/atm/waccm containing all the necessary boundary and initial datasets.

  1. Create a directory separate from the original distribution code (which resides in cam1). You may give any name you wish to this directory, but it is customary to name it after the particular calculation, or “case”, being run. In what follows we assume that the case directory has been named tst_waccm1b.
  2. Within tst_waccm1b create a new directory named src and copy into it any WACCM1b subroutines that you wish to modify; make your changes to these copies (not to the original source!)
  3. Copy into test_waccm1b the file cam1/models/atm/cam/bld/run-ibm.csh. This is a sample script that can be used to build and run WACCM1b on the IBM platforms at NCAR.
  4. Edit the script run-ibm.csh as follows:

·      Change the account_no information as appropriate.

·      Search for camroot and replace that setting with the complete path to the cam1 directory (where the original source code has been stored).

·      Search for CSMDATA and replace that setting with the complete path to the inputdata directory (where the input datasets have been stored).

·      The script will find the necessary datasets and initial files. It will also compile the code, storing the executable in directory $blddir; output logs and files are written to $rundir. These two directories are subdirectories of $wrkdir, the working directory for the model run; they are defined automatically by the script depending on the setting of $wrkdir. If the model is not being run at NCAR, the user will need to change the default setting of $wrkdir, as appropriate. If the model is run at NCAR, the default setting may be used.

·      The model is run by executing run-ibm.csh. Note that the script needs to know at configuration time where the netCDF library is located. If running at NCAR, the following two settings need to be added to the script, or, preferably, to one’s login shell:

setenv LIB_NETCDF /usr/local/lib32/r4i4

setenv INC_NETCDF /usr/local/include

·      The script creates a namelist file in $rundir. The necessary namelist variables to run the model as it is distributed are assigned default values when the model is built. Specific namelist variables that are relevant to the WACCM1b user are:

a.     cftgcm contains the path and filename of the TIME-GCM boundary dataset to be used by Fomichev’s parameterization and to calculate shortwave heating rates. The default value of this namelist variable is $CSMDATA/atm/waccm/waccm.12mos.nc.

b.     itgcmcyc is a numeric flag that tells the model what use should be made of the file in cftgcm :

itgcmcyc = 1, twelve once-a-month samples are expected. Interpolation at the current model time step is carried out. The dataset is recycled every year. This is the default setting.

itgcmcyc = 0, it is assumed that the same dataset is used at every model time step.

itgcmcyc = 2, only two datasets are provided and temporal interpolation takes place between those two.

c.     ncdata contains the name of the initial file for the atmospheric component of the model. Default value is ./cam2.i.1959-07-01-00000.nc.

d.     finidat contains the name of the initial file for the land component of the model. Default value is ./lsmi_19590701_00000.nc. In place of the initial file name, finidat accepts the string ‘arbitrary initialization’ which forces the model to create its own surface vegetation characteristics, based on the horizontal resolution. Datasets necessary for arbitrary initialization are provided in the inputdata directory.

If the model is not being run at NCAR, the instructions given above can still serve as a guide for arranging the model code, initial and boundary datasets, etc. The build/run script run-ibm.csh can also be adapted to the user’s specific circumstances.

 

References

Acker, T.L., L. Buja, J.M. Rosinski, and J.E. Truesdale. User’s Guide to NCAR CCM3. NCAR/TN-421+IA, 1996.

Banks, P. M., and G. Kockarts, Aeronomy, part B, 355 pp., Academic, San Diego, Calif., 1973.

Dickinson, R. E., E. C. Ridley, and R. G. Roble, Meridional circulation in the thermosphere, I, Equinox conditions, J. Atmos. Sci., 32, 1737–1754, 1975.

Fomichev, V. I., and J. P. Blanchet, Development of the New CCC/GCM longwave radiation model for extension into the middle atmosphere, Atmos. Ocean, 33, 513–531, 1995.

Fomichev, V. I., J. P. Blanchet, and D. S. Turner, Matrix parameterization of the 15µm CO2 band cooling in the middle and upper atmosphere for variable CO2 concentration, J. Geophys. Res., 103, 11,505–11,528, 1998.

Garcia, R.R., F. Stordal, S. Solomon, and J.T. Kiehl, A new numerical model of the middle atmosphere. 1. Dynamics and transport of tropospheric source gases. J. Geophys. Res., 97, 12,967-12,991, 1992.

Garcia, R.R. and S. Solomon, A new numerical model of the middle atmosphere. 2. Ozone and related species, J. Geophys. Res., 103, 12,937-12,951, 1994.

Lindzen, R. S., Turbulence and stress due to gravity wave and tidal breakdown, J. Geophys. Res., 86, 9701–9714, 1981.

Kiehl, J.T., J.J. Hack, G.B. Bonan, B.A. Boville, B.P. Briegleb, D.L. Williamson, P.J. Rasch. Description of the NCAR Community Climate Model (CCM3). NCAR/TN-420+STR, 1996.

Roble, R. G., and E. C. Ridley, A Thermosphere-Ionosphere-Mesosphere-Electrodynamics General Circulation Model (TIME-GCM): Equinox solar cycle minimum simulations (30–500km), Geophys. Res. Lett., 21, 417–420, 1994.

Williamson, D. L., Climate simulations with a spectral, semi-Lagrangian model with linear grids, in Numerical Methods in Atmospheric and Ocean Modelling: The André J. Robert Memorial Volume, edited by C. Lin, R. Laprise, and H. Ritchie, pp. 279–292, Can. Meteorol. and Oceanogr. Soc., Ottawa, Ont., Canada, 1997.

Williamson, D. L., and J. G. Olson, Climate simulations with a semi-Lagrangian version of the NCAR Community Climate Model, Mon. Weather Rev., 122, 1594–1610, 1994.