Dr. David Straus, Center for Ocean-Land-Atmosphere Studies, 4041 Powder Mill Road, Suite 302, Calverton, Maryland 20705-3106; Phone: +1-301-595-7000 or +1-301-902-1255; Fax: +1-301-595-9793; e-mail: straus@cola.iges.org; World Wide Web URL: http://grads.iges.org/home.html.
COLA COLA1.1 (R40 L18 ) 1993
The COLA model is derived from the NMC Medium-Range Forecast (MRF) model (cf. NMC Development Division 1988
[1]), but with substantial modifications in the treatment of vertical diffusion, radiation and cloud-radiative interactions, surface characteristics and fluxes, and land-surface processes.
Key documentation of the basic model framework is provided by Kinter et al. (1988)
[2], with subsequent modifications described by Sato et al. (1989a
[3],b
[4]), Xue et al. (1991)
[5], and Hou (1991)
[6].
Spectral (spherical harmonic basis functions) with transformation to a Gaussian grid for calculation of nonlinear quantities and physics.
Spectral rhomboidal 40 (R40), roughly equivalent to 1.8 x 2.8 degrees latitude-longitude.
Surface to 10 hPa. For a surface pressure of 1000 hPa, the lowest atmospheric level is at a pressure of about 995 hPa.
Finite-difference sigma coordinates.
There are 18 unevenly spaced sigma levels. For a surface pressure of 1000 hPa, 5 levels are below 800 hPa and 4 levels are above 200 hPa.
The AMIP simulation was run on a Cray 2 computer using a single processor in CTSS and UNICOS environments.
For the AMIP experiment, about 18 minutes Cray 2 computation time per simulation day.
For the AMIP experiment, the model atmospheric state is initialized from the NMC analysis for 1 January 1979, with nonlinear normal mode initialization also performed (cf. Machenauer 1977)
[33]. January soil moisture and snow cover/depth are obtained from GFDL climatologies.
Time integration is by a leapfrog semi-implicit scheme with an Asselin (1972)
[7] frequency filter. A time step of 12 minutes is used for dynamics and physics, except for full calculation of atmospheric radiation, which is done hourly for the shortwave fluxes and every 3 hours for the longwave fluxes. An implicit scheme with explicit coefficients also is used to eliminate numerical oscillation while integrating the coupled heat and mass exchanges between the surface and the atmospheric boundary layer (cf. Sato et al. 1989a)
[3].
Mean silhouette orography is determined for each Gaussian grid box (see Orography). Negative atmospheric moisture values arising from the model's spectral truncation are filled by resetting these to zero.
For the AMIP simulation, the model history is written every 6 hours.
Primitive-equation dynamics are expressed in terms of vorticity, divergence, surface pressure, specific humidity, and virtual temperature.
- Fourth-order (del^4) horizontal diffusion is applied to the vorticity, divergence, specific humidity, and virtual temperature. (The coefficient of diffusion for the divergence is 0.61 x 10^16 m^4/s, while it is 0.81 x 10^16 m^4/s for the other fields.) For the specific humidity and virtual temperature, the del^2 correction necessary to account for diffusion on constant pressure (rather than constant sigma) surfaces is also applied. This correction includes a priori specification of estimates for global-mean specific humidity and temperature.
- Stability-dependent vertical diffusion with Mellor and Yamada (1982)
[8] level-2 turbulence closure is used in the planetary boundary layer and free atmosphere. To obtain the eddy diffusion coefficients, a prognostic equation is solved for the turbulent kinetic energy, with other second-order moments being calculated diagnostically.
Gravity-wave drag is simulated as described by Kirtman et al. (1993)
[9] and Alpert et al. (1988)
[10]. The parameterization includes determination of the momentum flux due to gravity waves at the surface, as well as at higher levels. The gravity-wave drag (stress) is given by the convergence of the vertical momentum flux. The surface stress is calculated as a nonlinear function of both the surface wind speed and the local Froude number, following Pierrehumbert (1987)
[11]. Vertical variations in the momentum flux occur when the local value of the wave-modified Richardson number becomes less than 0.25 and the stress vanishes (cf. Eliassen and Palm 1961
[12]), or when wave breaking occurs (the local Froude number becomes critical); in the latter case the momentum flux is reduced according to the wave saturation hypothesis of Lindzen (1981)
[13]. See also Orography.
The solar constant is the AMIP-prescribed value of 1365 W/(m^2). Both seasonal and diurnal cycles in solar forcing are simulated.
The carbon dioxide concentration is the AMIP-prescribed value of 345 ppm. Seasonal zonal profiles of ozone are prescribed from GFDL climatologies, with daily values obtained by linear interpolation. Radiative effects of water vapor, but not those of aerosols, also are included (see Radiation).
- The radiation code follows Harshvardhan et al. (1987)
[14]. The shortwave scheme is based on the method of Lacis and Hansen (1974)
[15]. Six absorption bands are considered, one for ozone in the ultraviolet (wavelengths < 0.35 micron) and visible (wavelengths 0.50 to 0.70 micron) spectral ranges, and five near-infrared bands (wavelengths 0.70 to 4.0 microns). At the surface, solar radiative fluxes are separated into four components: the ultraviolet and visible direct and diffuse beams, and the near-infrared direct and diffuse beams (see Surface Characteristics). For clear-sky conditions, a combined surface-atmosphere treatment is employed for Rayleigh scattering. For cloudy skies, multiple scattering effects are treated by a delta-Eddington approach (cf. Joseph et al. 1976
[16]). In this case optical depth is estimated from cloud temperature and pressure thickness after Harshvardhan et al. (1989)
[17], while the asymmetry factor and single-scattering albedo are prescribed.
- The longwave scheme follows the broadband transmission approach of Chou (1984)
[18] for water vapor, that of Chou and Peng (1983)
[19] for carbon dioxide, and that of Rodgers (1968)
[20] for ozone. The treatment of continuum absorption by water vapor follows Roberts et al. (1976)
[21]. Water vapor absorption is calculated in spectral domains corresponding to two band centers (0-3.4 x 10^4 m^-1 and 1.38 x 10^5-1.90 x 10^5 m^-1) and the associated band wings (in the range 3.40 x 10^4-3.00 x 10^5 m^-1). Carbon dioxide absorption is treated similarly (with band-center domain 6.80 x 10^4-7.20 x 10^4 m^-1 and band-wing domains in the range 5.40 x 10^4-8.00 x 10^4 m^-1). Ozone absorption is calculated in the interval 9.80 x 10^4-1.10 x 10^5 m^-1. For cloudy-sky conditions, longwave emissivity is a function of the optical thickness of the cloud layer (see above). Cumuloform clouds are treated as fully overlapped in the vertical, and stratiform clouds as randomly overlapped (see Cloud Formation).
- Penetrative convection is simulated following Kuo (1965)[22] with modifications as described by Sela (1980)[23]. Convection occurs in the presence of large-scale moisture convergence accompanied by a moist unstable lapse rate under moderately high relative humidity conditions. The vertical integral of the moisture convergence determines the total moisture available for moistening vs heating (through precipitation formation) the environment. If the moisture convergence in the first several lowest layers of a vertical column exceeds a critical threshold (2 x 10^-8/sec), a moist adiabat is computed assuming the bottom layer is saturated, and using the preliminary pressure and temperature prediction of the model.
- An unstable subcolumn is then defined which extends from the bottom layer to the first layer for which a moist adiabatically lifted air parcel is not warmer than the environment. Within this subcolumn, the departures of the temperature and specific humidity of a saturated parcel from the respective environmental profiles in each layer determine the fraction of the total available moisture contributed to latent heat release vs moistening of that layer; the temperature and humidity profiles are revised accordingly. In addition, if the revised temperature profile exceeds a dry adiabatic lapse rate, a dry convective adjustment is performed and the moisture in the column is vertically redistributed to reflect the adjusted temperature profile. Cf. Sela (1980)
[23] for further details.
- Following Tiedtke (1983)
[24], simulation of shallow (nonprecipitating) convection is parameterized as an extension of the vertical diffusion scheme (see Diffusion).
- Cloud formation is simulated following the diagnostic method of Slingo (1987)
[25]. Two basic cloud types--cumuloform and stratiform--are represented. The height of cumuloform cloud is determined by the level of nonbuoyancy for moist adiabatic ascent in the model's convective scheme (see Convection). The cumuloform cloud fraction (not exceeding 0.8) is estimated from the scaled 3-hour mean convective precipitation rate (see Precipitation). In the case of convection penetrating above the 400 hPa level, the cumuloform cloud is capped by a cirrus anvil.
- Up to three separated layers of stratiform cloud are allowed in predefined domains (high, middle, and low). Clouds associated with fronts/tropical disturbances have fractional extent determined by a quadratic function of the difference between the local relative humidity and a threshold value of 80 percent. In the low-cloud domain, the fractional cloudiness is reduced in regions of moist subsidence, while stratus cloud forms if a temperature inversion is present under dry layers.
- Precipitation is produced both from large-scale condensation and from the convective scheme (see Convection). The large-scale precipitation algorithm compares the predicted specific humidity with a modified saturation value that is a function of temperature and pressure of a vertical layer (cf. Sela 1980)
[23]. If the predicted humidity exceeds this threshold value, condensation occurs and the predicted temperature field is adjusted to account for the associated latent heat release.
- To prevent convective precipitation when the environment is unstable but relatively dry, the falling condensed water evaporates as it acts to progressively saturate lower layers. Large-scale precipitation may similarly evaporate. In both cases all precipitation that penetrates the bottom atmospheric layer is allowed to fall to the surface. See also Snow Cover.
The PBL is typically represented by the first six vertical levels, which correspond to pressures of 995, 981, 960, 920, 857, and 779 hPa for a surface pressure of 1000 hPa. See also Diffusion, Surface Characteristics, and Surface Fluxes.
Silhouette orography is derived from the U.S. Navy terrain height data at a resolution of 10-minute arc (cf. Joseph 1980
[26]). These data are scanned to obtain the maximum height for each longitude and latitude encompassed by each model Gaussian grid box. The mean of the combined sets of maxima is then assigned as the silhouette height for each grid box.
AMIP monthly sea surface temperature fields are prescribed, with daily values determined by linear interpolation.
AMIP monthly sea ice extents are prescribed. Sea ice covers an entire Gaussian grid box with a fixed thickness of 2 m, and the temperature below the ice is assumed to be 271.16 K. Snow does not accumulate on sea ice. The surface sensible and latent heat fluxes are determined from an energy balance calculation that includes heat conduction through the sea ice (see Surface Fluxes).
If the temperature at the lowest atmospheric level is less than 271.16 K, precipitation falls as snow, but its accumulation is accounted for only on land. Snow cover affects both the surface albedo and the heat transfer/capacity of the soil. For purposes of calculating the surface albedo, the fractional coverage of a grid box by snow is a function of its depth and of the maximum snow cover for a grid box. Snow mass is determined prognostically from a budget equation that accounts for accumulation and melting. Snowmelt contributes to soil moisture, and sublimation of snow to the surface evaporation. See also Surface Characteristics, Surface Fluxes, and Land Surface Processes.
- Roughness lengths over oceans are determined from the surface wind stress after the method of Charnock (1955)
[27]. The roughness length over sea ice is a uniform 1 x 10^-4 m. Over land, the 12 vegetation/surface types of the Simple Biosphere (SiB) model of Sellers et al. (1986)
[28] and associated monthly varying roughness lengths are specified from data of Dorman and Sellers (1989)
[29].
- Over oceans the surface albedo depends on zenith angle, but not spectral interval (cf. Payne 1972
[30]). The albedo of sea ice is a constant 0.50. Surface albedos of vegetated land are prescribed after data of Dorman and Sellers (1989)
[29], and vary monthly according to seasonal changes in vegetation. The land albedo is specified separately for visible (0.0-0.70 micron) and near-infrared (0.70-4.0 microns) spectral intervals, and is also a function of solar zenith angle. The changes in land albedo associated with partial snow cover (including effects of multiple reflections between snow and the vegetation canopy) are parameterized as described by Xue et al. (1991)
[5].
- Surface longwave emissivity is everywhere prescribed to be unity (blackbody emission).
- Surface solar absorption is determined from the surface albedos, and longwave emission from the Planck equation with uniform emissivity of 1.0 (see Surface Characteristics).
- In the lowest atmospheric layer, surface turbulent eddy fluxes of momentum, heat, and moisture follow Monin-Obukhov similarity theory. To avoid an iterative solution for the surface fluxes, the associated drag and transfer coefficients are approximated as analytical functions of the surface characteristics and bulk Richardson number. Over the oceans, the equations formulated by Miyakoda and Sirutis (1986)
[31], expressed as bulk formulae, are used to compute surface fluxes. Over land, stability-dependent drag and transfer coefficients (expressed as aerodynamic and surface resistances) are determined after Xue et al. (1991)
[5].
- Surface evaporation is at the potential rate over oceans, snow, and ice. Over land, the surface moisture flux includes both evapotranspiration via vegetation root uptake (including the effects of bulk stomatal resistance) and direct evaporation from the vegetation canopy and from bare soil (see Land Surface Processes).
- Above the constant-flux surface layer, diffusion of momentum, heat, and moisture are predicted by the Mellor and Yamada (1982)
[8] level-2 turbulence closure scheme (see Diffusion).
- Land-surface processes are simulated following the Xue et al. (1991)
[5] modification of the SiB model of Sellers et al. (1986)
[28]. Within the single-story vegetation canopy, evapotranspiration from dry leaves includes detailed modeling of stomatal and canopy resistances; direct evaporation from the wet canopy and from bare soil is also treated (see Surface Fluxes). Precipitation interception by the canopy is simulated, and its infiltration into the ground is limited to less than the hydraulic conductivity of the soil.
- Soil temperature is determined in two layers by the force-restore method of Deardorff (1978)
[32]. Soil moisture, which is predicted from diffusion equations in three layers, is increased by infiltrated precipitation and snowmelt, and is depleted by evapotranspiration, direct evaporation, and drainage. Both surface runoff and deep runoff from gravitational drainage are simulated. See also Surface Characteristics and Surface Fluxes.
Go to COLA References
Return to COLA Table of Contents
Last update April 19, 1996. For further information, contact: Tom Phillips (
phillips@tworks.llnl.gov)
LLNL Disclaimers
UCRL-ID-116384