Modeling forest plantations for carbon uptake with the LPJmL dynamic global vegetation model

We present an extension of the dynamic global vegetation model, Lund–Potsdam–Jena Managed Land (LPJmL), to simulate planted forests intended for carbon (C) sequestration. We implemented three functional types to simulate plantation trees in temperate, tropical, and boreal climates. The parameters of these functional types were optimized to fit target growth curves (TGCs). These curves represent the evolution of stemwood C over time in typical productive plantations and were derived by combining field observations and LPJmL estimates for equivalent natural forests. While the calibrated model underestimates stemwood C growth rates compared to the TGCs, it represents substantial improvement over using natural forests to represent afforestation. Based on a simulation experiment in which we compared global natural forest versus global forest plantation, we found that forest plantations allow for much larger C uptake rates on the timescale of 100 years, with a maximum difference of a factor of 1.9, around 54 years. In subsequent simulations for an ambitious but realistic scenario in which 650 Mha (14 % of global managed land, 4.5 % of global land surface) are converted to forest over 85 years, we found that natural forests take up 37 PgC versus 48 PgC for forest plantations. Comparing these results to estimations of C sequestration required to achieve the 2 ^∘C climate target, we conclude that afforestation can offer a substantial contribution to climate mitigation. Full evaluation of afforestation as a climate change mitigation strategy requires an integrated assessment which considers all relevant aspects, including costs, biodiversity, and trade-offs with other land-use types. Our extended version of LPJmL can contribute to such an assessment by providing improved estimates of C uptake rates by forest plantations.

1 Introduction

It is increasingly clear that the stringent climate targets of the Paris Agreement cannot be achieved without negative emissions, i.e., net removal of carbon (C) from the atmosphere later during the 21st century, to compensate for emissions in the first half of the century (Gasser et al., 2015; Rogelj et al., 2018). Of the many proposed techniques to achieve C uptake, the two options currently most discussed for large-scale implementation are bioenergy in combination with carbon capture and storage and afforestation (Williamson, 2016). Both approaches will require considerable amounts of land and thus compete with other land-use functions, for example, food production and biodiversity. While bioenergy is receiving considerable attention (van Vuuren et al., 2013), less consideration has been given to afforestation as a tool for land-based mitigation. C uptake occurs when natural vegetation is allowed to grow back on former croplands and pasture. While deliberately taking cropland or pasture out of production may involve costs, the direct management costs of natural regrowth are negligible. The carbon uptake rate of such natural regrowth, however, will usually achieve only a fraction of the potential C uptake rate at short timescales. Considerably higher C uptake rates are possible by planting forests (Paquette and Messier, 2010). Assisting regrowth by planting trees can substantially boost growth rates compared to natural forests because initial stages of primary succession (with herbaceous or shrub vegetation) are skipped and because fast-growing tree species can be selected. Moreover, trees are usually planted as saplings, cultivated under controlled conditions, which improves the chance of successful establishment compared to development from seeds (Gladstone and Thomas Ledig, 1990).

Assessing the potential of land-based approaches for climate mitigation requires reliable estimates of C sequestration rates. Process-based models, such as dynamic global vegetation models (DGVMs), are a crucial tool for providing these estimates. DGVMs simulate carbon stocks and fluxes based on mechanistic descriptions of underlying processes, such as photosynthesis and organic matter decomposition in relation to environmental conditions. However, since the focus of DGVM development has traditionally been on natural ecosystems, very few of these models have an explicit representation of planted forests. Therefore, previous modeling studies on large-scale afforestation represented afforestation as natural regrowth (Krause et al., 2017), in some cases applying corrections to account for higher growth rates (Humpenoder et al., 2014; van Minnen et al., 2008).

In this paper, we present an updated version of the Lund–Potsdam–Jena Managed Land (LPJmL) DGVM (Bondeau et al., 2007; Schaphoff et al., 2013), modified to explicitly represent afforestation. Three new plant functional types have been implemented in order to represent planted forests in temperate, tropical, and boreal regions. The parameters of these plantation types were estimated based on observations of stemwood carbon from real-world forest plantations. Using this new LPJmL version, we present a global assessment of potential carbon sequestration rates in forest plantations and compare these to rates achieved by letting forests grow back naturally.

2 Methods

2.1 The LPJmL dynamic global vegetation model

LPJmL is a global process-based model simulating vegetation dynamics and fluxes of carbon and water in the vegetation and soil of terrestrial ecosystems (Bondeau et al., 2007; Schaphoff et al., 2013; Sitch et al., 2003), including agricultural land and biomass plantations for bioenergy production (Beringer et al., 2011). The model runs primarily on a daily time step, except for C allocation, vegetation dynamics, and disturbances for natural vegetation and biomass plantations, which are resolved annually. Forcing consists of monthly climate variables (air temperature, precipitation, cloud fraction, and number of wet days per month) – which are interpolated to daily values (Gerten et al., 2004) – and annual atmospheric CO₂ concentrations. Using a combination of plant physiological relations, generalized empirically established functions, and plant trait parameters, LPJmL simulates processes such as photosynthesis, plant growth, maintenance and regeneration losses, fire disturbance, soil moisture dynamics, runoff, evapotranspiration, irrigation, and vegetation structure (Schaphoff et al., 2013). Natural vegetation is represented as a number of plant functional types (PFTs): aggregated vegetation classes representing variation in leaf type (broadleaf, needleleaf), phenology (summergreen, evergreen, raingreen), and climate preference (boreal, temperate, tropical). Most model parameters related to vegetation are defined separately for each PFT. The model simulates the occurrence of each PFT based on bioclimatic limits and competition with other PFTs for resources. Agricultural ecosystems are handled in a separate module and are represented by a range of crop functional types (Bondeau et al., 2007). Additionally, one herbaceous and two woody PFTs are implemented to simulate short-rotation bioenergy plantations (Beringer et al., 2011). Area fractions specifying allocation to different land-use types are part of the model input. Finally, the model can simulate river discharge, surface water reservoirs, and several types of irrigation. LPJmL has been coupled to the IMAGE integrated assessment model, serving as the land surface component (Müller et al., 2016; Stehfest et al., 2014).

In all simulations for this study, the model was forced by semi-constant monthly climate input, representative of the period 1980–2010. This dataset was derived by repeating a cycle of detrended time series for this period, taken from the Climate Research Unit (CRU) TS3.23 global gridded (0.5^∘ × 0.5^∘) climate dataset (Harris et al., 2014). For simplicity, we chose to ignore the effect of atmospheric CO₂ concentration change at this stage; hence, this variable was held fixed at the mean global value for 1980–2010 (362.4 ppmv). Fire disturbance was not considered. Further information on the model input and configuration is given in subsequent sections.

2.2 Forest plantations

LPJmL was extended to represent forest plantations. Specifically, a new land-use type was added, as well as three functional types to represent trees in temperate, tropical, and boreal plantations. These types – referred to as forest plantation functional types (FPFTs) – are derived from the natural PFTs (temperate broadleaved summergreen tree, tropical broadleaved evergreen tree, and boreal needleleaf evergreen tree, respectively). The occurrence of the FPFTs is subject to the same establishment and mortality rules used for natural PFTs. However, the bioclimatic limits are set such that they do not overlap; hence, co-occurrence of different FPFTs in a single grid cell is rare, occurring only when climate fluctuates near a boundary between two types.

Structurally, the implementation largely follows that of the woody bioenergy plantations implemented in LPJmL (Beringer et al., 2011), which in turn are based on equivalent natural PFTs. Contrary to bioenergy trees, forest plantations are not automatically clear cut after a fixed rotation period, but a fraction of the plantation fraction may be harvested, specified as model input. However, for the purpose of this study, harvest was set to zero. Forest plantation PFTs also differ from other PFTs with regard to establishment of new trees. A fixed initial planting density (P_init) was introduced, which determines the number of trees per unit area at planting. After planting, establishment of new trees occurs similarly to natural PFTs: at fixed maximum rate, downscaled according to an exponentially declining function of foliar projective cover. Generally, stand density will decrease after plantation due to self-thinning, implemented according to Reineke's rule (Reineke, 1933), which relates stem diameter to crown area. When the area-fraction of forest plantations in a given grid cell increases over time, establishment is determined as a combination of P_init and the standard establishment rate, weighted according to the old forest plantation fraction and the fraction added.

2.3 Calibration

2.3.1 General setup

To obtain realistic growth rates, we calibrated several FPFT-specific parameters, based on published observed growth data for forest plantations. Ideally, calibration of dynamic vegetation models should be performed using detailed observations for a given site. However, this requires a large amount of data, both for model input and to compare to model output to assess performance. While much data on growth of forest plantations have been published, the number of forest plantation sites for which calibration data as well as data for model input are available for sufficiently long time periods is not enough to derive globally applicable parameter sets. Therefore, we chose a different approach. Rather than aiming to reproduce site-level observations, we calibrated the model in order to obtain desired mean biome-level behavior for each of the three FPFTs. For every iteration in the calibration, the model was run for a selection of 100 grid cells from the spatial domain of the FPFT being calibrated. Subsequently, model output for the relevant variables was aggregated over all grid cells and compared to observed values to determine model performance.

Within the spatial domain for a given FPFT, many grid cells exist where growth is marginal due to unfavorable climate and/or soil properties. The observations used in the calibration are not representative of these locations, since forest plantations from which data have been retrieved can be assumed to represent locations where productivity is sufficient for economic profitability. Therefore, rather than choosing grid cells randomly, the selection was limited to locations for which LPJmL simulates relatively high productivity. This was done based on results from a 300-year simulation with only natural vegetation, in the same setup as the one used in the calibration (see Sect. 2.3). For each FPFT, 100 cells were selected for which the simulated stemwood C storage of the corresponding natural PFT (see Sect. 2.2) exceeds the 70th percentile over the complete domain where this PFT is dominant, i.e., has highest foliar projective cover (Fig. 1). During the calibration, LPJmL was run only for these cells, with land-use type set to forest plantations.

Figure 1 Location of the grid cells included in calibration simulations (100 per FPFT). The map shows simulated stemwood C (kgC m⁻²) for a simulation with potential natural vegetation.

2.3.2 Observations

Target growth curves

Time series of stand-level stemwood C were collected from various sources in the literature. We required observations in the form of time series for a sufficiently long period to assess the growth behavior on timescales relevant to this study – at least 50 years. Due to limited data availability (see Sect. 2.3.1), a rigorous data-selection procedure was not possible; hence, the observations were collected in an ad hoc fashion. For the tropical FPFT, we used data from Brown et al. (1986), who derived time series of stemwood biomass for several species and species groups for tropical forest plantations. For the temperate and boreal FPFTs, no such compilations were available; hence, we used datasets for typical plantation species for wood production. Observations for natural poplar (Populus × euramericana) forests were taken from Cannell (1982) for the calibration of the temperate FPFT. For the boreal FPFT, we used data for Scots pine (Pinus sylvestris) plantations from Vanninen et al. (1996). Outliers in the observations were removed using Hampel filtering (Pearson, 2002). The data are depicted in Fig. 2.

Figure 2 Target growth curves for stemwood C, associated observations, and LPJmL output for natural vegetation in the cells selected for calibration.

Since most forest plantations are grown for timber production, they are harvested approximately at the optimal rotation length for maximum wood production, which is well before the trees reach maturity. Hence, growth data for higher tree ages are scarce. Calibrating LPJmL against these observations alone would result in excessive weight on the earlier part of the curve, leading to unpredictable results for the later part. Therefore, we did not use the observations directly in the LPJmL calibration but used them to derive growth curves representing the typical growing behavior of productive plantations for each FPFT. We refer to these as the target growth curves (TGCs). The general structure of the TGCs is given by the Chapman–Richards function, which is widely used to model forest growth (e.g., Von Gadow and Hui, 1999). It defines the stemwood C (C_SWC) at time t as

$$\begin{array}{}\text{(1)}& {\mathrm{C}}_{\mathrm{SWC}}\left(t\right)={\mathrm{C}}_{\mathrm{SWC},\max }{\left(\mathrm{1}-e^{-kt}\right)}^p,\end{array}$$

where C_SWC,max is the maximum C_SWC, k is the growth rate, and p is a shape parameter determining the timing of maximum growth. These parameters were estimated using a Bayesian non-linear regression approach. The scarcity of high-age observations was solved by constraining C_SWC,max using a prior distribution based on LPJmL output for the corresponding natural PFTs from the 300-year simulation used to select the calibration grid cells (see Sect. 2.3.1). Specifically, for each FPFT, we used the mean simulated stemwood C of the last 10 simulation years, averaged over the 100 calibration cells as a representative value.

The parameters were estimated using Markov chain Monte Carlo (MCMC) sampling. The samples with highest posterior density, together with the variances over the marginal posterior distributions, were used in the LPJmL calibration. Further details are given in the Supplement.

Additional constraints

Initial tests showed that parameter sets derived by calibration with the TGCs alone result in unrealistically high values of net primary production (NPP), leading to similarly high litter fluxes and soil carbon storage. This was traced to a higher carbon use efficiency (CUE) – the ratio of NPP to gross primary productivity – and a lower vegetation carbon turnover time (τ_vegC; vegetation-C-to-NPP ratio) compared to the natural PFT counterparts. Therefore, to assure realistic carbon fluxes and storage, we implemented additional constraints for these variables in the calibration. These constraints are based on LPJmL simulations for the natural PFT counterparts of the FPFTs, similar to the maximum stemwood C of the target growth curves (see Table 2 and Fig. 4).

Additionally, it was found that certain parameter sets, while leading to acceptable mean results, cause simulated trees for certain cells to die-off repeatedly at regular intervals. In order to avoid this, we modified the calibration such that a penalty was added to the cost function when this occurs.

Table 1 LPJmL parameters included in the calibration. “Prior mode” refers to the most probable value indicated by the prior distribution.

2.3.3 Parameter estimation

In the calibration, 15 parameters were estimated, separately for each FPFT (Table 1). The calibration was performed on a transformed scale (logit for α_a; log for all other parameters), in view of the lower bound at zero (and upper bound at 1 for α_a). We applied a Bayesian cost function, including informative prior distributions. Priors express belief about reasonable parameter values before the calibration in the form of probability distributions and help to avoid unrealistic values, particularly for parameters that have little influence on the relevant model output. The priors were chosen such that their central tendency reflects existing parameter values for the corresponding natural PFTs, with a relatively wide variance to avoid overly strong influence on the calibration. Full specification of the priors is given in the Supplement (Sect. S1).

Similar to the parameters, all observations were transformed in the calibration (logit for CUE; log for all other observations). For the calibration simulations, LPJmL was started from zero vegetation and soil C and run for a period of 300 years, sufficient for the vegetation C to reach equilibrium with reasonable parameter values. LPJmL simulates heartwood and sapwood C pools but does not distinguish between stem, branches, and coarse roots. For the purpose of the calibration, we assumed that all heartwood and 66 % of the sapwood are located aboveground (Müller et al., 2016), and 84 % of aboveground wood is located in the stem (which is representative of mature trees; Pretzsch, 2010).

Table 2 Observations and corresponding fits for the 100 included grid cells included in the calibration. Observed values correspond to the mode of the likelihood function.

After simulation, the Chapman–Richards function was fitted to the time series of simulated stemwood C for the 100 grid cells (using non-linear least squares) to derive FPFT-mean estimates of C_SWC,max, k, and p based on LPJmL predictions. Carbon use efficiency and vegetation turnover time were determined for the last 10 years of the simulation, averaged over the 100 grid cells. The observations were subsequently compared to the corresponding observations to determine log-likelihood, and combined with log prior density to determine the overall cost C(θ) for the given parameter set θ. Further details are given in Sect. S2.

The optimal parameter set (with minimal value of C) was derived using the GENOUD algorithm (Mebane Jr. and Sekhon, 2011), which combines a genetic algorithm with a gradient search approach. This algorithm has previously been applied to calibrate LPJmL (Forkel et al., 2014). An additional description is given in Sect. S2.

2.4 Global simulations

After calibration, several global simulations were performed. First, in order to assess sequestration potential of afforestation, a simulation was run in the same setup as the one used for the calibration, i.e., starting with zero vegetation and soil C and with land fully allocated to forest plantations and running for 300 years so that the vegetation C pool can reach equilibrium. Additionally, a simulation with land fully allocated to natural vegetation was performed to compare natural regrowth and afforestation as land-based mitigation options.

Second, we applied the model for an ambitious scenario of large-scale afforestation, assuming that from 2015 onwards approximately 14 % of global managed land is (corresponding to 650 Mha or 4.5 % of global land surface) gradually replaced by forest plantations over the course of 85 years. This afforestation area is in line with the average land area used for land-based mitigation (both bio-energy and afforestation) in 1.5^∘ mitigation scenarios in integrated assessment models (Doelman et al., 2019b; Rogelj et al., 2018). To bring soil C to reasonable values, the simulation was initialized by two spin-up phases: (1) a 1000-year phase with natural vegetation only until 1901, and (2) a phase from 1900 to 2015 with transient cropland and pasture fractions, based on the HYDE dataset (Klein Goldewijk et al., 2010). From 2015, forest plantation area was increased and crop and pasture area was decreased, balancing each other so that total area of managed (i.e., non-natural) land remained constant. From 2100, the simulation was continued for another 50 years with constant land use. For this analysis, two complementary simulations were performed. First, a simulation where fractions of natural vegetation were increased instead of forest plantations, and second, a “baseline” simulation where land-use fractions were held constant in time from 2015. Figures S2 and S3 depict the development of land-use fractions for the three scenarios.

3 Results

3.1 Target growth curves

Figure 2 depicts the stemwood C observations for the three FPFTs, LPJmL simulations for the corresponding natural PFTs, and the TGCs resulting from the fitting procedure. The values of the maximum stemwood C (C_SWC,max), growth rate (k), and shape parameter (p) and their marginal variances are given in Table 2 (see also Fig. S1 in the Supplement). As expected, the observations show substantially higher growth rates than simulated for the natural PFTs, as well as an earlier timing of maximum growth. The TGCs represent a compromise between the observations and the C_SWC,max for natural PFTs, predicted by LPJmL: the initial high growth rate is representative of the observations, while C_SWC,max is closer to that of the simulated natural PFTs and notably lower than the level indicated by the observations.

The tropical FPFT has substantially higher C_SWC,max, approximately twice as high as the other two FPFTs. With respect to the relative growth rate k, however, the tropical FPFT is comparable to the boreal FPFT. The temperate FPFT approaches its maximum stemwood approximately 4 times faster. The boreal TGC has the highest value of p, resulting in a later timing of maximum growth.

3.2 LPJmL calibration

The parameter estimates resulting from the calibration are shown in Fig. 3, together with the range of the prior distributions. Most estimates are within interquartile range of the priors, but for several parameters the calibration resulted in relatively strong changes, in particular k_mort1, which controls mortality due to low growth efficiency. Specifically for the tropical FPFT, the estimates also clearly deviate from the prior for E_max, k_allom1, lr_max, and τ_sapwood.

Figure 3 Prior distributions and estimated values of the FPFT parameters estimated in the calibration. The box plots indicate the 5th and 95th percentiles (whiskers), the median (red line), and 25th and 75th percentiles (box) or the priors. The final parameter estimate is indicated by the asterisk (^*). See Table 1 for an explanation of the parameters.

The ranges of the observed variables are depicted in Fig. 4, together with the LPJmL predictions for the calibration grid cells, based on the optimized parameter sets. The parameter C_SWC,max is fit well by the model for all FPFTs. However, k is clearly underestimated for all three FPFTs, compared to the observed ranges. This affects the simulated growth curve for stemwood C, as shown in Fig. 5. In the LPJmL simulation, vegetation needs a longer time to reach its maximum stemwood biomass than the target growth curve. Nevertheless, the growth rate based on the optimized parameters represents a substantial improvement compared to the natural PFT counterparts. The CUE and vegetation turnover time (τ_vegC) are also reasonably well fitted.

Figure 4 Ranges of the observations used in the calibration and LPJmL estimates after calibration. The box plots indicate the 5th and 95th percentiles (whiskers), the median (red line), and 25th and 75th percentiles (box) of the likelihood function. The fitted value is indicated by the asterisk (^*). C_SWC,max: maximum stemwood C, k: growth rate; p: shape factor, CUE: carbon use efficiency (NPP-to-GPP ratio); τ_vegC: vegetation C turnover time (vegetation-C-to-NPP ratio).

Figure 5 Predicted stemwood C for 100 calibration grid cells of each FPFT based on the optimal parameter sets. Note the different scales of the y axes.

3.3 Global simulations

3.3.1 Global afforestation and natural regrowth

Figure 6 depicts the predicted spatial distribution of forest plantation functional types for a global simulation experiment with land fully allocated to forest plantations. The total area for the temperate, tropical, and boreal plantation types is 2472 Mha (10¹⁰ m²), 6242, and 3094 Mha, respectively, corresponding to 17 %, 43 %, and 22 % of global land surface. In 2579 Mha (18 %), no tree growth is simulated due to a too-cold or too-dry climate. Note that in many non-marginal regions, tree growth may still be very low due to unfavorable conditions – the depicted distribution simply results from the bioclimatic limits set in the model. The distribution the FPFTs corresponds roughly to zones C, A, and D of the first level of the Köppen climate classification (Peel et al., 2007). Since there is no type for forest plantations in arid climates, the three FPFTs extend also into desert regions.

Figure 6 Spatial distribution of the different forest plantation functional types resulting from the bioclimatic limits. In marginal regions, no trees are simulated, but grass may be present.

The global vegetation C stock over time is depicted in Fig. 7a (see also Fig. S2). Tropical plantations contribute most to C storage due to their larger area and higher productivity. Comparison with the simulation where all land is allocated to natural vegetation shows considerably faster C uptake for forest plantations (Fig. 7b), with a maximum difference of 308 PgC (193 %) after 54 years. After 300 years, global vegetation C is 102 PgC higher (112 %) for afforestation simulation. Soil and litter C storage is also proportionally higher for forest plantations. Note that the soil and litter C uptake rate is extremely high due to the fact that the simulation was started with zero C. In reality, soil C will already be present before land-use change and uptake will be much slower, possibly even negative, depending on previous land use.

Figure 7 Global total ecosystem C over time for simulations with global forest plantations or global natural vegetation. (a) Vegetation C storage per biome for the global forest plantation simulation only; (b) C storage per compartment for both simulations.

The potential for C uptake is illustrated by the mean annual increment (MAI) of vegetation C since the start of plantation (Fig. 8). There are remarkable differences between the two simulations. After an initial similar increase, MAI sharply drops after approximately 10 years for natural regrowth, while for afforestation MAI keeps rising until approximately 30 years. The behavior for natural regrowth can be explained by vegetation succession, leading to a shift from grasses to trees. This succession does not occur for forest plantations, where trees start growing immediately, resulting in a substantially higher MAI in the early part of the simulation. From spatial differences in MAI after 50 years, it is evident that tropical regions contribute most to this difference.

Figure 8 Mean ecosystem sequestration rate (mean annual increment, MAI), determined as total C storage divided by time since the start of LPJmL simulations with only forest plantations or only natural vegetation. (a) Global total versus time; (b) difference between (afforestation minus natural regrowth) after 50 years.

3.3.2 Transient afforestation and natural regrowth

Figure 9 depicts results of the global simulation scenarios with gradual increase in forest, applying either afforestation or natural regrowth. Since changes in C storage – particularly for the soil – result also from land-use changes before 2015, we focus on the difference in global C stocks compared to the baseline simulation with constant land use from 2015. Until 2015, all three simulations have very similar results, but small differences arise from the stochastic generation of daily precipitation. Gradual afforestation of 650 Mha of land between 2015 and 2100 results in 19, 48, and 75 PgC additional C storage by 2065, 2100, and 2150, respectively, versus 16, 37, and 61 PgC for natural regrowth. Most of the difference between the two simulations is due to vegetation C, but from 2100 the difference for soil C grows and would ultimately dominate, had the simulation been continued after 2150. Global C sequestration rate peaks between 2090 and 2100 at approximately 0.91 and 0.68 PgC yr⁻¹ for afforestation and natural regrowth, with average rates of 0.25 and 0.19 PgC yr⁻¹ until 2100. Around 2130, global sequestration rates are higher for the natural regrowth simulation because land use remains constant from 2100, allowing natural ecosystems to “catch up”.

Figure 9 Results of the simulations for gradual afforestation and natural regrowth. Both graphs show differences relative to the baseline simulation with constant land use. (a) Global C storage; (b) global C sequestration rate, smoothed using a 10-year moving average window.

4 Discussion

4.1 LPJmL calibration

4.1.1 Parameters changes

Compared to the prior distributions – which are largely based on values for corresponding natural PFTs – the calibration resulted in a substantial shift for several parameters. We will discuss the more notable changes. First, k_mort1 is substantially higher compared to the prior mode for all FPFTs. This parameter controls tree mortality related to low growth efficiency, which is defined as the ratio of the annual net biomass increment to leaf area. A high value of k_mort1 results in higher mortality under unfavorable conditions. The increase of this parameter can be explained by the fact that the target growth curves have substantially higher growth rates than the natural PFT equivalents, while maximum biomass is approximately the same. The apparent conflict between these two constraints can in part be resolved by increasing first-order mortality. A higher value for k_mort1 for forest plantation trees is not necessarily unrealistic since it is likely that fast-growing tree species have low tolerance for low growth conditions (Pacala et al., 1996).

The parameter k_allom3, which relates tree stem diameter to tree height, has also substantially increased for all FPFTs. Higher values of this parameter mean higher trees for the same diameter, resulting in higher maximum biomass per tree. Again, this is in agreement with field observations, which have shown a positive relationship between tree growth rate and this parameter (Martínez Cano et al., 2019)

The maximum leaf-to-root mass ratio, lr_max, is also high compared to the prior, particularly for the tropical and boreal FPFTs. This causes higher allocation of C to leaves, compared to roots, which positively affects growth rate via leaf area index and absorbed photosynthetically active radiation. Conversely, in LPJmL, lower root biomass does not reduce growth since there is no link between root biomass and water uptake. Hence, higher values of lr_max unequivocally lead to higher productivity in the model.

For the tropical FPFT, E_max, τ_sapwood, and k_allom1 have increased. E_max is the maximum water transport capacity and controls the transpiration rate. τ_sapwood is the turnover time of sapwood – higher values result in more sapwood biomass, which allows for larger leaf area. Finally, k_allom1 relates crown area to stem diameter. The generally stronger shifts for the tropical FPFT compared to the temperate and boreal ones is explained by the lower uncertainties of parameters k and p of the target growth curve (Figs. 4 and S1), which is in turn caused by the larger number of stemwood C observations (Fig. 2).

4.1.2 Fit to observations

The calibration resulted in good fits with respect to most observations, with the exception of the growth rate parameter of the target growth curves. Despite substantial improvement compared to the corresponding natural PFTs, this parameter is underestimated for all three FPFTs. As a result, predicted initial C uptake rates are lower than implied by the stemwood C observations, possibly underestimating the potential efficacy of forest plantations for climate mitigation.

As discussed in Sect. 2.3.2, we incorporated data into the calibration to constrain NPP-to-GPP ratio and vegetation C turnover time to values similar to those of the corresponding natural PFTs. Earlier calibrations, in which these constraints were not included, yielded a substantially better fit to the growth rate but with unrealistically high litter fluxes, which points to a trade-off between the fit to these observations. From a mass-balance perspective, this result is explicable: fast growth requires high NPP, which will result in high litter fluxes once vegetation reaches equilibrium biomass. This is exacerbated by the fact that we constrained maximum stemwood C (C_SWC,max) at levels close to those of the corresponding natural PFTs. The fact that fast growth results in very high litter fluxes when trees reach equilibrium relates to the fact that LPJmL does not represent certain mechanisms that lead to declining productivity with age (Zaehle et al., 2006). In reality, NPP reduction with age has been frequently observed (Ryan et al., 1997). Multiple reasons for this phenomenon have been proposed, but the leading hypothesis is that hydraulic resistance increases with tree height due to longer distance between soil and leaves (Ryan and Yoder, 1997). This results in lower photosynthesis rates and gross productivity (GPP). Since LPJmL does not include such mechanisms, it is mostly representative of mature forests. Incorporation of a more realistic representation of age dependence of forest growth rate is likely to improve the fit to observations (Zaehle et al., 2006).

4.2 Global simulations

4.2.1 Global afforestation versus natural regrowth

Despite the underestimated growth rates, our results show C can be sequestered substantially faster by forest plantations compared to natural regrowth (Fig. 7), particularly in the first 50 years following land conversion. The largest potential for plantations lies in tropical regions, which is not surprising, given that the maximum biomass of tropical FPFTs is more than twice as high compared to the temperate and boreal FPFTs. In addition to faster C sequestration, LPJmL also predicts a 12 % higher equilibrium global vegetation C pool for forest plantations, despite the fact that the FPFTs were calibrated to produce a value of C_SWC,max comparable to that of natural equivalent PFTs (Fig. 5). This contradiction is in part explained by a larger productivity of forest plantations in less productive regions, which were not included in the calibration grid cells. Another reason is a larger spatial extent of the tropical FPFT compared to that of the natural tropical PFTs.

Combined soil and litter C is also higher for forest plantations after 300 years (Fig. 7b), but its proportion to total ecosystem C (60 %) is globally almost identical to that of the natural vegetation, due to the constraints on NPP-to-GPP ratio and vegetation C turnover time included in the calibration (see Sect. 2.3.2). It is difficult to compare these results to observations for real-world plantations since studies on this topic have generally compared natural forests to tree plantations for production of wood or other products, where the effects of harvest and other management on soil C are likely considerable (Guo and Gifford, 2002; van Straaten et al., 2015). Such effects are not relevant for plantations intended for C sequestration.

4.2.2 Gradual afforestation versus natural regrowth

According to our projections, gradual conversion of 650 Mha managed land to natural forest between 2015 and 2100 results in additional C uptake of 16 and 37 PgC by 2065 and 2100, respectively. If these lands are converted to forest plantations, the estimated C uptake is 19 and 48 PgC, i.e., 19 % and 30 % higher. These should be seen as conservative estimates, in view of the underestimated growth rates resulting from the calibration. To put these numbers into perspective, we compare them to results of Gasser et al. (2015), who estimated the negative emissions needed to limit global warming to 2 ^∘C for a range of scenarios in which both the start time and the rate of reduction of greenhouse gases varied. In their most favorable scenario (energy and industry emission reduction starting in 2015 at a rate of 5 % per year), they estimated an average cumulative negative emission of 25–100 PgC is needed by 2100, compared to 450–800 PgC in the most unfavorable scenario (energy and industry emission reduction starting in 2030 at 1 % per year). Hence, large-scale forest plantations can offer a substantial contribution to climate mitigation but will likely not be sufficient.

4.3 Comparison to previous work

The results of the simulations for transient afforestation and natural regrowth compare well to results of previous studies on potential C sequestration rates of forest plantations and natural regrowth. For example, using the IMAGE integrated assessment model, van Minnen et al. (2008) performed a simulation experiment based on the Intergovernmental Panel on Climate Change (IPCC) Special Report on Emissions Scenarios (SRES) A1B scenario where 831 Mha of agricultural land is converted to permanent forest plantations between 2000 and 2100, taking into consideration land demand for food production and other uses. They estimated an additional 93 PgC can be sequestered but mostly after 2050, when land becomes gradually available due to decreasing population and increasing agricultural efficiency.

Humpenoder et al. (2014) presented a much more ambitious afforestation scenario, in which 2773 Mha of land is converted to forest plantations. The authors used maximum C storage for natural vegetation predicted by LPJmL but corrected sequestration rates using stylized growth curves for plantations in different climate regions. They estimate an additional C uptake of 192 PgC after 80 years. Roughly converting our estimate to the same land area yields a similar result (205 PgC). This similarity is not surprising, given that we used the same model, and our FPFTs were calibrated to produce the same maximum biomass as the natural PFT equivalents.

Potential sequestration rates by natural regrowth were studied by Krause et al. (2017), using the dynamic global vegetation model LPJ-GUESS. In two scenarios, derived by IMAGE and the agricultural land-use model MAgPIE, 1119 and 914 Mha were converted to natural lands, resulting in a predicted additional C uptake of 76 and 55 PgC, respectively, between 2000 and 2099. This compares well with our estimates for natural regrowth.

4.4 Model limitations

In our implementation of planted forests, the diversity of plantation tree species is reduced to three functional types with fixed properties. While the functional diversity of plantation tree species is not as vast as that of natural forests – especially in the context of C sequestration – the predictions would likely improve from implementation of additional FPFTs, particularly for the tropical biome. The model currently predicts a relatively large C storage for dry tropical zones compared to natural regrowth, which may not be fully realistic, given water limitations. Addition of a dry tropical FPFT would allow for a more accurate assessment of C sequestration in these regions.

This study does not consider the effects of climate change and CO₂ concentration on productivity of forest plantations. Although there is still considerable uncertainty regarding this topic, accounting for CO₂ fertilization will likely increase the C sequestration potential (Schimel et al., 2015) in both natural and managed forests. However, in order to properly assess this, it is important to take into account nutrient limitation to productivity as well (Norby et al., 2010).

We also did not consider possible management options that may improve C uptake rates. In particular, regular thinning can result in substantially higher C uptake rates (van Minnen et al., 2008). The model supports harvesting, but this feature was not used in this study. However, continual thinning would result in export of nutrients from the ecosystem, which would ultimately slow down growth rates, unless plantations are fertilized. Thus, representing regular harvest in LPJmL would also require representation of nutrient limitation.

4.5 Considerations beyond C uptake

Evaluation of afforestation and natural regrowth as strategies for climate change mitigation involves a range of considerations other than carbon sequestration. First, converting agricultural land to forest involves a number of costs. For both natural and planted forests, this includes price for acquiring land, while specifically for the latter costs related to establishing and maintaining the plantation are relevant (e.g., land preparation, planting of seedlings). The costs per unit C sequestered will rise with increasing area of (planted) forest, mainly due to competition for land (Doelman et al., 2019a).

Second, the positive effects of carbon uptake of changing land cover to forest can be offset due biophysical changes in the surface energy budget, related to changes in albedo, evapotranspiration, and surface roughness (Perugini et al., 2017). This may result in a net warming effect, regionally, and possibly globally, depending on the extent of land-cover change.

Third, the reduction of cropland and pasture might also have a negative impact on food security due to increased competition for land (Hasegawa et al., 2018). In order to maintain food production for the growing population, strong intensification of the agricultural sector would be required. Locally, this will result in a range of negative effects on the environment, due to higher application of fertilizers and plant protection products, as well as water extraction for irrigation (Smith et al., 2013). Furthermore, in terms of climate change mitigation, agricultural intensification will likely partially offset the benefits of afforestation and regrowth, e.g., due to higher N₂O emissions from fertilizers (Burney et al., 2010).

Finally, biodiversity is a particularly important aspect to consider, given that plantation forests have usually substantially lower species richness than primary or secondary forests (Barlow et al., 2007). A more balanced solution may be a compromise between biodiversity and C sequestration by establishing a mixture of native and plantation species, or plantation forest with a native undergrowth (Barlow et al., 2007; Bremer and Farley, 2010).

5 Conclusions

To our knowledge, the extension of LPJmL presented here represents the first model of forest plantations for C sequestration as part of a DGVM for global-scale applications. Although calibration of the model still resulted in underestimated growth rates compared to observations of stemwood C, this represents an improvement over previous approaches. According to our simulations, conversion of 650 Mha of land to forest over 85 years results in an additional C uptake of 48 PgC for forest plantation versus 37 PgC for natural regrowth, with greatest potential in the tropics. We conclude that large-scale afforestation can offer a substantial contribution to C uptake, particularly on a timescale of approximately 50–100 years. Evaluating afforestation as a strategy for climate change mitigation requires consideration of all relevant aspects in a comprehensive assessment. Our model can contribute to such an evaluation by providing improved estimates of C uptake rates.