Physically-based reduction function to model unsaturated flow associated with plant transpiration

. extracted by the roots due to transpiration, this increases soil suction and, hence, soil shear strength. Transpiration occurs in two different regimes, energy-limited and the water-limited regimes respectively. These two regimes are reflected in the two branches of the transpiration reduction function used to model the hydraulic boundary conditions for vegetated ground. The water-limited branch accounts for the reduced transmissivity of the soil-root system when the degree of saturation and, hence, the hydraulic conductivity declines. The water-limited branch of existing reduction functions (e.g., Feddes function) is defined in purely phenomenological fashion with parameters that have no clear link with the complex interaction between soil hydraulic properties and root architecture. A paradigm shift can be achieved through physically-based reduction functions. These require analytical closed-form solutions of radial water flow at the soil-root interface that, in turn require introducing simplifying assumptions, i.e., steady-state flow and a simplified hydraulic conductivity function. This paper explores the implications of these assumptions by i) benchmarking the water-limited branch of the reduction function derived analytically against the one derived numerically for more realistic hydraulic behaviour and ii) assessing the steady-state assumption.


Introduction
Nature Based Solutions, NBS, are increasingly being used as engineering solutions [1][2][3]. They work with and enhance the positive benefits of natural systems to adapt and promote resilience to climate change. The use of vegetation to improve the stability of natural and engineered slopes is one engineering NBS. One effect of vegetation is to reinforce slopes ydrologicall by generating suction by the removal of soil water via transpiration. To improve upon this stabilising technique, it is key to develop transpiration models that account for the hydraulic characteristics of the soil and plant (belowand above-ground). In this way, modelling can guide the choice of the plant functional traits. Evapotranspiration can occur under two different regimes respectively. Energy limited evapotranspiration, also referred to as potential evapotranspiration, occurs when water is made available by the soil-plant system. Evapotranspiration is driven by the evaporative demand of the atmosphere, i.e., the energy supplied via solar radiation required to convert liquid water into vapour water (latent heat of evaporation), air relative humidity, and wind speed.
When the degree of saturation and, hence, the hydraulic conductivity of the soil declines, the soil-plant * Corresponding author: eve.roberts-self@strath.ac.uk system is not able to accommodate the evaporative demand of the atmosphere and evapotranspirationinduced outward water flux reduces. Evapotranspiration A very common approach to model macroscopic water uptake by vegetation is to consider actual transpiration, AET, as the product of the potential (energylimited) evapotranspiration, PET, times a reduction factor , assumed to be a function of the suction in the root zone, sbulk: Under optimal soil water conditions, root water extraction rate is equal to the maximum evapotranspiration rate, PET ( =1). Under non-optimal conditions, i.e., the soil is either too dry or too wet, transpiration is reduced by means of the factor <1). Feddes et al. [4] assumed that the reduction factor is a function of soil suction in the root zone as presented in Fig. 1. The evapotranspiration is assumed to be equal to zero for suction lower than s1 and above the wilting point s4; the transpiration is maximum ( =1) between s2 and s3, with the latter corresponding to the suction in the soil above which plant uptake starts to be limited. The suction s3 marks the transition from the energy-limited (potential) evapotranspiration to the water-limited evapotranspiration and is the most critical parameter of the Feddes function [5]. The approach proposed by Feddes et al. [4] to model the reduction factor is widely used in geotechnical applications [5][6][7][8][9][10][11]. This approach is convenient because it only requires information about the suction in the root zone without the need to address the complex interaction between the soil, the plant, and the atmosphere.
However, this simplicity is only apparent because the complexity of such an interaction is hidden in the choice of the parameters of the Feddes function. It is indeed the water potential in the leaf and not the suction s3 in the soil that dictates the transition from the energy-limited to the water-limited regime, due to mechanisms of stomatal closure triggered when the water potential in the leaf falls below a threshold [12].
The problem of the choice of the Feddes parameters is reflected in the very wide range of parameters adopted in the literature for s3 [13]. When the Feddes function is used in geotechnical applications, the parameter s3 is generally borrowed from the agricultural literature. This approach may be questionable as the parameters derived for crop species and often loosely compacted organic agricultural soils may significantly differ from non-crop species in often densely-compacted soils that are typically encountered in geotechnical applications [14].
A paradigm shift in modelling of the effect of transpiration is the development of physically-based reduction functions based on the concept of the Soil-Plant-Atmosphere Continuum (SPAC) (Fig. 2). Soil, plant, and atmosphere taken together form a physically integrated, dynamic system, where water flow takes place from regions of higher water potential in the bulk soil to regions of lower water potential in the atmosphere. This flow involves water movement in the bulk soil toward the roots, absorption into the roots, transport through the roots to the xylem and through the xylem to the leaves, evaporation in the intercellular air spaces of the leaves, vapour diffusion through stomatal openings to the boundary air layer in contact with the leaf surface, where the vapour is finally transported to the external atmosphere. Each component of the SPAC system represents an in-series resistance to the water flow and its response is controlled by its hydraulic properties.
The most critical component of the SPAC is represented by the soil-root interface, which involves radial flow towards individual roots [15]. To develop a physically-based closed-form reduction functions and, in particular, the water-limited branch of reduction function, the radial water flow needs to be modelled analytically. In turn, this requires introducing simplifying assumptions, i.e., steady-state flow and simplified hydraulic conductivity functions. This paper explores the implications of these assumptions and assess the robustness of an analytical closed-form equation used to model water flow at the soil root-interface. The analysis of radial water flow towards cylindrical root has been widely addressed in the literature [16,17]. The novelty of this work consists in the development of a closed-form solution that can be successfully used as a basis for an analytical expression of the water-limited branch of a transpiration reduction function. The cylindrical single root model (Fig. 3) is considered here to analyse the soil-root system [18,19]. Only passive pathways (water passes through roots without plant actively dedicating its resources into transporting) is considered here [20,21]. The radial flow equation is written as follows: where r is the radial coordinate in the reference system which has an origin at the centre of the root, k is the hydraulic conductivity, is the volumetric water content, and h is the pore-water pressure head. The latter is given by: ( 2) where uw is the pore-water pressure, s is the matric suction, and w is the unit weight of water. Under steady-state conditions, the water flow equation becomes: Let us consider the case where the hydraulic conductivity is modelled via an exponential function: where ksat is the hydraulic conductivity at zero pressure head and is a soil parameter. With an exponential hydraulic conductivity function, the solution of the steady-state equation can be derived as follows: (5) where hbulk and hroot are the pressure head in the bulk soil and at the root respectively, root is the length of the root segment, and Q is the rate of the flow through the single root. The flow rate through a single root is linked to the evapotranspiration ET as follows: (6) where N is the number of roots per unit horizontal area. It is interesting to observe that the flow rate through a single root is bounded by a maximum value that can be derived analytically by considering the limit for hroot - This limiting flow rate Qmax depends on the pressure head in the bulk soil, hbulk. Under steady-state conditions, this limiting flow rate controls the transition from energylimited to the water-limited regime:

Geometry of flow domain
The geometry of the radial flow domain is characterised by rroot, rbulk, and root. These geometrical parameters can be linked to root architecture parameters including: Under the assumption that all roots are vertical and parallel, arranged to form a regular square mesh, and assuming that they are spaced at twice the outer radius rbulk of cylindrical flow domain, we can write: The root architecture parameters considered in this exercise together with the geometry of the cylindrical water flow domain derived from Eqs. (9) are shown in Table 1.

Soil hydraulic properties
Two soil types were considered for the analyses, a coarsegrained and a fine-grained soil, having air-entry suction sAE equal to 6 kPa and 24 kPa respectively. These were characterised in terms of water retention and hydraulic conductivity functions using the van Genuchten model [22] to enable water flow analysis under steady-state and transient conditions. The two soils were also characterised in terms of hydraulic conductivity using the exponential function given by Eq. (4) to evaluate the analytical solution given by Eq. (5).

Van Genuchten model (VG)
The water retention function of the two soils considered was modelled through the following equation: where s and r are the saturated and residual water content respectively, and a, n, and m, are soil parameters. The water retention parameters considered for the finegrained and coarse-grained soils are shown in Table 2 and the resulting water retention functions in Fig. 4. The hydraulic conductivity function was modelled as follows: (11) where e is the effective degree of saturation given by:

Exponential model (EXP)
The parameter of the exponential hydraulic conductivity function given by Eq. (4) was derived by fitting the exponential function against the van Genuchten function as shown in Fig. 5 Table 3). 3 Results

Steady-state limiting flow rate
The existence of a limiting flow rate Q through the root that emerges from the analytical solution implementing the exponential hydraulic conductivity function as shown in Eq. (7) was investigated numerically for the case of the more realistic three-parameter van Genucthen hydraulic conductivity functions (Fig. 4). Fig. 6 shows the flow rate Q against the suction differential across the cylindrical water flow domain (Fig.  3) for different values of suction in the bulk soil, sbulk. There appears to be a limiting flow rate that that can be transmitted through the soil-root interface, which depends on sbulk. This confirms that the limited flow rate, Q, as per Eq. (7) is not an artefact of the exponential hydraulic conductivity function. The wide range of suction differential between root and soil is not unusual, it has been observed numerically [16] and experimentally [23].     Fig. 7 shows the evapotranspiration rate ETlim as limited by the soil-root interface versus the suction in the bulk soil, sbulk. ETlim was obtained by deriving the maximum flow rate through the root Qmax, either analytically for the exponential function as per Eq. (7) or numerically for the van Genuchten function as shown in Fig. 6, and then deriving the limiting ET via Eq. (6). Fig. 6Fig. 7 also shows an evapotranspiration of 5 mm/day taken as the potential evapotranspiration, PET. The intersection between the PET and the limiting ETlim represents conceptually the suction s3 that marks the transition from energy-limited to water-limited evapotranspiration in the Feddes function as shown in Fig.  1. The portion of the limiting ETlim beyond the suction, s3, represents the water-limited branch of the Feddes function. When considering the exponential hydraulic conductivity function, the water-limited branch of the reduction function decays very rapidly as compared to the water-limited branch derived by considering a more realistic van Genuchten function for the hydraulic conductivity. This stems from the very rapid decline of the exponential hydraulic conductivity function as emerges when plotting the hydraulic conductivity on a log-scale (Fig. 8). The exponential function declines very rapidly as suction increases and this returns an unrealistic steep water-limited branch of the reduction function as shown in Fig. 7. The exponential hydraulic conductivity is clearly not adequate to represent the water-limited branch of the reduction function. The interest in exploring the exponential function lied on its integrability (as opposed to the van Genuchten function) in turn leading to a closed-form expression for the water-limited branch of the reduction function. Different integrable hydraulic conductivity functions should be explored in the future.

Steady-state versus transient flow
A closed-form solution for the radial flow at the soil-root interface, which is required to derive a closed-form expression for the water-limited branch of the evapotranspiration reduction function, can only be derived for steady-state conditions. The question arises about whether such an assumption is appropriate. To this end, a series of transient water flow analyses were performed (Fig. 9).
Starting from an initial condition of uniform suction across the cylindrical water flow domain, a sinusoidal signal with 1-day period was applied at the outer boundary of the cylindrical flow domain (suction in the bulk, sbulk) while assuming no flux at the root (inner boundary of the cylindrical flow domain). The suction at the root (sroot) was therefore monitored and benchmarked against the sbulk. At low initial suction (Fig. 9a), the hydraulic conductivity of the soil is relatively high and there is no phase shift between the signal imposed at the outer A phase shift is observed if the initial suction is increased and, hence, the baseline hydraulic conductivity is decreased. As expected, the higher the initial suction, the lower the baseline hydraulic conductivity, and the higher is the observed phase shift. The phase shift was equal to 0.1 days for the baseline suction sini=500 kPa (Fig. 9b) and 0.2 days for the baseline suction sini=1000 kPa (Fig. 9c). Because of the reduced radial dimension of the water flow domain (rbulk -rroot = 3.36 mm in this exercise), the phase shift appears to be significantly smaller than the 1-day period of the sinusoidal signal even  for the case where suction is relatively high and, hence, the baseline hydraulic conductivity is relatively small.

Conclusion
This paper has explored the implications of steady-state flow and a simplified hydraulic conductivity function used to develop a closed-form solution for the radial water flow at the soil-root interface, in tun to develop a physically-based closed-form water-limited regime branch of the evapotranspiration reduction function.
The analytical solution based on the exponential function has been benchmarked again the solution derived from numerical simulations using a more realistic van Genuchten hydraulic conductivity function. It has been shown that the exponential function correctly captures the hydraulic response of the soil-root interface from a qualitative standpoint. However, it returns a water-limited branch of the reduction function that declines very rapidly and is therefore not realistic. Future research will therefore explore other more realistic but integrable hydraulic conductivity functions. Transient water flow analyses have shown that the steady-state assumption appears to be acceptable. Suction propagates relatively rapidly through the cylindrical volume representing the soil root-interface due to its millimetric radial dimension.