About Us | Help Videos | Contact Us | Subscriptions

Vadose Zone Journal - Article



This article in VZJ

  1. unlockOPEN ACCESS
    Received: Dec 19, 2016
    Accepted: Mar 08, 2017
    Published: May 18, 2017

    * Corresponding author(s): felicien.meunier@uclouvain.be
Request Permissions


Measuring and Modeling Hydraulic Lift of Lolium multiflorum Using Stable Water Isotopes

  1. Félicien Meunier *a,
  2. Youri Rothfussbc,
  3. Thierry Bariacc,
  4. Philippe Bironc,
  5. Patricia Richardc,
  6. Jean-Louis Durandd,
  7. Valentin Couvreura,
  8. Jan Vanderborghtb and
  9. Mathieu Javauxab
  1. a Earth and Life Institute, Univ. Catholique de Louvain, B-1348 Louvain-la-Neuve, Belgium
    b Institute of Bio- and Geosciences, IBG-3 Agrosphere, Forschungszentrum Jülich GmbH, D-52425 Jülich, Germany
    c UMR 7618 iEES-Paris (UPMC, CNRS, INRA, IRD, Paris-Diderot, UPEC), F-78850 Thiverval-Grignon, France
    d UR P3F (INRA), F-86600 Lusignan, France


This study tested a method to quantify and locate hydraulic lift (HL, defined as the passive upward water flow from wetter to dryer soil zones through the plant root system) by combining an experiment using the stable water isotope 1H218O as a tracer with a soil–plant water flow model. Our methodology consisted in (i) establishing the initial conditions for HL in a large rhizobox planted with Italian ryegrass (Lolium multiflorum Lam.), (ii) labeling water in the deepest soil layer with an 18O-enriched solution, (iii) monitoring the water O isotopic composition in soil layers to find out changes in the upper layers that would reflect redistribution of 18O-enriched water from the bottom layers by the roots, and (iv) comparing the observed soil water O isotopic composition to simulation results of a three-dimensional model of water flow and isotope transport in the soil–root system. Our main findings were that (i) the depth and strength of the observed changes in soil water O isotopic composition could be well reproduced with a modeling approach (RMSE = 0.2‰, i.e., equivalent to the precision of the isotopic measurements), (ii) the corresponding water volume involved in HL was estimated to account for 19% of the plant transpiration of the following day, i.e., 0.45 mm of water, and was in agreement with the observed soil water content changes, and (iii) the magnitude of the simulated HL was sensitive to both plant and soil hydraulic properties.


    DaS, days after seeding; HL, hydraulic lift; RLD, root length density; RWU, root water uptake

Hydraulic lift (HL) is defined as the upward water flow from wetter to drier soil layers through passive root transport (Richards and Caldwell, 1987) and has been reported for numerous plant species (Neumann and Cardon, 2012). In principle, HL takes place when the root water potential is higher in upper soil layers than the soil water potential containing these roots, while in deeper layers the soil water potential is higher than the root water potential contained in this layer. Hydraulic lift is a particular case of hydraulic redistribution, which may happen in any soil direction, e.g., from shallower to deeper layers or in the horizontal direction (Burgess et al., 1998).

The quantitative and qualitative implications of HL on the plant water balance, nutrient mobilization, and competition for space and nutrients between plants are still unclear among the scientific community (Burgess, 2011). From a hydrological perspective, the absolute volume of hydraulically lifted water has been debated; it has been reported to vary by nearly two orders of magnitude depending on the species and environmental conditions (representing between 2 and 143% of the plant transpiration, with a mean of 30%) (Neumann and Cardon, 2012). Discriminating between soil and root hydraulic redistribution is also very challenging (Katul and Siqueira, 2010). This leads to divergence among researchers on the importance of HL: in some studies, it was estimated to be a negligible term of the water budget (e.g., Meinzer et al., 2004, Scholz et al., 2010), while according to other studies it might significantly impact evapotranspiration and net ecosystem C exchange at the stand and catchment scales (Domec et al., 2010; Jackson et al., 2000). The ecological interest in HL might be not only to improve water uptake from soils but also to mobilize nutrients and enhance rhizosphere microbial communities (e.g., see reviews by Neumann and Cardon, 2012; Prieto et al., 2012).

The stable water isotope 1H218O is a tracer of ecosystem processes (Rothfuss and Javaux, 2016; Javaux et al., 2016; Sprenger et al., 2016; Werner et al., 2012). Because plant root water uptake (RWU) does not lead to fractionation of stable water O isotopes for a wide variety of plants (e.g., Barnard et al., 2006), when the soil water natural abundance ratio (expressed in ‰ on the international accepted scale and referred to in the following as δ18O; Gonfiantini, 1978) varies within the root zone, the stem xylem δ18O is a valuable tool for determining the origin of water extracted by the root system (e.g., Asbjornsen et al., 2007; Dawson and Ehleringer, 1991; McCole and Stern, 2007). When RWU and release do not lead to fractionation, stem xylem water δ18O of plants with different rooting depths can also be used to identify the injection of the isotopic signal from deeper soil layers into shallower ones. Caldwell and Richards (1989), for example, inferred from measurements of stem xylem water δ18O that isotopically labeled water injected in deep roots of sagebrush (Artemisia tridentata Nutt.) and released through HL was extracted by shallow roots of Agropyron desertorum (Fisch. ex Link) Schult. Hydraulic lift is often inferred from diurnal variations in the soil surface water potential (Dawson, 1993; Ludwig et al., 2003) and water content (Brooks et al., 2002; Warren et al., 2005). Only a very few direct observations of changes in soil water isotopic composition have been made (e.g., Zegada-Lizarazu and Iijima, 2004). Because no accurate isotopic “anomaly” in the soil water isotopic profile was reported, the depths where HL occurs and the water volumes that are involved could not be determined accurately.

Hydraulic lift can be explicitly simulated with existing models of the soil–plant–atmosphere continuum that describe the water flows in one (Couvreur et al., 2014a; Volpe et al., 2013; Guswa, 2012) or three dimensions (e.g., Javaux et al., 2008) as a passive process based on soil water potential gradients. Recently, Couvreur et al. (2012) developed a three-dimensional solution of RWU based on water flow equations in roots. They demonstrated that the RWU was a sum of two terms—one proportional to the root distribution and the other a function of the soil water potential distribution, as well as the root distribution. At locations where this second term exceeds the first one, RWU is negative and HL occurs. This implies that the existence and intensity of such a phenomenon depend on both environmental conditions and plant properties: low transpiration, high soil water potential variability, and high root conductance are the three conditions for HL (see Appendix A for more details). Interestingly, several researchers have developed similar root water uptake models (Javaux et al., 2013), confirming that this is an appropriate way to model HL (Volpe et al., 2013). Such models and their predictions can help in building an experiment in which the environmental prerequisites are met and quantifying the water and isotope flow in the soil–root system continuum.

The general objective of this study was to develop a methodology to quantify HL in the soil–plant system and investigate and understand the factors controlling the magnitude and location of the water efflux. This objective was addressed by combining experiment and modeling tools. First, we set up a unique controlled experiment in which the stable water isotope 1H218O was supplied to the deepest soil layers of a rhizobox and was monitored in soil and plants. In addition, an existing soil–root model (R-SWMS, Javaux et al., 2008), which includes all relevant processes to simulate RWU and HL according to basic physical principles, was extended to the transport of the 1H218O isotope. This model was used to test whether the observed distributions of 1H218O could be reproduced using meaningful parameters and to assess the water volumes that were involved in HL.

Material and Methods

Experimental Setup

A macro-rhizobox with transparent polycarbonate sides (dimensions length by depth by height = 100 by 20 by 160 cm, Fig. 1a) was filled with soil (Dystric Cambisol, particle size distribution 15% sand, 65% silt, and 20% clay, sieved to 2 mm) in 10-cm layers to reach a bulk density similar to the field, i.e., 1.42 g cm−3 (Rothfuss et al., 2010). It was placed in a greenhouse (Lusignan, France, 46.44°N 0.13°E) and continuously weighed (KE1500, Mettler-Toledo, 1500-kg range, 20-g resolution) to monitor the soil water balance and calculate water effluxes (i.e., evapotranspiration). Underneath the rhizobox and in contact with the soil bottom surface, a water reservoir (10-cm height) filled with gravel acted as a water table and allowed the discontinuous supply of water from the bottom (see Fig. 1 for illustration). The soil water content θ (cm3 cm−3) and pressure heads hsoil (hPa) (a list of the symbols used in the present study as well as the dimensions, units of the variables, and parameters and reference or optimized values are given in Appendix B) were monitored within the soil at six depths (z = −5, −10, −25, −60, −90, and −130 cm) with water content reflectometers (CS 616, Campbell Scientific, ?2.5% precision), micro-psychrometers (Soils– PST55, WESCOR), and tensiometers (SMS 2000 with pressure transducers SKT850, SDEC).

Fig. 1.
Fig. 1.

(a) Soil macro-rhizobox experimental setup with Lolium multiflorum cover and different profiles of soil probes and (b) sketch of the macro-rhizotron with height H, length L, and width l with a representation of a single plant root system.


In the greenhouse, atmospheric temperature and relative humidity were measured but not regulated. The amount and isotopic composition of water applied to the soil from the reservoir were controlled during the experiments. The rhizobox soil was initially saturated with input water (δIW = −8.1‰) from the bottom reservoir. Italian ryegrass was seeded (seeding density = 0.03 kg m−2) when the soil surface reached field capacity and emerged 12 d later.

The water isotopic tracer experiment was performed 116 d after seeding (DaS 116). At that time, there was a strong soil water potential gradient between deep (100–140 cm) and surficial (0–30 cm) soil layers as measured by the tensiometers and micro-psychrometers, and roots were observed to reach the water reservoir (see below). At 2:00 PM, the water reservoir was labeled with an 18O-enriched water solution (δLW = 400‰) and put in contact with the soil outlet. To assess the impact of this tracer, soil water isotopic composition profiles were determined at two times: on DaS 116 at 3:15 PM and on DaS 117 at 5:20 AM. To do so, horizontal soil samples were extracted with a 1.2-cm-diameter auger (corresponding to a soil sample volume range of ∼0.6–1.7 cm3) from the surface down to the 130 cm-depth at 5- to 15-cm increments and then subjected to isotopic composition analyses (see below).

At 14 times between 2:00 AM and 5:20 AM on DaS 117, three tillers were sampled homogeneously from the vegetation’s cover and placed together in a single vial. All dead material and the sheath of the oldest living leaf around each tiller were discarded to remove any organ that could have been transpiring (Durand et al., 2007).

After the second soil sampling for isotopic analyses, additional soil samples were taken between the polycarbonate sides using the same 1.2-cm-diameter auger to determine the root length density (RLD) profile. Each depth was sampled two or three times and combined in a single vial for each sampling depth. From each sample, soil particles were washed off and roots were manually placed over a 0.2-mm-mesh filter and dried at 60°C for 48 h. The RLD was determined from the obtained root dry masses using the specific root length for Italian ryegrass (181 m g−1) as determined by Gonzalez-Dugo et al. (2005).

Isotopic Analysis

Water was removed from the soil and tiller samples by cryogenic vacuum distillation (−196°C) and sealed in closed glass containers. Recovery was >98.5% and a correction, ranging from 0 to 0.4‰, was applied according to the distillation yield by a Rayleigh-type law (Araguás-Araguás et al., 1995). The table water O isotope composition was measured by isotope ratio mass spectrometry (Aquaprep Micromass Isoprime) after CO2–H2O equilibration (Gilson X222 system) using the Epstein and Mayeda (1953) technique.


A fully mechanistic, three-dimensional model of the soil–plant system, R-SWMS, was used to analyze the isotopic data set and quantify the HL. This model is based on a physical representation of the water flow in the soil and in the root system (Javaux et al., 2008). To ease the optimization procedure, we first made use of the root–rhizosphere module of this model, thereby neglecting bulk soil water flow. Based on the assumption that bulk water flow was minimal during a short period, we used this fast module with a global optimization algorithm to adjust the root and rhizosphere hydraulic parameters. The full R-SWMS, accounting for water flow and solute transport in the bulk soil, was then used to verify these assumptions and to fine tune our model parameters. Note that to generate a realistic three-dimensional root system architecture for R-SWMS and the root–rhizosphere module, the root growth model Root Typ was used.

Root–Rhizosphere Model

The root–rhizosphere model solves water flow from the bulk soil to the root collar. In this module, bulk water flow and 1H218O transport are neglected in the soil, and the bulk water potential and water isotopic compositions are model boundary conditions.

The water flow equation in the root system is solved with the finite difference approach of Doussan et al. (1998). The root system is divided into N segments leading to N linear equations that are solved for the root water potential given known boundary conditions at the root collar and soil–root interfaces. Once solved, these equations allow calculation of the radial water flows in the ith root segment, Qi (cm3 d−1):where kr,i (cm hPa−1 d−1) is the radial conductivity of the root (i.e., between the root surface and xylem of the ith root segment), Ai = 2πrili (cm2) is the surface of the root segment of radius ri (cm) and length li (cm), and Hsr,i (hPa) and Hx,i (hPa) are the soil–root interface and the root xylem water potential of the segment, respectively. The radial flow is defined as positive toward the root. Here and in the following we distinguish between the total water potential H and the matric potential h (hPa). The latter is simply the former without the gravitational component. The subscripts associated with the water potentials locate them. It is worth noting that in this study we used the nomenclature and the units of soil and root conductivity or conductance of Javaux et al. (2013). To account for the soil hydraulic resistance in the rhizosphere, the Doussan et al. (1998) hydraulic architecture was extended by an extra soil shell around the roots with a hydraulic conductivity ksoil (cm2 hPa−1 d−1) through which radial flow from the bulk soil to the soil–root interface takes place. The equivalent conductance between the bulk soil around and the xylem of the ith root element, Ki,eq (cm3 hPa−1 d−1), relates to the flow rate to the segment Qi and the water potential difference between the xylem and the bulk soil, Hsoil,i, aswhere Ki,eq = kr,iksoilliriB/(Bksoil + rikr,i (cm3 hPa−1 d−1) is the equivalent conductance between the soil node and the xylem of the ith root element and B (dimensionless) is a geometric parameter. Appendix C describes the calculation for Ki,eq and analyzes the limiting cases of the equivalent soil–root conductances. In this approach, soil water flow is not explicitly modeled but soil hydraulic conductivity is accounted for through ksoil and its dependence on the soil water potential. The variable ksoil is indeed a function of the soil water matric potential h (hPa) through the hydraulic conductivity curve k(h):where hsoil and hsr,i are the matric potential of the bulk soil and of the soil–root interface, respectively, and φsr,i and φsoil are the matric flux potentials of the soil–root interface and the soil, respectively (de Jong van Lier et al., 2008). The matric flux potential is defined as the integral of the soil conductivity curve. Figure 2a illustrates one of the numerous segments of the whole root system and its surrounding soil. The total xylem potential of this segment Hx,i (black node), the soil–root total potential Hsr,i (brown node), and the bulk soil potential total Hsoil (blue node) are represented as well as the soil and root radial conductances Ksoil and Kr,i (see Appendix C). It is worth mentioning that several studies reported different root conductivities for in- and outflows, but this aspect was not implemented in our modeling exercise. Other researchers have shown that shallow and deep roots may have different root hydraulic conductivities (Wan et al., 1994), but this is explicitly accounted for in our model through root ages and types (see below).

Fig. 2.
Fig. 2.

(a) Layout of a root segment of radius ri and length Ii in the root–rhizosphere model, where water flows through the rhizosphere conductance (Ksoil of characteristic length lr) and the root radial conductance (Kr,i) to reach the root xylem vessels according to the gradient of potential between the bulk soil (Hsoil), the soil–root interface (Hsr,i), and the xylem potential (Hx,i), and (b) the same root segment situated in a soil element of volume Vj in R-SWMS, where Hx,i is still the root xylem potential and Sj is the sink term of the soil element.


R-SWMS: Transient Soil–Root Water Flow Model

In a second step, we used the full transient soil–root water flow model R-SWMS to account for water flow in the soil in addition to water flow in the root and the rhizosphere. The R-SWMS model (Javaux et al., 2008) solves the coupled problem of water flow in the soil–root continuum by iteratively solving soil and root water flux equations to convergence at each time step. A convection–dispersion equation is solved for describing solute movement in the soil, and the soil water movement is modeled by the Richardson (1922) equation with a sink term to account for RWU. The sink term Sj (d−1) for each jth soil element (having a volume Vj [cm3]) in this equation is calculated by summing up the water radial fluxes for each ith root segment contained in this particular soil element:

Calculation of Qi is obtained via the root–rhizosphere module discussed above. Transport of isotopes in the soil is modeled via a classical convection–diffusion equation with passive uptake only (Schroeder et al., 2012). The model thus explicitly accounts for soil water flow, isotopic dispersion in the soil, and root water uptake and release. Figure 2b shows the ith root segment, one of those composing the whole hydraulic architecture, in the surrounding jth soil element, the element sink term Sj, and its volume Vj.

Isotope Uptake Modeling

The prognostic variable for the computation of isotopic flow is the O isotopic liquid concentration C (kg m−3), which can be obtained from the O isotopic liquid composition δ (‰) (Braud et al., 2005) at any time t (d):where ρ (kg m−3) is the volumetric mass density of water, Rref (dimensionless) is the Vienna Standard Mean Ocean Water (VSMOW) O isotopic ratio [Rref = (2005.2 ± 0.45)10−6], and Mw and Mi (kg mol−1) are the molar masses of 1H216O and 1H218O water, respectively. Considering perfect mixing inside the root xylem and no root capacitance (no water stored in root compartments), we define an equivalent O stable isotopic concentration inside the root system, CTi(t) (kg m−3) as a sink-weighted average O isotopic concentration of absorbed water:

where the local sink term is defined positive if water is taken up. When HL occurs, the outgoing flow is supposed to have the equivalent stable O isotopic concentration assuming an instantaneous mixing of all the water sources within the plant. In this study, we modeled the stable water isotope as a “standard” non-excluded solute—we assumed that root uptake does not lead to fractionation and we did not consider evaporation.

Model–Data Comparison

Root System Architecture

The root system architecture needed by the models was inferred by inverse modeling of the measured RLD with a root growth model, Root Typ (Pagès et al., 2004). Lolium multiflorum root growth parameters (maximum length and initial velocity) of both primary and lateral roots as well as the intermodal distance between two successive lateral roots were optimized to fit the measured RLD data. Root Typ was thus used to generate the realistic root system architecture that is needed as input for the root–rhizosphere module and R-SWMS.

Boundary Conditions for Soil and Root

Weight loss was considered due to transpiration solely (this assumption is justified thanks to observations; see below). The actual transpiration rate per plant, Tact (cm3 d−1), was obtained by dividing the weight loss by the number of plants, which was derived from the plant density observed 12 d after sowing. The transpiration rate was used as a boundary condition for the root system (a flux boundary condition) at the root collar in both modeling steps. The soil water content profile measured at the beginning of the HL experiment (DaS 116, 3:15 PM) was used as the initial condition. The bulk potential used as boundary conditions for the root–rhizosphere model came from the measured tensiometer or psychrometer profiles. No flow and a pressure condition were used as boundary conditions for R-SWMS at the top and bottom of the soil, respectively. The pressure values used were obtained from the tensiometers at the very bottom of the soil profile.

The time step in the root–rhizosphere model was fixed to 10−2 d. During a time step, the boundary conditions, i.e., the transpiration rate at the collar and the bulk water potentials in the soil, were assumed to be constant. The time step is adaptive in R-SWMS but was always between 10−5 and 10−2 d.

Soil Parameterization

The rhizobox was simulated as a bilayer medium formed by the soil and the reservoir underneath. Experimental retention curve data were derived from simultaneous measurements of the soil water content and potential provided by collocated water reflectometers, micro-psychrometers, and tensiometers (see Fig. 3; Table 1). The van Genuchten parameters were obtained by fitting the water retention function to the observed soil water content and soil water matric potential pairs (see Fig. 3). A Burdine-type condition was used to link the shape parameters m and n (m = 1 − 2/n; Burdine, 1963) of the water retention curve.

Fig. 3.
Fig. 3.

Soil water content vs. soil water potential as measured by the probes (blue symbols) or as fitted by the van Genuchten parameters given in the first line of Table 1.


View Full Table | Close Full ViewTable 1.

Soil (fitted to simultaneous measurements of water content and water potential) and reservoir (fixed) van Genuchten (1980) and Brooks and Corey (1964) retention and conductivity curve parameters, including saturated and residual volumetric water content (θs and θr, respectively), shape parameters α and n, saturated hydraulic conductivity (kS), and the conductivity function λ. The soil saturated conductivity and the parameter λ are obtained after optimization of the isotopic composition changes.

Source θs θr α n kS λ
cm3 cm−3 cm−1 cm d−1
Soil 0.4 0.044 0.0285 2.29 0.26 6.67
Reservoir 0.45 0.4 0.000344 1.42 1500 0.5

The retention and conductivity function parameters of the reservoir were imposed so that it stayed always saturated and much more (i.e., three orders of magnitude) conductive than the soil during the experiment. Table 1 summarizes the hydraulic parameters used in the simulations for both the soil and the reservoir.

Because the RWU model simulates one single plant while many distinct plants were present in the rhizobox, a continuous soil domain was used. It mimics the presence of plant neighbors by considering the vertical boundaries to be periodic for soil and root water fluxes and the root system architecture. In particular, outgoing roots reappear at the opposite face of the domain in the x and y directions (Couvreur et al., 2012).

Optimization Scheme

Missing information on root and soil hydraulic conductivities were obtained by inverse modeling from the experimental data set using, in a first step, a global single-objective optimization of the root–rhizosphere model that neglects soil water flow (see above). The parameters provided by the global optimization algorithm were then further locally optimized, taking into account the soil water movement using R-SWMS. The optimizations were achieved in MATLAB with the Multistart heuristic algorithm OQNLP (Optimal Methods Inc.) for the global optimization and with a local algorithm (Nelder–Mead simplex direct search) for the local one. Both root and soil hydraulic parameters were optimized. The optimized soil properties were the soil saturated conductivity (kS [cm2 d−1 hPa−1]) and the shape parameter of the Brooks and Corey (1964) conductivity function (λ [dimensionless]), whereas the other soil hydraulic parameters were derived directly from the retention curve (see Appendix D). Root hydraulic conductivity functions consisted of a monotonic (increasing for axial conductivities and decreasing for radial conductivities) piecewise linear curve function of the root age. To reduce the number of parameters to be optimized, the following simplifying assumptions were made: (i) root conductivity depended on the age and order of the root segment only (see, e.g., Bramley et al., 2009; Sanderson, 1983; Zwieniecki et al., 2002); (ii) the transition ages between successive conductivity plateaus were considered as fixed; and (iii) the primary roots were not responsible for root water uptake. The root and soil hydraulic property functions are shown in Fig. 4a. The following terminology for the root hydraulic properties was used: x and r subscripts stand for axial and radial conductivities, respectively, Roman numerals I and II describe the root order, i.e., primary and secondary (i.e., lateral) roots and, finally, a and b refer to the age of the root section (“young” and “old,” respectively).

Fig. 4.
Fig. 4.

(a) Soil and root hydraulic conductivity functions: the soil conductivity (ksoil) is a function of the soil water potential (Hsoil) and is parametrized with the conductivity function λ and the saturated conductivity kS. The root radial (kr) and axial (kx) conductivities depend on the root age (a = young, b = old) and order (I = primary, II = laterals); and (b) model optimization scheme: the models integrate information from different sources (such as the measured transpiration rate (Tact) and the initial isotopic profile (δmeas) at different depths z) between the initial (t0) and the final (tf) times to calculate objective functions (OF1 and OF2) according to the simulated isotopic profiles (δsim).


Two objective functions were used. The first objective function, OF1 (‰), was defined as the RMSE of measured (δmeas [‰]) vs. simulated (δsim [‰]) isotopic composition profiles:

where P is the soil and root parameter set, and N and Np are the numbers of observations and parameters, respectively.

The second objective function, OF2 (‰), is the maximum deviation value of the absolute difference between measured and simulated isotopic composition:

These two indicators allowed (i) estimating whether the anomaly in the isotopic profile was correctly represented, and (ii) discriminating a null model (i.e., no change in the isotopic profile).

Soil parameters obtained from infiltration measurements (Rothfuss et al., 2012) and root segment conductivity functions of a wheat (Triticum aestivum L.) plant (Couvreur et al., 2014a) were used as initial guesses for the soil and root hydraulic parameters, respectively. Figure 4b presents the two-step optimization scheme that was used. First, a global optimization with the root–rhizosphere model was performed using soil parameters obtained from infiltration measurements (Rothfuss et al., 2012) and the root segment conductivity functions of wheat (Couvreur et al., 2014b) as initial parameter estimates. Subsequently, the optimized parameters from the first optimization were used as the initial parameters for the local optimization and the transient model R-SWMS. The root–rhizosphere model globally optimized the soil and root hydraulic parameters without considering any bulk soil water movement, which makes the model fast. R-SWMS was then used with the optimized parameter set to take into account water flow and solute transport in the bulk soil, requiring additional information about soil and solute. The first parameter set guess for R-SWMS was the optimal parameter set provided by the root–rhizosphere model. Both models used a root system architecture generated by Root Typ and both optimizations used the objective functions OF1 and OF2 for assessing the simulation fitness. To evaluate the uniqueness of the solution and the sensitivity of the simulations to soil and/or root hydraulic parameters, a local sensitivity analysis was also performed by varying the parameters by ±20% around the found global optimum.


Experimental Results

The experimental results are summarized in Fig. 5. The first subplot represents the water content profile as measured by the reflectometers at the beginning and end of the experiment. The upper 80 cm are highlighted in the frame. The development of the vegetation and the observed surface volumetric water content (θ < 0.125 cm3 cm−3 corresponding to low soil potential, see Fig. 3) led to the assumption that the rhizobox water losses were due to transpiration flux solely (i.e., evapotranspiration = transpiration). On DaS 116, a strong water potential gradient existed between the soil surface (hsoil ? 105) and the bottom (hsoil ? 102) of the rhizobox soil. In the topsoil, the changes in water content were absent (15- and 20-cm depths), slightly negative (30- and 35-cm depths, maximum −0.2%), or slightly positive (10- and 25-cm depths, maximum 0.8%).

Fig. 5.
Fig. 5.

Experimental results: initial (116 d after seeding [DaS], 3:15 PM, red circles) and final (117 DaS, 5:20 AM, blue triangles) (a) soil volumetric water contents θ, (b) isotopic composition profiles (the isotopic anomaly is highlighted while no significant water content change could be observed), and (c) actual transpiration Tact as measured by weight loss during the experiment.


Figure 5b represents the measured isotopic profiles. The observed surface isotopic enrichment was mainly obtained before seeding, when the soil was bare and soil water content was still high and thus the soil was actively evaporating. As the vegetation layer progressively covered the soil and surface water content decreased, water losses were then mainly due to non-fractionating root water uptake. As a consequence, the differences observed between the two isotopic profiles following labeling of the reservoir water could only be due to the spatial heterogeneity of soil water potentials entailing water redistribution through capillary flow and HL. Measured soil water δ18O spanned between −8.3 and 182.5‰ on DaS 116 and between −8.2 and 188.6‰ on DaS 117 (the two lowest points in the profiles represent the δ18O values of the labeled water in the reservoir). The δmeas profiles were comparable with the exception of the 20- and 25-cm depths, where enrichments of 2.0 and 1.3‰ were observed on DaS 117 with respect to the values measured at the same depths on DaS 116 (see detail plot in Fig. 5b). These differences were ≥10 times larger than the isotopic measurements’ precision (i.e., 0.2‰). At around the 1.35-m depth, significant differences between both measured δmeas profiles were also observed due to the bottom supply of isotopically labeled water under pressure (the gravimetric water content measurements show significant differences of around 4% at these depths). Finally, Fig. 5c presents the actual transpiration rate of the rhizotron as measured by the weight losses by the scale from the beginning to the end of the experiment. The transpiration (even at night) was relatively high except for two short episodes of almost null transpiration. Consequently, during this experiment, we covered a wide range of plant actual transpiration rates.

The highest RLD was measured at the soil surface (38.40 cm cm−3) and decreased exponentially to the bottom of the soil profiles (0.02 cm cm−3), with a mean value of 3.5 cm cm−3 as shown in Fig. 6a. The RLD profile could be well fitted (r2 = 0.92, RMSE = 0.68 cm cm−3) with the ryegrass root system architecture generated by Root Typ. The optimal root system architecture is represented in Fig. 6b in its periodic domain. It is the same root system as in Fig. 1b in a shrunk (continuous, see above) domain. This root system architecture was used in both models (the root–rhizosphere model and R-SWMS).

Fig. 6.
Fig. 6.

(a) Measured (blue dots) and simulated (black dots) root length density (RLD) and (b) optimized root system architecture in the periodic soil domain.


Modeling Results

In Fig. 7a, simulation results obtained by running the R-SWMS model with the locally optimized plant and soil parameters are compared with data from DaS 116, 3:15 PM, and DaS 117, 5:20 AM. The soil water δ18O anomalies observed at 20 cm (−5.3‰) and 25 cm (−6.2‰) were well reproduced by the model (i.e., −5.5 and −6.3‰, respectively). The RMSE is comparable to the precision of the isotopic measurements (OF1 = 0.21‰) and corresponds to a correlation coefficient of 0.98. These good results are the consequence of the good calibration procedure. For both of the impacted layers, a change in θ from 0 to 0.02 cm3 cm−3 was modeled, in agreement with the observations made in situ by the water reflectometers, taking into account their measurement’s accuracy and representative measuring soil volume. The simulated volume of water involved in the HL process corresponded to 19% of the transpiration of the next day, i.e., to 0.45 mm of water. Figure 7b details the water radial flow entering (positive) or leaving (negative) the root as a function of the vertical position of the segment during a low transpiration episode. The root system takes up water in the middle and deep soil layers and releases it in shallow layers but not at the very top. The difference between the uptake and the release corresponds to the plant actual transpiration.

Fig. 7.
Fig. 7.

(a) Best-fit soil water isotopic composition profile (black dashed line) compared with initial (116 d after seeding [DaS], red symbols) and final profiles (117 DaS, blue symbols) and (b) root radial flow (Qi) profile during the first low transpiration rate episode.


Figure 8 summarizes the sensitivity analysis results. Each single subplot represents the OF1 obtained when modifying two plant and/or soil parameters concomitantly around their optimal values. The white crosses indicate where the optimal parameters are in each subplot In this contour plot, blue corresponds to low values of the OF (or a good model-to-data fit, i.e., OF close to its minimal value of 0.2‰), while red zones are made of high OF values (bad model-to-data fit, i.e., misrepresented peak of HL). When no HL was simulated, the OF reached a value of 0.5‰. When either the location or the intensity of the HL was misrepresented, the OF took values above 0.2‰ and could exceed 0.5‰. It was observed that the simulation was insensitive to the young xylem conductivity of the primaries (kxIa); in other words, the OF was not impacted by a change in the kxIa value. This was due to the age of the total root system: because the experiment was conducted 4 mo after plants were sown, the percentage of young primary roots among the root system was very low. Moreover, we assumed (see above) that these primary root segments had a zero radial conductance and did not take part in the uptake or release of water. However, because the kxIIa of these young laterals were low and systematically lower than the primary axial conductivity in the simulations (i.e., kxIIa < kxIa), the young primary axial conductivities never limited water flow. In the bottom-right plot, a strong correlation between kS and λ was observed, illustrating the antithetical effects of the saturated conductivity and the shape parameters in the Brooks and Corey (1964) conductivity curve formulation. A larger kS makes the unsaturated conductivity increase across the entire range of matric potentials, whereas a larger λ makes the conductivity decrease with decreasing potentials or water contents. Note that the two last columns show opposite trends due to the contradictory effect of these two soil parameters. Finally, the root axial conductivity of the old segments (primary and laterals) as well as the root radial conductivity of the old laterals controlled the water volumes flowing from the bottom to the top soil layers through the root system and thus the δ18O anomaly. In this 4-mo-old root system, water mostly entered in the deep soil layers and exited the root system in the upper layers through old lateral segments and was transported by the vessels of primaries and laterals. If one of the respective conductances was too low, the isotopic anomaly could not be reproduced in terms of absolute values because not enough 1H218O was brought to the top layers. To the contrary, if these conductances were too high, the root system would have released too much water and 1H218O to the topsoil.

Fig. 8.
Fig. 8.

Root–rhizosphere model sensitivity analysis of the soil and root hydraulic conductivity functions. The optimal parameter set is shown by the white crosses. The root radial (kr) and axial (kx) conductivities depend on the root age (a = young, b = old) and order (I = primary, II = laterals). The soil conductivity is a function of the shape parameter λ and the saturated conductivity ks. See Fig. 3a for details. The local sensitivity analysis was performed by varying the parameters by ±20% around the found global optimum.


Figure 9 shows a comparison between measurements and best-fit simulation (solid line, root–rhizosphere model) of the xylem sap water isotopic composition in the tiller, δTi (‰), for the late simulation period. The general downward trend could satisfactorily be reproduced by the model (RMSE = 15.4‰) without any optimization of model parameters to fit the tiller water isotopic composition data, which was remarkable. However, the model was not able to capture the variations exhibited by the observed data. Incidentally, it is likely that no other model parameterization would be able to reproduce such variations because there were no model input variables that exhibited such high dynamics. These variations could well be due to the measurement protocol, not coming from isotopic analytical precision but rather from the spatial variability of root length density and uptake between plants. In particular with the doped skewed soil isotopic composition profile, it could be expected that little changes in the root hydraulic architecture could generate big changes in terms of tiller concentration. The light envelope contains 95% while the dark envelope contains 50% of the 1000 simulation runs with a modified parameter data set. The parameters were always kept in such a range that it led to satisfactory OF1 values for the simulated vs. modeled soil water isotopic profiles (<0.3‰).

Fig. 9.
Fig. 9.

Best-fit equivalent isotopic composition and uncertainties (the light and dark blue envelopes contain 95 and 50%, respectively, of the simulations run with a set of parameters whose values lead to a good fit of the delta isotopic profile) compared with measured tiller water isotopic composition δTi (red circles) on 117 d after seeding (DaS).


Figure 10a shows the relation between the actual transpiration rate and the actual release rate by the plant roots as simulated by R-SWMS; the absolute value of the instantaneous water release was negatively correlated with the transpiration rate with a high correlation coefficient: r2 = 0.92.

Fig. 10.
Fig. 10.

(a) Instantaneous water release rate as a function of the actual transpiration rate (Tact) during the experiment and (b) change in the standard uptake fraction when incorporating soil resistance as a function of soil depth.


A comparison of the local uptake normalized by the total uptake with and without soil resistance shows the main differences that occurred in the topsoil (Fig. 10b), with the soil resistance considerably reducing the release. The normalized uptake of all other segments naturally increased because of the normalization.


The three-dimensional soil–root model used in this study was perfectly able to reproduce the peak observed in the isotopic composition profile. The optimized parameters needed to obtain this result are in the range of observed values in the literature (see Appendix B), and HL as modeled in this study is consequently a likely mechanism explaining the experimental results. As mentioned above, we did not consider any change in root hydraulic conductivity for in- and outflow, but the root properties changed according to the root type and age.

The estimated volume of water moving upward during the HL episode represented 19% of the transpiration of the next day (according to the RWU model). This is in agreement with the results of other studies from the literature (Neumann and Cardon, 2012). For instance, Couvreur et al. (2014b) found that HL contributed to between 11 and 20% of the daily transpired water for a maize (Zea mays L.) plant in a simulation study. Caldwell and Richards (1989) and Brooksbank et al. (2011) found similar proportions for sagebrush and Eucalyptus kochii Maiden & Blakely (i.e., 34 and 27%, respectively). It is worth noting that the water efflux could not be observed through the water sensors because the HL fluxes were lower than the probe precision level while far from being negligible.

High values of transpiration at night (e.g., DaS 116, 11 PM–DaS 117, 1 AM) illustrated the particular atmospheric conditions prevailing in the glasshouse during the period of soil and plant sampling. Night transpiration has been shown to be an important mechanism (for a complete discussion, see Caird et al., 2006). Relative humidity reached unity only on a few occasions (e.g., DaS 117, 1–4 AM, data not shown), as the air temperature was still high (e.g., >15°C until DaS 117, 1 AM; data not shown), which did not prevent night transpiration.

Stable water isotopes are powerful tools for partitioning raw fluxes (e.g., transpiration) into net fluxes (e.g., RWU profile). However, the methodology used to measure the isotopic composition of plant and soil waters suffers from major drawbacks, which ultimately limits its potential. The present study illustrates some of these limitations: while θ and hsoil profiles could be measured in a nondestructive manner at high temporal resolution, the soil water δ18O profile could only be determined destructively on two occasions during the experiment. Even though this was partly due to limited soil availability inside the rhizobox (many soil profiles could have been sampled under field conditions where soil is not limiting), soil water δ18O profiles could not provide insights into the HL temporal dynamics. New methods (unavailable at the time when the experiment was designed and performed) for monitoring soil water isotopic compositions have been developed in recent years that propose to overcome these limitations (Gaj et al., 2016; Rothfuss et al., 2013; Volkmann and Weiler, 2014). They open avenues toward a dynamic disentangling of ecosystem fluxes and identification of water sources and sinks in the soil–plant–atmosphere continuum system (Gangi et al., 2015; Rothfuss et al., 2015; Volkmann et al., 2016).

Another limitation of this study is the relatively poor knowledge of the three-dimensional hydraulic architecture of the root system used in the R-SWMS model. We fitted a root architectural model (Root Typ, Pagès et al., 2004) to the measured RLD profile. However, different Root Typ parameters could have fitted the RLD profiles equally well and led to distinct hydraulic architectures. The optimized soil and plant hydraulic parameters accurately reproduced the experiment, but the solution is not unique. A root system with more laterals in the topsoil (offset by a shorter primary length) would fit the RLD profile as well, and the optimization would lead to a reduced radial hydraulic conductivity for these segments.

The soil water potential was very low in the topsoil (the soil there was very dry, see shallow layers of Fig. 2); consequently, the soil hydraulic conductivity was small and became the main resistance against water flow between the root xylem and the soil. Adding a soil resistance to the root architecture of the root system consequently had an important effect on the vertical locations of water uptake and release. The low soil water content also explains why the experiment is sensitive to soil hydraulic properties. If kS is too low or λ too high, the soil resistance is too large in the topsoil (there the soil conductivity greatly limits the flow) and HL could not be simulated in the top layers. Because this upper soil layer was hardly conductive, HL appeared deeper in the soil profile where soil resistances were lower. Two main reasons explained why the HL occurred around 20 cm. First, the soil resistance was lower there than in the very top of the soil, while the gradient between the root xylem potential and bulk soil was still high. Second, numerous lateral roots could be found in the soil layers, while below (especially between 30–60 cm) the root density was drastically reduced.


Debate about the magnitude of HL is still ongoing and has important impacts for plant nutrient and water uptake, among others. We designed an experimental setup to monitor HL fluxes during a short time period in a planted soil column using a 1H218O labeling method. An isotopic enrichment of the soil water sampled at the 20- and 25-cm depths was observed after one night of the experiment while the soil water content at the same depths did not change dramatically. This result led to the conclusion that a HL process was observed: isotopically enriched water from the labeled reservoir located at the bottom of the rhizobox and from the rhizobox labeled bottom layers was released by the root system into the shallow layers.

Combining these experimental results with a root water uptake model, we could:

  • Confirm that hydraulic lift was a likely mechanism to explain the peak of isotopically enriched water observed around the 20-cm depth because root and soil parameters optimized via inverse modeling (r2 = 0.98 and RMSE = 0.21‰, equivalent to the accuracy of the isotopic measurements) are in the range of magnitude of expected values;

  • Quantify the water efflux and demonstrate that although no change in soil moisture could be directly measured by soil water probes, the water efflux was not negligible but lower than the probe precision level;

  • Explain which root (root radial conductivity) and soil (soil hydraulic conductivity) factors control the magnitude of the efflux;

  • Show that such experiments can provide valuable estimates of hardy measurable variables (like water fluxes) and root–soil hydraulic properties.

We therefore suggest that further quantitative investigation of hydraulic lift should focus on the hydraulic resistances of the soil–plant system, particularly at the soil–root interface. The need for additional information about both the soil system and hydraulic architecture is also clear from this modeling exercise. If such an experiment is used to retrieve plant hydraulic properties, (i) a good knowledge of the hydraulic conductivity curve of the soil medium (especially under dry conditions) as well as (ii) regular observations of root system growth for measuring the age distribution of the roots would be mandatory. If these measurements are included in the methodology, we expect that the estimates of plant hydraulic conductivities will be more accurate.


The experimental part of this study was financed by the French national scientific program EC2CO/CITRIX. Félicien Meunier is supported by the “Fonds National de la Recherche Scientifique” (FNRS) of Belgium as a Research Fellow and is grateful to this organization for its financial support. Valentin Couvreur is supported by a post-doctoral grant on the PAI MARS P7/29 project.”

Appendix A

Hydraulic Lift in a Simple Root Water Uptake Model

The model of Couvreur et al. (2012) describes water flow in a root system hydraulic architecture under a non-limiting xylem hypothesis. It provides an appropriate framework to understand root water flow under heterogeneous soil conditions (such as the ones imposed on the plants in this experiment). Couvreur et al. (2012) demonstrated analytically that the water uptake of a root segment under heterogeneous soil conditions can be written as the sum of two terms: the “standard uptake” or water uptake under uniform soil water potential conditions and the “compensatory uptake” that occurs only under heterogeneous soil water potential distribution. They correspond, respectively, to the two terms on the right-hand side of the following equation:

where Qi (cm3 d−1) is the water uptake of the ith segment of the root system (defined positive toward the root), Tact (cm3 d−1) is the actual transpiration rate, SUFj (dimensionless) is the standard uptake fraction (defined as the fraction of water taken up by the ith root segment under uniform soil water potentials and is consequently between 0 and 1), Kcomp (cm3 hPa−1 d−1) is the root system compensatory water uptake conductance, Hsr,i (hPa) is the potential at the soil–root interface, and Heq (hPa) is the equivalent soil water potential “sensed” by the plant, defined as a SUF-weighted soil–root interface potential (mathematically, Heq = ΣjSUFjHsr,j). In Eq. [A1], Tact and Hsr,i depend on environmental conditions that can be measured in situ, while Kcomp and SUFi are intrinsic root system properties that can be calculated from the root system architecture and root hydraulic properties. Hydraulic lift occurs when Qi takes negative values. Equation [A1] yields:

It can be observed that two environmental factors promote HL: a low transpiration rate (i.e., small Tact) and large variations in the soil water potential with depth (i.e., Hsr,i very different from Heq, e.g., for low water content in the topsoil and high water content in the layers where the deepest roots are located). To observe HL, the root system must also have a large Kcomp and active root segments (i.e., SUFi > 0) in both dry and wet soil regions, i.e., the soil and root radial conductivity must be high enough to enable root water efflux. Note that we can still observe HL even when the actual transpiration is not zero; then simply the right-hand side of Eq. [A2] must be larger than the actual transpiration rate.

Appendix B

Symbol Description Dimension Unit Values
Ai (ri, li) surface of the ith cylindrical root segment (of radius ri and length li) L2 (L and L) cm2 (cm and cm) 0.3 (0.05,1)
C soil water O stable isotopic concentration M L−3 kg m−3 1.75–2.5
CTi root water equivalent O stable isotopic concentration M L−3 kg m−3 1.89
Hsoil, Hsr,i, Heq soil, soil–root interface, and root system equivalent water potentials P hPa −100,000–0
i, j indices relative to root segment and soil element 1 1–3124, 1–12,160
kr,i radial conductivity of the ith root segment L P−1 T−1 cm hPa−1 d−1 0–0.0005
ksoil soil conductivity L2 P−1 T−1 cm2 hPa−1 d−1 10−12–10
Ki,eq equivalent conductance between the soil node and xylem of the ith root element L3 P−1 T−1 cm3 hPa−1 d−1 10−7–10−1
Kcomp compensatory water uptake conductance L3 P−1 T−1 cm3 hPa−1 d−1 0.0302
krI, krII radial conductivity of primary and secondary root orders L P−1 T−1 cm hPa−1 d−1 0, 0.00005–0.0005
kxI, kxII axial conductivity of primary and secondary root orders L4 P−1 T−1 cm4 hPa−1 d−1 0–0.35, 0–0.00056
Mi, Mw molar masses of 1H216O and 1H218O M L−3 g mol−1 16, 18
N, Np number of observations of isotopic compositions 1 25
Qi, Qsoil radial water flux in the ith root segment, soil water flux L3 T−1 cm3 d−1 −0.05–0.05, −0.05–0.05
Rref Vienna standard mean ocean water (VSMOW) O isotopic ratio 1 2005.2 × 10−6
RLD root length density L L−3 cm cm−3 0.02–40
Sj(Vj) local sink term (in the jth soil element with volume j) T−1 (L3) d−1 (cm3) −1–2.6 (1)
SUFi standard uptake fraction corresponding to the ith root segment 1 0–0.005
Tact plant actual transpiration rate L3 T−1 cm3 d−1 0.05–4
δIW, δLW, δTi stable soil input water, labeled water, and tiller O isotopic composition 1 −8.1–400
δmeas, δsim stable measured and simulated soil water oxygen O composition 1 −8.1–400
θ, θs, θr soil, saturated, and residual volumetric water content L3 L−3 cm3 cm−3 0.08–0.44, 0.4, 0.044
α, n, and m van Genuchten shape parameters L−1, 1, 1 cm−1, –, – 0.0285, 2.29, 0.13
Φsr,i, Φsoil matric flux potential L2 T−1 cm2 d−1 10−12–10
Typical value.
Optimized value.

Appendix C

Radial Flow from or to a Root through a Soil Resistance

The mass balance under steady state imposes that the radial flow from root surface to xylem Qi can be equalized with the water flow from bulk soil to root surface Qsoil (cm3 d−1), given by the axisymmetrical Darcy law (Schröder et al., 2009):where B (dimensionless) is a geometric parameter that depends on the distance between the soil node and the soil–root interface fixed to the mean distance between roots in the simulations (de Jong van Lier et al., 2008) and Hsoil is the potential of the bulk soil linked to the root segment. Note that the products ksoilliB and kr,iAi have the dimensions of conductance (cm3 hPa−1 d−1) and are determined by Ksoil and Kr,i, respectively. Under flux and water potential type boundary conditions at the inner (ri, the root radius) and outer edges (rbulk) representing the bulk soil, and under steady-rate behavior, i.e., (∂θ/∂t) null, the general solution for the radial flux (qr, cm d−1) without the assumption of uniform K(h) is given by (Schröder et al., 2009):whereWhen it is explained in terms of flow, Jr (cm3 d−1), we haveWhen the outer flux qbulk is 0, it becomeswhere the conductance term Ksoil is in cubic centimeters per day per hectopascal.When this soil conductance is associated in series with a root radial resistance, the equivalent resistance isTwo limiting cases can be analyzed. First, if the soil resistance is much higher than the root resistance, the equivalent conductance isThis last expression may be simplified when the root density is high. Then the bulk soil radius gets closer and closer to the root radius and the equivalent conductance becomesbecause the geometric term tends toward 1. If it is the root resistance that much limits the flow, then the equivalent resistance is close to the root resistance:and is independent of the soil resistance.




Be the first to comment.

Please log in to post a comment.
*Society members, certified professionals, and authors are permitted to comment.