the PMIP models (Bonfils et al. 1998)
PMIP Documentation for CCSR1
Center for Climate System Research:
Model CCSR/NIES 5.4 02 T21/L20 1995
Dr. Ayako Abe-Ouchi, Center for Climate System Research (CCSR), University
of Tokyo, 4-6-1 Komaba, Meguro-ku, Tokyo, 153 Japan, Phone: +81-3-5453-3955;
Fax: +81-3-5453-3964; e-mail: email@example.com
WWW URL: http://climate3.ccsr.u-tokyo.ac.jp/home.html
CCSR/NIES AGCM 5.4 02 (T21 L20) 1995
Model Identification for PMIP
0fix, 6fix, 21fix
Number of days in each month: 30 30 30 30 30 30 30 30 30 30 30 30
Model CCSR/NIES AGCM (T21 L20) 1995 is based on a simple global atmospheric
model first developed at the University of Tokyo (cf. Numaguti 1993), and
further refined as a collaboration between CCSR and the National Institute
of Environmental Studies (NIES). It is intended for use as a community
The model is identical to the latest AMIP model exept for different
initial conditions and the Earth's orbital parameters.
Numaguti and others (1996) Development of an Atmospheric General Circulation
Model.(prepared for submittion. available upon request)
A summary of model features including fundamental equations is provided
by Numaguti et al. (1995). The spectral formulation of atmospheric dynamics
follows closely Bourke (1988). The radiation scheme is described by Nakajima
and Tanaka (1986) and Nakajima et al. (1996). The convective parameterization
is based on the work of Arakawa and Schubert (1974) and Moorthi and Suarez
(1992). Cloud formation is treated prognostically after the method of Le
Treut and Li (1991). Gravity-wave drag is parameterized as in McFarlane
(1987). The planetary boundary layer (PBL) is simulated by the turbulence
closure scheme of Mellor and Yamada (1974, 1982)[10,11]. The representation
of surface fluxes follows the approach of Louis (1979), with inclusio of
adjustments recommended by Miller et al. (1992) for low winds over the
Spectral (spherical harmonic basis function) with transformation to a Gaussian
grid for calculation of nonlinear quantities and some physics.
Spectral triangular 21 (T21), roughly equivalent to a 5.6 x 5.6 degree
latitude/longitude grid. dim_longitude*dim_latitude: 64*32
Surface to 8 hPa. For a surface pressure of 1000 hPa, the lowest atmospheric
level is at a pressure of about 995 hPa.
Sigma coordinates with discretization following the vertical differencing
scheme of Arakawa and Suarez (1983) that conserves global mass integrals
of potential temperature and total energy for frictionless adiabatic flow
There are 20 unevenly spaced sigma levels. For a surface pressure of 1000
hPa, 5 levels are below 800 hPa and 8 levels are above 200 hPa. The model
uses a sigma coordinate in the vertical as: 1: 0.99500 2: 0.97999 3: 0.94995
4: 0.89988 5: 0.82977 6: 0.74468 7: 0.64954 8: 0.54946 9: 0.45447 10: 0.36948
11: 0.29450 12: 0.22953 13: 0.17457 14: 0.12440 15: 0.08468 16: 0.05980
17: 0.04493 18: 0.03491 19: 0.02488 20: 0.00830
The PMIP simulation was run on a HITAC S-3800 computer using a single processor
in the VOS3 operational environment.
For the PMIP experiment, about 0.3 minutes of HITAC S-3800 computation
time per simulated day.
Initial conditions of model atmosphere, atmospheric state, soil moisture,
and snow cover/depth were obtained from a AMIP run.
Time Integration Scheme(s)
Semi-implicit leapfrog time integration with an Asselin (1972) time filter.
The time step length is 40 minutes. Shortwave and longwave radiative fluxes
are recalculated every 3 hours, but with the longwave fluxes assumed constant
over the 3-hour interval, while the shortwave fluxes are assumed to vary
as the cosine of the solar zenith angle.
Orography is smoothed (see Orography). Spurious negative atmospheric moisture
values are filled by borrowing from the vertical level immediately below,
subject to the constraint of conservation of global moisture.
For the PMIP simulation, the model history is written once per 24-hour
Primitive equation dynamics are expressed in terms of vorticity and divergence,temperature,
specific humidity, cloud liquid water, and surface pressure.
Eighth-order linear (Ñ8) horizontal
diffusion is applied to vorticity, divergence, temperature, specific humidity,
and cloud liquid water on constant sigma surfaces. Stability-dependent
vertical diffusion of momentum, heat, and moisture in the planetary boundary
layer (PBL) as well as in the free atmosphere follows the Mellor and Yamada
(1974, 1982)[10,11 ] level-2 turbulence closure scheme. The eddy diffusion
coefficient is diagnostically determined as a function of a Richardson
number modified to include the effects of condensation. The diffusion coefficient
also depends on the vertical wind shear and on the square of an eddy mixing
length with an asymptotic value of 300 m. Cf. Numagati et al. (1995) for
further details. See also Planetary Boundary Layer and Surface Fluxes.
Orographic gravity-wave drag is parameterized after McFarlane (1987). Deceleration
of the resolved flow by dissipation of orographically excited gravity waves
is a function of the rate at which the parameterized vertical component
of the gravity-wave momentum flux decreases in magnitude with height. This
momentum-flux term is the product of local air density, the component of
the local wind in the direction of that at the near-surface reference level,
and a displacement amplitude. At the surface, this amplitude is specified
in terms of the mesoscale orographic variance, and in the free atmosphere
by linear theory, but it is bounded everywhere by wave saturation values.
See also Orography.
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, 280
and 200 ppm for 0fix, 6 fix and 21fix run, respectively. Radiative effects
of water vapor, oxygen, ozone, nitrous oxide (0.3 ppm, globally uniform),
and methane (1.7 ppm, globally uniform) are included. Monthly zonal ozone
profiles are specified from data of Keating and Young (1985) and Dütsch
(1978), and they are linearly interpolated for intermediate time points.
Although the model is able to treat radiative effects of aerosols, they
are not included for the PMIP simulation. See also Radiation.
The radiative transfer scheme is based on the two-stream discrete ordinate
method (DOM) and the k-distribution method described in detail by Nakajima
et al. (1996)).The radiative fluxes at the interfaces of each vertical
layer are calculated considering solar incidence, absorption, emission,
and scattering by gases and clouds, with the flux calculations being done
in 18 wavelength regions. Band absorption by water vapor, carbon dioxide,
ozone, nitrous oxide, and methane is considered in from 1 to 6 subchannels
for each wavelength region. Continuum absorption by water vapor, oxygen,
and ozone also is included. Rayleigh scattering by gases and absorption
by clouds are considered as well. See also Chemistry.
The radiative flux in each wavelength region is calculated as a sum
of the products of the fluxes over the subchannels and their respective
k-distribution weights, where each subchannel's fluxes are calculated by
the two-stream DOM. The optical depth of each subchannel is estimated as
the sum of the optical thicknesses of band absorption and of continuum
absorption by gases. The transmissivity, reflectivity, and source function
in each layer then are calculated as functions of optical depth, single-scattering
albedo, asymmetry factor, cutoff factor, Planck function, solar incidence,
and solar zenith angle. At each layer interface, fluxes are computed by
the adding technique.
In the presence of clouds, radiative fluxes are weighted according to
the convective and large-scale cloud fractions of each grid box. The fluxes
are computed by treating clouds as a mixture of scattering and absorbing
water and ice particles in the shortwave; shortwave optical properties
and longwave in the emissivity are functions of optical depth. Fluxes therefore
depend on the prognostic liquid water content (LWC) as well as on the fraction
of ice cloud (see Cloud Formation). Radiative transfer in large-scale and
convective clouds is treated separately, assuming random and full overlap,
respectively, in the vertical. Cf. Numaguti et al. (1995) for further details.
Penetrative and shallow cumulus convection are simulated by the Relaxed
Arakawa-Schubert (RAS) scheme of Moorthi and Suarez (1992) , a modification
of the Arakawa and Schubert (1974) parameterization. The RAS scheme predicts
mass fluxes from a spectrum of clouds that have different entrainment/detrainment
rates and levels of neutral buoyancy (i.e., different cloud-top heights).
The thermodynamic properties of the convective clouds are determined from
an entraining plume model and the vertical profile of cloud liquid water
(see Cloud Formation) is calculated from the difference between the adiabatic
total water mixing ratio (a function of the grid-scale specific humidity)
and the saturated specific humidity at the same level, given a prescribed
vertical profile of precipitation.
The predicted convective mass fluxes are used to solve budget equations
that determine the impact of convection on the grid-scale fields of temperature
(through latent heating and compensating subsidence) and moisture (through
precipitation and detrainment). The vertical mass flux at the base of each
cloud type is predicted from the cloud work function A, defined as the
integral over the cloud depth of the product of the mass flux (with a linear
vertical profile assumed) and the buoyancy (proportional to the difference
between the cloud virtual temperature and that of the grid-scale environment
at the same height). Because a nonzero cloud-base mass flux implies a positive-definite
work function, the former is determined assuming that the work function
vanishes in a specified time scale T > that is longer than the convective
time step t of 80 minutes.
In the RAS scheme, the new cloud-base mass flux at each time step is
estimated by the method of virtual displacement: the amount of grid-scale
warming and drying expected from a unit mass flux is calculated, and a
new cloud work function A' is determined; the new cloud-base mass flux
M' then is derived from a simple proportionality relation. The grid-scale
mass flux is obtained by summing over the contributions from the spectrum
of cloud types. The profile of the mass flux associated with convective
downdrafts is also simulated from a fixed fraction of the evaporation of
convective precipitation (see Precipitation). Updated values of convective
cloud fraction and convective liquid water content (LWC) for the grid box
also are determined from the grid-scale mass flux (see Cloud Formation
and Radiation). Cf. Numaguti et al. (1995) for further details.
The convective cloud fraction in a grid box is estimated as proportional
to the grid-scale convective mass flux. The grid-scale liquid water content
(LWC) at a given height due to convective cloud is determined by a sum
over the cloud-type spectrum of the products of LWC and mass flux for each
cloud type (see Convection) Large-scale (stratiform) cloud formation is
determined from prognostic cloud liquid water content (LWC) following Le
Treut and Li (1991) . The stratiform LWC follows a conservation equation
involving rates of large-scale water vapor condensation, evaporation of
cloud droplets, and the transformation of small droplets to large precipitating
drops (see Precipitation). The stratiform LWC (including ice content) also
determines the large-scale cloud fraction (see below) and cloud optical
properties (see Radiation).
The fraction of stratiform cloud C in any layer is determined from the
probability that the total cloud water (liquid plus vapor) is above the
saturated value, where a uniform probability distribution with prescribed
standard deviation is assumed. (For purposes of the radiation calculations,
the square root of C is taken as the cloud fraction). At each time step,
new values of LWC and vapor are determined by iteration, subject to conservation
of moist internal energy. The portion of C that is ice cloud is assumed
to vary as a linear function of the temperature departure below the freezing
point 273.15 K, with all of C being ice cloud if the temperature is <
258.18 K. Cf. Numaguti et al. (1995) and Le Treut and Li (1991) for further
The autoconversion of cloud liquid water into precipitation is estimated
from the prognostic liquid water content (LWC) divided by a characteristic
precipitation time scale which is an exponential function of temperature
(see Cloud Formation). Precipitation conversion is distinguished for liquid
vs ice particles. Snow is assumed to fall when the local wet-bulb temperature
is less than the freezing point of 273.15 K, with melting of falling snow
occurring if the wet-bulb temperature exceeds this value. See also Snow
Cover. Falling liquid precipitation evaporates proportional to the difference
between the saturated and ambient specific humidities and inversely proportional
to the terminal fall velocity (cf. Kessler (1969)). Falling ice and snow
melts if the ambient wet-bulb temperature exceeds the freezing point (273.15
K); evaporation may follow, as for falling liquid precipitation.
Planetary Boundary Layer
The Mellor and Yamada (1974, 1982) [10,11] level-2 turbulence closure scheme
represents the effects of the PBL. The scheme is used to determine vertical
diffusion coefficients for momentum, heat, and moisture from the product
of the squared mixing length (whose asymptotic value is 300 m), the vertical
wind shear, and a Richardson number that is modifed to include the effect
of condensation on turbulent fluxes. (A diffusion coefficient is never
allowed to fall below 0.15 m2/s.) Cf. Numaguti et al. (1995)
for further details. See also Surface Fluxes.
Raw orography is obtained from the ETOPO5 dataset (cf. NOAA/NGDC, 1989)
at a resolution of 5 x 5 minutes for 0fix and 6fix run. It is modified
after PMIP recommendation for 21fix run using the difference of the orographic
difference between present and 21ka. Orographic variances required for
the gravity-wave drag scheme are obtained from the same dataset. Orography
is smoothed by first expanding the grid point data in spectral space, then
filtering according to the formula [1-(n/N)4], where n is the
spectral wavenumber and N = 21 corresponds to the horizontal resolution
of the model. Finally, the smoothed spectral data is returned to the T21
Gaussian grid. See also Gravity-wave Drag.
The prescribed monthly climatological SSTs for 0fix and 6fix integrations
were made by averaging AMIP monthly sea surface temperature fields, with
daily values determined by linear interpolation.
For 21fix ? ? ? ?
For 0fix and 6fix run, monthly AMIP sea ice extents are prescribed. The
thickness of the ice can vary: the local thickness is determined from the
observed fractional coverage multiplied by a constant 1 m. The surface
temperature of the ice is predicted from a surface energy balance that
takes account of conduction heating from the ocean below. The temperature
of the underlying ocean is assumed to be 273.15 K, the freezing point of
the sea ice. Snow may accumulate on sea ice, but modifies only the thermal
conductivity of the ice. See also Surface Fluxes and Snow Cover.
For 21fix run, sea ice edge of CLIMAP data is used for Feb. and Aug.,
the shape of sea ice edge for other months is determined by considering
summer SST over both hemisphere. Grid with lower summer SST is freeze earlier
and melt later.
Precipitation falling on a surface with skin temperature < 273.15 K
accumulates as snow, and a snowpack melts if the skin temperature exceeds
this value. Fractional coverage of a grid box is determined by the ratio
of the local snow mass to a critical threshold of 200 kg/(m2).
Sublimation of snow contributes to the surface evaporative flux (see Surface
Fluxes), and snowmelt augments soil moisture and runoff (see Land Surface
Processes). Snow cover alters the evaporation efficiency and permeability
of moisture, as well as the albedo, roughness, and thermal properties of
the surface (see Surface Characteristics).
The surface is classified according to the 32 vegetation types of Matthews
(1983), but with only the locally dominant type specified for each grid
box. The stomatal resistance of the vegetation is a prescribed spatially
uniform value, but is set to zero in desert areas. Over ice surfaces, the
roughness length is a constant 1 x 10-3 m. Over ocean, the roughness
length is a function of the surface momentum flux, following the formulation
of Miller et al. (1992). Over land, roughness lengths are assigned according
to vegetation type following Takeuchi and Kondo (1981). IN areas with snow
cover, the roughness length is decreased proportional to the square root
of the fractional snow cover. The roughness length for calculation of surface
momentum fluxes is 10 times the corresponding value for heat and moisture
Over ice surfaces, the albedo is a constant 0.7 (unaffected by snow
accumulation). Over ocean, the albedo depends on sun angle and the optical
thickness of the atmosphere. The albedos of the land surface are specified
according to vegetation type from the data of Matthews (1983). For snow-covered
land, the albedo increases over that of the background proportional to
the square root of the fractional snow cover. Longwave emissivity is everywhere
specified to be 1.0 (i.e., blackbody emission). See also Surface Fluxes
and Land Surface Processes.
Solar absorption at the surface is determined from the albedo, and longwave
emission from the Planck equation with prescribed emissivities (see Surface
Characteristics). The representation of turbulent surface fluxes of momentum,
heat, and moisture follows Monin-Obukhov similarity theory as expressed
by the bulk formulae of Louis (1979). The requisite wind, temperature,
and humidity values are taken to be those at the lowest atmospheric level
(see Vertical Domain).
The associated drag/transfer coefficients are functions of the surface
roughness (see Surface Characteristics) and vertical stability expressed
as a function of a modified Richardson number (see Planetary Boundary Layer).
The effect of free convective motion is incorporated into the surface wind
speed following Miller et al. (1992), and the surface wind speed also is
not allowed to fall below 4 m/s.
For calculation of the moisture flux over ocean, ice, and snow-covered
surfaces, the evaporation efficiency beta is unity. Over partially snow-covered
grid boxes, beta increases as the square root of the snow fraction (see
Snow Cover). The evaporation efficiency over vegetation is limited by the
specified stomatal resistance. Cf. Numaguti et al. (1995) for further details.
See also Land Surface Processes).
Land Surface Processes
The skin temperature of soil and land ice is predicted by a heat diffusion
equation that is discretized in 3 layers with a zero-flux lower boundary
condition; heat capacity and conductivity are spatially uniform values.
Surface snow is treated as part of the uppermost soil layer, and thus modifies
its heat content, as well as the heat conduction to lower layers.
Soil liquid moisture is predicted in a single layer according to the
"bucket" formulation of Manabe et al. (1965)) The moisture field capacity
is a spatially uniform 0.15 m, with surface runoff occurring if the predicted
soil moisture exceeds this value. Snowmelt contributes to soil moisture,
but if snow covers a grid box completely, the permeability of the soil
to falling liquid precipitation becomes zero. For partial snow cover, the
permeability decreases proportional to increasing snow fraction (see Snow
Soil moisture is depleted by surface evaporation; the evaporation efficiency
beta (see Surface Fluxes) is not determined solely by the ratio of soil
moisture to its saturation value, but is limited by the specified stomatal
resistance of the vegetation. Other effects of vegetation, such as the
interception of precipitation by the canopy and its subsequent reevaporation,
are not included. See also Surface Characteristics and Surface Fluxes.
Last update November 9, 1998. For further information, contact: Céline