The InTEC model is a process-based biogeochemical C-budget model and considers all major C cycle components. This model adopts a distinct approach to simulate C components by combining (a) Farquhar’s leaf-level biochemical model , (b) a soil biochemical model CENTURY modified to include forest-specific C pools such as coarse roots and woody detritus, and (c) a set of empirical NPP and age relationships derived from forest growth and yield data.

82

●
Earth System Subject Earth Surface System Synthesis

The stand age distribution A(y, i ) can be determined using a Weibull distribution:

where q is the total fire and insect occurrence frequency, is the gamma function, and s is the shape parameter (Kasischke et al., 1995). We further assume that only mature forests are harvested. At y=0 (i.e. for areas disturbed and yet to regenerate), A(0, i ) was the difference between the disturbed and planted areas in previous n years, where n is the average age that a forest may need to start to regenerate. For example, n=5 for Canada’s forests, with a range of 1–10 years (Bunce, 1989). For each subsequent year, A(y, i ) is calculated by increasing the age by one year for not disturbed forest areas, returning the age to zero for newly disturbed forest areas, and entering age one for previously disturbed forest areas which at that time begin to regenerate. In this study, all areas are treated equally, and no spatial detail is involved. With the calculated A(y, i ) and F_{NPP}(y), the overall effect of disturbances on NPP is then given by

For C_{3} plants, which include all Canada’s forests species, the instantaneous photosynthesis rate of a single leaf p is limited by the minimum of the two values (Farquhar et al., 1980; Bonan,1995; Luo et al., 1996):

where p_{1} and p_{2} are leaf gross photosynthesis limited by electron transport and rubisco activity, respectively. The values of J, V_{m}, c_{i}, , and k_{co} in Eq. (4) are given by (Sellers et al., 1992; Bonan,1995):

J_{m}=[J_{m25}(N_{l}/N_{lmax})a_{jm}^{(Ta-25)/10}]/[1+exp((85.4T_{a}-3147.7)/(T_{a}+273))].The meanings of

other terms are listed in Notation. Eq. (5) shows that p is affected by climatic variables (i.e. T_{a}, S, and h_{r}) and atmospheric variables (i.e. c_{a} and N deposition). The effect of N deposition on p is incorporated through N_{l}/N_{lmax} (see details in Section 3).

The area-averaged annual gross photosynthesis rate of a forest region in year i, P(i ), is then given by integrating p for all leaves (x) over the whole forest region (y) and time periods during the year (t)

There are many ways to carry out this integration for all leaves in a stand (Norman, 1993). One way, simple yet effective, is to stratify a canopy into sunlit and shaded leaves (Norman, 1993), since T_{a}, c_{a}, and h_{r} are more or less the same for all leaves in a canopy because canopy air is often well mixed during daytime. This stratification is essential because irradiation changes greatly for different leaves depending on their positions relative to the sun, resulting in different J values for different leaves in the canopy. With this stratification and canopy radiation models (Black et al., 1991; Chen et al., 1999a), we calculate the instantaneous canopy photosynthesis rate, p_{can}, by the minimum of

where L_{sun} is the sunlit leaf area index, given by, whereis the gap fraction at the view zenith angle , calculated as (Nilson, 1971; Chen et al.,1997). The shaded leaf area index, Lshad, thus equals the difference of total leaf area index L_{t }and L_{sun}. Assuming no significant changes in morphology of leaves occurred since the industrialization, we can calculate the change in area-averaged L_{t} in year i by

The values of J for sunlit leaves and shaded leaves, J_{sun} and J_{shad}, are calculated by replacing S in Eq. (5) with S_{sun} and S_{shad} respectively, where S_{sun} and S_{shad} are calculated following Black et al. (1991) and Chen et al. (1999a).

Assuming f_{p} is the fraction of canopy photosynthesis limited by P_{can1}, we calculate the canopy photosynthesis rate over the time period by

With p_{can}, we reduce Eq. (6) to

While it is theoretically possible to calculate P(i) for each year since the industrialization, such an operation is practically limited by data availability. An alternative is to calculate P(i) only for a recent year for which quality data are available, and to determine P(i) in other years using a relationship between the interannual relative change in P(i),(dP(i)/[P(i)di]), and the external forcing factors. Differentiating Eq. (10) gives

where term I represents the effect on dP(i)/di caused by changes in p_{can}(y, t), while term II and III represent, respectively, the effects caused by changes in forest cover area and growing season length (l_{g}). The value of dp_{can}(y, t) is given by differentiating Eq. (9):

where L_{1}, L_{2}, L_{T1,1}, L_{T1,2}, L_{T1,3}, L_{N1}, L_{L1,1}, L_{L1,2},L_{T2,1}, L_{T2,2}, L_{T2,3}, L_{T2,4}, L_{N2}, and L_{L2} are coefficients for the effects of CO_{2} fertilization, climate variability, N availability, and leaf area changes (see Appendix A). Due to the lack of historical data about changes in h_{r}, S_{sun}, and S_{shad} since industrialization, we omit their impacts in this study. All these L terms (shortened as L_{x}) and pcan vary diurnally and seasonally as well as between locations. Due to the covariance between L_{x} and pcan, L_{x} cannot be factored out of the two-dimension integration in Eq. (11). To carry out this integration, detailed data for L_{x} and p_{can} are required. In reality, this is not feasible, especially for the long historical periods in this study. To avoid this difficulty, we use a 3-step spatial and temporal scaling algorithm: (1) To replace the integration by a discrete summation; (2) To estimate the discrete summation using the correlation coefficient, r, between the two variables L_{x} and p_{can}, namely,

where n is the number of data points in terms of both time periods and spatial locations, is the spatial and temporal assemble average of L_{x} in year i, and s is the standard deviation; and (3) To introduce a conversion coefficient,,that gives

, whereis calculated using annual mean values of climate,

N availability, and CO2 concentration. Using this algorithm, we express term I as follows:

Term II in Eq. (11) is the effect of the change in forest area on the total photosynthesis. Since we consider only the existing forests, and LUC effects are outside the scope of this study. The changes in forest area due to disturbances are considered in Section 1.

Term III in Eq. (11) is the result of growing season length changes, i.e.

This term can be very important at high latitudes where the growing season is short and air temperature increased at a higher rate than that at low and middle latitudes (Frolking, 1997; Chen et al., 1999b). Inserting Eqs. (13) and (14) into Eq. (11), we calculate the interannual relative change in P(i ) by

So far, we have considered only gross photosynthesis rate P(i ). NPP is only about 25–60% of P(i ), dependent on plant species, because a large part of P(i ) is consumed by autotrophic respiration (Ryan et al., 1997). Yet, the ratio of NPP to P(i ) is conservative with climate change and N status (Ryan et al., 1997), so that