2.1 Existing model description

JULES is a community land surface model that can be run stand-alone (as described here) or used as the land surface component of the Met Office's Earth system models (Collins et al., 2011). JULES is described in Best et al. (2011) and Clark et al. (2011). JULES calculates the surface energy and water fluxes, along with gross and net primary productivity, on a half-hourly or hourly time step. The net primary productivity (NPP) for each PFT is accumulated during each time step, to be later used for calculating changes in vegetation structure and coverage in TRIFFID, the dynamic global vegetation model built into JULES. TRIFFID is called at the end of a user-defined number of days (typically 1 or 10 d), and the accumulated NPP is allocated between “growth” and “spreading”. The former is used for increasing leaf area index (LAI) and canopy height, while the latter is used to allow PFTs to take up more space in a grid cell. Competition for space is determined based on PFT heights: the tallest plants get first access to space in a grid cell, but may not be able to compete if their NPP is too low.

In JULES, crops are represented in one of two ways. Major food crops such as wheat, maize and soya are represented by the JULES-crop module (Osborne et al., 2015). However, JULES-crop is suitable only for annual seed crops, and is not compatible with TRIFFID and the wider carbon cycle representation within JULES. Therefore, the TRIFFID-crop module was developed to represent crops within the carbon cycle and climate system. When the TRIFFID-crop option is enabled within JULES, multiple types of agricultural land are represented separately. The user defines the fraction of each grid cell dedicated to food crops, pasture and bioenergy. The fractions can vary in time with new values prescribed annually or less frequently. Each of these crop area types forms a separate “land class” for which specific PFTs are allocated. TRIFFID-crop requires height-based competition (Harper et al., 2018b), which allows for a flexible number of PFTs. Each PFT is assigned to only one land class and competes only with PFTs of the same land class, within the defined fraction. Any land within the fraction that cannot be filled by the assigned PFTs is occupied by bare soil. Multiple identically parametrised PFTs may be used if the same type of plant (e.g. C₃ grass) is desired in multiple land classes (e.g. natural, food crops and pasture). TRIFFID-crop also introduces harvesting of biomass from crop areas, described in Sect. 2.2.1 as “continuous harvest”. JULES-BE describes a set of options within JULES, building upon the TRIFFID-crop functionality to enable periodic harvesting and assisted expansion of the bioenergy PFT area.

2.2 Harvesting regimes

Two methods of representing crop harvest are used. A new TRIFFID parameter, harvest_type (Table 3), may be set to 0, 1 or 2 for each PFT. A value of 0 represents no harvest; the two harvest types are described below.

2.2.2 Periodic harvest (type 2)

At defined intervals, specified in days by the user, the PFT is reduced to a short height, also specified by the user (see Table 2 for a list of parameters). New values for wood (woodC), leaf (leafC) and root (rootC) biomass are calculated based on this height, per Eqs. (46), (56)–(58) and (60) given by Clark et al. (2011) and reproduced in the Supplement. The difference between old and new above-ground carbon is allocated to the harvest pool (Eq. 3), whereas the change in root (below-ground) carbon is added to the plant litter flux (lit_c), as given in Eq. (4). A time coefficient (Δt) is used to convert stocks to fluxes.

$$\begin{array}{}\text{(3)}& {\displaystyle }{\mathrm{harvest}}_t={\displaystyle \frac{\left({\mathrm{leafC}}_{t-\mathrm{1}}+{\mathrm{woodC}}_{t-\mathrm{1}}\right)-\left({\mathrm{leafC}}_t+{\mathrm{woodC}}_t\right)}{\mathrm{\Delta }t}}\text{(4)}& {\displaystyle }{{\mathrm{lit}}_{\mathrm{c}}}_t={{\mathrm{lit}}_{\mathrm{c}}}_t+{\displaystyle \frac{({\mathrm{rootC}}_{t-\mathrm{1}}-{\mathrm{rootC}}_t)}{\mathrm{\Delta }t}}\end{array}$$

Since the model describes a constant perfect correlation between PFT height and balanced-growth LAI, minimum LAI must also be set low enough to accommodate the prescribed harvest_ht (Table 2). The PFT then begins to regrow again from its new shorter height.

Table 2TRIFFID parameters required for JULES-BE. An explanation of the use of these parameters is given in Sect. 2.2.

Download Print Version | Download XLSX

2.3 Assisted expansion

This section describes new functionality which directs the model to simulate planting of new agricultural areas. In the existing scheme, when the fractional area of a land class increases, the new area is covered by bare soil, until the existing vegetation expands into it. Expansion of PFTs in the absence of competition follows Eq. (5). Equation (5) is a simplified version of Eq. (52) in Clark et al. (2011), assuming that only one PFT is assigned to the land class, the PFT occupies at least 1 % of the total grid cell and the plant has already reached its maximum height. Cveg represents the PFT's biomass density, and g_area is a constant parameter representing total mortality.

This arrangement represents competition and growth in natural landscapes, but where land is dedicated to a specific purpose such as bioenergy crops, it is less realistic to represent it as such; it is equivalent to humans clearing an area of land for cropping but then neglecting to plant anything.

Where the agricultural areas consist of ordinary C₃ and C₄ grasses, this does not pose much of a problem since Cveg is usually small relative to NPP during the growing season; therefore, $\frac{\mathrm{NPP}}{\mathrm{Cveg}}$ can attain sufficient size to allow the grass to increase its area. The problem is more significant in the case of high-density lignocellulosic bioenergy grasses, in which NPP may be 1–3 times that of an ordinary grass but Cveg is 5–10 times larger. Annual harvesting also reduces the capacity of crop grasses to increase their area, since more of their NPP is dedicated to increasing their height (i.e. one of the assumptions of Eq. (5) does not hold for much of the year).

Therefore, in order to represent the establishment of new agricultural areas, without sacrificing the benefits of dynamic vegetation, i.e. that plants can die off where the environment is unsuitable, a new planting mechanism has been implemented. This mechanism, activated using the switch l_ag_expand globally and the ag_expand switch on individual PFTs (Table 2), alters the value of Δfrac returned by TRIFFID. Land class fractions may change once per year, whereas TRIFFID (where plant competition and fractional allocation takes place) is run once per simulation day. At each grid cell, the current land class fraction is compared to the value used at the last TRIFFID call. Where the land class fraction has increased, the assisted expansion function is activated. Δfrac is calculated as it would have been without land-use change (Δfrac_na in Eq. (6), which could be positive or negative), but then the value of the increase (Δfrac_ag) is added to it. This is equivalent to assuming that agricultural expansion is accompanied by planting new crops. Δfrac is then added to the previous PFT fraction. If two or more PFTs (for which assisted expansion is enabled) share the same land class, the new area is divided equally between them (NPFT_ag). This process is also illustrated in Fig. 1.

Figure 1Schematic of the agricultural expansion functionality in JULES-BE. A full description of the process is provided in Sect. 2.3. The area marked with ∗ represents the change in the plant functional type (PFT) area that would occur without the change in crop area (Δfrac_na in Eq. 6). The area marked with † represents the newly available agricultural area (Δfrac_ag in Eq. 6), which is immediately populated with the crop PFT where assisted expansion is enabled, or left bare where it is disabled.

Download

2.4 New PFT parametrisation

A new bioenergy PFT was developed representing Miscanthus, a perennial grass of particular interest in the bioenergy literature due to its robust growth and low input requirements (Heaton et al., 2008; Zub and Brancourt-Hulmel, 2010; McCalmont et al., 2017). An earlier representation of Miscanthus in JULES (Hughes et al., 2010) focused on realistic representation of height and LAI, and estimated yields based on NPP. In the new method of periodic harvesting, above-ground biomass (AGB) is the most important factor determining yields, and therefore this aspect was emphasised in the development of this PFT (Fig. 2d; Fig. S1 in the Supplement).

Figure 2(a, b) Modelled leaf area index (LAI) and height of Miscanthus, compared against observations at Lincolnshire, UK, for the period 2010–2013. (c) Relationship between height and LAI, model compared against observations at Lincolnshire, UK. (d) Relationship between height and above-ground biomass (AGB), with the generic equation from the model compared against observations from the UK (Christian et al., 2008), Poland (Jeżowski et al., 2011) and Italy (Cosentino et al., 2007).

Download

In the current version of JULES, around 90 PFT parameters and 13 TRIFFID parameters govern a PFT's response to its environment, although they are not all used at once because many parameters are only required by specific configurations. The Miscanthus PFT presented here was developed based on a generic C₄ grass in the 9-PFT JULES scheme (Harper et al., 2016), with 14 parameters redefined specifically for this study. Table 3 gives an overview of the main features of the Miscanthus PFT. A full list of parameters and their relevance in JULES is given in the Supplement (see also Harper et al., 2018b for further information about PFT parameter selection).

Table 3PFT parameters distinguishing the Miscanthus PFT used in this study. “C₄ Grass” parameters are taken from Harper et al. (2018b); “Miscanthus (Hughes)” are given by Hughes et al. (2010). A full list of parameters, with definitions and explanations for values used, is given in the Supplement. “–” indicates parameters that were not yet introduced in the older version of JULES used by Hughes et al. (2010). Parameters described as “Allometry” were determined via an iterative process to improve the relationships between above-ground biomass, leaf area index (LAI) and height, as described in the Supplement. Parameters described as “BETYdb” were taken from observations in the Biofuel Ecophysiological Traits and Yields database (LeBauer et al., 2018). “GPP calibration” indicates tlow was determined via an iterative process to improve the fit of modelled gross primary productivity (GPP) to the flux data obtained from the Lincolnshire site. “Litter calibration” indicates g_leaf_0 was determined via an iterative process to approximate the observed ratio of leaf litter to yield (Amougou et al., 2012). Details of these calculations are provided in the Supplement.

Download Print Version | Download XLSX

JULES-BE can represent any type of plant as a bioenergy crop. In addition to perennial grasses, short-rotation coppicing (SRC) with willow or poplar can be simulated, or softwood or hardwood trees for forestry (Table 1). This study introduces examples of tree types grown for biomass or bioenergy in Sect. 3.4, using two poplar PFTs developed for JULES by Oliver et al. (2015).

2.5 Methods of evaluation

Simulations were carried out to evaluate and illustrate the new functionality in JULES-BE. These simulations were all based on the JULES-ES configuration, a set of options designed for best representation of carbon cycle and climate dynamics over decadal to centennial timescales. All simulations began with initial conditions from a spin-up to equilibrium, then included a transient spin-up period prior to the main run.

3.1 Lincolnshire site

Model results from the Lincolnshire site are shown in Figs. 2a–c and 3, compared against observational data from the site. The seasonal cycle of growth through to harvest in mid-February is illustrated by the seasonal fluctuation of height and LAI (Fig. 2a–b). The observations show more year-to-year variation in peak seasonal height and LAI than the model. The modelled peak heights (2.4–2.55 m during 2010–2012) and LAIs (2.75–2.9) are also generally lower than those observed (height: 2.8–3.1 m; LAI: 3.1–4.1), although observed height and LAI tended to decline after their peaks to values closer to those produced by the model. The modelled crop also increased in height and LAI earlier in the season compared to observations. The correlation between observed and modelled GPP at this site is excellent (R=0.956; Fig. 3).

Figure 3Modelled gross primary productivity of Miscanthus, compared against observations at Lincolnshire, UK, for the period 2008–2012.

Download

The mean modelled yield was 6.0±0.5 t C ha⁻¹ yr⁻¹, equivalent to a dry matter (DM) yield of 12.4±1.1 t DM ha⁻¹ yr⁻¹ assuming 48 % carbon in dry biomass (Baxter et al., 2014). This significantly exceeds the observed yields of 7.6±1.6 t DM ha⁻¹ yr⁻¹ at this site (Robertson et al., 2017), though it sits squarely within the range of yields observed in the UK (12.4±5.9 t DM ha⁻¹ yr⁻¹; 11 studies compiled by Li et al., 2018a).

3.2 Modelled Miscanthus yields against observations

A comparison of yields was conducted between the JULES-BE model results and observed Miscanthus yields compiled from the literature by Li et al. (2018a). The results of this comparison are given in Figs. 4 and 5. Across all sites and years, observed yields were much more variable, with a mean ± SD of 12.5±9 t DM ha⁻¹ yr⁻¹ (n=981), compared to 14.3±7 t DM ha⁻¹ yr⁻¹ for the modelled yields (n=1360). In a few cases, yields up to 51 t DM ha⁻¹ yr⁻¹ were observed, exceeding the maximum modelled yield of 37 t DM ha⁻¹ yr⁻¹; but more significantly, low yields of less than 4 t DM ha⁻¹ yr⁻¹ were much more common in the observations (Fig. 5b).

Figure 4Comparison of modelled Miscanthus yields against observations from Li et al. (2018a).

Figure 5Modelled yields compared to observed yields of Miscanthus collated by Li et al. (2018a). In (a), the error bars give the range of data at each half-degree grid cell. The observed range (horizontal error bars) accounts for variation between sites (different sites may exist within a grid cell area), years and fertiliser or irrigation treatment; the modelled range (vertical error bars) reflects interannual variability only. Panel (b) shows the range of values by relative frequency. This figure may be compared to Fig. 3e–f from Li et al. (2018b).

Download

The modelled yields showed a consistent positive correlation with both mean annual precipitation (R=0.752) and mean annual temperature (R=0.718; Fig. S2). For wider comparison, Fig. 6 shows simulated yields of Miscanthus across the global land surface. For the observed yields, the correlation with precipitation was much weaker (R=0.094), and while correlation with mean annual temperature was weak overall (R=0.252), yield appears to peak around 14–15 ^∘C and decline with higher temperatures (Fig. S2). This difference between modelled and observed results is clearly illustrated in the southern United States, where modelled yields are as much as 20 t DM ha⁻¹ yr⁻¹ higher than observations (Fig. 4). These observations were of Miscanthus × giganteus, a cultivar that produces very high yields in temperate climates but appears less well-adapted to high temperatures (Fedenko et al., 2013). Other perennial grasses may be more appropriate for hot climates. The model PFT would benefit from some further tuning to better represent properties such as stomatal conductance and photosynthetic temperature response, particularly the tupp and vsl parameters to better calibrate the relationship between leaf temperature and maximum rate of carboxylation of Rubisco (Vcmax; Sect. S1 in the Supplement).

Figure 6Modelled yields of the Miscanthus PFT, averaged over 2010–2019.

Figure 6 shows modelled yields for the whole Earth area, averaged over 2010–2019, in order to show the general spatial pattern of productivity of Miscanthus. Yields of 8–20 t DM ha⁻¹ yr⁻¹ are typical for most temperate climates, increasing to a maximum of about 35 t DM ha⁻¹ yr⁻¹ in the humid tropics. Yields are positively correlated with both temperature and precipitation (Fig. S2). This may help to contextualise the yields shown in Fig. 4.

3.3 Assisted expansion, global and future yields

To assess the impact of the assisted expansion feature on simulated global Miscanthus crop area, Fig. 7 shows total Miscanthus crop area in the RCP2.6-SSP2 scenario (van Vuuren et al., 2017). This scenario features a rapid increase in bioenergy crop area (“Available area”; black) from 29 Mha in 2025 to 282 Mha in 2045. “Natural expansion” (green) represents the Miscanthus PFT parametrised as discussed here, without using the new agricultural expansion functionality. In this scenario, Miscanthus occupies 13 Mha of the bioenergy crop area in 2025, increasing to 104 Mha in 2045 – leaving 178 Mha as bare soil. In 2035, only 31 Mha, or 25 % of the bioenergy crop area, is occupied by Miscanthus. With “Assisted expansion” (blue), the Miscanthus PFT occupies a consistently larger proportion of the available area throughout this period of rapid increase. In 2035, the PFT covers 119 Mha, 96 % of the available area. The proportion of area covered begins to decline after 2040, as the bioenergy production area shifts from the tropics into temperate biomes which are somewhat less favourable for growth in this representation of Miscanthus. The difference in crop area between the old and new expansion methods declines toward the end of the simulation, as the crop area begins to stabilise and the two simulations begin to converge.

Figure 7Modelled area of Miscanthus under RCP2.6-SSP2, showing the effect of the agricultural expansion functionality. “Assisted expansion” shows the model run using the assisted expansion function, compared to “Natural expansion” in which this function is disabled. “Available area” shows the total available area for bioenergy crops under this scenario.

Download

In Fig. 8, the total global Miscanthus yield is shown, using the “assisted expansion” method shown in Fig. 7. The bioenergy crop yield supplied in the IMAGE model is shown for reference (Huppmann et al., 2018; Doelman et al., 2018; Daioglou et al., 2019). Following the rapid increase in bioenergy crop area, from 2040 to 2099, bioenergy crop yields remain fairly steady in JULES-BE at 4.3 Gt DM yr⁻¹ globally, compared to 3.3 Gt DM yr⁻¹ in IMAGE over the same period. IMAGE uses a management factor when projecting energy yields, assuming that yields are currently used inefficiently (typical values are 60 % in 2020) but that improvements to crop breeding and management will increase yields to 120 %–140 % of physical potential by 2100 (Stehfest et al., 2014). This accounts for a portion of the gap in the early years of this scenario, which closes between the two models by the 2090s. The Miscanthus PFT in JULES-BE probably overestimates yields in hot climates (Fig. 4); as such, the yields projected by IMAGE may be more reliable. This scenario, and the comparison between JULES-BE and IMAGE, will be explored in greater detail in a future publication.

Figure 8Bioenergy crop yield from Miscanthus under the land-use scenario RCP2.6-SSP2, compared to equivalent bioenergy yield in IMAGE (Huppmann et al., 2018).

Download

3.4 Demonstrations of forestry and short-rotation coppicing

Figure 9 shows illustrative simulations of short-rotation coppicing and managed forestry using JULES-BE. Over the 20 harvest cycles of poplar SRC, the yield was 2.4±0.3 t C ha⁻¹ yr⁻¹ (P. Nigra) and 2.2±0.5 t C ha⁻¹ yr⁻¹ (P. × euramericana). This falls within the range observed by Sabbatini et al. (2016) over the 2011–2012 growing seasons (3.1±1.5 t C ha⁻¹ yr⁻¹) at the IT-CA1 site (growing Populus × canadensis on a 2-year coppicing rotation). The site received some supplemental irrigation during dry spells, which is not represented in the model; this may account for some underestimation of yields. Residue harvesting based on wood litter produced generally small yields of 0.15–0.25 t C ha⁻¹ yr⁻¹ in addition to forest carbon stock accumulation of 45–80 t C ha⁻¹ over the 60-year period. For rotation forestry, the yield over the 40-year rotation was 41 t C ha⁻¹ for broadleaf and 69 t C ha⁻¹ for needleleaf, equivalent to 1.0 and 1.7 t C ha⁻¹ yr⁻¹ respectively. This is higher than the average productivity for European forests (around 0.8 t C ha⁻¹ yr⁻¹, assuming 250 kg C m⁻³ of harvested roundwood; Payn et al., 2015), but lower than recent estimates from France for Douglas fir of 3.1 t C ha⁻¹ yr⁻¹ following a 40-year rotation (Bréda and Brunette, 2019). These examples show that with appropriate tuning and validation of the PFT and harvest parameters, JULES-BE could be used to facilitate decision-making on questions such as species selection, harvesting regime, harvest frequency and timing.

Figure 9Illustrated cumulative harvests and vegetation regrowth over a 60-year period, showing the application of JULES-BE to short-rotation coppicing (a, b); permanent forestry with harvest of wood litter (c, d), and rotation forestry (e, f).

Download

4.1 Main findings and limitations

The modelled yields of Miscanthus were broadly consistent with observations from sites in the USA and Europe, but showed much less variability. A major reason for this is that the harvest frequency is fixed in the model, with no option for irregular frequency or for harvests to be skipped. For example, in practice Miscanthus is generally allowed 1–2 years after planting to establish before being harvested annually, followed by 1–2 years of low yields. The largest yields generally occur during years 4–10 and decline thereafter, with a typical rotation length of 20 years (Zub and Brancourt-Hulmel, 2010). In the model, there is no representation of a plant's age, so it is not possible to establish an age-dependent harvest regime. Another reason for reduced variability in yields is that the root system reverts to the same small size after each harvest (Eq. 4), dropping its surplus biomass into the soil C pool. In reality, a relatively small proportion of root biomass is shed at harvest, and the mature plant gets a regrowth benefit from an established root system. The model currently relies on a fixed relationship between above-ground height and root biomass, and breaking this link would create other problems in the model relating to PFT scaling. Future versions of JULES will use the Reduced Ecosystem Demography (RED) approach, which represents separate mass classes within a PFT (Moore et al., 2018). Alternatively, an approach could be implemented similar to that of Black et al. (2012), in which three PFTs are used to represent different age classes of sugarcane, although this would not be compatible with dynamic vegetation. Given these difficulties, and the fact that JULES is a global model, accurate average yields with reduced variability compared to observations is likely to be an acceptable compromise for most applications of JULES-BE.

The Miscanthus PFT has not been tested with other advanced modules within TRIFFID, such as nitrogen cycling or layered soil carbon (Burke et al., 2017), and will likely require additional updating and tuning of parameters to yield useful results with other functions. Since nitrogen content is recorded for the harvested biomass, with appropriate tuning JULES-BE could also be used to quantify nitrogen loss from bioenergy crop ecosystems due to harvesting.

4.2 Further work

An example of rotation forestry has been shown in Fig. 9 for a single point. To represent forestry on a country, regional or global scale, further development of the model is required. The harvest frequency and timing are currently fixed for each PFT, meaning that all grid cells are harvested at the same time. Over a large number of grid cells, this would not be realistic and would produce undesirable hydrological and climatic effects. Further improvements to the model could enable the user to stagger the timing of harvesting. Allowing harvest frequency to vary regionally would better represent rotation forestry and increase yield by enabling the user to choose a regionally appropriate harvest frequency (shorter for more productive regions). Allowing harvest day-of-year to vary regionally would improve global-scale assessment of any bioenergy crop, since harvest timing is dependent on local climatology and affects local land surface properties, such as roughness length, albedo and transpiration rate, which in turn affect the climate. This functionality may be best implemented by allowing these variables to be user-prescribed for each grid cell. However, providing these data may be burdensome for the user, and some predictive algorithms based on climatology and growth, built into the model, may be more appropriate.

The algorithms for competition between PFTs within TRIFFID can potentially be used to determine the most suitable type of bioenergy crop in each grid cell. However, some modifications would need to be made to the existing code. In the simulations presented in this study, TRIFFID competition was enabled, allowing the bioenergy PFT to adjust its area to scale with its productivity – for example, allowing a crop to die back in response to an unsuitable environment. The current competition scheme is not useful for allowing different types of bioenergy PFTs to compete with each other within a grid cell, since it is based on height. This favours plants that can gain height easily, rather than shorter species with greater biomass density. A competition scheme based on above-ground biomass rather than height would be the first modification to make. This could help select between species within a harvesting regime, e.g. helping determine the best perennial grass for annual harvesting, or the best tree species for short-rotation coppicing. However, this may not necessarily select for the highest-yielding plant, because above-ground biomass is only a good proxy for yield at the end of the growing season, and competition for area is invoked at every iteration of TRIFFID (once per day in these simulations). Also, this development would not be useful for mixing PFTs with different harvest frequency or harvest day-of-year, since it would continue to bias competition towards the PFT that has been harvested less recently. The best solution would be to reapportion the bioenergy crop area between PFTs once per harvest cycle, based on the previous cycle's yield, but that would be a complex development given the existing model structure. Ultimately however, a yield-based competition scheme would still ignore the biophysical, economic and environmental factors that influence choice of crop type. As such, JULES-BE may always be more useful for informing these land-use decisions based on its output, rather than integrating these decisions into the existing model.