Abstract
The estimation of root water uptake and water flow in plants is crucial to quantify transpiration and hence the water exchange between land surface and atmosphere. In particular the soil water extraction by plant roots which provides the water supply of plants is a highly dynamic and non-linear process interacting with soil transport processes that are mainly determined by the natural soil variability at different scales. To better consider this root-soil interaction we extended and further developed a finite element tree hydro-dynamics model based on the one-dimensional (1D) porous media equation. This is achieved by including in addition to the explicit three-dimensional (3D) architectural representation of the tree crown a corresponding 3D characterisation of the root system. This 1D xylem water flow model was then coupled to a soil water flow model derived also from the 1D porous media equation. We apply the new model to conduct sensitivity analysis of root water uptake and transpiration dynamics and compare the results to simulation results obtained by using a 3D model of soil water flow and root water uptake. Based on data from lysimeter experiments with young European beech trees (Fagus silvatica L.) is shown, that the model is able to correctly describe transpiration and soil water flow. In conclusion, compared to a fully 3D model the 1D porous media approach provides a computationally efficient alternative, able to reproduce the main mechanisms of plant hydro-dynamics including root water uptake from soil.
Similar content being viewed by others
Introduction
An approach to model plant water relations has developed during the last two decades based on two main concepts: the cohesion-tension theory (Tyree and Zimmermann 2002) and the electrical circuit analogy applied for simulating water transport in plants using resistances, capacitances, water potentials and flow (Cruiziat et al. 2002). By applying these two concepts, the description of tree hydraulic architecture has made a strong improvement towards a more realistic and comprehensive vision of tree water relationships (Cruiziat et al. 2002). Only recently, this approach was further improved by substituting the electrical circuit analogy by a porous media description for sap flow in wood that is based on Darcy’s law and includes a water capacity term to better account for the dynamic behaviour of the hydraulic storage in trees in a way that mass conservation is considered, e.g. by Bohrer et al. (2005), see also Arbogast et al. (1993), Früh and Kurth (1999), Kumagai (2001) and Aumann and Ford (2002). This avoids calculations of negative water contents and unlimited water withdrawal from the tree as occurring with applications of the electrical circuit model (Chuang et al. 2006). By using a three-dimensional (3D) representation of the tree structure and the physically based representation of tree hydro-dynamics, questions can now be adressed of how the variability in crown architecture (i.e. inter-specific, age or ecosystem dependent structure) would lead to different transpiration responses (Bohrer et al. 2005). In particular by considering the stomatal response to the corresponding branch water potential, we can then simulate the stomatal closure, the corresponding change in leaf water conductance and the resulting actual transpiration rate. Moreover, considering a loss of conductivity with decreasing water potential allows the representation of vulnerability to cavitation at critical highly negative water potential values (Hölttä et al. 2005; Cochard et al. 1999).
Since water transport from the soil through the plant into the atmosphere takes place in a soil-plant-air continuum that is interconnected by a continuous film of water, modelling of plant water transport has to consider both the water exchange at the leaf-air interface and the water flux at the soil-root interface. Therefore, also root water uptake including root hydraulic architecture has to be described to integrate the tree hydro-dynamic model into ecosystem models (Doussan et al. 2003), a step still to be taken to replace current representations of plant water transport, e.g. big leaf or resistor-capacitor approaches (Bohrer et al. 2005) that are often coupled to one-dimensional (1D) effective root water uptake models (Cowan 1965; Gardner 1960; Nimah and Hanks 1973; Feddes et al. 1978; Campbell 1985). New root water uptake models are available that explicitly describe root architecture and related soil-plant processes in three dimensions by explicitly considering the 3D distribution of the uptake (Clausnitzer and Hopmans 1994; Somma et al. 1998; Vrugt et al. 2001) and also considering water flow in the root system (Doussan et al. 2006; Javaux et al. 2008). They can be based on root growth models that allow the integration of a great diversity of environmental conditions and their impact on root system development (Dunbabin et al. 2002).
Despite the more realistic assumptions for predicting soil-root interactions, a disadvantage of the 3D models is their high level of complexity and consequently their high computational demand. Therefore, to describe and simulate transpiration at the stand level a less complex but still realistic approach is needed. By extending the hydro-dynamics model for aboveground parts of trees of Bohrer et al. (2005), we developed a simple root water uptake and transpiration model that couples a 1D soil water flow model with a 1D plant xylem water flow model. Notwithstanding its 1D character the xylem water flow model takes 3D plant architecture into account and can consider different root properties for each root node and different root types. Moreover, by specifying different properties of the soil directly surrounding each root, to a certain extent also horizontal variability might be represented.
Model development
In our model both the water flow within the plant and the water flow in the soil are described by applying the 1D porous medium or Richards equation (Richards 1931). This is based on the assumption that similar to the soil matrix, also the plant xylem can be conceived as a porous medium. The xylem is seen as a bundle of small parallel capillary pores, which are filled with water and air. A capillary can conduct water only if it is completely saturated. With increasing suction head more and more pores cavitate and the ability of the xylem to conduct water decreases. To simulate the xylem water flow and the water uptake from the soil by the roots we need a representation of the plant architecture as model input to define the flow domain and the interfaces between the plant and its environment consisting of the soil and the atmosphere.
Description of plant architecture
The plant architecture is described as a combination of stem, branches, leaves, gross and fine roots including branching roots of several orders (Fig. 1a). The area of leaves on an outer branch is prescribed by a relative leaf area distribution. Stem, branches and roots are divided into slices that are approximated by finite cylinders of different lengths and diameters (Fig. 1a). Each of these cylinders consists of an inner cylinder representing the heartwood, an outer hollow cylinder for the sapwood or xylem and further outer hollow cylinders that include the cambium and the phloem (Fig. 1c). It is assumed that the diameters of these cylinders are maximal diameters in the sense that they do not increase anymore due to an increase in water content of the plant and that a loss of plant water can lead to smaller diameters (Fig. 1d). The flow domain of the xylem is further abstracted and represented by a 3D graph consisting of 1D elements each standing for an inner or hollow cylinder of xylem tissue (Fig. 1b).
Water flow within the plant
The basic state variable to describe water flow in plants is the hydraulic or water potential of the xylem which quantifies the energy status of the liquid phase. The hydraulic potential can be defined as the sum of the xylem potential or xylem hydrostatic potential and the gravitational potential assuming that changes of the osmotic potential are not important for long distance water flow in the xylem (Früh and Kurth 1999). Hence it is assumed that the only component of the xylem potential comes from the energy to move water isothermally and reversibly into or out of the porous matrix represented by the xylem tissue. Because this potential is determined by the porous matrix of the xylem in the following it is called xylem matric potential in analogy to the soil matric potential (Passioura 1980). If the xylem as the water conducting tissue of the plant is conceived as a porous medium (Siau 1984), water flow can be described by Darcy’s law, the relation between the water potential gradient of the xylem and the mass flow of water (Früh and Kurth 1999; Chuang et al. 2006). Following further Früh and Kurth (1999) we assume homogeneity of the xylem hydraulic characteristics, the capacitance or water retention and the axial hydraulic conductivity for each subunit of the plant xylem as represented by the xylem cylinders in Fig. 1a. By further assuming that the xylem volumetric water contents are defined in relation to fixed volumes of the xylem cylinders, we can apply the principle of mass conservation as expressed mathematically by the 1D continuity equation (Früh and Kurth 1999). For the sake of simplicity, in a first step, we neglect water exchange between xylem, phloem, cambium and heartwood and confine the water flow model mainly to the longitudinal direction of the stem and branch axes. By adding a sink-source term to represent transpiration and in our case also root water uptake, we get the following 1D Richards equation as outlined in detail by Früh and Kurth (1999) and Chuang et al. (2006). This equation can then be solved on the domain given by the graph, which represents the plant xylem porous medium by straight lines not only for stem and branches, but also below-ground for roots (Fig. 1b) extending in this way the approach of Bohrer et al. (2005):
where for all of the tree elements represented by straight lines θ x denotes the volumetric water content (m3 m − 3) of the xylem as function of the xylem matric potential ψ x (mm). The xylem matric potential is given on a weight basis, i.e. it is expressed as hydraulic head (mm). t (s) denotes the time, z (mm) the axial length of the element and k x (ψ x ) represents the xylem hydraulic conductivity (mm s − 1) defined as function of the xylem matric potential using the maximal xylem cross sectional area as reference surface of the water flux. The vertical position is given by the height above (positive upward) or the depth below the soil surface (negative downward). For a branch or root element α x is the zenith angle (–) and S x (s − 1) the sink or source term, which refers to roots for root water uptake and outer branches for transpiration. In a further model development also the radial exchange of water between the xylem and the phloem, the cambium and the heartwood of the tree might be added to this term.
Hydraulic characteristics of the xylem
Due to the form of the mass conservation principle involved in the formulation of the Richards’ equation, xylem water contents are defined in reference to fixed xylem volumes, i.e. in our case to the maximal xylem volumes. Therefore, if one wants to apply Richards equation and also consider xylem desaturation which keeps the xylem vessels still saturated but leads to significantly smaller xylem volumes, then this desaturation has yet to be described in relation to the maximal xylem volumes. In this case the usual form of soil water retention curves, that are nearly constant between zero matric potential and air entry value, cannot be applied to describe the xylem water–xylem matric potential relationship on a maximal sapwood volume basis. Hence to account for changes in xylem volume, in view of measurements for European beech (Fagus silvatica L.) branches by Oertli (1993), see Fig. 3a, we assume a linearly decreasing xylem water retention curve between zero matric potential at maximal saturation and air entry value. Thus we extend the xylem water retention curve of Chuang et al. (2006) describing for each tree element e of stem, branch or root the dependency of the volumetric xylem water content on the xylem matric potential at a certain height or depth z (mm) by
where the maximal volumetric water content of the xylem is assumed to be equal to the xylem sapwood porosity ϵ x (m3 m − 3) given for each tree element by the fraction of the maximal volume of xylem water at saturation to the maximal total volume of the xylem sapwood in the tree element or due to the assumed cylindrical geometry of the xylem elements by the fraction of the maximal water-filled xylem sapwood area s xw, max ,e (m2) and the maximal xylem sapwood area s x, max ,e (m2) of the xylem element e. In a similar way we define the xylem water content at air entry θ x,a (m3 m − 3) by use of the corresponding water-filled cross sectional xylem sapwood area s xw,a (m2) at the xylem matric potential when cavitation begins to occur:
Assuming that tension changes the volume of water conducting sapwood only (Perämäki et al. 2001) and neglecting corresponding changes in phloem water content, which are considered to be comparatively small, then, in combination with the xylem water retention curve, the xylem water content θ x,a at air entry can be derived using the relation between xylem matric potential and relative volume change of the xylem element (Steppe and Lemeur 2007) expressed by
where E (mm) denotes the elastic modulus of the xylem and V e the volume of the xylem element. By evaluating this relation, which is valid in the range of xylem water contents between ϵ x and θ x,a when the vessels are saturated, we get
and hence a value for the xylem water content at air entry can be determined for a given value of the elastic modulus. Moreover, inserting Eq. 5 into Eq. 2 we get
and from this the actual xylem cross-sectional area s x,e (m2):
We then can calculate the actual xylem diameter exploiting the assumed cylindrical volume geometry of the xylem element e (Fig. 1d).
The elastic modulus E (mm), the porosity ϵ x (m3 m − 3), the air entry value a (mm) and the exponent λ (–) that determine the xylem water retention curve, are model input parameters and assumed to be fixed values for the whole tree.
From the xylem water retention curve we can derive the xylem hydraulic conductivity k x based on the law of Hagen and Poiseuille for the mass flow rate of water in a cylindrical pipe. According to Burdine (1953) by considering capillary bundles of such pipes, the relative conductivity K(θ x ) (–) of the volumetric water flux in the capillary bundle as defined by Eq. 11 can be calculated by the following integral equation (Hoffmann-Riem et al. 1999):
Here also other functional forms such as Mualem’s equation (Mualem 1976) might be applicable. To account for the tortuosity of soil pores usually p = 2 is assumed, but for the sake of simplicity and lack of appropriate data to estimate an adequate value of p we assume p = 0, although tortuosity of water films in xylem vessels might occur during embolism. The evaluation of Burdine’s integral equation for θ x ≤ θ x,a and the application of the Hagen-Poiseuille law to consider the effect of an increase ΔR of the effective mean hydraulic radius of the capillary bundle for θ x > θ x,a leads to the following function of the xylem water content to represent the xylem hydraulic conductivity:
For the saturated case, i.e. for θ x > θ x,a we consider the increase of the mean effective hydraulic radius ΔR of the capillary bundle with increasing radii of the water-filled tubes resulting in xylem volume increase, and get for the xylem radii r = r(θ x ), r max = r(ϵ x ), r a = r(θ x,a) (mm) at the respective xylem water contents θ x , ϵ x and θ x,a and the mean value θ m = (ϵ x − θ x,a)/2:
To account for the change of cross-sections between tree elements, for each tree element at depth or height z the hydraulic conductivity is calculated in relation to the water-filled cross sectional area of the xylem sapwood (Bohrer et al. 2005), normalised by the maximal water-filled cross sectional xylem sapwood area of the stem element at stem-base, i.e. at height zero s xw, max (0) (m2):
where k max (mm s − 1) denotes the maximal hydraulic conductivity of the xylem. k max is an additional input parameter, which for the sake of simplicity is assumed to have fixed values for either shoot and root. For young European beech trees we apply in the following different values for k max (see Table 1).
Water flow in the soil
To simulate vertical soil water flow we apply again the 1D Richards equation, an approach frequently used to model soil-plant-atmosphere systems (Priesack and Gayler 2009; Simunek et al. 2008; Kroes and van Dam 2003):
where θ s is the volumetric soil water content (m3 m − 3) as a function of the soil matric potential ψ s (mm), t (s) denotes time and z (mm) soil depth (here positive downward). k s (ψ s ) (mm s − 1) is the soil hydraulic conductivity which is given as a function of the soil matric potential. Both soil hydraulic property functions θ s and k s are described by parameterisations according to van Genuchten (1980). The last term S w (s − 1) represents the sink term of root water uptake (per unit soil depth). Typical boundary conditions of the flow equation are a mixed boundary condition at the soil surface to represent atmospheric conditions leading either to infiltration, evaporation or water ponding and a free drainage condition at the bottom boundary (Simunek et al. 2008; Kroes and van Dam 2003).
Soil water uptake and transpiration
Water flow between the soil-root interface and the root xylem is defined by the 1D soil-root volumetric flux j r,e (mm3 s − 1) which describes the volume of water exchanged between a root element and a soil layer per unit of time
where k r (s − 1) is the radial conductivity of the root through the root cell structure between root xylem and soil and s r,e(z) (mm2) the outer surface of the root element e within the soil layer at depth z (mm). ψ s (z) and ψ x (z) denote the matric potentials in soil at the root surface and in the xylem (mm) at soil depth z. Instead of applying only the radial conductivity k r , we alternatively consider the rhizosphere radial conductivity k rs (s − 1), defined by
the geometric mean of the radial conductivity k r of the root and the soil hydraulic conductivity k s at the soil matric potential near the roots divided by the radial thickness l rs (mm) of the rhizosphere soil around the roots assumed to be 2.0 mm (Hinsinger et al. 2005). To account for the impact of the strong root-soil interaction in the rhizosphere on the water status of the soil directly surrounding the roots, we take instead of the soil water matric potential ψ s a weighted mean \(\hat{\psi}_s\) between xylem and soil matric potentials
where γ 1 (–) and γ 2 (–) are weighting constants.
The calculated flux j r,e into or out of the cylindrical root element e is then used to calculate the below ground water uptake sink or source terms S x and S w of the Richards equations of xylem and soil water flow. The below ground xylem sink or source term S x (mm3 mm − 3 s − 1) is obtained by relating the fluxes to the corresponding xylem cylinder volumes V e = s x, max ,e(z) l e of the root elements of maximal xylem sapwood area s x, max ,e (mm2) and element length l e (mm):
The water fluxes between soil and root cylinders located in a certain numerical soil layer i are related to the corresponding soil volume ΔV i given by the thickness of the soil layer Δz i and the unit area of a square meter and are summed up to give the sink or source term S w,i (mm3 mm − 3 s − 1) of the soil layer i. If the axis of a root cylinder element e intersects with more than one of the numerical soil layers, the water flux is divided between these layers proportionally to the fractions f e,i (–) of the cylinder axis intersections with the particular layers:
Similar to water uptake, also transpiration, the loss of water to the atmosphere, is described by a sink term. At the outer branches, where the leaves are located, the corresponding water flux is prescribed by the atmospheric evaporative demand and by the stomatal hydraulic conductivity, which similar to the xylem hydraulic conductivity is assumed to depend on the xylem matric potential of the branch element (Bohrer et al. 2005). For simplicity the atmospheric water demand is prescribed by a constant flux condition or calculated from the potential evapotranspiration rate (mm s − 1) estimated by the FAO Penman-Monteith method (Allen 2000) and is distributed according to the leaf area. Therefore, in contrast to the water exchange with the soil, the exchange with the atmosphere is prescribed and no direct feedback of transpiration to the micro-climate at the leaf level is considered.
The prescribed flux condition or the potential transpiration T pot (mm s − 1) of the total tree is then distributed to single branch elements according to the relative leaf area distribution LAD e (z) (m2 m − 2) at each element e, given by the leaf area of all leaves that are connected to the branch element divided by the total leaf area of the tree. Hence the potential transpiration T pot,lf,e (mm s − 1) at the branch element e can be given by:
The stomatal hydraulic conductivity k x,stom (mm s − 1) is modelled in dependence of the xylem matric potential in the leaf ψ x,lf (mm) by
following Bohrer et al. (2005) with a stomatal maximal hydraulic conductivity k max ,stom = 1.0 mm s − 1 and stomatal hydraulic conductivity parameters b (mm) and c (–) and by additionally assuming ψ x,lf = ψ x (z ep ) at the end point z ep of the branch element where the leaf is anchored. The actual prescribed transpiration rate T act,lf,e (mm s − 1) is then calculated from the corresponding potential transpiration rate T pot,lf,e (mm s − 1) at the branch element by
for which the maximal evaporation rate from the leaf q max ,lf (mm s − 1) at the actual leaf xylem potential ψ x,lf is given by
where \(\psi_{0,lf}= - 3.0\ \cdot 10^5\) mm denotes a potential that corresponds to a minimal potential in the leaf, which can be reached under dry air conditions. Δz ep is the length of the element belonging to the end point where the leaf is anchored.
Finally, the above ground xylem sink term S x,e (mm3 mm − 3 s − 1) which represents transpiration of the outer branch element e is calculated by relating the volumetric water flux out of the xylem given by T act,lf,e s x, max ,e (mm3 s − 1) to the xylem volume V e (mm3) of the branch element:
Numerical implementation
The soil and tree water flow equations are coupled via the below ground water uptake sink terms and each solved by iterations using a standard Newton–Picard fixed point iteration scheme (Celia et al. 1990; Priesack 2006) until tree and soil water contents and potentials converge below a threshold value (Huang et al. 1996). Although both flow equations may be coupled by a further fix-point iteration via the water uptake term to account for the non-linearity of the systems, similar to Tseng et al. (1995) we use a single-pass solution scheme without iteration between the two subsystems, since for the considered simulation scenarios the additional iteration does not lead to significantly different results. The solution procedures of both flow equations 1 and 12 are based on a finite element discretisation following the approach of the HYDRUS-1D model (Priesack 2006; Simunek et al. 2008). Since the solution method of the Richards equation on a graph as represented in Fig. 1b has not been documented in more detail, we describe the applied finite element method in Appendix A.
Statistical criteria
Three statistical criteria are used to assess the adequateness of simulations. The root mean square error RMSE as defined by Mayer and Butler (1993) is used to assess over the whole period of the lysimeter experiment the deviation between predicted, P i , and observed values, O i , proportioned against the mean observed value \(\bar{O}\):
where n is the number of observations during the experiment. The lower the value of RMSE, the better is the agreement between simulation and observation. The Modelling efficiency ME (Willmott 1982) is used to assess the model adequacy in describing the observed variability.
where P j and O j are the predicted and observed values and \(\bar{O}\) and \(\bar{P}\) are the corresponding mean values. A value of ME < 0 indicates that \(P_j = \bar{O}\) would be a better estimation for observed variability than the simulation results. If simulated values are completely in accordance with observed values, ME equals 1.
Finally, a necessary criterion to evaluate the accuracy of a numerical scheme is given by the global mass balance error \(MB = \sum_{m=1}^{\rm end} MB(t^m)\) at the end of the simulation, where at each time step m we have for the water balance of the soil-plant system:
and where \(\theta_l^m\) denotes the volumetric water content at time step m of the tree at element l = e or of the soil at node l = i, ΔV i is the volume of the soil layer i and \(Q_0^j\) and \(Q_k^j\) denote the volumetric water fluxes (mm3 s − 1) across the boundaries of the soil profile.
Simulation scenarios and parameterisation
The new model of water flow in soil-plant systems is applied to three different scenarios mainly to illustrate that important features of water flow in plants and soils can be simulated. In a first scenario, the drying soil scenario, we simulate soil water uptake and transpiration under conditions of a prescribed constant transpiration rate without refill of soil water. By this scenario the effects of stomatal conductivity, xylem conductivity, radial root conductivity and soil conductivity on water uptake and transpiration are considered. In particular the onset of cavitation is examined. In a second scenario, the hydraulic lift scenario, we demonstrate the basic ability of the model to simulate transport of water through the roots from wetter to dryer soil regions. In a third scenario, the lysimeter scenario, we test if the model can simulate daily and seasonal water balances based on measured input parameters. In this case for the xylem water flow model only parameters of the tree architecture are changed: the stem diameter and the distributions of leaf and root area index based on measurements from the lysimeter experiment.
In the drying soil scenario, we simulate soil water uptake and transpiration also for comparison with a 3D root water uptake modelling approach. In the scenario, the water is taken up by the roots from a homogeneous soil column of 0.4 m depth and of 0.1 m by 0.1 m soil surface, which is assumed to be initially in hydrostatic equilibrium with an aquifer located 3.0 m below the soil surface (Javaux et al. 2008). The soil boundary conditions are a no flux condition at the bottom of the soil profile and no rainfall or evaporation at the soil surface. For reasons of comparison with the numerical experiment of Javaux et al. (2008) we use the same parameters to study the effect of different soil types and to analyse the sensitivity of different soil and xylem hydraulic conductivity values (Table 1). For the reference case input parameter values are taken from the second column of Table 1, and for the sensitivity analysis only the conductivity parameters, k sat for the soil, k max for the xylem and k r for root water uptake, were changed as indicated in the third column. This perturbation of parameter values corresponds to a realistic variation in observed conductivity values for soil or xylem (Javaux et al. 2008; Doussan et al. 1998). In both cases for the above-ground tree a constant potential transpiration rate of 1.0 mm d − 1 was assumed and prescribed along the tree at the ends of stem or branch elements according to the assumed leaf area distribution. Related to the 0.01 m2 soil column surface this corresponds to a volumetric water flux rate of 1.0 ·104 mm3 d − 1 from the soil to the atmosphere.
In the hydraulic lift scenario, instead of the no flux of the first scenario we consider a Dirichlet boundary condition which prescribes saturated conditions at the lower boundary of the soil column now assumed to be 1.5 m deep. Furthermore, we now prescribe a daily sinusoidal cycle of the potential transpiration according to (Childs and Hanks 1975; Priesack 2006), but keep the daily rate at 1.0 mm d − 1. All other conditions and parameters are kept the same as in the loam reference case of the first scenario.
In the lysimeter scenario, the soil-plant systems are represented by 2–3 years old European beech trees growing in cylindrical lysimeters of 2.0 m depth and 1.0 m2 surface area. Here the soil boundary conditions are a seepage flow condition at the bottom of the soil profile and an atmospheric flow condition at the top soil prescribing either infiltration, evaporation or water ponding. The third scenario is based on a lysimeter study about effects of elevated ozone and below-ground pathogen infection on juvenile European beech (Winkler et al. 2009). From this study we take hourly and daily water balance data of one of the control lysimeters, i.e. lysimeter S3, that has seen neither an ozone nor a pathogen treatment. The input data that characterize the hydraulic properties of the lysimeter soil profile are derived from measurements given in Gayler et al. (2009) and summarized by Table 2. The stem diameter d stem is 15.8 mm, the xylem conductivity parameter values are the reference values given in Table 1.
In all three scenarios we use a simplified plant architecture of a young European beech tree (Fig. 2) obtained using digital photographs and photogrammetric evaluations. In the second scenario we extend a central root to the bottom of the soil column to facilitate water uptake from deeper soil layers and in the third scenario the stem diameter is increased keeping the diameter fractions between stem and branch elements the same, whereas the root system is now adapted to represent the measured root biomass distribution.
We estimate the maximum water-filled cross sectional xylem sapwood area s x, max of each tree element e from the diameter d e (mm) of the tree element according to Aranda et al. (2005). For young European beech trees we fit the new parameterisation of the xylem water retention curve given by Eq. 2 to the measured retention data of Oertli (1993) for shoots of European beech. This fit results in an elastic modulus value of E = 1.5 ·106 mm and an air entry value of a = − 3.0 ·105 mm with exponent λ = 0.854 for each tree element (Fig. 3). Alternatively to the xylem water retention curve given by Eq. 2 we also apply and fit a van Genuchten parameterisation (van Genuchten 1980) obtaining the parameters α = 2.8 ·10 − 6 mm − 1, n = 2.77, θ sat = ε x and θ res = 0 (Fig. 3).
Results
Numerical accuracy and convergence
All numerical simulations were conducted on an ordinary personal computer with a single core 1.5 GHz AMD processor. CPU time is given in seconds (Table 3). During each simulation the water mass balances are recorded to compare for every individual time step the input-output balances of the whole tree respectively of the soil profile with the corresponding internal mass change in the tree or soil. For the simulation of the soil drying scenario using the reference parameters the global mass balance error of the soil-plant system is MB = 1.3 ·10 − 3%, for the hydraulic lift scenario MB = 3.1 ·10 − 3% and for the lysimeter scenario MB = 3.2 ·10 − 4%. Table 3 summarizes results obtained for the xylem water flow simulation in the drying soil scenario using the reference parameters. By relaxation of the iteration tolerance δ tol it is seen, that already at the criterion of δ tol = 0.01 a rather low global mass balance error is reached without strongly increasing CPU time, iteration and time step numbers. In contrast, the CPU time is considerably increased by increasing the number of tree elements up to 482 or 903 by prescribing a maximal axial element length of only 20 or 10 mm instead of 100 mm as in the reference case. This increase indicates the first order accuracy in space of the applied numerical solver.
Drying soil scenario
Figures 4 and 5 show the simulated course of transpiration and root water uptake with related plant water storage for different parameterisations of the stomatal hydraulic conductivity curve. Because initially, at simulation start, tree water content is not in equilibrium with soil water and atmospheric conditions, the tree first looses water to the atmosphere before the hydraulic gradients of xylem water potentials get strong enough and force the roots to take up the available water from the soil. Since a constant maximum water flux is prescribed at the tree leaves and no water is added to the soil, the water flow in the soil-plant system quickly reaches an almost steady state as long as the soil provides enough water to fulfil the prescribed flux. After some time, when the soil gets dry, the root water uptake rate decreases due to decreased soil hydraulic conductivities and low soil water availability. This low water availability induces low xylem water potentials and consequently low stomatal hydraulic conductivities such that the transpiration rate is reduced. The transpiration rate is strongly determined by the assumed parameterisation of the stomatal hydraulic conductivity curve which reflects the importance of stomatal control on xylem water flow. For the lower value of the stomatal hydraulic conductivity parameter b the stomatal hydraulic conductivity decreases already at higher xylem potentials and hence protects the tree from water losses earlier during the drying experiment (Fig. 4). If we then assume that the closure of the stomata would lead to zero stomatal hydraulic conductivity the plant could almost keep its water content (Fig. 4b, solid line) without reaching the air entry value. But if we assume a maximal possible reduction of stomatal hydraulic conductivity to only 10% of its maximal value, then the plant internal water is almost completely transpired (Fig. 4b, dotted and dashed lines) and cavitation occurs. Generally, a higher b value leads to a higher stomatal hydraulic conductivity also at lower xylem potentials. Therefore higher water losses of tree water storage take place, which results in critically low tree water contents. In our example in three cases xylem water content falls below the air entry value θ x,a = 0.8, so that cavitation is induced (Fig. 4b). In the cases where the stomatal hydraulic conductivity remains above a non-zero limit, the transpiration gets finally limited when the plant internal water is almost exhausted (Fig. 4a).
Due to the evaporative water flux from the leaves to the atmosphere, the simulated xylem matric potentials show a strong vertical gradient along the tree height (Fig. 2) having lower potentials at the tree crown, where the leaf area is highest and higher potentials towards the roots and within the roots which take up the water from the soil.
The perturbation of conductivity values gives only very small differences between the simulations of the reference case and the cases with one-tenth- or ten-fold maximal xylem conductivity (Fig. 5). This small deviation is similar to results obtained by Javaux et al. (2008) when using the same reference conductivity values for root water uptake simulations (Table 1). Additionally, no difference can be detected if we neglect the second term of the relative conductivity curve in the case θ x > θ x,a, i.e. if we neglect the impact of an increased mean radius of the water-filled xylem vessel elements. Moreover, for the small trees the exact form of the conductivity curve, which is difficult to evaluate considering the high variability of the measured conductivity data (Fig. 3), may play a less important role, if only the decrease of the curve with decreasing xylem matric potential is appropriately described. Because of the low sensitivity of maximal xylem conductivity variation on the simulated actual transpiration in case of the considered small trees, the relative xylem conductivity curve given by Eq. 9 could be further simplified in the saturated range above the air entry value by omitting the second term of the right hand side without significantly changing our simulation results of xylem water flow dynamics (results not shown). This low sensitivity also means, that the impact of the elastic change on the diameter and hence the conductivity of water-filled tubes of the capillary bundle representing the xylem can be neglected in our model description of xylem water flow in small European beech trees. In this case the low sensitivity of the axial xylem conductivities seems to be caused by their relatively high values compared to the low stomatal and radial root conductivity values.
When compared to the reference case in the cases of increased soil or radial root conductivity the begin of limited water supply for root water uptake is 2–5 days delayed (Fig. 5a). Applying the radial rhizosphere conductivity the delay of the onset of limitation of transpiration is 9–10 days (Fig. 5b). In the sequel, the higher conductivities lead to earlier soil water exhaustion and reduced transpiration rates below 0.2 mm d − 1.
To represent effects of horizontal soil variability we introduce a possibility to account for different soil properties in the direct surroundings of different roots or root elements. This representation is achieved by considering the radial rhizosphere conductivity defined by Eq. 14. The root water uptake rates for the three soil types applying either the constant radial root conductivity k r or the rhizosphere conductivity k rs are shown in Fig. 5c and d. Generally, we observe, that if we apply the radial rhizosphere conductivity k rs instead of the radial root conductivity k r in Eq. 13, the simulated root water extraction decreases earlier and remains slightly higher when the soil gets drier (Fig. 5b and d). For k r the root water uptake limitation appears first for the clay, then for the loam and eventually for the clay loam (Fig. 5c). For k rs the picture changes in case of the loam soil, since the time-span until transpiration gets limited is now shortest for this soil type (Fig. 5d). In both cases the clay loam soil supplies the transpiration water demand by far the best, up to 15 d longer than the clay and loam soils.
The simulated depth distribution of root water uptake is shown by Fig. 6. Under non limiting soil water conditions the root demand is met by available soil water and the depth distribution of the uptake rate almost coincides with the root area index distribution. Only in the case of the increased radial root conductivity available soil water gets already scarce in the region of high root area index, i.e. m2 root surface area m − 2 soil surface area (Fig. 6a). In this region water uptake by roots begins to be limited after 15 d, and is already very low after 30 d (Fig. 6b and c). Then water is increasingly taken up from regions of lower root area index, where soil water is not completely exhausted and still can satisfy almost half of the prescribed demand (Fig. 6d).
In the cases with increased soil and radial root conductivity, similar to the 3D root water uptake simulations of Javaux et al. (2008), the limitation of root water uptake is significantly delayed by 2–5 days (Fig. 5a). Furthermore, the finding that the clay loam best sustains the evaporative demand of the plant, about 10–15 d longer than the other soil types could be confirmed by our simulations (Fig. 5c and d). But, in contrast to the 3D water uptake simulations under the second collar boundary condition (CBC 2) by Javaux et al. (2008), in our 1D simulations water stress occurred earlier for clay than for loam soil, if the radial root conductivity was used (Fig. 5c). Only if the radial rhizosphere conductivity is applied, the same order of soil types concerning the appearance of water stress is found, i.e. transpiration is reduced first for loam, then for clay, and eventually for clay loam soil (Fig. 5d).
Hydraulic lift scenario
Plant roots do not only take up water but can also release water to dry soil (van Bavel and Baker 1985; Richards and Caldwell 1987). In this way the plant root system may transport water from moist into drier, often upper soil horizons and hydraulic lift occurs. Assuming non polar flow across root membranes we are interested to see if our proposed soil-plant water flow model is in principle capable to simulate such a hydraulic lift. Under the conditions of the second scenario with periodic daily transpiration cycle the simulated water uptake rates show a diurnal variation (Fig. 7). During the day when potential transpiration rates are above zero, water is taken up mainly from the upper soil at 0.0–0.4 m depth (Fig. 7a), where most of the root biomass is located. During the night, when no transpiration occurs, deep roots take up water from the bottom of the soil column, where the soil is nearly water saturated, replenish xylem water and additionally release water to the drier soil near the top of the soil column also at 0.0–0.4 m depth (Fig. 7a). Simulations considering a time span of more than 240 days show that water flow in the soil-plant system gradually reaches a quasi steady state situation, where water is taken up from the lower part of the soil column and redistributed to the upper soil during the night (Fig. 7b), whereas during the day mainly water from the upper soil layer is taken up by the roots to accommodate transpiration demand (Fig. 7c). Moreover, compared to a simulation where we assume no water flux from the roots to the soil, the daily transpiration with hydraulic lift is about 5% higher during the first 100 days.
Lysimeter scenario
A more realistic scenario is given by the lysimeter scenario. Based on the measured daily precipitation, the estimation of daily potential evapotranspiration by the Penman–Monteith dual crop coefficient method (Allen 2000; Loos et al. 2007; Gayler et al. 2009) and the partitioning of potential evaporation and potential transpiration according to Allen (2000), the daily actual evapotranspiration, the daily changes in water storage and the daily percolation rates were simulated and compared to the corresponding measurements (Fig. 8). Results show that the model is able to reproduce the actual evapotranspiration and the changes in water storage quite well, with a model efficiency ME of 0.901 and a root mean square error RMSE of 0.298 for the evapotranspiration, and with the values ME = 0.977 and RMSE = 0.586 for the water storage change. The simulated dates and sizes of the percolation peaks differ more strongly from the measured data giving the values ME = 0.845 and RMSE = 1.812, in particular, the simulated percolation peaks seem to occur always somewhat later than the observed peaks.
Figure 9 shows the simulated diurnal course of transpiration, stem diameter, water flux at the root collar and the rate of change of xylem water content during three consecutive days of the lysimeter experiment. Transpiration and water flux at the root collar strictly follow the prescribed sinusoidal potential evapotranspiration. After the onset of transpiration, the rate of change of tree water content decreases due to the onset of transpiration but soon increases because root water uptake begins to refill plant water capacity, before it decreases again indicating the end of water replenishment of the xylem water by soil water taken up during the night. Correspondingly, the xylem diameter changes with a maximal diameter difference of 0.194 mm between day and night.
Discussion
A new model approach is introduced to describe water flow in the soil-plant system. The model extends previous approaches that were developed to describe hydro-dynamics in above-ground plants (Früh and Kurth 1999; Aumann and Ford 2002; Bohrer et al. 2005) by additionally including the hydraulic root system architecture in combination with a soil water flow model. In this way a complete water flow model for the whole soil-plant system is obtained, which then vice versa extends existing root water uptake modelling approaches (Doussan et al. 1998; Javaux et al. 2008) by including the above-ground water dynamics of the plant. Similar to the soil water flow model, which is based on the continuum approach as described by the Richards’ equation, also the plant water flow is described by this porous media equation following Bohrer et al. (2005). This description is in contrast to more discrete approaches that consider the plant as built up of storage compartments and describe the water flow by water exchange processes between these compartments (Perämäki et al. 2001; Steppe and Lemeur 2007). But, besides the advantages of a continuous mathematical formulation representing mass conservation (Früh and Kurth 1999), due to the similarity between the porous media approaches to describe water flow in the soil and the xylem, it seems also to be conceptually easier to derive a consistent model of hydro-dynamics of soil-plant systems based on the continuum approach.
In a first analysis, the new model was applied to investigate model behaviour and model sensitivity by scenarios with young European beech trees that have already been analysed in two previous studies (Javaux et al. 2008; Gayler et al. 2009). Model tests and numerical results show, that the model is successfully implemented and can efficiently simulate water flow in the soil-plant system based on the a simple representation of plant architecture. Mass balance errors are low and in an acceptable range. But since the global mass balance is a necessary but not a sufficient criterion for having a correct solution (Celia et al. 1990), the model needs further analysis, which may be based on comparison with exact solutions (Lehmann and Ackerer 1998) or as in our case by evaluating simulations of practical examples using experimental data.
In the case of the 3D simulation by Javaux et al. (2008) the effect of the steeper slope of the loam hydraulic properties generates a heterogeneous distribution of soil water potentials in the horizontal direction with drier soil near the roots. This distribution leads to a larger conductivity drop and produces an earlier water stress (Javaux et al. 2008). In the case of our 1D simulation this effect is mimicked by amplifying the direct impact of the lower xylem potential on the soil matric potentials directly near the roots resulting in a lower radial hydraulic conductivity between root axis and soil, which is modelled by Eq. 14 based on the definition of the radial rhizosphere conductivity k rs . If the distribution of radial conductivities of the roots is adapted in this way, the main conclusions about effects of soil type and conductivity values remain, although the root architecture was different in the drying scenario simulations we compared. But amplitudes or local values of considered variables are different due to the 1D simulation of soil water flow. Therefore, a full data comparison cannot be performed, even if locally the water uptake by the root system may be of similar magnitude.
Besides the model sensitivity of the radial root or rhizosphere conductivity, the most significant model sensitivity on simulated soil water uptake and transpiration is due to the stomatal hydraulic conductivity. This model behaviour is mainly caused by eventually rather low conductivity values, reflecting the significant role of both conductivities in controlling tree water flow which has been shown in several experimental studies (Steudle 2000; Lemoine et al. 2002; Sperry et al. 2003). Moreover, the water has to flow through nonxylary pathways in the root across the root cortex and endodermis and in the leaf through the mesophyll. These flow paths have low hydraulic conductivities and can account to a high proportion of the overall root respectively shoot hydraulic resistance (Sperry et al. 2003). In particular, the radial root conductivity may be strongly changed by aquaporin activity (Ehlert et al. 2009), which itself may be controlled by metabolism or be triggered by environmental factors (Steudle 2000). Also changes of leaf hydraulic conductivity may be mediated by plasma membrane aquaporins (Cochard et al. 2007).
Hydraulic lift has been shown to occur in more than 60 species (Espeleta et al. 2004), “but there is no fundamental reason why it should not be more common as long as active root systems are spanning a gradient in soil water potential and that the resistance to water loss from roots is low” (Caldwell et al. 1998). However, for European beech trees it seems to be an open question whether the radial hydraulic root conductivity from inside the roots towards the soil is large enough to permit considerable water flow to dry soil (Rewald 2008). Therefore our simulations of hydraulic lift only show that in principle the model is able to describe this kind of water redistribution from wet to dry soil under almost optimal conditions for hydraulic lift, i.e. assuming water saturation at the bottom of the soil profile with a resulting simulated leaf xylem matric potential ranging from −40,000 to −70,000 mm. Nevertheless, the hydraulic lift could be significant in mitigating plant water stress and need to be accounted for in models describing root water uptake and the onset of plant water stress in a water limited ecosystem (Siqueira et al. 2008).
The lysimeter scenario simulations confirm the applicability of the new model to describe daily water balance dynamics in soil-plant systems. Also the diurnal course of transpiration seems to be realistically described. In particular, the simulated diurnal changes of xylem diameter up to 0.194 mm (Fig. 9) are within the range of xylem or stem diameter changes that were observed for young European beech and Norway spruce (Pinus abies L.) trees (Steppe and Lemeur 2007; Perämäki et al. 2001). This agreement results from the fact that based on the parameterisation of the xylem water retention curve obtained from the experimental curve of Oertli (1993), we get a value for the elastic modulus E which is similar to values found by Steppe and Lemeur (2007) and Perämäki et al. (2001). However, this value might have been underestimated, since we assume that changes of xylem water potential correspond only to changes in the xylem water content neglecting related changes in water contents of phloem, cambium and heartwood.
Since the new plant water flow model focuses on the water flow in the sapwood tissue in the longitudinal direction of the stem and branch axes, the model still needs to be extended to describe water flow in the phloem and to account for the radial exchanges of water between xylem, heartwood, phloem and cambium. In particular the active uploading and unloading of sugars according to Münch’s hypothesis (Münch 1928) needs to be considered, because near the sources of sugars considerable water flow from the xylem to the phloem can occur and vice versa close to sugar sinks significant amounts of water can flow from the phloem to the xylem. Locally this water exchange can strongly influence the xylem water flow dynamics which also can lead to xylem diameter changes (Sevanto et al. 2002).
At the present state of development the model is supposed to be already useful to simulate the impact of differences in plant architecture and soil properties on transpiration and water uptake under various climatic conditions, in particular considering differences of leaf area and root area distribution or of root and shoot branching. The model may be also helpful to describe the impact of different distributions of conductivity values along the tree architecture in particular of larger trees to study their role in protecting the xylem against cavitation, while maintaining transpirational water flow. This description may lead to a better understanding of observed differences in xylem structures and functions between deep and shallow roots, roots, stem and branches (Nadezhdina 2010; North 2004) and of the relevance of nonxylary pathways for water uptake and transpiration (Sperry et al. 2003). The model may also be applied to better describe differences of stand water budgets between mixed forest stands of different species composition.
Conclusion
The newly proposed model of water flow in soil-plant systems presented in this study combines approaches to simulate water flow in the above-ground plant based on the porous media equation with root water uptake and soil water flow models also based on the porous media equation. In this way the model concept exploits the similarity between water flow in the soil and the xylem of plants. It could be shown that the model can simulate water balance dynamics in soil-plant systems with young European beech trees including some main known features of water flow in plants and soil such as the diurnal cycle of plant water content and xylem diameter change, stomatal control of transpiration, cavitation due to water stress, sensitivity of root water uptake to radial root conductivity, hydraulic lift, interaction of hydraulic soil properties, root distribution and plant water availability. Due to the simplified representations of plant architecture and distribution of soil properties that lead to 1D representations of water flow domains, the model is numerically very efficient and may therefore be also applicable to simulate water flow in canopies accounting for spatial and temporal variability of canopy structures based on rather simple representations of individual plant architecture.
References
Allen RG (2000) Using the FAO-56 dual crop coefficient method over an irrigated region as part of an evapotranspiration intercomparison study. J Hydrol 229:27–41
Aranda I, Gil L, Pardos J (2005) Seasonal changes in apparent hydraulic conductance and their implications for water use of European beech (Fagus sylvatica L.) and sessile oak (Quercus petraea (Matt.) Liebl) in South Europe. Plant Ecol 179:155–167
Arbogast T, Obeyesekere M, Wheeler MF (1993) Numerical methods for the simulation of flow in root-soil systems. SIAM J Numer Anal 30:1677–1702
Aumann CA, Ford DE (2002) Modeling tree water flow as an unsaturated flow through a porous medium. J Theor Biol 219:415–429
Bohrer G, Mourad H, Laursen TA, Drewry D, Avissar R, Poggi D, Oren R, Katul GG (2005) Finite element tree crown hydrodynamics model (FETCH) using porous media flow within branching elements: a new representation of tree hydrodynamics. Water Resour Res 41:W11404. doi:10.1029/2005WR004181
Burdine N (1953) Relative permeability calculations from pore size distribution data. Trans AIME 198:71–77
Caldwell MM, Dawson TE, Richards JH (1998) Hydraulic lift: consequences of water efflux from the roots of plants. Oecologia 113:151–161
Campbell G (1985) Soil physics with BASIC, Transport Models for Soil-Plant Systems, Developments in Soil Science, vol 14. Elsevier, New York
Celia M, Bouloutas E, Zarba R (1990) A general mass-conservative numerical solution for the unsaturated flow equation. Water Resour Res 26:1483–1469
Childs S, Hanks R (1975) Model of soil salinity effects on crop growth. Soil Sci Soc Am Proc 39:617–622
Chuang YL, Oren R, Bertozzi AL, Phillips N, Katul GG (2006) The porous media model for the hydraulic system of a conifer tree: Linking sap flux data to transpiration rate. Ecol Model 191:447–468
Clausnitzer V, Hopmans J (1994) Simultaneous modeling of transient three-dimensional root growth and soil water flow. Plant Soil 164:299–314
Cochard H, Lemoine D, Dreyer E (1999) The effects of acclimation to sunlight on the xylem vulnerability to embolism in Fagus sylvatica L. Plant Cell Environ 22:101–108
Cochard H, Venisse JS, Barigah TS, Brunel N, Herbette S, Guilliot A, Tyree MT, Sakr S (2007) Putative role of aquaporins in variable hydraulic conductance of leaves in response to light. Plant Physiol 143:122–133
Cowan I (1965) Transport of water in the soil-plant-atmosphere system. J Appl Ecol 2:221–239
Cruiziat P, Cochard H, Améglio T (2002) Hydraulic architecture of trees: main concepts and results. Ann For Sci 59:723–752
Doussan C, Vercambre G, Pagès L (1998) Modelling of the hydraulic architecture of root systems: an integrated approach to water absorption–distribution of axial and radial conductances in maize. Ann Bot 81:225–232
Doussan C, Pagès L, Pierret A (2003) Soil exploration and resource acquisition by plant roots: an architectural and modelling point of view. Agronomie 23:419–431
Doussan C, Pierret A, Garrigues E, Pagès L (2006) Water uptake by plant roots: II—modelling of water transfer in the soil root-system with explicit account of flow within the root system—comparison with experiments. Plant Soil 283:99–117
Dunbabin VM, Diggle AJ, Rengel Z (2002) Simulation of field data by a basic three-dimensional model of interactive root growth. Plant Soil 239:39–54
Ehlert C, Maurel C, Tardieu F, Simonneau T (2009) Aquaporin-mediated reduction in maize root hydraulic conductivity impacts cell turgor and leaf elongation even without changing transpiration. Plant Physiol 150:1093–1104
Engeln-Müllges G, Uhlig F (1996) Numerical algorithms with C. Springer, Berlin
Espeleta J, West J, Donovan L (2004) Species-specific patterns of hydraulic lift in co-occurring adult trees and grasses in a sandhill community. Oecologia 138:341–349
Feddes R, Kowalik P, Zaradny H (1978) Simulation of field water use and crop yield. Simulation monographs, Pudoc, Wageningen
Früh T, Kurth W (1999) The hydraulic system of trees: theoretical framework and numerical simulation. J Theor Biol 201:251–270
Gardner W (1960) Dynamic aspects of water availability to plants. Soil Sci 89:63–73
Gayler S, Klier C, Mueller C, Weis W, Winkler J, Priesack E (2009) Analysing the role of soil properties, initial biomass and ozone on observed plant growth variability in a lysimeter study. Plant Soil 323:125–141
Hinsinger P, Gobran GR, Gregory PJ, Wenzel WW (2005) Rhizosphere geometry and heterogeneity arising from root-mediated physical and chemical processes. New Phytol 168:293–303
Hoffmann-Riem H, van Genuchten M, Flühler H (1999) General model for the hydraulic conductivity of unsaturated soils. In: van Genuchten M, Leij F, Wu L (eds) Proceedings of the International Workshop on Characterization and Measurement of the Hydraulic Properties of Unsaturated Porous Media, vol 1. University of California, Riverside, pp 31–42
Hölttä T, Vesala T, Nikinmaa E, Perämäki M, Siivola E, Mencuccini M (2005) Field measurements of ultrasonic acoustic emissions and stem diameter variations. new insight into the relationship between xylem tensions and embolism. Tree Physiol 25:237–243
Huang K, Mohanty B, van Genuchten M (1996) A new convergence criterion for the modified picard iteration method to solve the variably-saturated flow equation. J Hydrol 178:69–91
Javaux M, Schröder T, Vanderborght J, Vereecken H (2008) Use of a three-dimensional detailed modeling approach for predicting root water uptake. Vadose Zone J 7:1079–1088
Kroes JG, van Dam J (2003) Reference manual swap 3.03. Tech. rep. Alterra-report 773, Alterra, Green World Research
Kumagai T (2001) Modeling water transportation and storage in sapwood—model development and validation. Agric For Meteorol 109:105–115
Lehmann F, Ackerer P (1998) Comparison of iterative methods for improved solutions of the fluid flow equation in partially saturated porous media. Transp Porous Media 31:275–292
Lemoine D, Cochard H, Granier A (2002) Within crown variation in hydraulic architecture in beech (Fagus sylvatica L.): evidence for a stomatal control of xylem embolism. Ann For Sci 59:19–27
Loos C, Gayler S, Priesack E (2007) Assessment of water balance simulations for large-scale weighing lysimeters. J Hydrol 335:259–270
Magnani F, Borghetti M (1995) Interpretation of seasonal changes of xylem embolism and plant hydraulic resistance in Fagus sylvatica. Plant Cell Environ 18:689–696
Mayer D, Butler D (1993) Statistical validation. Ecol Model 68:21–32
Mualem Y (1976) A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resour Res 12:513–522
Münch E (1928) Versuche über den Saftkreislauf. Dtsch Bot Ges 45:340–356
Nadezhdina N (2010) Integration of water transport pathways in a maple tree: responses of sap flow to branch severing. Ann For Sci 67:107
Nimah M, Hanks R (1973) Model for estimation of soil water, plant, and atmospheric interrelations: I. description and sensitivity. Soil Sci Soc Am Proc 37:522–527
North G (2004) A long drink of water: how xylem changes with depth. New Phytol 163:449–451
Oertli JJ (1993) Effect of cavitation on the status of water in plants. In: Borghetti M, Grace J, Raschi A (eds) Water transport in plants under climatic stress. Cambridge University Press, Cambridge, pp 27–40
Passioura JB (1980) The meaning of matric potential. J Exp Bot 31:1161–1169
Perämäki M, Nikinmaa E, Sevanto S, Ilvesniemi H, Siivola E, Hari P, Vesala T (2001) Tree stem diameter variations and transpiration in scots pine: an analysis using a dynamic sap flow model. Tree Physiol 21:889–897
Priesack E (2006) Expert-N dokumentation der modell-bibliothek. FAM Bericht 60, Hieronymus, München
Priesack E, Gayler S (2009) Agricultural crop models: concepts of resource acquisition and assimilate partitioning. In: Lüttge U, Beyschlag W, Büdel B, Francis D (eds) Progress in botany, vol 70. Springer, Berlin, pp 192–222
Rewald B (2008) Impact of climate-change induced drought on tree root hydraulic properties and competition belowground. PhD thesis, University of Göttingen, Germany
Richards JH, Caldwell MM (1987) Hydraulic lift: substantial nocturnal water transport between soil layers by Artemisia tridentata roots. Oecologia 73:486–489
Richards LA (1931) Capillary conduction of liquids through porous mediums. Physics 1:318–333
Sevanto S, Vesala T, Perämäki M, Nikinmaa E (2002) Time lags for xylem and stem diameter variations in a scots pine tree. Plant Cell Environ 25:1071–1077
Siau J (1984) Transport Processes in Wood. Springer, Berlin
Simunek J, Sejna M, Saito H, Sakai M, van Genuchten MT (2008) The HYDRUS-1D software package for simulating the one-dimensional movement of water, heat, and multiple solutes in variably-saturated media. Tech. rep., Department of Environmental Sciences, University of California, Riverside
Siqueira M, Katul G, Porporato A (2008) Onset of water stress, hysteresis in plant conductance, and hydraulic lift: scaling soil water dynamics from millimeters to meters. Water Resour Res 44:W01432. doi:10.1029/2007WR006094
Somma F, Hopmans J, Clausnitzer V (1998) Transient three-dimensional modeling of soil water and solute transport with simultaneous root growth, root water and nutrient uptake. Plant Soil 202:281–293
Sperry JS, Stiller V, Hacke UG (2003) Xylem hydraulics and the soil-plant-atmosphere continuum: Opportunities and unresolved issues. Agron J 95:1362–1370
Steppe K, Lemeur R (2007) Effects of ring-porous and diffuse-porous stem wood anatomy on the hydraulic parameters used in a water flow and storage model. Tree Physiol 27:43–52
Steudle E (2000) Water uptake by plant roots: an integration of views. Plant Soil 226:45–56
Tseng PH, Sciortino A, van Genuchten MT (1995) A partitioned solution procedure for simulating water flow in a variably saturated dual-porosity medium. Adv Water Resour 18:335–343
Tyree MT, Zimmermann MH (2002) Xylem structure and the ascent of sap, 2nd edn. Springer Series in Wood Sciences. Springer, Berlin
van Bavel CHM, Baker JM (1985) Water transfer by plant roots from wet to dry soil. Naturwissenschaften 72:606–607
van Genuchten M (1980) A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci Soc Am J 44:892–898
van Genuchten M (1982) A comparison of numerical solutions of the one-dimensional unsaturated-saturated flow and mass transport equations. Adv Water Resour 5:47–55
Vrugt JA, van Wijk MT, Hopmans JW, Simunek J (2001) One-, two-, and three-dimensional root water uptake functions for transient modeling. Water Resour Res 37:2457–2470
Willmott C (1982) Some comments on the evaluation of model performance. Bull Am Meteorol Soc 64:1309–1313
Winkler J, Lang H, Graf W, Reth S, Munch JC (2009) Experimental setup of field lysimeters for studying effects of elevated ozone and below-ground pathogen infection on a plant-soil-system of juvenile beech (Fagus sylvatica L.). Plant Soil 323:7–19
Acknowledgements
We are grateful to the Deutsche Forschungsgemeinschaft which funded this study within the frame of Forschergruppe 788 ‘Competitive mechanisms of water and nitrogen partitioning in beech-dominated deciduous forests’. We also want to thank an anonymous reviewer whose comments helped to considerably improve the manuscript and we thank Sebastian Bittner for his help during the revision of the manuscript.
Open Access
This article is distributed under the terms of the Creative Commons Attribution Noncommercial License which permits any noncommercial use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.
Author information
Authors and Affiliations
Corresponding author
Additional information
Responsible Editor: Alexia Stokes.
Appendix A
Appendix A
The finite element method is applied according to the Galerkin method to discretise the 1D Richards equation (van Genuchten 1982; Priesack 2006). For this purpose the 3D graph, which represents the tree is divided into finitely many 1D intervals [z i ,z i + 1] (1 ≤ i ≤ n), the finite elements i (cf. Fig. 1b). For each knot z i , (1 ≤ i ≤ n + 1) the piecewise linear basis function ϕ i is introduced, which is equal to 1 at knot z i and zero for all other knots. By this, every continuous, piecewise linear function f(z) on the graph can be represented as a superposition of basis functions ϕ i , where the coefficients of ϕ i are given by the values f(z i ) of the function f at knot z i . The piecewise linear approximation \(\widetilde{\psi(}t,z)\) of the unknown ψ = ψ(t, z) of the Richards equation therefore can be expressed by the following function:
where \(\widetilde{\psi_i}(t):=\widetilde{\psi(}t,z_i)\) is defined for all 1 ≤ i ≤ n + 1.
The energy method from the calculus of variations leads to the following condition of orthogonality for each of the n + 1 basis functions φ i , i.e. to the requirement that the following integrals vanish on the entire solution domain Ω represented by the graph:
Integration by parts and the insertion of \(\widetilde{\psi }=\) \(\sum_{i=1}^{n+1} \widetilde{\psi_i}\ \phi_i\) leads to
for all 1 ≤ i ≤ n + 1, where \(\partial \Omega\) denotes the boundary of the domain Ω and q b denotes the water flux across the boundary. Since we only consider binary, ternary up to m-ary branching, the basis functions φ i are in each case at most on m + 1 finite elements different from zero and it is advantageous to accomplish the integration element wise per element Ω e = [z e ,z e + 1], 1 ≤ e ≤ n:
where the summation ∑ e needs to be carried out only for those elements Ω e , that include the knot z i . To further evaluate the integrals (Eq. 29) for the term with time derivative the following approximation, also denoted as ‘mass lumping’, is taken as definition for all 1 ≤ i ≤ n + 1 (van Genuchten 1982):
Additionally it is assumed, that the hydraulic conductivity K and also the sink term S x are continuous, piecewise linear functions on Ω, i.e. can be represented similar to the function \(\widetilde{\psi }\) by basis functions ϕ i (1 ≤ i ≤ n + 1), with K i (t) = K(t,z i ), and S i (t) = S(t,z i ) representing the values at the knot i:
Using these representations the integrals can be explicitly calculated and, in case of an equidistant decomposition of Ω into the finite elements Ω e each of equal length Δz, the following generally non-linear equation system (in matrix notation) results:
for the vectors \({\boldsymbol\vartheta}\), \({\boldsymbol\psi}\) and d and for the matrices A and B, where we get in case of a m i -ary branching (1 ≤ m i ≤ m) at node z i for the m i + 1 intervals \([z_{p_i},z_i]\), \([z_i,z_{c_{i,1}}],\ldots,[z_i,z_{c_{i,m_i}}]\) that contain z i :
This nonlinear system of Eq. 33 is then solved in a standard way using a Newton–Picard iteration method (Celia et al. 1990) and elementary Gaussian elimination (Engeln-Müllges and Uhlig 1996).
Rights and permissions
Open Access This is an open access article distributed under the terms of the Creative Commons Attribution Noncommercial License (https://creativecommons.org/licenses/by-nc/2.0), which permits any noncommercial use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.
About this article
Cite this article
Janott, M., Gayler, S., Gessler, A. et al. A one-dimensional model of water flow in soil-plant systems based on plant architecture. Plant Soil 341, 233–256 (2011). https://doi.org/10.1007/s11104-010-0639-0
Received:
Accepted:
Published:
Issue Date:
DOI: https://doi.org/10.1007/s11104-010-0639-0