alexa The Water Isotopic Version of the Land-Surface Model ORCHIDEE: Implementation, Evaluation, Sensitivity to Hydrological Parameters
ISSN: 2157-7587
Hydrology: Current Research

Like us on:

Make the best use of Scientific Research and information from our 700+ peer reviewed, Open Access Journals that operates with the help of 50,000+ Editorial Board Members and esteemed reviewers and 1000+ Scientific associations in Medical, Clinical, Pharmaceutical, Engineering, Technology and Management Fields.
  • Research Article   
  • Hydrol Current Res 2016, Vol 7(4): 258
  • DOI: 10.4172/2157-7587.1000258

The Water Isotopic Version of the Land-Surface Model ORCHIDEE: Implementation, Evaluation, Sensitivity to Hydrological Parameters

Camille Risi1*, Jerome Ogée2, Sandrine Bony1, Thierry Bariac3, Naama Raz-Yaseef4,5, Lisa Wingate2, Jeffrey Welker6, Alexander Knohl7, Cathy Kurz-Besson8, Monique Leclerc9, Gengsheng Zhang9, Nina Buchmann10, Jiri Santrucek11,12, Marie Hronkova11,12, Teresa David13, Philippe Peylin14 and Francesca Guglielmo14
1LMD/IPSL, CNRS, UPMC, Paris, France
2INRA Ephyse, Villenave d’Ornon, France
3UMR 7618 Bioemco, CNRS-UPMC-AgroParisTech-ENS Ulm-INRA-IRD-PXII Campus AgroParisTech, Bâtiment EGER, Thiverval-Grignon, 78850, France
4Earth Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, USA
5Department of Environmental Sciences and Energy Research, Weizmann Institute of Science, PO Box 26, Rehovot 76100, Israel
6Biology Department and Environment and Natural Resources Institute, University of Alaska, Anchorage, AK 99510, USA
7Bioclimatology, Faculty of Forest Sciences and Forest Ecology, Georg-August University of Göttingen, 37077 Göttingen, Germany
8Instituto Dom Luiz, Centro de Geofísica IDL-FCUL, Lisboa, Portugal
9University of Georgia, Griffin, GA 30223, USA
10Institute of Agricultural Sciences, ETH Zurich, Zurich, Switzerland
11Biology Centre ASCR, Branisovska 31, Ceske Budejovice, Czech Republic
12University of South Bohemia, Faculty of Science, Branisovska 31, Ceske Budejovice, Czech Republic
13Instituto Nacional de Investigação Agrária e Veterinária, Quinta do Marquês, Portugal
14LSCE/IPSL, CNRS, UVSQ, Orme des Merisiers, Gif-sur-Yvette, France
*Corresponding Author: Camille Risi, LMD/IPSL, CNRS, UPMC, Paris, France, Tel: +33144276272, Email: [email protected]

Received Date: Aug 30, 2016 / Accepted Date: Sep 24, 2016 / Published Date: Oct 01, 2016


Land-Surface Models (LSMs) exhibit large spread and uncertainties in the way they partition precipitation into surface runoff, drainage, transpiration and bare soil evaporation. To explore to what extent water isotope measurements could help evaluate the simulation of the soil water budget in LSMs, water stable isotopes have been implemented in the ORCHIDEE (ORganizing Carbon and Hydrology In Dynamic EcosystEms: the land-surface model) LSM. This article presents this implementation and the evaluation of simulations both in a stand-alone mode and coupled with an atmospheric general circulation model. ORCHIDEE simulates reasonably well the isotopic composition of soil, stem and leaf water compared to local observations at ten measurement sites. When coupled to LMDZ (Laboratoire de Météorologie Dynamique-Zoom: the atmospheric model), it simulates well the isotopic composition of precipitation and river water compared to global observations. Sensitivity tests to LSM (Land-Surface Model) parameters are performed to identify processes whose representation by LSMs could be better evaluated using water isotopic measurements. We find that measured vertical variations in soil water isotopes could help evaluate the representation of infiltration pathways by multi-layer soil models. Measured water isotopes in rivers could help calibrate the partitioning of total runoff into surface runoff and drainage and the residence time scales in underground reservoirs. Finally, co-located isotope measurements in precipitation, vapor and soil water could help estimate the partitioning of infiltrating precipitation into bare soil evaporation.

Keywords: Water isotopes; Land-surface model; Global models; Soil water budget; Rain infiltration; Runoff; Evapo-transpiration partitioning


Land-surface models (LSMs) used in climate models exhibit a large spread in the way they partition radiative energy into sensible and latent heat [1,2] precipitation into evapo-transpiration and runoff [3-5], evapo-transpiration into transpiration and bare soil evaporation [6,7], and runoff into surface runoff and drainage [8-10]. This results in an large spread in the predicted response of surface temperature [11] and hydrological cycle [12,13] to climate change [11] or land use change [14,15]. Therefore, evaluating the accuracy of the partitioning of precipitation into surface runoff, drainage, transpiration and bare soil evaporation (hereafter called the soil water budget) in LSMs is crucial to improve our ability to predict future hydrological and climatic changes.

The evaluation of LSMs is hampered by the difficulty to measure over large areas the different terms of the soil water budget, notably the evapo-transpiration terms and the soil moisture storage [16,17]. Single point measurements of evapo-transpiration fluxes [18] and soil moisture [19] are routinely performed within international networks, but those measurements remain difficult to upscale to a climate model grid box due to the strong horizontal heterogeneity of the land surface [20,21]. Spatially-integrated data such as river runoff observations are very valuable to evaluate soil water budgets at the regional scale [22,23], but are insufficient to constrain the different terms of the water budget. Additional observations are therefore needed.

In this context, water isotope measurements have been suggested to help constrain the soil water budget [24,25], its variations with climate or land use change [26], and its representation by large-scale models [27,28]. For example, water stable isotope measurements in the different water pools of the soil-vegetation-atmosphere continuum have been used to quantify the relative contributions of transpiration and bare soil evaporation to evapo-transpiration [29-32], to infer plant source water depth [33], to assess the mass balance of lakes [34-36] or to investigate pathways from precipitation to river discharge [37- 40]. These isotope-based techniques generally require high frequency isotope measurements and are best suitable for intensive field campaigns at the local scale. At larger spatial and temporal scales, some attempts have been made to use regional gradients in precipitation water isotopes for partitioning evapo-transpiration into bare soilevaporation and transpiration [41-43].

To explore to what extent water isotope measurements could be used to evaluate and improve land surface parameterizations, water isotopes were implemented in the LSM ORCHIDEE (ORganizing Carbon and Hydrology In Dynamic EcosystEms [44,45]. This isotopic version of ORCHIDEE has already been used to explore how tree-ring cellulose records past climate variations [46] and to investigate the continental recycling and its isotopic signature in Western Africa [47] and at the global scale [48].

The first goal of this article is to evaluate the isotopic version of the ORCHIDEE model against recently-made-available new datasets combining water isotopes in precipitation, vapor, soil water and rivers. The second goal is to evaluate the isotopic version of the ORCHIDEE model when coupled to the atmospheric General Circulation Model (GCM) LMDZ (Laboratoire de Météorologie Dynamique Zoom [49]). The third goal is to perform sensitivity tests to LSM parameters to identify processes whose representation by LSMs could be better evaluated using water isotopic measurements.

After introducing notations and models in section 4, we present ORCHIDEE simulations in a stand-alone mode at measurement sites and global ORCHIDEE-LMDZ coupled simulations.

Notation and Models


Isotopic ratios hydrology-current-research in the different water pools are expressed in‰ relative to a standard: hydrology-current-research where Rsample and RSMOW are the isotopic ratios of the sample and of the Vienna Standard Mean Ocean Water (V-SMOW) respectively [50,51]. To first order, variations in δD are similar to those in δ18O but are 8 times larger. Deviation from this behavior can be associated with kinetic fractionation and is quantified by deuterium excess (d=δD-8.δ18O [50,52]). Hereafter, we note δ18Op, δ18Ov, δ18Os, δ18Ostem and δ18Oriver the δ18O of the precipitation, atmospheric vapor, soil, stem, river water respectively. The same subscripts apply for d.

The LMDZ model

LMDZ is the atmospheric GCM (General Circulation Model) of the IPSL (Institut Pierre Simon Laplace) climate model [53,54]. We use the LMDZ-version 4 model [49] which was used in the International Panel on CLimate Change’s Fourth Assessment Report simulations [55,56]. The resolution is 2.5 ° in latitude, 3.75 ° in longitude and 19 vertical levels. Each grid cell is divided into four sub-surfaces: ocean, land ice, sea ice and land (treated by ORCHIDEE) (Figure 1a). All parameterizations, including ORCHIDEE, are called every 30 min. The implementation of water stable isotopes is similar to that in other GCMs [57,58] and has been described in [59,60]. LMDZ captures reasonably well the spatial and seasonal variations of the isotopic composition in precipitation [60] and water vapor [61].


Figure 1: The four sub-surfaces in the LMDZ GCM: land, ocean, sea ice and land ice. Their relative fraction in each grid box is prescribed. The sea surface temperature of the ocean is prescribed, and interactively calculated for sea-ice and land-ice. Over land, the land-Surface model (LSM) ORCHIDEE calculates interactively the surface temperature and outgoing water fluxes. b) Water fluxes and pools represented in the ORCHIDEE LSM. Water pools are the soil water in the superficial (qsg) and bottom (qsb) layers, the water intercepted by the canopy (qw) and the snow pack (qsnow). Fluxes onto the land surface are the total rain (P) and snow (S), and possibly dew or frost. As some rain is intercepted by the canopy, only throughfall rain (Ps) arrives at the soil surface. Evaporation fluxes are the evaporation of intercepted water (Ew), transpiration by the vegetation (T), bare soil evaporation (E) and snow sublimation (Es). Snow melt may be transferred from the snow pack to the soil (M). Water from rainfall, melt (and possibly dew) exceeding the soil capacity is converted to surface runoff (R) and drainage (D). The routing model then transfers surface runoff and drainage to streams.

The ORCHIDEE (ORganizing Carbon and Hydrology In Dynamic EcosystEms: the land-surface model) model

The ORCHIDEE model is the LSM component of the IPSL climate model. It merges three separate modules: (1) SECHIBA (Schématisation des EChanges Hydriques a l’Interface entre la Biosphère et l’Atmosphère [44,62]) that simulates land-atmosphere water and energy exchanges, (2) STOMATE (Saclay-Toulouse-Orsay Model for the Analysis of Terrestrial Ecosystems [45]) that simulates vegetation phenology and biochemical transfers ; and (3) LPJ (Lund- Postdam-Jena [63]) that simulates the vegetation dynamics. Water stable isotopes were implemented in SECHIBA, and we use prescribed land cover maps so that the two other modules could be de-activated.

Each grid box is divided into up to 13 land cover types: bare soil, tropical broad-leaved ever-green, tropical broad-leaved rain-green, temperate needle-leaf ever-green, temperate broad-leaved ever-green, temperate broad-leaved summer-green, boreal needle-leaf ever-green, boreal broad-leaved summer-green, boreal needle-leaf summer-green, C3 grass, C4 grass, C3 agriculture and C4 agriculture. Water and energy budgets are computed for each land cover type.

Figure 1b illustrates how ORCHIDEE (ORganizing Carbon and Hydrology In Dynamic EcosystEms: the land-surface model) represents the surface water budget. Rainfall is partitioned into interception by the canopy and through-fall rain. Through-fall rain, snow melt, dew and frost fill the soil. The soil is represented by two water reservoirs: a superficial and a bottom one [64,65]. Taken together, the two reservoirs have a water holding capacity of 300 mm and a depth of 2 m. Soil water undergoes transpiration by vegetation, bare soil evaporation or runoff. Transpiration and evaporation rates depend on soil moisture to represent water stress in dry conditions. Runoff occurs when the soil water content exceeds the soil holding capacity and is partitioned into 95% drainage and 5% surface runoff [66]. Snowfall fills a singlelayer snow reservoir, where snow undergoes sublimation or melt. By comparison, when not coupled to ORCHIDEE, the simple bucketlike LSM in LMDZ makes no distinction neither between bare soil evaporation and transpiration nor between surface runoff and drainage [67].

Surface runoff and drainage are routed to the coastlines by a water routing model [68]. Surface runoff is stored in a fast ground water reservoir which feeds the stream reservoir with residence time of 3 days. Drainage is stored in a slow ground water reservoir which feeds the stream reservoir with residence time of 25 days. The water in the stream reservoir is routed to the coastlines with a residence time of 0.24 days.

Implementation of water stable isotopes in ORCHIDEE

We represent isotopic processes in a similar fashion as other isotopeenabled LSMs [69-73]. Some details of the isotopic implementation are described in Risi [74]. In absence of fractionation, water stable isotopes hydrology-current-research are passively transferred between the different water reservoirs. We assume that surface runoff has the isotopic composition of the rainfall and snow melt that reach the soil surface. Drainage has the isotopic composition of soil water [24]. We calculate the isotopic composition of bare soil evaporation or of evaporation of water intercepted by the canopy using the Craig and Gordon equation [75] (Equation 3). We neglect isotopic fractionation during snow sublimation (Equation 2). We consider isotopic fractionation at the leaf surface (Equation 5) but we assume that transpiration has the isotopic composition of the soil water extracted by the roots (Equation 2).

In the control coupled simulation, we assume that the isotopic composition of soil water is homogeneous vertically and equals the weighted average of the two soil layers. However, transpiration, bare soil evaporation, surface runoff and drainage draw water from different soil water reservoirs whose isotopic composition is distinct [76-78]. Therefore, we also implemented a representation of the vertical profile of the soil water isotopic composition.

Stand-alone ORCHIDEE Simulations at MIBA (Moisture in Biosphere and Atmosphere: Network for Water Isotopes in Soil, Stem and Leaf Water) and Carbo- Europe Measurement Sites

First, we performed simulations using ORCHIDEE as a standalone model at ten sites. Using isotopic measurements in soil, stem and leaf water, simulations are evaluated at each site at the monthly scale. Sensitivity tests to evapo-transpiration partitioning and soil infiltration processes are performed.

Measurements used for evaluation

To first order the composition of all land surface water pools is driven by that in the precipitation [79]. Therefore, a rigorous evaluation of an isotope-enabled LSM requires to evaluate the difference between the composition in each water pool and that in the precipitation. Besides, to better isolate isotopic biases, we need a realistic atmospheric forcing. We tried to select sites where (1) isotope were measured in different water pools of the soil-plant-atmosphere continuum, during at least a full seasonal cycle and (2) meteorological variables were monitored at a frequency high enough (30 minutes) to ensure robust forcing for our model and (3) water vapor and precipitation were monitored to provide isotopic forcing for the LSM. Only two sites satisfy these conditions: Le Bray and Yatir. Relaxing some of these conditions, we got a more a representative set of ten sites representing diverse climate conditions (Table 1 and Figure 2).


Figure 2: Location of the ten stations used in this study for single-point model-data comparison. The background represents the annual-mean precipitation from GPCP (Global Precipitation Climatology Project) to illustrate the diversity of climate regimes covered by the ten stations. Each station is described in more detail in Table 1.

Site name Country Location Network Years Reference
Le Bray France 44.70° N, 0.77° W MIBA, Carbo- Euroe 2007-2008 [87
Yatir Israel 31.33° N, 35.0° E MIBA, Carbo-Euroe 2004-2005 [89,104]
Morgan-Monroe United States 39.32° N, 86.42° W MIBA-US 2005-2006 [167,172]
Donaldson Forest United States 29.8° N, 82.163° W MIBA-US 2005-2006 [91,169]
Anchorage United States 61.2° N, 149.82° W MIBA-US 2005-2006  -
Mitra Portugal 38.5° N, 8.00° W Carbo-Euroe 2001-2002 [171]
Bily Kriz Czech Republic 49.5° N, 18.53° E MIBA, Carbo-Euroe 2005 [92]
Brloh Czech Republic 49.80° N, 14.66° E MIBA 2004-2010 [93]
Hainich Germany 50.97° N, 13.57° E Carbo-Euroe 2001-2002 [170]
Tharandt Germany 51.08° N, 10.47° E Carbo-Euroe 2001-2002  -

Table 1: Information on the 10 sites used in this study: geographical location, network the sites are part of, years during which the istopic measurements were made and are used in this study, reference.

Description of the ten sites: The ten sites belong to two kinds of observational networks: MIBA (Moisture Isotopes in the Biosphere and Atmosphere [80-82] or Carbo-Europe [83,84].

Le Bray site, in South-Eastern France, joined the MIBA and GNIP (Global Network for Isotopes in Precipitation) network in 2007. It is an even-aged Maritime pine forest with C3 grass understory that has been the subject of many eco-physiological studies since 1994, notably as part of the Carbo-Europe flux network [85]. In 2007 and 2008, samples in precipitation, soil surface, needles, twigs and atmospheric vapor were collected every month and analyzed for δ18O following the MIBA protocol [82,86]. This site was also the subject of intensive campaigns where soil water isotope profiles were collected between 1993 and 1997, and in 2007 [87].

The Yatir site, in Israel, is a semi-arid Aleppo pine forest. It is an afforestation growing on the edge of the desert, with mean-annual precipitation of 280 mm [88,89]. It has also been the subject of many eco-physiological studies as part of the Carbo-Europe flux network [89] and joined the MIBA network in 2004. It. In 2004-2005, samples of soil water at different depth, stems and needles were collected following the MIBA protocol. The water vapor isotopic composition has been monitored daily at the nearby Rehovot site (31.9 ° N, 34.65 ° E, [90]) and is used to construct the water vapor isotopic composition forcing. We must keep in mind however that although only 66 km from Yatir, Rehovot is much closer to the sea and is more humid than Yatir. The precipitation isotopic composition has been monitored monthly at the nearby GNIP station Beit Dagan (32 ° N, 34.82 ° E) and is used to construct the precipitation isotopic composition forcing.

The Morgan-Monroe State Forest, Donaldson Forest and Anchorage sites are part of the MIBA-US (MIBA-United States) network and are located in Indiana, in Florida and in Alaska respectively (Table 1). Sampling took place in 2005 and 2006 according to the MIBA protocols. The Donaldson Forest site, which jointed the MIBA-US network in 2005, is located at the AmeriFlux Donaldson site near Gainesville, Florida, USA. The site is flat with an elevation of about 50 m. It was covered by a forest of managed slash pine plantation, with an uneven understory composed mainly of saw palmetto, wax myrtle and Carolina jasmine [91]. The leaf area index was measured during a campaign in 2003 and estimated at 2.85. We use this value in our simulations.

The Mitra, Bily Kriz, Brloh, Hainich and Tharandt sites are part of the Carbo-Europe project. Hainich and Tharandt are located in Germany. The experimental site of Herdade da Mitra (230 m altitude, nearby Évora in southern Portugal) is characterized by a Mediterranean mesothermic humid climate with hot and dry summers. It is a managed agroforestry system characterized by an open evergreen woodland sparsely covered with Quercus suber L. and Q. ilex rotundifolia trees (30 trees/ha), with an understorey mainly composed of Cistus shrubs, and winter-spring C3 annuals. The isotopic samplings of leaves, twigs, soil, precipitation and groundwater were performed on a seasonal to monthly basis. All samples where extracted and analyzed at the Paul Scherrer Institute (Switzerland).

Bily Kriz and Brloh are both located on the Czech Republic. Bily Kriz is an experimental site in Moravian–Silesian Beskydy Mountains (936 m a.s.l.) with detailed records of environmental conditions [92]. It is dominated by Norway spruce forest. It joined the MIBA project in the season 2005. Brloh is a South Bohemian site in the Protected Landscape Area Blanskýles (630 m a.s.l.). It is dominated by deciduous beech forest and was used as MIBA sampling site from 2004 to 2010 [93].

Isotopic measurements: Samples of soil water, stems and leaves were collected at the monthly scale. The MIBA and MIBA-US protocols recommend sampling the first 5-10 cm excluding litter and the Carbo-Europe protocol recommends sampling the first 5 cm [84], but in practice the soil water sampling depth varies from site to site. At some sites, soil water was sampled down to 1 m. For evaluating the seasonal evolution of soil water δ18O , we focus on soil samples collected in the first 15 cm only. Observed full soil water δ18O profiles were used only at Le Bray and Yatir for evaluating the shape of simulated soil water δ18O profiles.

Carbo-Europe samples were extracted and analyzed at the Department of Environmental Sciences and Energy Research, Weizmann Institute of Science, Israel. MIBA-US samples were extracted and analyzed at the Center for Stable Isotope Biogeochemistry of the University of California, Berkeley. Analytical errors for δ18O in soil, stem and leaf water vary from 0.1‰ to 0.2‰ depending on the sites and involved stable isotope laboratory (Table 2).

Site name Biome Dominant Species Annual-mean temperature (°C) Annual-mean precipitation
Elevation (m)
Le Bray Temperate coniferous forest Maritime pine 12 1022 60
Yatir semi-arid forest Aleppo pine 15.3 270 650
Morgan-Monroe Temperate deciduous forest Liriodendron tulipifera 12.4 1094 275
Donaldson Forest Tropical pine plantation Pinus palustris 21.7 1330 50
Anchorage Boreal coniferous forest Picea glauca 2.3 408 35
Mitra Mediteranean forest Sparse holm oak trees with patches of cork trees 13.9 480 230
Bily Kriz Temperate coniferous forest Pine forest 3.4 1024 936
Brloh Temperate deciduous forest Beech forest 7.6 832 630
Hainich Temperate deciduous forest Fagus Sylvatica 8 800 440
Tharandt Temperate deciduous forest Pine forest 8.1 1000 380

Table 2: Vegetation and climtological information on the 10 sites used in this study: biome, dominant species, annual-mean temperature and precipitation, elevation.

Meteorological, turbulent fluxes and soil moisture measurements: At most of the sites, meteorological parameters (radiation, air temperature and humidity, soil temperature and moisture) are continuously measured and are used to construct the meteorological forcing for ORCHIDEE.

Fluxes of latent and sensible energy are measured using the eddy co-variance technique and are used for evaluating the hydrological simulation. Gaps are filled using ERA-Interim reanalyses [94]. Soil moisture observations are available at most sites.

Simulation set-up

To evaluate in detail the isotope composition of different water pools, stand-alone ORCHIDEE simulations on the ten MIBA and Carbo-Europe sites were performed. We prescribe the vegetation type and properties and the bare soil fraction based on local knowledge at each site (Table 3).

Site name Prescribe vegetation in ORCHIDEE Meteoro-logical forcing Isotopic forcing for precipitation and vapor local, GNIP, USNIP or Carbo-Europe stations used to calculate isotopic forcing
Le Bray 70% temperate needleleaf evergreen (LAI=0.4), 30% C3 grass (LAI=0.4) obs obs_iso Le Bray local data for both precipitation and water vapor
Yatir 100% temperate needleleaf evergreen (LAI=4) obs obs_iso Rehovot for water vapor and Beit Dagan GNIP station for precipitation
Morgan-Monroe 100% temperate broad-leaved summergreen (LAI=4.5) obs_ERA NIP_LMDZ USNIP_IN22, USNIP_KY03
Donaldson Forest 100% temperate needleleaf evergreen (LAI=2.85) obs_ERA NIP_LMDZ USNIP_FL14, USNIP_FL99
Anchorage 40% boreal needle-leaved evergreen (LAI=4), 60% boreal broad-leaved summergreen (LAI=4.5) ERA NIP_LMDZ Bethel, USNIP_SOGR_10, USNIP_CA45
Mitra 50% temperate broad-leaved evergreen (LAI=2), 50% C3 grass (LAI=0.4) obs_ERA NIP_LMDZ Beja, Faro, Penhas, Mitra, Portoallegre
Bily Kriz 100% temperate needleleaf evergreen (LAI=7.5) obs_ERA NIP_LMDZ Vienna, Podersdorf, Apetlon, Liptovsky, Krakow
Brloh 100% temperate broad-leaved summergreen (LAI=4.5) ERA NIP_LMDZ Leipzig, Hohhohensaas, Regensburg, Vienna, Petzenkirchen
Hainich 80% temperate broad-leaved summergreen (LAI=4.5), 20% C3 grass (LAI=0.4) obs_ERA NIP_LMDZ Leipzig, Hohhohensaas, Braunschweig, BadSalzuflen, Wuerzburg, Wasserkuppe
Tharandt 80% temperate needleleaf evergreen (LAI=4), 20% C3 grass (LAI=0.4) obs_ERA NIP_LMDZ Leipzig, Berlin, Hohhohensaas, Regensburg

Table 3: Information on the offline simulations performed on the 10 sites listed in Table 1: meteorological forcing (6 hourly observations of temperature, humidity, winds, precipitation and radiative fluxes), isotopic forcing (monthly isotopic composition of the precipitation and near-surface water vapor), and prescribed vegetation type and LAI (leaf area index) properties. We give proportions (in %) of the total vegetated area, excluding bare soil. For example, if a given vegetation type covers 100% of the vegetated area and the bare soil fraction is 30%, then the vegetation type covers only 70% of the total area. Three kinds of meteorological forcing are possible: meteorological observations only (obs), meteorological observations filled with ERA-Interim for missing variables (obs_ERA) or ERA-Interim (ERA). Two kinds of isotopic forcing are possible: isotopic composition of precipitation and water vapor observed on the site (obs_iso), or interpolation between GNIP, USNIP or Carbo-Europe stations using the LMDZ atmospheric general circulation model. In the former case, the datasets used for prescribing the water vapor and precipitation isotopic composition forcing are mentionned. In the latter case, GNIP, USNIP or Carbo-Europe stations used to construct the interpolated precipitation isotopic composition forcing are listed.

ORCHIDEE offline simulations require as forcing several meteorological variables: near-surface temperature, humidity and winds, surface pressure, precipitation, downward longwave and shortwave radiation fluxes. At Le Bray and Yatir, we use local meteorological measurements available at hourly time scale. At other sites, we use local meteorological measurements when available and combine them with ERA-Interim reanalyses at 6-hourly time scale for missing variables. At other sites, no nearby meteorological measurements are available and only ERA-Interim reanalyses [94] are used (Table 3).

At each site, we run the model three times over the first year of isotopic measurement (e.g., 2007 at Le Bray). These three years are discarded as spin-up. Then we run the model over the full period of isotopic measurements (e.g., 2007-2008 at Le Bray). We checked that at all sites, the seasonal distribution of δ18O , which is the slowest variable to spin-up, is identical between the last year of spin-up and the following year.

We force ORCHIDEE with monthly isotopic composition of precipitation and near-surface water vapor. Since we evaluate the results at the monthly time scale, we assume that monthly isotopic forcing is sufficient. At Le Bray and Yatir, monthly observations of isotopic composition of precipitation and near-surface water vapor are available to construct the forcing. Unfortunately, these observations are not available on the other sites. Therefore, we create isotopic forcing using isotopic measurements in the precipitation performed on nearby GNIP or USNIP stations. To interpolate between the nearby stations, we take into account spatial gradients and altitude effects by exploiting outputs from an LMDZ simulation.

Model-data comparison methods

Simulated isotopic composition in soil, stem and leaf water: The soil profile option is activated in all our stand-alone ORCHIDEE simulations. We compare the soil water samples collected in the first 15 cm of the soil (in the first 5-10 cm at many sites) to the soil water composition simulated in the uppermost layer.

The observed composition of stem water is compared to the simulated composition of the transpiration flux.

When comparing observed and simulated composition of leaf water, the Peclet effect, which mixes stomatal water with xylem water (Equation 8), is deactivated. Neglecting the Peclet effect may lead to overestimate of δ18Oleaf values.

Impact of the temporal sampling: Over the ten sites, samples were collected during specific days and hours. This temporal sampling may induce artifacts when comparing observations to monthly-mean simulated ORCHIDEE values. For soil and stem water, the effect of temporal sampling can be neglected because simulated soil and stem water composition vary at a very low frequency. For leaf water however, there are large diurnal variations [95]. For example, if leaf water is sampled every day at noon when δ18Oleaf is maximum, then observed δ18Oleaf will be more enriched than monthly-mean δ18Oleaf . The exact sampling time is available for Le Bray site only, where we will estimate the effect of temporal sampling.

Spatial heterogeneities: We are aware of the scale mismatch between punctual in-situ measurements and an LSM designed for large scales (a typical GCM grid box is more than 100 km wide). However, for soil moisture it has been shown that local measurements represent a combination of small scale (10-100 m) variability [20,21] and a large-scale (100-1000 km) signal [96] that a large-scale model should capture [97]. The sampling protocol allows us to evaluate the spatial heterogeneities. For example at Le Bray, two samples were systematically taken a few meters apart, allowing us to calculate the difference between these two samples. On average over all months, the difference between the two samples is 3.5‰ for δ18Os, 4.8‰ for δ18Ostem and 1.3‰ for δ18Oleaf . At Yatir, samples were taken several days every month, allowing us to calculate a standard deviation between the different samples for every month. On average of all months, the standard deviation is 0.9‰ for δ18Os, 0.4‰ for δ18Ostem and 1.2‰ for δ18Oleaf . These error bars need to be kept in mind when assessing modeldata agreement.

Soil moisture: Soil moisture have a different physical meaning in observations and model. Soil moisture is measured as volumetric soil water content (SWC) and expressed in %. In ORCHIDEE, the soil moisture is expressed in mm and cannot be easily converted to volumetric soil water content: the maximum soil water holding capacity of 300 mm and soil depth of 2 m are arbitrary choices and do not reflect realistic values at all sites. In LSMs, soil moisture is more an index than an actual soil moisture content [3]. In this version of ORCHIDEE in particular, it is an index to compute soil water stress, but it was not meant to be compared with soil water content measurements. Therefore, to compare soil moisture between model and observations, we normalize values to ensure that they remains between 0 and 1. The observed normalized SWC is calculated as hydrology-current-research where SWCmin and SWCmax are the minimum and maximum observed values of monthly SWC at each site. Similarly, simulated normalized SWC is calculated as hydrology-current-research where SWCmin and SWCmax are the minimum and maximum simulated values of monthly SWC at each site.

Evaluation at measurement sites

In this section, we evaluate the simulated isotopic composition in different water reservoirs of the soil-vegetation-atmosphere continuum at the seasonal scale.

Hydrological simulation: Before evaluating the isotopic composition of the different water reservoirs, we check whether the simulations are reasonable from a hydrological point of view. ORCHIDEE captures reasonably well the magnitude and seasonality of the latent and sensible heat fluxes at most sites (Figures 3 and 4). At Le Bray for example, the correlation between monthly values of evapotranspiration is 0.98 and simulated and observed annual mean evapotranspiration rates are 2.4 mm/d and 2.0 mm/d respectively. However, the model tends to overestimate the latent heat flux at the expense of the sensible heat flux at several sites. This is especially the case at the dry sites Mitra and Yatir: the observed evapo-transpiration is at its maximum in spring and then declines in summer due to soil water stress. ORCHIDEE underestimates the effect of soil water stress on evapo-transpiration and maintains the evapo-transpiration too strong throughout the summer.


Figure 3: Evaluation of hydrological and isotopic variables simulated by ORCHIDEE on different MIBA or Carbo-Europe sites. a, d, g, j, m: latent (green) and sensible (red) heat fluxes observed locally when available (circles), simulated in the ERA-Interim reanalyses (stars) and simulated by ORCHIDEE (lines). b, e, h, k, n: normalized soil moisture content (SWC, without unit) observed locally (circles) and simulated by ORCHIDEE (lines). c, f, i, l, o: δ18O of the surface soil (brown) and stems (green) simulated by ORCHIDEE in the control offline simulations (thin curves) and observed (circles). Observed δ18O in precipitation (thick dashed red) and vapor (thick dashed blue) used as forcing are also shown. a-c: Le Bray, d-f: Yatir, g-i: Morgan-Monroe, j-l: Donaldson Forest, m-o: Anchorage. The normalized SWC (soil water content) is calculated.


Figure 4: Same as Figure 3 but for Mitra (a-c), Bily Kriz (d-f), Brloh (g-i), Hainich (j-l: Donaldson Forest), and Tharandt (m-o)

The soil moisture seasonality is very well simulated at all sites where data is available (Figures 3 and 4), except for a two-month offset at Yatir (Figure 3f).

Water isotopes in the soil water: The evaluation of the isotopic composition of soil water is crucial before using ORCHIDEE to investigate the sensitivity to the evapo-transpiration partitioning or to infiltration processes, or in the future to simulate the isotopic composition of paleo-proxies such as speleothems [98].

In observations, at all sites, δ18Os remains close to δ18Op, within the relatively large month-to-month noise and spatial heterogeneities (Figures 3 and 4) At most sites (Le Bray, Donaldson Forest, Anchorage, Bily Kriz and Hainich), observed δ18Os exhibits no clear seasonal variations distinguishable from month-to-month noise. At Morgan- Monroe and Mitra, and to a lesser extent at Brloh and Tharandt, δ18Os progressively increases throughout the spring, summer and early fall, by up to 5‰ at Morgan-Monroe. The increase in δ18Os in spring can be due to the increase in δ18Op. The increase in δ18Os in late summer and early fall, while δ18Op starts to decrease, is probably due to the enriching effect of bare soil evaporation. At Yatir, δ18Os increases by 10‰ from January to June, probably due to the strong evaporative enrichment on this dry site. Then, the δ18Os starts to decline again in July. This could be due to the diffusion of depleted atmospheric water vapor in the very dry soil.

ORCHIDEE captures the order of magnitude of annual-mean δ18Os on most sites, and captures the fact that it remains close to δ18Op. ORCHIDEE captures the typical δ18Os seasonality, with an increase in δ18Os in spring-summer at Morgan-Monroe, Donaldson Forest, Mitra and Bily Kriz. However, the sites with a spring-summer enrichment in ORCHIDEE are not necessarily those with a spring-summer enrichment in observations. This means that ORCHIDEE misses what controls the inter-site variations in the amplitude of the δ18Os seasonality. The seasonality is not well simulated at Yatir. This could be due to the missed seasonality in soil moisture and evapo-transpiration. This could be due also to the fact that at Yatir ORCHIDEE underestimates the proportion of bare soil evaporation to total evapo-transpiration: less than 10% in ORCHIDEE versus 38% observed [89], which could explain why the spring enrichment is underestimated. Besides, ORCHIDEE does not represent the diffusion of water vapor in the soil, which could explain why the observed δ18Os decrease at Yatir in fall is missed.

When comparing the different sites, annual-mean δ18Os follows annual-mean δ18Op, with an inter-site correlation of 0.99 in observations. Therefore, it is easy for ORCHIDEE to capture the intersite variations in annual-mean δ18Os. A more stringent test is whether ORCHIDEE is able to capture the inter-site variations in annual-mean δ18Os - δ18Op. This is the case, with a correlation of 0.85 (Figure 5a) between ORCHIDEE and observations. In ORCHIDEE (and probably in observations), spatial variations in δ18Os - δ18Op are associated with the relative importance of bare soil evaporation.


Figure 5: a) Relationship between simulated and observed annual-mean δ18O in the soil water (red), stem water (blue) and leaf water (green), to which the precipitation-weighted annual-mean precipitation δ18O is subtracted. In the case of perfect model-data agreement, markers should fall on the y=x line. b) Relationship between the annual-mean δ18O in the soil water and in stem water, to which the precipitation-weighted annual-mean precipitation δ18O is subtracted, for both ORCHIDEE (magenta) and observations (cyan). When soil and stem water share the same δ18O, they fall on the y=x line.

Water isotopes in the stem water: In observations, observed δ18Ostem exhibits no seasonal variations distinguishable from month-tomonth noise (Figures 3 and 4). At Le Bray, Yatir, Mitra, Brloh, Hainich, observed δ18Ostem is more depleted than the surface soil water. It likely corresponds to the δ18O values in deeper soil layers, suggesting that the rooting system is quite deep. For example, at Mitra, the root system reaches least 6 m deep, and could at some places reach as deep as 13 m where it could use depleted ground water. At Donaldson Forest, Morgan-Monroe, Anchorage and Tharandt, δ18Ostem is very close to δ18Os, maybe reflecting small vertical variations in isotopic composition within the soil or shallow root profiles.

At Bily Kriz, observed δ18Ostem is surprisingly more enriched than surface soil water. Several hypotheses could explain this result: (1) the surface soil water could be depleted by dew or frost at this mountainous, foggy site; (2) spruce has shallow roots and therefore sample soil water that is not so depleted; (3) the twigs that were sampled were relatively young so that evaporation from their surface could have occurred when they were still at tree; (4) twigs were sampled in sun-exposed part of the spruce crowns during sunny conditions, which could favor some evaporative enrichment. Additional measurements show a lower Deuterium excess in the stem water compared to the soil water, supporting evaporative enrichment of stems.

ORCHIDEE captures the fact that δ18Ostem is nearly uniform throughout the year. As for soil water, it is easy for ORCHIDEE to capture the inter-site variations in annual-mean δ18Ostem (inter-site correlation between ORCHIDEE and observations of 0.90). ORCHIDEE is able to capture some of the inter-site variations in annual-mean δ18Ostem - δ18Op, with a inter-site correlation between ORCHIDEE and observations of 0.60. However, ORCHIDEE simulates δ18Ostem values that are very close to δ18Os values (Figure 5b). It is not able to capture δ18Ostem values that are either more enriched or more depleted than δ18Os. This could be due to the fact that ORCHIDEE underestimates vertical variations in soil isotopic composition. Also, ORCHIDEE is not designed to represent deep ground water sources or photosynthesizing twigs.

Vertical profiles of soil water isotope composition: At Le Bray, we compare our offline simulation for 2007 with soil profiles collected from 1993 to 1997 and in 2007 (Figure 6a-6b). The year mismatch adds a source of uncertainty to the comparison. In summer (profiles of August 1993 and September 1997), the data exhibits an isotopic enrichment at the soil surface of about 2.5‰ compared to the soil at 1 m depth (Figure 6a), likely due to surface evaporation [99]. Then, by the end of September 1994, the surface becomes depleted, likely due to the input of depleted rainfall. Previously enriched water remains between 20 and 60 cm below the ground, suggesting an infiltration through piston-flow [100]. ORCHIDEE predicts the summer isotopic enrichment at the surface, but slightly later in the season (maximum in September rather than August) and underestimates it compared to the data (1.5‰ enrichment compared to 2.5‰ observed, Figure 6b). The model also captures the surface depletion observed after the summer, as well as the imprint of the previous summer enrichment at depth. However, ORCHIDEE simulates the surface depletion in December, whereas the surface depletion can be observed sooner in the data, at the end of September 1994.


Figure 6: Vertical profiles of soil δ18O measured (a,c) and simulated by ORCHIDEE for the control offline simulations (b,d) on the Bray site (a,b) and the Yatir sites (b,d). Beware that the y-scales for observations and simulations are different. This is because the representation of the soil water content is very rudimentary in the ORCHIDEE model, preventing any quantitative comparison of measured and simulated soil depth. The horizontal black dashed line represents the bottom of the observed profiles. Model outputs are sampled at the same time as the data. For the Yatir sites, frequent soil sampling for the same year allowed us plot representative bi-monthly averages for both measured and simulated profiles. This could not be the case for Le Bray. Some soil profiles were observed at Le Bray in 2007, but we do not show them because they are limited to the top 24 cm of the soil only.

At Yatir, observed profiles exhibit a strong isotopic enrichment from deep to shallow soil layers in May-June by up to 10‰ (Figure 6c). As for Le Bray, the model captures but underestimates this isotopic enrichment in spring and summer by about 3‰ (Figure 6d). This discrepancy could be the result of underestimated bare soil evaporation. Observed profiles also feature a depletion at the surface in winter that the model does not reproduce. This depletion could be due to back-diffusion of depleted vapor in dry soils [99,101-103], a process that is not represented in ORCHIDEE but likely to be significant in this region. Soil evaporation fluxes measured with a soil chamber at Yatir shows that when soils are dry, there is adsorption of vapor from the atmosphere to the dry soil pores before sunrise and after sunset [104].

Water isotopes in leaf water: It is important to evaluate the simulation of the isotopic composition of leaf water by ORCHIDEE if we want to use this model in the future for the simulation of paleoclimate proxies such tree-ring cellulose [105,106], for the simulation of the isotopic composition of atmospheric CO2 which may be used to partition CO2 fluxes into respiration from vegetation and soil [107,108] or for the simulation of the isotopic composition of atmospheric O2 which may be used to infer biological productivity [109,110].

In the observations, δ18Oleaf exhibits a large temporal variability reflecting a response to changes in environmental conditions (e.g., relative humidity and the isotopic composition of atmospheric water vapor). At all sites except at Yatir, δ18Oleaf is most enriched in summer than in winter, by up to 15‰ (Figures 3 and 4). This is because the evaporative enrichment is maximum in summer due to drier and warmer conditions.

ORCHIDEE captures the maximum enrichment in summer. However, ORCHIDEE underestimates the annual-mean δ18Oleaf at most sites (Figure 5). This could be due to the fact that most leaf samples were collected during the day, when the evaporative enrichment is at its maximum, while for ORCHIDEE we plot the daily-mean δ18Oleaf . At Le Bray, if we sample the simulated δ18Oleaf during the correct days and hours, simulated δ18Oleaf increases by 4‰ in winter and by 10‰ in summer. Such an effect can thus quantitatively explain the model-data mismatch. After taking this effect into account, simulated δ18Oleaf may even become more enriched than observed. This is the case at Le Bray, especially in summer. The overestimation of summer δ18Oleaf could be due to neglecting diffusion in leaves or non-steady state effects.

Again, Yatir is a particular case. Minimum δ18Oleaf occurs in spring-summer while the soil evaporative enrichment is maximum. In arid regions and seasons, leaves may close stomata during the most stressful periods of the day, inhibiting transpiration, and thus retain the depleted isotopic signal associated with the moister conditions of the morning [111,112]. ORCHIDEE does not represent this process and thus simulates too enriched δ18Oleaf .

Summary: Overall, ORCHIDEE is able to reproduce the main features of the seasonal and vertical variations in soil water isotope content, and seasonal variations in stem and leaf water content. Discrepancies can be explained by some sampling protocols, by shortcomings in the hydrological simulation or by neglected processes in ORCHIDEE (e.g., fractionation in the vapor phase).

The strong spatial heterogeneity of the land surface at small scales does not prevent ORCHIDEE from performing reasonably well. This suggests that in spite of some small-scale spatial heterogeneities at each site, local isotope measurements contain large-scale information and are relevant for the evaluation of large-scale LSMs.

Sensitivity analysis

Sensitivity to evapo-transpiration partitioning: Several studies have attempted to partition evapo-transpiration into the transpiration and bare soil evaporation terms at the local scale [29-31,113]. Estimating E/ET, where E is the bare soil evaporation and ET is the evapo-transpiration, requires measuring the isotopic composition of soil water, stem water and of the evapo-transpiration flux. The isotopic composition of the evapo-transpiration can be estimated through “Keeling plots” approach [114], but this is costly [29] and the assumptions underlying this approach are not always valid [115].

Considering a simple soil water budget at steady state and with vertically-uniform isotopic distribution, we show that although estimating E/ET requires measuring the isotopic composition of the evapo-transpiration flux, estimating E/I (where I is the precipitation that infiltrates into the soil) requires measuring temperature, relative humidity (h) and the isotopic composition of the soil water (δ18Os), water vapor (δ18Ov) and precipitation (δ18Op) only. Such variables are available from several MIBA and Carbo-Europe sites. More specifically, E/I is proportional to δ18Op - δ18Os:

hydrology-current-research (1)

where αeq and αk are the equilibrium and kinetic fractionation coefficients respectively.

Below, we show that this equation can apply to annual-mean quantities, neglecting effects associated with daily or monthly covariations between different variables. We investigate to what extent this equation allows us to estimate the magnitude of E/I at local sites.

At the Yatir site, all the necessary data for equation 1 is available. An independent study has estimated E/I =38% [89]. Using annually averaged observed values (δ18Op=-5.1‰ and δ18Os=-3.7‰ in the the surface soil), we obtain E/I =46%. However, in ORCHIDEE, the annually averaged surface δ18Os is 0.8 lower when sampled at the same days as in the data. When correcting for this bias, we obtain E/I =28%. Observed E/I lies between these two estimates. This shows the applicability of this estimation method, keeping in mind that estimating E/I is the most accurate where E/I is lower.

When we perform sensitivity tests to ORCHIDEE parameters at the various sites, the main factor controlling δ18Os is the E/I fraction. This is illustrated as an example at Le Bray and Mitra sites (Figure 7). Sensitivity tests to parameters as diverse as the rooting depth or the stomatal resistance lead to changes in δ18Os - δ18Op and in E/I that are very well correlated, as qualitatively predicted by equation 13. This means that whatever the reason for a change in E/I, the effect on δ18Os - δ18Op is very robust.


Figure 7: Isotopic difference between soil water and precipitation (δ18Os- δ18Op) as a function of E/I (fraction of the infiltrated water that evaporates at the bare soil surface), for different sensitivity tests in ORCHIDEE. a) at Le Bray and b) at Mitra. All values are annual means. The horizontal dashed line represents the observed values for δ18Os- δ18Op. The orange dashed line shows the best linear fit between the different sensitivity tests.

Quantitatively, the slope of δ18Os - δ18Op as a function of E/I among the ORCHIDEE tests is of 0.78‰/% (r=0.94, n=6) at Le Bray and of 0.25‰/% (r=0.999, n=5) at Mitra, compared to about 0.25-0.3‰/% predicted by equation 13. The agreement is thus very good at Mitra. The better agreement at Mitra is because it is a dry site where E/I varies greatly depending on sensitivity tests. In contrast, Le Bray is a moist site where E/I values remains small for all the sensitivity tests, so numerous effects other than E/I and neglected in equation 13 can impact δ18Os - δ18Op.

To summarize, local observations of δ18Os - δ18Op could help constrain the simulation of E/I in models. This would be useful since the evapo-transpiration partitioning has a strong impact on how an LSMs represents land-atmosphere interactions [116].

Sensitivity to soil infiltration processes: Partitioning between evapo-transpiration, surface runoff and drainage depends critically on how precipitation water infiltrates the soil [5,8,10], which is a key uncertainty even in multi-layer soil models where infiltration processes are represented explicitly [62]. It has been suggested that observed isotopic profiles could help understand infiltration processes at the local scale [100]. The capacity of ORCHIDEE to simulate soil profile allows us to investigate whether measured isotope profiles in the soil could help evaluate the representation of these processes also in largescale LSMs.

With this aim, we performed sensitivity tests at Le Bray. The simulated profiles are sensitive to vertical water fluxes in the soil. When the diffusivity of water in the soil column is decreased by a factor 10 from 0.1 to 0.01 compared to the control simulation, the deep soil layer becomes more depleted by about 0.7‰ (Figure 8) and the isotopic gradient from soil bottom to top becomes 30% steeper in summer, because the enriched soil water diffuses slower through the soil column.


Figure 8: Sensitivity of simulated δ18Os profiles to the parameterization of infiltration processes in the soil at Le Bray. July (a) and December (b) are shown for three different parameterizations in offline simulations: control simulation (solid red), a simulation in which the soil water diffusivity was divided by 10 (dashed blue) and a simulation is which the water infiltrates the soil uniformly in the vertical (crude representation of preferential pathways, dash-dotted green) rather than in a piston-like way as is the case for other simulations.

Simulated profiles are also sensitive to the way precipitation infiltrates the soil. When precipitation is added only to the top layer (piston-flow infiltration) the summer enrichment is reduced by mixing of the surface soil water with rainfall, and it propagates more easily to lower layers during fall and winter. Conversely, when rainfall is evenly spread throughout the soil column (a crude representation of preferential pathway infiltration), the surface enrichment is slightly more pronounced and the deep soil water is more depleted by up to 0.8‰ in winter (Figure 8). However, the observed surface depletion occurs in February with preferential pathways, compared to December in the piston-like in infiltration. The quick surface depletion observed after the summer suggests that infiltration is dominated by the pistonlike mechanisms.

To summarize, we show that vertical and seasonal variations of δ18Os are very sensitive to infiltration processes, and are a powerful tool to evaluate the representation of these processes in LSMs.

Global-scale Simulations Using the Coupled LMDZORCHIDEE Model

Simulation set-up

To compare with global datasets, we performed LMDZ-ORCHIDEE coupled simulations. In all our experiments, LMDZ three-dimensional fields of horizontal winds are nudged towards ECMWF (European Center for Medium range Weather Forecast) reanalyses [117]. This ensures a realistic simulation of the large-scale atmospheric circulation and allows us to perform a day-to-day comparison with field campaign data [60,118]. At each time step, the simulated horizontal wind field hydrology-current-research is relaxed towards the reanalysis following this equation:


where hydrology-current-research is the reanalysis horizontal wind field,hydrology-current-research is the effect of all simulated dynamical and physical processes on hydrology-current-research , and τ is a time constant set to 1 h in our simulations [119].

To compare with global datasets, LMDZ-ORCHIDEE simulations are performed for the year 2006, chosen arbitrarily. We are not interested in inter-annual variations and focus on signals that are much larger. To ensure that the water balance is closed at the annual scale, we performed iteratively 10 times the year 2006 as spin-up. In these simulations, the Peclet and non-steady state effects are de-activated.

To compare with field campaign observations in 2002 and 2005, we use simulations performed for these specific years, initialized from the 2006 simulation. In these simulations, we test activating or deactivating the Peclet effect.

In all LMDZ-ORCHIDEE simulations, canopy-interception was de-activated (consistent with simulations that our modeling group performed for the Fourth Assessment Report).

Evaluation of water isotopes in leaf water at the diel scale during campaign cases

Daily data from field campaigns: Two field campaigns are used to evaluate the representation of δ18Oleaf diurnal variability. The first campaign covers six diurnal cycles in May and July 2002 in a grassland prairie in Kansas (39.20 ° N 96.58 ° W [120]). The second campaign covers four diurnal cycles in June 2005 in a pine plantation in Hartheim, Germany (7.93 ° N, 7.60 ° E [121]).

Because meteorological and isotopic forcing are not available for the entire year, we prefer to compare these measurements with LMDZORCHIDEE simulations. At both sites, the simulated δ18Ov and δ18Ostem are consistent with those observed (model-data mean difference lower than 1.4‰ in Kansas and 0.4‰ at Hartheim), allowing us to focus on the evaluation of leaf processes.

Evaluation results: At the Kansas grassland site, δ18Oleaf exhibits a diel cycle with an amplitude of about 10‰ [120]. LMDZ-ORCHIDEE captures this diel variability, both in terms of phasing and amplitude (Figure 9). The model systematically overestimates δ18Oleaf by about 4‰, in spite of the underestimation of the stem water by 1.4‰ on average. This may be due to a bias in the simulated relative humidity (LMDZ is on average 13% too dry at the surface, which translates into an expected enrichment bias of 3.9‰ on the leaf water assuming steady state based on Equation 7) or to uncertainties in the kinetic fractionation during leaf water evaporation.


Figure 9: δ18O of stem and grass leaves measured during two series of 3 diurnal cycles in May and July 2002 over the plains of Kansas [120] and simulated by LMDZ-ORCHIDEE for the same year in the grid box containing the observation site. δ18O of vapor (blue), pine leaves (pink and red) and stems (green) measured during four diurnal cycles in June 2005 in Hartheim, Germany [121] and simulated by LMDZ-ORCHIDEE for the same year in the grid box containing the observation site. Simulated values are dashed, observed values solid. Two kinds of leaves were sampled during this campaign: one-year-old leaves (solid pink) and current-year leaves (solid brown). Two leaf water diagnostics were computed for in LMDZ-ORCHIDEE: stationary state at the evaporative site (dashed red, equation 7) or non-stationary state in the lamina, taking into account the Peclet effect (dashed brown, equation 9, using an effective length scale of 25 mm).

At the Hartheim pine plantation, δ18Oleaf is on average 8‰ more depleted for current-year needles than for 1-year-old needles. Also, the observed diel amplitude is weaker for current-year needles (5 to 8‰) than for 1-year-old needles (10 to 15‰). These observations are consistent with a longer diffusion length for current-year needles (15 cm) than for 1-year-old needles (5 cm) [121] and with a larger transpiration rate, leading to a stronger Peclet effect. When neglecting Peclet and non-steady state effects, ORCHIDEE simulates an average δ18Oleaf close to that of 1-year-old needles, consistent with the small diffusion length and evaporation rate of these leaves. ORCHIDEE captures the phasing of the diurnal cycle, but underestimates the diel amplitude by about 4‰. This is probably due to the underestimate of the simulated diel amplitude of relative humidity by 20%. Accounting for Peclet and non-steady state effects strongly reduces both the average δ18Oleaf and its diel amplitude (Figure 9), in closer agreement with current-year needles.

To summarize, ORCHIDEE simulates well the leaf water isotopic composition. The leaf water isotope calculation based on Craig et al. [75] simulates the right phasing and amplitude for leaves that have short diffusive lengths or low transpiration rates. Non-steady state and diffusion effects need to be considered in other cases. By activating or de-activating these effects, ORCHIDEE can simulate all cases.

Evaluation of water isotopes in precipitation

Precipitation datasets: To evaluate the spatial distribution of precipitation isotopic composition simulated by the LMDZORCHIDEE coupled model, we use data from the Global Network for Isotopes in Precipitation (GNIP [122]), further complemented by data from Antarctica [123] and Greenland [124]. We also use this network to construct isotopic forcing at sites where the precipitation was not sampled, complemented with the USNIP (United States Network for Isotopes in Precipitation [125]) network.

Evaluation results: At the global scale, the LMDZ-ORCHIDEE coupled model reproduces the annual mean distribution in δ18Op and dp observed by the GNIP network reasonably well (Figure 10), with correlations of 0.98 and 0.46 and Root Mean Square Errors (RMSE) of 3.3‰ and 3.5‰ respectively.


Figure 10: a) Annual mean δ18Op from GNIP [122], Antarctica [123] and Greenland [124] data. The data is gridded over a coarse 7.5 × 6.5° grid for visualization purposes. b) Same as a) but for annual mean dp. c) Annual mean δ18Op simulated by coupled LMDZ-ORCHIDEE model for the control simulation. d) same as c) but for annual mean dp.

This good model-data agreement can be obtained even when we deactivate ORCHIDEE. When we use LMDZ in a stand-alone mode, in which the isotope fractionation at the land surface is neglected [60], the model-data agreement is as good as when we use LMDZ-ORCHIDEE. Therefore, fractionating processes at the land surface have a second order effect on precipitation isotopic composition, consistent with [28,71-73].

To quantify in more detail the effect of fractionation at the land surface, we performed additional coupled simulations with LMDZORCHIDEE. We compare the control simulation described above (ctrl) to a simulation in which fractionation at the land surface was deactivated (nofrac) (Figure 11). In nofrac, the composition of bare soil evaporation equals that of soil water. Even when restricting the analysis to continental regions, the spatial correlations between the ctrl and nofrac simulations are 0.999 and 0.95 for δ18Op and dp respectively, and the root mean square differences are 0.27‰ and 1.1‰ for δ18Op and dp respectively. This confirms that fractionation at the land surface has a second-order effect on precipitation isotopic composition compared to the strong impact of atmospheric processes.


Figure 11: a) Annual-mean δ18Op in the ctrl simulation (LMDZ-ORCHIDEE) minus annual mean δ18Op in the nofrac simulation (LMDZ-ORCHIDEE in which the isotopic fractionation was de-activated during bare soil evaporation). This shows the effect of isotopic fractionation at the soil surface on δ18Op. b) Same as a) but for dp.

However, to second order, a detailed representation of fractionation at the land surface lead to a slight improvement in the simulation of δ18Op and to a significant improvement in that of dp . In ctrl, δ18Op is lower by up to 1.5‰ and dp higher by up to 5‰ than in nofrac over boreal continental regions such as Siberia, Canada and central Asia, consistent with the expected effect of fractionation at surface evaporation [42]. Taking into account fractionation at the land surface leads to a better agreement with the GNIP data over these regions, where δ18Op is overestimated by about 4‰ and dp underestimated by 4 to 7‰ when neglecting fractionation at the land surface. The effect of fractionation is maximal over these boreal regions because (1) the fraction of bare soil evaporation is maximal, (2) a significant proportion of evaporativelyenriched soil water is lost by drainage and (3) a larger proportion of the moisture comes from land surface recycling [48,126,127]. Similar results were obtained with other models [128].

To summarize, LMDZ-ORCHIDEE simulates well the spatial distribution of precipitation isotopic composition, but this distribution is not a very stringent test for the representation of land surface processes in ORCHIDEE. In the next section, we argue that the distribution of river isotopic composition is a more stringent test.

Evaluation of water isotopes in river water: Large rivers integrate a wide range of hydrological processes at the scale of GCM grid boxes [22,23,129-131]. Here we evaluate the isotopic composition of river water simulated by ORCHIDEE using data collected by the Global Network for isotopes in Rivers (GNIR [132,133]).

Observed annual mean δ18Oriver follows to first order the isotopic composition of precipitation [79], and is thus also well simulated by LMDZ-ORCHIDEE (Figure 12a and 12b), with a spatial correlation between measured and simulated δ18Oriver of 0.80 and a RMSE of 3.2‰ over the 149 LMDZ grid boxes containing data. Regionally however, the δ18O difference between precipitation and river water (δ18Oriver - δ18Op) can be substantial and provides a stronger constraint for the model. Over South America, Europe and some parts of the US, the river water is typically 1‰ to 4‰ more depleted than the precipitation (Figure 12a), because precipitation contributes more to rivers during seasons when it is the most depleted [134]. In contrast, over central Asia or northern America, river water is more enriched than precipitation, due to evaporative enrichment of soil water [79,134,135]. This is further confirmed by a simulation where fractionation at the land surface was neglected (not shown), for which the river water is in global average 5‰ more depleted.


Figure 12: (a) Annual mean δ18Op in rivers (δ18Oriver) measured from the GNIR database. (b) Simulated by LMDZ-ORCHIDEE for the control simulation. (c) Annual mean δ18Oriver18Op observed from the GNIR and GNIP databases. (d) Simulated by LMDZ-ORCHIDEE for the control simulation. On sub-plots d and f the United States, where the GNIR network is the densest, are enlarged for better readability.

ORCHIDEE reproduces moderately well the magnitude and patterns of δ18Oriver - δ18Op, with a spatial correlation of 0.39 and a RMSE of 2.7‰ over the 22 LMDZ grid boxes that contain δ18Oriver observations. It simulates the negative values over the western US, Europe and South America and the positive value over Mongolia. However, the model does not capture the positive δ18Oriver - δ18Op in Eastern US, though positive values are simulated further North. This suggests that such a diagnostic may help identify biases in the representation of the soil water budget, as discussed in the following section.

Sensitivity to the representation of pathways from precipitation to rivers

At the local scale, water isotopes have already been used to partition river discharge peaks into the contributions from recent rainfall and soil water [37-39]. Given the property of rivers to integrate hydrological processes at the basin scales [22,23,129-131], we now explore to what extent δ18Oriver could help evaluate pathways from precipitation to rivers in LSMs. We illustrate this using seasonal variations in δ18Oriver on two well established GNIR and GNIP stations in Vienna (Danube river) and Manaus (the Amazon) (Figure 13). The seasonal cycle in δ18Oriver is attenuated compared to that in δ18Op, and δ18Oriver lags δ18Op (by 5 month at Vienna and 1-3 months at Manaus).


Figure 13: Seasonal variations in δ18Op (a,b) and δ18Oriver (c,d) observed (solid black) and simulated for the control LMDZ-ORCHIDEE simulation (dashed black) for (a,c) the Danube river in Vienna and (b,d) the Amazon river in the Manaus region (average over the 8 ° S-3 ° S-56 ° W 63 ° W domain). Also shown are δ18Oriver for simulations where the total runoff is partitioned into surface runoff only without drainage (dash-dotted blue) and where we multiplied by two the time residence in the reservoir collecting drainage in the routing scheme (dash-dotted red). Beware that the y-scale is different on the two sites. The difference in the annual-mean values between the two sites reflect the difference in the annual-mean δ18Op.

LMDZ-ORCHIDEE (control simulation) simulates qualitatively well the amplitude and the phasing observed in δ18Op and δ18Oriver . To understand better what determines the attenuation and lag of the seasonality in δ18Oriver compared to that in δ18Op, we perform sensitivity tests to ORCHIDEE parameters. Parameters tested include the partitioning of excess rainfall into surface runoff and drainage and the residence time scale of different reservoirs (slow, fast and stream) in the routing scheme. River discharge is extremely sensitive to these parameters [136].

If all the runoff occurs as surface runoff (Figure 13), then the seasonal cycle of δ18Oriver is similar to that of δ18Op. This shows that the attenuation and lag of the seasonality in δ18Oriver compared to that in δ18Op are caused by the storage of water into the slow reservoir, which accumulates drainage water.

When the residence time scale of the slow reservoir is multiplied by 2 (i.e., the water from the slow reservoir is poured twice faster into the streams, Figure 12), the simulated lag of δ18Oriver at Vienna increases from 4 to 5 months (in closer agreement with the data). In contrast, the seasonal cycle in δ18Oriver is not sensitive to residence time scales in the stream and fast reservoirs, which are too short to have any impact at the seasonal scale.

To summarize, ORCHIDEE performs well in simulating the seasonal variations in δ18Oriver . In turn, δ18Oriver observations could help estimate the proportion of surface runoff versus drainage and calibrate empirical residence time constants in the routing scheme, offering a mean to enhance model performance.

Evapo-transpiration partitioning

In this section, we generalize at the global scale our results on evapo-transpiration partitioning estimates.

We apply equation 1 to annual-mean outputs from a LMDZORCHIDEE simulation. We compare E/I estimated from Equation 1 to E/I directly simulated by LMDZ-ORCHIDEE. The spatial pattern of E/I is remarkably well estimated by Equation 1 (Figure 14). The equation captures the maximum over the Sahara, Southern South America, Australia, central Asia, Siberia and Northern America. The isotope-derived spatial distribution of E/I correlates well with the simulated distribution (r=0.91). Average errors are lower than 50% of the standard deviation at the global scale. This confirms that covariation between the different variables at sub-annual time scales has a negligible effect, so that the equation can be applied to annual-mean quantities. Generally, E/I estimates are best where E/I is relatively small.


Figure 14: a) annual mean E/I (proportion of infiltrating water recycled back to the atmosphere as bare soil evaporation) simulated by LMDZ-ORCHIDEE for the control simulation. b) E/I estimated from water isotopes measurements. We perform the estimations only on grid points where the denominators in the equation are different from 0 and where the soil water contents and the water fluxes whose compositions we need are strictly positive. Grid points where estimations cannot be performed are left white.

To test the effect of the assumption that the soil water isotopic composition is vertically constant, we applied Equation 1 using δ18Os - δ18Op from a simulation with soil profiles activated. This assumption is a significant source of uncertainty on estimating E/I (Table 4). We also analyzed the effect of potential measurement errors in δ18Os, δ18Op, δ18Ov temperature or relative humidity on the E/I reconstruction. Results are relatively insensitive to small errors in these measurements (Table 4). However, results are sensitive to the choice of the n exponent in the calculation of the kinetic fractionation αk (Table 4): knowing the n exponent with an accuracy of 0.07 (e.g., estimated n ranges from 0.63 to 0.70) is necessary to estimate E/I with an absolute precision of 2%.

Absolute or relative error RMS absolute error on rE/I RMS relative error on rE/I,
when rE/I> 4% (37% of total land aread)
soil profiles 12% 50%
ΔT =1°C 0.2% 1%
Δrh =1% 0.5% 1%
Δδp=1 3% 35%
Δδv=1 1% 8%
Δδs=1 5% 49%
Δn = 0.5 14% 52%

Table 4: Uncertainties in the estimation of E/I related to measurement errors and assumptions necessary in the simple conceptual model. Values give absolute (in ratio) and relative variations (in %) in estimated E/I when temperature T is modified by 1 ° C (line 4), when relative humidity rh is modified by 1% (line 5), when δ18O v , δ18O p and δ18O s are modified by 1, when n in the kinetic fractionation is varied from 0.5 to 1, and when the soil δ18O is not homogeneous vertically. The resulting variations in estimated E/I are averaged over all land grid points where the estimation could be performed.

Finally, estimating E/I using equation 1 bears additional sources of uncertainty in that we cannot estimate using the ORCHIDEE model. These are related to all processes that ORCHIDEE does not simulate. For example, ORCHIDEE underestimates or mis-represents the vertical isotopic gradients in soil water at some sites and does not represent the effect of water vapor diffusion in the soil. These effects may disturb the proportionality between E/I and δ18Os - δ18Op in practical applications.

To summarize, co-located isotope measurements in precipitation, vapor and soil water could provide an accurate constrain on the proportion of bare soil evaporation to precipitation infiltration.

Conclusion and Perspectives

The ORCHIDEE LSM, in which we have implemented water stable isotopes, reproduces the isotopic compositions of the different water pools of the land surface reasonably well compared to local data from MIBA and Carbo-Europe and to global observations from the GNIP and GNIR networks. Despite the scale mismatch between local measurements and a GCM grid box, and despite the strong spatial heterogeneity in the land surface, the capacity of ORCHIDEE to reproduce the seasonal and vertical variations in the soil isotope composition suggests that even local measurements can yield relevant information to evaluate LSMs at the large scale.

We show that the simulated isotope soil profiles are sensitive to infiltration pathways and diffusion rates in the soil. The spatial and seasonal distribution of the isotope composition of rivers is sensitive to the partitioning of total runoff into surface runoff and drainage and to the residence time scales in underground reservoirs. The isotopic composition of soil water is strongly tied to the fraction of infiltrated water that evaporates through the bare soil. These sensitivity tests suggest that isotope measurements, combined with more conventional measurements, could help evaluate the parameterization of infiltration processes, runoff parameterizations and the representation of surface water budgets in LSMs.

Evaluating an isotopic LSM requires co-located observations of the isotope composition in precipitation, vapor and soil at least at the monthly scale. However, such co-located measurements are still very scarce, and most MIBA and Carbo-Europe sites are missing one of the components. Therefore, for LSM evaluation purpose, we advocate for the development of co-located isotope measurements in the different water pools at each site, together with meteorological variables. Our results suggest that isotope measurements are spatially relatively well representative and that even monthly values are already valuable to identify model bias or to estimate soil water budgets. Therefore, in the perspective of LSM evaluation, if a compromise should be made with sampling frequency and spatial coverage, we favor co-located measurements of all the different water pools at the monthly scale on a few sites representative of different climatic conditions, rather than multiplying sites where water pools are not all sampled. Additionally, at each observation site, collecting different soil samples a few meters apart is helpful to check that they are spatial representative. In the future, development in laser technology [137,138] will allow the generalization of water vapor isotope monitoring at the different sampling sites, which has long been a very tedious activity [90].

From the modeling point of view, kinetic fractionation processes during bare soil evaporation are a source of uncertainty, and a better understanding and quantification of this fractionation is necessary [103,139]. In addition, the accuracy of isotopic simulations by LSM is expected to improve as the representation of hydrological processes improves. In particular, given the importance of vertical water exchanges for the isotopic simulation, implementing water isotopes in a multi-layer hydrological parameterization with sufficient vertical resolution [69] is crucial. In the future, we plan to implement water isotopes in the latest version of ORCHIDEE, which is multi-layer and more sophisticated [140-142]. Finally, latest findings largely based on water isotopic measurements suggest that different water pools coexist within a soil column and that evaporation, transpiration, runoff and drainage tap from these different pools [77,143,144]. These effects are not yet represented explicitly in global LSMs. These effects were mainly evidenced based on isotope measurements, and in turn, their representation expected to significantly impact isotopic simulations. Such feedbacks between isotopic research and hydrological parameterization improvements should lead to LSM improvements in the future. With this in mind, LSM inter-comparison projects would strongly benefit from including water isotopes as part of their diagnostics, in the lines of iPILSP (isotope counterpart of the Project for Intercomparison of Land-surface Parameterization Schemes [27]).

Representation of isotope fractionation during evaporation from land surface water pools

Processes for which we neglect fractionation: Snow sublimation is associated with a slight fractionation due to exchanges between snow and vapor in snow pores [115,145,146]. However, we assume that these effects are small enough to be neglected, as in other GCMs [58].

Water uptake by roots has been shown to be a non-fractionating process [147,148], but fractionation at the leaf surface during transpiration impacts the composition of transpired fluxes at scales shorter than daily [95,137]. As the application of ORCHIDEE in the context of our study focuses mainly on time scales of a month or longer, we assume here that the transpiration and stem water have the composition of soil water extracted by the roots.

Evaporation from bare soils and canopy-intercepted water: We represent isotope fractionation during evaporation of soil and canopyintercepted water using the model of Craig [75]: at any time t, the isotopic composition of evaporation RE is given by:

hydrology-current-research (2)

where Rl and Rv are the isotopic compositions of liquid water at the evaporative site and of water vapor respectively, h is the relative humidity normalized to surface temperature, αeq is the isotopic fractionation during liquid-vapor equilibrium [149] and αk is the kinetic fractionation during water vapor diffusion. The kinetic fractionation during soil evaporation is still very uncertain [103,150]. We use the very widespread formulation of [99,151]:

hydrology-current-research (3)

where D and Di are the molecular diffusivities of light and heavy water vapor in air, respectively, and n is an exponent that depends on the flow regime (0.5, 0.67 and 1 for turbulent, laminar and stagnant regimes respectively) but remains difficult to estimate [103,150]. In this study, we take n =0.67 for both evaporation of soil and canopyintercepted water, corresponding to moist conditions in the case of soils [99]. However, we also tried 0.5 and 1.0 to estimate the range of uncertainty related to this parameter. The isotopic composition of precipitation is only slightly sensitive to the formulation of the kinetic fractionation: when n varies from 0.5 to 1, significant changes in δ18Op and dp are restricted to areas where bare soil covers more than 70%. Even in those case, changes in δ18Op and dp never exceed 2‰ and 7‰ respectively. The impact is slightly stronger on soils. Varying n from 0.5 to 1 leads to δ18Os variations of 2‰ in offline simulations on the Bray site, of the order of the observed average difference between two samples collected on the same day (2.2‰). In coupled simulations, the impact on δ18Os and ds reaches 8‰ and 20‰ respectively on very arid regions such as the Sahara.

To calculate the temporal mean isotopic composition of evaporation over the time step Δt,hydrology-current-research, we assume Rv and h are constant throughout each time step. On the other hand, we allow the isotopic ratio of liquid water to vary over the simulation time step Δt following [151]. While assuming constant Rl is a valid assumption for models with very short time steps [152], it is not the case in ORCHIDEE (Δt =30min). We then calculate hydrology-current-research as:

hydrology-current-research (4)

where Rl0 is the initial isotopic ratio of liquid water, f is the remaining liquid fraction in the water reservoir affected by isotopic enrichment, and β and γ are parameters defined by Stewart [151]:




For canopy-intercepted water, the water reservoir is sufficiently small to assume that the water reservoir affected by isotopic enrichment is the total canopy-intercepted water. For soil evaporation on the other hand, we assume that the depth of the water reservoir affected by isotopic enrichment equals the average distance traveled by water molecules in the soil:

hydrology-current-research (5)

where KD is the effective self-diffusivity of liquid water in the soil column. Neglecting the dispersion term, KD is given by Munnich et al. [100,147,151-153]:

hydrology-current-research (6)

where hydrology-current-research is the molecular liquid water selfdiffusivity [154,155], τ is the soil tortuosity and θl is the volumetric soil water content. In the control simulation, we assume θl.τ =0.1 leading to L =0.67 mm. This choice is consistent with a τ of 0.67 [151] and an average θl of about 15%. At the Bray, measurements along profiles show θl varying from about 5 to 30%. Since these values are difficult to constrain observationally and very variable spatially and temporally, sensitivity tests to θl .τ are performed and described. We neglect the vapor phase in the soil and associated fractionation and diffusion processes [153].

Dew formation: We assume fractionation during dew and frost formation following a Rayleigh distillation of the vapor in the lowest 10 hPa (~ 80 m) of the atmosphere. Since the atmospheric water vapor condenses in small proportion during frost and dew, this choice of the depth of atmosphere involved in the condensation has almost no impact on the composition of the dew and frost formed. Following common practice, we use equilibrium fractionation coefficient from Merlivat et al. [148,156,157] and the kinetic fractionation formation of [158] with λ =0.004, whose choice has very little impact on the results.

Leaf water evaporation: At isotopic steady state, the composition of water transpired by the vegetation is equal to that of the soil water extracted by the roots. In default simulations, we assume that isotopic steady state for plant water is established at any time and we diagnose the composition of the leaf water at the evaporation site, SS e R , by inverting the Craig and Gordon equation [75]:

hydrology-current-research (7)

where Rs and Rv are the isotopic ratio in soil water and water vapor respectively, h is the relative humidity normalized to surface temperature, αeq is the isotopic fractionation during liquid-vapor equilibrium [148] and αk is the kinetic fractionation during water vapor diffusion. We take the same kinetic fractionation formulation as for the soil evaporation [150], with n =0.67 [31,69]. Leaf water compositions are significantly sensitive to parameter n, with variations of the order of 10‰ as n varies from 0.5 to 1. We assume that the leaf temperature used to calculate αeq is equal to the soil temperature, but results are very little sensitive to this assumption.

The isotopic composition of leaf water has been the subject of many observational and numerical modeling studies [86,159-161]. Several studies have shown that the composition of the leaves is affected by mixing with xylem water and by non-stationary effects [161,162]. Nonsteady state effects are also incorporated in ORCHIDEE following [159]. The isotopic ratio in the leaf mesophyll RLSS is the result of the mixing between leaf water at the evaporative site and xylem water (Peclet effect):

hydrology-current-research (8)

where f is a coefficient decreasing as the Peclet effect increases:


and p is the Peclet parameter [120,160]:


E is the transpiration rate per leaf area, Leff is the effective diffusion length and W is the leaf water content per leaf volume (assumed equal to 103 kg/m3, order of magnitude in [121]). The Peclet number P can be tuned by changing Leff , that depends on leaf geometry and drought intensity (e.g., 7 to 12 mm in Cuntz et al. [161], 50 to 150 mm in Barnard et al. [121]). We take Leff =8 mm to optimize our simulation on Hartheim.

For some simulations, we account for the effect of water storage in leaves (leading to some memory in the leaf water isotopic composition) following Dongmann [163]. Assuming that W is constant, we calculate the leaf lamina composition RL as Farquhar [159]:

hydrology-current-research (9)



and g is the sum of the total (stomatic and boundary layer) conductances. The isotopic composition of transpiration is then calculated so as to conserve isotope mass.

Representation of the vertical distribution of soil water isotopic composition

Principle: In control simulations, we assume that the isotopic composition of soil water is homogeneous vertically and equals the weighted average of the two soil layers. In addition, to test this assumption, we implemented a representation of the vertical distribution of the soil water isotopic composition: the soil water is spread vertically between several layers. The first layer contains a water height hydrology-current-research, where KD is the diffusivity of water molecules in water and Δt is the time step of the simulation, and the other layers contain a water height resol.L. The parameter resol can be tuned to find a compromise between vertical resolution and computational time. Layers are created from the top to bottom until all layers are full with water except the deepest one that contains the remaining soil water. For example, with L =0.67 mm, up to 16 layers can thus be created if the soil is saturated. Bare soil evaporation is extracted from the first layer. Transpiration is extracted from the different layers following a root extraction profile that reflects the sensitivity of transpiration to soil moisture [164]. Drainage takes water from the deepest layer. In the control simulation, rain and snow melt are added to the first layer (piston-like flow). In a sensitivity test, that can also be homogeneously distributed in the different layers, to crudely represent preferential pathways through fractures or pores in the soil.

At each time step, the soil water isotopic composition in each layer is re-calculated by taking into account the sources and sinks for each layer and ensuring that each layer remains full except the deepest one. Isotopic diffusion between adjacent layers is applied at each time step (Equation 6). The water budget of the total soil remains exactly the same as without vertical discretization.

Evaluation for an idealized case: The module representing vertical distribution of water isotopes in the soil is first evaluated for an idealized case when it is not yet embedded into ORCHIDEE.

First, we use a case in which the soil column evaporates at its top and is permanently refilled at the bottom by a water with δ18O of -8‰ [152]. The soil remains saturated, and we focus on the steady state reached after a few hundreds of days [152]. An analytical solution is available for this case [100,165]. The analytical solution and a much more sophisticated model of soil water isotopes (MuSICA [166]) yield very similar results (Figure 15a): the bottom of the soil is at -8‰ while the top of the soil is enriched up to 15‰. The soil module of ORCHIDEE is able to reproduce these results when the value of θl .τ is set to be very low (0.001) and when the vertical resolution is sufficiently high (layers of 0.75 mm). Whatever the value for θl .τ, ORCHIDEE results become less sensitive to the vertical discretization when layers are thinner than about 2 mm.


Figure 15: Vertical profile of soil water δ18O in idealized cases described by Braud [152]. a) The soil column evaporates at its top and is permanently refilled at the bottom by a water with δ18O =-8 ‰ b) The soil column is evaporated progressively until its soil water content is only 20%. See appendix 8.2 for more details. Simulations using the soil profile module of the isotopic version of ORCHIDEE (colors) with different parameters and vertical resolution are compared with the more sophisticated MuSICA and SiSPAT models and with an analytical solution. For θl.τ =0.005, the vertical resolution for ORCHIDEE is 0.15 mm for the first layer and 0.75 mm (resol=5), 1.5 mm (resol=10), 3 mm (resol=20) or 6 mm (resol=40) for the other layers. For θl.τ =0.01, the vertical resolution for ORCHIDEE is 0.21 mm for the first layer and 2.12 mm (resol=10) or 4.24 mm (resol=20) for the other layers. For θl.τ =0.1, the vertical resolution for ORCHIDEE is 0.67 mm for the first layer and from 1.34 mm (resol=2) to 3.35 mm (resol=5) for the other layers.

Second, we use a case in which the soil column, initially with a soil water of -8‰, evaporates at its top until the soil water content is only 20% [99]. The atmosphere has a relative humidity of 20% and a vapor δ18O of -15‰. The sophisticated models MuSICA and SiSPAT [152] feature a typical evaporative enrichment profile, with δ18O increasing from its initial value of -8‰ at the bottom to a maximum δ18O of 13‰ about 10 mm below the surface (Figure 15b). In the uppermost 10 mm, there is a slight depletion due to diffusion of water vapor into the soil column [101]. ORCHIDEE is not able to reproduce this vertical profile. First, since diffusion of water vapor in the soil is neglected, it is not able to simulate the depletion near the surface. Second, since θl .τ is temporally and vertically constant in ORCHIDEE, it is not able to adapt to the drying of the soil. In the sophisticated model, as the soil dries, the soil water content θl decrease, thus inhibiting vertical mixing of soil water and favoring strong isotopic gradients. In contrast in ORCHIDEE, θl .τ remains constant at a value representative of a moister soil, thus favoring vertical mixing of soil water and leading to a nearly uniform enrichment with depth [167-170].

To summarize, our representation of isotopic vertical profiles in ORCHIDEE is probably most suited when soil moisture remains high and does not vary too strongly.

Calculation of isotopic forcing from LMDZ outputs and nearby GNIP or USNIP stations

When precipitation and water vapor isotopic observations are not available at a given site, we create isotopic forcing using isotopic measurements in the precipitation performed on nearby GNIP (Global Network for Isotopes in Precipitation [122]) or USNIP (United States Network for Isotopes in Precipitation [125]) precipitation stations. To interpolate between the nearby stations, taking into account spatial gradients and altitude effects, we use outputs from an LMDZ simulation.

Let’s assume there are n GNIP or USNIP stations around the site of interest (MIBA or Carbo-Europe). The isotopic composition of precipitation at the site of interest and for a given month, δp,site, is calculated as:




and where Di is the geographical distance between the site of interest and the GNIP or USNIP station, δp,lmdz (s) is the precipitation isotopic composition simulated by LMDZ in the grid box containing the site s, δp,lmdz (i) is the precipitation isotopic composition simulated by LMDZ in the grid box containing the GNIP or USNIP station, δp,NIP (i) is the precipitation isotopic composition observed at the GNIP or USNIP station, Zsite is the altitude of the site of interest, Zlmdz (s) is the altitude of the LMDZ grid box containing the site of interest and as is the slope of the isotopic composition as a function of altitude simulated by LMDZ in the grid boxes containing and surrounding the site of interest [171]. The first term on the right hand side corresponds to the raw LMDZ output for the site of interest. The second term allows us to correct for the altitude effect. Since LMDZ is run at a 2.5 ° latitude * 3.75 ° longitude resolution, we cannot expect the average grid box size to be representative of the local altitude at the site. The third term allows us to correct for possible biases in LMDZ compared to GNIP and USNIP observations. Table 3 lists the GNIP and USNIP stations used to construct the forcing at each site of interest.

To calculate the isotopic composition of the water vapor, we assume that although LMDZ might have biases for simulating the absolute values of precipitation and water vapor composition, it simulates properly the precipitation-vapor difference [47,60]. Therefore, the isotopic composition of water vapor at the site of interest, δv,site, is calculated as:


where δv,lmdz (s) is the isotopic composition of water vapor simulated by LMDZ in the grid box containing the site of interest.

A simple equation to relate the soil water isotopic composition to the surface soil water budget

To explore how the isotopic composition of soil water can help estimate terms of the soil water budget, we derive here a very simple theoretical framework.

We assume that the water mass balance is:

hydrology-current-research (10)

where P is the precipitation, R the surface runoff, E is the bare soil evaporation, T the transpiration and D the drainage. Similarly, the isotopic mass balance is:

hydrology-current-research (11)

Where Rp, RE, RT, RD and RR are the isotopic ratios of incoming water at the soil surface, bare soil evaporation, transpiration, drainage and surface runoff respectively.

We assume that the bare soil evaporation isotope ratio depends on that of the soil (Rs) following the Craig [75] relationship (Equation 2) and that the transpiration composition is equal to that of the soil (RT=Rs) implying little vertical variations in soil water isotope ratios [172]. We assume that the isotopic composition of surface runoff is that of the incoming water (RR=Rp) and that the isotopic composition of drainage is that of the soil water (RD=Rs). In doing so, we neglect again vertical isotope variations in the soil and the temporal co-variation between Rs , D and T. Combining equations for the mass balance of water (Equation 11) and of water isotopes (Equation 10) then yields:

hydrology-current-research (12)

where I=P-R represents the incoming water that infiltrates into the soil. E/I represents the proportion of the infiltrated water which is evaporated at the soil surface.

The composition of the bare soil evaporation flux, RE, is a function of Rs following the Craig [75] formulation (Equation 2). Replacing RE by its function of Rs in Equation 12 allows us to deduce E/I :

hydrology-current-research (13)

Therefore, E/I is a function of the isotopic difference between the soil water and the precipitation water, which is easy to observe on instrumented sites such as MIBA or Carbo-Europe sites.


We thank a reviewer for his thorough review and detailed comments. LMDZ and ORCHIDEE simulations were performed on IDRIS machines to which access was granted by GENCI under project 0292. We thank Katia Laval for fruitful discussion and comments on an earlier version of this manuscript. We thank Matthias Cuntz for discussions. We thank Arthur Gessler and Romain Barnard for providing their data from Hartheim, and thank Chun-Ta Lai for providing his data from the Kansas prairie. We thank Danilo Dragoni, Kim Novick and Rich Phillips for providing information and data on the Morgan-Monroe site. We thank Marion Devaux, Cathy Lambrot (Inra-Ephyse, France), Rolf Siegwolf (Paul Scherrer Institute, Switzerland), Glyn Jones and Howard Griffiths (University of Cambridge, UK) for sampling and analysis of the isotopic data on the Bray and Mitra sites. We thank Eyal Rotenberg and Jean-Marc Bonnefond for providing the meteorological forcing over Yatir and the Bray respectively. We thank Dan Yakir for the isotopic and meteorological data collection in Yatir, his role in the MIBA initiative and comments on the manuscript. Part of the work was done while Camille Risi was a post-doc advised by David Noone, who I thank as well. This work benefited from financial support of the LEFE project MISSTERRE. Cathy Kurz-Besson was supported by the Fundação para a Ciência e Tecnologia (PTDC/AAG-REC/7046/2014). Lisa Wingate was supported by a Marie Curie Career Development Fellowship, thus some of the research leading to these results has received funding from the [European Community’s] Seventh Framework Programme ([FP7/2007-2013] under grant agreement n ° [237582]. The research was supported partly by the Czech Science Foundation project to JS (14-12262S) and by the Czech research infrastructure for systems biology C4SYS project (LM2015055).


Citation: Risi C, Ogée J, Bony S, Bariac T, Raz-Yaseef N, et al. (2016) The Water Isotopic Version of the Land-Surface Model ORCHIDEE: Implementation, Evaluation, Sensitivity to Hydrological Parameters. Hydrol Current Res 7: 258. Doi: 10.4172/2157-7587.1000258

Copyright: © 2016 Risi C, et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Select your language of interest to view the total content in your interested language

Post Your Comment Citation
Share This Article
Recommended Conferences

European Conference on Aquaculture and Fisheries

Stockholm, Sweden

11th International conference on Fisheries &Aquaculture

Vancouver, Canada

12th Global Summit on Aquaculture & Fisheries

Sydney, Australia
Article Usage
  • Total views: 10062
  • [From(publication date): 12-2016 - Feb 28, 2020]
  • Breakdown by view type
  • HTML page views: 9892
  • PDF downloads: 170
Share This Article