20 th century changes in carbon isotopes and water-use efficiency : tree-ring-based evaluation of the CLM 4 . 5 and LPX-Bern models

Measurements of the stable carbon isotope ratio (δ13C) on annual tree rings offer new opportunities to evaluate mechanisms of variations in photosynthesis and stomatal conductance under changing CO2 and climate conditions, especially in conjunction with process-based biogeochemical model simulations. The isotopic discrimination is indicative of the ratio between the CO2 partial pressure in the intercellular cavities and the atmosphere (ci/ca) and of the ratio of assimilation to stomatal conductance, termed intrinsic wateruse efficiency (iWUE). We performed isotope-enabled simulations over the industrial period with the land biosphere module (CLM4.5) of the Community Earth System Model and the Land Surface Processes and Exchanges (LPX-Bern) dynamic global vegetation model. Results for C3 tree species show good agreement with a global compilation of δ13C measurements on leaves, though modeled 13C discrimination by C3 trees is smaller in arid regions than measured. A compilation of 76 tree-ring records, mainly from Europe, boreal Asia, and western North America, suggests on average small 20th century changes in isotopic discrimination and in ci/ca and an increase in iWUE of about 27 % since 1900. LPXBern results match these century-scale reconstructions, supporting the idea that the physiology of stomata has evolved to optimize trade-offs between carbon gain by assimilation and water loss by transpiration. In contrast, CLM4.5 simulates an increase in discrimination and in turn a change in iWUE that is almost twice as large as that revealed by the tree-ring data. Factorial simulations show that these changes are mainly in response to rising atmospheric CO2. The results suggest that the downregulation of ci/ca and of photosynthesis by nitrogen limitation is possibly too strong in the standard setup of CLM4.5 or that there may be problems associated with the implementation of conductance, assimilation, and related adjustment processes on long-term environmental changes.

unequivocal testimony to the input of isotopically light fossil and terrestrial carbon by human activities (Keeling et al., 1979;Francey et al., 1999). The δ 13 C data representing atmospheric air are used to quantify the global ocean and land carbon sources and sinks (Keeling et al., 1989;Joos and Bruno, 1998;Trudinger et al., 2002;Bauska et al., 2015). Furthermore, δ 13 C observations allow identification of the imprint of fossil-fuel carbon in atmospheric air to quantify regional-tolocal-scale land carbon sources and sinks (Torn et al., 2011;Vardag et al., 2016), or to evaluate air-sea transfer velocity parameterizations (Krakauer et al., 2006). The δ 13 C data from the modern ocean are applied to infer the oceanic uptake of anthropogenic carbon (Heimann and Maier-Reimer, 1996;Gruber et al., 1999;Sonnerup and Quay, 2012;Becker et al., 2016), while paleoproxy δ 13 C data from ocean sediments and ice cores permit inference of land carbon changes between the last glacial maximum and the current warm period (Shackleton, 1977;Ciais et al., 2012;Peterson et al., 2014). Paleo-δ 13 C data are also used to trace water mass, circulation and biological productivity changes on glacialinterglacial timescales and during past abrupt events Schmittner and Somes, 2016), to disentangle processes of past glacial-interglacial carbon-cycle changes Schneider et al., 2013;Eggleston et al., 2016), and of ancient climate events (Kennett and Stott, 1991;Korte and Kozur, 2010).
Box models, ocean-and land-only models, and Earth system models of intermediate complexity (Siegenthaler and Joos, 1992;Aranibar et al., 2006;Lai et al., 2006;Tschumi et al., 2011;Holden et al., 2013;Schmittner et al., 2013) have been traditionally evaluated by δ 13 C data and used for the interpretation of δ 13 C observations. Yet, despite such potential, 13 C has only been implemented recently in comprehensive Earth system models and their subcomponents (Tagliabue and Bopp, 2008;Oleson et al., 2013;Jahn et al., 2015). Now carbon isotopes have been implemented in the ocean component of the Community Earth System Model (CESM) (Jahn et al., 2015). In this paper we present the implementation of δ 13 C in the CESM land module.
Isotopic discrimination of plants following the C3 photosynthesis pathway depends on CO 2 assimilation and stomatal conductance (Farquhar et al., 1982;Lloyd and Farquhar, 1994), which themselves depend on the availability of nitrogen, water, and light, as well as species-specific leaf traits. Discrimination and variations thereof are thus indicative of the extent to which carbon assimilation by plants, fueling plant growth, is limited by factors such as drought and nitrogen limitation. In other words, quantification of isotopic discrimination changes over time and permits the evaluation of responses of stomatal conductance and assimilation to environmental variation and extremes. Environmental changes on the policy-relevant timescale of global warming include increasing atmospheric CO 2 , climate change, and increasing nitrogen deposition. Tree-ring δ 13 C records capture the influence of local climate variability with measured δ 13 C vari-ations used to reconstruct, for example, fluctuations in temperature (Treydte et al., 2009;Sidorova et al., 2013), precipitation (Schubert and Jahren, 2011), or cloud cover Young et al., 2012).
The intrinsic water-use efficiency (iWUE), defined as the ratio between assimilation and stomatal conductance, is closely related to 13 C discrimination. Rising atmospheric CO 2 concentrations can have a fertilizing effect on plants, which in turn potentially increases iWUE (Keenan et al., 2013;Saurer et al., 2014) -at least as long as plant growth is not limited by other factors such as nitrogen or phosphorus supply (Reich et al., 2014;Yang et al., 2016) or water stress (Walker et al., 2015).
There is a rich body of literature on changes in 13 C isotopic discrimination and iWUE and on observational evidence from δ 13 C tree-ring records McCarroll and Loader, 2004;Fichtler et al., 2010;Leonardi et al., 2012;Churakova , Sidorova, S;Lévesque et al., 2014;Liu et al., 2014b;Saurer et al., 2014;Hartl-Meier et al., 2015;Voelker et al., 2016), free-air CO 2 enrichment (FACE)-type experiments (Battipaglia et al., 2013;Klein et al., 2016), δ 13 C site measurements (Pataki et al., 2003;Bowling et al., 2014), and δ 13 C paleodata (Voelker et al., 2016). These data generally suggest small to moderate decadal-to-century-scale changes in discrimination that correspond to a 20th century increase in iWUE On the order of 20 % and physiological control towards a constant ratio of the partial pressure of CO 2 within the leaf's substomatal cavities to the CO 2 pressure outside the leaf (c i /c a ; Saurer et al., 2004;Leonardi et al., 2012;Frank et al., 2015) and, more generally, a pattern of stomatal optimization towards minimizing water loss per unit carbon assimilated (Voelker et al., 2016).
The δ 13 C data are used to evaluate global and local models of plant growth, carbon cycling, and land-atmosphere isotopic fluxes for atmospheric carbon balancing (Scholze et al., 2003(Scholze et al., , 2008Suits et al., 2005;Danis et al., 2012;van der Velde et al., 2014). The δ 13 C data from leaf material  are used by Prentice et al. (2014) and Wang et al. (2016) to develop the representation of assimilation in land biosphere models following an optimization principle to balance carbon gain by assimilation and water loss.
The goal of this study is to present the implementation of δ 13 C in the land component, Community Land Model version 4.5 (CLM4.5), of CESM and Land Surface Processes and Exchanges (LPX-Bern) and to discuss the model performance for δ 13 C on the global scale. This is a step towards fully coupled isotope-enabled CESM applications and complements recent advances in simulating marine carbon isotopes. We compare the CLM4.5 and LPX-Bern results to a dataset of δ 13 C measurements on modern leaf material , a comprehensive compilation of century-scale δ 13 C tree-ring records, and results from the isotope-enabled LPX-Bern dynamic global vegetation model. We discuss spatial and century-scale trends in isotopic discrimination and iWUE of the two models in light of observational evidence, and the models' contrasting implementations of stomatal conductance and the balance between carbon assimilation and water loss.
2 Methods 2.1 CLM4.5 We use the CLM4.5 (Oleson et al., 2013), the land component of the Community Earth System Model version 1.2.0 (CESM1.2; Hurrell et al., 2013). The implementation of 13 C is outlined below (Sect. 2.1.1). A comprehensive description of the implementation of carbon isotopes in CLM4.5 is given in the technical description of Oleson et al. (2013), with further details elsewhere for 13 C Duarte et al., 2016) and 14 C . In addition to the land, carbon isotopes are also implemented in the ocean model of CESM1.2 (Jahn et al., 2015).
CLM4.5 features fully prognostic terrestrial carbon and nitrogen cycling which comprises all vegetation, litter, and soil organic matter pools (Oleson et al., 2013). Each grid cell is composed of multiple, independently represented land-use classes. Each class has its own set of plant functional types (PFTs), snow and soil columns.
Vegetation is comprised of 15 different PFTs, which are classified into three different phenological groups: evergreen, seasonal deciduous, and stress deciduous. A total of 14 of these PFTs follow the C3 photosynthetic pathway (11 tree, 2 grasses and crops) and 1 follows the C4 path (warm grasses). Altogether, 20 carbon and 19 nitrogen pools per PFT represent carbon (C) and nitrogen (N) in vegetation. C and N are tracked for leaf, live stem, dead stem, live coarse root, dead coarse root, and fine root pools and corresponding storage pools representing, respectively, short-term and long-term storage of nonstructural carbohydrates and labile nitrogen.
Decomposition of fresh litter material (including C and N) into progressively more recalcitrant forms of soil organic matter is represented as a cascade of transformations between decomposing coarse woody debris, litter, and soil organic matter pools. Depending on the C : N ratios of involved pools and the amount of carbon lost by respiration, each transformation can generate either a source or a sink of new mineral nitrogen.
Steps that result in an uptake of mineral nitrogen (e.g., immobilization fluxes) are subject to rate limitation, depending on the availability of mineral nitrogen and the sum of nitrogen demands from immobilization, photosynthesis, nitrification, and denitrification. If mineral N is less than the sum of these demands, fluxes are downregulated in order to match N supply. We note that the "Relative Demand" downregulation in CLM4.5 has recently been shown to be inaccurate in several tropical (Zhu et al., 2016b), tundra (Zhu et al., 2016a), and grassland (Zhu et al., 2017) systems. In addition to the cycling of nitrogen within the plant-litter-soil organic mat-ter system, CLM represents external sources, including atmospheric deposition and biological nitrogen fixation. CLM also represents other N sinks not included in this budgeting, including leaching and losses due to fire.
Photosynthesis in C3 and C4 plants is based on Farquhar et al. (1980) and Collatz et al. (1992), respectively. The maximum rate of carboxylation at 25 • C varies with an assumed static foliage nitrogen concentration and specific leaf area and is a PFT-specific parameter. It is assumed that leaf nitrogen and sunlight decrease exponentially with cumulative leaf area index from the canopy top to canopy bottom. Accordingly, the carboxylation rate and other photosynthesis parameters decrease exponentially within the canopy. Leaflevel photosynthesis is scaled to the canopy level by integration over the total leaf area. This is done separately for sunlit and shaded leaves and by considering the exponential scaling.
The allocation of carbon and nitrogen is determined in the following steps: first, gross primary productivity (GPP) is calculated under the assumption of unlimited nitrogen supply. Then, the maintenance respiration demand is subtracted from GPP. Following this, the actual nitrogen supply is compared against the nitrogen demand and GPP, if necessary, is downregulated. Finally, the available carbon is either utilized for plant growth and growth respiration or stored for growth in the subsequent years. Ghimire et al. (2016) recently presented an improved scheme which avoids this instantaneous downregulation by N limitation.
2.1.1 Carbon isotope discrimination during photosynthesis in CLM4.5 Isotopic ratios are usually reported as deviation from a standard material: where R sample and R std denote the 13 C / 12 C molar ratios of the sample and the standard material, respectively. Isotopic fractionation factors, α, are here defined as the ratio of the carbon isotope ratios in reactant to product (Farquhar et al., 1989). (We note that α is also defined by the ratio of product ratio to reactant ratio in the literature.) The fractionation factor for photosynthetic assimilation of CO 2 from canopy air is then as follows: where R air and R GPP denote the 13 C / 12 C molar ratios in the canopy air and in the resulting GPP flux incorporated in the plant material. A α value larger than unity results in a discrimination against the heavier isotope and therefore to a depletion of 13 C in GPP and plant material compared to air. Discrimination is also expressed by i , the deviation of α from unity and here multiplied by 1000 for conformity with the δ 13 C notation ( i = (α − 1) · 1000).
K. M. Keller et al.: Land carbon isotopes in CESM1 Photosynthesis in CLM4.5 and, embedded in this process, photosynthetic discrimination, are implemented in two steps following Farquhar et al. (1989). In step 1, given that no enzymatic fractionation is of relevance, the diffusion of CO 2 across the leaf boundary layer and into the stomata is associated with a kinetic isotope effect of a = 4.4. During step 2, enzymatic fixation, the effect on C3 plants is b = 27. These two steps are additive and result in the leaf-level fractionation factor; however, note that in the case of C4 plants only step 1 is of relevance. The CAM photosynthetic pathway is not considered in the model. The leaf-level fractionation factors (α psn ) of C3 and C4 plants are defined as follows: C4 plants: C3 plants: where c * i and c a represent the intercellular and atmospheric concentration of CO 2 (mol mol −1 ), respectively. This results in an isotopic discrimination between assimilated plant material (δ 13 C plant ) and atmospheric CO 2 (δ 13 c a ), expressed in per mill units (Farquhar et al., 1989) for C3 plants as follows: Isotopic discrimination can be approximated by i = a + (b −a)c * i /c a (Farquhar et al., 1989). Fractionation factors for all other land biosphere fluxes are set to unity and, thus, no further discrimination occurs in the land model.
The kinetic isotope effect during CO 2 fixation is constrained by c * i which, in turn, depends on the net carbon assimilation during photosynthesis (a n ; µmol m −2 s −1 ). The asterisk in c * i denotes the consideration of nitrogen downregulation in the photosynthesis calculation (for details, see Oleson et al., 2013), which is implemented as follows: where f dreg is the downscaling factor due to nitrogen limitation, g b is leaf boundary-layer conductance, and g s is leaf stomatal conductance for water (µmol m −2 s −1 ). In CLM4.5, assimilation is calculated before nitrogen limitation is considered. Thus, c i in the photosynthesis module is different from c * i . The value of c i follows from Eq. (6) with f dreg = 0. In other words, plants photosynthesize carbon at a potential, nitrogen-unlimited rate. Photosynthesis a n and, in turn, stomatal conductance g s are not directly affected by nitrogen limitation. The flow of water and carbon through the stomata is simulated to occur by molecular diffusion, and a ratio of 1.6 applies between the diffusivity of water and CO 2 , while a diffusivity ratio of 1.4 is assumed for the leaf boundary layer (Oleson et al., 2013). Stomatal conductance, g s , is itself linearly related to the product of leaf net photosynthesis (a n ), the relative humidity (RH) at the leaf surface, and the inverse of the CO 2 concentration at the leaf surface, c s (Ball et al., 1987;Collatz et al., 1991;Oleson et al., 2013): where b is a minimal conductance as a function of soil water availability, and m a constant parameter, defining the slope between a n and g s . For C 3 plants m is set to 9, and for C 4 plants it is set to 4 (Oleson et al., 2013). Both assimilation and stomatal conductance are downregulated by a drought stress factor which depends on soil moisture availability and can vary between 0 and 1. Nitrogen limitation and water stress can therefore result in a downregulation of c * i /c a and thus to a reduction in the carbon isotope discrimination of C3 plants in CLM4.5.
Equation (6) represents the constraint that CO 2 consumption by photosynthesis (A = a n · (1 − f dreg )) equals the transport of CO 2 through the leaf boundary layer and the stomata into the stomatal cavity. The latter is the product of the concentration gradient between the atmosphere and the stomatal cavity (c a − c * i ) and an overall conductance for CO 2 (g CO 2 = (g b g s )/(1.4g s + 1.6g b ) (µmol m −2 s −1 ). Thus, we can also write the following: In the following sections, for simplicity and for consistency with the isotopic literature and the LPX description, we will omit the superscript of c * i and denote c * i also by c i .

LPX-Bern 1.3
The LPX-Bern 1.3 model Stocker et al., 2013;Saurer et al., 2014;Ruosch et al., 2016;Keel et al., 2016) is based on the Lund-Potsdam-Jena (LPJ) dynamic global vegetation model (Sitch et al., 2003; see also Joos et al., 2001;Gerber et al., 2004;Strassmann et al., 2008). LPX combines process-based, large-scale representations of terrestrial vegetation dynamics, the dynamics of terrestrial carbon and nitrogen stocks and fluxes, and land-atmosphere exchanges of water, carbon dioxide, methane, nitrous oxide, and carbon and water isotopes in a modular framework. Following the assimilation of multiple experimental constraints, such as net primary productivity (NPP) (Olson et al., 2013) and the seasonal fraction of absorbed photosynthetically active radiation (Gobron et al., 2006), several key model parameters were updated. Note that no δ 13 C observational data, e.g., from tree-rings or atmospheric samples, were used as a constraint in the assimilation. A list of the updated parameters and their values can be found in Table 1.
Each grid cell in LPX is subdivided into different land-use classes (areas with natural mineral soils, peatland, other wet- Fraction of direct oxidation of leaf turnover on cropland 0.9 0.9094 lands, cropland, pasture, urban). LPX simulates the distribution of PFTs based on bioclimatic limits for plant growth and regeneration, and plant-specific parameters that govern plant competition for light, water (Sitch et al., 2003), and nitrogen (Xu-Ri and Prentice, 2008;Xu-Ri et al., 2012). Here, seven generic tree PFTs all follow the C3 photosynthetic pathway together with C3 grasses, and C4 grasses are considered on natural nonpeatland areas. The two PFTs on cropland and pastures have the same properties as the C3 and C4 grasses on natural land and grow depending on climatic conditions. On peatland, two PFTs (C3 graminoids and Sphagnum moss) grow. Seven carbon and nitrogen pools per PFT represent leaves, sapwood, heartwood, fine roots, aboveground leaf litter, aboveground woody litter, and belowground litter. Separate soil organic carbon and nitrogen pools receive input from litter of all PFTs. N input by biological N fixation is implied by maintaining prescribed C to N ratios in litter and soil pools. Inorganic soil nitrogen pools receive input from litter and soil decomposition and through atmospheric N deposition and are subject to plant uptake, leakage, and gaseous losses (Xu-Ri and Prentice, 2008). Photosynthesis and stomatal control in LPX is described by Haxeltine and Prentice (1996b), as summarized elsewhere (Keel et al., 2016;Sitch et al., 2003). The equations for water supply from soil and transpiration, assimilation, and canopy conductance are solved simultaneously by varying the ratio c i /c a , also termed λ, to yield self-consistent values for these properties.
Total daytime net photosynthesis, A dt , is modeled following Collatz et al. (1991Collatz et al. ( , 1992, which is a Farquhar model (Farquhar et al., 1980) generalized for global modeling purposes (for details, see Haxeltine and Prentice, 1996a). A dt is a function of daily radiation, temperature, day length, and atmospheric CO 2 partial pressure and λ. A dt is computed from a formulation which gives a gradual transition between light-limited and Rubisco-limited rate of photosynthesis. The amount of photosynthetically active radiation absorbed by the entire canopy increases with the modeled leaf area index following Beer's law and is used to compute the light-limited photosynthesis rate. The N content and Rubisco activity of leaves are assumed to vary seasonally and with canopy position in such a way as to maximize A dt .
Canopy conductance, g c , is linked to A dt through where g min is a PFT-specific minimum canopy conductance. Daily transpiration of water is calculated for each PFT as the minimum of a plant-and soil-limited supply function (E supply ) and the demand for transpiration (E demand ). E supply is the product of root-weighted soil moisture availability and a maximum water supply rate that is equal for all PFTs. E demand is calculated following the empirical relation between evaporation efficiency and surface conductance of Monteith (1995): where E eq is the equilibrium evaporation rate, dependent on temperature and radiation, g m and α m are empirical parameters, g c is the canopy conductance, and φ is the fraction of present foliage area to ground area (i.e., projected leaf area). Equation (10) is solved for E demand using the non-waterstressed potential canopy conductance which is calculated using Eq. (9) and a fixed ratio λ to compute A dt . The value of λ is set equal to 0.8 following Sitch et al. (2003) to approximate non-water-stressed conditions and as a starting value for the iterative computation of carbon assimilation and transpiration under water shortage. In the case of water-stressed conditions when E demand exceeds E supply , canopy conductance and photosynthesis are downregulated; E demand is set to E supply and Eq. (10) is solved for g c . Knowing g c and c a , λ is varied in the photosynthesis module until the following relationship is satisfied: NPP is downregulated on the daily model time step if N demand exceeds N availability from the inorganic soil nitrogen pools. Daily NPP is integrated over a year and allocated annually to vegetation. Thus, in contrast to CESM, there is no immediate feedback of nitrogen limitation on isotopic discrimination on a daily timescale, but there is a long-term feedback by annual changes in vegetation structure and, in turn, photosynthesis and carbon assimilation.

Carbon isotope discrimination during photosynthesis in LPX
The isotopic discrimination during assimilation is calculated on a daily time step following Farquhar et al. (1989). Thus, the same discrimination formulations as in CESM (Eqs. 3 and 4) are used in LPX. We recall, however, that the intercellular CO 2 concentration is, unlike in CLM4.5, not downregulated by nitrogen limitation in LPX. For Sphagnum mosses, discrimination during assimilation is fixed to −30 ‰. As in CESM, no further discrimination is assumed in LPX. In earlier LPX applications, δ 13 C was implemented following Scholze et al. (2003) and discrimination modeled following Lloyd and Farquhar (1994). Here, we adjusted the discrimination formulations and used instead the simpler formulations of Farquhar et al. (1989) for consistency with CLM4.5 and the computation of iWUE as described below. We note that our conclusions do not depend on this choice. LPX yields similar 20th century changes in discrimination and iWUE for both formulations. However, simulated δ 13 C of carbon assimilated by C3 trees is generally less negative (by about 2 ‰) when applying the Lloyd and Farquhar (1994) formulation, and agreement with leaf δ 13 C is less favorable.

Spin-up and transient simulations and model forcings
CLM: Starting with empty pools, the CLM model is brought into equilibrium with the following steps: (1) accelerated decomposition of soil organic matter with perpetual 1850 Common Era (CE) forcing (1000 model years), (2) normal decomposition with perpetual 1850 CE forcing (500 model years), (3) transition to 1900 CE conditions (100 model years), and (4) transient simulation over the 20th century . During steps (1) and (2), with both default model options, atmospheric δ 13 C is held on a constant value of −6 ‰. The first step, based on the accelerated decomposition technique by Thornton and Rosenbloom (2005), accelerates the equilibration of the soil carbon pools. During step (3) the model is adjusted to 1900 CE conditions by prescribing transient atmospheric pCO 2 and δ 13 C (years 1800-1900) together with CRU NCEP climate forcing data (repeated years 1901-1920Viovy, 2011).
Step (4), the 20th century simulation, is run with both transient atmospheric pCO 2 and δ 13 C and climatic forcing data (years 1901Viovy, 2011). In both steps (3) and (4), atmospheric CO 2 is prescribed following Law Dome data (Etheridge et al., 1996;MacFarling Meure et al., 2006) and atmospheric samples; prescribed atmospheric δ 13 C is a combination of data from Rubino et al. (2013) and, from 1993 on, White et al. (2015) ( Fig. 1). Land-use area and change is prescribed following Hurtt et al. (2006). LPX: The LPX model is forced with CRU TS3.23 climate data (Harris et al., 2014), and the same atmospheric CO 2 and δ 13 C data is prescribed as in CLM4.5. Land-use maps are prescribed following Hurtt et al. (2006). After the spin-up procedure, the model is adjusted to 1900 CE conditions by applying transient atmospheric pCO 2 and δ 13 C (years 1765-1900) together with recycled climate data (years 1901-1930). Then, the 20th century is simulated with transient atmospheric CO 2 and δ 13 C and climatic forcing data (years 1901-2007).
Factorial runs were performed in addition to the standard simulation. Each setup is the same as for the 20th century simulations, although with one driving factor held constant: (i) cCLIM, "constant climate" (repeated years 1901-1920 for CLM4.5 and 1901-1930 for LPX); (ii) cCO 2 , constant atmospheric CO 2 concentrations, only done with LPX; (iii) cN-DEP, constant atmospheric nitrogen deposition, only done with LPX; and (iv) cLU, constant land use, only done with LPX. The difference between the factorial simulation and standard simulation is attributed to the driving factor that was held constant in the factorial run. An interaction term is computed from the difference between the change simulated in the standard simulations and the sum of all attributed changes.
The two models are forced with two slightly different climate datasets as CLM4.5 is run on a subhourly time step, while a daily time step is used in LPX. CRU NCEP, used to force CLM4.5, is a combination of the CRU TS3.2 monthly climatology (resolution of 0.5 • × 0.5 • ) and the NCEP reanalysis product with a time step of 6 h (2.5 • × 2.5 • ). The results of the factorial runs, presented later, suggest a generally small influence of this difference in forcings on simulated discrimination and iWUE. We further investigated the potential impact of using two different climate data products. The CRU NCEP data are interpolated on the 1 • × 1 • LPX grid and integrated to monthly values for use in LPX. LPX is run with the modified CRU NCEP and, as usual, with CRU-TS3.22 in the standard model setup where climate is changing transiently and the factorial setup with cCLIM. Differ- ence in simulated changes in discrimination and iWUE are small for both cases. In conclusion, the large model differences, presented in the result section, are not caused by differences in climate input data.

Observations
Two different observational datasets are used to evaluate the model simulations. The first dataset is a compilation of δ 13 C time series measured on tree rings (for an overview, see Table A1). The time series cover at least the years 1900-1985 (and 2005 if possible) and comprise a wide range of different tree species and locations. Since absolute values are highly variable, depending on factors such as, for example, species, tree age, or location (McCarroll and Loader, 2004), we refrain from a direct comparison of absolute values and focus on the changes over the course of the century. The second dataset is a compilation of δ 13 C measurements on leaves . The dataset is comprised of data from C3, C4, and CAM plants; however, to ensure comparability between models and observations, only data of C3 plants are used.
The observational δ 13 C data are compared with 10-year averages of available model output to remove interannual variability from the model data. For CESM, data are compared to the simulated grid-cell average of the δ 13 C signature of the "live stem" pool. This pool has a very fast turnover time (0.7 per year; Oleson et al., 2013) and no post photosynthetic fractionation occurs in the model during the allocation of assimilated carbon. Modeled δ 13 C in the live stem pool thus represents on average the δ 13 C signature of the leaves of C3 trees. For LPX the δ 13 C signature of daily GPP, GPP (t,PFT) , is averaged for all tree PFTs within a grid cell and the land-use class "mineral soil": t is time and the sum is over 10 years. PFT is the index for PFTs and the sum is over all tree PFTs (without including grasses).

K. M. Keller et al.: Land carbon isotopes in CESM1
The evolution of δ 13 C in (C3) trees is a consequence of both the changes in atmospheric δ 13 C as well as changes in fractionation, related to changes in c i /c a (Eq. 4) and thus to a combination of physiological adjustments of trees and ecosystems to CO 2 , climate, and other environmental changes. The atmospheric δ 13 CO 2 record shows a centuryscale decrease in δ 13 C caused by the input of isotopically depleted CO 2 from deforestation and fossil-fuel burning. This century-scale decrease, known as the Suess effect (e.g., Tans et al., 1979), is a precisely known external forcing. To account for this effect, we focus on the isotopic discrimination i between assimilated plant material and atmospheric CO 2 (Eq. 5).
Changes in discrimination of C3 plants are related to the ratio of photosynthesis to the conductivity of carbon through the stomata, two relevant tree physiological parameters regulating carbon and water fluxes. The ratio of net daytime photosynthesis to stomatal conductance for water vapor (A/g H 2 O ) is also known as iWUE and discussed in the literature (e.g., Scheidegger et al., 2000;Saurer et al., 2004Saurer et al., , 2014. Using Eqs. (5) and (6) for CLM4.5, Eq. (11) for LPX, the definition of iWUE, and a ratio of 1.6 between the conductance of water and CO 2 , we obtain the following: In the following, the symbol is used to denote a temporal change (not to be confused with the symbol for discrimination i ). For example, iWUE is the difference between iWUE at time t 2 and time t 1 . We approximate iWUE to better understand how iWUE depends on discrimination, changes in discrimination, and atmospheric CO 2 . The CO 2 concentration c a (t 2 ) is expressed as the sum of c a (t 1 ) and the change in CO 2 is expressed as c a . In analogy, we write d(t 2 )=d(d 1 )+ d for the difference between the atmospheric and plant δ 13 C (d = δ 13 c a − δ 13 C plant ). The term 1+ δ 13 C plant 1000 is denoted as k. Variations in k are of order 1 ‰ and have a small influence on iWUE. Then, the following applies for the difference in iWUE ( iWUE): where q denotes the term 1.6(b − a) and is, with b = 27 and a = 4.4, equal to 36.16. The approximation given in Eq. (14) holds well within 0.1 %. Equation (14) shows that iWUE is linearly related to the isotopic discrimination at time t 1 , i (t 1 ) and the change in discrimination between t 1 and t 2 , i . In other words, the change in iWUE depends both on the change in discrimination and on the initial magnitude of the discrimination. The sensitivity in iWUE is about 7 times larger for a unit change in i compared to a unit change in i (t 1 ) when assuming a change in CO 2 from 300 to 350 ppm as reconstructed for the period 1900 to 1980. Under increasing CO 2 , discrimination and thus iWUE may change. The larger the decrease in discrimination, the larger the increase in iWUE, and vice versa.
The relative change in iWUE may be expressed, again very accurately, as follows: The relative change in iWUE scales in proportion to the relative change in atmospheric CO 2 in the absence of a change in discrimination. It also scales with changes in discrimination. For typical values, a unit change in discrimination changes the relative change in iWUE by 50 to 100 %. Further, the relative change in iWUE depends on the difference between b(= 27) and the isotopic discrimination i (t 1 )(≈ 20). This difference is typically 7 and so a unit change in i (t 1 ) affects the second term in Eq. (15) by about 14 % and the relative change in iWUE typically between 0 and 10 %. We note that Eq. (15) is nonlinear and so the numerical values given here are illustrative.
Equation (13) is evaluated for LPX and CESM by using the annual-mean δ 13 C signature of C3 plants for CESM (live stem) and LPX (GPP) within each pixel. We note that in CLM4.5 leaf boundary conductance is also explicitly modeled, and g H 2 O in Eq. (13) thus reflects the ratio of photosynthesis to the total conductivity through the leaf boundary layer and the stomata of water, whereas in LPX g H 2 O reflects the resistance by the stomata only. In CLM4.5, the ratio between the conductivity for H 2 O and CO 2 may vary between 1.4(g s g b ) and 1.6(g b g s ) (Eq. 6). Here, we have assumed for simplicity a factor of 1.6 for CLM4.5, LPX, and observational data. This simplification introduces an uncertainty of up to 11 % in computed trends for CLM4.5.
Turning to the tree-ring data, there are isotopic offsets between assimilated material and cellulose. The leaflevel model for fractionation is more representative for the bulk matter rather than a specific chemical compound such as cellulose. We correct for this isotope offset between cellulose and total organic matter (δ 13 C plant-corrected = δ 13 C plant cellulose − offset) following Saurer et al. (2014) and using an offset of 1.2‰ for all tree species. This correction has a relatively small influence: as discussed above, a change in δ 13 C plant by 1 ‰ affects the relative change in iWUE by about 10 %.
In the following, we will present and discuss annual or multiannual averages of the related variables δ 13 C, c i /c a , discrimination i , iWUE, and A/g. These values represent, similarly to Eq. (12), weighted averages. The time-varying carbon assimilation rate, or more generally, the assimilation  (mean 1996-2005) and change over the 20th century for both models and an observational dataset (Carvalhais et al., 2014 3 Results

Primary productivity and carbon pools
We first compare GPP, transpiration, vegetation, and soil carbon stocks of the two models (Table 2, Fig. 2), as well as model results to observation-based estimates of vegetation carbon ( Fig. 3; Carvalhais et al., 2014). Both models show the highest GPP per unit land area in the tropics and low GPP in arid regions. Overall, the GPP patterns are similar between the two models. However, the tropical maxima in GPP are about a third larger in CLM4.5 than in LPX. Simulated transpiration is similar for the two models.
Estimates of total vegetation carbon, including below and aboveground biomass, were derived by Carvalhais et al. (2014) from a collection of estimates for pantropical regions and for northern and temperate forests based on radar remote-sensing retrievals. The global carbon inventory of vegetation is overestimated by a factor of 2 by CLM4.5 and slightly underestimated by LPX compared to the estimate of Carvalhais et al. (2014, Table 2). CLM4.5 simulates too much vegetation carbon in northern mid-to-high latitudes and overestimates vegetation carbon stock by more than a factor of 2 in tropical forests (Fig. 3) (Negrón-Juárez et al., 2015;Koven et al., 2015a). Correspondingly, we expect a negative bias in the global-mean δ 13 C signature of vegetation in CLM4.5. LPX simulates vegetation carbon stocks in relatively good agreement with observation-based estimates (Fig. 3). Both models underestimate vegetation carbon in parts of southeastern Asia. Overall, the spatial correlation (r) between modeled and observation-based vegetation carbon is 0.83 for CLM4.5 and 0.85 for LPX.
The global carbon inventory in soils of LPX is 1700 GtC, while CLM4.5 simulates a soil C inventory that is twice as large at 3900 GtC (Table 2). This larger inventory mainly stems from the large carbon stocks simulated by CLM4.5 in peatland and permafrost regions (Koven et al., , 2015b and in northern mid-latitude regions (Fig. 2). CLM4.5 includes formulation for carbon storage in deep soils. The large soil carbon inventory of CLM4.5 is in contrast to that of its predecessor, CLM4 (Oleson et al., 2010), which severely underestimated the global soil carbon reservoir by 50 % and more (Todd-Brown et al., 2013;Anav et al., 2013;Tian et al., 2015), especially at the northern high latitudes (Foereid et al., 2014). LPX, in the version applied here, does not include formulations for carbon storage below 2 m. An evaluation of peat and permafrost carbon stocks simulated by LPX is provided by Spahni et al. (2013). Data on soil carbon have particularly high uncertainties, with notable discrepancies between different datasets (e.g., Carvalhais et al., 2014;Hugelius et al., 2014).
Simulated GPP and carbon stocks change over the industrial period in response to increasing atmospheric CO 2 , land use, climate change, and nitrogen deposition. Both models indicate an increase in GPP almost everywhere over the course of the 20th century (Fig. 2). Globally averaged GPP increases from 148 GtC yr −1 in 1900 to 168 GtC yr −1 in 2000 for CLM4.5 and from 132 to 149 GtC for LPX. In the case of CLM4.5 the increase is especially pronounced in the tropics, particularly in Indonesia (Fig. 2), whereas LPX simulates a large increase in GPP in Australia.
Simulated 20th century changes in annual-mean transpiration are generally small in both models, except in Australia for LPX and in parts of Latin America for CLM4.5. Generally, an increase in transpiration is found in boreal and temperate forest regions in both models. Transpiration is slightly decreasing in LPX and slightly increasing in CLM4.5 in most tropical forest regions.
For vegetation carbon, both models show large areas with decreasing C storage in the mid-latitudes (Europe, North America) and parts of the tropics, especially Southeast Asia, in response to deforestation (Fig. 2). Areas with a negative vegetation C balance are more widespread in LPX than CLM4.5 in tropical and subtropical Africa and South America. LPX simulates a relatively small increase in the vegetation stocks of remaining "natural" tropical forests, whereas CLM4.5 shows a considerable increase in vegetation C in the remaining tropical forests of Africa, America, and Indonesia. This large increase in CLM4.5 is related to the strong stimulation of GPP and the long overturning timescale of vegetation carbon in tropical forests. As a result, the decrease in globally integrated vegetation C storage is more than 2 times larger for LPX than for CLM4.5 (−55 to −21 GtC; Table 2, Fig. 1b).
For soil C (Fig. 2), CLM4.5 shows a general increase over the 20th century which is especially pronounced in the tropics. For LPX the picture is more heterogeneous, with increases (decreases) in South America, Africa, Australia, and the high latitudes (India and parts of the US and Eurasia). Globally, the models simulate an increase of 45 (CLM) and 19 GtC (LPX), respectively (Table 2, Fig. 1). Overall, global carbon storage in the land biosphere increases by 34 GtC in CLM4.5 and decreases by 35 GtC in LPX over the 20th century. Independent estimates of the difference between anthropogenic emissions and ocean uptake plus an increase in atmospheric CO 2 suggest an increase in land carbon storage by roughly 20 GtC (Le Quéré et al., 2015) over the 20th century.
3.2 Simulated δ 13 C of primary productivity, vegetation, and soil carbon: spatial distributions and 20th century trends Global maps of the δ 13 C signature of GPP, vegetation, and soil carbon in the year 2000 and related 20th century changes are shown in Fig. 4. We note that isotopic fractionation, typically on the order of 20 ‰ for C3 plants, is only 4.4 ‰ for C4 plants, growing in high-temperature, high-light-intensity regions. The grid-cell average δ 13 C signature of GPP therefore reflects to a large extent the occurrence or absence of C4 plants. Vegetation distribution is prescribed in CLM4.5 and simulated in LPX, and CAM plants are not represented in either model. Both models show less negative isotopic signatures in arid and semi-arid low-latitude regions. In comparison to CLM4.5, LPX simulates in the year 2000 much larger areas with δ 13 C larger than −16 ‰. These areas include central and western parts of the US, Argentina, large parts of southern Asia, and parts of southern Australia. This is related to a larger share of C4 plants in these regions in LPX than in CLM4.5.
Atmospheric δ 13 C decreases by about 1.3 ‰ over the 20th century (Fig. 1a). Thus, without a change in isotopic discrimination, we would expect an equal decrease in the GPP signature. Accordingly, CLM4.5 and LPX show a decrease in δ 13 C of GPP in most areas (Fig. 4).
In LPX, large 20th century increases in the δ 13 C signature of GPP of up to 5 ‰ are simulated in many low latitude and mid-latitude areas, indicating an expansion of C4 plants under warming and rising CO 2 as well as the influence of land use. Processes such as the conversion of forests into C4 pastures in the tropics can substantially affect the isotopic signature of total terrestrial vegetation at a location (Townsend et al., 2002;Kaplan et al., 2002). In CLM4.5, the influence of the negative atmospheric δ 13 C trend is offset mainly by changes in discrimination in the eastern US, Europe, and eastern Asia, resulting in slightly positive δ 13 C changes in these areas (Fig. 4). Note that, in contrast to LPX, crops are treated as C3 plants in CLM4.5 (Oleson et al., 2013). This explains the land-use-related, differing trends between models in regions such as North and South America, Africa, and Southeast Asia. Distribution and 20th century change in δ 13 C of vegetation and soil are very similar to the distribution of δ 13 C of GPP. However, temporal trends are generally smaller in the soil carbon pools than for GPP.
The average δ 13 C signature of the globally averaged pool of vegetation carbon is more negative than that of soil and leaf carbon (Fig. 1). This reflects the share of organic matter assimilated by the C4 path versus that by the C3 path in these globally averaged pools. Generally, the share of C4-derived organic matter is much smaller than the contribution by C3derived matter. Correspondingly, the average δ 13 C signature of globally averaged C pools is strongly negative and closer to that of C3 plants. The contribution of C4-derived organic matter is very small for the globally averaged vegetation pool . The estimates are based on decadal means (1896-1905 and 1996-2005).
and absent for the stem carbon pool, but noticeable for the leaf and soil carbon pools. Further, the globally averaged leaf δ 13 C signature decreases in both models over the 20th century, mainly in response to changes in the atmospheric δ 13 C source. Trends are less developed in the other global pools due to their longer residence times. There is a substantial offset between the models for globally averaged soil and vegetation pools, with CLM4.5 being around 5 ‰ more negative (Table 2). In addition to the already mentioned widespread, comparably less negative δ 13 C signatures simulated by LPX, this can be attributed to the spatial distribution of carbon pools in CLM4.5. As is evident from Fig. 3, CLM4.5 simulates higher vegetation and soil carbon stocks in the high latitudes and tropical forest regions than LPX. Since these regions are characterized by C3 plant cover and thus comparably negative signatures (see Fig. 4), the result is a shift of the global averages towards more negative δ 13 C in CLM4.5 compared to LPX.

Modeled versus measured leaf δ 13 C
In this section, we evaluate how well modeled δ 13 C, and thus isotopic discrimination, of C3 tree species compares with a global compilation of δ 13 C measurements on leaf material from C3 plants Cornwell et al. (2016, Figs. 5, 6). Unlike the previous analyses, this comparison is not (or is hardly) affected by land use as the comparison is only for C3 trees on unmanaged lands. However, a few caveats apply. We compare grid-cell averages (with GPP or mass as weights; see Sect. 2.4) from all simulated C3 plants with site measurements for individual species. Difference in δ 13 C are reported for different species. In addition, differences within a species, even growing at the same site, can be as large or even larger as those between species (McCarroll and Loader, 2004;, and references therein). In addition, δ 13 C in the canopy air may deviate from the prescribed atmospheric δ 13 C, representative of tropospheric background air. Deviations may arise due to varying additions of isotopically depleted respired carbon and of carbon from fossil sources to local air. The bar plots in Figs. 5 and 6 indicate the mean δ 13 C of the leaf data and of modeled grid-cell averages for C3 plants at the measurement locations for observations, CLM4.5, and LPX, respectively. On average, δ 13 C of all 344 leaf samples is −27.54 ‰ compared to a corresponding average of −26.29 ‰ for CLM4.5 and −26.14 ‰ for LPX. This suggests a discrimination of 19.5, 18.3 and 18.1 ‰ relative to the atmospheric δ 13 C signature of approx. −8 ‰, respectively. According to Eq. (4), these values correspond to a ratio between internal and atmospheric partial pressure, c i /c a , of 0.67 for the leaf data and of 0.63 and 0.61 for CLM4.5 and LPX, respectively. The bias of the models is particularly large in arid regions and larger for LPX than for CLM4.5, whereas good agreement between models and leaf data is found in the remaining regions, e.g., Europe or Alaska (Fig. 6).
The δ 13 C of C3 plants from CLM4.5 and LPX show similar spatial gradients. Both models reasonably capture the observations, though spatial correlations between data and models are low (CLM: r = 0.36; LPX: r = 0.34). Root mean square errors (RMSEs) between models and data are 2.42‰ (CLM) and 3.22 ‰ (LPX). The largest discrimination and most negative δ 13 C is simulated in high northern latitude regions and in tropical moist forests, with the lowest discrimination in arid regions of inner Asia, southwestern North America, Patagonia, southern Africa, and parts of Australia. Low discrimination and high (less negative) δ 13 C correspond to low stomatal conductance, low internal leaf CO 2 partial pressure and therefore reduced assimilation. However, these high δ 13 C values are not supported by the leaf measurements, when assuming that the leaf samples are regionally representative (Fig. 6). Maxima in δ 13 C are a few per mill higher for LPX than for CLM4.5 in these regions. This suggests a stronger downregulation of stomatal conductance by water stress in LPX than in CLM4.5 in these regions.
In summary, δ 13 C and, by implication, discrimination of C3 trees is reasonably represented in both models, when considering the caveats discussed at the beginning of this section and the simplicity of the isotopic model approach. However, both models tend to overestimate δ 13 C in comparison with the leaf data in many arid regions.

δ 13 C -century-scale trend
Next, simulated century-scale trends in the discrimination, i (Eq. 5, see also methods), of C3 tree species are discussed and compared to the corresponding trends derived from 76  tree-ring δ 13 C records (Figs. 7 to 10 and Table A1). Again, grid-cell averages for C3 trees are compared to site data. Similar caveats to those discussed for the comparison between leaf measurements and modeled δ 13 C apply for this trend analysis. The term "average" refers again to an average over only those grid cells with measurements. Out of 76 sites, 74 are located in three regions (Fig. 9), with more than half of the sites (51) in Europe, and with only 7 records in Asia and 16 records in western North America. For these two regions, the comparison of regional model and data averages is therefore hampered by the scarcity of site data and the complex topography along the Rocky Mountains, while the "global" average is biased towards Europe.

CLM4.5 LPX-Bern
On average, discrimination inferred from the tree-ring records varies within a few tens of a per mill in the first and second half of the 20th century (Fig. 7). There is a transition to less discrimination in the 1940s, with no clear trend before and afterwards. Similar to the measurements, LPX-Bern yields no long-term trend in the discrimination of C3 trees. In contrast, discrimination of C3 trees simulated by CLM4.5 shows a trend towards smaller values over the 20th century. The discrepancy between tree-ring data and CLM4.5 results grows towards the present and is larger than 1 ‰ after year 2000.
Next, we analyze spatial and regional changes in the discrimination of C3 trees by comparing the differences in i between the decades 1980and 1900-1909.5 and iWUE (%) as measured on tree rings (black) and modeled by CLM4.5 (blue) and LPX (red). Observations are represented by the average over all measurements (see Table A1 simulates a decreasing trend in the discrimination i of C3 trees almost everywhere where C3 trees grow. In other words, the ratio of internal to atmospheric CO 2 partial pressure, c i /c a , decreases over the 20th century. The simulated decrease in i for CLM4.5 is larger than suggested by the tree-ring data. The average change in discrimination ((1980-1989) − (1900-1909)) is −0.36 in the tree-ring data and −0.94 in CLM4.5. The observations show on average a slight positive change in i (+0.18) in the North American region, a slight negative change (−0.13) in the Asian region, and a relatively large decrease in i in Europe (−0.54) (Fig. 9). The three corresponding regional averages for CLM4.5 are more negative by between 0.52 and 0.79. LPX simulates generally small (positive or negative) changes in i of C3 trees over the last century. Exceptions are areas at the margin of C3 tree-covered regions such as the southern limit of the boreal tree belt that show a large positive change (Fig. 8). For grid cells with measurements, the average change in i for LPX is very small with +0.06 (Fig. 10). For the three selected regions, LPX simulates changes in i that are between 0.45 to 0.54 more positive than suggested by the tree-ring records (Fig. 9). We note that the agreement between tree-ring data and LPX is better for previous and later decades (Fig. 7). The correlations r between models and observationderived changes in i are relatively low, with values of 0.40 (CLM) and 0.27 (LPX). Here, possible explanations are (i) the strong spatial heterogeneity of both measurements and model simulations combined with a limited number of obser-vational data and the overall small changes in discrimination of, for example, only −0.36 for the average of tree-ring data or (ii) that the models do not properly represent photosynthesis and stomatal conductance. The RMSEs of changes in i are 0.93 (CLM) and 1.07 (LPX).
The factorial simulations with LPX reveal that average changes in i attributable to individual drivers are small. Thus, the relatively small changes in discrimination simulated by LPX are not the result of offsetting influences of different drivers. In detail, a somewhat complex interaction of the driving factors in shaping the change in i and nonnegligible interaction between climate and CO 2 is found. Keeping nitrogen deposition or land use constant generally has a negligible influence on the change in i (Fig. 10). On average, a positive change in i is attributed to the CO 2 change for all sites (+0.03), the European sites (+0.05), and the Asian sites (+0.09) and a negative change for the North American sites (−0.18). Thus, the increase in CO 2 does not cause a general downregulation of c i /c a in LPX. The influence of climate is small on average over all sites (+0.07) and for the European sites (−0.08). At the North American sites and the few Asian sites, climate change tends to increase discrimination by about 0.7 and 0.4, respectively.
For CLM4.5, the change in i attributed to climate change is small, except for the North American sites. There, a positive influence (+0.55) similar to that found for LPX is inferred. The other drivers together cause on average a change in discrimination between −0.79 and −1.16 in the three re-CLM4.5   Figure 9. Same as Fig. 8 but for three selected regions. The bar plots show the regional average change in discrimination of C3 trees as calculated from tree-ring δ 13 C data (gray) and from results of the standard simulations of CLM4.5 (filled blue) and LPX (filled red) as well as from factorial runs (pattern). Individual driving factors (climate, CO 2 , N-deposition, and land use) were kept constant in the factorial runs, as explained in the main text and indicated by the legend. Model averages are computed from the mean grid-cell discrimination of C3 trees and for all model grid cells where tree-ring estimates are available. The number of available tree-ring records is indicated.

LPX-Bern
gions. This suggests that the simulated decrease in discrimination and in c i /c a in the CLM4.5 runs is mainly linked to increasing CO 2 and a corresponding downregulation of c i /c a (Eq. 6). The downregulation in CLM4.5 is larger than suggested by the tree-ring records. moisture simulated by CLM and LPX are also small in large parts of Europe and Asia, where most of our tree-ring data are located. Despite these small changes in soil moisture, large differences in discrimination changes are found between the two models in these simulations and regions (Fig. 9). This suggests that the primary reason for the modelmodel difference in simulated 20th century changes in discrimination is rooted in the different parameterizations of the photosynthesis-conductance coupling.

δ 13 C -century-scale trend in iWUE
Of particular importance are changes in stomatal conductance and photosynthetic carbon assimilation since they are the primary controls of plant-atmosphere water and CO 2 exchange, carbon assimilation, and, ultimately, tree growth (Lambers et al., 2008). Changes in the ratio of these two processes, also termed iWUE (see Eq. 13), thus reflect a change in assimilation, in conductance, or in both (Scheidegger et al., 2000;Saurer et al., 2004Saurer et al., , 2014. Similar to in the previous sections, we compare relative changes in iWUE (Eq. 15) from the tree-ring records with model results (Figs. 10, 11, 12). Changes are evaluated between the decade 1990and decade 1900-1909.
As shown by Eq. (15), the relative change in iWUE increases equally to the relative change in atmospheric CO 2 in the absence of changes in discrimination. This is equal to 23 % for the change from 1900-1909 to 1990-1999. A negative change in discrimination, as simulated for CLM4.5, contributes to a further increase in iWUE.
In the standard setup (see also maps), both models show a substantial increase in iWUE over the century of approx. 52 % (CLM) and 25 % (LPX). The globally averaged iWUE derived from δ 13 C measurements show a value of 28 %, which concurs with previous studies reporting increases in iWUE of 27.8 and approx. 30 % for European forests Frank et al., 2015) and 17 % for high-latitude forests (Trahan and Schubert, 2016). While LPX is in good agreement with the available observations (both globally and regionally), CLM4.5 appears to overestimate 20th century changes in iWUE by about a factor of 2 in Europe and the Asian region.
The drivers behind modeled changes in iWUE can be investigated based on the sensitivity experiments already discussed in the previous section for discrimination. As expected, all simulations show similar changes in iWUE in LPX. The small changes in discrimination simulated by LPX do not strongly affect iWUE, and the relative change in iWUE of LPX (24.6 %) is close to the relative change in atmospheric CO 2 (23 %). Simulations with and without climate change with CLM4.5 yield similar results, pointing again to a downregulation of c i /c a under rising CO 2 by the implemented nitrogen limitation.

Discussion
This study presents the implementation of δ 13 C in two global land biosphere models, CLM4.5 and LPX-Bern 1.3, and a global compilation of 20th century tree-ring δ 13 C records. Results from isotope-enabled simulations over the industrial period are investigated for regional carbon stocks and carbon isotope signatures (Figs. 1 to 4) and for the global carbon and 13 C budgets (Table 2) of the land biosphere. The performance of the two models in terms of modern isotopic discrimination of C3 plants is evaluated by comparing model results with a global compilation of δ 13 C measurements on leaf material (Figs. 5, 6). The main focus of the study is on the analysis of 20th century changes in isotopic discrimination, i , in the ratio of CO 2 concentrations between the stomatal cavity and the atmosphere, c i /c a , and in the ratio between assimilation and stomatal conductance, A/g, termed intrinsic water-use efficiency iWUE. The comparison of model results with observational estimates from the δ 13 C tree-ring records permits us to evaluate the performance of the two models with respect to linkages between photosynthesis and stomatal conductance on the leaf level (Figs. 7 to 12).
The "minimal" discrimination model of Farquhar et al. (1989) was implemented in CLM4.5 and LPX-Bern. This model assumes fixed discrimination for C4 plants and that discrimination by C3 plants is directly proportional to c i /c a . Genetic species-specific variations in δ 13 C are not considered here, though they can be considerable (Marshall et al., CLM4.5 2008; Yang et al., 2015). The influence of rooting depth, water transport systems, root-to-leaf distance, leaf morphology, and irradiance (sunlit versus shaded) on discrimination is neglected. However, these factors may affect discrimination. Reported isotopic differences within a species, even between those growing at the same site, can be as large as or even larger than those between species (McCarroll and Loader, 2004;. We also neglect fractionation between different pools within plants and soils (Wingate et al., 2010;Brüggemann et al., 2011). Despite these limitations, we provide a comprehensive framework for a model-data comparison and derive new insights into patterns of physiological plant adaptation to 20th century environmental changes.

Isotopic signatures, pools, and fluxes of the global land biosphere
A high spatial correlation of simulated vegetation carbon with observational estimates is found for both models. However, CLM4.5 overestimates carbon stocks in vegetation by a factor of 2 or more in many regions compared to the data of Carvalhais et al. (2014), while data-model agreement is reasonable for LPX (Fig. 3). Global soil carbon stock is 1700 GtC in LPX and much smaller than the 3900 GtC simulated by CLM4.5. Larger soil carbon inventories in CLM4.5 than LPX are mainly simulated in northern mid-latitudes and high latitudes (Fig. 2). Only the top 2 m of soils are considered in LPX, while CLM4.5 includes deeper soil layers.
Both models feature most negative δ 13 C in high-latitude and tropical-forest regions (Fig. 4) and isotopically heavier signatures in semi-arid and arid regions, where C4 plants are typically abundant. The global average assimilationweighted δ 13 C signature of GPP is −20.6 ‰ in LPX and −24.9 ‰ in CLM4.5 in the year 2000. The associated discrimination is 12.4 ‰ in LPX and 16.7 ‰ in CLM4.5. The CLM4.5 estimate is within the range (15.7 to 18.1 ‰) of previous studies, as summarized in Suits et al. (2005) and Scholze et al. (2008). This may point to a too-large abundance of C4 plants in LPX. Differences in soil and vegetation carbon and in C3 and C4 plant abundance also affect the global isotopic budgets of the two models, as the relative importance of C4-derived plant and soil material is less in CLM4.5 than in LPX. On global average, the mean δ 13 C signature of all land biosphere carbon is 4.6 ‰ more negative in CLM4.5 than in LPX (Table 2).
Turning to 20th century changes, both models simulate a widespread and strong increase in GPP (Fig. 2), mainly in response to CO 2 fertilization under increasing atmospheric CO 2 . On the other hand, 20th century changes in transpiration are generally small in both models (Fig. 2). The two models bracket independent estimates of the net change in global land biosphere carbon stocks over the last century (20 GtC), with CLM4.5 showing slightly larger uptake (35 GtC) and LPX-Bern 1.3 showing a release of carbon (35 GtC) instead of an uptake.
Land use is treated differently in the two models. In CLM4.5, C3 and C4 species are replaced by a C3 crop when land is converted to cropland. This results in a negative 20th century change in δ 13 C of GPP of several per mill (gridcell average) in tropical and subtropical regions affected by land use (Fig. 4). In contrast, in LPX, C3 and C4 crops have identical PFT parameters to C3 and C4 grasses and grow in competition on pasture and cropland. Thus, conversion of tropical and subtropical forests (C3 trees) to cropland or pasture results in an increase of C4 grasses and a decrease of C3 plants. Correspondingly, a positive 20th century change in δ 13 C of up to 5 ‰ is simulated in regions undergoing landuse changes (Fig. 4). The implication is that for atmospheric δ 13 C budget analyses (Keeling et al., 1989;Ciais et al., 1995;Joos and Bruno, 1998;Scholze et al., 2008;van der Velde et al., 2013), the correct crop type (C3 versus C4) should be specified in the models.

δ 13 C in leaves of C3 plants: modern distribution
The mean and spatial gradients in discrimination of C3 trees simulated by CLM4.5 and LPX are evaluated using the leaf δ 13 C data of Cornwell et al. (2016). Both models reasonably represent the observation-based distribution in discrimination of C3 trees, though modeled discrimination is on average too small compared to the measurements, particularly in arid regions (Figs. 5, 6). The low discrimination in arid regions may be due to a mismatch in scale between local site conditions and grid-cell averages. This could be due to the fact that trees at sites with relatively good growing conditions were selected for the δ 13 C analysis, while modeled trees experience grid-cell average soil water conditions. In applications with the predecessor of LPX Scholze et al., 2008), δ 13 C was implemented following Scholze et al. (2003) with discrimination modeled following Lloyd and Farquhar (1994). Here, we adjusted the discrimination formulations and instead used the simpler formulations of Farquhar et al. (1989) for consistency with CLM4.5 and with the computation of iWUE from tree-ring δ 13 C measurements. We note that our conclusions do not depend on this choice. LPX yields similar 20th century changes in discrimination and iWUE for both formulations. However, simulated δ 13 C of carbon assimilated by C3 trees is generally less negative (by about 2 ‰) when applying the Lloyd and Farquhar (1994) formulation, and agreement with leaf δ 13 C is less favorable.
The difference between the two implementations arises from additional processes considered in the approach of Lloyd and Farquhar (1994) and from different parameter choice. Here, fractionations associated with the diffusion of CO 2 from the stomatal cavity to the cell wall, entrance of CO 2 in solution at the cell wall, and transport within the cell as well as photorespiration are neglected. The discriminations by the first three processes are set to be constant in the earlier implementation, following Scholze et al. (2003), neglecting a minor temperature dependency associated with the carboxylation by PEP-C (Lloyd and Farquhar, 1994). The discrimination associated with respiration varies with atmospheric CO 2 and the photocompensation point, but is small. In addition, b, the discrimination during photosynthetic fixation, is set to 27 in this study and to 27.5 in the earlier implementation. Overall, the consideration of the additional processes and the difference in b result in an approximately constant offset between the two implementations. This highlights that absolute values of discrimination depend on uncertain model parameters. Agreement or disagreement between leaf and model data (or tree-ring and model data) concerning absolute levels of δ 13 C (Figs. 5, 6) may not be interpreted as an indication of the performance of the conductancephotosynthesis module.

20th century changes in carbon isotopes and water-use efficiency of C3 plants
A set of 76 20th century δ 13 C tree-ring chronologies, mainly from Europe, boreal Asia, and western North America, was compiled (Table A1). The δ 13 C tree-ring data show on average little to no change in isotopic discrimination (Fig. 7). It remains unclear whether there is an overall small decrease in discrimination over the 20th century, given the still limited number of records and the large variability of the averaged and individual δ 13 C records. Small or no changes in discrimination of C3 plants imply that c i /c a remained roughly unchanged over the 20th century. It also implies that the change in iWUE is approximately equal to the relative change in atmospheric CO 2 (Eq. 15). This change is about 25 % over the 20th century and 43 % since preindustrial times. Recall, however, that isotopic signatures and the related variables (c i /c a , i , iWUE, A/g) as considered here represent annual or multiannual averages that have been weighted by C assimilation or allocation.
Discrimination i , c i /c a , and iWUE (A/g) as well as their changes are closely related. These variables hold, in the framework applied in this study (e.g., Eqs. 4, 5, 13, 15), basically the same information. For iWUE this information is transformed by the known atmospheric CO 2 evolution. As shown by Eq. (15), iWUE increases in proportion to CO 2 in the absence of a change in discrimination. A substantial increase in discrimination of more than 2 ‰ would be required to offset this change. It is therefore no surprise that Silva and Horwath (2013) find an increase in iWUE when using randomized δ 13 C records in a Monte Carlo analysis instead of actual tree-ring data. Their finding does not, however, invalidate the usefulness of iWUE as a physiological meaningful interpretation of δ 13 C tree-ring records, particularly as many studies have shown that nonrandomized δ 13 C measurements contain well-known environmental signals ranging from the Suess effect to interannual to centennial changes in climate.
We applied the tree-ring δ 13 C data to test whether the current CLM4.5 and LPX-Bern implementation for δ 13 C, stomatal conductance, and photosynthesis as well as related CO 2 fertilization, nitrogen limitation, and water-use efficiency mechanisms are consistent with the tree-ring records. Simulated changes by LPX agree on average excellently with the tree-ring data (Fig. 7). LPX shows little change in the discrimination by C3 trees and the evolution of the average change in iWUE closely matches the tree-ring reconstruction. In contrast to the data, CLM4.5 results show a steadily decreasing trend in discrimination by C3 trees, roughly in parallel with rising CO 2 concentrations (Figs. 1, 7). The average decrease in discrimination by C3 trees for the grid cells with data amounts to about 1.5 ‰ by the end of the simulation, and the increase in iWUE is about 2 times the change indicated by the tree-ring data. The decrease in discrimination corresponds to a decrease in c i /c a (Eq. 4). In other words, the leaf internal partial pressure c i is downregulated too strongly in CLM4.5, at least in the regions with tree-ring data.
Factorial simulations demonstrate that these overestimated trends are dominated by the response of CLM4.5 to increasing atmospheric CO 2 , and not by the response to changes in climate. Factorial and sensitivity simulations also show that differences in trends between CLM4.5 and LPX are not related to trends in soil moisture nor to (small) differences in the applied climate forcing data. The downregulation in c i /c a in CLM4.5 therefore does not reflect a response to drought or worsening climatic conditions for CO 2 assimilation. These findings suggest that the main reason for the model-model and model-tree-ring data difference in the trends of discrimination and iWUE is rooted in the parameterization of the photosynthesis-conductance coupling.
It is difficult to trace origins of data-model mismatch in complex models such as CLM4.5. There are two aspects in the current implementation of photosynthesis and conductance in CLM4.5 that may be problematic (Bonan et al., 2014;Ghimire et al., 2016). First, nitrogen downregulation of photosynthesis and conductance occurs on the subhourly model time step, although it is not plausible that leaf structures adjust so quickly. Second, the relationship between assimilation of carbon and transpiration of water (Ball et al., 1987) is prescribed with time-invariant and globally constant parameters.
Regarding nitrogen limitation, photosynthesis for isotopic discrimination is downregulated immediately, on the subhourly time step of the model, by limited nitrogen availability in CLM4.5 (Eq. 6). This can lead to a depression of assimilation in the isotope calculation during times of high assimilation. As explained in the method section, photosynthesis of carbon and stomatal conductance for water is computed, unlike for discrimination, without nitrogen limitation. Instead, the allocation of GPP to carbon pools is downreg-ulated under nitrogen limitation. This downregulation is on annual average generally less than 20 % in forested regions, and 20th century changes in downregulation are generally small (within ±2 %) in areas with tree-ring data. Raczka et al. (2016) applied a recalibrated version of CLM4.5 at a single site (Niwot Ridge, US) in simulations with and without nitrogen limitation. These authors show that nitrogen downregulation strongly affects the change in average c i /c a over the 20th century. At their site, changes in c i /c a and in discrimination were positive over the industrial period and the century-scale change in discrimination was smaller when nitrogen limitation was not active. The observation-estimated seasonal cycle in discrimination at Niwot Ridge is better reproduced without nitrogen limitation than with limitation. Indeed, Raczka et al. (2016) suggest that downregulation of assimilation and c i /c a by nitrogen limitation (Eq. 6) may be too strong in CLM4.5. Ghimire et al. (2016) implemented an alternative formulation for nitrogen limitation in CLM4.5, where nitrogen availability affects the maximum rate of photosynthesis (Vcmax) on slow timescales through changes in the leaf carbon-to-nitrogen ratio. This reduces the global bias in GPP, leaf area index, and biomass and improved water-use efficiency predictions compared to the original CLM4.5 formulation.
Another possibility for the model-data mismatch is that the photosynthesis formulation in CLM4.5 may not adequately represent the relationship between stomatal conductance and assimilation and related adjustment processes to long-term environmental changes. CLM4.5 employs the experimentally well-verified and widely used Ball-Berry equation (Eq. 7) applying globally uniform, time invariant slope parameters, m, for C3 and C4 plants. Thus, any potential adjustment of m to changes in environmental conditions, including the century-scale increase in atmospheric CO 2 or changes in water stress are not considered. It is currently unclear whether such adjustment processes occur and the current understanding of the underlying physiological mechanisms of stomatal responses is incomplete (e.g. Miner et al., 2017). The value of m is set to nine in CLM4.5 for all C3 plants world-wide. This value is within the observational range, but species differences (e.g., evergreen versus deciduous) are neglected by using a single value. The review of Miner et al. (2017) yields mean values for m of 9.8, 8.7, and 6.8 for angiosperm evergreen, angiosperm deciduous and for gymnosperm trees, respectively (their Fig. 1). Aranibar et al. (2006) inferred values of m around 10 to 12 using observed foliar 13 C values in pine.
Equations (7) and (6), defining the photosynthesisconductance coupling are evaluated on the subhourly time step of the model and the variables A, RH, CO 2 , and f dreg vary daily and seasonally, while we investigate here decadalto-century-scale trends. This hampers a simple interpretation of the two equations, in particular in terms of expected temporal changes in iWUE and whether expected century-scale changes in discrimination may be positive or negative. As shown in Fig. 8, simulated changes in discrimination are positive in small areas of South America and southeastern Australia, and thus opposite to the world-wide trend simulated by CLM4.5. Sato et al. (2015) investigated the influence of the use of vapor pressure deficit (VPD) versus relative humidity (RH) in Ball-Berry-type stomatal conductance formulations. About half of the investigated models apply RH as a driving variable (Eq. 7), despite the fact that VPD is considered the more relevant controlling factor. The global warming simulations reveal an increase in VPD and little change in RH. Their results suggest that the increase in VPD under global warming leads to a stronger downregulation of stomatal conductance (g s ) and c i /c a for the VPD formulations compared to the RH formulation. This implies that replacing RH by VPD in the Ball-Berry equation in CLM4.5 may, without further adjustments, even increase the disagreement between modeled and reconstructed changes in discrimination and in iWUE. Duarte et al. (2016) and Raczka et al. (2016) recalibrated CLM4.5 to match site-specific conditions in a conifer forest in the northwestern US and at Niwot Ridge, Colorado. These studies find a substantial 20th century increase in discrimination for their site-calibrated versions. One does expect from Eq. (7) that a change in the slope parameter m will result in a proportional change in the absolute magnitude of iWUE and a corresponding change in discrimination. Duarte et al. (2016) reduced m by a third in their single-site version. In turn, simulated discrimination was altered by about 2.5 ‰, roughly corresponding to the expected 33 % change in iWUE. In addition, the 20th century change in isotopic discrimination is reduced. In brief, model structure in CLM4.5 permits both positive and negative changes in discrimination under rising CO 2 . The slope parameter m, potentially affecting trends in discrimination, is set as time invariant and globally uniform for all C3 plants. The relationship between A and g s is assumed to be linear for constant CO 2 and relative humidity.
In LPX, the shape of the relationship between g and A is not as tightly prescribed as in CLM4.5. Rather, A and c i /c a are optimized for given environmental conditions as described in the method section. This optimization leads to small 20th century changes in simulated isotopic discrimination in agreement with the observational evidence. The parameterization used in LPX is for the daily model time step and may not be readily transferred to models with a subhourly time step.
Alternative conductance-photosynthesis formulations may be preferable, compared to the Ball-Berry relation as used in CLM4.5. The Ball-Berry relation is viewed as consistent with an optimization (Cowan and Farquhar, 1977) of stomatal conductance towards maintaining a constant water-use efficiency by optimizing carbon gain per unit water loss. However, a number of alternative formulations are found in the literature. These differ in potentially important ways. For example, Medlyn et al. (2011) suggest that the slope parameter increases with increasing temperature; further, the inverse of the square root of vapor pressure deficit is used instead of relative humidity in their relationship. Prentice et al. (2014) and Wang et al. (2016) rely on an optimality hypothesis considering the costs of maintaining both water flow and photosynthetic capacity. These authors use information on the spatial gradients in δ 13 C from stable isotope measurements on leaf material  to develop their photosynthesisconductance relation. The slope parameter proposed by Prentice et al. (2014) depends on a number of variables, e.g., the temperature-dependent Michaelis-Menten coefficient for Rubisco-limited photosynthesis, viscosity of water, or the water potential difference between soil and leaf or foliage height. Hence, the slope parameter varies with environmental conditions in this approach. Bonan et al. (2014) implemented different photosynthesis-conductance modules within a CLM4.5 model version featuring a multilayer canopy and compared results with leaf analyses and eddy covariance fluxes at six forest sites. The continuous soil-plant-atmosphere module, optimizing carbon gain per unit water loss and considering hydraulic safety, performs similarly to or better than the Ball-Berry formulations used in CLM4.5. A better performance is achieved under soil-moisture-stressed conditions in particular.
It is a task for the future to explore whether the implementation of such optimizing modules will lead to a more realistic simulation of the spatiotemporal changes in isotopic discrimination. Taken together, agreement between CLMmodeled and observational trends in discrimination may be improved in the future by adjusting model parameters such as the slope parameter m, replacing formulations for nitrogen limitation, e.g., following Ghimire et al. (2016), and by implementing photosynthesis-conductance routines that adhere to an optimization principle as proposed in the literature (Bonan et al., 2014;Prentice et al., 2014).

Consistency with previous studies
The small changes in discrimination, an increase in iWUE proportional to the atmospheric CO 2 increase, and approximately constant c i /c a over the 20th century as reconstructed by our tree-ring compilation and simulated by LPX-Bern is consistent with most, but not all, studies. As reported by Voelker et al. (2016), studies that are consistent with the notion of constant c i /c a include early work by Wong et al. (1979) measuring A and g of leaves in a chamber, and a range of tree-ring δ 13 C studies (Saurer et al., 2004;Ward et al., 2005;Bonal et al., 2011;Franks et al., 2013), as well as a metaanalysis of FACE experiments (Ainsworth and Long, 2005). The European δ 13 C tree-ring records analyzed by Frank et al. (2015) and Saurer et al. (2014) also point to a moderate control towards a constant c i /c a ratio. Leonardi et al. (2012) conclude that the temporal variation in δ 13 C in their long-term isotope tree-ring chronologies for 53 sites worldwide supports the hypothesis of an active plant mechanism that maintains a constant ratio between intercellular and ambient CO 2 concentrations. Lévesque et al. (2014) reported, from their δ 13 C tree-ring data, an increase in iWUE over the last 50 years in the range of 8 to 29 % for xeric and mesic sites in the Alps and Switzerland. At their sites, drought-induced stomatal closure has reduced transpiration at the cost of reduced carbon uptake and growth. Churakova (Sidorova) report different iWUE strategies with almost constant c i /c a for European larch since the 1990s and continuously increasing iWUE for mountain pine trees since the 1980s in the Swiss National Park. Liu et al. (2014a) find moderate changes in c i /c a and an increase in iWUE by 36 % in a riparian forest in northwestern China from 1920 to 2012. Similarly, Peñuelas et al. (2011), analyzing changes in treering δ 13 C and growth at 47 sites worldwide, inferred little change in discrimination, c i /c a , and an increase in iWUE of 20.5 % from 1960 to 2000. Out of 35 studies with growth data, 18 show an increase in growth with time, while the others show no growth changes or negative growth changes. Voelker et al. (2016) analyzed studies of δ 13 C and photosynthetic discrimination in woody angiosperms and gymnosperms that grew across a range of CO 2 spanning at least 100 ppm and combining paleodata, tree-ring records, and FACE-type experiments. They conclude that woody plants respond to increasing CO 2 by regulating leaf gas exchange along a continuum of c a − c i and c i /c a that minimizes water loss for a given amount of C gain and therefore increasingly minimizes the likelihood of exposure to drought stress. Summarizing 5 years of results from the Basel FACE experiment, Klein et al. (2016) find that iWUE increased in their experiment by 38 % at the needle level, as a result of higher assimilation at constant conductance. Interestingly, Klein et al. (2016) could not identify an increase in plant carbon stocks corresponding to the increase in assimilation, and the fate of the additionally assimilated carbon remains unclear. The 38 % increase in iWUE is comparable to the increase from ambient (400 ppm) to elevated (550 ppm) CO 2 concentrations and in line with observations reported by De Kauwe et al. (2013) as well as a metaanalysis of FACE experiments (Ainsworth and Long, 2005). Streit et al. (2014) measured needle gas exchange and analyzed δ 13 C of needles and tree rings of Larix decidua and Pinus mugo after 9 years of free air CO 2 enrichment. Both species showed an increase in net photosynthesis and, in agreement with Keel et al. (2006), small changes or no changes in conductance under elevated CO 2 and a change in iWUE roughly in proportion to the increase in CO 2 . Elevated CO 2 induced increased basal area growth in L. decidua, but not in P. mugo. Neither nitrogen limitation, end-product limitation, nor, in agreement with other FACE studies (Bader et al., 2010;Liberloo et al., 2007), downregulation of maximal photosynthetic rate was found.
In contrast, Battipaglia et al. (2013) report an increase in iWUE between 50 and 90 % for the ORNL, DUKE, and POP-EUROFACE Face sites and that c i was likely maintained constant and 13 C discrimination reduced under elevated CO 2 . We note that this finding, based on 13 C treering isotope measurements, differs remarkably from results by De Kauwe et al. (2013). De Kauwe et al. (2013) suggest a proportional increase in iWUE at the Duke and ORNL FACE sites based on in situ foliage gas-exchange measurements (see their Fig. 6). This difference might be caused by the comparison of leaf-level gas-exchange and stem-level data. Integration times may differ between the two data streams and stem level data are influenced by growth rate and storage compounds. A scenario of constant c i under increasing CO 2 concentrations is also suggested by Keenan et al. (2013), who report a strong increase in water-use efficiency summarizing results from eddy-covariance measurements at about 20 temperate and boreal forest sites in the Northern Hemisphere. However, these studies may be affected by local conditions which are not representative for larger scales and may also be influenced by temporal sampling biases, as the record length of FACE and eddy-covariance measurements is limited. Indeed, the individual tree-ring records in Fig. 7 show considerable site-to-site variability as well as decadal variability. In addition, factors other than c i /c a , e.g., vapor pressure deficit, may influence trends at eddy covariance sites (Knauer et al., 2017). The suggestion of a constant c i under increasing CO 2 is challenged by Knauer et al. (2017). These authors find that the ecosystem trends reported by Keenan et al. (2013) and a scenario of a constant c i , as also suggested by Battipaglia et al. (2013), are in conflict with observed large-scale trends in continental discharge, evapotranspiration, and the seasonal CO 2 exchange. Rather, the comparison of observational data and model outcome by Knauer et al. (2017) support the finding of a physiological regulation towards a constant c i /c a under rising atmospheric CO 2 .

Conclusions
We compiled 20th century δ 13 C tree-ring records for 73 sites (Table A1). The records are used to reconstruct changes in the ratio of CO 2 concentration within the stomatal cavity of C3 trees and the atmosphere, c i /c a , and to reconstruct changes in the intrinsic water-use efficiency, iWUE, denoting the ratio of photosynthesis to conductance (A/g). The treering results suggest, on average, constant c i /c a over the 20th century and iWUE to have increased in proportion to atmospheric CO 2 . The simulated changes in discrimination, c i /c a , and iWUE are consistent with the reconstructions for LPX-Bern 1.3, whereas CLM4.5 shows larger trends and overestimates 20th century changes in iWUE by almost a factor of 2 (Fig. 7). This suggests that it is desirable to adjust the implementation of photosynthesis and conductance in CLM4.5 towards a better agreement with observation-derived centuryscale trends in 13 C discrimination and intrinsic water-use efficiency. In conclusion, this study demonstrates that existing information on the magnitude and trends of stable carbon isotopes permits the evaluation of global land carbon cycle models. In particular, δ 13 C data provide insight into the fundamental relationship between carbon assimilation and stomatal conductance controlling the flow of CO 2 and water. Treering records are useful as they retrospectively cover decadalto-century timescales which are not accessible by laboratory or field experiments and the relevant instrumental records, and which are directly relevant for the response of plants to the decadal-to-century-scale rise in CO 2 and to global warming.
Data availability. Model output data and the averaged tree ring data shown in Fig. 7 Leavitt et al. (2007)