# A vertically discretised canopy description for ORCHIDEE (SVN r2290) and the modifications to the energy, water and carbon fluxes

K. Naudts, J. Ryder, M. J. Mcgrath, J. Otto, Y. Chen, A. Valade, Valentin

Bellassen, G. Berhongaray, G. Bönisch, M. Campioli, et al.

► To cite this version:

K. Naudts, J. Ryder, M. J. Mcgrath, J. Otto, Y. Chen, et al.. A vertically discretised canopy description for ORCHIDEE (SVN r2290) and the modifications to the energy, water and carbon fluxes. Geoscientific Model Development, 2015, 8 (7), pp.2035. 10.5194/gmd-8-2035-2015. hal-01196479

HAL Id: hal-01196479

https://hal.science/hal-01196479

Submitted on 27 May 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

<figure><img src="image_3.png" /></figure>

# A vertically discretised canopy description for ORCHIDEE (SVN r2290) and the modifications to the energy, water and carbon fluxes

K. Naudts¹,¹⁴, J. Ryder¹, M. J. McGrath¹, J. Otto¹,¹⁰, Y. Chen¹, A. Valade¹, V. Bellasen², G. Berhongaray³, G. Bönisch⁴, M. Campioli³, J. Ghattas¹, T. De Groote³,¹¹, V. Haverd⁵, J. Kattge⁴, N. MacBean¹, F. Maignan¹, P. Merilä⁶, J. Penuelas⁷,¹², P. Peylin¹, B. Pinty⁸, H. Pretzsch⁹, E. D. Schulze⁴, D. Solyga¹,¹³, N. Vuichard¹, Y. Yan³, and S. Luyssaert¹

¹LSCE, IPSL, CEA-CNRS-UVSQ, 91191 Gif-sur-Yvette, France

²INRA, 21079 Dijon, France

³University of Antwerp, 2610 Wilrijk, Belgium

⁴MPI-Biogeochemistry, Jena, Germany

⁵CSIRO-Ocean and Atmosphere Flagship, 2600 Canberra, Australia

⁶METLA, Oulu, Finland

⁷CSIC, Global Ecology Unit CREAF-CSIC-UAB, Cerdanyola del Valles, Spain

⁸European Commission, Joint Research Centre, Ispra, Italy

⁹TUM, Munich, Germany

¹⁰Helmholtz-Zentrum Geesthacht, Climate Service Center 2.0, Hamburg, Germany

¹¹VITO, 2400 Mol, Belgium

¹²CREAF, Cerdanyola del Vallès, Spain

¹³CGG, 91341 Massy, France

¹⁴MPI-Meteorology, Hamburg, Germany

Correspondence to: K. Naudts (kim.naudts@lsce.ipsl.fr)

Received: 29 October 2014 – Published in Geosci. Model Dev. Discuss.: 05 December 2014

Revised: 04 May 2015 – Accepted: 22 May 2015 – Published: 13 July 2015

**Abstract.** Since 70% of global forests are managed and forests impact the global carbon cycle and the energy exchange with the overlying atmosphere, forest management has the potential to mitigate climate change. Yet, none of the land-surface models used in Earth system models, and therefore none of today's predictions of future climate, accounts for the interactions between climate and forest management. We addressed this gap in modelling capability by developing and parametrising a version of the ORCHIDEE land-surface model to simulate the biogeochemical and bio-physical effects of forest management. The most significant changes between the new branch called ORCHIDEE-CAN (SVN r2290) and the trunk version of ORCHIDEE (SVN r2243) are the allometric-based allocation of carbon to leaf, root, wood, fruit and reserve pools; the transmittance, absorbance and reflectance of radiation within the canopy; and

the vertical discretisation of the energy budget calculations. In addition, conceptual changes were introduced towards a better process representation for the interaction of radiation with snow, the hydraulic architecture of plants, the representation of forest management and a numerical solution for the photosynthesis formalism of Farquhar, von Caemmerer and Berry. For consistency reasons, these changes were extensively linked throughout the code. Parametrisation was revisited after introducing 12 new parameter sets that represent specific tree species or genera rather than a group of often distantly related or even unrelated species, as is the case in widely used plant functional types. Performance of the new model was compared against the trunk and validated against independent spatially explicit data for basal area, tree height, canopy structure, gross primary production (GPP), albedo and evapotranspiration over Europe. For all tested variables,
ORCHIDEE-CAN outperformed the trunk regarding its ability to reproduce large-scale spatial patterns as well as their inter-annual variability over Europe. Depending on the data stream, ORCHIDEE-CAN had a 67 to 92% chance to reproduce the spatial and temporal variability of the validation data.

# 1 Introduction

Forests play a particularly important role in the global carbon cycle. Forests store almost 50% of the terrestrial organic carbon and 90% of vegetation biomass (Dixon et al., 1994; Pan et al., 2011). Globally, 70% of the forest is managed and the importance of management is still increasing both in relative and absolute terms. In densely populated regions, such as Europe, almost all forest is intensively managed by humans. Recently, forest management has become a top priority on the agenda of political negotiations to mitigate climate change (Kyoto Protocol, http://unfccc.int/resource/docs/convkp/kpeng.pdf). Because forest plantations may remove CO₂ from the atmosphere, if used for energy production, harvested timber is a substitute for fossil fuel. Forest management thus has great potential for mitigating climate change, which was recognised in the United Nations Framework Convention on Climate Change and the Kyoto Protocol.

Forests not only influence the global carbon cycle, but they also dramatically affect the water vapour and energy fluxes exchanged with the overlying atmosphere. It has been shown, for example, that the evapotranspiration of young plantations can be so great that the streamflow of neighbouring creeks is reduced by 50% (Jackson et al., 2005). Modelling studies on the impact of forest plantations in regions that are snow-covered in winter suggest that because of their reflectance (the so-called albedo), forest could increase regional temperature by up to four degrees (Betts, 2000; Bala et al., 2007; Davin et al., 2007; Zhao and Jackson, 2014). Management-related changes in the albedo, energy balance and water cycle of forests (Amiro et al., 2006a, b) are of the same magnitude as the differences between forests, grasslands and croplands (Luyssaert et al., 2014). Moreover, changes in the water vapour and the energy exchange may offset the cooling effect obtained by managing forests as stronger sinks for atmospheric CO₂ (Pielke et al., 2002). Despite the key implications of forest management on the carbon-energy-water exchange, there have been no integrated studies on the effects of forest management on the Earth's climate.

Earth system models are the most advanced tools for predicting future climate (Bonan, 2008). These models represent the interactions between the atmosphere and the surface beneath, with the surface formalised as a combination of open oceans, sea ice and land. For land, five classes are distinguished: glacier, lake, wetland, urban and vegetated. Vegetation is typically represented by different plant functional types. ORCHIDEE is the land-surface component of

the IPSL (Institut Pierre Simon Laplace) Earth system model. Hence, by design, the ORCHIDEE model can be run coupled to the LMDz global circulation model. In this coupled set-up, the atmospheric conditions affect the land surface and the land surface, in turn, affects the atmospheric conditions. Coupled land-atmosphere models thus offer the possibility to quantify both the climatic effects of changes in the land surface and the effects of climate change on the land surface. The most advanced land-surface models used, for instance, in Earth system models to predict climate changes (see the recent CMIP5 exercise), account for changes in vegetation cover but consider forests to be mature and ageless, e.g. JSBACH (Reick et al., 2013), CLM (Stöckli et al., 2008), MOSES (Cox et al., 1999), ORCHIDEE (Krinner et al., 2005) and LPJ-DVGM (Bonan et al., 2003). At present, none of the predictions of future climate thus accounts for the essential interactions between forest management and climate. This gap in modelling capability provides the motivation for further development of the ORCHIDEE land-surface model to realistically simulate both the biophysical and biogeochemical effects of forest management on the climate. The ORCHIDEE-CAN (short for ORCHIDEE-CANOPY) branch of the land-surface model was specifically developed to quantify the climatic effects of forest management.

The aim of this study is to describe the model developments and parametrisation within ORCHIDEE-CAN and to evaluate its performance. ORCHIDEE-CAN is validated against structural, biophysical and biogeochemical data on the European scale. To allow comparison with the standard version of ORCHIDEE, ORCHIDEE-CAN was run with a single-layer energy budget. A more detailed description and evaluation of the new multi-layer energy budget and multi-level radiative transfer scheme is given by Ryder et al. (2014), Chen et al. (2015) and McGrath et al. (2015b). A new forest management reconstruction, which is needed to drive forest management in ORCHIDEE-CAN, is presented in McGrath et al. (2015a), and the interactions between forest management and the new albedo scheme have been discussed by Otto et al. (2014).

# 2 Model overview

## 2.1 The starting point: ORCHIDEE SVN r2243

The land-surface model used for this study, ORCHIDEE, is based on two different modules (Krinner et al., 2005, their Fig. 2). The first module describes the fast processes such as the soil water budget and the exchanges of energy, water and CO₂ through photosynthesis between the atmosphere and the biosphere (Ducoudré et al., 1993; de Rosnay and Polcher, 1998). The second module simulates the carbon dynamics of the terrestrial biosphere and essentially represents processes such as maintenance and growth respiration, carbon allocation, litter decomposition, soil carbon dynamics and phenol-
ogy (Viovy and de Noblet-Ducoudré, 1997). The trunk version of ORCHIDEE describes global vegetation by 13 meta-classes (MTCs) with a specific parameter set (one for bare soil, eight for forests, two for grasslands and two for croplands). Each MTC can be divided into a user-defined number of plant functional types (PFTs) which can be characterised by at least one parameter value that differs from the parameter settings of the MTC. Parameters that are not given at the PFT level are assigned the default value for the MTC to which the PFT belongs. By default, none of the parameters is specified at the PFT level; hence, MTCs and PFTs are the same for the standard ORCHIDEE-trunk version. A concise description of the main processes in the ORCHIDEE-trunk version and a short motivation to change these modules in ORCHIDEE-CAN is given in Table 1.

Before running simulations, it is necessary to bring the soil carbon pools into equilibrium due to their slow fill rates, an approach known as model spin-up (Thornton and Rosenbloom, 2005; Xia et al., 2012). For a long time, spin-ups have been performed by brute force, i.e. running the model iteratively over a sufficiently long period which allows even the slowest carbon pool to reach equilibrium. This naïve approach is reliable but slow (in the case of ORCHIDEE it takes 3000 simulation years) and thus comes with a large computational demand, often exceeding the computational cost of the simulation itself. Alternative spin-up methods calling only parts of the model, e.g. subsequent cycles of 10 years of photosynthesis only followed by 100 year cycles of soil processes only, have been used for ORCHIDEE to reduce the computational cost in the past. These approaches, however, tend to lead to instabilities in litter and carbon pools. In recent years, semi-analytical methods have been proposed as a cost-effective solution to the spin-up issue (Martin et al., 2007; Lardy et al., 2011; Xia et al., 2012). A matrix-sequence method has been implemented in ORCHIDEE following the approach used by the PaSim model (Lardy et al., 2011). The semi-analytical spin-up implemented in ORCHIDEE relies on algebraic methods to solve a linear system of equations describing the seven carbon pools separately for each PFT. Convergence of the method and thus equilibrium of the carbon pools is assumed to be reached when the variation of the passive carbon pool (which is the slowest) drops below a pre-defined threshold. The net biome production (NBP) is used as a second diagnostic criterion to confirm equilibrium of the carbon pools. In order to optimise computing resources, the semi-analytical spin-up will stop before the end of the run once the convergence criteria are met. ORCHIDEE's implementation of the semi-analytical spin-up has been validated on regional and global scales against a naïve spin-up, and has been found to converge 12 to 20 times faster. The largest gains were realised in the tropics and the smallest gains in boreal climate (not shown).

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Schematic overview of the changes in ORCHIDEE-CAN. For the trunk the most important processes and connections are indicated in black, while the processes and connections that were added or changed in ORCHIDEE-CAN are indicated in red. Numbered arrows are discussed in Sect. 2.2.</figcaption></figure>

## 2.2 Modifications between ORCHIDEE SVN r2243 and ORCHIDEE-CAN SVN r2290

One major overarching change in the ORCHIDEE-CAN branch is the increase in internal consistency within the model by adding connections between the different processes (Fig. 1, red arrows). A more specific novelty is the introduction of circumference classes within forest PFTs, based on the work of Bellassen et al. (2010). For the temperate and boreal zone, tree height and crown diameter are calculated from allometric relationships of tree diameter that were parametrised based on the French, Spanish, Swedish and German forest inventory data and the observational data from Pretzsch (2009). The circumference classes thus allow calculation of the social position of trees within the canopy, which justifies applying an intra-tree competition rule (Deleuze et al., 2004) to account for the fact that trees with a dominant position in the canopy are more likely to intercept light than suppressed trees, and, therefore, contribute more to the stand level photosynthesis and biomass growth. To respect the competition rule of Deleuze et al. (2004), a new allocation scheme was developed based on the pipe model theory (Shinozaki et al., 1964) and its implementation by Sitch et al. (2003). The scheme allocates carbon to different biomass pools (leaves, fine roots, and sapwood) while
<table><caption>**Table 1.** Concise description of the modules in the standard ORCHIDEE version with the motivation to change the modules in ORCHIDEE-CAN.</caption><thead><tr><th>Module</th><th>Description</th><th>Motivation for change</th></tr></thead><tbody><tr><td>Albedo</td><td>For each PFT the total albedo for the grid square is computed as a weighted average of the vegetation albedo, the soil albedo, and the snow albedo.</td><td>The scheme overlooks the effect of vegetation shading bare soil for sparse canopies and gives the ground in all PFTs the same reflectance properties as bare soil.</td></tr><tr><td>Soil hydrology</td><td>Vertical water flow in the soil is based on the Fokker-Planck equation that resolves water diffusion in non-saturated conditions from the Richards equation (Richards, 1931). The 2 m soil column consists of 11 moisture layers with an exponentially increasing depth (D'Orgeval et al., 2008).</td><td>No change</td></tr><tr><td>Soil temperature</td><td>The soil temperature is computed according to the Fourier equation using a finite difference implicit scheme with seven numerical nodes unevenly distributed between 0 and 5.5 m (Hourdin, 1992).</td><td>No change</td></tr><tr><td>Energy budget</td><td>The coupled energy balance scheme, and its exchange with the atmosphere, is based on that of Dufresne and Ghattas (2009). The surface is described as a single layer that includes both the soil surface and any vegetation.</td><td>A big leaf approach does not account for within canopy transport of carbon, water and energy. Further, it is inconsistent with the current multi-layer photosynthesis approach and the new multi-layer albedo approach.</td></tr><tr><td>Photosynthesis</td><td>C3 and C4 photosynthesis is calculated following Farquhar et al. (1980) and Collatz et al. (1992), respectively. Photosynthesis assigns artificial LAI levels to calculate the carbon assimilation of the canopy. These levels allow for a saturation of photosynthesis with LAI, but have no physical meaning.</td><td>The scheme uses a simple Beer law transmission of light to each level, which is inconsistent with the new albedo scheme.</td></tr><tr><td>Autotrophic respiration</td><td>Autotrophic respiration distinguishes maintenance and growth respiration. Maintenance respiration occurs in living plant compartments and is a function of temperature, biomass and, the prescribed carbon/nitrogen ratio of each tissue (Ruimy et al., 1996). A prescribed fraction of 28 % of the photosynthates allocated to growth is used in growth respiration (McCree, 1974). The remaining assimilates are distributed among the various plant organs using an allocation scheme based on resource limitations (see allocation).</td><td>No change</td></tr><tr><td>Carbon allocation</td><td>Carbon is allocated to the plant following resource limitations Friedlingstein et al. (1999). Plants allocate carbon to their different tissues in response to external limitations of water, light and nitrogen availability. When the ratios of these limitations are out of bounds, prescribed allocation factors are used.</td><td>The resource limitation approach requires capping LAI at a predefined value. Due to this cap, the allocation rules are most often not applied, reducing the scheme to prescribing allocation.</td></tr><tr><td>Phenology</td><td>At the end of each day, the model checks whether the conditions for leaf onset are satisfied. The PFT-specific conditions are based on long- and short-term warmth and/or moisture conditions (Botta et al., 2000).</td><td>No change</td></tr><tr><td>Mortality and turnover</td><td>All biomass pools have a turnover time. Living biomass is transferred to the litter pool; litter is decomposed or transferred to the soil pool.</td><td>This approach is not capable of modelling stand dimensions.</td></tr><tr><td>Soil and litter carbon and heterotrophic respiration</td><td>Following (Parton et al., 1989), prescribed fractions of the different plant components go to the metabolic and structural litter pools following senescence, turnover or mortality. The decay of metabolic and structural litter is controlled by temperature and soil or litter humidity. For structural litter, its lignin content also influences the decay rate.</td><td>No change</td></tr><tr><td>Forest management</td><td>An explicit distribution of individual trees (Bellassen et al., 2010) is the basis for a process-based simulation of mortality. The aboveground stand-scale wood increment is distributed on a yearly time step among individual trees according to the rule of (Deleuze et al., 2004): the basal area of each individual tree grows proportionally to its circumference.</td><td>The concept of the original implementation were retained, however, the implementation was adjusted for consistency with the new allocation scheme and to have a larger diversity of management strategies.</td></tr></tbody></table>
respecting the differences in longevity and hydraulic conduc-
tivity between the pools. In addition to the biomass of the
different pools, leaf area index (LAI), crown volume, crown
density, stem diameter, stem height and stand density are cal-
culated and now depend on accumulated growth. The new
scheme allows for the removal of the parameter that caps the
maximum LAI (Table 1).

The calculation of tree dimensions (e.g. sapwood area and tree height) that respect the pipe theory supports making use of the hydraulic architecture of plants to calculate the plant water supply (Fig. 1, arrow 1), which is the amount of water a plant can transport from the soil to its stomata. The representation of the plant hydraulic architecture is based on the scheme of Hickler et al. (2006). The water supply is calculated as the ratio of the pressure difference between soil and leaves, and the total hydraulic resistance of the roots, leaves and sapwood, where the sapwood resistance is increased when cavitation occurs. Species-specific parameter values were compiled from the literature. As the scheme makes use of the soil water potential, it requires the use of the 11-layer hydrology scheme of de Rosnay (2002) (Table 1). When transpiration based on energy supply exceeds transpiration based on the water supply, the latter restricts stomatal conductance directly, which is a physiologically more realistic representation of drought stress than the reduction of the carboxylation capacity (Flexas et al., 2006) done in the standard version of ORCHIDEE (further also referred to as the “trunk” version). In line with this approach, the drought stress factor used to trigger phenology and senescence is now calculated as the ratio between the transpiration based on water supply and transpiration based on atmospheric demand (Fig. 1, arrow 2).

The new allocation scheme also drastically changed the way forests are represented in the ORCHIDEE-CAN branch. Although the exact location of the canopies in the stand is not known, individual tree canopies are now spherical elements with their horizontal location following a Poisson distribution across the stand. Each PFT contains a user-defined number of model trees, each one corresponding to a circumference class. Model trees are replicated to give realistic stand densities. Following tree growth, canopy dimensions and stand density are updated (Fig. 1, arrow 3). This formulation results in a dynamic canopy structure that is exploited in other parts of the model, i.e. precipitation interception, transpiration, energy budget calculations, a radiation scheme (Fig. 1, arrow 4) and absorbed light for photosynthesis (Fig. 1, arrow 5). In the trunk version these processes are driven by the big-leaf canopy assumption. The introduction of an explicit canopy structure is thought to be a key development with respect to the objectives of the ORCHIDEE-CAN branch, i.e. quantifying the biogeochemical and biophysical effects of forest management on atmospheric climate.

The radiation transfer scheme at the land surface benefits
from the introduction of canopy structure. The trunk version
of ORCHIDEE prescribes the vegetation albedo solely as a

function of LAI. In the ORCHIDEE-CAN branch each tree
canopy is assumed to be composed of uniformly distributed
single scatterers. Following the assumption of a Poisson dis-
tribution of the trees on the land surface, the model of Haverd
et al. (2012) calculates the transmission probability of light to
any given vertical point in the forest. This transmission prob-
ability is then used to calculate an effective LAI, which is a
statistical description of the vertical distribution of leaf mass
that accounts for stand density and horizontal tree distribu-
tion. The complexity and computational costs are largely re-
duced by using the effective LAI in combination with the 1-
D two-stream radiation transfer model of Pinty et al. (2006)
rather than resolving a full 3-D canopy model. By using the
effective LAI, the 1-D model reproduces the radiative fluxes
of the 3-D model. The approach of the two-stream radia-
tion transfer model was extended for a multi-layer canopy
(McGrath et al., 2015b) to be consistent with the multi-layer
energy budget and to better account for non-linearities in
the photosynthesis model. The scattering parameters and the
background albedo (i.e. the albedo of the surface below the
dominant tree canopy) for the two-stream radiation transfer
model were extracted from the Joint Research Centre Two-
stream Inversion Package (JRC-TIP) remote sensing product
(Sect. 4.7). This approach produces fluxes of the light ab-
sorbed, transmitted, and reflected by the canopy at vertically
discretised levels, which are then used for the energy bud-
get (Fig. 1, arrow 6) and photosynthesis calculations (Fig. 1,
arrow 5).

The canopy radiative transfer scheme of Pinty et al. (2006) separates the calculation of the fluxes resulting from downwelling direct and diffuse light, with different scattering parameters available for near-infrared (NIR) and visible (VIS) light sources. The snow albedo scheme in the trunk does not distinguish between these two short-wave bands. Therefore, the snow scheme of the Biosphere-Atmosphere Transfer Scheme (BATS) for the Community Climate Model (Dickinson et al., 1986) was incorporated into the ORCHIDEE-CAN branch, since it distinguishes between the NIR and VIS radiation. The radiation scheme of Pinty et al. (2006) requires snow to be put on the soil below the tree canopy instead of on the canopy itself. The calculation of the snow coverage of a PFT therefore had to be revised according to the scheme of Yang et al. (1997), which allows for snow to completely cover the ground at depths greater than 0.2 m. The parameter values of Yang et al. (1997) were used in the ORCHIDEE-CAN branch.

The ORCHIDEE-CAN branch differs from any other land-
surface model by the inclusion of a newly developed multi-
layer energy budget. There are now subcanopy wind, tem-
perature, humidity, long-wave radiation and aerodynamic re-
sistance profiles, in addition to a check of energy closure
at all levels. The energy budget represents an implementa-
tion of some of the characteristics of detailed single-site, it-
erative canopy models (e.g. Baldocchi, 1988; Ogee et al.,
2003) within a system that is coupled implicitly to the at-
mosphere. As an enhancement to the trunk version of ORCHIDEE (Table 1), the new approach also generates a leaf temperature, using a vegetation profile and a vertical short-wave and long-wave radiation distribution scheme (Ryder et al., 2014), which will be fully available when parametrisation of the scheme has been completed across test sites corresponding to the species within the model (Chen et al., 2015). As with the trunk version, the new energy budget is calculated implicitly (Polcher et al., 1998; Best et al., 2004). An implicit solution is a linear solution in which the surface temperature and fluxes are calculated in terms of the atmospheric input at the same time step, whereas an explicit solution uses atmospheric input from the previous time step to calculate the surface temperature and fluxes. Although it is less straightforward to derive, the implicit solution is more computationally efficient and stable, which allows the model to be run over a time step of 15 min when coupled to the LMDz atmospheric model – much longer than would be the case for an explicit model. Parameters were derived by optimising the model against the observations from short-term field campaigns. The new scheme may also be reduced to the existing single layer case, so as to provide a means of comparison and compatibility with the ORCHIDEE-trunk version.

The combined use of the new energy budget and the hydraulic architecture of plants required changes to the calculation of the stomatal conductance and photosynthesis (Fig. 1, arrow 7). When water supply limits transpiration, stomatal conductance is reduced and photosynthesis needs to be recalculated. Given that photosynthesis is among the computational bottlenecks of the model, the semi-analytical procedure as available in previous trunk versions (r2031 and further) is replaced by an adjusted implementation of the analytical photosynthesis scheme of Yin and Struik (2009), which is also implemented in the latest ORCHIDEE-trunk version. In addition to an analytical solution for photosynthesis, the scheme includes a modified Arrhenius function for the temperature dependence that accounts for a decrease in carboxylation capacity ($k_V_{cmax}$) and electron transport capacity ($k_{Jmax}$; see Table 2 for variable explanations) at high temperatures and a temperature-dependent $k_{Jmax}/V_{cmax}$ ratio (Kattge and Knorr, 2007). The temperature response of $k_V_{cmax}$ and $k_{Jmax}$ was parametrised with values from reanalysed data in the literature (Kattge and Knorr, 2007), whereas $k_V_{cmax}$ and $k_{Jmax}$ at a reference temperature of 25 °C were derived from observed species-specific values in the TRY database (Kattge et al., 2011). As the amount of absorbed light varies with height (or canopy depth), the absorbed light computed from the albedo routines is now directly used in the photosynthesis scheme, resulting in full consistency between the top of the canopy albedo and absorption. This new approach replaces the old scheme which used multiple levels based on the leaf area index, not the physical height.

ORCHIDEE-CAN incorporates a systematic mass balance closure for carbon cycling to ensure that carbon is not getting

created or destroyed during the simulation. Hence, budget closure is now consistently checked for water, carbon and energy throughout the model.

The trunk uses 13 PFTs to represent vegetation globally: one PFT for bare soil, eight for forests, two for grasslands, and two for croplands. The ORCHIDEE-CAN branch makes use of the externalisation of the PFT-dependent parameters by adding 12 parameter sets that represent the main European tree species. Species parameters were extracted from a wide range of sources including original observations, large databases, primary research and remote sensing products (Sect. 4). The use of age classes is introduced through externalisation of the PFT parameters as well. Age classes are used during land cover change and forest management to simulate the regrowth of a forest. Following a land cover change, biomass and soil carbon pools (but not soil water columns) are either merged or split to represent the various outcomes of a land cover change. The number of age classes is user defined. Contrary to typical age classes, the boundaries are determined by the tree diameter rather than the age of the trees.

Finally, the forest management strategies in the ORCHIDEE-CAN branch were refined from the original forest management (FM) branch (Bellassen et al., 2010). Self-thinning was activated for all forests regardless of human management, contrary to the original FM branch. The new default management strategy thus has no human intervention but includes self-thinning, which replaces the fixed 40 year turnover time for woody biomass. Three management strategies with human intervention have been implemented: (1) “high stands”, in which human intervention is restricted to thinning operations based on stand density and diameter, with occasional clear-cuts. Aboveground stems are harvested during operations, while branches and belowground biomass are left to litter; (2) “coppices” involve two kinds of cuts. The first coppice cut is based on stem diameter and the aboveground woody biomass is harvested, whereas the belowground biomass is left living. From this belowground biomass, new shoots sprout, which increases the number of aboveground stems. In subsequent cuts the number of shoots is not increased, although all aboveground wood biomass is still harvested; and (3) “short rotation coppices”, where rotation periods are based on age and are generally very short (3–6 years). The different management strategies can occur with or without litter raking, which reduces the litter pools and has a long-term effect on soil carbon (Gimmi et al., 2012). All management types are parametrised based on forest inventory data, yield tables and guidelines for forest management. The inclusion of forest management resulted in two additional carbon pools, branches and coarse roots (i.e. aboveground and belowground woody biomass) and therefore required an extension to the semi-analytical spin-up method (Sect. 2.1). The semi-analytical spin-up is now run for nine C pools.
<table><caption>**Table 2.** Variable description. Variables were grouped as follows: *F* = flux, *f* = fraction, *M* = pool, *m* = modulator, *d* = stand dimension, *T* = temperature, *p* = pressure, *R* = resistance, *q* = humidity, *g* = function.</caption><thead><tr><th>Symbol in text</th><th>Unit</th><th>Symbol in ORCHIDEE-CAN</th><th>Description</th></tr></thead><tbody><tr><td>F<sub>rm</sub></td><td>gC m<sup>-2</sup> s<sup>-1</sup></td><td>resp_maint</td><td>Maintenance respiration</td></tr><tr><td>F<sub>rg</sub></td><td>gC m<sup>-2</sup> s<sup>-1</sup></td><td>resp_growth</td><td>Growth respiration</td></tr><tr><td>F<sub>LW,i</sub></td><td>W m<sup>2</sup></td><td>r_lw</td><td>Long-wave radiation incident at vegetation level i</td></tr><tr><td>F<sub>SW,i</sub></td><td>W m<sup>2</sup></td><td>r_sw</td><td>Short-wave radiation incident at vegetation level i</td></tr><tr><td>F<sub>Trs</sub></td><td>m s<sup>-1</sup></td><td>Transpirsupply</td><td>Amount of water that a tree can get up from the soil to its leaves for transpiration</td></tr><tr><td>T<sub>a,i</sub></td><td>K</td><td>temp_atmos_pres</td><td>Atmospheric temperature at the "present" and "next" time step, respectively, at level i</td></tr><tr><td>T<sub>L,i</sub></td><td>K</td><td>temp_atmos_next</td><td>Leaf temperature at level i</td></tr><tr><td>q<sub>a,i</sub></td><td>kg kg<sup>-1</sup></td><td>temp_leaf_pres</td><td>Leaf-specific humidity at level i</td></tr><tr><td>q<sub>L,i</sub></td><td>kg kg<sup>-1</sup></td><td>temp_leaf_pres</td><td>Leaf-specific humidity at level i</td></tr><tr><td>M<sub>l</sub></td><td>gC plant<sup>-1</sup></td><td>Cl</td><td>Leaf mass of an individual plant</td></tr><tr><td>M<sub>s</sub></td><td>gC plant<sup>-1</sup></td><td>Cs</td><td>Sapwood mass of an individual plant</td></tr><tr><td>M<sub>h</sub></td><td>gC plant<sup>-1</sup></td><td>Ch</td><td>Heartwood mass of an individual plant</td></tr><tr><td>M<sub>r</sub></td><td>gC plant<sup>-1</sup></td><td>Cr</td><td>Root mass of an individual plant</td></tr><tr><td>M<sub>linc</sub></td><td>gC plant<sup>-1</sup></td><td>Cl_inc</td><td>Increment in leaf mass of an individual plant</td></tr><tr><td>M<sub>sinc</sub></td><td>gC plant<sup>-1</sup></td><td>Cs_inc</td><td>Increment in sapwood mass of an individual plant</td></tr><tr><td>M<sub>rinc</sub></td><td>gC plant<sup>-1</sup></td><td>Cr_inc</td><td>Increment in root mass of an individual plant</td></tr><tr><td>M<sub>totinc</sub></td><td>gC</td><td>b_inc_tot</td><td>Total biomass increment</td></tr><tr><td>M<sub>inc</sub></td><td>gC plant<sup>-1</sup></td><td>b_inc</td><td>Increment in plant biomass of an individual plant</td></tr><tr><td>M<sub>swc</sub></td><td>m<sup>3</sup> m<sup>-3</sup></td><td>swc</td><td>Volumetric soil water content</td></tr><tr><td>m<sub>w</sub></td><td>—</td><td>wstress_fac</td><td>Modulator for water stress as experienced by the plants</td></tr><tr><td>m<sub>ψ</sub></td><td>MPa</td><td>psi_soil_tune</td><td>Modulator to account for resistance in the soil-root interface</td></tr><tr><td>m<sub>Ndeath</sub></td><td>—</td><td>scale_factor</td><td>Normalisation factor for mortality</td></tr><tr><td>m<sub>LAIcorr</sub></td><td>—</td><td>laiCorrection_factor</td><td>Adjustable parameter in the calculation of gap probabilities of grasses and crops</td></tr><tr><td>d<sub>h</sub></td><td>m</td><td>height</td><td>Plant height</td></tr><tr><td>d<sub>l</sub></td><td>m<sup>-2</sup></td><td>—</td><td>One-sided leaf area of an individual plant</td></tr><tr><td>d<sub>s</sub></td><td>m<sup>-2</sup></td><td>—</td><td>Sapwood area of an individual plant</td></tr><tr><td>d<sub>hinc</sub></td><td>m</td><td>delta_height</td><td>Height increment</td></tr><tr><td>d<sub>dhb</sub></td><td>m</td><td>dia</td><td>Plant diameter</td></tr><tr><td>d<sub>ba</sub></td><td>m<sup>2</sup> plant<sup>-1</sup></td><td>ba</td><td>Basal area</td></tr><tr><td>d<sub>bainc</sub></td><td>m<sup>2</sup> plant<sup>-1</sup></td><td>delta_ba</td><td>Basal area increment</td></tr><tr><td>d<sub>circ</sub></td><td>m</td><td>circ</td><td>Stem circumference of an individual plant</td></tr><tr><td>d<sub>ind</sub></td><td>trees</td><td>n_circ_class</td><td>Number of trees in diameter class l</td></tr><tr><td>d<sub>c</sub></td><td>m<sup>2</sup></td><td>crown_shadow_h</td><td>Projected area of an opaque tree crown</td></tr><tr><td>d<sub>csa</sub></td><td>m<sup>2</sup></td><td>csa_sap</td><td>Projected crown surface area</td></tr><tr><td>d<sub>LAI</sub></td><td>m<sup>2</sup> leaf m<sup>-2</sup> ground</td><td>—</td><td>Leaf area index</td></tr><tr><td>d<sub>LAl efficiency</sub></td><td>—</td><td>laieff</td><td>Effective leaf area index</td></tr><tr><td>d<sub>LAl above</sub></td><td>—</td><td>lai_sum</td><td>Sum of the LAI of all levels above the current level</td></tr><tr><td>d<sub>A,i</sub></td><td>m<sup>2</sup></td><td>—</td><td>Cross-sectional area of vegetation level i</td></tr><tr><td>d<sub>hl,i</sub></td><td>m</td><td>delta_h</td><td>Vegetation height of level i</td></tr><tr><td>d<sub>V,i</sub></td><td>m<sup>3</sup></td><td>—</td><td>Volume of vegetation level i</td></tr><tr><td>d<sub>rd</sub></td><td>—</td><td>root_dens</td><td>Root density</td></tr><tr><td>d<sub>λ</sub></td><td>ind m<sup>2</sup></td><td>—</td><td>Inverse of the individual plant density</td></tr><tr><td>p<sub>delta</sub></td><td>MPa</td><td>delta_P</td><td>Pressure difference between leaves and soil</td></tr><tr><td>p<sub>ψsr</sub></td><td>MPa</td><td>psi_soilroot</td><td>Bulk soil water potential in the rooting zone</td></tr><tr><td>p<sub>ψs</sub></td><td>MPa</td><td>psi_soil</td><td>Soil water potential for each soil layer</td></tr><tr><td>R<sub>r</sub></td><td>MPa s m<sup>-3</sup></td><td>R_root</td><td>Hydraulic resistance of roots</td></tr><tr><td>R<sub>sap</sub></td><td>MPa s m<sup>-3</sup></td><td>R_sap</td><td>Hydraulic resistance of sapwood</td></tr><tr><td>R<sub>l</sub></td><td>MPa s m<sup>-3</sup></td><td>R_leaf</td><td>Hydraulic resistance of leaves</td></tr><tr><td>R<sub>temp</sub></td><td>MPa s m<sup>-3</sup></td><td>—</td><td>Hydraulic resistance of roots, sapwood or leaves adjusted for temperature</td></tr><tr><td>R<sub>a,i</sub></td><td>s m<sup>-1</sup></td><td>big_r</td><td>Aerodynamic resistance of vegetation at level i in the canopy</td></tr><tr><td>R<sub>s,i</sub></td><td>s m<sup>-1</sup></td><td>big_r_prime</td><td>Sum of the stomatal and leaf boundary layer resistance terms for latent heat</td></tr></tbody></table>
<table><caption>**Table 2. Continued.**</caption><thead><tr><th>Symbol in text</th><th>Unit</th><th>Symbol in ORCHIDEE-CAN</th><th>Description</th></tr></thead><tbody><tr><td>f<sub>Pwc</sub></td><td>–</td><td>Pwc_h</td><td>Porosity of a tree crown</td></tr><tr><td>f<sub>trees</sub></td><td>–</td><td>PgapL</td><td>Gap probability for trees</td></tr><tr><td>f<sup>gC</sup><sub>Pgap</sub></td><td>–</td><td>PgapL</td><td>Gap probability for grasses and crops</td></tr><tr><td>f<sup>bs</sup><sub>Pgap</sub></td><td>–</td><td>PgapL</td><td>Gap probability for bare soil</td></tr><tr><td>f<sup>icir</sup><sub>death</sub></td><td>–</td><td>mortality</td><td>Mortality fraction per circumference class</td></tr><tr><td>f<sub>KF</sub></td><td>–</td><td>KF</td><td>Leaf allocation factor</td></tr><tr><td>f<sub>LF</sub></td><td>–</td><td>LF</td><td>Root allocation factor</td></tr><tr><td>f<sub>γ</sub></td><td>–</td><td>gamma</td><td>Slope of the intra-specific competition</td></tr><tr><td>f<sub>s</sub></td><td>m</td><td>s</td><td>Slope of linearised relationship between height and basal area</td></tr><tr><td>f<sub>rl</sub></td><td>–</td><td>leaf_reflectance</td><td>Reflectance of a single leaf</td></tr><tr><td>f<sub>ul</sub></td><td>–</td><td>leaf_transmittance</td><td>Transmittance of a single leaf</td></tr><tr><td>f<sub>Rbgd</sub></td><td>–</td><td>bdg_reflectance</td><td>Reflectance of the ground beneath the canopy</td></tr><tr><td>f<sup>Rk</sup><sub>Coll,veg</sub></td><td>–</td><td>Collim_alb_BB, Isotrop_alb_BB</td><td>Reflected fraction of light to the atmosphere which has collided with canopy elements, separated for direct and diffuse sources, respectively</td></tr><tr><td>f<sup>Rk</sup><sub>UnColl, bgd</sub></td><td>–</td><td>Collim_alb_BC, Isotrop_alb_BC</td><td>Reflected fraction of light to the atmosphere which has not collided with any canopy elements, separated for direct and diffuse sources, respectively</td></tr><tr><td>f<sup>T</sup><sub>UnColl, veg</sub></td><td>–</td><td>Collim_Tran_Uncoll</td><td>Transmitted fraction of light to the ground which has not collided with any canopy elements</td></tr><tr><td>f<sup>Rk</sup><sub>Coll, bgd,1</sub></td><td>–</td><td>–</td><td>Reflected fraction of light which has struck the background a single time and has collided with vegetation</td></tr><tr><td>f<sup>Rk</sup><sub>Coll, bgd,n</sub></td><td>–</td><td>–</td><td>Reflected fraction of light which has struck the background multiple times and has collided with vegetation</td></tr><tr><td>z</td><td>m</td><td>z_array</td><td>Height above the soil</td></tr><tr><td>θ<sub>z</sub></td><td>radians</td><td>solar_angle</td><td>Solar zenith angle</td></tr><tr><td>θ<sub>μ</sub></td><td>radians</td><td>–</td><td>Cosine of the solar zenith angle</td></tr><tr><td>g<sub>G</sub></td><td>–</td><td>–</td><td>Leaf orientation function</td></tr><tr><td>g<sub>σ</sub></td><td>–</td><td>sigmas</td><td>Cut-off circumference of the intra-specific competition, calculated as a function of k<sub>ncirc</sub></td></tr></tbody></table>

# 3 Description of the developments

## 3.1 Allocation

Following bud burst, photosynthesis produces carbon that is added to the labile carbon pool. Labile carbon is used to sustain the maintenance respiration flux ($F_{rm}$), which is the carbon cost to keep existing tissue alive (Amthor, 1984). Maintenance respiration for the whole plant is calculated by summing maintenance respiration of the different plant compartments, which is a function of the nitrogen concentration of the tissue following the Beer–Lambert law and subtracted from the whole-plant labile pool (up to a maximum of 80% of the labile pool).

The remaining labile carbon pool is split into an active and a non-active pool. The size of the active pool is calculated as a function of plant phenology and temperature and was formalised following Ryan (1991), Sitch et al. (2003) and Zaehle and Friend (2010). The remaining non-active pool is used to restore the labile and carbohydrate reserve pools according to the rules proposed in Zaehle and Friend (2010). The labile pool is limited to 1% of the plant biomass or 10 times the actual daily photosynthesis. Any excess carbon is transferred to the non-respiring carbohydrate reserve pool. The carbohydrate reserve pool is capped to reflect limited

starch accumulation in plants, but carbon can move freely between the two reserve pools. After accounting for growth respiration ($F_{rg}$), i.e. the cost for producing new tissue excluding the carbon required to build the tissue itself (Amthor, 1984), the total allocatable C used for plant growth is obtained ($M_{totinc}$).

New biomass is allocated to leaves, roots, sapwood, heartwood, and fruits. Allocation to leaves, roots and wood respects the pipe model theory (Shinozaki et al., 1964) and thus assumes that producing one unit of leaf mass requires a proportional amount of sapwood to transport water from the roots to the leaves as well as a proportional fraction of roots to take up the water from the soil. The different biomass pools have different turnover times, and therefore at the end of the daily time step, the actual biomass components may no longer respect the allometric relationships. Consequently, at the start of the time step carbon is first allocated to restore the allometric relationships before the remaining carbon is allocated in the manner described below. The scaling parameter between leaf and sapwood mass is derived from:

$$d_1 = k_{ls} \times m_w \times d_s \quad (1)$$

where $d_1$ is the one-sided leaf area of an individual plant, $d_s$ is the sapwood cross-section area of an individual plant, $k_{ls}$ is a parameter linking leaf area to sapwood cross-section area,
and $m_w$ is the water stress as defined in Sect. 3.2. Alternatively, leaf area can be written as a function of leaf mass ($M_l$) and the specific leaf area ($k_{sla}$):

$$d_l = M_l \times k_{sla}. \qquad (2)$$

Sapwood mass $M_s$ can be calculated from the sapwood cross-section area $d_s$ as follows:

$$M_s = d_s \times d_h \times k_{\rho s}, \qquad (3)$$

where $d_h$ is the tree height and $k_{\rho s}$ is the sapwood density. Following substitution of Eqs. (2) and (3) into Eq. (1), leaf mass can be written as a function of sapwood mass:

$$M_l = (M_s \times f_{KF}) / d_h, \qquad (4)$$

where,

$$f_{KF} = (k_{ls} \times m_w) / (k_{sla} \times k_{\rho s}), \qquad (5)$$

where $k_{ls}$ is calculated as a function of the gap fraction as supported by site-level observations (Simonin et al., 2006):

$$k_{ls} = k_{lsmin} + f_{pgap, trees} \times (k_{lsmax} - k_{lsmin}). \qquad (6)$$

$k_{lsmin}$ is the minimum observed leaf area to sapwood area ratio, $k_{lsmax}$ is the maximum observed leaf area to sapwood area ratio and $f_{pgap, trees}$ is the actual gap fraction. By using the gap fraction as a control of $k_{ls}$ more carbon will be allocated to the leaves until canopy closure is reached.

Following Magnani et al. (2000), sapwood mass and root mass ($M_r$) are related as follows:

$$M_s = k_{sar} \times d_h \times M_r, \qquad (7)$$

where the parameter $k_{sar}$ is calculated according to Magnani et al. (2000) (their Eq. 17):

$$k_{sar} = \sqrt{(k_{rcon}/k_{scon}) \times (k_{\tau s}/k_{\tau r}) \times k_{\rho s}}, \qquad (8)$$

where $k_{rcon}$ is the hydraulic conductivity of roots, $k_{scon}$ is the hydraulic conductivity of sapwood, $k_{\tau s}$ is the longevity of sapwood and $k_{\tau r}$ is the root longevity. Following substitution of Eq. (4) into Eq. (7) and some rearrangement, leaf mass can be written as a function of root mass:

$$M_l = f_{LF} \times M_r, \qquad (9)$$

where,

$$f_{LF} = k_{sar} \times f_{KF}. \qquad (10)$$

Parameter values used in Eqs. (1) to (9), i.e. $k_{lsmax}$, $k_{lsmin}$, $k_{sar}$, $k_{sla}$, $k_{\rho s}$, $k_{rcon}$, $k_{scon}$, $k_{\tau s}$ and $k_{\tau r}$, are based on literature review (Tables S1, S2 and S3 in the Supplement). The allo-metric relationships between the plant components and the hydraulic architecture of the plant (Sect. 3.2) are both based on the pipe model theory; hence, both the allocation and the

hydraulic architecture module use the same parameter values for root and sapwood conductivity.

In this version of ORCHIDEE, forests are modelled to have $k_{ncirc}$ circumference classes with $d_{ind}$ identical trees in each one. Hence, the allocatable biomass ($M_{totinc}$) needs to be distributed across $l$ diameter classes:

$$M_{totinc} = \sum_{l} (l)[d_{ind(l)} \times M_{inc(l)}], \qquad (11)$$

where $M_{inc(l)}$ is the biomass that can be allocated to diameter class $l$. Mass conservation thus requires:

$$M_{inc(l)} = M_{linc(l)} + M_{rinc(l)} + M_{sinc(l)}, \qquad (12)$$

where $M_{linc(l)}$, $M_{rinc(l)}$ and, $M_{sinc(l)}$ are the increase in leaf, root and wood biomass for a tree in diameter class $l$, respectively. Equations (4) and (9) can be rewritten as

$$(M_{l(l)} + M_{linc(l)}) / (M_{s(l)} + M_{sinc(l)}) = f_{KF} / (d_{h(l)} + d_{hinc(l)}) \qquad (13)$$

$$(M_{l(l)} + M_{linc(l)}) = (M_{r(l)} + M_{rinc(l)}) \times f_{LF} \qquad (14)$$

An allometric relationship is used to describe the relationship between tree height and basal area (Pretzsch, 2009):

$$d_{h(l)} = k_{\alpha 1} \times (4/\pi \times d_{ba(l)})^{(k_{\beta 1}/2)}. \qquad (15)$$

The change in height is then calculated as

$$d_{hinc(l)} = [k_{\alpha 1} \times (4/\pi \times (d_{ba(l)} + d_{bainc(l)}))^{(k_{\beta 1}/2)}] - d_{h(l)}, \qquad (16)$$

where $d_{ba(l)}$ and $d_{bainc(l)}$ are the basal area and its increment, respectively. $k_{\alpha 1}$ and $k_{\beta 1}$ are allometic constants relating tree diameter and height. The distribution of C across the $l$ diameter classes depends on the basal area of the model tree within each diameter class. Trees with a large basal area are assigned more carbon for wood allocation than trees with a small basal area, according to the method of Deleuze et al. (2004).

$$d_{bainc(l)} = f_{\gamma} \times (d_{circ(l)} - k_{m} \times g_{\sigma} + \sqrt{(k_{m} \times g_{\sigma} + d_{circ(l)})^{2} - (4 \times g_{\sigma} \times d_{circ(l)})}) / 2, \qquad (17)$$

where $k_m$ is a parameter, $f_{\gamma}$ and $g_{\sigma}$ are calculated from parameters and $d_{circ(l)}$ is the circumference of the model tree in diameter class $l$. $g_{\sigma}$ is a function of the diameter distribution of the stand at a given time step.

Equations (10) to (16) need to be simultaneously solved. An iterative scheme was avoided by linearising Eq. (15), which was found to be an acceptable numerical approximation as allocation is calculated at a daily time step, and hence the changes in height are small and the relationship is locally linear:

$$d_{hinc(l)} = d_{bainc(l)} / f_s, \qquad (18)$$
where $f_s$ is the slope of the locally linearised Eq. (15) and is calculated as

$$f_s = k_{\text{step}} / (k_{\alpha 1} \times (4/\pi \cdot (d_{\text{ba}} + k_{\text{step}}))^{(k_{\beta 1}/2)} - k_{\alpha 1} \times (4/\pi \times d_{\text{ba}})^{(k_{\beta 1}/2)}). \quad (19)$$

Equations (10), (11), (12), (13), (14), (16), (17) and (18) are then solved for $f_γ$. $f_γ$ distributes photosynthates across the different diameter classes and as such controls the intra-species competition within a stand. $f_γ$ thus depends on the total allocatable carbon and needs to be optimised at every time step. Once $f_γ$ has been calculated, $M_{\text{linc}(l)}$, $M_{\text{rinc}(l)}$ and $M_{\text{sinc}(l)}$ can be calculated.

## 3.2 Hydraulic architecture

The representation of the impact of soil moisture stress on water, carbon and energy fluxes has been identified as one of the major uncertainties in land-surface models (De Kauwe et al., 2013). Neither the empirical functions nor the soil moisture stress functions, which are commonly used in land-surface models, fully capture stomatal closure and limitation of C uptake during drought stress (Bonan et al., 2014; Verhoef and Egea, 2014). Therefore, we replaced the soil moisture stress function which limits C assimilation through a constraint on $k_{\text{Vcmax}}$ in the ORCHIDEE trunk, with a constraint based on the amount of water plants can transport from the soil to their leaves.

The model calculates plant water supply according to the implementation of hydraulic architecture by Hickler et al. (2006). Plant water supply is the amount of water the plant can transport from the soil to its stomata, accounting for the resistances to water transport in the roots, sapwood and leaves. If transpiration rate exceeds plant water supply, the stomatal conductance is reduced until equilibrium is reached.

The water flow from the soil to the leaves is driven by a gradient of decreasing water potential. Using Darcy's law (Slatyer, 1967; Whitehead, 1998), the supply of water for transpiration through stomata can be described as

$$F_{\text{Trs}} = p_{\text{delta}} / (R_r + R_{\text{sap}} + R_l), \quad (20)$$

where $p_{\text{delta}}$ is the pressure difference between the soil and the leaves; and $R_r$, $R_{\text{sap}}$ and $R_l$ are the hydraulic resistances of fine roots, sapwood and leaves, respectively. $p_{\text{delta}}$ is calculated following Whitehead (1998):

$$p_{\text{delta}} = p_{\psi \text{sr}} - k_{\psi 1} - (d_{\text{h}} \times k_{\psi \text{w}} \times k_{\text{g}}) \quad (21)$$

where $k_{\psi 1}$ is a PFT-specific minimal leaf water potential, which means that plants are assumed to maximise water uptake by lowering their $k_{\psi 1}$ to the minimum, if transpiration exceeds $F_{\text{Trs}}$ (Tyree and Sperry, 1989). The product of $d_h$, $k_{\psi \text{w}}$ and $k_g$ accounts for the loss in water potential by lifting a mass of water from the soil to the place of transpiration at height $d_h$, $k_{\psi \text{w}}$ is the density of water, and $k_g$ is the gravitational constant. The soil water potential in the rooting zone

$(p_{\psi \text{sr}})$ was calculated by adding a modulator $(m_{\psi})$ to the bulk soil water potential, which was calculated as the sum of the soil water potential in each soil layer weighted by the relative share of roots ($d_{\text{rd}}$) in the individual soil layer:

$$p_{\psi \text{sr}} = \sum_{l} [p_{\psi s} \times d_{\text{rd}}] + m_{\psi}. \quad (22)$$

The soil water potential for each layer $p_{\psi \text{sl}}$ is calculated from soil water content according to Van Genuchten (1980).

$$p_{\psi s(l)} = \frac{1}{k_{\text{av}}} \left( \left( \frac{M_{\text{swc}} - k_{\text{swcr}}}{k_{\text{swcs}} - k_{\text{swcr}}} \right)^{-1/k_{\text{mv}}} - 1 \right)^{1/k_{\text{nv}}}, \quad (23)$$

where $M_{\text{swc}}$ is the volumetric soil water content, $k_{\text{swcr}}$ and $k_{\text{swcs}}$ are respectively the residual and saturated soil water content and $k_{\text{av}}$, $k_{\text{mv}}$ and $k_{\text{nv}}$ are parameters.

Root resistance is related to the root mass and thus can be expressed as (Weatherley, 1982):

$$R_r = \frac{1}{(k_{\text{rcon}} \times M_r)}, \quad (24)$$

where $k_{\text{rcon}}$ is the fine root hydraulic conductivity per unit biomass. Sapwood resistance is calculated according to Magnani et al. (2000):

$$R_{\text{sap}} = \frac{d_{\text{h}}}{(d_{\text{s}} \times k_{\text{scon}})}, \quad (25)$$

where $k_{\text{scon}}$ is the sapwood-specific conductivity, which is decreased when cavitation occurs. The loss of conductance as a result of cavitation is a function of $p_{\psi \text{sr}}$ and was implemented by using an s-shaped vulnerability curve

$$k_{\text{scon}} = k_{\text{scon}} \times e^{(-p_{\psi \text{rs}}/k_{\psi 50})^{k_c}}, \quad (26)$$

where $k_{\psi 50}$ is the $p_{\psi \text{sr}}$ that causes 50% loss of conductance; and $k_c$ is a shape parameter.

$R_l$ is related to the specific leaf conductivity per unit leaf area ($k_l$) and the leaf area index:

$$R_l = \frac{1}{(k_{\text{lcon}} \times d_{\text{lAv}})}, \quad (27)$$

The response of water viscosity to low temperatures increases the resistance (Cochard et al., 2000). The relationship is described as

$$R_{\text{temp}} = \frac{R}{(k_{\alpha 1v} + k_{\alpha 2v} \times T)}, \quad (28)$$

where $k_{\alpha 1v}$ and $k_{\alpha 2v}$ are empirical parameters (Cochard et al., 2000), $R_{\text{temp}}$ is the temperature adjusted $R_l$, $R_{\text{sap}}$ or $R_r$, $T$ is air temperature for $R_l$ and $R_{\text{sap}}$ and $T$ is soil temperature for $R_r$.

If, for any time step, the transpiration calculated by the energy budget exceeds the amount of water the plant can transport from the soil to its stomata, transpiration is limited to
the plant water supply. As the transpiration is now reduced, the initial calculations of the energy budget and photosynthesis, solely based on atmospheric information, are no longer valid. As a result the energy budget and photosynthesis must be recalculated for the time step in question. For this recalculation, stomatal conductance at the canopy level is calculated such that transpiration equals the amount of water the plant can transport. Owing to the feedback between stomatal conductance, leaf surface temperature and transpiration, this calculation may require up to 10 iterations to converge, using a stationary iterative method. When the multi-layer energy budget is reduced to its single-layer implementation, however, canopy level stomatal conductance is decomposed to obtain the stomatal conductance at each canopy layer assuming that each layer is equally restricted by drought stress. Finally, the restricted stomatal conductance is used to calculate CO₂ assimilation rate according to the photosynthesis model by Farquhar, von Caemmerer and Berry (Sect. 3.6).

## 3.3 Canopy structure

Stand structure controls the amount of light that penetrates to a given depth in the canopy. For example, the amount of light reaching the forest floor will be higher for a stand with few mature trees compared to many young trees, even if both stands have the same leaf area index. Where a big leaf approach assumes a homogeneous block-shaped canopy (as in the trunk version of ORCHIDEE) and can therefore rely on the law of Beer and Lambert, a geometric approach is required to calculate light penetration through structured canopies. Light penetration needs to be simulated to calculate albedo (Sect. 3.4), photosynthesis (Sect. 3.6), partitioning of energy fluxes (Sect. 3.5) and the amount of light reaching the forest floor (see for example Sect. 3.1). The gap fraction, which is the basic information in calculating light penetration at different depths in the canopy, is calculated following the approach presented by Haverd et al. (2012) and formalised in their semi-analytical model. Rather than a spatially explicit approach, Haverd et al. (2012) follow a statistical approach which reduces the memory requirements for the simulations and limits the space requirements for storing the model output files.

The model of Haverd et al. (2012) represents the canopy by a statistical height distribution with varying crown sizes and stem diameters for each height class. The crown canopies are treated as spheroids containing homogeneously distributed single scatterers. Although this $f_{Pgap}$ model can explicitly include trunks, we made the decision to exclude them, as the spectral parameters for our radiation model (Sect. 3.4) are extracted from remote sensing data (Sect. 4.8) without distinguishing between leafy and woody masses. This gives the gap probability for trees as a function of height (z) and solar zenith angle ($\theta_z$):

$$f_{Pgap}^{\text{trees}}(\theta_z, z) = e^{-d_{\lambda} \times d_c(\theta_z, z) \times (1 - f_{Pwc}(\theta_z, z))}, \quad (29)$$

where $d_{\lambda}$ is the inverse of the tree density, $d_c$ is the projected crown area (for an opaque canopy), and $f_{Pwc}$ is the mean crown porosity. The overbar depicts the mean over the tree distribution as a function of tree height or, in our case, the mean over the l circumference classes. Following minor adaptations, the implementation of Haverd and Lovell (Haverd et al., 2012) was incorporated into ORCHIDEE-CAN. As there also exist crops, grasses, and bare soil in the model, $f_{Pgap}$ was adjusted for these situations as well. For grasses and crops, the same formulation is used:

$$f_{Pgap}^{\text{gc}}(\theta_z, z) = e^{-0.5 \times d_{\text{LAlabove}} \times m_{\text{LAlcorr}} / \cos(\theta_z)}, \quad (30)$$

where $d_{\text{LAlabove}}$ is the total amount of LAI above height z, and $m_{\text{LAlcorr}}$ is a correction factor to account for the fact that grasses and crops are treated as homogeneous blocks of vegetation with no internal structure, and is often referred to as a clumping factor. Here it is treated as a tunable parameter and therefore the term “correction factor” was used. For bare soil, there is no vegetation to intercept radiation, and therefore $f_{Pgap}^{\text{bs}}(\theta_z, z)$ is always unity.

## 3.4 Multi-layer two-way radiation scheme for tall canopies

Species-specific radiation absorbance, reflectance and transmittance by the forest canopy were calculated from a radiation transfer model (Pinty et al., 2006) which was parametrised by satellite-derived species-specific scattering values (Sect. 4.8). Given the complexity of radiation transfer, it remains challenging to accurately simulate radiation transfer through structurally and optically complex vegetation canopies without using explicit 3-D models. The applied 1-D model belongs to the family of two-stream models (Meador and Weaver, 1980) and thus calculates transmittance, absorbance and reflectance of both the incoming and outgoing radiation. The calculation of the reflectance at the top of the canopy due to a collimated source (i.e. the Sun) is divided into three components:

1. scattering of radiation between the vegetated elements with a black background

$$f_{\text{Coll, veg}}^{f_R} = f(\theta_{\text{mu}}, f_{\text{rl}}, f_{\text{tl}}, g_G, d_{\text{LAleff}}) \quad (31)$$

2. scattering of radiation by the background with a black canopy

$$f_{\text{UnColl, bgd}}^{f_R} = f_{\text{Rbgd}} \times e^{(-d_{\text{LAleff}}/(2 \times \theta_{\text{mu}}))} \times f_{\text{UnColl, veg}}^T \quad (32)$$

3. multiple scattering of radiation between the canopy and the background

$$f_{\text{Coll, bgd}}^{f_R} = f_{\text{Rbgd}} \times \left[ f_{\text{Coll, bgd},1}^{f_R} + f_{\text{Coll, bgd},n}^{f_R} \right] \quad (33)$$
Term (1) is widely used in cloud reflectance calculations, and depends on the cosine of the solar zenith angle ($\theta_{mu}$), the reflectance and transmittance of the single leaves ($f_{rl}$ and $f_{l}$, respectively), the leaf orientation function ($g_{G}$), and the effective LAI ($d_{LAleff}$). The exact definition of this term is given in Eq. (B2) in Pinty et al. (2006). In term (2), $f_{Rbgd}$ is the reflectance of the ground beneath the canopy and $f_{\text{UnColl, veg}}^{T}$ is the transmitted fraction of light to the ground which has not collided with any canopy elements. In term (3), $f_{\text{Coll, bgd,1}}^{fR}$ is the fraction of light which has struck vegetation and collided with the background a single time, while $f_{\text{Coll, bgd, n}}^{fR}$ is the fraction which has collided multiple times (n) with the background.

The sum of the three components results in the canopy albedo (Pinty et al., 2006). Similar equations can be derived for light originating from diffuse sources (e.g. clouds and other atmospheric scattering) (Pinty et al., 2006). Implementations of the calculations of the canopy fluxes for a single level are available from the JRC, and these implementations were used as the basis of the routines put into ORCHIDEE-CAN for both the single- and multi-level cases (McGrath et al., 2015b). This implementation relies on the use of the effective LAI, which is the LAI that needs to be used in a 1-D process representation to obtain the same reflectance, absorbance and transmittance as would be obtained by a 3-D-canopy representation (Pinty, 2004). In this study, the effective LAI was calculated by first computing the canopy gap probability, i.e. the probability that light is transmitted to a specified height in the canopy at a given solar angle. The gap probability is then converted into the effective LAI by passing it as an input to the inverted Beer-Lambert law (with an extinction coefficient of 0.5 to ensure compatibility with the two-steam inversion of Pinty et al., 2011a).

$$d_{\text{LAleff}} = -2.0 \times \cos(\theta_z) \times \log(f_{\text{Pgap}}), \quad (34)$$

where $f_{\text{Pgap}}$ can be $f_{\text{Pgap}}^{\text{trees}}$, $f_{\text{Pgap}}^{\text{gc}}$, $f_{\text{Pgap}}^{\text{bs}}$. Following the introduction of multi-layer photosynthesis and energy budget submodels, the approach proposed by Pinty (2004) had to be adjusted such that it could be applied for every level for which absorbance needs to be known to calculate photosynthesis (Sect. 3.6) and reflectance needs to be known to calculate the net short-wave radiation (Sect. 3.5). The multi-layer approach basically applies the 1-D two-stream canopy radiation transfer model by Pinty et al. (2006) to each canopy level where the light transmitted by the overlaying level becomes the input for the lower level.

As the multi-level approach is built around the solution of the one-level scheme for each canopy level, no new equations are introduced. The method can be summarised by the following algorithm for which the details are given in McGrath et al. (2015b). First, three fluxes are calculated for each level independently: the fraction of light transmitted through the layer without striking vegetation, the fraction of light reflected after striking vegetation, and the fraction of

light transmitted through the layer after striking vegetation. These three fluxes represent the only possible fate of light (any light not taking one of these paths must be absorbed for energy conservation). Next, an iterative approach is invoked which follows the path of a single photon entering the top level. Based on the solutions for each single level, probabilities can be calculated that the photon will be transmitted to a lower level or reflected to a higher level. Any fraction which is reflected upwards from the top level is added to the total canopy albedo and is not considered further. The fraction which is transmitted through the top level enters the next highest level, and again the single-level solutions determine where this light goes. Any fraction reflected upwards is considered in the next iteration as part of the light entering the upper level. The steps continue until the bottom canopy level is reached. Here, any fraction which is transmitted into the soil is removed from consideration and added to the total transmittance through the canopy. The algorithm then proceeds to the above canopy level. Now the transmitted fluxes are moving in the upwards direction towards to the sky, while reflected fluxes are moving towards the ground. The code continues towards the top level, taking as input from below both the flux reflected by downwelling light from the level below the current level and the flux transmitted from the lower level by upwelling light. After each iteration (moving from the top of the canopy to the bottom and back to the top), the total amount of light considered active has been reduced by light escaping to the sky or being absorbed by the canopy or ground. Eventually, this “active” light falls below a pre-defined threshold, and the calculation is considered to be converged.

Due to the iterative procedure, energy is not strictly conserved, although we have attempted to choose a threshold which minimises this loss. The multi-level albedo calculation is currently the most expensive part of the model, due to the iterations and the fact that it must be performed over all canopy levels (currently set to 10), grid points, and PFTs at every physical time step. Levels with no LAI are no less expensive to compute, either, although we have arranged our canopy levels to make sure no levels are empty in most cases.

## 3.5 Multi-layer energy budget

The present generation of land-surface models have difficulties in reproducing consistently the energy balances that are observed in field studies (Pitman et al., 2009; Jiménez et al., 2011; de Noblet-Ducoudré et al., 2012). The ORCHIDEE-CAN branch implemented an energy budget scheme that represents more than one canopy layer to simulate the effects of scalar gradients within the canopy for determining more accurately the net sensible and latent heat fluxes that are passed to the atmosphere. As outlined in Polcher et al. (1998), the use of an implicit solution for coupling between the atmospheric model and the surface layer model is the only way to keep profiles of temperature and humidity synchronised
across the two models when the coupled model is run over large time steps (e.g. of 30 min). The difference between explicit and implicit schemes is that an explicit scheme will calculate each value of the variable (e.g. temperature and humidity) at the current time step entirely in terms of values from the previous time step. An implicit scheme requires the solution of equations written only in terms of those at the current time step.

The modelling approach formalises three constraints that ensure energy conservation. The three equations that describe the main interactions are the following.

1. The energy balance at each layer is the sum of incoming and outgoing fluxes of latent and sensible heat and of short-wave and long-wave radiation:

$$k_{\text{lhc},i} k_{\rho_{v},i} \frac{\delta T_{\text{L},i}}{\delta t} = \left( -k_{\text{shc}} k_{\rho_{a}} \frac{(T_{\text{L},i} - T_{\text{a},i})}{R_{\text{a},i}} - k_{\lambda,\text{LE}} \rho_{a} \frac{q_{\text{L},i} - q_{\text{a},i}}{R_{\text{s},i}} + F_{\text{SW},i} + F_{\text{LW},i} \right) \left( \frac{1}{\Delta d_{\text{hl},i}} \right), \quad (35)$$

where $F_{\text{LW},i}$ is the sum total of long-wave radiation, that is, the net long-wave radiation absorbed into layer i, and $F_{\text{SW},i}$ is the net absorbed short-wave radiation as calculated by the radiation scheme in Sect. 3.4. $k_{\text{shc}}$ is the specific heat capacity of air. The source sensible heat flux from the leaf at level i is the difference between the leaf temperature ($T_{\text{L},i}$) and the atmospheric temperature at the same level ($T_{\text{a},i}$), divided by $R_{\text{a},i}$, which is the leaf resistance to sensible heat flux (a combination of stomatal and boundary layer resistance). Similarly, the source latent heat flux from the leaf at level i is the difference between the saturated humidity in the leaf ($q_{\text{L},i}$) and that in the atmosphere at level i ($q_{\text{a},i}$), divided by $R_{\text{s},i}$, which is the leaf resistance to latent heat flux. $R_{\text{a},i}$ is calculated based upon the leaf boundary layer resistance, and is described in the present model according to Baldocchi (1988). $R_{\text{s},i}$ is an abbreviation for the sum of the stomatal and leaf boundary layer resistance terms for latent heat.

2. The sensible heat flux between the vegetation (“the leaf”) and the surrounding atmosphere at each level, and between adjacent atmospheric levels above and below, is provided by the following expression:

$$\frac{\delta T_{\text{a},i}}{\delta t} \Delta d_{\text{V},i} = k_{\text{k},i} \frac{\delta^2 T_{\text{a},i}}{\delta z^2} \Delta d_{\text{A},i} - \left( \frac{T_{\text{L},i} - T_{\text{a},i}}{R_{\text{a},i}} \right) \left( \frac{1}{\Delta d_{\text{hl},i}} \right) \Delta d_{\text{V},i}, \quad (36)$$

where z denotes the height above the soil surface. We have re-written the scalar conservation equation, as applied to canopies, in terms of the sensible heat flux, temperature and source sensible heat from the vegetation at each layer.

3. The latent heat flux between the vegetation and surrounding atmosphere at each level, and between adjacent atmospheric levels above and below is described in a form that is analogous to Eq. (36), above:

$$\frac{\delta q_{\text{a},i}}{\delta t} \Delta d_{\text{V},i} = k_{\text{k},i} \frac{\delta^2 q_{\text{a},i}}{\delta z^2} \Delta d_{\text{A},i} - \left( \frac{q_{\text{L},i} - q_{\text{a},i}}{R_{\text{s},i}} \right) \left( \frac{1}{\Delta d_{\text{hl},i}} \right) \Delta d_{\text{V},i} \quad (37)$$

In addition to these three basic equations, various terms had to be parametrised. The 1-D second-order closure model of Massman and Weil (1999) was used to simulate the vertical transport coefficients $k_{\text{k},i}$ within the canopy while accounting for the vertical and horizontal distribution of LAI (Sect. 3.3). This set of equations was then written in an implicit form and solved by induction. More details on the implicit multi-layer energy budget and a complete mathematical documentation are given in Ryder et al. (2014).

To complete the energy budget calculations, the multi-layer 1-D canopy radiation transfer model (Sect. 3.4) was used to calculate the net short-wave radiation at each canopy layer. Furthermore, the canopy radiation scheme makes use of the Longwave Radiation Transfer Matrix (LRTM) (Gu, 1988; Gu et al., 1999). This approach separates the calculation of the radiation distribution completely from the implicit expression. Instead, a single source term for the long-wave radiation is added at each level. This means that the distribution of LW radiation is now explicit (i.e. makes use of information only from the “previous” and not the “current” time step), but the changes within the time step were small enough not to affect the overall stability of the model. However, an advantage of the approach is that it accounts for a higher order of reflections from adjacent levels than the single order assumed in the process above.

## 3.6 Analytical solution for photosynthesis

The photosynthesis model by Farquhar, von Caemmerer and Berry (Farquhar et al., 1980) predicts net photosynthesis of C3 plants as the minimum of the Rubisco-limited rate of CO₂ assimilation and the electron transport-limited rate of CO₂ assimilation (Farquhar et al., 1980). The ORCHIDEE-CAN branch calculates net photosynthesis following an analytical algorithm as described by Yin and Struik (2009). In addition, the C4 photosynthesis is calculated by an equivalent version of the Farquhar, von Caemmerer and Berry model that was extended to account for noncyclic electron transport (Yin and
Struik, 2009). A detailed derivation of the analytical solution of the Farquhar, von Caemmerer and Berry model is given in Yin and Struik (2009). Although the exclusion of mesophyll conductance from the photosynthesis model could lead to an underestimation of the CO₂ fertilization effect in Earth system models (Sun et al., 2014), mesophyll conductance was not included in ORCHIDEE-CAN, to maintain compatibility between the model formulation and its parametrisation. Because values of $k_{Vcmax}$ and $k_{Jmax}$ differ between different formulations of the photosynthesis model (Kattge and Knorr, 2007; Medlyn et al., 2002) and the parametrisation that was used in ORCHIDEE-CAN did not include mesophyll conductance, it was also not accounted for in the model formulation. The analytical photosynthesis model implemented in ORCHIDEE-CAN could be easily extended to include mesophyll conductance, but that would require reparameterising the photosynthesis model.

Owing to the canopy structure simulated in this model version and the layering of the canopy, the amount of absorbed light now varies with canopy depth. This new approach replaces the old scheme which uses multiple levels based on the leaf area index, not the physical height within the canopy. Photosynthesis is now calculated at each vertically resolved canopy level independently, using the total amount of absorbed light calculated by the radiation transfer scheme, which means that radiation transfer inside the canopy and photosynthesis are now fully consistent. In the new photosynthesis scheme, photosynthesis thus indirectly depends on canopy structure.

## 3.7 Forest management and natural mortality

Although forest management has developed a wide range of locally appropriate and species-specific strategies (Pretzsch, 2009), the nature of large-scale land-surface models such as ORCHIDEE-CAN requires only a limited number of contrasting strategies that are expected to be relevant on the spatial scale (e.g. 50 × 50 km) of global and regional modelling studies. Four management strategies were implemented based on their expected impact on biogeochemical and biophysical processes.

1. In unmanaged stands self-thinning drives stand dynamics and continues until too few trees are left on site. Subsequently, a stand replacing disturbance moves all standing biomass into the appropriate litter pools and a new stand is established.

2. High stand management is characterised by regular thinning and a final harvest cut. Thinning is decided on the basis of the deviation between the actual and potential stand density for any given diameter. This approach relates to the so-called relative density index (Fortin et al., 2012), the land use disturbance index (Luyssaert et al., 2011) or hemeroby and naturalness approaches (Schall and Ammer, 2013). Exceeding a threshold di-

ameter results in a clear cut and the stand is replanted in the next year. For both thinning and harvest, leaves, roots and belowground wood are transferred to the appropriate litter pools, whereas the aboveground woody biomass is removed from the site and stored in a product pool. Trees with a diameter below a species-specific threshold are stored in a short-lived product pool which mimics wood uses for fuel, paper and cardboard. Trees with larger dimensions are moved to medium- and long-lived product pools which mimic, for example, particle boards and timber usages, respectively.

3. Coppicing of the aboveground biomass is decided on stem diameter. At harvest, the root system is left intact and, in between coppicing, no wood is harvested. Note that at present it is not possible to simulate coppicing-with-standards in ORCHIDEE-CAN.

4. In ORCHIDEE-CAN, stands under short rotation management are limited to poplar (Populus spp.) and willow (Salix spp.) forests. Stands are harvested at a prescribed age. Following a set number of harvest cycles, the root system is uprooted and the whole stand is replanted.

Different age classes are distinguished to better account for the structural diversity and its possible effects on the element, energy and water fluxes. A clear hierarchy was established for the mortality processes regarding the actual killing of trees (i.e. move their biomass to the litter or harvest pools). All of the processes determine first how much biomass they would remove in the absence of all the other processes. Afterwards, the killing is arranged in the most realistic way possible. A clear-cut event has the highest priority, followed by human thinning and finally natural mortality including self-thinning. If, for example, a forest is scheduled to be clear-cut, the entire forest biomass is subjected to the rules of the clear-cut and no other mortality occurs in that time step.

In addition to forest management and natural prescribed mortality, a variety of changes have been made to processes involving vegetation mortality. A whole PFT within a grid cell is now killed if, at the end of the day, the labile pool is empty and there is no carbon available in the leaf or carbohydrate reserve pool to refill it. In this situation, it will be impossible for the plant to assimilate new carbon from the atmosphere as it will not be able to grow new leaves and thus initiate plant recovery. Furthermore, a forest can die if the density falls below a certain prescribed value. In the next time step a new young forest will be prescribed.

If a forest is thinned, it is assumed that the weakest trees will be thinned, and therefore human thinning reduces or even eliminates the natural mortality for that time step. Natural mortality still happens on a daily time step, while human-induced mortality happens only at the end of the year. Self-thinning, as described below, takes priority over environmental mortality, which is the mortality of individuals by intersects, lightening, wind, drought, frost and heart rot. Envi-
ronmental mortality is calculated by multiplying the stand biomass by an assumed mortality fraction of 1/k<sub>resid</sub>. Where self-thinning is less than this assumed environmental mortality, self-thinning is complemented by additional mortality to reach the set environmental mortality. Where self-thinning mortality exceeds the set environmental mortality, simulated self-thinning is assumed to include environmental mortality. The fire module that is available for the trunk but not for the ORCHIDEE-CAN branch simulates stand replacing fires rather than individual-tree-based mortality due to lightening. The approach implemented in the ORCHIDEE-CAN branch could therefore be extended with models that simulate stand replacing mortality from fire, insects and storms.

The use of circumference classes adds a good deal of realism and flexibility to the ORCHIDEE-CAN simulations, but it also raises additional questions. For example, which trees should be targeted by which mortality? Given that self-thinning reflects the outcome of continuous resource competition, the largest trees are expected to be most successful when competing for resources, and therefore we assume that the smallest trees die first to reduce the stand density. Conversely, larger trees are more likely to die because of environmental stress factors, being more prone to cavitation, wind damage, lightening, and, heart rot. Therefore, we select more older trees to die from environmental mortality. While doing this also trees in the other circumference classes were killed based on the following recursive definition (cf. Bellassen et al., 2010):

$$f_{\text{death}}^{\text{icir}} = \frac{f_{\text{death}}^{\text{icir}-1} \times k_{\text{ddf}}^{1-(k_{\text{ncirc}}-1)}}{m_{\text{Ndeath}}}, \quad (38)$$

where k<sub>ddf</sub> is the death distribution factor, which is the factor by which the smallest and largest circumference classes differ (e.g. k<sub>ddf</sub> = 10 means that the largest circumference class will lose 10 times as much biomass as the smallest as a result of the mortality), m<sub>Ndeath</sub> is a normalisation factor so that the sum of f<sub>death</sub><sup>icir</sup> is unity, and f<sub>death</sub><sup>1</sup> is set equal to unity before normalisation. As the stands are very close to even-aged, we set the factor k<sub>ddf</sub> to be equal to 1. This means the same number of trees is killed in each circumference class. If, for some reason, there is not enough biomass in a given class to satisfy this distribution, the extra biomass is taken from the next smallest class (in case the smallest class does not have enough, it is taken from the largest class).

Related to mortality is the question of the circumference class distribution. As mentioned above, trees in different circumference classes are preferentially killed by different processes. If the simulation is long enough (or if the mortality is aggressive enough), eventually the number of trees in some circumference classes may become 0. This would reduce the numerical resolution of the allocation scheme. When only one circumference remains populated, the scheme effectively loses its meaning, as all the newly produced biomass is now being allocated to the only remaining circumference class. In

order to maintain the same level of detail through the simulation, the distribution of all the circumference classes is recalculated at the end of each day. A normalised target distribution is specified as an input parameter (an exponential distribution is currently used), and this distribution is scaled to produce a target distribution for the current number of individuals. All of the current individuals are placed in these new classes until the target distribution is satisfied. The target distribution now contains, however, trees of multiple sizes, so we need to average them to find the new model tree for each class. By changing the size of the model tree in each class, we are able to preserve the total biomass of the stand as well as the total number of individuals. Note that the boundaries of each diameter class are recalculated at each time step; this approach is a numerically efficient alternative to fixing the boundaries of each diameter class with a varying distribution.

## 4 Description of the parametrisation

The ORCHIDEE-CAN branch was specifically developed to quantify the climate effects of forest management over Europe. Although the developments are sufficiently general to be applied outside of Europe, the model was initially parametrised for the boreal, temperate and Mediterranean climate zones and validation focused on Europe. Parametrisation of the tropical zone is subject of a follow-up study. The parametrisation of the model, including parameter optimisation and tuning, consisted of five major steps:

1. Parameters related to carbon allocation (Sect. 4.2), forest management and mortality (Sect. 4.3), hydraulic architecture (Sect. 4.4, canopy structure (Sect. 4.5)), photosynthesis (Sect. 4.6), and canopy radiation transfer (Sect. 4.7), and for which observations exist at the species level (Sect. 4.1), were extracted from a wide range of sources (Tables S1–S5). Using the extracted species-level parameter values in ORCHIDEE without further processing avoids hidden model tuning and largely reduces the likelihood that simulation results will be biased by hidden calibration owing to a poor taxonomic definition of PFTs (Scheiter et al., 2013).

2. The phenology-related parameters of the deciduous MTCs were optimised by MacBean et al. (2015), using MODIS-derived NDVI data normalised to model fAPAR over the 2000–2008 time period.

3. The modulator (m<sub>ψ</sub>) which accounts for processes in the soil-plant continuum that are currently not modelled, was manually tuned against species distribution maps (Sect. 4.4).

4. The coefficient for maintenance respiration was optimised making use of Bayesian calibration (Sect. 4.8)
against a compilation of 100+ observations of biomass production efficiency.

5. The leaf to sapwood area ratio was manually tuned (Sect. 4.9) to match 100+ site-level gross primary production (GPP) and LAI observations recorded over Europe.

## 4.1 Introducing 12 new PFTs

Similarly to the ORCHIDEE trunk, the ORCHIDEE-CAN branch distinguishes 13 metaclasses (MTC) for vegetation. Outside Europe the original MTC classification of ORCHIDEE was kept, while inside Europe 12 new parameter sets representing the main European tree species were added. The default vegetation distribution map in ORCHIDEE, i.e. Olson et al. (1983), was replaced by an up-to-date global MTC map which has been produced using the ESA CCI ECV Land Cover map (http://www.esa-landcover-cci.org/) (Poulter et al., 2015). The mapping from land cover to MTC basically followed Poulter et al. (2011), although Table 5 (the “cross-walking” table) has been updated following discussions with the LC-CCI team at Universite Catholique de Louvain. For the European domain, the global MTC distribution was overlaid by a tree species distribution map (Brus et al., 2012).

This study focusses on tree species with a coverage of more than 2% in Europe, yielding seven species groups covering in total 78.8% of the European forest area: *Betula* sp., *Fagus sylvatica*, *Pinus sylvestris*, *Picea* sp., *Pinus pinaster*, *Quercus ilex* and a group combining *Quercus robur* and *Quercus petraea*. For *Pinus sylvestris*, *Picea* sp. and *Betula* sp. An additional distinction between boreal and temperate forest was made for the species map and parametrisation: trees located in Norway, Sweden and Finland were considered boreal, while trees growing at lower latitudes were categorised as temperate. Given the potential role of tree species of the Salicacea genus in short rotation coppice management, a separate PFT was parametrised for *Populus* sp. Furthermore, to improve the parametrisation of the MTC of boreal needleaved deciduous forest, observations from *Larix* sp. were included when possible.

For these 12 forest species, 12 new PFTs were created, with each PFT belonging to a single MTC (Tables S2, S3 and S4). Almost 79% of the European forest was parametrised at the species level. The remaining 21% was reclassified into four residual groups, i.e. a temperate and boreal needleleaf evergreen and a temperate and boreal broadleaved residual group. For use outside Europe, the original MTC classification of ORCHIDEE was kept. The parameters of the residual groups and MTCs are the mean of the parameters of the species-level PFTs that are in the MTC, with the exception of albedo parameters that could be extracted from remote-sensing products. Finally, separate PFTs were introduced for boreal grasses and croplands, which allowed for a boreal

parametrisation of phenology, senescence and growth. This approach, which distinguishes a total of 28 PFTs, allows a higher taxonomic resolution over Europe, better defines forest types compared to the more general MTC approach and facilitates the use of observations to derive parameters.

## 4.2 Allocation

The allocation scheme relies on the leaf to sapwood area ratio (Sect. 4.9) and the relationship between diameter and height. Following a logarithmic transformation of the more than 150 000 data points from the national forest inventory data of Spain, France, Germany and Sweden, the two parameters (i.e. *k*<sub>α1</sub> and *k*<sub>β1</sub>) describing the relationship between diameter and height (Eq. 15) were fitted at the species level making use of a least square regression. Parameter values for MTCs were derived by grouping the species into MTCs and fitting the parameters. Data sources and parameter estimates are presented in Tables S2 and S3.

## 4.3 Forest management and mortality

Forest management and tree mortality are controlled by (Sect. 3.7): (1) maximum tree diameter (no symbolic notation; called largest_tree_diam in ORCHIDEE-CAN), (2) minimum stand density (no symbolic notation; called ntrees_dia_profit in ORCHIDEE-CAN), (3) environmental mortality (no symbolic notation; called residence_time in ORCHIDEE-CAN), (4) self-thinning (*k*<sub>α2</sub> and *k*<sub>β2</sub>) and, (5) anthropogenic thinning (no symbolic notation; called alpha_RDI_upper, alpha_RDI_lower, beta_RDI_upper and beta_RDI_lower in ORCHIDEE-CAN) where the parameters depend on the management strategy.

Maximum tree diameter was extracted from the French, Swedish, German and Spanish forest inventories as the observed 50% quantile for diameter at breast height. The 50% quantile rather than the observed maximum was used to account for the fact that large-scale land-surface models are expected to reproduce large-scale patterns rather than local extremes. Minimum stand density was estimated as the expected stand density for the maximum tree diameter for a stand under self-thinning. Although both criteria are related to each other through the observed self-thinning relationship (see below), the minimum number of trees is used to decide when unmanaged forests should be replaced, whereas both the maximum diameter and the minimum number are used for managed sites as criteria to initiate a clear cut. Parameters for anthropogenic thinning are based on the national forest inventory data and checked against the JRC database of species-specific yield tables. Parameter values are presented in Table S5. Resource competition between trees in the same stand has been reported to result in the so-called self-thinning relationship that relates the number of individuals within a stand to the stand biomass (Reineke, 1933; Kira et al., 1953;
Yoda et al., 1963):

$$ (M_s + M_h) \times k_{\rho s} = k_{\alpha} \times (d_{\text{ind}})^{-k_{\beta}}, \quad (39) $$

where $k_α$ and $k_β$ are the constants of the self-thinning relationship. Furthermore, stem volume can be written as a function of tree diameter ($d_{\text{dbh}}$), tree height and stem form factor ($k_{\alpha'}$) to account for the fact that the stem shape is not a perfect cylinder:

$$ (M_s + M_h) \times k_{\rho s} = k_{\alpha'} \times (d_{\text{dbh}})^2 \times d_h. \quad (40) $$

Following the allometric relationship given in Eq. (15), tree height can be written as a function of tree diameter. Hence, the self-thinning relationship can be re-written to relate stand diameter to stand density:

$$ d_{\text{dbh}} = k_{\alpha 2} \times (d_{\text{ind}})^{-k_{\beta 2}}, \quad (41) $$

where, $k_{\beta 2}$ relates to $k_{\beta 1}$ (as in Eq. 15) as follows:

$$ k_{\beta 2} = -3/2 \times (2 + k_{\beta 1}) \quad (42) $$

$k_{\alpha 1}$ and $k_{\beta 1}$ were estimated by fitting Eq. (15) to observed diameter and height of individual trees from NFI of Sweden, Germany, France and Spain. $k_{\beta 2}$ was calculated from Eq. (42) and $k_{\alpha 2}$ was estimated by fitting Eq. (41) to observations of the quadratic mean stand diameter and stand density from NFI data.

## 4.4 Hydraulic architecture

Initial choices of parameters for this scheme were based on the values and parameter sources listed by Hickler et al. (2006). All data sources were revisited and the search was extended to obtain values at the PFT rather than MTC level. Given that plant hydrology is rather well studied, observed parameters were available for most of the species. Data sources are listed in Table S1, whereas the parameter values are shown in Table S3. Our implementation of hydraulic architecture required the introduction of a tuning parameter ($m_ψ$) to account for processes that are currently absent in the scheme, e.g. plant water storage and soil-root resistance. A process-based description of these processes (i.e. Sperry et al., 1998; Steppe et al., 2006) is being tested and should reduce the effect of the tuning parameter and eventually allow its removal from the model.

For the time being, the modulator $m_ψ$ was tuned manually against the species distribution map to obtain a match between the simulated and observed species distributions. When the modulator is set to zero, all PFTs experience excessive water stress resulting in large-scale plant mortality. The modulator was increased until the prescribed vegetation distribution which was based on remote-sensing observations (Sect. 4.1), survived where it was prescribed. To this aim, the model was run for 50 years, forced with v5.2 of the CRU-NCEP climatology for Europe (Climatic Research

Unit, University of East Anglia). Note that the values of the modulator depend on the climate data that are used to force the model. Similarly the modulators may need to be retuned when ORCHIDEE-CAN is coupled to an atmospheric model.

## 4.5 Canopy structure

The relationship between diameter and projected crown surface area follows the model proposed by Pretzsch (2009):

$$ d_{\text{csa}} = k_{\text{ap}} \times d_{\text{dbh}}^{k_{\text{bp}}} \quad (43) $$

with parameters estimated using the data set presented in Pretzsch and Dieler (2012). This data set contains diameter and projected crown surface area observations for over 37 000 individual trees in Europe covering almost 30 species. Following logarithmic transformation of the observations a linear least square regression was used to fit species-specific parameter values. Parameter values are shown in Table S2. Parameter values for MTCs were derived by grouping the species into MTCs and fitting the parameters. No observations were available for the boreal zone and temperate evergreen deciduous species. For the boreal species, a subset of the temperate observations (*Pinus sylvestris*, *Picea abies* and *Betula pendula*) was used, i.e. the relationship between $d_{\text{csa}}$ and $d_{\text{dbh}}$ was fitted to all available data for *Pinus sylvestris*. Next, all observations with a $d_{\text{csa}}$ that falls below the predicted $d_{\text{csa}}$ were selected as considered to represent a boreal subset. Given the importance of snow pressure on crown structure, selecting observations with sub average $d_{\text{csa}}$ is justifiable as a first approximation. Subsequently, the parameters were fitted to this subset of data. For *Quercus ilex* no data were available and parameters were tuned such that the crown diameter was 0.85 m less than the tree height.

## 4.6 Analytical solution for photosynthesis

Three originally MTC-specific photosynthetic parameters ($k_{\text{Vcmax}}$, $k_{\text{Jmax}}$ and $k_{\text{sla}}$) were derived at the species level by obtaining weighted site means for each species from the TRY global leaf trait database (Kattge et al., 2011) and additionally from Medlyn et al. (2002). Only $k_{\text{Vcmax}}$ and $k_{\text{Jmax}}$ standardised to a common formulation and parametrisation of the photosynthesis model by (Farquhar et al., 1980) were used. Most $k_{\text{Vcmax}}$ and $k_{\text{Jmax}}$ values in the TRY database had already been standardised to a reference temperature of 25 °C (Kattge and Knorr, 2007). Subsequently, a species-specific $k_{\text{Jmax, opt}}/k_{\text{Vcmax, opt}}$ ratio was calculated from the records which included both $k_{\text{Vcmax, opt}}$ and $k_{\text{Jmax, opt}}$ measurements. From this ratio, which was within a range of 1.91–2.47 for each species, $k_{\text{Jmax, opt}}$ was calculated for records which originally only included $k_{\text{Vcmax}}$. Only geo-referenced observations within Europe were used and the distinction between boreal and temperate forest was made similar to the species map. Depending on the species this resulted in 5
to 183 observations for *k*<sub>sla</sub> and 11 to 173 observations for *k*<sub>Vcmax opt</sub> and *k*<sub>Jmax opt</sub>. From these observations species-specific means were calculated, weighted for differences in the number of observations per site. The parameter values are shown in Table S3.

## 4.7 Multi-layer two-way radiation scheme for tall canopies

The radiation transfer scheme makes use of parameters describing leaf and background properties, i.e. leaf single scattering and preferred scattering direction (for both visible (VIS) and near-infrared (NIR) wavelengths) and the so-called background albedo or the albedo of the surface below the dominant tree canopy (VIS and NIR). All parameters were taken from the Joint Research Centre Two-stream Inversion Package (JRC-TIP) (Pinty et al., 2011a, b). This is a software package (Pinty et al., 2007) which inverts a two-stream model (Pinty et al., 2006) to best fit the MODIS broadband visible and near-infrared white sky surface albedo from 2001 to 2010 at 1 km resolution (Pinty et al., 2011a). The inverse procedure implemented in the JRC-TIP is shown to be robust, reliable, and compliant with large-scale processing requirements (Pinty et al., 2011a). Furthermore, this package ensures the physical consistency between sets of observations, the two-stream model parameters, and radiation fluxes.

Only parameter values for which the posterior standard deviation of the probability density functions were significantly smaller than the prior standard deviation were selected from the JRC-TIP optimisation (Pinty et al., 2011a), since this condition ensures statistically significant values. Species- and MTC-specific values were derived from JRC-TIP by performing a multiple regression. This methods determines, in an objective way, how the fractions of each MTC or species explain the JRC-TIP parameter. The multiple regression was performed separately for the six parameters: the single scattering of leaves (for both VIS and NIR), the scattering direction of leaves (VIS and NIR) and the background albedo (VIS and NIR). Each JRC-TIP parameter was used as the dependent variable and the independent variables consisted of the fractions of each MTC (Poulter et al., 2015) or species (Brus et al., 2012). These fractions were used to find a linear function that best predicted each JRC-TIP parameter. The corresponding slope of a regression of each MTC or species fraction gives the MTC or species dependent JRC-TIP value. The multiple regression was performed without an intercept. To avoid pollution by the seasonal cycle, the multiple regression was applied only for the pixels of the Northern Hemisphere. Only pixels that were less than 10% covered by non-vegetative fractions where selected for the analysis and only significant results following an *F* test and positive *r*² values were selected. The derived parameter values are shown in Table S4.

## 4.8 Maintenance respiration

Both the trunk and ORCHIDEE-CAN branch reduce the definition of net primary production to biomass production; hence, carbon leaching from the roots, volatile organic emissions from the leaves, dissolved and particulate carbon losses through water fluxes and carbon subsidies to mycorrhizae are not accounted for in the model. These fluxes are (incorrectly) accounted for in the modelled autotrophic respiration. Modelled autotrophic respiration should therefore be considered an effective rather than a true value. For this reason, the basal rate of autotrophic respiration was optimised against 126 site observations of the biomass production efficiency (*k*<sub>cmaint</sub>) calculated as the ratio between annual biomass production and annual photosynthesis (Vicca et al., 2012; Campioli et al., 2015), using a Bayesian optimisation scheme. The scheme, for which more details are given in Santaren et al. (2007), uses a standard variational method based on the iterative minimisation of a cost function that measures both the model data misfit and the parameter deviations from prior knowledge (Tarantola, 2005).

The simulations that were used in the Bayesian optimisation prescribed a 20 m tall vegetation for temperate tree species, a 15 m tall vegetation for boreal tree species and a 10 m tall vegetation for Mediterranean tree species as its initial condition. This approach reduced the need for several decades of simulations to a single year to grow a mature forest. In total, the simulations were run for 10 years and covered the European domain. The first year was discarded and the ratio between modelled GPP and NPP was averaged over the remaining 9 years. Prior to the optimisation, the observations were averaged for agricultural PFTs (0.57), and deciduous (0.44) and evergreen (0.53) forest PFTs; the observed uncertainty was 0.03. The parameter values were set to range between 0.0032 and 0.160. The optimisation converged within 11 iterations and the optimised parameter values are shown in Table S2.

It remains untested how well the simulated effective autotrophic respiration represents the (rarely) observed autotrophic respiration. Note that in the cases of both the trunk and the ORCHIDEE-CAN branch of ORCHIDEE, a match between effective and observed autotrophic respiration should not be interpreted as evidence of desired model behaviour because several components of net primary production are not modelled yet.

After the optimisation of the maintenance respiration coefficient (*k*<sub>cmaint</sub>), the model simulates reasonable biomass production efficiency for a unit of photosynthesis. Hence, the final step of the parametrisation focussed on optimising the leaf area, as this is one of the main drivers of photosynthesis.

## 4.9 Sapwood to leaf area ratio

The vegetation structure simulated by the ORCHIDEE-CAN branch is sensitive to the value of *k*<sub>ls</sub> which describes the ratio
between the leaf and sapwood area of an individual tree. The available observations show a wide range within and across forest species. Dependencies of $k_{ls}$ on tree height (McDowell et al., 2002; Novick et al., 2009), tree diameter following stand thinning (Simonin et al., 2006) and $CO_2$ (Pataki et al., 2006) have been reported. Most observations, however, come from experiments where time was substituted by space which hampers teasing apart the sources of variability. Given the variation and uncertainty in the observations and the model sensitivity to this parameter, we manually tuned its value within the observed range, to match European-wide observations of leaf area index as recorded in the Database of Global Forest Ecosystem Structure and Function Luyssaert et al. (2007).

This database was used to calculate a mean and maximum observed leaf area index at the species level for the temperate and boreal region. Initially 20 year long European-wide simulations were used to simulate leaf area index of a species, when the large-scale leaf area index approached the mean target value and did not exceed the maximum value, the simulations were extended to reach 100 years for checking the temporal evolution of leaf area index. We deliberately optimised the sapwood to leaf area ratio ($k_{ls}$) by making use of stand-level data to reduce circularity with the model validation (see below).

Limited tests over a period of 100 years in a Scots pine forest at 51–52° N, 13–14° E (Fig. S1 in the Supplement) suggested that optimising $k_{cmaint}$ and $k_{ls}$ had the largest effect on the maximum LAI, which decreased by almost 17% after optimisation compared to a simulation with prior parameter values. Mean annual GPP, mean annual transpiration and basal area decreased by, respectively, 6, 6 and 7% compared to a simulation with prior parameter values (Fig. S1).

## 5 Validation

ORCHIDEE-CAN is designed as the land-surface model to be coupled to the LMDz atmospheric model. As such, future applications of ORCHIDEE-CAN are expected to be regional to global in the spatial domain and to span several years in the temporal domain. Given its anticipated uses, the ability of the model to reproduce large-scale spatial patterns as well as their inter-annual variability is essential. The first applications of the model, both offline and coupled to the atmosphere, will focus on Europe. The validation, therefore, reports performance indices both over Europe as over eight separate regions within Europe (Bellprat et al., 2012). These eight regions, which partially overlap, are defined after Bellprat et al. (2012). Furthermore, the performance indices are calculated for winter, spring, summer and autumn, and thus allow one to evaluate the capacity of the model to reproduce observed annual cycles.

In addition to the root mean square error, a land performance index (LPI) based on the principles laid out for the

Climate Performance Index (Murphy et al., 2004, their SI) was also calculated. LPI normalises the root of the squared differences between the simulations and observations by the observed spatial and temporal variance. The LPI was used to estimate the likelihood that the simulated variable belongs to the same population as the observed variable, defined as $exp(-0.5LPI^2)$. An LPI equal to 1 indicates that the model correctly reproduces the mean observed value and implies a likelihood of 61% (Murphy et al., 2004) that the simulations and observations come from the same population. Similarly, an LPI of 2 reduces this likelihood to 13%. An LPI of less than 0.32 has a likelihood of more than 95% and therefore indicates a statistically significant result.

While developing ORCHIDEE-CAN, the numerical approaches that added functionality to the code were selected on the basis of their performance at the site level (see below). Rather than running the same site-level tests for our implementation, we performed a complementary large-scale validation. The strength of our approach lies not in the details, as is the case for site-level validation, but in its width by simultaneously testing model performance for structural variables such as basal area (de Rigo et al., 2014), canopy structure (Pinty et al., 2011a) and canopy height (Simard et al., 2011), biogeochemical fluxes such as GPP (Jung et al., 2008), biophysical fluxes such as albedo (Schaaf et al., 2002) and fluxes at the interface of biogeochemistry and biophysics such as evapotranspiration (Jung et al., 2008). The selection of variables was limited by the availability of spatially explicit data-derived products for Europe.

For the validation, both the trunk and ORCHIDEE-CAN branch were run from 1850 to 1900 using CRU-NCEP climate forcing from 1901 to 1950 at 0.5 degree resolution. From 1901 until 2012, the corresponding CRU-NCEP forcing data for each year were used. Both versions used the 11 layer soil hydrology, the single-layer energy budget and the same land cover map (Poulter et al., 2015). Given that no European-wide, spatially explicit and data-derived products were found for the validation of the net carbon flux, there was no need for a carbon spin-up. For the ORCHIDEE-CAN branch, the observed tree height and basal area were compared against the simulation values at the end of 2010 (the trunk does not simulate these variables). For both the trunk and the ORCHIDEE-CAN branch, the observed GPP, evapotranspiration, effective LAI and VIS and NIR albedos were compared against monthly means between 2001 and 2010.

### 5.1 Species versus PFTs

In ORCHIDEE-CAN the PFT concept was refined by parametrising the main European tree species groups (Sect. 4.1). To evaluate the effect of the species parametrisation, we performed a companion simulation for the configuration described above, but at the MTC level. Model performance was barely affected by the use of the MTC parameters,
compared to the simulation with the species parameters (see
Fig. S2 for RMSE scores).

5.2 Allocation

In ORCHIDEE-CAN, functional relationships which vary by
species and light stress are used to allocate carbon among
the fine roots, foliage and sapwood. The allocation scheme
largely follows Zaehle and Friend (2010), who in turn was
inspired by Sitch et al. (2003). Approaches simulating allo-
cation based on functional relationships were found to out-
compete allocation schemes based on constant fractions or
resource limitation (De Kauwe et al., 2014). The ability of
these schemes to reproduce foliage, fine root and sapwood
reported in large observational data sets (for example, Luys-
saert et al., 2007) demonstrates that these schemes capture
the main observed features (Zaehle and Friend, 2010). In
addition, allocation schemes making use of functional rela-
tionships were also capable of simulating the observed ef-
fect of elevated CO2 on two mature forest ecosystems (De
Kauwe et al., 2014). Despite these successes, the schemes
were reported to be sensitive to their parametrisation. Dif-
ferences in parameters were reported to result in substan-
tial differences in the simulated allocation. The parameters
for the functional relationships used in ORCHIDEE-CAN
are given in Table S2. The main conceptual difference be-
tween the allocation scheme by Zaehle and Friend (2010) and
ORCHIDEE-CAN is that the latter was designed to simulate
one or more diameter classes.

Given that photosynthesis is still calculated at the stand
level (and thus not at the tree level) the allocation rule of
Deleuze et al. (2004) was integrated in the functional allo-
cation scheme to account for light and resource competition
within a stand. Where the functional relationships are used
to simulate carbon allocation within an individual tree of a
given diameter, the rule of Deleuze et al. (2004) allocates
carbon across the different diameter classes. The allocation
rule which models the radial increment for individual trees
in pure even-aged stands was successfully tested for Nor-
way spruce and Douglas fir stands in France (Deleuze et al.,
2004). A similar approach for modelling radial increment has
already been implemented in a version close to the trunk of
ORCHIDEE (Bellassen et al., 2010) and was able to suc-
cessfully simulate stand characteristics such as height, basal
area and stand diameter (Bellassen et al., 2011). This previ-
ous implementation differs from the current implementation
in its time resolution (which is now daily instead of yearly),
its analytical solution and the underlying allocation scheme
(which is now based on functional relationships instead of
resource limitation).

The aforementioned studies performed a detailed valida-
tion of the two approaches dealing with carbon allocation,
which were combined in ORCHIDEE-CAN. Complemen-
tary to these studies, we performed a European-wide vali-
dation of our implementation and parametrisation of these

well-tested schemes against a remote-sensing-based map of
tree height (Simard et al., 2011), upscaled eddy-covariance
observations for GPP (Jung et al., 2008) and a map of basal
area based on national forest inventory data (de Rigo et al.,
2014). The model's ability to reproduce GPP is thought to
reflect its capacity to simulate the foliage biomass, a correct
simulation of height reflects the model's capacity to simulate
aboveground woody biomass, and its capacity to reproduce
observed basal areas suggests that the interaction of stand
density and individual tree diameter are well captured.

The new implementation and parametrisation of the
within-tree and within-stand allocation schemes were found
to have a 91, 68 and 72% chance that the simulations will
reproduce the observations for GPP, tree height and basal
area for Europe, respectively (Table 3). Given that basal
area and height are not available from the trunk version
of ORCHIDEE, we could not compare the performance of
model versions in this respect. With respect to GPP, the
ORCHIDEE-CAN branch was found to outperform the trunk
by 12% and thus increased the likelihood that ORCHIDEE-
CAN is an unbiased simulator of the spatial and temporal
variability of GPP from 79 to 91%. Improved performance
of the ORCHIDEE-CAN branch compared to the trunk is
observed for all regions in summer where the RMSE of GPP
was halved from 2.5–5 to 1–2 gC m⁻² day⁻¹ (Figs. 2, 3 and
4).

Although part of the high likelihood could be due to the
fact that the observed GPP was upscaled making use of sim-
ilar climatologies being used as the forcings of the mod-
els, this circularity could neither have contributed to the im-
proved performance between the trunk and the ORCHIDEE-
CAN branch nor to the decrease in RMSE. The improve-
ments are thought to be due to structural changes to the
model such as allocation, hydraulic architecture and canopy
structure as well as to the use of more consistent parametri-
sation.

5.3 Plant water supply

Our implementation of plant hydraulic architecture was
largely based on the scheme of Hickler et al. (2006), which
was tested globally and at site level. Global simulation re-
sults for actual evapotranspiration were found to reproduce
available data (Baumgartner and Reichel, 1975; Henning,
1989). At the site level, the model agreed well with the mag-
nitude and seasonality of eddy-covariance measurements of
actual evapotranspiration for 15 European forest sites (EU-
ROFLUX), with a tendency to slightly overestimate actual
evapotranspiration for 6 sites (Hickler et al., 2006).

The maximum amount of water that can be transported
by a tree relies on the hydraulic architecture of the tree and
therefore on the capacity of the model to simulate tree and
stand dimensions as well as on the model's capacity to sim-
ulate soil water content. As an additional test, our imple-
mentation of the model was compared against the upscaled
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Root mean square error of ORCHIDEE-CAN for gross primary production, evapotranspiration, visible and near-infra-red albedo, effective leaf area index, basal area and height for different regions and periods (DJF: December–February, MAM: March–May, JJA: June–August, SON: September–November). The gray-scale of the symbols indicates the number of pixels included in the calculation. The transition from green to white indicates an RMSE of 100 %.</figcaption></figure>

eddy-covariance measurements for GPP and actual evapo-
transpiration (Jung et al., 2008). The capacity to jointly re-
produce GPP and actual evapotranspiration is an indicator
that the model successfully reproduces the coupling between
CO2 and water exchange. Model validation showed 91 and
87 % chance (compared to 79 and 45 % for the trunk) that
ORCHIDEE-CAN reproduces the upscaled GPP and actual
evapotranspiration data (Table 3, Fig. 4). The RMSE for
actual evapotranspiration during summer dropped well be-
low 1 mm day⁻¹ for most regions (Fig. 2), whereas it never
dropped below 1 mm day⁻¹ for the trunk (Fig. 3).

5.4 Canopy structure

The canopy structure model by Haverd et al. (2012) was
previously validated against ground-based LIDAR data for
several test sites with varying density, structural complexity,
layering and clumping (Lovell et al., 2012). Model-derived
canopy gap probabilities compared with observations using
a one-sample *t* test were significant for 11 out of 12 test sites.
We considered this result to be a sufficient proof to use this
canopy structure model in the ORCHIDEE-CAN branch and

added to its validation by comparing the simulated canopy
structure model over Europe against a remote-sensing-based
map of tree height (Simard et al., 2011) and the JRC-TIP
effective LAI product (Pinty et al., 2011a). The effective
LAI value expresses the capability of the canopy to inter-
cept direct radiation, and is thus associated with the proba-
bility distribution function of the canopy gaps (Haverd et al.,
2012). Thus the effective LAI contains information about the
forest structure and leaf distribution of the canopy. In the
ORCHIDEE-CAN branch, canopy structure is used to cal-
culate the albedo, roughness length, absorbed light for pho-
tosynthesis and leaf area that is coupled to the atmosphere
for e.g. transpiration and interception of precipitation.

The ORCHIDEE-CAN branch is the first branch of OR-
CHIDEE that makes use of an effective LAI to calculate
the interaction between the canopy and the atmosphere. The
LPI and RMSE of the branch, therefore, cannot be compared
against the trunk. Overall, the combined implementation of
the allocation scheme and the canopy structure model shows
a 67 % chance to reproduce the satellite-based estimates for
effective LAI. Surprisingly, effective LAI is better simulated
in spring and autumn when dynamics within the canopy are
<figure><img src="image_6.png" /></figure>

**Figure 3.** Root mean square error of ORCHIDEE trunk for gross primary production, evapotranspiration and visible and near-infrared albedo for different regions and periods (DJF: December–February; MAM: March–May; JJA: June–August; SON: September–November). The grey scale of the symbols indicates the number of pixels included in the calculation. The transition from green to white indicates an RMSE of 100%.

<figure><img src="image_7.png" /></figure>

**Figure 4.** Comparison between observations and simulations of ORCHIDEE-CAN for gross primary production and basal area over Europe. Gross primary production represents the mean for June–August between 2001–2010 and basal area is the value at the end of 2010.

substantial due to leaf on-set and senescence. For the periods when the effective LAI is expected to be most stable, i.e. summer and winter, LPI approached and frequently exceeded 1 (data not shown). Part of this shortcoming may be due to the lack of shrubs in the land cover classification. In the model, shrublands are replaced by forest and/or grasslands, likely resulting in differences between the observed and simulated canopy structure. This lapse also appears in the RMSE of effective LAI (RMSE higher than 0.8, Fig. 2).

## 5.5 Top of the canopy albedo

The radiation transfer model (Pinty et al., 2006) has been validated extensively against realistic complex three-dimensional canopy scenarios (Pinty et al., 2006) and as part of the RAdiation transfer Model Intercomparison (RAMI) project. The 1-D canopy radiation transfer model by Pinty et al. (2006) was demonstrated to accurately simulate both the amplitude and the angular variations of all radiant fluxes with respect to the solar zenith angle (Widlowski et al., 2011). In addition, the radiation transfer model and its effective values extracted from the JRC-TIP data set were successfully applied to a single forest site (Pinty et al., 2011c).

Previously we reported on the capacity of the radiation transfer model to simulate the effects of forest management on albedo (Otto et al., 2014). For the latter, forest properties were prescribed and the radiation transfer model was vali-
<table><caption>**Table 3.** Likelihood that the simulated variable comes from the same population as the data. The ORCHIDEE-trunk version does not include effective LAI, basal area and height. Note that the likelihood of Europe cannot be derived from the values of the other regions due to the overlap between regions.</caption><thead><tr><th rowspan="2"></th><th colspan="8">ORCHIDEE-CAN</th><th colspan="8">ORCHIDEE-TRUNK</th></tr><tr><th>GPP</th><th>EVAPO</th><th>ALB_NIR</th><th>ALB_VIS</th><th>EFFLAI</th><th>BA</th><th>HEIGHT</th><th>GPP</th><th>EVAPO</th><th>ALBEDO</th><th>EFFLAI</th><th>BA</th><th>HEIGHT</th><th>GPP</th><th>EVAPO</th><th>ALBEDO</th><th>EFFLAI</th><th>BA</th><th>HEIGHT</th></tr></thead><tbody><tr><td>British Isles</td><td>0.91</td><td>0.87</td><td>0.78</td><td>0.45</td><td>0.55</td><td>0.47</td><td>0.13</td><td>0.91</td><td>0.49</td><td>0.74</td><td>0.04</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Iberian Peninsula</td><td>0.80</td><td>0.80</td><td>0.73</td><td>0.65</td><td>0.60</td><td>0.09</td><td>0.66</td><td>0.65</td><td>0.37</td><td>0.25</td><td>0.04</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>France</td><td>0.86</td><td>0.90</td><td>0.92</td><td>0.46</td><td>0.60</td><td>0.66</td><td>0.60</td><td>0.69</td><td>0.46</td><td>0.75</td><td>0.02</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Mid-Europe</td><td>0.92</td><td>0.93</td><td>0.88</td><td>0.86</td><td>0.68</td><td>0.80</td><td>0.76</td><td>0.81</td><td>0.48</td><td>0.64</td><td>0.46</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Scandinavia</td><td>0.92</td><td>0.83</td><td>0.47</td><td>0.91</td><td>0.59</td><td>0.62</td><td>0.24</td><td>0.81</td><td>0.31</td><td>0.55</td><td>0.65</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Alps</td><td>0.92</td><td>0.86</td><td>0.46</td><td>0.83</td><td>0.68</td><td>0.80</td><td>0.47</td><td>0.77</td><td>0.52</td><td>0.25</td><td>0.52</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Mediterranean</td><td>0.84</td><td>0.77</td><td>0.77</td><td>0.80</td><td>0.65</td><td>0.51</td><td>0.72</td><td>0.54</td><td>0.45</td><td>0.43</td><td>0.45</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Eastern Europe</td><td>0.93</td><td>0.94</td><td>0.70</td><td>0.93</td><td>0.73</td><td>0.71</td><td>0.76</td><td>0.84</td><td>0.52</td><td>0.51</td><td>0.75</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Europe</td><td>0.91</td><td>0.87</td><td>0.71</td><td>0.92</td><td>0.67</td><td>0.72</td><td>0.68</td><td>0.79</td><td>0.45</td><td>0.61</td><td>0.69</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td></tr></tbody></table>

dated against top-of-the-canopy albedo data from five observational sites. Differences in the spatial scales between the observed and simulated albedo values were accounted for by presenting the mean June albedo during 2001–2010 (Otto et al., 2014). The simulated summertime canopy albedo falls within the range of observation. However, there occurs a slight overestimation in the near-infrared wavelength band compared to the single site measurement. Overly high near-infrared single scattering albedo values for pine, as obtained from the JRC-TIP product, are the most likely cause. The observed deviation is not due to a shortcoming in the model itself, but reflects the difficulties the JRC-TIP has with optimising parameter values in the absence of field observations in the specific case of sparse canopies (Otto et al., 2014).

For the spatial validation we use the white-sky albedo (VIS and NIR) from Moderate Resolution Imaging Spectroradiometer (MODIS, Schaaf et al., 2002) at 0.5° resolution (distributed in netCDF format by the Integrated Climate Data Center (ICDC, http://icdc.zmaw.de) University of Hamburg, Hamburg, Germany). Over large spatial and temporal domains the ORCHIDEE-CAN branch reproduces the observed VIS and NIR albedo and its variability; LPI for the albedo in the visible light is especially satisfying with a likelihood of 92% for the simulations to come from the same population as the observations (Table 3). This high overall performance index, however, hides performance issues over Scandinavia and the Alps during the snow season. The RMSE for VIS and NIR albedo without snow lies around 0.05, whereas during the snow season the RMSE increases to 0.20 (VIS) and 0.18 (NIR) over these regions (Fig. 2). When the ORCHIDEE-CAN branch is coupled to an atmospheric model, however, these deviations will only have a minor effect on the climate, owing to low incoming radiation during most of the snow season, especially in Scandinavia.

Previous validation of the radiation transfer model showed that the largest discrepancies were occurring in the near-infrared domain with a snow-covered background (Pinty et al., 2006). With the exception of the snow-covered season, the new albedo scheme, which relies on the simulated canopy structure, resulted in a substantial improvement of 0.05–0.15 compared to the trunk for the RMSE in both the VIS and NIR range in Scandinavia and the Alps (Figs. 2 and 3). The European LPI-based likelihood that our model simulations come from the same populations as the MODIS albedo increased by a remarkable 11 and 23% for, respectively, NIR and VIS albedo (from 61 and 69% for the trunk to 72 and 92% for the ORCHIDEE-CAN, Table 3).

Given that the parametrisation of the canopy radiation transfer model used in ORCHIDEE-CAN relies on MODIS, the high likelihood may not come as a surprise. However, our implementation of the radiation transfer model also relies on the simulated absorbed light, simulated GPP, simulated allocation and simulated canopy structure (which depends on mortality and forest management). In the absence of all these processes our canopy radiation transfer model is expected
to reproduce the MODIS data with a probability of 100%. Hence, the likelihood of 72 and 92% (for NIR and VIS, respectively) could also be interpreted as a verification of the aforementioned calculations; all calculations that determine the canopy structure reduce the reproducibility of the data by only 8–28% (100 to 72 or 92%).

## 5.6 Energy fluxes

The multi-layer scheme is in the process of a detailed evaluation across a range of test conditions (Ryder et al., 2014), and further validation across a range of sites is ongoing. The scheme is able to produce within-canopy temperature and humidity profiles, and successfully simulates the in-canopy radiation distribution, as well as the separation of the canopy from the soil surface. However, in order to preserve a measure of continuity with previous evaluations of the model, the multi-layer solution is here set to single-layer operation mode, which includes the effects of hydraulic limitation (Sect. 3.2) and canopy structure (Sect. 3.3) on the energy budget.

The single-layer set-up of the multi-layer solution makes use of an improved albedo estimation and is therefore expected to better simulate the net radiation that needs to be re-distributed in the canopy. This has been confirmed at a single site with a sparse canopy (Ryder et al., 2014). Furthermore, the improvements in actual evapotranspiration in addition to the low RMSE (Fig. 2) are expected to be propagated in the performance of the energy budget.

## 5.7 Forest management strategies

Model comparison has previously demonstrated that explicitly treating thinning processes is essential to reproduce local and large-scale biomass observations (Wolf et al., 2011). This finding justifies the implementation of generic approaches to forest management despite the difficulties associated with defining and quantifying forest management and its intensity (Schall and Ammer, 2013). Although the use of so-called naturalness indices, in which the current state of the forest in referenced against the potential state of the forest, has been criticised because of difficulties in defining the potential state of the forest (Schall and Ammer, 2013), such approaches were demonstrated to correctly rank different management strategies according to their intensity (Luyssaert et al., 2011).

Naturalness indices making use of only diameter and stand density or the so-called relative density index (RDI) have been previously implemented at the stand level (Fortin et al., 2012) as well as in large-scale models (Bellassen et al., 2010). This approach was shown to successfully reproduce the biomass changes during the life cycle of a forest (Bellassen et al., 2011; Fortin et al., 2012). The implementation of a forestry model based on the relative density index was reported to perform better than simple statistical models for

<figure><img src="image_8.png" /><figcaption>**Figure 5.** Impact of the different forest management strategies on an oak forest for unmanged (green), high stand (orange) and coppice (blue) compared to a Poplar short rotation coppicing (red) at 48° N, 2° E. The simulation was run without spin-up to better visualise carbon build-up in the coarse woody debris (C.W.D.) pool. Simulation cycled of a single year (1990) of climate data to minimise the inter-annual variability due to climatic year-to-year variability</figcaption></figure>

stand-level variables such as stand density, basal area, standing volume and height (Bellassen et al., 2011). Although the performance of the model was reported as less satisfying for tree-level variables, the approach is nevertheless considered reliable for modelling the effects of forest management on biomass stocks of forests across a range of scales from plot to country (Bellassen et al., 2011).

In the absence of forest management, ORCHIDEE-CAN simulates that the stands develop into tall canopy (Fig. 5a), with a high biomass (Fig. 5b), a substantial dead wood and litter pool (Fig. 5c) and no harvest (Fig. 5d). High stand management reduces the height, standing biomass and litter pools (Fig. 5a–c) but produces biomass for harvest (Fig. 5d). Under coppicing, the reduction in forest age is reflected in a shorter canopy and lower biomass and litter pools (Fig. 5a–c) compared to high stand management. The harvest is more evenly spread in time but falls below the harvest generated by high stand management (Fig. 5d). Given the shorter rotations, canopy height, standing biomass and litter pools are lower for short rotation coppicing with poplar and willow compared to all other management strategies applied on oak forest (Fig. 5a–c). Short rotation coppice was harvested every 3 years resulting in a quasi-continuous supply of woody biomass (Fig. 5d).

The forestry model implemented in ORCHIDEE-CAN is based on the RDI approach by Bellassen et al. (2010). We complemented earlier validation of such an approach over France (Bellassen et al., 2011) by a new European-wide validation for basal area. On the European scale we verified the simulated basal area and height against observed basal area
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Root mean square error (RMSE) of tree diameter for different species (shown as different markers) for different regions over France (shown as A to K). Open triangle, *Pinus sylvestris*; open circle, *Pinus pinaster*; open square, *Picea* Sp.; filled diamond, *Quercus ilex*/*suber*; filled triangle, *Betula* Sp.; filled circle, *Fagus sylvatica*; filled square, *Quercus robur*/*petraea*.</figcaption></figure>

from national forest inventories (de Rigo et al., 2014) and height from remote sensing (Simard et al., 2011). With an RMSE of 3–7 for height and 7–15 for BA, and a chance of, respectively, 68 and 72% to reproduce the data on the European scale (Table 3), our model is capable of correctly simulating the mean height and basal area but fails to capture much of the spatial variability (Fig. 4; temporal variability was not considered because the data products were only available for one time period).

Furthermore, we evaluated basal area and tree diameter at the species level for 11 regions over France, which represents a finer spatial scale than targeted by the model developments and their parametrisation. The data were extracted from the French forest inventory between 2005 and 2010 and we used the same simulations as for the European validation in the previous paragraph. We selected pixels included in the French inventory data and for both simulations and observations we calculated a moving average for the diameter and basal area per age class to then calculated the RMSE (Fig. 6). To account for intrinsic species differences in diameter and basal area, we normalised the RMSE. The normalised RMSE was lower than 30% of the mean tree diameter or mean basal area for each region for *Betula* sp., *Pinus pinaster* and *Quercus ilex*. For *Fagus sylvatica*, *Pinus sylvestris*, *Picea* sp. and *Quercus robur*/*petraea* the normalised RMSE of diameter

and basal area exceeded 50% for one to four regions for tree diameter and basal area (not shown).

The inability to fully capture the observed spatial variability in the simulation could be due to the simulation protocol that started in 1850 with 2 to 3 m tall trees all over Europe. A longer simulation accounting for the major historical changes in forest management such as the reforestation in the 1700s following an all time low in the European forest cover, the start of high stand management at the expense of coppicing in the early 1800s, and the reforestation programs following World War II (Farrell et al., 2000) is expected to improve the spatial variability in tree height and basal area. Regional deviations such as those observed on the Iberian Peninsula or over the entire Mediterranean (thus including part of the Iberian Peninsula) may be due to the lack of shrubs in the land cover map and parametrisation of the ORCHIDEE-CAN branch. Therefore the models simulates a higher stand density and higher basal area for regions where in reality shrubs occur (Fig. 4).

The parametrisation of the forestry module strongly depends on the national forest inventories from Spain, France, Germany and Sweden. Therefore verification against the same data contains little information about the model quality. Nevertheless, no time-dependent relationships were used in the ORCHIDEE-CAN branch; thus the model's capacity to reproduce the relationship between basal area and stand age, diameter and stand age or wood volume and stand age could be considered a largely independent test of the model quality. These tests were performed over eight bioclimatic regions of France and the ORCHIDEE-CAN branch was found to largely capture the time dependencies of basal area, diameter and wood volume (not shown).

## 6 Conclusions

ORCHIDEE-CAN (SVN r2290) differs from the trunk version of ORCHIDEE (SVN r2243) by the allometric-based allocation of carbon to leaf, root, wood, fruit and reserve pools; the transmittance, absorbance and reflectance of radiation within the canopy; and the vertical discretisation of the energy budget calculations. Conceptual changes towards a better process representation were made for the interaction of radiation with snow, the hydraulic architecture of plants, the representation of forest management and a numerical solution for the photosynthesis formalism of Farquhar, von Caemmerer and Berry. Furthermore, these changes were extensively linked throughout the code to improve the consistency of the model. By making use of observation-based parameters, the physiological realism of the model was improved and significant reparametrisation was done by introducing 12 new parameter sets that represent specific tree species or genera rather than a group of phylogenetically often unrelated species, as is the case in widely used plant functional types (PFTs). As PFTs have no meaning outside the
scientific community, the species-level parametrisation of the ORCHIDEE-CAN branch can deliver actionable information to decision-makers and forest owners on the implications of management strategies for the climate.

Model performance was tested against spatially explicit or upscaled data for basal area, tree height, canopy structure, GPP, albedo and evapotranspiration over Europe. The tested data streams represented biogeochemical fluxes, biophysical fluxes and forest management related vegetation characteristics. Enhanced process representation in ORCHIDEE-CAN compared to the trunk version, was found to increase model performance regarding its ability to reproduce large-scale spatial patterns of all tested data streams as well as their inter-annual variability over Europe. Although this validation approach gives us confidence in the large-scale performance of the model over Europe, additional validation is recommended for other regional applications or higher resolution studies.

## Code availability

The code and the run environment are open source (http://forge.ipsl.jussieu.fr/orchidee). Nevertheless readers interested in running ORCHIDEE-CAN are encouraged to contact the corresponding author for full details and latest bug fixes.

**The Supplement related to this article is available online at doi:10.5194/gmd-8-2035-2015-supplement.**

*Author contributions.* Developed and parametrised the ORCHIDEE-CAN model: Kim Naudts, James Ryder, Matthew J. McGrath, Juliane Otto, Sebastiaan Luyssaert, Nicolas Vuichard, Didier Solyga.

Kim Naudts, James Ryder, Matthew J. McGrath, Juliane Otto and Sebastiaan Luyssaert equally contributed to model development and parametrisation.

Evaluated the performance of the ORCHIDEE-CAN model: Kim Naudts, James Ryder, Juliane Otto, Matthew J. McGrath, Sebastiaan Luyssaert, Aude Valade, Yiying Chen, Fabienne Maignan.

Contributed Fortran code: Vanessa Haverd (canopy gaps), Bernard Pinty (albedo), Valentin Bellassen (forestry).

Provided/shared observational data sets or tools for model parametrisation: Hans Pretzsch, Päivi Merilä, Jens Kattge, Gerhard Bönisch, Matteo Campioli, Josep Penuelas, Detlef Schulze, Toon De Groote, Gonzalo Berhongaray, Yuan Yan, Philippe Peylin.

Developed driver data: Natasha MacBean.

Maintained and developed the run environment: Josefine Ghattas.

*Acknowledgements.* J. Ryder, Y. Chen, M. J. McGrath, J. Otto, K. Naudts and S. Luyssaert were funded through ERC starting grant 242564 (DOFOCO), and AV was funded through ADEME (BiCaFF). ESA ECV land cover also supported this work. The research leading to these results has received funding from the European Community's Seventh Framework Programme (FP7/2007–2013) under the Grant Agreement no. 284181-TREES4FUTURE. The authors would like to thank Daniele de Rigo for providing a basal area map for Europe.

Edited by: T. Kato

## References

Amiro, B., Barr, A., Black, T., Iwashita, H., Kljun, N., Mccaughey, J., Mogenstern, K., Murayama, S., Nesic, Z., and Orchansky, A.: Carbon, energy and water fluxes at mature and disturbed forest sites, Saskatchewan, Canada, Agr. Forest Meteorol., 136, 237–251, doi:10.1016/j.agrformet.2004.11.012, 2006a.

Amiro, B., Orchansky, A., Barr, A., Black, T., Chambers, S., Chapin III, F., Goulden, M., Litvak, M., Liu, H., McCaughey, J., McMillan, A., and Randerson, J.: The effect of post-fire stand age on the boreal forest energy balance, Agr. Forest Meteorol., 140, 41–50, doi:10.1016/j.agrformet.2006.02.014, 2006b.

Amthor, J. S.: The role of maintenance respiration in plant growth, Plant Cell Environ., 7, 561–569, doi:10.1111/1365-3040.ep11591833, 1984.

Bala, G., Caldeira, K., Wickett, M., Phillips, T. J., Lobell, D. B., Delire, C., and Mirin, A.: Combined climate and carbon-cycle effects of large-scale deforestation, Proc. Natl. Aca. Sci. USA, 104, 6550–6555, doi:10.1073/pnas.0608998104, 2007.

Baldocchi, D.: A multi-layer model for estimating sulfur dioxide deposition to a deciduous oak forest canopy, Atmos. Environ., 22, 869–884, doi:10.1016/0004-6981(88)90264-8, 1988.

Ball, J. T., Woodrow, I. E., and Berry, J. A.: A model predicting stomatal conductance and its contribution to the control of photosynthesis under different environmental conditions, in: Progress in Photosynthesis Research, edited by: Biggins, J. and Nijhoff, M., 221–224, Martinus-Nijhoff Publishers, Dordrecht, Netherlands, 1987.

Baumgartner, A. and Reichel, E.: Die Weltwasserbilanz, R. Oldenburg Verlag, Munich, 1975.

Bellassen, V., Le Maire, G., Dhôte, J., Ciais, P., and Viovy, N.: Modelling forest management within a global vegetation model – Part 1: Model structure and general behaviour, Ecol. Model., 221, 2458–2474, doi:10.1016/j.ecolmodel.2010.07.008, 2010.

Bellassen, V., le Maire, G., Guin, O., Dhôte, J., Ciais, P., and Viovy, N.: Modelling forest management within a global vegetation model – Part 2: Model validation from a tree to a continental scale, Ecol. Model., 222, 57–75, doi:10.1016/j.ecolmodel.2010.08.038, 2011.

Bellprat, O., Kotlarski, S., Lüthi, D., and Schär, C.: Objective calibration of regional climate models, J. Geophys. Res., 117, D23115, doi:10.1029/2012JD018262, 2012.

Best, M. J., Beljaars, A., Polcher, J., and Viterbo, P.: A Proposed Structure for Coupling Tiled Surfaces with the Planetary Boundary Layer, J. Hydrometeorol., 5, 1271–1278, doi:10.1175/JHM-382.1, 2004.
Betts, R. A.: Offset of the potential carbon sink from boreal forestation by decreases in surface albedo., Nature, 408, 187–90, doi:10.1038/35041545, 2000.

Bonan, G. B.: Forests and climate change: forcings, feedbacks, and the climate benefits of forests., Science (New York, N.Y.), 320, 1444–1449, doi:10.1126/science.1155121, 2008.

Bonan, G. B., Levis, S., Sitch, S., Vertenstein, M., and Oleson, K. W.: A dynamic global vegetation model for use with climate models: concepts and description of simulated vegetation dynamics, Glob. Change Biol., 9, 1543–1566, doi:10.1046/j.1365-2486.2003.00681.x, 2003.

Bonan, G. B., Williams, M., Fisher, R. A., and Oleson, K. W.: Modeling stomatal conductance in the earth system: linking leaf water-use efficiency and water transport along the soil–plant–atmosphere continuum, Geosci. Model Dev., 7, 2193–2222, doi:10.5194/gmd-7-2193-2014, 2014.

Botta, A., Viovy, N., Ciais, P., Friedlingstein, P., and Monfray, P.: A global prognostic scheme of leaf onset using satellite data, Glob. Change Biol., 6, 709–725, doi:10.1046/j.1365-2486.2000.00362.x, 2000.

Brus, D., Hengeveld, G., Walvoort, J., Goedhart, P., Heidema, A., Nabuurs, G., and Gunia, K.: Statistical mapping of tree species over Europe, Eur. J. Forest Res., 131, 145–157, 2012.

Campioli, M., Vicca, S., Luyssaert, S., Bilcke, J., Ceschia, E., Chapin III, F., Ciais, P., Fernández-Martínez, M., Malhi, Y., Obersteiner, M., Olefeldt, D., Papale, D., Piao, S., Peñuelas, J., Sullivan, P., Wang, X., Zenone, T., and Janssens, I.: Management improves the efficiency of biomass production of global terrestrial ecosystems, Nat. Geosci., submitted, 2015.

Chen, Y., Ryder, J., Naudts, K., Bastriko, V., Mcgrath, M. J., Otto, J., Launiainen, S., Ogée, J., Elbers, J. A., Foken, T., Tiedemann, F., Heinesch, B., Black, A., Haverd, V., Loustau, D., Gorsel, E. V., Knohl, A., Moors, E., Vessala, T., Ottlé, C., Pelin, P., Polcher, J., and Luyssaert, S.: Improving energy partitioning and the nighttime energy balance by implementation of a multilayer energy budget in ORCHIDEE-CAN, in preparation, 2015.

Cochard, H., Martin, R., Gross, P., and Bogeat-Triboulot, M.: Temperature effects on hydraulic conductance and water relations of Quercus robur L., J. Experiment. Botany, 51, 1255–1259, 2000.

Collatz, G., Ribas-Carbo, M., and Berry, J.: Coupled photosynthesis-stomatal conductance model for leaves of C4 plants, Aust. J. Plant Physiol., 19, 519–538, 1992.

Cox, P. M., Betts, R. A., Bunton, C. B., Essery, R. L. H., Rown-tree, P. R., and Smith, J.: The impact of new land surface physics on the GCM simulation of climate and climate sensitivity, Clim. Dynam., 15, 183–203, doi:10.1007/s003820050276, 1999.

Davin, E. L., de Noblet-Ducoudré, N., and Friedlingstein, P.: Impact of land cover change on surface climate: Relevance of the radiative forcing concept, Geophys. Res. Lett., 34, L13702, doi:10.1029/2007GL029678, 2007.

De Kauwe, M. G., Medlyn, B. E., Zaehle, S., Walker, A. P., Di- etze, M. C., Hickler, T., Jain, A. K., Luo, Y., Parton, W. J., Prentice, I. C., Smith, B., Thornton, P. E., Wang, S., Wang, Y.-P., Wårlind, D., Weng, E., Crous, K. Y., Ellsworth, D. S., Hanson, P. J., Seok Kim, H., Warren, J. M., Oren, R., and Norby, R. J.: Forest water use and water use efficiency at elevated CO2: a model-data intercomparison at two contrasting temper- ature forest FACE sites., Glob. Change Biol., 19, 1759–1779, doi:10.1111/gcb.12164, 2013.

De Kauwe, M. G., Medlyn, B. E., Zaehle, S., Walker, A. P., Di- etze, M. C., Wang, Y.-P., Luo, Y., Jain, A. K., El-Masri, B., Hick- ler, T., Wårlind, D., Weng, E., Parton, W. J., Thornton, P. E., Wang, S., Prentice, I. C., Asao, S., Smith, B., McCarthy, H. R., Iversen, C. M., Hanson, P. J., Warren, J. M., Oren, R., and Norby, R. J.: Where does the carbon go? A model-data intercompari- son of vegetation carbon allocation and turnover processes at two temperate forest free-air CO2 enrichment sites, New Phytologist, 203, 883–899, doi:10.1111/nph.12847, 2014.

de Noblet-Ducoudré, N., Boisier, J.-P., Pitman, A., Bonan, G. B., Brovkin, V., Cruz, F., Delire, C., Gayler, V., van den Hurk, B. J. J. M., Lawrence, P. J., van der Molen, M. K., Müller, C., Reick, C. H., Strengers, B. J., and Voldoire, A.: Determining Robust Impacts of Land-Use-Induced Land Cover Changes on Surface Climate over North America and Eurasia: Results from the First Set of LUCID Experiments, J. Climate, 25, 3261–3281, doi:10.1175/JCLI-D-11-00338.1, 2012.

de Rigo, D., Caudullo, G., Busetto, L., and San Miguel, J.: Supporting EFSA assessment of the EU environmental suitability for exotic forestry pests: Final Report, Tech. rep., EFSA Supporting publications, 2014.

de Rosnay, P.: Impact of a physically based soil water flow and soil-plant interaction representation for modeling large- scale land surface processes, J. Geophys. Res., 107, 4118, doi:10.1029/2001JD000634, 2002.

de Rosnay, P. and Polcher, J.: Modelling root water uptake in a com- plex land surface scheme coupled to a GCM, Hydrol. Earth Syst. Sci., 2, 239–255, doi:10.5194/hess-2-239-1998, 1998.

Deleuze, C., Pain, O., Dhôte, J.-F., and Hervé, J.-C.: A flexible radial increment model for individual trees in pure even-aged stands, Ann. Forest Sci., 61, 327–335, doi:10.1051/forest:2004026, 2004.

Dickinson, R., Henderson-Sellers, A., Kennedy, P., and Wilson, M.: Biosphere-Atmosphere Transfer Scheme (BATS) for the NCAR Community Climate Model, Tech. Rep. December, 1986.

Dixon, R. K., Solomon, A. M., Brown, S., Houghton, R. A., Trexier, M. C., and Wisniewski, J.: Carbon pools and flux of global forest ecosystems, Science, 263, 185–190, doi:10.1126/science.263.5144.185, 1994.

d'Orgeval, T., Polcher, J., and de Rosnay, P.: Sensitivity of the West African hydrological cycle in ORCHIDEE to in- filtration processes, Hydrol. Earth Syst. Sci., 12, 1387–1401, doi:10.5194/hess-12-1387-2008, 2008.

Ducoudré, N. I., Laval, K., and Perrier, A.: SECHIBA, a new set of parametrizations of the hydrologic exchanges at the land- atmosphere interface within the LMD atmospheric general cir- culation model, J. Climate, 6, 248–273, 1993.

Dufresne, J. and Ghattas, J.: Description du schéma de la couche limite turbulente et l'interface avec la surface planetaire dans LMDZ, 2009.

Farquhar, G. D., von Caemmerer, S., and Berry, J. A.: A biochem- ical model of photosynthetic CO2 assimilation in leaves of C3 species, Planta, 149, 78–90, doi:10.1007/BF00386231, 1980.

Farrell, E. P., Führer, E., Ryan, D., Andersson, F., Hüttl, R., and Piussi, P.: European forest ecosystems: building the future on the legacy of the past, Forest Ecol. Manage., 132, 5–20, doi:10.1016/S0378-1127(00)00375-3, 2000.

Flexas, J., Bota, J., Galmés, J., Medrano, H., and Ribas- Carbó, M.: Keeping a positive carbon balance under ad-
verse conditions: responses of photosynthesis and respi ra tion to water stress, Physiologia Plantarum, 127, 343–352, doi:10.1111/j.1399-3054.2006.00621.x, 2006.

Fortin, M., Ningre, F., Robert, N., and Mothe, F.: Quantifying
the impact of forest management on the carbon balance of
the forest-wood product chain: A case study applied to even-
aged oak stands in France, Forest Ecol. Manage., 279, 176–188,
doi:10.1016/j.foreco.2012.05.031, 2012.

Friedlingstein, P., Joel, G., Field, C. B., and Fung, I. Y.: Toward an allocation scheme for global terrestrial carbon models, Glob. Change Biol., 5, 755–770, doi:10.1046/j.1365-2486.1999.00269.x, 1999.

Gimmi, U., Poulter, B., Wolf, A., Portner, H., Weber, P., and Bürgi, M.: Soil carbon pools in Swiss forests show legacy effects from historic forest litter raking, Landscape Ecology, 28, 835–846, doi:10.1007/s10980-012-9778-4, 2012.

Gu, L.: Longwave radiative transfer in plant canopies, Ph.D. thesis, Universit of Virginia, 1988.

Gu, L., Shugart, H. H., Fuentes, J. D., Black, T., and Shewchuk, S. R.: Micrometeorology, biophysical exchanges and NEE decomposition in a two-story boreal forest – development and test of an integrated model, Agr. Forest Meteorol., 94, 123–148, doi:10.1016/S0168-1923(99)00006-4, 1999.

Haverd, V., Lovell, J., Cuntz, M., Jupp, D., Newnham, G., and Sea, W.: The Canopy Semi-analytic Pgap And Radiative Transfer (CanSPART) model: Formulation and application, Agr. Forest Meteorol., 160, 14–35, doi:10.1016/j.agrformet.2012.01.018, 2012.

Henning, D.: Atlas of the surface heat balance of the continents, Gebrüder Bornträger, Berlin, Stuttgart, 1989.

Hickler, T., Prentice, I. C., Smith, B., Sykes, M. T., and Zaehle, S.: Implementing plant hydraulic architecture within the LPJ Dynamic Global Vegetation Model, Global Ecol. Biogeogr., 15, 567–577, doi:10.1111/j.1466-822X.2006.00254.x, 2006.

Hourdin, F.: A new representation of the absorption by the CO 2 15-µm band for a Martian general circulation model, J. Geophys. Res., 97, 18319, doi:10.1029/92JE01985, 1992.

Jackson, R. B., Jobbágy, E. G., Avissar, R., Roy, S. B., Barrett, D. J., Cook, C. W., Farley, K. A., le Maitre, D. C., McCarl, B. A., and Murray, B. C.: Trading water for carbon with biological carbon sequestration, Science, 310, 1944–1947, doi:10.1126/science.1119282, 2005.

Jiménez, C., Prigent, C., Mueller, B., Seneviratne, S. I., Mc-
Cabe, M. F., Wood, E. F., Rossow, W. B., Balsamo, G., Betts,
A. K., Dirmeyer, P. A., Fisher, J. B., Jung, M., Kanamitsu,
M., Reichle, R. H., Reichstein, M., Rodell, M., Sheffield, J.,
Tu, K., and Wang, K.: Global intercomparison of 12 land
surface heat flux estimates, J. Geophys. Res., 116, D02102,
doi:10.1029/2010JD014545, 2011.

Jung, M., Verstraete, M., Gobron, N., Reichstein, M., Papale, D., Bondeau, A., Robustelli, M., and Pinty, B.: Diagnostic assessment of European gross primary production, Glob. Change Biol., 14, 2349–2364, doi:10.1111/j.1365-2486.2008.01647.x, 2008.

Kattge, J. and Knorr, W.: Temperature acclimation in a biochem-
ical model of photosynthesis: a reanalysis of data from 36
species, Plant Cell Environ., 30, 1176–1190, doi:10.1111/j.1365-
3040.2007.01690.x, 2007.

Kattge, J., Díaz, S., Lavorel, S., Prentice, I. C., Leadley, P., Bönisch, G., Garnier, E., Westoby, M., Reich, P. B., Wright, I. J., Cor-

nelissen, J. H. C., Violle, C., Harrison, S. P., Van Bodegom, P. M., Reichstein, M., Enquist, B. J., Soudzilovskaia, N. A., Ackery, D. D., Anand, M., Atkin, O., Bahn, M., Baker, T. R., Baldocchi, D., Bekker, R., Blanco, C. C., Blonder, B., Bond, W. J., Bradstock, R., Bunker, D. E., Casanoves, F., Cavender-Bares, J., Chambers, J. Q., Chapin III, F. S., Chave, J., Coomes, D., Cornwell, W. K., Craine, J. M., Dobrin, B. H., Duarte, L., Durka, W., Elser, J., Esser, G., Estiarte, M., Fagan, W. F., Fang, J., Fernández-Méndez, F., Fidelis, A., Finegan, B., Flores, O., Ford, H., Frank, D., Freschet, G. T., Fyllas, N. M., Gallagher, R. V., Green, W. A., Gutierrez, A. G., Hickler, T., Higgins, S. I., Hodgson, J. G., Jalili, A., Jansen, S., Joly, C. A., Kerkhoff, A. J., Kirkup, D., Kitajima, K., Kleyer, M., Klotz, S., Knops, J. M. H., Kramer, K., Kühn, I., Kurokawa, H., Laughlin, D., Lee, T. D., Leishman, M., Lens, F., Lenz, T., Lewis, S. L., Lloyd, J., Llusià, J., Louault, F., Ma, S., Mahecha, M. D., Manning, P., Mas-sad, T., Medlyn, B. E., Messier, J., Moles, A. T., Müller, S. C., Nadrowski, K., Naeem, S., Niinemets, U., Nöllert, S., Nüske, A., Ogaya, R., Oleksyn, J., Onipchenko, V. G., Onoda, Y., Ordoñez, J., Overbeck, G., Ozinga, W. A., Patiño, S., Paula, S., Pausas, J. G., Peñuelas, J., Phillips, O. L., Pillar, V., Poorter, H., Poorter, L., Poschlod, P., Prinzing, A., Proulx, R., Rammig, A., Reinsch, S., Reu, B., Sack, L., Salgado-Negret, B., Sardans, J., Shiodera, S., Shipley, B., Siefert, A., Sosinski, E., Soussana, J.-F., Swaine, E., Swenson, N., Thompson, K., Thornton, P., Waldram, M., Weiher, E., White, M., White, S., Wright, S. J., Yguel, B., Zaehle, S., Zanne, A. E., and Wirth, C.: TRY – a global database of plant traits, Glob. Change Biol., 17, 2905–2935, doi:10.1111/j.1365-2486.2011.02451.x, 2011.

Kira, T., Ogawa, H., and Sakazaki, N.: Intraspecific competition among higher plants. I. Competition-yield-density interrelationship in regularly dispersed populations, Journal of the Institute of Polytechnics (Osaka University), 4, 1–16, 1953.

Krinner, G., Nicolas, V., de Noblet-Ducoudre, N., Ogée, J., Polcher, J., Friedlingstein, P., Ciais, P., Sitch, S., and Prentice, I.: A dynamic global vegetation model for studies of the coupled atmosphere-biosphere system, Global Biogeochem. Cy., 19, GB1015, doi:10.1029/2003GB002199, 2005.

Lardy, R., Bellocchi, G., and Soussana, J.-F.: A new method to de-
termine soil organic carbon equilibrium, Environ. Model. Softw.,
1759–1763, doi:10.1016/j.envsoft.2011.05.016, 2011.

Lovell, J., Haverd, V., Jupp, D., and Newnham, G.: The Canopy Semi-analytic Pgap And Radiative Transfer (CanSPART) model: Validation using ground based lidar, Agr. Forest Meteorol., 158-159, 1–12, doi:10.1016/j.agrformet.2012.01.020, 2012.

Luyssaert, S., Inglima, I., Jung, M., Richardson, A. D., Reich-
stein, M., Papale, D., Piao, S. L., Schulze, E. D., Wingate,
L., Matteucci, G., Aragao, L., Aubinet, M., Beer, C., Bern-
hofer, C., Black, K. G., Bonal, D., Bonnefond, J. M., Cham-
bers, J., Ciais, P., Cook, B., Davis, K. J., Dolman, A. J., Gie-
len, B., Goulden, M., Grace, J., Granier, A., Grelle, A., Griffis,
T., Grünwald, T., Guidolotti, G., Hanson, P. J., Harding, R.,
Hollinger, D. Y., Hutyra, L. R., Kolari, P., Kruijt, B., Kutsch,
W., Lagergren, F., Laurila, T., Law, B., Le Maire, G., Lindroth,
A., Loustau, D., Malhi, Y., Mateus, J., Migliavacca, M., Mis-
son, L., Montagnani, L., Moncrieff, J., Moors, E., Munger, J. W.,
Nikinmaa, E., Ollinger, S. V., Pita, G., Rebmann, C., Roup-
sard, O., Saigusa, N., Sanz, M. J., Seufert, G., Sierra, C., Smith,
M. L., Tang, J., Valentini, R., Vesala, T., and Janssens, I. A.:
ter (MISR) albedo products, J. Geophys. Res., 112, D10116, doi:10.1029/2006JD008105, 2007.

Pinty, B., Andredakis, I., Clerici, M., Kaminski, T., Taberner, M., Verstraete, M. M., Gobron, N., Plummer, S., and Widlowski, J.-L.: Exploiting the MODIS albedos with the Two-stream Inversion Package (JRC-TIP): 1. Effective leaf area index, vegetation, and soil properties, J. Geophys. Res., 116, D09105, doi:10.1029/2010JD015372, 2011a.

Pinty, B., Clerici, M., Andredakis, I., Kaminski, T., Taberner, M., Verstraete, M. M., Gobron, N., Plummer, S., and Widlowski, J.-L.: Exploiting the MODIS albedos with the Two-stream Inversion Package (JRC-TIP): 2. Fractions of transmitted and absorbed fluxes in the vegetation and soil layers, J. Geophys. Res., 116, D09106, doi:10.1029/2010JD015373, 2011b.

Pinty, B., Jung, M., Kaminski, T., Lavergne, T., Mund, M., Plummer, S., Thomas, E., and Widlowski, J.-L.: Evaluation of the JRC-TIP 0.01° products over a mid-latitude deciduous forest site, Remote Sensing of Environment, 115, 3567–3581, doi:10.1016/j.rse.2011.08.018, 2011c.

Pitman, A. J., de Noblet-Ducoudré, N., Cruz, F. T., Davin, E. L., Bonan, G. B., Brovkin, V., Claussen, M., Delire, C., Ganzeveld, L., Gayler, V., van den Hurk, B. J. J. M., Lawrence, P. J., van der Molen, M. K., Müller, C., Reick, C. H., Seneviratne, S. I., Strengers, B. J., and Voldoire, A.: Uncertainties in climate responses to past land cover change: First results from the LU-CID intercomparison study, Geophys. Res. Lett., 36, L14814, doi:10.1029/2009GL039076, 2009.

Polcher, J., McAvaney, B., Viterbo, P., Gaertner, M.-A., Hahmann, A., Mahfouf, J.-F., Noilhan, J., Phillips, T., Pitman, A., Schlosser, C., Schulz, J.-P., Timbal, B., Verseghy, D., and Xue, Y.: A proposal for a general interface between land surface schemes and general circulation models, Glob. Planet. Change, 19, 261–276, doi:10.1016/S0921-8181(98)00052-6, 1998.

Poulter, B., Ciais, P., Hodson, E., Lischke, H., Maignan, F., Plummer, S., and Zimmermann, N. E.: Plant functional type mapping for earth system models, Geosci. Model Dev., 4, 993–1010, doi:10.5194/gmd-4-993-2011, 2011.

Poulter, B., MacBean, N., Hartley, A., Khlystova, I., Arino, O., Betts, R., Bontemps, S., Boettcher, M., Brockmann, C., De-fourny, P., Hagemann, S., Herold, M., Kirches, G., Lamarche, C., Lederer, D., Ottlé, C., Peters, M., and Peylin, P.: Plant functional type classification for Earth System Models: results from the European Space Agency’s Land Cover Climate Change Initiative, Geosci. Model Dev. Discuss., 8, 429–462, doi:10.5194/gmdd-8-429-2015, 2015.

Pretzsch, H.: Forest dynamics, growth and yield, Springer-Verlag, Berlin, 2009.

Pretzsch, H. and Dieler, J.: Evidence of variant intra- and interspecific scaling of tree crown structure and relevance for allometric theory, Oecologia, 169, 637–49, doi:10.1007/s00442-011-2240-5, 2012.

Reick, C. H., Raddatz, T., Brovkin, V., and Gayler, V.: Representation of natural and anthropogenic land cover change in MPI-ESM, J. Adv. Model. Earth Syst., 5, 459–482, doi:10.1002/jame.20022, 2013.

Reineke, L.: Perfecting a stand-density index for even-aged forests, J. Agr. Res., 46, 627–638, 1933.

Richards, L. A.: Capillary conduction of liquids through porous mediums, Physics, 1, p. 318, doi:10.1063/1.1745010, 1931.

Ruimy, A., Dedieu, G., and Saugier, B.: TURC: A diagnostic model of continental gross primary productivity and net primary productivity, Global Biogeochem. Cy., 10, 269–285, doi:10.1029/96GB00349, 1996.

Ryan, M.: The effects of climate change on plant respiration, Ecol. Appl., 1, 157–167, 1991.

Ryder, J., Polcher, J., Peylin, P., Ottlé, C., Chen, Y., van Gorsel, E., Haverd, V., McGrath, M. J., Naudts, K., Otto, J., Valade, A., and Luyssaert, S.: A multi-layer land surface energy budget model for implicit coupling with global atmospheric simulations, Geosci. Model Dev. Discuss., 7, 8649–8701, doi:10.5194/gmdd-7-8649-2014, 2014.

Santaren, D., Peylin, P., Viovy, N., and Ciais, P.: Optimizing a process-based ecosystem model with eddy-covariance flux measurements: A pine forest in southern France, Global Biogeochem. Cy., 21, GB2013, doi:10.1029/2006GB002834, 2007.

Schaaf, C. B., Gao, F., Strahler, A. H., Lucht, W., Li, X., Tsang, T., Strugnell, N. C., Zhang, X., Jin, Y., Muller, J.-P., Lewis, P., Barnsley, M., Hobson, P., Disney, M., Roberts, G., Dunderdale, M., Doll, C., D'Entremont, R. P., Hu, B., Liang, S., Privette, J. L., and Roy, D.: First operational BRDF, albedo nadir reflectance products from MODIS, Remote Sens. Environ., 83, 135–148, doi:10.1016/S0034-4257(02)00091-3, 2002.

Schall, P. and Ammer, C.: Can land use intensity be reliably quantified by using a single self-thinning relationship?, Ecol. Appl., 23, 675–677, doi:10.1890/12-0847.1, 2013.

Scheiter, S., Langan, L., and Higgins, S. I.: Next-generation dynamic global vegetation models: learning from community ecology, The New Phytologist, 198, 957–69, doi:10.1111/nph.12210, 2013.

Shinozaki, K., Yoda, K., Hozumi, K., and Kira, T.: A quantitative analysis of plant form-the pipe model theory: I. Basic analyses, Jpn. J. Ecol., 14, 97–105, 1964.

Simard, M., Pinto, N., Fisher, J. B., and Baccini, A.: Mapping forest canopy height globally with spaceborne lidar, J. Geophys. Res., 116, G04021, doi:10.1029/2011JG001708, 2011.

Simonin, K., Kolb, T. E., Montes-Helu, M., and Koch, G. W.: Restoration thinning and influence of tree size and leaf area to sapwood area ratio on water relations of Pinus ponderosa, Tree Physiol., 26, 493–503, doi:10.1093/treephys/26.4.493, 2006.

Sitch, S., Smith, B., Prentice, I. C., Arneth, A., Bondeau, A., Cramer, W., Kaplan, J. O., Levis, S., Lucht, W., Sykes, M. T., Thonicke, K., and Venevsky, S.: Evaluation of ecosystem dynamics, plant geography and terrestrial carbon cycling in the LPJ dynamic global vegetation model, Glob. Change Biol., 9, 161–185, doi:10.1046/j.1365-2486.2003.00569.x, 2003.

Slatyer, R.: Plant-Water Relationships, vol. 158, Academic Press, New York, 1967.

Sperry, J. S., Adler, F. R., Campbell, G. S., and Comstock, J. P.: Limitation of plant water use by rhizosphere and xylem conductance: results from a model, Plant Cell Environ., 21, 347–359, doi:10.1046/j.1365-3040.1998.00287.x, 1998.

Steppe, K., De Pauw, D. J. W., Lemeur, R., and Vanrolleghem, P. A.: A mathematical model linking tree sap flow dynamics to daily stem diameter fluctuations and radial stem growth, Tree Physiol., 26, 257–273, doi:10.1093/treephys/26.3.257, 2006.

Stöckli, R., Lawrence, D. M., Niu, G.-Y., Oleson, K. W., Thornton, P. E., Yang, Z.-L., Bonan, G. B., Denning, A. S., and Running, S. W.: Use of FLUXNET in the Community
Land Model development, J. Geophys. Res., 113, G01025, doi:10.1029/2007JG000562, 2008.

Sun, Y., Gu, L., Dickinson, R. E., Norby, R. J., Pallardy, S. G., and Hoffman, F. M.: Impact of mesophyll diffusion on estimated global land CO2 fertilization, Proc. Natl. Aca. Sci., 111, 15774–15779, doi:10.1073/pnas.1418075111, 2014.

Tarantola, A.: Inverse Problem Theory and Methods for Model Parameter Estimation, SIAM, Philadelpia, 2005.

Thornton, P. E. and Rosenbloom, N. A.: Ecosystem model spin-up: Estimating steady state conditions in a coupled terrestrial carbon and nitrogen cycle model, Ecol. Model., 189, 25–48, doi:10.1016/j.ecolmodel.2005.04.008, 2005.

Tyree, M. T. and Sperry, J. S.: Vulnerability of xylem to cavitation and embolism, Ann. Rev. Plant Physiol. Molecular Biol., 40, 19–38, 1989.

Van Genuchten, M.: A Closed-form Equation for Predicting the Hydraulic Conductivity of Unsaturated Soils, Soil Sci. Soc. Am., 44, 892–898, 1980.

Verhoef, A. and Egea, G.: Modeling plant transpiration under limited soil water: Comparison of different plant and soil hydraulic parameterizations and preliminary implications for their use in land surface models, Agr. Forest Meteorol., 191, 22–32, doi:10.1016/j.agrformet.2014.02.009, 2014.

Vicca, S., Luyssaert, S., Peñuelas, J., Campioli, M., Chapin, F. S., Ciais, P., Heinemeyer, A., Högberg, P., Kutsch, W. L., Law, B. E., Malhi, Y., Papale, D., Piao, S. L., Reichstein, M., Schulze, E. D., and Janssens, I.: Fertile forests produce biomass more efficiently., Ecol. Lett., 15, 520–526, doi:10.1111/j.1461-0248.2012.01775.x, 2012.

Viovy, N. and de Noblet-Ducoudré, N.: Coupling water and carbon cycle in the biosphere, Sci. Geol. Bull., 50, 109–121, 1997.

Weatherley, P. E.: Water Uptake and Flow in Roots, in: Physiological Plant Ecology II, Springer Berlin Heidelberg, 79–109, doi:10.1007/978-3-642-68150-9_4, 1982.

Whitehead, D.: Regulation of stomatal conductance and transpiration in forest canopies, Tree Physiol., 18, 633–644, doi:10.1093/treephys/18.8.9.633, 1998.

Widlowski, J.-L., Pinty, B., Clerici, M., Dai, Y., De Kauwe, M., de Ridder, K., Kallel, A., Kobayashi, H., Lavergne, T., Ni-Meister, W., Olchev, A., Quaife, T., Wang, S., Yang, W., Yang, Y., and Yuan, H.: RAMI4PILPS: An intercomparison of formulations for the partitioning of solar radiation in land surface models, J. Geophys. Res., 116, G02019, doi:10.1029/2010JG001511, 2011.

Wolf, A., Ciais, P., Bellassen, V., Delbart, N., Field, C. B., and Berry, J. A.: Forest biomass allometry in global land surface models, Global Biogeochem. Cy., 25, 1546–1556, doi:10.1029/2010GB003917, 2011.

Xia, J. Y., Luo, Y. Q., Wang, Y.-P., Weng, E. S., and Hararuk, O.: A semi-analytical solution to accelerate spin-up of a coupled carbon and nitrogen land model to steady state, Geosci. Model Dev., 5, 1259–1271, doi:10.5194/gmd-5-1259-2012, 2012.

Yang, Z., Dickinson, R., Robock, A., and Vinnikov, K.: Validation of the snow submodel of the biosphere – atmosphere transfer scheme with Russian snow cover and meteorological observational data, J. Climate, 10, 353–373, 1997.

Yin, X. and Struik, P.: C3 and C4 photosynthesis models: An overview from the perspective of crop modelling, NJAS – Wageningen Journal of Life Sciences, 57, 27–38, doi:10.1016/j.njas.2009.07.001, 2009.

Yoda, K., Kira, T., Ogawa, H., and Hozumi, K.: Self-thinning in overcrowded pure stands under cultivated and natural conditions, J. Inst. Polytech. (Osaka University), 14, 107–129, 1963.

Zaehle, S. and Friend, A. D.: Carbon and nitrogen cycle dynamics in the O-CN land surface model: 1. Model description, site-scale evaluation, and sensitivity to parameter estimates, Global Biogeochem. Cy., 24, GB1005, doi:10.1029/2009GB003521, 2010.

Zhao, K. and Jackson, R. B.: Biophysical forcings of land-use changes from potential forestry activities in North America, Ecol. Monogr., 84, 329–353, doi:10.1890/12-1705.1, 2014.