The rate at which isoprene is emitted by a forest depends on an array of environmental variables, the forest’s biomass, and its species composition. At present it is unclear whether errors in canopy-scale and process-level isoprene emission models are due to inadequacies in leaf-to-canopy integration theory or the imperfect assessment of the isoprene-emitting biomass in the flux footprint. To address this issue, an isoprene emission model (CANVEG) was tested over a uniform aspen stand and a mixed-species, broad-leaved forest.
The isoprene emission model consists of coupled micrometeorological and physiological modules. The micrometeorological module computes leaf and soil energy exchange, turbulent diffusion, scalar concentration profiles, and radiative transfer through the canopy. Environmental variables that are computed by the micrometeorological module, in turn, drive physiological modules that calculate leaf photosynthesis, stomatal conductance, transpiration and leaf, bole and soil/root respiration, and rates of isoprene emission.
The isoprene emission model accurately predicted the diurnal variation of isoprene emission rates over the boreal aspen stand, as compared with micrometeorological flux measurements. The model’s ability to simulate isoprene emission rates over the mixed temperate forest, on the other hand, depended strongly upon the amount of isoprene-emitting biomass, which, in a mixed-species forest, is a function of the wind direction and the horizontal dimensions of the flux footprint. When information on the spatial distribution of biomass and the flux footprint probability distribution function were included, the CANVEG model produced values of isoprene emission that compared well with micrometeorological measurements. The authors conclude that a mass and energy exchange model, which couples flows of carbon, water, and nutrients, can be a reliable tool for integrating leaf-scale, isoprene emission algorithms to the canopy dimension over dissimilar vegetation types as long as the vegetation is characterized appropriately.
Plants release numerous hydrocarbon compounds into the atmosphere (Fehsenfield et al. 1992; Lerdau et al. 1997). The hydrocarbon compound, isoprene (C5H8), can be very abundant in the air over deciduous forests. Isoprene is of interest to atmospheric chemists because it contributes to the formation of the radiatively active, oxidant gas, ozone, and its reaction with the hydroxyl radical reduces the oxidative capacity of the troposphere (Fehsenfield et al. 1992; Lerdau et al. 1997). The emission of isoprene from forests is also of interest to physiological ecologists studying the carbon cycle. Precursors for isoprene synthesis are derived from the carbon fixation pathway in certain plant genera, such as aspen (Populus) and oak (Quercus) (Sharkey et al. 1991; Monson et al. 1992). Hence, carbon lost as isoprene can be a significant fraction of photosynthetic carbon fixation. Under normal growing conditions plants release 1%–2% of carbon fixed by photosynthesis to the atmosphere as isoprene (Monson and Fall 1989; Sharkey et al. 1991;Harley et al. 1995). Under moderate water deficits this percentage can exceed 10%, as a result of foliage warming and reductions in net photosynthesis (Fang et al. 1996).
The measurement and modeling of isoprene emission rates from forests is still in its infancy. Consequently, tests between field data and model calculations have yielded inconsistent and inconclusive results. For example, recent estimates of isoprene emission rates differed with measurements by 30%–200% (Baldocchi et al. 1995; Lamb et al. 1996; Fuentes et al. 1996). One major source of modeling uncertainty involves the relatively small pool of data that has been used to test isoprene emission models. Before 1992, few direct measurements of canopy-scale isoprene flux studies were conducted or reported in the literature. At present, fewer than a dozen independent field campaigns exist in the literature (Table 1). We also note that the majority of these studies were conducted for only a few weeks, over mixed stands of isoprene emitters, and used new and evolving chemical sensors and measurement methods. Larger datasets from reliable measurement techniques are needed for model testing because individual field measurements are subject to significant sampling errors due to the stochastic nature of the atmosphere (Moncrieff et al. 1996). For instance, elementary statistics calculations indicate that at least 25 samples are needed to detect the difference between two means, within 5% of the population mean, when the coefficient of variation equals 25%. Longer field studies also enable the investigator to study the behavior of a forest across a wider range of environmental conditions, thus providing a more robust test of a biogenic emission model.
Another source of modeling uncertainty stems from the interrelationship between the horizontal flux footprint, detected by a flux measurement system, and the horizontal distribution of isoprene-emitting biomass. Lamb et al. (1996), for example, show that the mean isoprene-emitting biomass of a mixed temperate forest ranges between 110 and 220 g m−2, within a horizontal “footprint” that extended 500 m upwind of the flux measurement tower, a zone wherein 95% of the isoprene flux emanated. This cited span of the emitting biomass depended on whether one computed the arithmetic mean or a distance-weighted average.
At present, it is unclear whether errors in isoprene emission models are due to statistical uncertainty arising from small datasets, inaccurate measurement methods, theoretical limitations of leaf-to-canopy integration models, or from the inadequate specification of the isoprene-emitting biomass. Circumstantial evidence exists to suggest that isoprene fluxes can be integrated from the leaf-to-canopy dimensions with fidelity. Isoprene emission algorithms have been tested and verified at the leaf scale (Guenther et al. 1991, 1993; Harley et al. 1995, 1997). Additionally, it has been demonstrated that leaf-to-canopy integration theory, which will be used in this exercise, can produce realistic and field-validated estimates of water vapor, sensible heat, and CO2 exchange over a broad-leaved forest (Baldocchi and Harley 1995) and an array of crops (Baldocchi and Meyers 1998).
The goal of this paper is to narrow the uncertainties about modeling isoprene emission over forests. To do so, calculations of isoprene emission rates are examined using two new and contrasting sets of test data. One dataset comes from a homogeneous, boreal aspen stand and the other comes from a temperate, mixed-species, oak–maple forest.
The homogeneous forest stand consists of trembling aspen (Populus tremuloides Michx.), an effective isoprene emitter (Fuentes et al. 1996). The forest is situated in the southern range of the boreal forest of Canada (53°63′N, 106°20′W), near Prince Albert, Saskatchewan; these measurements were made under the auspices of the Boreal Ecosystem Atmosphere Study (BOREAS) (Sellers et al. 1995). The aspen forest is 21 m in height and its leaf area index is 2.3. The understory consisted of hazel (Corylus cornuta Marsh), a nonisoprene emitter. The hazel is 2 m tall and its leaf area index is 3.3. The uniform forest extended several kilometers in each direction on flat terrain. These conditions made the site suitable for eddy covariance and flux–gradient measurements. Additional information on the site and canopy characteristics is reported in Black et al. (1996) and Blanken et al. (1997).
Mass and energy fluxes and meteorological variables were measured during the summer of 1994. The eddy-covariance method was used to measure flux densities of CO2, water vapor, and sensible heat. This method involved the implementation of a three-dimensional sonic anemometer and a fast-response infrared gas analyzer (Black et al. 1996; Blanken et al. 1997). These instruments were situated 17 m above the aspen canopy. The flux–gradient method was used to evaluate isoprene flux densities. Eddy exchange coefficients used to compute isoprene flux densities were determined from the ratio between water vapor gradient and flux density measurements. Concentrations of isoprene, water vapor, temperature, and carbon dioxide were measured at 27 and 39 m above the ground. Concentration gradients were computed as finite differences between the concentrations measured at the cited heights.
Isoprene concentrations were determined on continuously sampled air. A gas chromatograph that was coupled to a flame ionization detector was used to make the concentration measurements (Fuentes et al. 1996).
Isoprene emission rates from leaves of trembling aspen were measured to determine leaf-scale model parameters. A glass cuvette system was used for this purpose (Fuentes and Gillespie 1991). The cuvette system controlled temperature within ±2°C of ambient temperature and used natural light and humidity.
The second field experiment was performed near Oak Ridge, Tennessee (35°57.534′N, 84°17.260′W) during the summer of 1996. The forest was a mixed-species stand, vegetated with oak, hickory, maple, and pine. The stand was over 50 years old, its mean height is 24 m, and its leaf area index is 4.9. At this site, oaks (Quercus alba, Q. prinus, Q. velutina) and black gum (Nyssa sylvatica) are strong isoprene emitters and sweet gum (Liquidambar styraciflua) is a weak emitter (Guenther et al. 1996). Maple (Acer rubrum), tulip poplar (Liriodendron tulipifera), hickory (Carya spp.), dogwood (Cornus florida), and beech (Fagus spp.) also inhabit this forest, but they are nonemitters (Guenther et al. 1996).
Flux densities of isoprene between the canopy and the atmosphere were measured using the eddy-covariance and relaxed eddy accumulation (REA) methods (Bowling et al. 1998). The eddy-covariance method employed a three-dimensional sonic anemometer (model SAT-211/3K, Applied Technology, Boulder, CO) and an ozone-induced chemiluminesence instrument (FIS, Hills Scientific, Boulder, CO). The eddy-covariance method was used also to measure flux densities of CO2, water vapor, and heat (see Baldocchi and Harley 1995; Baldocchi 1997a).
The REA system drew air through Teflon tubing to a potassium iodide ozone trap and filter. The air then passed through a manifold. Solenoid valves directed the airstream into separate reservoirs, according to whether the sample was associated with an up- or downdraft. Air samples were collected in 8-L Tedlar bags. The isoprene concentration of air was measured with a gas chromatograph equipped with a flame ionization detector (details are reported in Bowling et al. 1998).
Rates of isoprene emission from leaves were studied in a prior experiment with a cuvette-based gas exchange system. These data are reported in a recent paper by Harley et al. (1997).
Isoprene emission rates from leaves
Isoprene emission rates from leaves were computed with an algorithm, derived by Guenther et al. (1991) and Guenther et al. (1993). The isoprene emission rates, predicted by this algorithm, are a function of several factors, light energy absorbed by a leaf and its temperature:
where SI is a species specific, standardized emission factor (nmol m−2 s−1) at a specified leaf temperature (Tl) and photon flux density (Qp). The functions, f(Qp) and f(Tl), are non dimensional and adjust the standard emission rate according to variations in photon flux density and temperature, respectively.
Specific leaf content, leaf nitrogen content scale, and photosynthetic capacity vary with depth into a canopy (Ellsworth and Reich 1993; Harley and Baldocchi 1995). Since isoprene emission rates scale with photosynthetic capacity scale, they vary with specific leaf weight (Harley et al. 1995, 1997) and leaf nitrogen (Monson et al. 1994; Litvak et al. 1996), too. To accommodate for this behavior, we reduced Sl with depth in the canopy according to an exponential function of leaf area index.
The adjustment factor algorithm for photon flux density is
where α [m2 s−1 μmol−1 (quanta)] and CL are empirical constants. Based on isoprene quantum yield measurements made on aspen leaves (Monson and Fall 1989), α equals 0.0025. The value is close to the standard value of 0.0027 that is used in the model of Guenther et al. (1991) and Guenther et al. (1993), which was derived from measurements on eucalyptus, sweet gum, aspen, and velvet bean. For comparison, Harley et al. (1997) report that α equals 0.0029 and CL equals 1.05 for Quercus alba leaves exposed to shade. For sun leaves, the respective values of α and CL are 0.0013 and 1.25.
Here, R is the universal gas constant (8.314 J K−1 mol−1), CT1 and CT2 are coefficients, and Topt is the optimum temperature and Tk is the leaf temperature in kelvins. The scaling factor equals 1 when temperature equals the reference value of 303 K. For Quercus alba, typical values for CT1, CT2, and Topt are 78 000 J mol−1, 380 000 J mol−1, and 315.8 K, respectively (Harley et al. 1997). Tests indicate that Guenther et al.’s leaf algorithm can account for 90% of observed variability of isoprene flux densities during the day, and it can predict emission rates within 35% (Guenther et al. 1993).
The CANVEG model was used to calculate fluxes of isoprene, carbon dioxide, water vapor, and sensible heat between the canopy and the atmosphere (Baldocchi and Harley 1995; Baldocchi and Meyers 1998). The model consists of coupled micrometeorological and physiological modules. The micrometeorological model computes leaf and soil energy exchange, turbulent diffusion, scalar concentration profiles, and radiative transfer through the canopy. Environmental variables, computed with the micrometeorological module, in turn, drive the physiological models that compute leaf photosynthesis, stomatal conductance, transpiration, leaf, bole and soil/root respiration, and rates of isoprene emission. The model is driven by a minimal number of external variables that were measured above the forest. Environmental inputs include incident photosynthetic photon flux density (Qp), air temperature, wind speed, relative humidity, and CO2 concentration. For the aspen study, these data came from days 207, 215, 216, 219, and 243 during 1994. Data from the oak forest study were acquired between days 227 and 235 during 1996. Plant structural variables include leaf area index, leaf angle orientation, a leaf clumping factor, and canopy height. Model parameters used for the calculations are listed in Table 2. A thorough description of the CANVEG model is presented in appendix A.
If one expects to model isoprene emissions correctly, one must be able to simulate the amount of light absorbed by leaves and their temperature. This exercise involves the modeling of the leaf energy balance through the domain of the canopy. Specifically, it involves simulation of photon transport through the canopy, absorption of photons by leaves, and the partitioning of this absorbed energy into latent and sensible heat fluxes. Details on a test of the CANVEG model for the temperate broad-leaved forest under study have been reported recently (Baldocchi and Harley 1995; Baldocchi 1997a). Hence, we summarize the results briefly, and refer the reader to the cited papers, rather than repeat that material here. Overall, the CANVEG model was able to simulate the diurnal patterns and magnitudes of net radiation (Rn), latent (LE) and sensible (H) heat, and CO2 fluxes densities with acceptable fidelity and accuracy. For instance, measured and modeled values of LE were nearly identical (daily means were 91 ± 22 and 92 ± 19 W m−2, respectively) when the canopy had ample soil moisture (Baldocchi 1997a). Comparisons between measured and modeled mean daily means of Rn agreed within 12% (their mean daily values were 119 ± 37 and 105 ± 34 W m−2, respectively). Finally, mean values of canopy CO2 exchange compared moderately well with modeled values, during the daylight period, when eddy-covariance measurements are less biased. For this comparison, mean measured CO2 flux densities during the day equaled 13.0 ± 1.37 μmol m−2 s−1 and the computed CO2 flux densities equaled 11.4 ± 1.38 μmol m−2 s−1.
A test of the ability of the CANVEG model to simulate Rn and LE and H heat fluxes densities over the aspen stand is shown in Fig. 1. To obtain statistically representative results for testing the model, the data were binned by hour and then averaged over the duration of the study. The model was able to reproduce the mean diurnal shape and magnitudes of measured Rn very well. The mean, daily average value of Rn measured over the canopy was 131 ± 47 W m−2, while the computed value was 139 ± 43 W m−2. With respect to LE, the model simulated the diurnal pattern and the mean, daily average value with acceptable fidelity; the mean daily value of LE measured over the canopy was 84 ± 18 W m−2, while the computed mean was 87 ± 14 W m−2. On a temporal basis, the means and standard errors associated with the model often overlapped those associated with the measurements. The model’s ability to simulate daily mean values of sensible heat flux density was only satisfactory (measured, 28.6 ± 12.4; calculated, 33.8 ± 23.6 W m−2), too, but simulation of the diurnal pattern of sensible heat flux was biased. The model underestimated H at night and overestimated it during the day, with differences on the order of 50–100 W m−2. Our tests show that the overestimation of H was compensated by an underestimation of soil heat flux density (G), and vice versa.
Stomatal conductance was computed by coupling transpiration and photosynthesis algorithms. Consequently, it is instructive to examine the ability of the CANVEG model to simulate CO2 exchange over the aspen/hazel stand. Figure 2 shows that the model was able to simulate CO2 flux densities well during the daylight period. Model calculations, on the other hand, may have overestimated the magnitude of nocturnal respiration, but this fact is uncertain.
Soil respiration rates were computed with an empirical algorithm that was determined on site by Black et al. (1996), so we have some confidence in that algorithm. Biases between nocturnal measurements and model calculations would have been derived from two other sources. One source of error is associated with the estimation of soil temperature. Although, the soil heat flux model was driven by measurements of soil temperature at 0.10 m, there is a potential to compute the near-surface soil temperature incorrectly since model calculations of G underestimated measured values. Eddy-covariance measurements of CO2 can be biased at night too. Often, an eddy-covariance system will not measure the CO2 efflux respired by soil, leaves, and boles if CO2 is stored in the underlying airspace or soil or if it drains out of the control volume under study (Lee et al. 1996). Under this circumstance, the neglect of errors associated with measurement of this storage term can cause a sizable bias between measured and modeled CO2 fluxes (Baldocchi and Meyers 1998).
Figure 3 shows a comparison between model calculations and measurements of isoprene emission rates over the aspen stand. Model calculations reproduced the pronounced diurnal pattern well. In general, both measurements and calculations of isoprene fluxes peaked around midday, achieving magnitudes of near 13 nmol (isoprene) m−2 s−1. Measured and calculated isoprene fluxes were also similar when integrated over the day. The daily mean, measured isoprene flux was 4.60 ± 1.01 nmol (isoprene) m−2 s−1, as compared to the calculated flux, which equaled 4.68 ± 1.03 nmol (isoprene) m−2 s−1. Favorable model agreement may be an artifact of the model’s ability to calculate the net radiation balance well since isoprene emission rates are a very strong function of the photosynthetic photon flux density.
In contrast to the aspen stand, isoprene emission rates measured over the temperate oak forest are more than threefold greater than the values observed over the boreal aspen stand (Fig. 4); maximum emission rates over the oak forest approached 55 nmol (isoprene) m−2 s−1. Data in Fig. 4 also show that modeled values of isoprene emission rates, over the mixed-species temperate forest stand, depend upon the chosen biomass factor. The model calculations can equal field measurements or overestimate them by 30%–40%. The highest calculations are associated with the biomass value that was representative of the average isoprene biomass of the entire forest (220 g m−2). Since the stand consists of a mixture of isoprene-emitting and nonemitting species (Guenther et al. 1996), it is more proper to compute a distance-weighted value of the isoprene-emitting biomass (Lamb et al. 1996).
The isoprene-emitting biomass (bI), sensed by a micrometeorological flux measurement system, along the wind-blown axis (x) is defined, theoretically, as
where p(x) is the source probability density function of material emitted by the surface. Mathematically, p(x) is a statistical measure of the probability that a fluid element released at a point will be measured at a given position downwind from the ground source (see Schmid 1994; Horst and Weil 1992); the integration of p(x) between zero and infinity equals one. Equation (4) was evaluated using biomass data from Lamb et al. (1996) and a flux footprint model [Baldocchi (1997b); see appendix B for additional information]).
The distance-weighted, isoprene-emitting biomass was calculated to be 146 g m−2 (Fig. 5). This value exceeds the figure of 110 g m−2 reported in Lamb et al. (1996), which was computed with an analytical footprint model. Use of this updated amount of isoprene-emitting biomass yields better agreement between measured and calculated isoprene emission rates over the mixed-species, temperate forest. The numerical footprint model has an advantage over an analytical model because it can accommodate the distinct fields of turbulence that exist within and above a forest.
With a tested canopy isoprene model at hand we use it in this section to address some relevant scientific questions. One issue involves how much carbon is lost from a forest as isoprene? Figure 6 shows theoretical estimates of the fraction of carbon gained by assimilation versus the amount of carbon lost via isoprene emission. A substantially higher high fraction of assimilated carbon is lost from the oak stand than from the aspen stand. The warmer air temperatures experienced at the temperate forest site (24.2° vs 18°C) account for a large portion of this difference.
The coordinated change of air temperature and available sunlight has a marked and nonlinear impact on isoprene emission rates from an aspen leaf (Monson and Fall 1989). How do these factors interact to influence canopy-scale isoprene fluxes? The highest emission rates occur under high light (Qp greater than 500 μmol m−2 s−1) and air temperatures exceeding 30°C (Fig. 7). This theoretical response is consistent with many field observations (e.g., Guenther et al. 1996).
The opening and closing of stomata affect photosynthesis and transpiration. How do variations in the parameterization of stomatal conductance affect isoprene emission rates? At the leaf scale, changes in stomatal conductance have no direct impact on isoprene emission rates (Sharkey et al. 1991; Fall and Monson 1992). However, variations in the stomatal conductance and its multiplication factor, m, can affect theoretical calculations of isoprene fluxes over aspen (Fig. 8). This response is indirect and is associated with a leaf’s energy balance. Higher values of isoprene emission rates are associated with lower m values because stomatal closure reduces rates of transpiration and causes the leaves to be warmer.
At this point it is also appropriate to comment on the choice of m used to model trace gas fluxes over the aspen stand. In this exercise, evaporation from the aspen stand was assumed to be emanating from an ensemble of freely transpiring leaves. In past exercises, this assumption has been found to be valid for crops and broad-leaved forests (Baldocchi 1997a; Baldocchi and Meyers 1998) when they have ample soil moisture. Under these ideal conditions the stomatal conductance coefficient (m) possesses a value near 10. For the situation under study, the stomatal conductance coefficient of the aspen leaves was about 20% lower than this optimal value. The diminished value of the stomatal conductance coefficient, m, may be chronic and attributed to the impact of tree age and low growth rates on hydraulic conductivity (Falge et al. 1996). Had either forest experienced severe water deficits, a more elaborate model would have been needed to simulate stomatal conductance and transpiration, for under stress conditions there is hysteresis between the demand and supply of water for transpiration. Under such conditions a capacitance term, due to storage of water in the boles, must be considered (Williams et al. 1996).
A trace gas exchange model (CANVEG), composed of micrometeorological and physiological modules, was applied to examine how realistic isoprene emission rates could be simulated over two contrasting forests, a uniform aspen and a mixed-species temperate forest. For the homogeneous aspen stand, modeled and measured isoprene emission rates agreed quite closely (within 2%). In contrast, estimates of isoprene emission over the mixed temperate forest were a strong function of the isoprene-emitting biomass, which in turn was a function of the wind direction and footprint detected at the tower reference point. Linking estimates of the flux footprint distribution function and biomass surveys yielded comparable values of isoprene emission between micrometeorological measurements and model computations.
The use of a coupled carbon–water–nutrient mass and energy exchange model proved to be a reliable driver for leaf-scale isoprene emission models for two dissimilar vegetation types since we were able to simulate most of the energy fluxes between the vegetation and atmosphere with reasonable fidelity. We encourage future investigations to explore the validity of using this approach for routine calculations of trace gas fluxes. An advantages of this approach is the systematic variation model parameters across biomes (see Baldocchi and Meyers 1998).
The senior author was supported by a contract from the U.S. Environmental Protection Agency, through Washington State University (Brian Lamb, PI), and the NOAA Office of Global Programs, which supported his participation in the BOREAS project. The University of Colorado team received financial support by the U.S. Environmental Protection Agency Graduate Fellowship Program, a National Science Foundation Research Training Grant, and the Department of Energy (DOE) NIGEC program.
The Canadian NOx/VOC Management Plan supported the field project to obtain isoprene fluxes at the aspen forest, which was conducted by a team of scientists from the Canadian Atmospheric Environment Service. We thank Drs. Harold Neumann, Gerry den Hartog, Daniel Wang, and Mr. John Deary for their efforts to acquire this field data. We also thank Dr. Joe Berry, Carnegie Institution of Washington, and Dr. Betsy Middleton, NASA, for information on the stomatal conductance and photosynthetic properties of aspen.
Colleagues from the National Center of Atmospheric Research and Washington State University have worked with us over the years to help characterize the structure and architecture of the Oak Ridge forest, with regards to isoprene studies. We specifically thank Drs. Alex Guenther, Lee Klinger, and Brian Lamb for their assistance and access to their data. Finally, we thank Tom Pierce, NOAA, EPA, RTP, and the referees for their constructive comments.
Funds from the DOE-NIGEC program provided the travel support for D. R. B., A. T., and R. K. M. and is gratefully acknowledged.
Overview of the CANVEG, a Canopy Trace Gas Flux Model
We used a one-dimensional, multilayer biosphere–atmosphere gas exchange model (CANVEG) to compute water vapor, CO2, and isoprene flux densities. The model consists of coupled micrometeorological and ecophysiological modules. The micrometeorological modules compute leaf and soil energy exchange, turbulent diffusion, scalar concentration profiles, and radiative transfer through the canopy. Environmental variables, computed with the micrometeorological module, in turn, drive the physiological modules that compute leaf photosynthesis, stomatal conductance, transpiration, and leaf, bole, and soil/root respiration, as well as isoprene emission rates. We discuss the salient aspects of the model system below.
The conservation budget for a passive scalar provides the foundation for computing scalar fluxes and their local ambient concentrations. If a canopy is horizontally homogeneous and environmental conditions are steady, the scalar conservation equation can be expressed as an equality between the change, with height, of the vertical turbulent flux (F) and the diffusive source/sink strength, S(C, z):
The diffusive source/sink strength of a scalar in a unit volume of leaves is proportional to the concentration gradient normal to individual leaves, the surface area of individual leaves, and the number of leaves in the volume. The diffusive source strength can be expressed in the form of a resistance–analog relationship (Meyers and Paw U 1987):
where a(z) is the leaf area density (m2 m−3), [C(z) − Ci] is the potential difference of scalar concentration or heat content between air outside the laminar boundary layer of leaves and the air within the stomatal cavity (mol mol−1), rb is the boundary layer resistance to molecular diffusion (mol−1 m2 s), and rs is the stomatal resistance (mol−1 m2 s).
The light environment on sunlit and shaded leaves in a forest is very distinct and the response of many biophysical processes to that light is highly nonlinear. Hence, we cannot evaluate Eq. (A2) on the basis of the mean light environment. Rather, we must evaluate the strengths of sources and sinks at a particular level in the canopy on the basis of the radiation balance on the sunlit (psun) and shaded (pshade) leaf fractions:
The transfer of photons through the canopy must be simulated to evaluate the probability of sunlit and shaded leaves, as well as photosynthesis, stomatal conductance, and leaf and soil energy balances. The radiative transfer model was derived from probabilistic theory (Norman 1979). The radiative transfer model assumes that foliage is randomly distributed in space and the sun is a point source. In this case the probability of beam penetration is calculated using a Poisson distribution:
where L is leaf area index, β is the solar elevation angle, and G is the foliage orientation function. Here, G represents the direction cosine between the sun and the mean leaf normal. For the ideal case, in which leaves have a spherical angle distribution, G is constant and equals one-half.
Native vegetation has clumped foliage. In these circumstances, the Poisson probability density function is inadequate for computing probabilities of photon transmission through vegetation. Instead, the Markov model can be employed (Myneni et al. 1989), where the probability of beam penetration is
Ω is a clumping factor and ranges between zero and one.
Mathematically, the probability of sunflecks is equal to the derivative of P0 with respect to L times the average cosine of the leaf–sun angle (Gutschick 1991):
The integration of Eq. (A6) with respect to leaf area yields the sunlit leaf area.
The scattering of light was computed for the visible and near-infrared wavebands using the slab, “adding” approach of Norman (1979). Downward-directed diffuse light (Di), at layer i for a given waveband, is a function of the downward directed diffuse radiation that was transmitted from the upper layer (TnDi+1) and the upward-directed radiation that is reflected (RUi):
Similarly, the upward-directed diffuse sunlight (Ui+1) is a function of that radiation that is transmitted through the layer and the reflected downward radiation:
Since we are considering the transmission and reflectance of sunlight through layers of vegetation (Δf) rather than from the surface of planes of leaves, the layer transmission and reflectance coefficients are defined as
The subscripts u and l refer to the upper and lower sides of the leaves, ρ is the leaf reflectance, and τ is the leaf transmittance for the specific waveband.
The interdependence between sources and sinks [S(C, z)] and scalar concentrations [C(z)] requires the use of a turbulent diffusion model. We used a Lagrangian random-walk model to evaluate turbulent diffusion in the vertical directions, above and within a plant canopy, using a Lagrangian framework (see Thomson 1987;Raupach 1988). Vertical displacement of fluid parcels was computed as a function of time:
where wL is the Lagrangian vertical velocity (m s−1) and dt is the differential time increment. Incremental changes in vertical velocity were computed with the Langevin equation, an algorithm that is weighted by a deterministic forcing (which is a function of the fluid parcel’s previous velocity) and a random forcing term (Thomson 1987):
The coefficients a(z, t, w) and b(z, t, w) are nonlinear functions of w and are defined to account for inhomogeneous turbulence. The term dξ defines a Gaussian random forcing with a mean of zero and a variance of dt.
The terms a(z, t, wL) and b(z, t, wL) are derived from the budget equation for the Eulerian probability density function of wL (the Fokker–Planck equation) (Thomson 1987). For the one-dimensional case, where velocity fluctuations and gradients occur only in the vertical direction, Eq. (A14) becomes
With regard to the Lagrangian model, concentration differences between an arbitrary level (Cz) and a reference level (Cr) (located above a plant canopy) were computed by summing the contributions of material diffusing to or from different layers in the canopy (Raupach 1988). Numerically, this relation is expressed as
The dispersion matrix (Dzj) has units of seconds per meter and was calculated using the random-walk algorithm of Thomson (1987). A discussion on how this random-walk model is implemented in our canopy micrometeorology model is presented elsewhere (Baldocchi 1992; Baldocchi and Harley 1995).
Leaf boundary layer resistances for molecular compounds were computed using flat plate theory. In principle such resistances, under forced convection, are a function of a leaf’s length scale (l), molecular diffusivity (d), and the Sherwood number, Sh:
In this form rb has units of seconds per meter, but it can be converted into units of moles per square meter per second with gas law conversion factors to be consistent with the inputs of the stomatal conductance and photosynthesis algorithms.
To simulate wind speed within the vegetation, which is required to assess Sh, we applied the logarithmic wind law above the canopy and the exponential wind profile within it. Under free convection, we computed the Sherwood number as a function of the Grasshof number and the leaf to air temperature difference.
Leaf temperature was calculated to determine enzymatic rates associated with carboxylation, electron transport, and respiration and to evaluate transpiration, sensible heat fluxes, and infrared emission. A quadratic equation, defining the difference between leaf and air temperature (ΔT), was derived from the leaf energy balance relationship so an analytical solution could be used to compute leaf temperature (Paw U 1987):
The coefficients are defined as
where ɛ is emissivity, σ is the Stefan–Boltzmann constant, ρa is air density, λ is the latent heat of vaporization, Tk is absolute temperature (K), gs is the stomatal conductance (m s−1), gh is the aerodynamic conductance for sensible heat transfer (m s−1), Cp is the specific heat of air, Q is absorbed energy (incoming short- and longwave radiation minus reflected shortwave radiation; W m−2), mυ and ma are the molecular weights of vapor and dry air (g mol−1), P is pressure (kPa), es is saturated vapor pressure (kPa), and ea is the ambient vapor pressure (kPa). The leaf energy balance can also be used to derive a quadratic equation for latent heat exchange (λE; W m−2):
The coefficients for Eq. (A19) are
Soil constitutes the lowest boundary of a canopy-scale, water vapor, CO2, and trace gas exchange model. Flux densities of convective and conductive heat transfer and evaporation at the soil/litter boundary and soil temperature profiles were computed using a 10-layer numerical soil heat transfer model (Campbell 1985). Surface energy fluxes were computed using an analytical solution to a surface’s energy balance. Soil evaporation was computed using an algorithm:
where Rsoil is the resistance of the soil to evaporation, ϕ is the relative humidity of the soil matrix, qa is the mixing ratio of the air, and qsat is the saturated mixing ratio. For these calculations we assumed Rsoil was 816 s m−1 and that the vapor pressure at the evaporating site was 85% of the saturation vapor pressure [ϕqsat(T)].
Ecophysiology: Photosynthesis, respiration, and stomatal conductance
The biochemical equations for the carbon exchange processes are taken from Farquhar et al. (1980). Leaf photosynthesis (A) is a function of the carboxylation (Vc), oxygenation (Vo, photorespiration), and dark respiration (Rd) rates of CO2 exchange between the leaf and the atmosphere (all have units of μmol m−2 s−1):
The term Vc − 0.5Vo is expressed by Farquhar et al. (1980) as
Adopting the minimum value between Wc, the rate of carboxylation when ribulose bisphosphate (RuBP) is saturated, and Wj, the carboxylation rate when RuBP regeneration is limited by electron transport, assesses Eq. (A22). The variable Γ is the CO2 compensation point in the absence of dark respiration (mol mol−1) and Ci is the intercellular CO2 concentration (mol mol−1). Evaluating Eq. (A22), in terms of Ci, implicitly assumes that the mesophyll conductance is infinite.
If Wc is minimal, then
In this case VCmax is the maximum carboxylation rate when RuBP is saturated and Ko and Kc are the Michaelis–Menten coefficients for O2 and CO2. If Wj is minimal, then
where J is the potential rate of electron transport (μmol e− m−2 s−1). Here, J is evaluated as a function of incident—photosynthetic photon flux density (I):
The variable α is the quantum yield (mol e− mol−1 quanta) and Jmax is the maximum rate of electron transport.
A simple conductance relation is employed to express Ci:
where Cs is the surface CO2 concentration and gs is stomatal conductance (mol m−2 s−1). Stomatal conductance was computed with the algorithm of Collatz et al. (1991), which couples it to leaf photosynthesis and relative humidity:
The coefficient m is a dimensionless slope, rh is relative humidity at the leaf surface, g0 is the zero intercept, and A (μmol m−2 s−1) is leaf photosynthesis. Finally, the system of equations and unknowns, for computing leaf photosynthesis, is closed by expressing the CO2 concentration at the leaf’s surface (Cs) in terms of the atmosphere’s CO2 concentration (Ca) and the conductance across the laminar boundary layer of a leaf (gb):
The variables, Ca and gb, are external inputs to the leaf biochemistry model and are determined from the micrometeorology of the canopy. Either numerical or analytical solutions for the coupled leaf photosynthesis–stomatal conductance model can be used to compute these fluxes.
The evaluation of some photosynthetic model parameters merits further comment. The coefficients for Jmax, VCmax, and Γ, KO2, KC, and Rd are strong, nonlinear functions of temperature (Harley and Tenhunen 1991). One temperature function used for Jmax and VCmax is
Ea is the activation energy, R is the universal gas constant, Tl is leaf temperature, and Topt is the optimum temperature. The terms ΔH and ΔS represent changes in enthalpy and entropy. The Arrenhius temperature function is used to describe temperature dependencies for Γ, KO2, KC, and Rd, with respect to a reference temperature (Tref):
Model inputs and parameters
An objective of this work is to examine how well a leaf-to-canopy integration model can simulate water vapor, sensible heat, and CO2 exchange rates over a range of environmental conditions using simple input variables and constrained parameters. Variable model inputs include photosynthetic photon flux density, air temperature, humidity, wind speed, and soil temperature at a deep reference point. The CO2 concentration is required for the photosynthesis model. For cases when CO2 data were not available we assumed it equaled 350 ppm. The key extrinsic plant input parameters are leaf area index and VCmax. Other model parameters are scaled to these two parameters, or they are representative of the vegetation’s functional type (e.g., roughness length, zero-plane displacement, canopy height).
The maximum carboxylation rate (Vcmax), the maximum rate of electron transport (Jmax), and dark respiration rate (Rd) are needed at a reference temperature to compute leaf photosynthesis and stomatal conductance. These data were derived from a prior study on oak leaves (Quercus alba). Information on the leaf photosynthesis model parameters at the oak–maple site is reported in Harley and Baldocchi (1995). Parameterization of the photosynthetic model used for the aspen was achieved with data collected by Dr. J. Berry (Carnegie Institution of Washington) and Dr. B. Middleton at the same boreal aspen stand (see Table 2).
Evaluation of the Lagrangian dispersion matrix requires information on the vertical variation in the standard deviation of vertical velocity (σw). Algorithms and parameters presented by Raupach (1988) were used to calculate the dispersion matrix for a given friction velocity. Dispersion matrices for other conditions were scaled to friction velocity.
For aerodynamic calculations, the roughness parameter was set at 10% of canopy height. For the deciduous forests, d was set at 85% of canopy height. An exponential relation was employed to calculate wind speeds within the canopy. The attenuation coefficient was set at 2.5.
The source probability density (or weighting) function used to evaluate the “flux footprint” was calculated using a numerical, random-walk Lagrangian model. Details of the model are reported in Baldocchi (1997a) and are based on theory developed by Thomson (1987) and Flesch and Wilson (1992). Flux footprints were computed by releasing a large ensemble of fluid elements (>5000) from the zero-plane displacement height of the forest, which is 85% of the canopy height. The trajectory of each fluid element was tracked and the number of elements that crossed specified heights, after a given travel distance, were counted. The sum of fluid elements captured, at each bin, was normalized by the source strength. Normalization ensured that the integral, with respect to distance (x) at each level, equaled one when integrated between zero and infinity.
Corresponding author address: Dennis Baldocchi, Atmospheric Turbulence and Diffusion Division/NOAA, 456 S. Illinois Ave., Oak Ridge, TN 37831.