Dynamics of Biotic Carbon Fluxes under Different Scenarios of Forest Area Changes

Abstract

This paper estimates the global biotic carbon fluxes into the atmosphere under various scenarios of changes in forest area for temperate/boreal and tropical zones. Forestry scenarios have been developed from both national forest land-cover inventories and various estimates for tree vegetation-cover areas derived from the Earth’s remote sensing data for recent decades. Three scenarios have been proposed for changes in the area of temperate/boreal forests: FAOSTAT (a growth scenario on the basis of inventories), CONST (a scenario of maintaining the present level), and LANDSAT (a scenario on the basis of satellite data from different sources). The trend of tropical deforestation is maintained for FAOSTAT and LANDSAT scenarios. Together with GEPL97 (the base scenario published earlier by the authors), these scenarios have been used to obtain model biotic carbon fluxes. A model of the global carbon cycle developed in MPEI and recently updated to take into account anthropogenic carbon emissions (both industrial and related to deforestation/land use) has been used to perform a thorough analysis of the plausible dynamics of biotic carbon sources and sinks. It has been shown that the magnitude and sign of biotic carbon fluxes into the atmosphere depends substantially on the adopted scenario—from preserving the biosphere as an effective sink of excess carbon from the atmosphere in the current century (GEPL97 and CONST) to its turning into an additional source of CO₂ (LANDSAT). The FAOSTAT scenario leads to the return of biotic carbon fluxes by 2100 to the modern level. The maximum difference in net carbon fluxes from the biosphere to atmosphere for the scenarios considered here is almost 2.5 Gt C/year and falls approximately on the middle of the current century. The difference between the corresponding CO₂ concentrations for these scenarios reaches 65 ppm by 2100.

UPGRADED VERSION OF THE GLOBAL CARBON CYCLE MODEL

The box-diffusion model constructed from the well-known model of the Bern group [34] describes the ocean as a surface (homogeneous) layer and deep layers where the carbon propagation satisfies the one-dimensional diffusion equation. The model also incorporates the upwelling of cold and deep water to the surface in the region of high latitudes, where the carbon exchange with the atmosphere is more intense. The key carbon reservoirs on the surface are soil humus and vegetation (especially forests). The natural exchange of atmospheric CO₂ with the biosphere is superimposed by anthropogenic emissions associated with intensive forest clearance and subsequent biomass burning. The model includes a five-component biosphere module (Fig. A1) representing the dynamics of carbon exchange by balance-type ordinary differential equations, which form a subsystem of the model. Its output parameter Q_b is important for the correct reproduction and prediction of changes in the atmospheric CO₂ concentration.

The notations and main equations of this model are as follows: let ν be the change in the carbon content N relative to the preindustrial level N^o (ν_а(t) in the atmosphere, ν_m(t) in the mixed layer, and ν_d(t, z), 0 ≤ z ≤ h_d, t > 0 in the deep ocean); then, we can write for the carbon balance in the reservoirs

$$\begin{gathered} \frac{{d{{\nu }_{a}}}}{{dt}} = \frac{{{{Q}_{a}}}}{{N_{a}^{o}}} + \frac{{{{Q}_{b}}}}{{N_{a}^{o}}} + {{k}_{{am}}}({{\xi }_{m}}{{\nu }_{m}} - {{\nu }_{a}}) + {{k}_{{ad}}}({{\xi }_{d}}\overline {{{\nu }_{d}}} - {{\nu }_{a}}), \\ \frac{{d{{\nu }_{m}}}}{{dt}} = {{k}_{{ma}}}({{\nu }_{a}} - {{\xi }_{m}}{{\nu }_{m}}) + \frac{1}{{(1 - {{a}_{c}}){{h}_{m}}}}K{{\left. {\frac{{\partial {{\nu }_{d}}}}{{\partial z}}} \right|}_{{z = 0}}}, \\ \frac{{\partial {{\nu }_{d}}}}{{\partial t}} = K\frac{{{{\partial }^{2}}{{\nu }_{d}}}}{{\partial {{z}^{2}}}} + {{k}_{{da}}}({{\nu }_{a}} - {{\xi }_{d}}{{\nu }_{d}}), \\ 0 \leqslant z \leqslant {{h}_{d}},\,\,\,\,t > 0, \\ \end{gathered} $$

(A1)

where t = 0 corresponds to the preindustrial state and is characterized by the stationary distribution of carbon. The initial conditions are

$${{\nu }_{a}}(0) = {{\nu }_{m}}(0) = {{\left. {{{\nu }_{d}}} \right|}_{{t = 0}}} = 0.$$

(A2)

For the diffusion equation, the following boundary conditions are given:

$${{\left. {{{\nu }_{d}}} \right|}_{{z = 0}}} = {{\nu }_{m}},\,\,\,\,{{\left. {\frac{{\partial {{\nu }_{d}}}}{{\partial z}}} \right|}_{{z = {{h}_{d}}}}} = 0.$$

(A3)

Here, k_am, k_ma, k_ad, k_da are exchange coefficients; ξ_m(t), ξ_d(t) are buffer coefficients; K is the eddy diffusivity in the deep ocean; ${{\bar {\nu }}_{d}} = \int_0^{{{h}_{d}}} {{{\nu }_{d}}\left( {t,z} \right)dz} {\text{,}}$Q_a(t) is the industrial carbon flux into the atmosphere; Q_b(t) is the net carbon flux from the biosphere to the atmosphere; and Q_a(t) = Q_b(t) = 0 for the stationary state.

The buffer coefficients ξ_m and ξ_d are found from empirical equations (A4):

$$\begin{gathered} {{\xi }_{m}} = \xi _{m}^{0} + 4.9{{\xi }_{m}}{{\nu }_{m}} - 0.1{{\left( {{{\xi }_{m}}{{\nu }_{m}}} \right)}^{2}}, \\ {{\xi }_{d}} = \xi _{d}^{0} + 4.9{{\xi }_{d}}{{\nu }_{d}} - 0.1{{\left( {{{\xi }_{d}}\overline {{{\nu }_{d}}} } \right)}^{2}}. \\ \end{gathered} $$

(A4)

The exchange coefficients are

$$\begin{gathered} {{k}_{{am}}} = \left( {1 - {{a}_{c}}} \right){{F_{{am}}^{0}} \mathord{\left/ {\vphantom {{F_{{am}}^{0}} {\left( {{{h}_{a}}{{C}_{0}}} \right)}}} \right. \kern-0em} {\left( {{{h}_{a}}{{C}_{0}}} \right)}}, \\ {{k}_{{ma}}} = {{F_{{am}}^{0}} \mathord{\left/ {\vphantom {{F_{{am}}^{0}} {\left( {{{h}_{m}}{{C}_{0}}} \right)}}} \right. \kern-0em} {\left( {{{h}_{m}}{{C}_{0}}} \right)}}, \\ {{k}_{{ad}}} = {{{{a}_{c}}F_{{ad}}^{0}} \mathord{\left/ {\vphantom {{{{a}_{c}}F_{{ad}}^{0}} {\left( {{{h}_{a}}{{C}_{0}}} \right)}}} \right. \kern-0em} {\left( {{{h}_{a}}{{C}_{0}}} \right)}}, \\ {{k}_{{da}}} = {{{{a}_{c}}F_{{ad}}^{0}} \mathord{\left/ {\vphantom {{{{a}_{c}}F_{{ad}}^{0}} {\left( {{{h}_{d}}{{C}_{0}}} \right)}}} \right. \kern-0em} {\left( {{{h}_{d}}{{C}_{0}}} \right)}}, \\ \end{gathered} $$

(A5)

where a_c is the relative area of deep water contact with the atmosphere in polar latitudes; $N_{a}^{0}$ is the preindustrial carbon content in the atmosphere; C₀ is the preindustrial concentration of inorganic carbon in surface waters of the ocean; A_oc is the ocean surface area; ${{h}_{a}} = {{N_{a}^{0}} \mathord{\left/ {\vphantom {{N_{a}^{0}} {{{C}_{0}}{{A}_{{oc}}}}}} \right. \kern-0em} {{{C}_{0}}{{A}_{{oc}}}}}$ is the thickness of the oceanic surface layer containing as much carbon as the preindustrial atmosphere; and $F_{{am}}^{0}$ and $F_{{ad}}^{0}$ are the stationary (preindustrial) fluxes between the atmosphere and the mixed layer, as well as between the atmosphere and the ocean in cold surface waters in polar latitudes, respectively. For the stationary state, the direct and reverse fluxes are equal in value. The current version of the model uses the relation $F_{{am}}^{0} = \lambda F_{{ad}}^{0}{\text{,}}$ where the coefficient λ is determined in the course of numerical experiments. The preindustrial flux from the atmosphere to the ocean is $F_{{as}}^{0} = \left( {1 - {{a}_{c}}} \right)$ × $F_{{am}}^{0} + {{a}_{c}}F_{{ad}}^{0}{\text{,}}$ and the exchange coefficients are calculated through $F_{{as}}^{0}$ using the formulas

$$\begin{gathered} {{k}_{{am}}} = \left( {1 - {{a}_{c}}} \right){{{{F_{{as}}^{0}} \mathord{\left/ {\vphantom {{F_{{as}}^{0}} {\left( {\lambda + \left( {1 - \lambda } \right){{a}_{c}}} \right)}}} \right. \kern-0em} {\left( {\lambda + \left( {1 - \lambda } \right){{a}_{c}}} \right)}}} \mathord{\left/ {\vphantom {{{{F_{{as}}^{0}} \mathord{\left/ {\vphantom {{F_{{as}}^{0}} {\left( {\lambda + \left( {1 - \lambda } \right){{a}_{c}}} \right)}}} \right. \kern-0em} {\left( {\lambda + \left( {1 - \lambda } \right){{a}_{c}}} \right)}}} {{{h}_{a}}{{C}_{0}}}}} \right. \kern-0em} {{{h}_{a}}{{C}_{0}}}}, \\ {{k}_{{ad}}} = {{{{{{a}_{c}}F_{{as}}^{0}} \mathord{\left/ {\vphantom {{{{a}_{c}}F_{{as}}^{0}} {\left( {\lambda + \left( {1 - \lambda } \right){{a}_{c}}} \right)}}} \right. \kern-0em} {\left( {\lambda + \left( {1 - \lambda } \right){{a}_{c}}} \right)}}} \mathord{\left/ {\vphantom {{{{{{a}_{c}}F_{{as}}^{0}} \mathord{\left/ {\vphantom {{{{a}_{c}}F_{{as}}^{0}} {\left( {\lambda + \left( {1 - \lambda } \right){{a}_{c}}} \right)}}} \right. \kern-0em} {\left( {\lambda + \left( {1 - \lambda } \right){{a}_{c}}} \right)}}} {{{h}_{a}}{{C}_{0}}}}} \right. \kern-0em} {{{h}_{a}}{{C}_{0}}}}. \\ \end{gathered} $$

(A6)

The value of $F_{{as}}^{0}$ is estimated from the condition of equilibrium for the preindustrial atmosphere-to-ocean flux of the radioactive carbon isotope ¹⁴C and its radioactive decay in the ocean. The respective value of the eddy diffusivity K is determined from the preindustrial profile of the ¹⁴C concentration in the deep ocean.

The partition of the biosphere into reservoirs, their carbon content, and characteristic turnover times were chosen based on the model described in [35]. The change in the carbon content N_i in the reservoirs of the biotic module is calculated from the system

$$\begin{gathered} \frac{{d{{N}_{W}}}}{{dt}} = {{F}_{{aW}}} - {{F}_{{WL}}} - \frac{3}{2}{{F}_{{def}}}\frac{{1 + \gamma {{\nu }_{a}}\left( {1 - {{{{\tau }_{W}}} \mathord{\left/ {\vphantom {{{{\tau }_{W}}} 2}} \right. \kern-0em} 2}} \right)}}{{1 + \gamma {{\nu }_{a}}\left( {{{t}_{{1980}}}} \right)}} + {{F}_{{aff}}}, \\ \frac{{d{{N}_{G}}}}{{dt}} = {{F}_{{aG}}} - {{F}_{{GL}}},\,\,\, \\ \,\frac{{d{{N}_{L}}}}{{dt}} = {{F}_{{GL}}} + {{F}_{{WL}}} - {{F}_{{La}}} - {{F}_{{LH}}}, \\ \frac{{d{{N}_{H}}}}{{dt}} = {{F}_{{LH}}} - {{F}_{{HS}}} - {{F}_{{Ha}}},\,\,\, \\ \,\frac{{d{{N}_{S}}}}{{dt}} = {{F}_{{HS}}} + \,\,\frac{1}{2}{{F}_{{def}}}\frac{{1 + \gamma {{\nu }_{a}}\left( {1 - {{{{\tau }_{W}}} \mathord{\left/ {\vphantom {{{{\tau }_{W}}} 2}} \right. \kern-0em} 2}} \right)}}{{1 + \gamma {{\nu }_{a}}\left( {{{t}_{{1980}}}} \right)}} - {{F}_{{Sa}}} - {{F}_{{soil}}}, \\ \end{gathered} $$

(A7)

where τ_i is the characteristic time of carbon turnover in reservoir i, F_ij is the flux from reservoir i to reservoir j, F_def is the carbon flux into the atmosphere from reservoir W due to deforestation corresponding to bioproductivity in 1980, and F_aff is the carbon flux into the atmosphere from reservoir W due to afforestation. The net carbon flux from the biosphere to the atmosphere is

$$\begin{gathered} {{Q}_{b}} = {{F}_{{ba}}} - {{F}_{{ab}}} = {{F}_{{La}}} + {{F}_{{Ha}}} + {{F}_{{Sa}}} + {{F}_{{def}}} \\ \times \,\,\frac{{1 + \gamma {{\nu }_{a}}\left( {1 - {{{{\tau }_{W}}} \mathord{\left/ {\vphantom {{{{\tau }_{W}}} 2}} \right. \kern-0em} 2}} \right)}}{{1 + \gamma {{\nu }_{a}}\left( {{{t}_{{1980}}}} \right)}} + {{F}_{{soil}}} - {{F}_{{aW}}} - {{F}_{{aG}}} - {{F}_{{aff}}}{\text{.}} \\ \end{gathered} $$

(A8)

Here,

$$\begin{gathered} {{F}_{{ba}}} = {{F}_{{La}}} + {{F}_{{Ha}}} + {{F}_{{Sa}}} \\ + \,\,{{F}_{{def}}}\frac{{1 + \gamma {{\nu }_{a}}\left( {1 - {{{{\tau }_{W}}} \mathord{\left/ {\vphantom {{{{\tau }_{W}}} 2}} \right. \kern-0em} 2}} \right)}}{{1 + \gamma {{\nu }_{a}}\left( {{{t}_{{1980}}}} \right)}} + {{F}_{{soil}}}, \\ \end{gathered} $$

(A9)

$${{F}_{{ab}}} = {{F}_{{aW}}} + {{F}_{{aG}}} + {{F}_{{aff}}}.$$

(A10)

The model was adjusted to refine the coefficient γ of bioproductivity growth for doubled atmospheric concentration of CO₂, as well as the coefficient that takes into account the change in the mass of long-living biota with changing temperature.

The model parameters and coefficients can be determined either from the results of direct measurements or by calibration according to the preindustrial distribution of radiocarbon ¹⁴C; the system of similar equations for ¹⁴C is solved simultaneously with the system for ¹²C. However, the relative area of deep water contact with atmosphere a_c cannot be estimated in this way, because the regulation of this area in the model allows the infinitely fast horizontal mixing to be compensated (according to [34], the value of a_c should be set in the range from 0.05 to 0.15). Therefore, the further adjustment of the model involves data on the atmospheric CO₂ content as well. The rate of change in atmospheric CO₂ concentration is reproduced by the model by fitting the coefficient a_c. In the previous version of the model, the choice of the pair a_c = 0.065 and the preindustrial carbon concentration in the atmosphere $P_{a}^{0}$ = 277 ppm failed to provide the best agreement between the simulation results and experimental data. In the current version, a_c = 0.15 and $P_{a}^{0}$ = 277 ppm, and the adjustment with respect to data on the atmospheric CO₂ concentration is performed by fitting the coefficient λ from the relation $F_{{am}}^{0} = \lambda F_{{ad}}^{0}$ and the coefficient γ of bioproductivity growth for doubled atmospheric concentration of CO₂. In addition, the fitted coefficient λ controls the amount of carbon sink into the ocean.