the PMIP models (Bonfils et al. 1998)
PMIP Documentation for CSIRO
Commonwealth Scientific and Industrial
Research Organization: Model CSIRO V4-7 R21L9 199?
Dr. Jozef Syktus, CSIRO Division of Atmospheric Research, Private Bag No
1, Aspendale, 3195 Victoria Australia, Phone: 61 3 9239 4548; Fax: 61 3
9239 4444; e-mail: Jozef.Syktus@dar.csiro.au
WWW URL: http://www.dar.csiro.au
CSIRO V4-7 R21L9 199?
Model Identification for PMIP
Number of days in each month: 31 28 31 30 31 30 31 31 30 31 30 31
The CSIRO model is derived from earlier two-level and four-level spectral
models based on the primitive equations expressed in conservative flux
form (cf. Gordon 1981 ,1993 ).
H. B. Gordon and S. P. O'Farrell, Transient Climate Change in the CSIRO
Coupled Model with Dynamic Sea Ice , 1997, Mon. Wea. Rev., 125 : 875-907
McGregor, J. L., H. B. Gordon, I. G. Watterson, M. R. Dix, and L. D.
Rotstayn, 1993, The CSIRO 9-level AGCM. CSIRO Dar Tech. Paper No. 26, CSIRO,
PMB1, Aspendale, Victoria 3195, Australia, 89 pp.
Watterson, I. G., Dix, M. R., Gordon, H. B., and McGregor, J. L. (1995).
The CSIRO nine-level atmospheric general circulation model and its equilibrium
present and doubled CO2 climates. Australian Meteorological Magazine, 44
Documentation of the present version of the CSIRO model is provided
by Gordon and O'Farell (1997).
Spectral (spherical harmonic basis functions) with transformation to a
Gaussian grid for calculation of nonlinear quantities and some physics.
The atmospheric moisture field is represented only in gridded form.
Spectral rhomboidal 21 (R21), roughly equivalent to 3.2 x 5.6 degrees latitude-longitude.
Surface to about 21 hPa. For a surface pressure of 1000 hPa, the lowest
atmospheric level is at about 979 hPa.
Finite-difference sigma coordinates.
There are 9 unevenly spaced sigma levels with the following levels: 0.0171
0.0802 0.1927 0.3381 0.5 0.6619 0.8073 0.9198 0.9829. For a surface pressure
of 1000 hPa, 3 levels are below 800 hPa and 3 levels are above 200 hPa.
The PMIP simulation was run on a Cray Y/MP computer using one processor
in the UNICOS environment.
For the PMIP experiment, about 35 seconds Cray Y/MP computation time per
Initial conditions of PMIP model for control simulation were taken from
previous integration of model with fixed SST run for few years. Initialization
of the atmosphere, soil moisture, and snow cover/depth for 1 January 1979
is from an earlier model simulation with climatological sea surface temperatures.
Time Integration Scheme(s)
A semi-implicit leapfrog time scheme with an Asselin (1972) frequency filter
is used for most calculations, with the momentum surface flux and vertical
diffusion above the surface computed by split backward implicit integration.
A time step of 30 minutes is used for dynamics and physics, except for
full calculations of all radiative fluxes and heating rates, which are
done every 2 hours.
Orography is truncated at the R21 resolution of the model (see Orography).
To counter the negative values of atmospheric moisture that may otherwise
develop, vertical transport of moisture is inhibited if the local water
vapor mixing ratio drops below 2 x 10-6 kg (water) per kg (air).
In addition, negative moisture values are removed by a proportional adjustment
method while conserving the global mean (cf. Royer 1986 ). Cf. McGregor
et al. 1993 for further details.
For the PMIP simulation, the model history is written every 6 hours.
Primitive-equation dynamics are expressed in conservative flux form (i.e.,
weighting vorticity, divergence, temperature, and specific humidity by
the prognostic surface pressure) as described by Gordon (1981) . Effects
of frictional heating are included in the temperature tendency equation,
and virtual temperature is used to compute geopotential height.
Linear second-order (Ñ2) horizontal
diffusion of the temperature, vorticity, divergence, and moisture fields
is computed via a split implicit time integration. Diffusion (with coefficient
106 m2/s) is applied to the upper-half of the rhomboid
for the spectral temperature, vorticity, and divergence. Diffusion of the
gridded moisture is calculated from a temporary spectral representation
of the moisture field without surface pressure weighting, and is applied
to the entire rhomboid with the diffusion coefficient halved. Diffusion
for temperature and moisture contains a first-order correction to constant
Vertical diffusion of momentum, heat, and moisture is parameterized
in terms of stability-dependent K-theory following Blackadar (1962) , with
the choice of asymptotic mixing length after Louis (1979) . The calculation
of vertical momentum diffusion is via backward implicit time differencing
(see Time Integration Scheme(s)).
Under conditions of vertical stability, orographic gravity-wave drag is
simulated after the method of Chouinard et al. (1986) . The drag at the
surface is dependent on sub-gridscale orographic variance (see Orography),
and it is parameterized by means of a "launching" height which is defined
to be twice the local standard deviation of the surface heights. Following
Palmer et al. (1986) , the maximum launching height is limited to 800 m
in order to prevent two-grid noise near steep mountains. At a particular
sigma level the frictional drag on the atmosphere from breaking gravity
waves depends on the projection of the wind on the surface wind and on
the Froude number, which in turn is a function of the launching height,
the atmospheric density, the Brunt-Vaisalla frequency, and the wind shear.
Gravity-wave drag is assumed to be zero above a critical level, which is
taken to be the top sigma level of the model (see Vertical Domain).
The solar constant is the AMIP-prescribed value of 1365 W/(m2).
The orbital parameters and seasonal insolation distribution are calculated
after PMIP recommendations. Both seasonal and diurnal cycles in solar forcing
The carbon dioxide concentration is the AMIP-prescribed value of 345 and
280ppm for 0fix and 6fix respectively. Ozone concentrations, specified
as a function of latitude and pressure, are interpolated from the Dopplick
(1974) seasonal climatology. Radiative effects of water vapor, but not
of aerosols, also are included (see Radiation).
The radiation code is after Fels (1985) , Fels and Schwarzkopf (1975 ,
1981 ), and Schwarzkopf and Fels (1991). The shortwave calculations are
based on a modified Lacis and Hansen (1974) approach. The shortwave spectrum
is divided into 9 bands, the first band covering the ultraviolet (wavelengths
0.1 to 0.4 micron) and visible (0.4 to 0.7 micron) spectral intervals,
while the other 8 bands are in the near-infrared (0.7 to 4.0 microns).
Rayleigh scattering by air molecules and absorption by ozone and water
vapor are treated in the first band. In the 8 near-infrared bands, variable
absorption by water vapor is included, and carbon dioxide absorption is
calculated after a modified Sasamori et al. (1972) method. Pressure corrections
and multiple reflections between clouds and the surface are treated, but
not the radiative effects of aerosols.
Longwave calculations follow the simplified exchange method of Fels
and Schwarzkopf (1975) and Schwarzkopf and Fels (1991) applied over seven
spectral bands (with wavenumber boundaries at 0, 4.0 x 104,
5.6 x 104, 8.0 x 104, 9.9 x 104, 1.07
x 105, 1.20 x 105, and 2.20 x 105 m-1).
Absorption by the vibrational and rotational lines of water vapor, carbon
dioxide, and ozone, as well as continuum absorption by water vapor are
treated, but some weak absorption bands of ozone and carbon dioxide are
neglected. Carbon dioxide transmission coefficients are calculated for
the actual temperature and pressure profile of each vertical column after
the interpolation method of Fels and Schwarzkopf (1981). Longwave ozone
and water vapor absorption (including temperature effects) are computed
by a random-band model.
Cloud optical properties are prescribed. In the visible, cloud absorptivity
is assumed to be zero; infrared absorptivity depends on cloud height, as
do the visible and infrared cloud reflectivities. The absorptivity and
reflectivity of clouds are also proportional to cloud amount. Longwave
emissivity is prescribed as unity (blackbody emission) for all clouds.
For purposes of calculating radiation and top-of-atmosphere cloud cover,
clouds in different vertical layers are assumed to be randomly overlapped.
Cf. McGregor et al. (1993) for further details. See also Cloud Formation.
After checking for supersaturation and attendant release of precipitation
(see Precipitation) and performing a dry convective adjustment if needed,
a modified Arakawa (1972) "soft" moist adjustment scheme predicts any subsequent
precipitation release and the redistribution of moisture and momentum that
may occur within sub-gridscale cumulus towers. These form when a layer
(other than the lowest layer) is moist unstable with respect to at least
one layer above, and when the relative humidity in the lowest unstable
layer is > 75 percent. It is assumed that a constant convective mass flux
effects a vertical redistribution of heat within the cumulus tower, such
that the moist instability at each level (the difference between the moist
static energy at cloud base and the saturation value at each level) decays
with an e-folding time of one hour. (The heating at cloud base is assumed
to be zero to ensure closure for the convective scheme.)
The attendant moisture redistribution and removal (see Precipitation)
results in drying of the environment if the ambient relative humidity is
at least 60 percent. The convective mass flux also transfers momentum upward
through the cumulus tower, and downward via the surrounding large-scale
descent. Shallow cumulus convection is parameterized by a modified Geleyn
(1987) scheme that operates as an extension of the vertical diffusion of
heat and momentum. The stability dependence of this diffusion is defined
by a modified moist bulk Richardson number. Cf. McGregor et al. (1993)
for further details.
Model description is to a large degree similar to AMIP database, however
there are few important changes in cloud formulation (based on NCAR CCM2
for stratiform clouds. If convective activity occurs in the previous timestep
(see Convection), convective cloud fraction is determined according to
a formula by Slingo (1987).
Large-scale cloud amounts are determined from a modified form of the
Rikus (1991) diagnostic, which is a quadratic function of the difference
between the relative humidity of a layer and critical humidities that depend
on cloud height (low, middle, and high cloud). Following Saito and Baba
(1988) , maximum cloud fractions are also specified for each cloud type
(0.70 for low cloud, 0.53 for middle cloud, and 0.50 for high cloud). Middle
and high clouds are restricted to single layers, while low clouds can occupy
Stability-dependent low cloud associated with temperature inversions
also may form if the relative humidity at cloud base is at least 60 percent.
The cloud fraction is a function of the intensity of the temperature inversion,
following Slingo (1987) . The overall fraction of low cloud is then set
to the largest value predicted by either this stability-dependent diagnostic
or by another operative mechanism (0.55 in the case of convective activity,
or 0.70 for large-scale condensation).
Precipitation forms as a result of supersaturation and/or the moist convective
adjustment process (see Convection). There is no subsequent evaporation
Planetary Boundary Layer
There are typically two atmospheric levels in the model PBL (whose top
is not explicitly determined, however). In order to validate the simulation
of surface atmospheric temperature against observations, a 2-meter (screen
height) temperature is calculated from the bulk Richardson number determined
for the lowest model layer (by applying the Monin-Obukhov assumption of
constant momentum and heat fluxes in the surface layer). For unstable conditions,
this requires an iterative solution. For purposes of computing surface
fluxes, however, atmospheric winds, temperatures, and humidities at the
first full atmospheric level above the surface (at sigma = 0.979) are used
(see Surface Fluxes). Cf. McGregor et al. (1993) .
Orography from the 1 x 1-degree data of Gates and Nelson (1975) is transformed
to spectral coefficients and truncated at the R21 resolution of the model.
Orographic variance data (supplied by the United Kingdom Meteorological
Office) are also used for the parameterization of gravity-wave drag (see
For control: PMIP datasets will be used for SSTs and sea-ice have been
prepared at PCMDI and were calculated by averaging the 10-year AMIP datasets
For 6 fix: SSTs and sea-ice prescribed at their present day value, as
in the control run.
For 0fix and 6fix run, a sea ice model is used. It is a fully dynamic model
described in Gordon and O'Farell 1997. Ocean currents are prescribed for
fixed SST in both cases. Sea Ice Mass derived from sea ice thicknes(m)
x concentration (fraction) x 900kg/m**3 (assumed density of sea ice).
SEAICE % is calculated from mean monthly concentration files.
Precipitation falls to the surface as snow if the temperature of the second
atmospheric vertical level above the surface (at sigma = 0.914) is below
0 degrees C. The latent heat is incorporated into the surface temperature
prognostic for non-ocean surfaces (see Sea Ice and Land Surface Processes).
Prognostic snow mass, with accumulation and melting over both land and
sea ice, is included, but the allowable snow depth is limited to 4 m. Snow
cover affects the heat capacity and conductivity of the land surface (see
Land Surface Processes), and sublimation of snow contributes to the surface
evaporative flux (see Surface Fluxes). Melting of snow, which contributes
to soil moisture, occurs when the surface (top soil layer or sea ice) temperature
is > 0 degrees C. Snow cover affects the albedo of the surface, but with
less impact if the snow is melting (see Surface Characteristics). No fractional
snow cover calculated in the model. Value of 3 mm of H2O was used as a
cutoff for snow/nosnow in calculating this variable from snow thickness
data. Snow mass. Maximum depth of snow is imposed in the model and is set
to 400 mm H2O (snow calculated in the CSIRO GCM as water equivalent). If
snow exceedes 400 mm it is turned to sea ice (over ice) and on land it
is assumed to vanish (glaciers flow etc), however GCM take into account
this into energy budget.
There are 12 vegetations types taken from Dorman and Sellers (1989). J.
Appl. Meteor., 28, 833-855. The roughness length over land is prescribed
as a monthly values from Dorman and Sellers (1989) dataset at each point.
However, the roughness length for other surface types is the same for heat
and momentum fluxes: over ice, the roughness is a uniform 0.001 m, while
over the ocean it is a function of surface wind stress, following Charnock
Surface albedo and surface charactersitics are prescribed as a monthly
values (and interpolated to a daily values) from Dorman and Sellers dataset.
Longwave emissivity is set to unity (blackbody emission) for all surface
types. Cf. McGregor et al. (1993) for further details.
Surface solar absorption is determined from surface albedos, and surface
longwave emission from the Planck equation with prescribed emissivity of
1.0 (see Surface Characteristics). Following Monin-Obukhov similarity theory
applicable to a constant-flux layer, surface turbulent eddy fluxes are
expressed as bulk formulae. The requisite atmospheric surface winds, potential
temperatures, and specific humidities are taken to be those at the first
vertical level above the surface (at sigma = 0.979). The effective ground
value of specific humidity needed for determination of the surface moisture
flux from the bulk formula is obtained as a fraction alpha of the saturated
humidity at the ground temperature, where alpha is a function of soil moisture
(see Land Surface Processes) or is prescribed to be 1 over oceans and ice.
The bulk drag and transfer coefficients are functions of roughness length
(see Surface Characteristics) and vertical stability (bulk Richardson number)
following Louis (1979) , with the same transfer coefficient used for the
heat and moisture fluxes. Over the oceans, the neutral transfer coefficient
for surface heat and moisture fluxes is a constant 8.5 x 10-4.
Land Surface Processes
There are 9 soil types is used from Zobler dataset. The soil uses following
parameters; saturated moisture content, willing moisture content, saturation
value of matrix potential, specific heat, density, albedo, emissivity,
and field capacity. Vegetation uses following parameters; min observed
conopy stomal resistance, fraction of vegetation, LAI, roughness lenght,
albedo and emissivity. SVAT model is after Kowalczyk et al. The model used
in PMIP experiment is exactly the same as in Gordon and O'Farell (1997),
except dynamic ocean is not used.
E. A. Kowalczyk and J. R. Garratt and P. B. Krummel, 1991, A Soil-Canopy
Scheme for Use in a Numerical Model of the Atmosphere - 1D Stand-Alone
Model, CSIRO Division of Atmospheric Research, Technical Paper, 23
E. A. Kowalczyk and J. R. Garratt and P. B. Krummel, 1994, Implementation
of a Soil-Canopy Scheme into the CSIRO GCMs, CSIRO Division of Atmospheric
Research, Technical Paper, 32.
Last update November 9, 1998. For further information, contact: Céline