the PMIP models (Bonfils et al. 1998)
PMIP documentation for UKMO
United Kingdom Meteorological Office:
Model UKMO HADAM2 (2.5x3.75 L19) 1997
Chris Hewitt, Hadley Centre for Climate Prediction and Research, United
Kingdom Meteorological Office,London Road,London Road, Bracknell, Berkshire
RG12 2SY, United Kingdom; Phone: +44 1344 85 4520; Fax:+44 1344 85 4898;
for fixed experiments:
UKMO HADAM2 (2.5x3.75 L19) 1997 or UKMO - Unified Model vn3.2 AGCM -
for computed experiments:
UKMO HADAM2b (2.5x3.75 L19) 1997 or UKMO - Unified Model vn3.4 Slab
GCM - 2x4L19
Model Identification for PMIP
0fix, 6fix, 0cal, 21cal
Number of days in each month: 30 30 30 30 30 30 30 30 30 30 30 30
The UKMO HADAM2 model part of the line of Unified Models (UM) intended
to provide a common framework for forecasting and climate applications
(cf. Cullen 1993 ). The dynamical formulations are those described by Bell
and Dickinson (1987) ; the physical parameterizations are substantially
modified from those of an earlier UKMO model documented by Slingo (1985)
The UKMO HADAM2b model is a coupled atmosphere-mixed layer ocean model.
The atmospheric model is described in :
Hewitt, C.D., and Mitchell, J.F.B., 1996: GCM simulations of the climate
of 6 kyr BP : Mean changes and interdecadal variability. J. Climate, 9,
Johns,T.C., R.E. Carnell, J.F. Crossley, J.M. Gregory, J.F.B. Mitchell,
C.A. Senior, S.F.B. Tett, and R.A. Wood, 1997. The second Hadley Centre
coupled ocean-atmosphere GCM : Model Description, Spinup and Validation.
Climate Dyn., 13, 103-134.
Computed experiments are described briefly in:
Hewitt, C.D., and Mitchell, J.F.B., 1997 : Radiative forcing and response
of a GCM to ice age boundary conditions : cloud feedback and climate sensitivity.
Climate Dynamics, 13, 821-834.
There is a Climate Research Technical Note by Tim Johns (1996), number
71 which describes the UKMO HADAM2 model in detail.
The differences between HADAM2b and HADAM2 are briefly described in
the Hewitt and Mitchell (1997) paper. The PMIP documentation (below) refers
Fourth-order finite differences on a B-grid (cf. Arakawa and Lamb 1977
, Bell and Dickinson 1987 ) in spherical polar coordinates. Mass-weighted
linear quantities are conserved, and second moments of advected quantities
are conserved under nondivergent flow.
2.5 x 3.75-degree latitude-longitude grid.
Surface to about 5 hPa; for a surface pressure of 1000 hPa, the lowest
atmospheric level is at about 997 hPa.
Finite differences in hybrid sigma-pressure coordinates after Simmons and
Strüfing (1981) . Mass and mass-weighted potential temperature and
moisture are conserved. See also Horizontal Representation.
The model uses a hydrid cor-ordinate in the vertical (the top 3 levels
are pressure layers, and the bottom 4 are terrain-following sigma layers,
with a smooth transition in between) with the following levels: 0.005,
0.015, 0.030, 0.057, 0.100, 0.149, 0.200, 0.250, 0.300, 0.355, 0.422, 0.505,
0.599, 0.699, 0.793, 0.870, 0.930, 0.975, 0.997.
There are 19 unevenly spaced hybrid levels. For a surface pressure of
1000 hPa, 4 levels are below 800 hPa and 7 levels are above 200 hPa.
The PMIP simulation was run on a Cray C90 computer in a UNICOS environment.
For the PMIP experiment, about 2 minutes Cray C90 computation time per
simulated day (about half this time being associated with output postprocessing).
For the PMIP experiment, the model atmosphere, soil moisture, and snow
cover/depth were initialized for 1 June from a previous model simulation.
Snow mass for areas of permanent land ice was initially set to 5 x 104
kg/(m2). The model was then integrated forward for months before
starting the fix experiments.
Time Integration Scheme(s)
Time integration proceeds mainly by a split-explicit scheme, where the
solution procedure is split into "adjustment" and "advection" phases. In
the adjustment phase, a forward-backward scheme that is second-order accurate
in space and time is applied. The pressure, temperature, and wind fields
are updated using the pressure gradient, the main part of the Coriolis
terms, and the vertical advection of potential temperature. In the advective
phase, a two-step Heun scheme is applied. A time step of 30 minutes (including
a 10-minute adjustment step) is used for integration of dynamics and physics,
except for full calculation of shortwave/longwave radiation once every
3 hours. In addition, an implicit scheme is used to compute turbulent vertical
fluxes of momentum, heat, and moisture in the planetary boundary layer
(PBL). Cf. Cullen et al. (1991) for further details. See also Diffusion,
Planetary Boundary Layer, and Surface Fluxes.
To prevent numerical instability, the orography is smoothed in high latitudes
(see Orography), and Fourier filtering is applied to mass-weighted velocity
and to increments of potential temperature and total moisture. Negative
values of atmospheric moisture are removed by summing the mass-weighted
positive values in each horizontal layer, and rescaling them to ensure
global moisture conservation after the negative values are reset to zero.
For the PMIP simulation, the model history is written once every 6 hours.
(All average quantities in the PMIP monthly-mean standard output data are
computed from samples taken at every 30-minute time step.)
Primitive-equation dynamics, formulated to ensure approximate energy conservation,
are expressed in terms of u and v winds, liquid/ice water potential temperature,
total water, and surface pressure (cf. White and Bromley 1988 ).
Linear conservative horizontal diffusion is applied at fourth-order (Ñ4)
to moisture, and at sixth-order (Ñ6)
to winds and to liquid water potential temperature (cf. Cullen et al. 1991).
Stability-dependent, second-order vertical diffusion of momentum and
of conserved cloud thermodynamic and water content variables (to include
the effects of cloud-water phase changes on turbulent mixing), operates
only in the PBL. The diffusion coefficients are functions of the vertical
wind shear (following mixing-length theory), as well as surface roughness
length and a bulk Richardson number that includes buoyancy parameters for
the cloud-conserved quantities (cf. Smith 1990b ). See also Cloud Formation,
Planetary Boundary Layer, and Surface Fluxes.
The parameterization of orographic gravity-wave drag follows a modified
Palmer et al. (1986) scheme, as described by Wilson and Swinbank (1989)
. The drag is given by the vertical divergence of the wave stress. Near
the surface, the stress is equal to the product of a representative mountain
wave number, the square of the wave amplitude (taken to be the subgrid-scale
orographic variance--see Orography), and the density, wind, and Brunt-Vaisalla
frequency evaluated in near-surface layers. At higher levels, the stress
is given by this surface value weighted by the projection of the local
wind on the surface wind. If this projection goes to zero, the stress is
also zero; otherwise, if the minimum Richardson number falls below 0.25,
the gravity wave is assumed to break. Above this critical level, the wave
is maintained at marginal stability, and a corresponding saturation amplitude
is used to compute the stress.
For 0fix, 6fix, 0cal and 21cal, 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 are simulated.The seasonal cycle of solar
insolation is based on a 360-day year (each month 30 days in length), with
the date of perihelion adjusted to minimize discrepancies (cf. Ingram 1993).
The carbon dioxide concentration are prescribed value of 323, 262 and 229
ppm for control, 6ka and 21ka, respectively. Zonally averaged monthly ozone
profiles are specified from the climatology of Keating et al. (1987) above
hybrid level 0.0225 (see Vertical Representation), and from satellite data
of McPeters et al. (1984) below this level. Radiative effects of water
vapor and clouds, but not those of aerosol, are also included (see Radiation).
The radiation schemes are as described by Ingram (1993)). Incoming insolation
is based on a 360-day year (see Solar Constant/Cycles). Shortwave computations
follow Slingo (1985)) , extended to include 4 spectral intervals (with
boundaries at 0.25, 0.69, 1.19, 2.38, and 4.0 microns) and the interactive
cloud optical properties of Slingo (1989) . Rayleigh scattering is represented
by reflection of 3 percent of the incoming insolation before any interaction
with the atmosphere. Shortwave absorption (by ozone, carbon dioxide, and
water vapor) is treated by use of look-up tables. The effects of pressure
broadening on carbon dioxide and water vapor absorption bands are included;
the overlap of these bands is assumed to be random within each spectral
Calculation of longwave fluxes follows the method of Slingo and Wilderspin
(1986) . Absorption by water vapor, carbon dioxide, and ozone is treated
in 6 spectral bands (with boundaries at 0.0, 4.0x104, 5.6x104,
8.0x104, 9.0x104, 1.1x105, and 1.2x105
m-1), using exponentials for the water vapor continuum and look-up
tables otherwise. The temperature dependence of e-type absorption in the
water vapor continuum is after Roberts et al. (1976) . Pathlengths of gaseous
absorbers are scaled to account for pressure-broadening effects(and for
carbon dioxide, also for temperature-broadening effects).
The convective cloud in each grid box is confined to a single tower
(with full vertical overlap). All absorption of shortwave radiation by
the convective cloud occurs in its top layer, with the remainder of the
beam passing unimpeded through lower layers. In each vertical column, the
shortwave radiation interacts only with 3 stratiform clouds defined by
vertical domain: low (levels 1 to 5), middle (levels 6 to 9), and high
(levels 10 to 18), all randomly overlapped. The cloud with the greatest
area in each domain interacts with the shortwave radiation, and it is assigned
the cloud water content of all the layer clouds in the domain. (Neither
the shortwave nor longwave scheme allows cloud in the top layer, however.)
In the longwave, similar to the formulation of Geleyn and Hollingsworth
(1979) , clouds in different layers are treated as fully overlapped if
there is cloud in all intervening layers, while clouds separated by cloud-free
layers are treated as randomly overlapped.
Sunlight reflected from a cloud is assumed to pass directly to space
without further cloud interactions, but with full gaseous absorption; the
surface albedo (see Surface Characteristics) is adjusted to account for
the increased absorption due to multiple reflections between clouds and
the surface. Cloud shortwave optical properties (optical depth, single-scattering
albedo, and asymmetry factor) are calculated from the cloud water path
(CWP), the effective radius of the drop-size distribution (7 microns for
water and 30 microns for ice), and the sun angle, following the Practical
Improved Flux Method of Zdunkowski et al. (1980) . Cloud longwave emissivity
is a negative exponential function of CWP, with absorption coefficient
65 m2/kg for ice clouds and 130 m2/kg for water clouds.
See also Cloud Formation and Precipitation.
Moist and dry convection are both simulated by the mass-flux scheme of
Gregory and Rowntree (1990)) that is based on the bulk cloud model of Yanai
et al. (1973) . Convection is initiated if a parcel in vertical layer k
has a minimum excess buoyancy beta that is retained in the next higher
level k + 1 when entrainment effects and latent heating are included. The
convective mass flux at cloud base is taken as proportional to the excess
buoyancy; the mass flux increases in the vertical for a buoyant parcel,
which entrains environmental air and detrains cloud air as it rises. Both
updrafts and downdrafts are represented, the latter by an inverted entraining
plume with initial mass flux related to that of the updraft, and with detrainment
occurring over the lowest 100 hPa of the model atmosphere.
When the parcel is no longer buoyant after being lifted from layer m
to layer m + 1, it is assumed that a portion of the convective plumes has
detrained in layer m so that the parcel in layer m + 1 has minimum buoyancy
b. Ascent continues until a layer n is reached at which an undiluted (without
entrainment) parcel originating from the lowest convectively active layer
k would have zero buoyancy, or until the convective mass flux falls below
a minimum value. Cf. Gregory (1990) and Gregory and Rowntree (1990) for
The convection scheme (see Convection) determines the vertical extent of
subgrid-scale convective cloud, which is treated as a single tower in each
grid box. The convective cloud base is taken as the lower boundary of the
first model layer at which saturation occurs, and the cloud top as the
upper boundary of the last buoyant layer. The fractional coverage of each
vertical column by convective cloud is a logarithmic function of the mass
of liquid water condensed per unit area between cloud bottom and top (cf.
Large-scale (stratiform) cloud is prognostically determined in a similar
fashion to that of Smith (1990a) . Cloud amount and water content are calculated
from the total moisture (vapor plus cloud water/ice) and the liquid/frozen
water temperature, which are conserved during changes of state of cloud
water (i.e., cloud condensation is reversible). In each grid box, these
cloud-conserved quantities are assumed to vary (because of unresolved atmospheric
fluctuations) according to a top-hat statistical distribution, with specified
standard deviation. The mean local cloud fraction is given by the part
of the grid box where the total moisture exceeds the saturation specific
humidity (defined over ice if the local temperature is < 273.15 K, and
over liquid water otherwise). See also Radiation for treatment of cloud-radiative
Large-scale precipitation forms in association with stratiform cloud (see
Cloud Formation). For purposes of precipitation formation and the radiation
calculations (see Radiation), the condensate is assumed to be liquid above
0 degrees C, and to be ice below -15 degrees C, with a liquid/ice fraction
obtained by quadratic-spline interpolation for intermediate temperatures.
The rate of conversion of cloud water into liquid precipitation is a nonlinear
function of the large-scale cloud fraction and the cloud-mean liquid water
content, following Sundqvist (1978 , 1981 ) and Golding (1986) . The precipitation
of ice is a nonlinear function of the cloud-mean ice content, as deduced
by Heymsfield (1977) . Liquid and frozen precipitation also form in subgrid-scale
convection (see Convection).
Evaporation/sublimation of falling liquid/frozen precipitation are modeled
after Kessler (1969) and Lin et al. (1983) . Frozen precipitation that
falls to the surface defines the snowfall rate (see Snow Cover). For purposes
of land hydrology, surface precipitation is assumed to be exponentially
distributed over each land grid box, with fractional coverage of 0.5 for
large-scale precipitation and 0.1 for convective precipitation. Cf. Smith
and Gregory (1990) and Dolman and Gregory (1992) for further details.
Planetary Boundary Layer
Conditions within the PBL are typically represented by the first 5 levels
above the surface (centered at about 997, 975, 930, 869, and 787 hPa for
a surface pressure of 1000 hPa), where turbulent diffusion of momentum
and cloud-conserved thermodynamic and moisture variables may occur (see
Cloud Formation and Diffusion). The PBL top is defined either by the highest
of these layers, or by the layer in which a modified bulk Richardson number
(that incorporates buoyancy parameters for the cloud-conserved variables)
exceeds a critical value of unity. Nonlocal mixing terms are included for
heat and moisture. See also Surface Characteristics, Surface Fluxes, and
Land Surface Processes.
For 0fix, 0cal and 6fix : Orography obtained from the U.S. Navy 10-minute
resolution dataset (cf. Joseph 1980 ) is grid-box averaged, and is further
smoothed with a 1-2-1 filter at latitudes poleward of 60 degrees. The orographic
variances required by the gravity-wave drag parameterization are obtained
from the same dataset (see Gravity-wave Drag).
For 21cal we use Peltier (1994).
For fixed SSTs experiments, the model used is HADAM2.
For 0fix: The SST and sea ice datasets come from the Met. Office observational
climatology MOHSST3 for 1951-1980.
For 6 fix: SSTs and sea-ice prescribed at their present day value, as
in the control run.
For computed SSTs experiments, the model was coupled to an ocean model
based on the Bryan-Cox primitive equation model (Bryan 1969a, Cox 1984),
with 20 depht levels in the vertical. (HADAM2b)
For 0cal: SSTs and sea ice were calculated using the ocean model.
For 21 cal: The 21 kyr BP experiments with computed SSTs and sea ice
should be performed in the same way as the CO2 experiments. Typically this
means that a coupled atmosphere-ocean model will be used with prescribed
present day ocean heat flux. This PMIP run uses Peltier ice sheets.
For 0fix and 6fix, sea ice extents are prescribed monthly. The ice may
occupy only a fraction of a grid box, and the effects of the remaining
ice leads are accounted for in the surface roughness length, shortwave
albedo and longwave emission, and turbulent eddy fluxes (see Surface Characteristics
and Surface Fluxes). The spatially variable sea ice thickness is prescribed
from climatological data. Snow falling on sea ice affects the surface albedo
(see Surface Characteristics), but not the ice thickness or thermodynamic
properties. Ice temperature is prognostically determined from a surface
energy balance (see Surface Fluxes) that includes a conduction heat flux
from the ocean below. Following Semtner (1976) , the conduction flux is
proportional to the difference between the surface temperature of the ice
and the subsurface ocean temperature (assumed to be fixed at the melting
temperature of sea ice, or -1.8 degrees C), and the conduction flux is
inversely proportional to the prescribed ice thickness. For computed SSTs
experiments, the model was coupled with a sea ice model. The thermodynamics
of the model is based on the zero-layer model of Stemtner (1976), and a
simple parameterization of sea ice dynamics based on Bryan (1969) is included
using the top layer ocean currents prescribed from values diagnosed in
a full coupled ocean-atmosphere model.
Surface snowfall is determined from the rate of frozen large-scale and
convective precipitation in the lowest vertical layer. (Snowfall, like
surface rainfall, is assumed to be distributed exponentially over each
land grid box--see Precipitation.) On land only, prognostic snow mass is
determined from a budget equation that accounts for accumulation, melting,
and sublimation. Snow cover affects the roughness and heat conduction of
the land surface, and it also alters the albedo of both land and sea ice
(see Surface Characteristics). Snow melts when the temperature of the top
soil/snow layer is > 0 degrees C, the snowmelt being limited by the total
heat content of this layer. Snowmelt augments soil moisture, and sublimation
of snow contributes to the surface evaporative flux over land. See also
Surface Fluxes and Land Surface Processes.
Surface types include land, ocean, sea ice, and permanent land ice. Sea
ice may occupy a fraction of a grid box (see Sea Ice). On land, 24 different
soil/vegetation types are specified from the 1 x 1-degree data of Wilson
and Henderson-Sellers (1985) . The effects of these surface types on surface
albedo (see below) and on surface thermodynamics and moisture (see Land
Surface Processes) are treated via parameters derived by Buckley and Warrilow
On each surface, roughness lengths are specified for momentum, for heat
and moisture, and for free convective turbulence (which applies in cases
of very light surface winds under unstable conditions). Over oceans, the
roughness length for momentum is a function of surface wind stress (cf.
Charnock 1955 ), but is constrained to be at least 10-4 m; the
roughness length for heat and moisture is a constant 10-4 m,
and it is 1.3 x 10-3 m for free convective turbulence. Over
sea ice, the roughness length is a constant 3 x 10-3 m, but
it is 0.10 m for that fraction of the grid box with ice leads (see Sea
Ice). Over land, the roughness length is a function of vegetation and small
surface irregularities; it is decreased as a linear function of snow cover,
but is at least 5 x 10-4 m. Cf. Smith (1990b) for further details.
The surface albedo of open ocean is a function of solar zenith angle.
The albedo of sea ice varies between 0.60 and 0.85 as a linear function
of the ice temperature above -5 degrees C, and it is also modified by snow
cover. Where there is partial coverage of a grid box by sea ice, the surface
albedo is given by the fractionally weighted albedos of sea ice and open
ocean. Surface albedos of snow-free land are specified according to climatological
vegetation and land use. Snow cover modifies the albedo of the land surface
depending (exponentially) on the depth (following Hansen et al. 1983 )
and (linearly) on the snow temperature above -2 degrees C (following Warrilow
et al. 1990 ), as well as on the background vegetation (following Buckley
and Warrilow 1988 ). Cf. Ingram (1993) for further details.
Longwave emissivity is unity (blackbody emission) for all surfaces.
Thermal emission from grid boxes with partial coverage by sea ice (see
Sea Ice) is calculated from the different surface temperatures of ice and
the open-ocean leads, weighted by the fractional coverage of each.
Surface solar absorption is determined from the surface albedos, and longwave
emission from the Planck equation with constant emissivity of 1.0 (see
Turbulent eddy fluxes are formulated as bulk formulae in a constant-flux
surface layer, following Monin-Obukhov similarity theory. The momentum
flux is expressed in terms of a surface stress, and the heat and moisture
fluxes in terms of cloud-conserved quantities (liquid/frozen water temperature
and total water content) to account for effects of phase changes on turbulent
exchanges (see Cloud Formation). These surface fluxes are solved by an
implicit numerical method.
The surface atmospheric variables required for the bulk formulae are
taken to be at the first level above the surface (at 997 hPa for a 1000
hPa surface pressure). (For diagnostic purposes, temperature and humidity
at 1.5 m and the wind at 10 m are also estimated from the constant-flux
assumption.) Following Louis (1979) , the drag/transfer coefficients in
the bulk formulae are functions of stability (expressed as a bulk Richardson
number) and roughness length, and the same transfer coefficient is used
for heat and moisture. In grid boxes with fractional sea ice, surface fluxes
are computed separately for the ice and lead fractions, but using mean
drag and transfer coefficients obtained from linearly weighting the coefficients
for ice and lead fractions (see Sea Ice).
The surface moisture flux also is a fraction beta of the local potential
evaporation for a saturated surface. Over oceans, snow and ice, and where
there is dew formation over land, beta is set to unity; otherwise, beta
is a function of soil moisture and vegetation (see Land Surface Processes).
Above the surface layer, momentum and cloud-conserved variables are
mixed vertically within the PBL by stability-dependent diffusion. For unstable
conditions, cloud-conserved temperature and water variables are transported
by a combination of nonlocal and local mixing. Cf. Smith (1990b, 1993 )
for further details. See also Diffusion and Planetary Boundary Layer.
Land Surface Processes
A vegetation canopy model (cf. Warrilow et al. 1986 and Shuttleworth 1988
) includes effects of moisture condensation, precipitation interception,
direct evaporation from wet leaves and from surface ponding, and evapotranspiration
via root uptake of soil moisture. The fractional coverage and water-storage
capacity of the canopy vary spatially by vegetation type (see Surface Characteristics),
with a small storage added to represent surface ponding. The canopy intercepts
a portion of the precipitation, which is exponentially distributed over
each grid box (see Precipitation). Throughfall of canopy condensate and
intercepted precipitation occurs in proportion to the degree of fullness
of the canopy; evaporation from the wet canopy occurs in the same proportion,
as a fraction of the local potential evaporation. Cf. Gregory and Smith
(1990) for further details.
Soil moisture is predicted from a single-layer model with spatially
nonuniform water-holding capacity; it is augmented by snowmelt, precipitation,
and the throughfall of canopy condensate. This moisture infiltrates the
soil at a rate depending on saturated soil hydraulic conductivity enhanced
by effects of root systems that vary spatially by vegetation type (see
Surface Characteristics). The noninfiltrated moisture is treated as surface
runoff. Subsurface runoff from gravitational drainage is also parameterized
as a function of spatially varying saturated soil hydraulic conductivity
and of the ratio of soil moisture to its saturated value (cf. Eagleson
1978 ). Soil moisture is depleted by evaporation at a fraction beta of
the local potential rate (see Surface Fluxes). The value of beta is obtained
by combining the canopy evaporation efficiency (see above) with a moisture
availability parameter that depends on the ratio of soil moisture to a
spatially varying critical value, and of the ratio of the stomatal resistance
to aerodynamic resistance (cf. Monteith 1965 ). The critical soil moisture
and stomatal resistance vary spatially by soil/vegetation type (see Surface
Characteristics). Cf. Gregory and Smith (1990) and Smith (1990b) for further
Soil temperature is predicted after Warrilow et al. (1986) from heat
conduction in four layers. The depth of the topmost soil layer is given
by the penetration of the diurnal wave, which depends on the spatially
varying soil heat conductivity/capacity (see Surface Characteristics).
The lower soil layers are, respectively, about 3.9, 14.1, and 44.7 times
the depth of this top layer. The top boundary condition for heat conduction
is the net downward surface energy balance (see Surface Fluxes), including
the latent heat of fusion for snowmelt (see Snow Cover); the bottom boundary
condition is zero heat flux. Heat insulation by snow (see Snow Cover) is
modeled by reducing the thermal conductivity between the top two soil layers;
however, subsurface moisture does not affect the thermodynamic properties
of the soil. Cf. Smith (1990b) for further details.
Last update November 9, 1998. For further information, contact: Céline