An optimality-based model of the coupled soil moisture and root dynamics

. The main processes determining soil moisture dynamics are inﬁltration, percolation, evaporation and root water uptake. Modelling soil moisture dynamics therefore requires an interdisciplinary approach that links hydrological, atmospheric and biological processes. Previous approaches treat either root water uptake rates or root distributions and transpiration rates as given, and calculate the soil moisture dynamics based on the theory of ﬂow in unsaturated media. The present study introduces a different approach to linking soil water and vegetation dynamics, based on vegetation optimality. Assuming that plants have evolved mechanisms that minimise costs related to the maintenance of the root system while meeting their demand for water, we develop a model that dynamically adjusts the vertical root distribution in the soil proﬁle to meet this objective. The model was used to compute the soil moisture dynamics, root water uptake and ﬁne root respiration in a tropical savanna over 12 months, and the results were compared with observations at the site and with a model based on a ﬁxed root distribution. The optimality-based model reproduced the main features of the observations such as a shift of roots from the shallow soil in the wet season to the deeper soil in the dry season and substantial root water uptake during the dry season. At the same time, simulated ﬁne root respiration rates never exceeded the upper envelope determined by the observed soil


Introduction
The weakest component of soil-vegetation-atmosphere transfer (SVAT) models is their link with the soil environment (Feddes et al., 2001). Therefore, improvements in our understanding and parameterisation of root water uptake are necessary to increase our confidence in the outputs of global circluation models (GCM) that depend on accurate estimates of vegetation water use. Typically, SVAT models within numerical weather prediction and climate models do a poor job at simulating soil moisture dynamics, which in turn means that the fluxes of water and heat to the atmosphere are also poorly represented. Improving this is therefore also a priority to advance weather forecasting (Giard and Bazile, 2000). Processbased models of soil moisture dynamics usually consider root water uptake by adding a sink term to Richards' equation. This sink term typically depends on the description of an "effective root-density distribution" of varying complexity and other variables (Varado et al., 2006). The choice of vertical distributions of roots in the soil profile can affect the vertical distribution of soil moisture and transpiration rates simulated by models (Feddes et al., 2001). Feddes et al. (2001) therefore suggest that detailed observations of root profiles within different biomes at the global scale may be necessary for improved parameterisation of root water uptake in models and ultimately for improving the predictions of GCMs. A number of such observations have been compiled by Jackson et al. (1997).
However, the necessary pooling of observations within and across sites masks the important spatial and temporal dynamics within a biome (Jackson et al., 1997). It is known that root distributions can be very dynamic, especially in savannas. For example, Chen et al. (2004) have shown for a tropical savanna that the active root distribution can change from being shallow to deeper root dominated within a couple of months. The observed seasonal variation of root abundance within the top meter of soil was one order of magnitude (Chen et al., 2002), presumably due to the dynamic nature of the tree-grass interactions. The dynamics of root systems are also reviewed in Schenk (2005). Prescribing static root profiles, or even empirically-based root growth algorithms, is particularly problematic if the model is intended to be used for the prediction of responses to long-term environmental change. Jackson et al. (2000) pointed out that the belowground changes associated with environmental change can have large impacts on biogeochemical cycles and neglecting this can lead to substantial errors in model outputs.
Models exploring the assumption that vegetation is optimally adapted to its environment are an alternative to models based on prescribed vegetation properties (Raupach, 2005). The optimal water use hypothesis (Cowan and Farquhar, 1977), for example, is useful for predicting the diurnal and day-to-day dynamics of water use if the photosynthetic properties of the vegetation and its monthly water use are known (Schymanski et al., 2008). Another study explored the assumption that natural vegetation self-organises in a way to maximise its net carbon profit (i.e. the difference between carbon acquired by photosynthesis and carbon spent on the maintenance of the organs involved in its uptake) . The model was able to reproduce the observed above-ground leaf area index and photosynthetic properties of a savanna vegetation given the observed water use and meteorological data, without prescribing any sitespecific vegetation properties during the wet season . However, this model did not consider the costs for the uptake and transport of water and therefore it was not able to predict the water use itself or the vegetation cover during the dry season, when the below-ground costs would have been more important than during the wet season . Therefore, an optimality-based model of root water uptake, quantifying the costs and benefits of the root system, is needed to take the vegetation optimality approach one step further and implement it for the whole plant system in coupled ecohydrological models.
The distinct advantage of optimality-based models is that they simulate the adaptation of vegetation to given environmental conditions and, in theory, do not rely on parameter tuning. Therefore, they are particularly well suited for predicting the long-term effects of environmental change, when it can be assumed that vegetation has adapted to the new conditions. Optimality approaches have been explored previously for modelling root water uptake (e.g. Heimann, 1996, 1998; van Wijk and Bouten, 2001;Laio et al. , 2006;Collins and Bras, 2007), but we are only aware of one that modelled a dynamically adapting root distribution (Kleidon and Heimann, 1996). However, the resulting model was deemed impractical for implementation into coupled biogeochemical or ecohydrological models due to its large computational demands.
The aim of this study is to present and test an optimalitybased model relating root water uptake to carbon costs that is simple enough to be implemented into a coupled ecohydrological model allowing for simultaneous optimisation of above-and below-ground vegetation. To test the model, we used the same data set as Schymanski et al. (2007Schymanski et al. ( , 2008, where the optimality approach has been applied previously to model the above-ground vegetation properties of a tropical savanna. This time, the root system was optimised to meet the observed canopy water demand while minimising the root maintenance costs. The observations available for comparison with model outputs were the dynamics of surface soil moisture, evapo-transpiration rates and below-ground respiration. As this study primarily focuses on fine root water uptake, a realistic simulation of transpiration and root respiration rates was taken as the main benchmark for model performance. For comparison with the conventional approach to modelling root water uptake, the same model was also run with a prescribed root profile, which is considered to be typical for a humid tropical savanna (Schenk and Jackson, 2002;Jackson et al., 1997). Model parameters were taken from the literature, or, where not available, reasonable values were prescribed without parameter tuning. Given that the optimalitybased model does not require any input about the root surface area on the site, we consider the model useful if it does not lead to a significantly worse correspondence between the model results and observations than the model based on a prescribed root profile.

Methods
In the following, we describe the soil water balance model used for calculating the soil water fluxes, the vegetation water balance model used for calculating root water uptake and the root optimisation algorithm separately. The soil water balance model calculates the soil moisture distribution within the soil profile, the position of the water table and the spontaneous flow of water in a conceptual lumped catchment, including infiltration during rainfall and ex-filtration into the channel. The vegetation water balance model simulates the water storage within the vegetation (e.g. water storage in tree trunks) and the root suction as a function of this water storage. Due to the explicit consideration of the water storage within the vegetation tissues, the model allows simulation of root water uptake during the night when no transpiration occurs and investigation of phenomena such as hydraulic redistribution by plant roots (Burgess et al., 1998;Meinzer et al., 2001). The purpose of the root optimisation algorithm was to simulate dynamic adaptation of the root system to the water availability in the soil profile and the water demand by the canopy.

Water balance model
To account for the transfer of water between the atmosphere, soil and the river channel, we initially followed the "Representative Elementary Watershed" (REW) approach formulated by Reggiani et al. (2000), but extended it to allow the calculation of the vertical distribution of water within the unsaturated zone and adjusted some of the closure relations to be consistent with our formulation. For simplicity, interactions between elementary watersheds and streamflow routing were neglected, so that all water reaching the channel was assumed to be instantaneous runoff.
In the REW model, an elementary watershed is subdivided into two soil layers, the saturated and unsaturated layer, both of variable thickness. The thickness of each zone and all fluxes are spatially averaged (divided by the catchment area), so that the model can be summarised with just a few variables (Reggiani et al., 2000): the average bedrock elevation from reference datum (z s , m), the average channel elevation from reference datum (z r , m), the average depth of the pedosphere (Z, m), the average thickness of the saturated zone (y s , m), the average thickness of the unsaturated zone (y u , m), the average saturation degree in the unsaturated zone (s u ), the unsaturated surface area fraction (ω u ) and the saturated surface area fraction contributing to seepage face and overland flow (ω o ). See Fig. 1 for a diagram of a simplified REW. For simplicity, we set the datum to coincide with the average bedrock elevation, so that z s =0 m. Any fluxes into and out of the saturated zone (Q u , Q sf or E T s , m 3 m −2 s −1 = m s −1 ) lead to changes in the thickness of both the saturated and the unsaturated zones, and also to changes in the unsaturated and saturated surface area fractions if y s > z r . The relationships between y s , y u , ω u and ω o depend on the geometry of the catchment and are given for linear hillslopes in this study, where ω u and y u are both calculated as a function of y s (for , while the variables on the right hand side denote water fluxes (precipitation (Q rain ), infiltration (Q inf ), infiltration excess runoff (Q iex ), soil evaporation from the saturated zone and the unsaturated zone (E ss and E su respectively), flow between saturated and unsaturated layer (Q u ) and outflow across the seepage face (Q sf ).
a derivation, see Appendix A.3.2.1 in Schymanski, 2007): and The saturated surface area fraction (ω o ) is the complement of ω u . The above relations can be used to calculate ω o , ω u and y u for any given values of Z, z r and y s in a linear hillslope.

Vertical subdivision of the unsaturated zone
In the original REW model, the unsaturated zone was treated as a lumped volume, but for the purpose of this study, the unsaturated zone was subdivided into several layers and soil moisture was calculated for each layer separately. Starting from the top of the unsaturated zone, we divided the unsaturated zone into soil layers of prescribed thickness δ yumin (in this study 0.1 m), until we reached the groundwater table. The soil layer adjacent to the water table was the only layer with differing thickness, which we set to a value between δ yumin and 2 × δ yumin to reflect the total thickness of the unsaturated zone (y u ). The number of soil layers in the unsaturated zone (n layers ) was equivalent to the lower integer value of y u /δ yumin .

Soil water fluxes
Water fluxes between different soil layers were calculated using a discretisation of the Buckingham-Darcy Equation (Radcliffe and Rasmussen, 2002), which is the 1-D equivalent of Richards' equation for steady flow. All fluxes were  Schymanski (2007) for their derivations). The flux between layer i and layer i+1 was calculated as: where h (m) is the "matric suction head" and K unsat (m s −1 ) the unsaturated hydraulic conductivity. The subscript i denotes the ith layer, Q i (m s −1 ) denotes the flux across the bottom boundary of layer i and δ yu,i (m) denotes the thickness of layer i. In the present work, h with units of hydraulic pressure head (m), is defined as positive and increases with decreasing soil saturation. Q i is defined as positive if water flows upwards and negative if it flows downwards.  In the saturated zone, the hydraulic conductivity is K sat (taken as 1.23×10 −5 m s −1 for sandy loam) and the matric suction head (h) is assumed to be 0 m, so that the flux across the boundary between the unsaturated and the saturated zone was written as: The above equations require the calculation of matric suction head (h i ) and unsaturated hydraulic conductivity (K unsat,i ) in each soil layer at each time step. These were obtained as a function of the soil saturation degree (s u,i ), from the widely used water retention model formulated by van Genuchten (1980): The parameters α vG (m −1 ), n vG and m vG have to be fitted to empirical soil water retention curves, and the saturated hydraulic conductivity (K sat m s −1 ) is also an empirical constant, specific for a given soil type. Standard values for different soil types can be found in the literature, where n vG and m vG are usually assumed to follow the relation (van Genuchten, 1980): The soil saturation degree (s u,i ) itself is a function of the volumetric water content (θ i ) and the empirical soil properties θ r and θ s (van Genuchten, 1980): The parameter values of the saturated hydraulic conductivity (K sat m s −1 ), α vg , n vG and θ r were taken from the software package Hydrus 1-D (Simunek et al., 2005) as typical values given for sandy loam (K sat =1.23×10 −5 m s −1 , α vg =7.5 m −1 , n vG =1.89, θ r =0.065 m 3 m −3 ), while θ s was adapted to reflect the highest measured soil water content in the data set (θ s =0.31 m 3 m −3 ). The parameterisation was the same for all soil layers.

Infiltration and runoff
Interception of rainfall by the vegetation was neglected in this study, so that all rainfall was assumed to hit the ground. Interception effects are expected to be relatively small at the study site, mainly because rainfall events tend to be large relative to the size of canopy storage. For application of the model to sites where interception is expected to play a larger role, it could be included by replacing the rainfall rate (Q rain , m s −1 ) in Eq. (9) by the throughfall rate after interception. In the present study, consideration of interception would have led to a slight delay of infiltration, rather than a substantial reduction of the simulated infiltration.
Infiltration was assumed to only occur into the unsaturated zone. The infiltration capacity was expressed by imagining an infinitely thin layer of water above the top soil layer and expressing the infiltration capacity as Q i for i=0, where K unsat,i is replaced by K sat , h i by 0 m, and δ yu,i by 0 m in Eq. (3). The rate of infiltration (Q inf , m s −1 ) was then formulated as the lesser of infiltration capacity and rainfall intensity (Q rain , m s −1 ): Rainfall exceeding Q inf was assumed to contribute to immediate runoff (Q out , m s −1 ). In the presence of a seepage face (i.e. when y s >z r ), flow across the seepage face also contributed to runoff. We calculated the seepage face flow (Q sf , m s −1 ) following Reggiani et al. (2000): The parameter γ 0 is the average slope angle of the seepage face (set to 0.033 radians (=1.9 degrees) in this study), while s is a typical horizontal length scale, which depends on the length of the hillslope. In the absence of a more rigorous treatment of this parameter, we followed the approach by Reggiani et al. (2000) and set the parameter value to s =10 m. It is obvious from Eq. (10) that a larger value of s would have the same effect as a smaller value of γ 0 , i.e. a reduction of Q sf .
The total runoff summed to: 2.1.4 Evaporative fluxes Figure 1 included two other fluxes: the soil evaporation from the unsaturated zone (E su , m s −1 ) and from the saturated zone (E ss , m s −1 ). For the layered model ( Fig. 2), we have to distinguish between soil evaporation, which occurs at the soil-air interface only, and transpiration by vegetation, which is linked to root water uptake (Q r,i ) from all layers within the rooting zone. Root water uptake is described in Sect. 2.2, so only equations for soil evaporation are presented here. Soil evaporation (E s , m s −1 ) was modelled using a fluxgradient approach: where W s denotes the mole fraction of water in the laminar air layer immediately above the soil, and W a the mole fraction of water in the atmosphere, while G soil (mol m −2 s −1 ) is the conductivity of the soil to water vapour fluxes. The molar weight of water (M w , 0.018 kg mol −1 ) and the density of water (ρ w , set to 1000 kg m −3 , irrespective of the temperature) were used to convert from molar units (mol m −2 s −1 ) to volumetric units of liquid water (m 3 m −2 s −1 = m s −1 ). The formulation is equivalent to common soil evaporation models (Lee and Pielke, 1992) if G soil is set to the product of the moisture transfer coefficient and wind speed. W s was calculated as the vapour pressure in the laminar layer immediately above the soil (p vs ) divided by air pressure. p vs was modelled as a function of the atmospheric vapour pressure, the saturation vapour pressure at the measured soil temperature, modelled volumetric soil moisture in the top soil layer and volumetric soil moisture at field capacity at the site (Lee and Pielke, 1992).
The value for soil moisture at field capacity for the site was set to 0.156, equivalent to the soil moisture at a matric suction head (h) of 1.0 m (Kelley, 2002). In the absence of data about near-soil wind speeds, the parameter G soil in Eq. (12) was set to the constant value of 0.03 mol m −2 s −1 , which led to reasonable soil evaporation rates on the site Schymanski, 2007). E s was calculated for the unsaturated and the saturated zones separately, and multiplied by their respective surface area fractions to obtain E su and E ss .

Conservation of mass
Changes in the state variables s u,i , δ yu,i , ω o , ω u , y s and y u due to water fluxes must satisfy conservation of mass. As implied above, we ignored density variations due to changes in temperature and expressed the mass of water per square meter of area by the volume of liquid water per square meter of area, so that the units of mass of water per unit area were given in m 3 m −2 =m, which is consistent with water flux units of m 3 m −2 s −1 =m s −1 .
The derivations of the following equations can be found in Appendix A.3.2.4 in Schymanski (2007). Downwards flux of water into the saturated zone results in its expansion into the unsaturated zone, while an upward flux results in its contraction and an increased unsaturated volume. The change in the thickness of the saturated zone (y s ) was thus expressed as a function of the fluxes in and out of the saturated zone and the saturation of the bottom layer of unsaturated zone (s u,nlayers ): where E ss (m s −1 ) is the soil evaporation rate from the saturated zone, Q nlayers (m s −1 ) is the flux across the bottom boundary of the unsaturated zone, Q sf (m s −1 ) is outflow across the seepage face and ε is the soil porosity (here taken as ε=θ s − θ r , m 3 m −3 ). The division by ε s u,nlayers (t) − 1 in Eq. 13 reflects the fact that the saturated zone can only expand into the air-filled volume of the unsaturated zone. The change in soil moisture for the top soil layer was written as: and for the soil layers between the top and the bottom layer: To calculate the change in the state variables for a finite time step from t 1 to t 2 , the above equations were solved at t 1 and then multiplied by the length of the time step. This gave y s and s u,i at time t 2 for all layers apart from the bottom layer of the unsaturated zone. The saturation degree in the bottom layer at time t 2 (s u,nlayers (t 2 )) was then calculated by difference from the fluxes into and out of the whole soil domain and the change in water storage. The water storage in each soil layer i (w s,i , m) was written as and the sum of fluxes in and out of the soil domain between time t 1 and t 2 was written as: where w c (m) denotes the change in the total water store of the soil domain per unit catchment area. The water store in the bottom soil layer of the unsaturated zone (w s,nlayers , m) at time t 2 was calculated as and the value of s u,nlayers at time t 2 was then obtained from w s,nlayers : The maximum length of each time step was restricted so that no state variable in the model could change by more than 10% in a single time step.

Vegetation water balance and root water uptake
Root water uptake was modelled using an electrical circuit analogy, where radial root resistivity and soil resistivity are in series in each soil layer (Hunt et al., 1991). Water uptake per unit root surface area in a soil layer (J r,i , m s −1 ) was thus written as: where r is root resistivity to water uptake per unit root surface area (taken as 1.02×10 8 s in this study), and s,i (s) is the resistivity to water flow towards the roots in the soil. The driving force for water uptake by roots is the difference between the forces holding the water in the soil (h i , m head) and the forces holding the water in the roots (h r,i , m head). Defining S Ar,i (m 2 m −2 ) as the root surface area per ground area in layer i, we can write the root water uptake rate per ground area in layer i (Q r,i , m s −1 ) as: The resistivity to water flow towards the roots in the soil ( s,i , s) was formulated as a function of the unsaturated hydraulic conductivity (K unsat,i , m s −1 ), root radius (r r , m) and root surface area density in soil layer i (S Adr,i , m 2 m −3 ): Equation (22) has the desired properties that s,i decreases with increasing K unsat,i and decreasing distance between roots (represented by the second term). A derivation of Eq. (22) is given in Appendix A.3.3.1 in Schymanski (2007).

Water storage and tissue balance pressure
In the above, root water uptake was modelled as a function of root and soil properties and the suction head difference between the soil and the inside of the roots. The suction head inside the roots is often considered to be linked to the suction head in leaves, which is caused by adhesive forces and is driven by transpiration. Thus, water transport from the soil to the leaves could occur passively, without the expenditure of energy other than for the maintenance of the plant tissues involved. A model to quantify the forces involved in such a passive process has been developed in Appendix A.3.3.2 in Schymanski (2007) and will only be summarised here. Based on work by Roderick and Canny (2005), who found a correlation between measurements of tissue balance pressure (P b , the pressure that has to be applied in order to force water out of the tissue) and the tissue water content (M q ), the relationship was written as: where P b (bar) is the tissue balance pressure, M qx and M q are the maximum and actual amount of water stored in plant tissues per unit catchment area respectively (kg m −2 ), M d (kg m −2 ) is the total mass of dry matter associated with living tissues per unit catchment area, and c 1 (750 bar) and c 2 (1 bar) are fitted to match the data presented by Roderick and Canny (2005). If the tissue balance pressure is assumed to represent the suction force exerted by the tissue, Eq. (23) implies that the suction force increases as the tissue water content decreases. However, P b can only increase until M q reaches a value of 0.9M qx , because any further decrease in water content is assumed to lead to tissue damage (see Appendix A.3.2.2 in Schymanski, 2007). In order to use the tissue balance pressure in plant organs above ground (P b ) as a driver for passive water uptake by roots in the model, P b was translated into the root suction head h r,i (m) by taking into account the hydrostatic head between roots and trunks: where h h,i (m) is the hydrostatic head difference between the soil surface and the depth of layer i, while c P bm =10.2 m bar −1 is a conversion coefficient to convert from units of P b (bar) to units of h r,i (m). The height of the canopy was not considered in the calculation of h r,i , as the model did not include any information about tree heights. While it is clear from the above that the value of M q is important for calculating water uptake rates by roots in the present model, it was also postulated that M q should not decrease below 90% of its maximum value, M qx . The rate of change in M q was written as a function of root water uptake and transpiration rate: where ρ w (1000 kg m −3 ) is the density of water, Q r,i (m s −1 ) is the water uptake rate by tree roots in soil layer i, i r is the deepest soil layer accessed by roots and E t (m s −1 ) is the transpiration rate. To hold M q above 0.9M qx , we prescribed root-induced stomatal closure whenever E t would otherwise exceed root water uptake at M q =0.9M qx . A large tree water storage capacity (M qx ) can act as a buffer for meeting peak foliage water demands that exceed root water uptake rates during the day. Ideally, water use by trees and grasses should be modelled separately because the grasses that dominate the fluxes during the wet season have a much smaller water storage capacity than trees, which can store water in their sapwood. For the purpose of the present study, however, it was not feasible to distinguish tree water use from grass water use, and hence trees and grasses were treated as one system with a common water storage capacity. The error was not expected to be very large, as the effect of the tissue water storage on root water uptake should be relatively small during the wet season, when soil moisture is high.
The total above-ground volume of sapwood at the study site was estimated to be 0.0032 m 3 per m 2 catchment area and observed mean values of sapwood density ranged between species from 0.81 to 0.94 g cm −3 (Cernusak et al., 2006). In rough terms, this gives an estimated 3.0 kg of sapwood dry matter per m 2 catchment area. To also account for the dry mass in tree leaves with a leaf area index of 0.6 and a specific leaf area of 5.5 m 2 kg −1 (Cernusak et al., 2006), we used 3.1 kg m −2 for M d in the model. Taking a typical value for Eucalyptus leaves (Roderick and Canny, 2005), we set M qx =M d in the present model, which also allows the maximum values of P b in Eq. (23).

Costs and benefits of the root system
We assumed that there would be a relationship between the vegetation's capacity to extract water from the soil and the amount of carbon that has to be invested in the root system. Unfortunately, such a general relationship was difficult to obtain from the literature, as respiration rates, root hydraulic properties and turnover rates are rarely measured on the same plants, while all of these parameters are highly variable, not only between species, but even within the same root over time (Steudle, 2000).
In the absence of a general empirical relationship between root costs and their water uptake capacity, we used measurements on citrus roots, for which observations of both respi-ration rates and hydraulic properties were available in the literature. Note that there is no evidence suggesting similarity between citrus roots and the roots of the savanna vegetation on the study site, but we hypothesise that the relationship between fine root water uptake capacity and fine root respiration would be roughly similar over a wide range of species. Huang and Eissenstat (2000) measured radial conductivity in different citrus species and found hydraulic conductivities of 1 to 3 µm s −1 MPa −1 per m 2 root area in first-order lateral roots (0.34 to 0.44 mm diameter). In second-order lateral roots (0.58 to 0.87 mm diameter), they found values of 0.2 to 0.75 µm s −1 MPa −1 m −2 . Using 1 µm s −1 MPa −1 m −2 as a typical value for radial root conductivity per root area and converting to units of hydraulic head, we set the root resistivity to water uptake ( r ) to 1.02×10 8 s. Bryla et al. (2001) gave values of root respiration for a single citrus fine root on a dry weight (DW) basis in the order of 10 nmol (g DW) −1 s −1 . The average fine root diameter of the measured roots was 0.6 mm. Eissenstat (1991) gave values for the dry mass to volume relationships of different citrus roots between 0.15 and 0.2 g cm −3 . Taking 0.17 g cm −3 as a typical value, 1 m 3 root volume would have a dry weight of 0.17×10 6 g. Consequently, the respiration rate for 1 m 3 of fine roots would be 0.0017 mol s −1 . Assuming cylindrical roots, we obtained root respiration per unit catchment area (R r , mol m −2 s −1 ) as a function of root radius (r r , m) and root surface area per unit ground area (S Ar , m 2 m −2 ): where, following the above, c Rr =0.0017 mol s −1 m −3 and r r =0.3×10 −3 m. The root surface area per unit ground area (S Ar , m 2 m −2 ) was modelled as the sum of the root surface area densities in all soil layers (S Adr,i , m 2 m −3 ) multiplied by the volumes of the respective soil layers per unit ground area:

Root optimisation
The canopy water demand, determined by the observed transpiration rates, has to be met by root water uptake, so that the optimisation problem for the root system is the minimisation of costs while meeting the water demand by the canopy. The optimisation of the root system was performed on a daily scale and involved two steps. The first step was to determine whether the actual root surface area was more or less than adequate to meet the water demand during the past day. This was performed by recording the minimum value of the tissue water store (M q , kg m −2 ) during the past day (M q min , kg m −2 ), which was then used to compute a coefficient of change for the root system (k r ): M q min can range between 0.9M qx and M qx . A value of M q min =0.9M qx suggests that root water uptake did not meet the canopy water demand, causing stomatal closure to prevent further depletion of M q , while M q min =M qx suggests that canopy water demand did not deplete M q at all. The use of the factors 0.95 and 0.05 in Eq. (28) results in k r ranging between 1 for the first case and −1 for the latter case. The second step was to determine the relative effectiveness of roots in different soil layers (k reff,i ) on the past day. This was performed by dividing the daily water uptake per unit root surface area in each soil layer (J rdaily,i , m d −1 ) by the daily water uptake per unit root surface area in the most effective soil layer.
The change in root surface area in each soil layer from day to day was computed as a function of k r and k reff,i : where G r max (m 2 per m 3 ) is the maximum daily root growth rate. The value of G r max determines the maximum fine root turnover rate and thereby the flexibility of the root system to respond to changes in soil moisture or canopy water demand. Based on trial and error, the value of G r max was set to 0.1 m 2 per m 3 soil volume. This allowed enough flexibility to satisfy the canopy water demand during soil drying while avoiding excessive day-to-day fluctuations of S Adr,i in response to the short-term variability of soil moisture. The optimisation procedure simulates a dynamic adaptation of the root system to the canopy water demand and the water supply in the soil, based on the memory of the past. It allows for a seasonal variation in the vertical root distribution in response to the seasonality of the vertical soil moisture distribution.

Prescribed root profile
For comparison with the conventional approach to modelling root water uptake, we ran the same model with a prescribed, static root distribution. The static root distribution used was the typical root distribution for humid tropical savannas given by Schenk and Jackson (2002), with a root area index of 43 m 2 m −2 as given for tropical grasslands/savannas by Jackson et al. (1997).

Study site
The site chosen for the present study is the Howard Springs eddy covariance site, which is located in the Northern Territory of Australia, 35 km South-East of Darwin, in the Howard River catchment (12 • 29 ′ 39.30 ′′ S, 131 • 9 ′ 8.58 ′′ E). The climate is sub-humid on an annual basis (1750 mm mean annual rainfall, 2300 mm mean annual class A pan evaporation), but with a very strong monsoonal seasonality. Approximately 95% of the 1750 mm mean annual rainfall is restricted to the wet season (December to March, inclusive), while the dry season (May to September) is characterised by virtually no rainfall and high atmospheric water demand . Air temperatures range between roughly 25 and 35 • C in the wet season and between 15 and 30 • C in the dry season.
The terrain at the study site is very flat, with slopes <1 • (Beringer et al., 2003. The surface of the lowland plains, where the study site is situated, is a late Tertiary depositional surface, with a sediment mantle that seldom reaches more than 30 to 40 m in depth. On the study site itself, the soil profile has been described as a red kandosol, with sandy loams and sandy clay loams in horizons A and B respectively and weathered laterite in the C horizon, below about 1.2 m (Kelley, 2002).
The site is situated between the Howard River (4.5 km to the West, around 20 m AHD (Australian Height Datum)), and a smaller river channel, (0.5 km to the East, around 30 m AHD). The terrain reaches a maximum elevation of roughly 40 m AHD between these two channels. In terms of the catchment conceptualisation in Fig. 1, we interpreted the catchment as having an average depth of the pedosphere (Z) of 15 m, and an average channel elevation (z r ) of 10 m from the reference datum, which was set to coincide with the average bedrock elevation, so that z s =0 m).
The vegetation has been classified as a Eucalypt open forest (Specht, 1981), with a mean canopy height of 15 m, where the overstorey has an estimated cover of 30-50% Schymanski et al., 2007) and is dominated by the evergreen Eucalyptus miniata and Eucalyptus tetrodonta. The dominant tree species contribute to 60-70% of the total basal area (i.e. the ground area covered by tree trunks) of this forest and are accompanied by some brevi-, semi-and fully deciduous tree species . The overstorey leaf area index (LAI) varies little seasonally, between roughly 0.6 during the dry season and 0.95 during the wet season . The understorey on the site is highly dynamic. During the dry season it is composed of small individuals of the tree species, some fully or partly deciduous shrubs and some perennial grasses with a total LAI of around 0.2, while during the wet season it is dominated by up to 2 m tall annual C 4 grasses of the genus Sarga sp. and reaches LAI values of 1.5 .
The root system of the vegetation on the site is mainly limited to the top 4-5 m of soil (Kelley, 2002), with single roots observed at depths of up to 9 m, but not in significant quantities (O'Grady, unpubl. data). For the present model, we assumed a constant rooting depth of 5 m after Kelley (2002).

Measurements
The data used for this study were the same as described in Schymanski et al. (2007Schymanski et al. ( , 2008. In summary, transpiration rates were obtained by subtracting estimated soil evaporation rates from the observed latent heat flux using the eddy covariance technique. Soil evaporation was hereby estimated using the same model as described in Eq. (12), but utilising measured soil moisture. The soil moisture was measured using time domain reflectometry (TDR) probes (CS615 probes, Campbell Scientific, Logan, UT, USA) at 10 cm depth, and soil temperature was obtained from an averaging soil thermocouple with sensors at 2 and 6 cm depth.

Dynamically optimised root profile
As described in Sect. 2.3, the optimality-based model increased or decreased the root surface area at the end of each day based on the daily minimum in the tissue water store (M q ) recorded for that day. If the water store has been drawn down too much, root surface area was increased, otherwise it was decreased. In each individual soil layer, the changes in root surface area were also dependent on this layer's effectiveness in root water uptake relative to all other soil layers.
This optimisation led to the effect that the simulated root distribution seemed to reflect the general direction of water flow. Starting with an equilibrium soil moisture distribution (i.e. no fluxes within the soil domain) and a uniform root distribution in the 5 m thick root zone, the modelled root system self-optimised to assume a distribution of root surface area (S Ar,i ) that was skewed towards the deeper soil after 40 days (Fig. 3). Over the 40 day period, S Ar increased from the initial 0.1 m 2 m −3 to 0.7 m 2 m −3 in the lowest root layer. During the dry season, when the root zone can only be recharged from the moister soil layers below, the roots were concentrated at the bottom of the profile. In the wet season, when the soil was wetted from the top, the simulated root distribution shifted towards the top soil layers, where also the majority of water uptake took place (Figs. 4 and 5). A high root surface area also remained in the bottom layer of the root zone, where water supply is lower, but steady. When the distribution of soil moisture was very heterogeneous in the soil profile, the model predicted temporal release of water by roots in the driest soil layers (not shown), but this effect was not evident at the daily time scale. Over 24 h, water uptake by roots was generally greater than water release in all soil layers, with some exceptions that did not exceed a net water release of 0.06 mm over 24 h into a layer 0.1 m thick. This demonstrates that the dynamic root optimisation, which continuously shifted roots from layers with little daily water uptake to layers with larger daily water uptake (see Sect. 2.3) led to an avoidance of water loss by the root system. Figure 7 demonstrates that the prescribed plant water storage capacity (M qx ) of 3.1 kg m −2 can allow for spikes in transpiration rates during the day and continuing root water uptake during parts of the night, which would not be possible in a model without a significant plant water storage capacity. However, the prescribed M qx of 3.1 kg m −2 did not have a large impact on the annual transpiration, which was simulated as 1092 mm, compared with 1089 mm simulated by the model if a negligible plant water storage capacity was prescribed. The slight discrepancy between the simulated annual transpiration and the observed annual transpiration of 1118 mm is caused by occasional limitations of transpiration by root water uptake, as shown in Fig. 7. These occurrences are limited to the period between September and November 2004, when the assumed initial soil moisture in the soil pro-file was depleted and prior to its replenishment by the first wet season rain falls. In the second modelled dry season, simulated root water uptake was never limiting for transpiration (Fig. 8a).
Modelled root respiration rates per unit ground area resulting from the dynamically optimised root surface area varied between roughly 0.1 and 0.7 µmol CO 2 m −2 s −1 (data not shown). This seems reasonable as the minimum observed soil respiration per unit ground area on the site was around 1.5 µmol CO 2 m −2 s −1 (Chen et al., 2002;Schymanski, 2007), so that the modelled root respiration never exceeded the bulk soil respiration estimated from observations. Note that soil respiration includes the CO 2 release by the soil due to the decomposition of organic matter and hence fine root respiration can only be smaller than soil respiration.
Modelled and observed surface soil moisture were very similar both in magnitude as well as dynamics, except for a moderate under-estimation of surface soil moisture by the model at the beginning of the model run (Fig. 6). The mean absolute error (MAE) of the difference between model and observations was 0.0156 m 3 m −3 , which is equivalent to about 10% of the wet season values. Major contributors to this error were an additional spike at the end of the modelled time series, which was not observed in the measurements and a generally faster onset of the simulated surface wetting compared with that observed during rainfalls, which was probably due to the neglected interception.

Prescribed, fixed root profile
The typical root distribution for a savanna found in the literature had an exponential decline of root surface area with depth (Schenk and Jackson, 2002) and a root area index of 43 (Jackson et al., 1997). Prescription of such a fixed root profile led to frequent water release by roots in the model (e.g. Fig. 9a), reaching values of up to 1.3 mm d −1 in a single soil layer. A comparison of the observed surface soil moisture time series with the one modelled using the fixed root distribution is given in Fig. 10. The modelled surface soil moisture generally decreased faster than that observed after rainfalls, leading to an under-estimation of surface soil moisture from the start of the wet season on and consequently an increased model error (MAE=0.02 m 3 m −3 ).
The modelled root water uptake failed to meet the ob-served water demand by the canopy on many occasions (Fig. 8b), leading to a reduced modelled annual transpiration of 1055 mm compared with the observed 1118 mm. The prescribed root area index of 43 m 2 m −2 led to a constant root respiration rate per unit ground area of 10.9 µmol m −2 s −1 , which exceeds the estimated dry season soil respiration rates of 1.5 µmol m −2 s −1 by an order of magnitude. Assuming a root area index of 10 m 2 m −2 , in comparison, led to an improvement in the match between modelled and observed surface soil moisture (MAE=0.0178 m 3 m −3 ), but an even higher reduction in modelled annual transpiration (953 mm, Fig. 8c).

Discussion
The aim of this study was to present and test an optimalitybased model relating root water uptake to carbon costs that is simple enough to be implemented into a coupled ecohydrological model. The presented root model is easy to couple with an aboveground vegetation model, because the root optimisation is performed dynamically from day to day, without the need for iterations. The optimality-based model also led to a reasonable reproduction of observed surface soil moisture dynamics and vegetation water use without the need for prescribing a root distribution function or any parameter tuning. In addition, the root respiration rates resulting from the dynamically optimised root distribution never exceeded the soil respiration rates estimated from observations. In contrast, running the same model with a fixed root distribution, obtained from the literature as a typical one for tropical savannas, led to an under-estimation of vegetation water use and an over-estimation of dry season root respiration by an order of magnitude. We found that arbitrary tuning of the root area index could reduce the over-estimation of root respiration, but would lead to a greater under-estimation of vegetation water use.
Note that the observed transpiration rates were used to prescribe the canopy water demand, and the simulated water use could only deviate from the observed when the root system was less than adequate to satisfy this demand. Hence, simulated water use could only be less than the observed, but should not deviate too much if the root water uptake model was realistic. At the same time, the observed soil respiration should always exceed the simulated root respiration, as soil respiration is the sum of root respiration and microbial respiration due to the decomposition of soil carbon. Hence, both the observed water use and the observed soil respiration represent two constraints for a realistic root water uptake model.
Based on these results, the model with a fixed root distribution has to be rejected in favour of the optimality-based model for modelling root water uptake and its associated carbon costs in the given savanna. It may appear trivial that a more complex model (dynamic root distribution) performs better than a simpler model (fixed root distribution) and the question arises whether the better performance of the optimality-based model is due to the optimisation or solely due to its larger complexity compared with the empirical model. The fact is that fine root dynamics can rarely be parameterised empirically, as appropriate observations are not available for most sites. The optimality-based model hence constitutes a major advance in the modelling of dynamic root systems, as it allows parameterising unobserved root system properties without the need for parameter tuning.
One of the reviewers of this paper suggested that the fine root dynamics simulated by the optimality-based model seem unrealistically high. However, recall that the study site has already been shown to have very dynamic roots (Chen et al., 2002(Chen et al., , 2004. The maximum possible fine root increment in the model is determined by the parameter G rmax in Eq. 30, which was set to 0.1 m 2 m −3 d −1 . In comparison, Eissenstat (1991) reported fine root growth rates of 0.07-0.15 cm cm −3 wk −1 for citrus fine roots with a specific root length of 2 cm mg −1 , as in our manuscript. This would translate to growth rates of 0.19-0.4 m 2 m −3 d −1 , which is more than twice the maximum value prescribed in the model. Eissenstat's values are given for roots growing into fresh soil and represent the behaviour after disturbance, not necessarily the response to natural variability in soil moisture. However, the values prove that roots have the capability to grow very fast, even faster than permitted in our model. The more important question is "under what conditions would plants have an advantage by rapidly adjusting their root distribution in the soil profile?". Our results suggest that a savanna with a highly seasonal soil moisture in the top soil represents such conditions. A recently published paper presents observations of fine root growth under natural conditions in the same vegetation type as the one at our study site (Janos et al., 2008). This paper shows an increase in fine root length in the top 1 m of soil from 10 to 120 m m −2 in 3 weeks, which is equivalent to an average increase in root surface area by 0.01 m 2 m −3 d −1 if the dominant fine root radius is assumed to be 0.3 mm. In comparison, the maximum fine root growth rate in our model occurred at the beginning of the simulation, when the root surface area increased from 0.1 to 0.7 m 2 m −3 in 40 days at the bottom of the root profile. This is equivalent to an average increase in root surface area by 0.015 m 2 m −3 d −1 . We conclude that Eissenstat's observations confirm that the maximum growth rate theoretically permitted in our model is realistic, while the observations by Janos et al. confirm that the maximum growth rate simulated by our model at the study site is realistic, too.

Effect of plant water storage
Besides avoiding the need for prescribing the abundance and distribution of roots in the soil profile, the model presented here has another innovative feature compared with conventional root water uptake models. Existing models assume that when stomata close, there is no movement of water through the plant system resulting in the abrupt shut-down of root water uptake. Recently, Amenu and Kumar (2007) formulated a root water uptake model that enables the simulation of simultaneous efflux of water from the plant to the soil in layers of dry soil and root water uptake in layers of wet soil at night, when stomata are assumed to be closed. However, observations of prolonged sap flow after the shut-down of canopy transpiration (e.g. Silberstein et al., 2001;Unsworth et al., 2004) and the finding that tree water storage can be important for tree water use suggest that there can be significant net water uptake at night when the plant water store is being filled up (e.g. Goldstein et al., 1998;Meinzer et al., 2003;. The present study explicitly accounts for the water storage capacity of living plant tissues (Roderick and Canny, 2005) and formulates the suction force exerted by the roots as a function of the amount of water stored in the plants. This results in net water uptake during part of the night (Fig. 7), and opens the way for investigating the costs and benefits of water storage tissues. Although the effect of the water store on the site investigated in the current study was relatively small, it is likely to become more important in catchments dominated by large trees (e.g. Meinzer et al., 2004b;Phillips et al., 2003;Unsworth et al., 2004;Waring and Running, 1978). To our knowledge, the presented model is the first one to allow consideration of such effects in hydrology.

Hydraulic redistribution
The present model also allows simulation of processes such as the uptake of water by roots in wet soil and simultaneous release of water by roots in dry soil layers ("hydraulic redistribution"). Hydraulic redistribution (HR) has been widely observed and could be seen as a passive process, which depends on the soil suction head and the root distribution within the soil column (Schulze et al., 1998;Burgess et al., 1998Burgess et al., , 2001Meinzer et al., 2004a;Hultine et al., 2004;Espeleta et al., 2004;Brooks et al., 2002). Species not showing HR are often those that shed roots in dry soil patches and thus avoid the loss of water from roots (Espeleta et al., 2004). HR could improve the uptake of nutrients from the surface soil, which would otherwise be inhibited by dryness (Burgess et al., 2001). On the other hand, HR could be an undesired "leak in the system", a view which is supported by the observation that root resistance to water release seems to be generally higher than root resistance to water uptake (Hunt et al., 1991). Although root resistance was assumed to be the same for water movement in both directions in the present study, the predicted water release by roots, when it occurred, was very small and hardly apparent at the daily scale if the root distribution was dynamically optimised. This is in line with field observations elsewhere, which showed that tree roots take up more water from shallow soil than they exude via hydraulic lift (Ludwig et al., 2004). In contrast, if the model was run with a fixed, prescribed root distribution, substantial water release by roots was predicted even at the daily scale, when soil layers with a high root abundance were relatively dry. Dawson (1993) documented the uptake of water by shallow-rooting plant species, that has previously been released near the surface by deep-rooted species. As the present study does not distinguish between tree and grass water use on the site, we could not investigate such behaviour, but the observation that the water released by shallow roots during the night is subsequently taken up during the day does not preclude the possibility that the roots releasing the water belong to different species than the roots taking it up later.
One could conclude that the model run based on root optimisation represents species whose roots avoid dry soils and hence do not show HR, while the model run with a fixed root distribution represents species whose fine roots survive even in dry soil and express pronounced HR. The latter corresponds to previous approaches to modelling HR (e.g. Amenu and Kumar, 2007;Personne et al., 2003), while modelling drought-avoiding roots in a hydrological model has not been done before to our knowledge. The better reproduction of observations by the optimisation-based model suggests that drought-avoiding roots dominate in the investigated savanna, but it does not allow conclusions about the generality of such a strategy. In fact, if root decay and construction were associated with an additional cost in the model, it might turn out to be more beneficial to have a less dynamic root distribution, leading to more frequent root water release. In addition, if the uptake of nutrients was considered as another objective function of fine roots (beside water uptake), it could turn out to be beneficial for the plants to increase the release of water into shallow, dry but well aerated soil layers to increase nutrient availability.

Caveats and need for further research
Although the results presented here appear very promising for the use of the optimality-based model, we wish to point out that the roots were only optimised for the uptake of water. The optimal root distributions may be different if nutrient uptake was made part of the objective function for the optimisation. In addition, variability in root resistivity to water uptake was neglected and the parameterisation of the root costs was very simplistic and based on observations in citrus roots only, although it is known that roots can be very versatile. For instance in Lotus japonicus, radial root hydraulic conductivity has been shown to vary 6-fold between 8.0×10 −9 m s −1 MPa −1 and 4.7×10 −8 m s −1 MPa −1 on a diurnal basis, due to its control by aquaporins (Henzler et al., 1999). Aquaporins can be thought of as another degree of freedom available to plants for the regulation of their water uptake. For example, they could open when tissue "suction" would lead to water uptake (during the day) and close to reduce the reverse effect at night. Alternatively, they could open in soil patches with high concentrations of particular nutrients and close in nutrient-poor patches, to use the transpiration stream for the selective uptake of nutrients. Perhaps, they could even discriminate salty water against fresher water in salinity-affected soils. However, despite their functional similarity to "Maxwell's demon" (Maxwell, 1871), which was a theoretical construct designed to highlight ways to violate the second law of thermodynamics, aquaporin regulation requires energy expenditure by plant cells, but quantitative data on these costs are not yet available. More research will be needed to understand all of the plants' degrees of freedom related to water uptake and the associated costs.
Given our limited understanding of the below-ground processes related to nutrient and water uptake, the question remains whether application of optimality assumptions or empirical parameterisation of the root system, based on incomplete observations, leads to less uncertainty in models. In vegetation with a highly dynamic root system, as the one investigated in this study, the optimality assumption clearly led to a better representation of root water uptake, respiration rates and surface soil moisture dynamics than prescribed root distributions. On the other hand, Jackson et al. (1997) found that many studies reported little seasonal variations in root biomass, with the exception of savannas. This is not surprising, as the dynamically changing root distribution should only be beneficial where water availability shifts seasonally between different soil depths. If water availability does not vary much or if there is no accessible water in deeper soil layers during the dry season, shifts in the root distribution would not be expected. Our finding that a static root distribution would not allow adequate water uptake to meet the observed canopy demand in the water-limited savanna is also consistent with the findings by Teuling et al. (2006), who found that land surface schemes with a static root distribution are likely to under-estimate root water uptake in water-limited conditions. Where the root system can be assumed to be reasonably static, fixed parameterisation based on empirical observations can avoid the uncertainty related to the correct param-eterisation of the objective functions and associated costs and benefits inherent in the optimality approach. However, the assumption of a static root system is clearly not reasonable for predictions of responses to long-term environmental change (Norby and Jackson, 2000), where, in our opinion, optimality approaches should be preferred over empirical parameterisations.

Conclusions
Optimality assumptions reduce the need for empirical model parameterisation, which is particularly important for modelling the adaptation of vegetation to environmental change. The present study introduced a model of root water uptake that simulates a dynamically optimising fine root surface area in the soil profile with the objective to meet the canopy water demand while minimising carbon expenditure for fine root maintenance. The simulation results obtained for a tropical savanna are consistent with observations in terms of total water use, surface soil moisture dynamics and soil respiration rates. Given that the results obtained using the dynamically optimising root surface area reproduced available observations even better than results based on an empirically prescribed root distribution, we conclude that the presented model is a useful tool to parameterise the costs and benefits of root water uptake, and allows consideration of belowground adaptation of vegetation to its environment. The model's independence from prescribed root distributions and its low computational demand could make it a powerful tool in conjunction with optimality-based above-ground models to simulate the effects of long-term environmental change on vegetation and the water balance.