NASA Ames Mars Climate Modeling Group

Brief Description

  • Models
  • »
  • 1-D Radiation Code
  • »
  • Brief Description

We describe here the radiation code used in our general circulation model (GCM).

This code is a stand-alone version that calculates instantaneous profiles of the solar and infrared fluxes given a set of basic input parameters.

The code can be easily interfaced to any GCM to provide radiative heating/cooling rates or, with some modification, it could also be adopted for 1-D time marching calculations to simulate a diurnal cycle.

The code treats both visible and infrared radiation and accounts for gaseous absorption by CO2 and water vapor, Rayleigh scattering, and scattering and absorption due to water ice and dust particles.

The code is based on a two-stream solution to the radiative transfer equation for plane parallel atmospheres. It can be configured for an arbitrary number of layers (N).

The user-supplied inputs are the:

  • » number of layers
  • » surface pressure
  • » tropopause pressure
  • » surface temperature
  • » atmospheric temperature and water vapor profiles
  • » solar zenith angle
  • » and layer opacities for dust and water ice clouds

The atmosphere is assumed to be 95% CO2 with layer masses distributed according to a sigma coordinate system. The model outputs the upward, downward, and net solar and infrared fluxes at layer boundaries, and heating rates at layer mid-points.

Two-Stream Solutions

Two-stream solutions yield the upward and downward flux at a given level in the atmosphere. They are obtained by solving a simplified set of coupled linear differential equations that come from integrating the radiative transfer equation separately over the upward- and downward-facing hemispheres.

In 2-Stream_Solutions.pdf we provide a detailed derivation of the 2-stream equations and how we solve them for both homogeneous and layered atmospheres. We derive their general form starting from the fundamental equation of radiative transfer and illustrate what they look like for two specific examples: the Eddington and Hemispheric mean approximations. Here we briefly summarize the key concepts. The reader is referred to 2-Stream_Solutions.pdf for further details.

Two-stream solutions are distinguished by the assumptions made regarding the angular dependence of the intensity field and phase function. Meador and Weaver (1980) have shown that all two-stream equations have the same general form and can be written as

Two-stream Equations

where F and F are the upward and downward fluxes, the gamma coefficients depend on the type of two-stream approximation employed, and the source terms are

Source Terms Equations

for solar radiation, and

Infrared Radiation Equation

for infrared radiation.

πFs is the solar flux at the top of the atmosphere, τ the extinction optical depth, μ0 the solar zenith angle, ω0 the single scattering albedo, and B(τ) is the Planck function.

For a single homogeneous layer, the solution to these equations is,

Single Homogeneous Layer Equation

where k1 and k2 are constants determined by the boundary conditions, Γ is a parameter that depends on the specific two stream solution employed (i.e., the coefficients), and the C terms are source functions which take specific forms depending on which part of the spectrum the solutions are applied (see 2-Stream_Solutions.pdf).

These equations can be extended to multiple vertically inhomogeneous layers by matching fluxes at layer boundaries. The solution then leads to an NxN matrix equation which we tridiagonalize for more efficient inversion.

Table 1 lists the distinguishing features of the three two-stream approximations available in our code along with their corresponding coefficients.

Table 1. 2-Stream approximations and their coefficients. Note that γ4 = 1-γ3

Table: Stream approximations

As shown in 2-Stream_Solutions.pdf, the accuracy of the Eddington 2-stream solution can be increased for cases where aerosols have strongly forward scattering peaks in the visible by scaling the asymmetry factor, single scattering albedo, and optical depth.

These solutions are referred to as Delta-Eddington solutions. The scalings are given by

Delra-Eddington Solutions

In our code, the default 2-stream approximation is the Delta-Eddington for the visible bands, and the hemispheric mean for the infrared bands. While the user can select between the Delta Eddington and Quadrature approximations in the visible, we only provide one choice in the infrared. There, the hemispheric mean approximation is the only approximation that conserves energy (see Toon et al., 1989). Hence, we use this approximation exclusively in the infrared.

Spectral Bands

The code treats visible and infrared wavelengths separately. The visible is divided into 7 spectral intervals from 0.40 - 4.5 μm, while the infrared is divided into 5 intervals from 4.5 - 1000 μm. Thus, there are 12 bands in total. Table 2 lists these intervals.

Table 2. GCM Spectral Intervals

GCM Spectral Intervals

* These fluxes are taken from the 1985 Wehrli Standard Extraterrestrial Solar Irradiance Spectrum. The full spectrum can be downloaded from The BB sun produces an integrated flux over these intervals of 1359 Wm-2 at 1AU compared to 1355 Wm-2 for the observed sun.

More bands are needed in the visible compared to the infrared because the wide spacing in spectral space of the weak and strong near infrared CO2 bands is not well resolved with broad visible bands by our correlated-k approach to calculating opacities (see below).

Gaseous Opacities

Gaseous opacities for binary mixtures of CO2 and water vapor are derived from correlated-k distributions. These distributions are calculated off-line using a line-by-line code that generates spectra with line parameters taken from the HITRAN database and


Figure 1. (Click to enlarge) Cumulative probability distribution of absorption coefficients (cm2 / molecule) for the 15 μm band for a temperature of 200 K and pressure of 6.1 mbars.

theoretical predictions. Line widths are adjusted to represent CO2 broadening. A Voigt profile is used at low pressures, and a Lorentz profile at high pressures. The line widths are extended at high pressures so as to include all significant absorption. The abundance of the deuterated species for H2O was adjusted for Mars conditions.

The absorption coefficients are sorted, reordered in g space (cumulative probability space), and stored in a look-up table on a 7 x 11 x 10 temperature (T), pressure (P), and water vapor volume mixing ratio (QH2O) grid with T = (50, 100, 150, …., 350K), P = (10-6, 10-5, …., 104) hPa, and QH2O = (10-7, 10-6, …., 10-1, 0.2, 0.3, 0.4). Actual values are linearly interpolated from this table (linear in log space for pressure).

An example of k-distributions produced in this manner for pure CO2 is shown in Figure 1.

Integration within each spectral interval is performed using a modified 16-point Gaussian quadrature. The modification places 8 gauss points between g = 0 - 0.95, and 8 between g = 0.95 - 1.0. The boundaries and weights of these intervals are listed in Table 3. This spacing was chosen since only a small fraction of the k-distribution in g space contributes to the sharp increase in k values at the high end of the distribution where absorption is dominated by line centers. However, as was the case for the solar wavelengths, when the lines are spread out and the band is too wide, even this split does not resolve that contribution to absorption.

A certain fraction (f) of each spectral interval can have negligible absorption. The line-by-line spectra-generating code keeps track of this fraction and we use it to modify the fluxes computed from the 2-stream solutions by multiplying them by (1-f). A seventeenth gauss point, our so-called “clear channel”, accounts for that fraction of the interval with negligible absorption. In practice, we add those gauss points with negligible opacity (τ < TLIMITS = 10-3 at visable wavelengths; τ < TLIMITI = 5 x 10-3 at infrared wavelengths) into the clear channel to avoid going through the full radiative transfer code. Rayleigh scattering cross-sections for CO2 are taken from Hansen and Travis (1974).

Aerosol Opacities

Aerosol opacities at a reference wavelength (0.67μm) are user supplied. However, the code assumes they have a specific set of optical properties. For dust, the present version uses single scattering properties (Qext, Qscat, and g) derived from a Mie code assuming a log-normal size distribution (cross-section weighted mean effective radius reff = 1.5 μm and standard deviation σ = 0.5), and the wavelength-dependent refractive indices of Wolff et al. (2009). For water ice particles we assume reff = 4.0 μm, σ = 0.1, and we use refractive indices given in Warren (1984). The wavelength dependence of the refractive indices for dust and ice, and their single scattering parameters are shown in Figures 2 & 3. Digital tables are provided here: Dust Refractive Indicies, Dust Scattering Properties, Water Ice Refractive Indicies, Water Ice Scattering Properties.

Table 3. Gauss weights for each spectral interval.

Interval and Weight

Figure 2. (Click to enlarge) Refractive indicies for dust (left) and water ice (right). Dust indicies are taken from Wolff (2009); water ice indicies from Warren (1984).

Figure 3. (Click to enlarge) Single scattering properties of dust (left) and water ice (right).

Calculating the Effective Scattering Parameters for the Two-Stream Code

The two-stream code requires the effective scattering parameters to solve for the fluxes. These parameters represent the combined effect of all radiatively active species. The effective opacity is the sum of all extinction (scattering + absorbing) opacities.

Effective opacity is the sum of all extinction opacities

while the effective single scattering albedo is

Effective single scattering albedo

while the effective asymmetry factor albedo is

Effective asymmetry factor albedo

*Note that the effective asymmetry parameter is weighted by the scattering opacities (not the extinction opacities), and that g = 0 for Rayleigh scattering.

Sample Results

Solar and infrared fluxes and heating rates for several representative cases are shown in Figure 4. Digital output from these figures can be found here (link). In these runs, the atmosphere is isothermal at 200K and the ground temperature is 205K.

Figure 4. (Click to enlarge) Sample output from the radiation code for different levels of dust in the atmosphere.