Calibration of the root water uptake spatial distribution of a young Melaleuca styphelioides, an Australian-native plant, by means of a large-scale apparatus experiment

The stability of slopes is greatly influenced by seasonal variations in pore water pressures (pwp) induced by rainfall infiltration and evapotranspiration processes. Despite that, the prediction of the hydrological effects of long-stem planting is often simplified or neglected because it is challenging to address. Its computation requires a proper definition of the plant root water uptake spatial distribution, which depends, in turn, on geometry and spatial root density. A well-suited case study in this field of application has been provided by a soil-filled embankment, close to an important traffic artery in Newcastle (Australia), which experienced shallow instability. The implementation of long-stem planting has been suggested as a remediation intervention. Based on this, an experimental study focusing on the effects of plant roots on the distribution of pwp in the site soil has been performed by means of a large-scale laboratory experiment on a 2-year-old native plant. Suction measurements were recorded within the vegetated soil mass under controlled boundary conditions and used to calibrate two different root spatial distributions in a seepage simulation. One is based on a flexible RWU spatial distribution function, and the other, specific for the plant RWU pattern, is simpler in its formulation and requires the definition of a lower number of parameters. A comparison between their performances in reproducing pwp distribution suggests that the second one is a better alternative. The methodological approach adopted has proven to be suitable for representing the hydraulic behaviour of a vegetated hillslope, to be eventually implemented in a proper stability assessment problem.


Introduction
Soil bioengineering refers to sustainable and low-impact solutions for soil slope stabilization, surface erosion prevention and environmental restoration of degraded areas by employing living or dead (wood) plant materials, often coupled with inert material (Gray and Leiser 1982). As indicated by its name, soil bioengineering operates in an interdisciplinary field, combining mechanical, biological and ecological concepts. Cost-effectiveness, the aesthetic of the interventions and sustainability are parts of the benefits that make natural bio-engineering techniques globally appreciated.
Among the case of application, long-stem planting has proven to be a valid protection measure against shallow mass movements in slopes, especially the ones triggered by rainfall events, through plant-root transpiration, rainfall interception by plant canopy (Wu et al. 2015) and the mechanical reinforcing effect given by roots permeating a soil volume (Wu 1976;Waldron 1977;Wu et al. 1979).
Living plants enhance soil mechanical behaviour through hydrological effects since soil strength parameters, which govern slope stability, are strongly affected by soil suction and water content. However, the influence of plant transpiration on soil strength, although still acknowledged, needs a more rigorous validation to be properly considered in many practical applications, especially in numerical simulations of local scale problems, such as slope failure assessment (Stokes et al. 2014).
A realistic root water uptake (RWU) spatial distribution and an effective root growth model, providing the RWU spatial pattern over time, are a fundamental component for reliable modelling of pore water pressure (pwp) distributions in the soil close to a high-stem plant. In particular, the RWU distribution is intrinsically related to information that is challenging to estimate, such as lateral spread, depth and architecture of a root system. Moreover, the parameters of the RWU models, that can be adopted to describe the spatial variability of RWU, are hardly measurable directly (Raats 1974;Vrugt et al. 2001), but they are calibrated in practice using observation data so that the model accurately simulates the water dynamics in the proximity of the plant.
In the scientific literature, many measurements of moisture distribution around a root bulb, obtained through soil water content sensors, can be found (Coelho and Or 1996;Vrugt et al. 2001;Zuo et al. 2006;Bufon et al. 2012, among others) for the calibration of the RWU spatial distribution. The use of these monitoring sensors has proven a valid alternative to destructive methods, such as root coring, or nondestructive methods as rhizotube or minirhizotron (Johnson et al. 2001;Rewald and Ephrath 2013), which also provided information on root length density (RLD) and root weight density (RWD).
The idea of using tensiometers to investigate phenomena related to the vadose zone is not new either, especially for agricultural purposes and for optimal water management of cultivated plants (Dirksen 1979;Bacci et al. 2003;Koumanov et al. 2006;Coelho et al. 2007 among others). Their application for slope stability issues has still been adopted in geotechnical research (Toll et al. 2011;Rocchi et al. 2018;Gragnano et al. 2020Gragnano et al. , 2021 but, in this field, no example of calibration of RWU models using water potential data exists to the knowledge of the authors. While water flow simulation within the vadose zone is performed as a routine in a one-, two-or three-dimensional space, the RWU is usually considered as only function of the vertical dimension. This is demonstrated by the large availability of one-dimensional (1D) uptake models in which the variability of RWU with depth is dealt with different approaches (Hoogland et al. 1981;Jarvis 1989;Gardner 1991;Berardi et al. 2022, among others). Using 1D RWU models is appropriate for crops characterized by a uniform RWU pattern and results an advantageous choice to reduce the number of model parameters. However, row crops, row high-stem plants and single-stand high-stem plants require the adoption of multidimensional RWU models due to an increased complexity in water uptake patterns and the description of the water uptake process (Vrugt et al. 2001). To these models correspond naturally a higher number of model parameters, most of them with no physical meaning, and thus difficult to properly estimate.
The present study describes the results of a large-scale laboratory experiment on a native Australian plant, Melaleuca styphelioides, whose influence on pwp distribution in soil has been investigated under different controlled boundary conditions. Water potential data were collected at different distances from the root bulb and used to calibrate and validate the two-dimensional (2D) RWU spatial distribution of the plant. Two different functions have been chosen for that purpose; one is the flexible function of Vrugt et al. (2001), potentially applicable to all RWU distribution patterns, and the other is user-defined, simpler in its formulation and particularly suited for this bulb architecture, and it requires the definition of a lower number of parameters. This study falls in the context of slope stability remediation by long-stem planting in Newcastle, NSW, Australia. The commercial code Hydrus 2D by PC-Progress has been used to simulate the water movement in the whole flow domain containing the plant (Simunek et al. 1999). The size of the described laboratory experiment offers the potential of controlled conditions while eliminating possible edge effects. The result is an accurate modelling of the evapotranspiration phenomenon for the investigated plant, with a good representativeness of the real site conditions.

Field site
A partially back-filled slope behind a retaining wall overlooking the Newcastle Inner City Bypass Cut 7, West Charlestown in Newcastle (NSW, Australia), has experienced shallow instabilities in March 2017 triggered by heavy rains (Fig. 1A, B and C). In attempting to solve this issue, several drains have been installed along the toe of the slope. The slope is inclined at 1.5 (H):1 (V) (i.e. ~ 34%) and then flattens out toward the retaining wall and the uppermost crest. The slope is constituted by a topsoil (A -Horizon, ~ 0.5 m thick) overlying a clayey layer (B-Horizon, 1 m thick) on top of a poor-quality rock. In a 10-year period, from 2007 to 2017, a decay of the vegetation has been observed (see Fig. 1F). Chemical analyses have shown that the soils characterizing the two horizons, in that particular location, are extremely acidic: horizon A has a pH ranging between 4.8 and 5.4 while for horizon B, it is between 4.9 and 5.0. The high level of hydrogen, which causes the soil acidity, has resulted in the availability of aluminium (between 7.8 and 84.6% for Horizon A and between 81.3 and 85.4% for Horizon B) which is toxic to plants (SESL AUSTRALIA Report 2019). To enhance slope stability, a new campaign of high stem planting has been proposed and designed. The Melaleuca styphelioides, a native eastern Australia plant, has been chosen to revegetate the slope ( Fig. 1D and E). It is an evergreen medium tree (height between 5 and 11 m and foliage spread of 2-3 m), tolerant to drought, frost, pollution and soil saturation. Because of current soil acidity, the soil volume that will host a long-stem plant needs to be dug to a depth of 1 m (for a total volume of removed soil equal to ~ 1 m 3 ) and substituted with a new mixture of soil. A blend of soil, topsoil, lime, sand and appropriate fertilizer has been proposed (SESL AUSTRALIA Report 2019). Lime needs to be added to the mixed soil to the rate of 0.04% of the mass of the dry soil in order to neutralize soil acidity and increase pH to ~ 6. The soil needs to be mixed with an additional 65%, in weight, of clean coarse river sand in order to increase the soil porosity and, in turn, the moisture available to the plant and its aeration state (Tuyishimire et al. 2022). The physical, hydraulic and retention properties of the in situ topsoil, of the improved soil blend and of the organic potting mix (which contains the plant roots) have been characterized by proper laboratory tests and are presented in the next section.

Topsoil (A-Horizon) and river sand
According to the Australian Standard AS 1726 (2017), Horizon A is an organic sandy clay at the boundary between medium and high plasticity, made of 73.5% of fine fraction (clay and silt) and 26.5% of coarse fraction (see Table 1; Fig. 2A). The organic matter content is high, around 4.6%. The grain size distribution has been obtained combining the standard sieving of the coarser fraction with the sedimentation of the finer fraction performed with an automated particle size analyser (SediGraph III by Micromeritics) (Micromeritics 2021a). The specific weight of the solid particles (G s ) has been determined using a gas pycnometer (AccuPyc II by Micromeritics) (Micromeritics 2021b). The saturated permeability has been determined by a falling head permeameter. The clean river sand, used as additional component of the A-Horizon soil, consists of 22.4% of gravel and 77.4% of sand (see Table 1; Fig. 2B). With a uniformity coefficient (C u ) of ~ 4.25 and a curvature coefficient (C c ) of ~ 1.81, the material is well-graded. Table 1 summarizes the main physical properties of the A-Horizon and of the river sand.

Improved mix of A-Horizon and river sand
As reported in "Field site", the improved soil has been obtained by mixing the site material (A-Horizon) with an additional 65%, in weight, of clean river sand. Table 2 reports a summary of the basic geotechnical properties of the mixture, and in Fig. 3, its grain size distribution is presented. A proctor compaction test was also conducted to identify the optimum moisture content (14.7%) and the maximum dry density (2.10 Mg/m 3 ) of this soil.
The water retention curve of the soil was measured using sample reconstituted at an initial void ratio of 0.8. The sample has been prepared to the target void ratio in fully saturated condition in an oedometer cell then cut and placed in a cylindrical airtight plastic container. The evaporation method (Gardner and Miklich 1962;Wind 1968) was then applied, where the sample is incrementally dried, and after Table 1 Main physical properties of the A-Horizon and of the river sand as the organic content, the soil specific gravity (G s ), the Atterberg limits, the Plasticity Index (PI), the saturated coefficient of per-meability (Ks), the grain size distribution, the uniformity coefficient (Cu) and the coefficient of curvature (Cc) The standard deviation of each physical property is reported in brackets while the number of tests performed for each soil property is reported in the 3rd and 5th column each period of moisture equilibration (~ 2 days), suction measurements are recorded using a high-capacity tensiometer (measurement range up to 1 MPa, accuracy ± 10 kPa) or a dewpoint potentiometer WP4-T (Decagon Device 2015).
To reproduce field conditions, soil samples were allowed  to shrink along the drying branch of the retention curve, and the volume change was tracked using the hand spray plaster method developed by Liu and Buzzi (2014). Figure 4 reports the SWRC in terms of gravimetric (LEFT) and volumetric water content (RIGHT) vs suction as obtained for the mixed soil. The saturated permeability of the improved soil has been determined in a Rowe Cell. The sample was prepared at a void ratio of 0.8 and at an initial water content of 5% (air dried soil). A constant pressure gradient was then applied between bottom (connected to a pressure controller device injecting water) and top (at atmospheric pressure and connected to an external precision balance) of the soil sample until a constant rate outflow (i.e. full saturation) was reached. A saturated hydraulic conductivity value of 2.93 × 10 −7 m/s was measured. Mercury intrusion porosimetry (MIP) was conducted to obtain the pore-size distribution of the improved mixed soil (see Fig. 5). The soil specimens were prepared by freeze drying, which is a well-accepted preparation technique. For matter of conciseness, the technique is not fully described here, and the reader can refer to Yuan et al. (2020) for more details. The improved soil shows a dual porosity with diameter of the dominant macropore and micropore of 0.07 mm and 0.00026 mm, respectively. The predominance of macropores over micropores ensures a high aeration of the soil (corresponding air entry value of ~ 2 kPa), an easier growth of roots and an access for plants to a larger quantity of free water and nutrients due to the intense microorganism activities hosted (Zheng et al. 2021).

Potting mix
All Melaleuca styphelioides plants used for this study were grown in an organic potting mix, and it is very difficult to separate this material from the root bulb without damaging the root system, especially the fine root system (i.e. roots ≤ 2 mm in diameter) which represents for the plant the principal pathway for water and nutrient absorption (Eissenstat 1992). For this reason, the potting mix has been kept when placing the plant in the laboratory experiment. It is also planned to do the same for the field planting. The potting mix is a blend of sphagnum peat moss, perlite, coir fibres, sand, limestone and fertilizers. Laboratory tests have been carried out on this material to evaluate its physical, hydraulic and retention properties. Table 3 reports the variation range of the main physical properties of the investigated soil while Fig. 6A the grain size distributions of two potting mix samples.
The drying branch of the SWRC has been determined similarly to what was done for the improved soil mix ("Improved mix of A-Horizon and river sand"), but suction measurements this time were performed at a constant void Fig. 4 The main drying branch of the Soil Water Retention Curve (SWRC) in terms of gravimetric water content (LEFT) and volumetric water content (RIGHT) vs suction of the improved soil (A-Horizon + river sand) ratio (~ 2). A loose and saturated potting mix is progressively dried and each time compacted to the target void ratio. Suction is then measured by a high-capacity tensiometer up to the cavitation pressure (1 MPa), then in a dewpoint potentiometer WP4 device. Figure 6B reports the experimental points of the main drying curve which has clearly a bimodal trend. Such soil pore size distribution is characteristic of highly organic soils (Burger and Shackelford 2001) as the one under analysis. The saturated permeability of the potting mix has been determined through a constant head permeability test and estimated equal to 8.42 × 10 −6 m/s.

Description of the apparatus
A specific apparatus has been designed and developed in order to investigate the spatial and temporal evolution of soil suction in close proximity of a 2-year-old Melaleuca styphelioides plant (see Fig. 7). The setup is composed of a 1 m 3 container filled with 720 mm of mixed soil overlying 130 mm of clean gravel. A non-woven fabric separates the two layers. The improved soil has an initial water content equal to ~ 10%, and it has been arranged in the containers in 50-mm-thick layers; no compaction has been performed during the setup, and the resulting void ratio is around 0.8 (average value between the denser lower layers and the looser upper ones). The plant is positioned in the centre of the container, under 150 mm of soil cover. The dimensions of the root bulb correspond to the ones of the pot where the plant grew, which was entirely filled by the root bulb at the time of plant installation. The underside of the root bulb is positioned at an approximate depth of 290 mm from the soil Results of mercury intrusion porosimetry of the improved soil (A-Horizon + river sand) prepared at a void ratio of 0.8 and moisture content of ~ 10% under ambient humidity: pore size distribution (left); cumulative intrusion curve (right) Table 3 Main physical properties of the potting mix, in which the plants have been grown in the nursery, as the organic content, the soil specific gravity (Gs), the saturated coefficient of permeability (Ks), the grain size distribution, the uniformity coefficient (Cu) and the coefficient of curvature (Cc) The standard deviation of each physical property is reported in brackets while the number of tests performed for each soil property is reported in the 3rd column surface and 430 mm from the top of the gravelly layer. A flexible plastic tube is placed along the plant stem to irrigate the bulb just after the replant. A water reservoir is connected at the base of the apparatus to control the water table level. A float valve installed at the water reservoir inlet maintains a constant water level, which is the hydraulic bottom boundary condition of the physical model. Fertilizer is poured into the water reservoir to supply the plant with proper nutrients. A plastic cover resting on the water reservoir is used to minimize evaporation. After an initial stage of irrigation of the plant from the tube positioned along the stem (three irrigation events), the groundwater table is imposed (in steps) to sustain plant life by raising capillary fringe, increasing the soil water content where the RWU occurs.
A plastic cover and insulating panels are placed at the top of the container to create an enclosed environment where relative humidity (RH) can be controlled. This is achieved by a device that circulates air (via a pump) from a water reservoir to the enclosed volume at the physical model top and back. Two RH sensors are placed to monitor RH above the soil surface and to guide RH apparatus functioning. The flow of air is adjusted until the RH in the enclosed space matches the target RH. Temperature is not controlled at the top boundary in this experiment, but the whole apparatus is located indoor, and the insulating panels reduced the daily temperature oscillations. The temperature and RH of the laboratory room, where the physical model is stored, are monitored by a sensor fixed to the plant stem and acquired by a data logger (see Fig. 7).
An indoor grown light (model Hi-Par Dynamic CMH 315w DE Control Kit) was used to provide a full spectrum light optimal for plant growth; from 10 to 12 h of light was provided per day. The lamp has been placed 60 cm from the plant top, providing a PPFD (photosynthetic photon flux density) of 856.61 μmol m 2 s −1 , equal to 60,819.31 lx (for a ceramic metal halide lamp) or 183 W/m 2 .
In order to monitor suction evolution in time, 14 TEROS 21 sensors by Meter group (Meter Group 2021) were installed in close proximity to the plant and in a control zone, unaffected by root bulb effects. TEROS 21 measures the dielectric permittivity of two ceramic disks in moisture equilibrium with the surrounding soil volume. The dielectric permittivity is then converted to matric potential and volumetric moisture content based on the moisture characteristic curve of the ceramic. The working range of the water potential is − 9 to − 100,000 kPa, with a resolution of 0.1 kPa. It is relevant to highlight that temperature changes can greatly affect TEROS 21 readings below − 500 kPa. The sensor measures temperature using a surface-mounted thermistor (range − 40 to 60 °C, resolution 0.1 °C; accuracy ± 1 °C). Even considering that TEROS 21 sensor estimates indirectly water potential and it is unable to measure pressures greater than − 9 kPa; on the other hand, it does not suffer, as tensiometers, from cavitation problems, and it does not require a user-defined calibration (it is provided with a manufacturer 5-points calibration curve). Figure 8 shows two cross-sections of the large-scale apparatus with the location of the TEROS 21 sensors with respect to the root bulb and surface. To monitor the influence of transpiration on suction data close to the plant (section A-A, Fig. 8), five sensors (from no. 1 to no. 5) were installed on the side of the root bulb at the same depth of 220 mm from the soil surface, horizontally spaced by 50 mm and with the closest to the bulb (no. 5) at a distance of 30 mm from it. Five sensors were installed at an incremental depth below the root bulb: the shallower (no. 10) at a distance of 30 mm, the remaining ones (no. 9, 8, 7, 6) vertically spaced by 50 mm. Another sensor (no. 11) was positioned at a radial distance, below the bulb, 80 mm from the root bulb underside (vertical distance) and 145 mm from the centreline axis (horizontal distance). A corner of the experiment (section B-B, Fig. 8) is separated from the remaining soil volume by means of two impervious boundaries to create a control zone unaffected by plant transpiration. Three sensors (nos. 12, 13, 14) were installed at an incremental depth (75 mm, 320 mm and 470 mm) from the soil surface in that control  Figure 9 shows the time-lapse of the large-scale apparatus setup through a series of photographs.

Experimental program
RH and water table height have been changed multiple times during the experiment in order to drive the spatial and temporal evolution of suction around the root bulb under different boundary conditions and to collect relevant information to calibrate and validate a RWU model of the investigated plant, during an overall period of 308 days. In Table 4, a summary of the six different experimental phases, with the applied boundary conditions, is reported. Phases 0 to 3 include the application of an increasing water level to the bottom boundary (from 135 to 432 mm from the base of the physical model) while phase 4 sees a reduction in water level with respect to phase 3 (from 432 to 312 mm). A variable RH was applied to the upper boundary ranging from ~ 99% for phase 1 (only transpiration effects on pwp distribution) to ~ 74-87% for the remaining phases (combination of evapotranspiration effects), while at phase 0, the RH was not controlled. These boundary conditions are equally applied to the soil volume containing the plant and to the control zone, to allow a comparison between the two sites.

Evapotranspiration data
In Fig. 10, the boundary conditions applied to the largescale apparatus are presented. In Figs. 11 and 12, the pwp and temperature evolutions are plotted with time, as recorded by sensors during the different testing phases. Each phase will be analysed specifically in the next sub-sections.

Phase 0
Boundary conditions Phase 0 (day 0 to day 30) was designed to be an equilibration phase, which follows the initial setting up of all the large-scale apparatus and the planting of the plant in the centre of the container. In this phase, the upper boundary was not covered, and the RH was not controlled (Table 4). At this stage, the container was filled with dry soil (initial suction around 800-900 kPa). The root bulb was irrigated via a plastic tube along the plant stem twice (at days 1 and 14) to keep the plant alive prior to the imposition of the water table. At day 15, the water table was set at 715 mm from the soil surface (5 mm above the gravelly layer, 195 mm below the lowest sensor, no. 6) (see Fig. 10A). Sensors in the control zone Sensor readings in the control zone (see Fig. 11A) show significant fluctuations with some values sometimes out of range. The noise is caused by the strong influence of temperature variation on suction measurements above 500 kPa, as reported also in TEROS 21 Technical Manual (Meter Group 2021). This influence is observable to a lesser extent also at lower suction values since total matric suction, temperature and RH are physically related by psychometric law (Fredlund and Rahardjo 1993). This explains the same trend easily observed comparing temperature (Fig. 11A) and suction values (Fig. 11B) in time (phases from 0 to 1A). As the cover depth increases and provides temperature insulation, these oscillations tend to reduce.
Sensors on the side of the root bulb After the irrigation at day 14, it has been observed that sensor no. 5 generates a suction increment of ~ 600 kPa in the following 17 days (see Fig. 13A), almost halfway to the wilting point, while the remaining sensors on the side of the bulb showed lower increments, below 130 kPa. Sensor no. 1 is still too far (230 mm) from the root bulb to be rapidly responsive by this wetting event and by the RWU, showing a slow equilibration with the local soil moisture. Fig. 9 Photographs of the largescale apparatus during the setup. A shows the installation of the non-woven fabric at the boundary between the gravelly layer and the first layer of improved soil (A-Horizon + river sand). B shows the installation of sensor no. 7 below the root bulb while C the installation of the sensors (from no. 1 to no. 5) on the side of the root bulb. D shows the planting of the young Melaleuca styphelioides and E the installation of the sensors TEROS 21 in the control zone Sensors below the root bulb Sensors no. 9 and no. 10, just below the root bulb, are affected by the wetting along the plant stem performed at day 14, and, in 17 days, they show a suction increment of ~ 250 and ~ 130 kPa (see Fig. 13B), respectively. The lower suction increment of sensor no. 10 with respect to sensor no. 5 can be due to a difference in root density and plant absorption between horizontal and vertical direction (unverified hypothesis) and/or to a major presence of water which ponds at the bottom interface between the highly permeable potting mix and the less permeable mixed soil. Sensor no. 7 starts recording a higher suction value at time zero (~ 1600 kPa, in Fig. 12C) compared to the other sensors (below 800 kPa), eventually due to a potential inhomogeneous initial water content, which is equalized by the water table imposition and the subsequent seepage.
Sensor positioned radially from the root bulb Sensor no. 11 shows the same trend of sensor no. 9 (see Fig. 12C), placed at the same depth under the bulb, but showing higher suction values due to a lower influence of irrigation events (water percolates mainly vertically in macropores for gravity effects).

Phase 1
Boundary conditions In phase 1, the RH has been initially set at 100%, to avoid evaporation from the soil surface (phase 1A), then let free to change (phase 1B) to observe the combined effect of evaporation and transpiration on suction distribution. The RH was maintained in phase 1A (day 31 to day 81) equal to a value of 99.5%, while in phase 1B (day 82 to day 152) to 73.8% (see Table 4; Fig. 10).
The water level was set at 200 mm from the bottom of the box in both phases. The plant was irrigated a third and last time at day 31, which marks the initial step of phase 1. No suction or temperature information was recorded between days 71 and 117 due to a logging issue.
Sensors in the control zone Suction readings of sensor no. 12 fell out the measurement range (water potential − 100 MPa) in the interval between 71 and 117 days, due to evaporation effects (Fig. 11A). The upward moisture migration, subsequent to water table imposition, causes a drop of suction values recorded by sensor no. 6 and no. 13.
Sensors on the side of the root bulb As still previously observed, sensors no. 5 no. 4, no. 3 and no. 2, positioned on the side of the root bulb (Fig. 12A), are affected by an immediate suction decrease as a consequence of plant irrigation, followed by an increment driven by plant transpiration (phase 1A). Sensor no. 5 shows a rapid increase in suction because of its vicinity to the root bulb (high RWU Fig. 11 The graph (A) shows the evolution of the soil suction with time recorded by the installed sensors in the control zone while (B) the temperature variations recorded by TEROS 21 sensors during the test contribution), while the remaining sensors only show a slight variation in data. Readings for sensor no. 1, oscillating in phase 1A in the range 400-600 kPa, seem to be affected mainly by the water content of the surrounding soil rather than the irrigation events and the transpiration contribution. In addition, sensor no. 1 is also significantly affected by temperature changes which cause oscillation of suction readings above 500 kPa, as discussed before, for limits of the specific sensor. Apart from sensor no. 5, suction seems to slightly increase with radial distance from the root bulb. At day 71, sensors no. 4 to no. 1 recorded 216, 285, 290 and 497 kPa of suction, respectively. This is due to a lower imbibition of the soil by plant irrigation moving away from the root bulb coupled to a limited or null absorption of water by RWU. Between days 71 (nearly the end of phase 1A) and 120 (phase 1B), sensors on the side of the root bulb (from no. 5 to no. 1) have recorded a suction increase in absolute value of 204, 27, 51, 62 and 136 kPa, respectively. It seems that, apart from sensor no. 5 which responds mainly to plant RWU, the remaining sensors are affected differently by evaporation from soil surface, albeit under the same soil cover. This could be explained by a different degree of moisture content of the superficial soil layer which is dryer close to the plant for RWU influence and wetter moving away from the bulb. The formation of a dry crust and the reduction of the coefficient of permeability according to the hydraulic conductivity function could create an unevenly evaporation flux through the soil surface, even at the same depth.
Sensors below the root bulb The rate of suction increment observed after the irrigation event at day 14 (phase 0) occurs again in correspondence of the third irrigation, at day 32. Sensor no. 10 shows an increment of 194 kPa, much lower than the one observed for sensor no. 5 showing 652 kPa of increment over the same period. The same hypotheses about this peculiar behaviour exposed in phase 0 could be applied again in phase 1A. However, it is interesting to compare sensor no. 4 and no. 9 ( Fig. 12A and C, respectively): the suction increment over the same period (from day 32 to 43) is lower on the side with respect to the bottom of the bulb (49 kPa vs 120 kPa). This could be due to an initial major availability of water underneath the root bulb (97 kPa of suction for sensor no. 9 vs 159 kPa for sensor no. 5 at day 32) for a fast downward percolation of the irrigated water and, as a result, greater water absorption rate in the vertical direction. Sensors below the root bulb, at the end of phase 1A, clearly show the influence of RWU on suction evolution. Moving vertically from the bulb to the bottom of the model, suction changes from 233 kPa (no. 10) to 219 kPa (no. 9) and to 59 kPa (no. 8). Readings of sensors no. 7 and no. 6 result out of the optimal measuring range (suction below 9 kPa) due to the upward seepage from the imposed water table. There is, indeed, a very little change in suction trend passing from phase 1A to 1B, i.e. allowing evaporation from soil surface. In fact, suction increment between day 71 and day 117 is equal to 37 kPa for sensor no. 10 and 20 kPa for sensor no. 9. No suction increments are recorded for the deeper sensors, meaning evaporation and transpiration have no tangible effect at these depths.
In order to understand the different impact of evaporation and transpiration over suction, sensor no. 10 is compared to sensor no. 13 and sensor no. 7 to sensor no. 14, as being no. 13 and no. 14 positioned at the same depths but in the control zone (compare sensors no. 10 and no. 7 in Fig. 12C with sensors no. 13 and no. 14 in Fig. 11A). For instance, in the period between 71 and 117 days, sensor no. 13 has shown a suction decrement of ~ 485 kPa due to the upward seepage gradient, while sensor no. 10 a suction increment of 37 kPa. Such significant difference can be attributed only to transpiration contribution. Sensors no. 14 and no. 7 show on the same period a similar suction variation (~ 1.5 kPa) confirming that at this depth, suction is not affected by RWU or by a significant evaporative flux.
Sensor positioned radially from the root bulb Sensor no. 11 shows an intermediate behaviour between sensor no. 10 and no. 9 (see Fig. 12C). Indeed, in the temporal period 71-117 days, it is affected by a suction increment of ~ 30 kPa, while sensors no. 10 and no. 9 by an increment of 37 and 20 kPa, respectively. It is likely that the lower RWU contribution in correspondence of sensor no. 11 (due to a greater distance from the bulb with respect to sensor no. 9) is balanced by a higher contribution of evaporation. As also observed in phase 0, evaporative flux seems to increase moving horizontally from the root bulb. Moreover, to suction values close to ~ 270 kPa, characterizing sensors no. 9 and no. 10 below the root bulb in the period 71-117 days, corresponds the inflection point of the bimodal SWRC of the potting mix. As reported by Simunek et al. (2012), the hydraulic conductivity function has a strong slowdown in correspondence of this point. To a reduction of the coefficient of permeability corresponds a simultaneous drop of the evaporative contribute.

Phase 2
Boundary conditions At the end of phase 1, sensor no. 5 is approaching 1400 kPa of suction. Due to the high risk to reach the wilting point of the plant, at day 152, the water table was raised to 355 mm from the base of the apparatus. This led to an increase in the water available for plant Fig. 12 Graphs A, C report the evolution of the soil suction with time recorded by the installed sensors on the side of the root bulb and below the root bulb, respectively. In graphs B, D, the temperature variations recorded by TEROS 21 sensors during the test are shown ◂ transpiration in the soil volume depleted of moisture in the previous phases. The RH maintained in the apparatus was 81.2% (standard deviation 3.71%) ( Table 4).

Sensors in the control zone
Sensor no. 13 only shows a suction reduction of 0.5 kPa/day suggesting that equilibrium requires a very long time at that point (Fig. 11A). Sensor no. 14 located in the control zone at a depth of 450 mm from soil surface shows no changes in suction readings during phase 2, evidence of the fact that equilibrium has been reached (see Fig. 11A).
Sensors on the side of the root bulb Sensor no. 5 shows a slow increase in suction which is the combined effect of evaporation and transpiration processes (see Fig. 12A). Suction readings of the remaining sensors, on the side of the plant, seem unchanged from phase 1B. Sensors seem affected by temperature change, especially sensors no. 5 and no. 1, which show the same fluctuations but wider, due to their higher suction values above 500 kPa.
Sensors below the root bulb Between days 162 and 213, sensors no. 10 and no. 9 have experienced a slight increase in suction (17 kPa and 2 kPa, respectively), while sensor no. 8 an overall reduction of 5 kPa due to the upward seepage. Sensors no. 7 and no. 6 are in hydraulic equilibrium (see Fig. 12C). Comparing sensors no. 13 and no. 10 (Figs. 11A and 12C), it is possible to observe that sensor no. 10 stands at higher suction values and shows opposite trend with respect to no. 13; in fact, suction is slightly decreasing at sensor no. 13 while increasing at sensor no. 10 due to plant transpiration activity.
Sensor positioned radially from the root bulb Sensor no. 11 is the most affected by temperature fluctuation because at highest suction value with respect to other sensors below the root bulb (see Fig. 11A and B). Due to temperature-induced oscillations of suction readings, it is difficult to interpret if suction is slightly increasing or oscillating around the same value.

Phases 3 and 4
Boundary conditions (phase 3) At day 220, the water table was raised of ~ 77 mm with respect to phase 2, while RH was kept the same adopted in the previous phase (see Table 4; Fig. 10).

Sensors on the side and below the root bulb (phase 3)
At day 232, 12 days after the water table was raised, a significant drop of suction was reported for all sensors at the side of and below the root bulb ( Fig. 12A and C). Around 30 days were necessary to reach equilibrium (232-262 days) with a final suction value of ~ 50/60 kPa.

Boundary conditions (phase 4)
In phase 4, the water table was lowered of 120 mm to avoid damage to the root system for prolonged anaerobic conditions. The RH was set to 100% (Table 4), eliminating the evaporative contribution through soil surface.

Sensors in the control zone (phase 4)
In phase 4, sensor no. 12 was unaffected by the boundary conditions change (suction values out of range) while the other sensors, underneath it, not influenced by plant transpiration, reached hydraulic equilibrium in the previous phase 3 (see Fig. 11A).
Sensors on the side of the root bulb (phase 4) During this phase, sensors no. 2, no. 3, no. 4 and no. 5 equalized around a value of ~ 25 kPa while sensor no. 1 remained unchanged from the equilibrium condition reached in the previous phase, as it is not influenced at any level by root water absorption (Fig. 12A). Important variations in temperature were observed during phase 4 (temperature oscillates between 16 and 24 °C) (see Fig. 10B), but suction readings in presence of a soil water content close to saturation seem unaffected by them. From this, it seems that when water content variation/water seepage in the soil is induced by a significant waterfront movement, temperature effects on suctions are negligible; otherwise, as long as water seepage is driven by a thermohydraulic process (e.g. evaporation, transpiration), suction and temperature variations are more strictly related.
Sensors below the root bulb (phase 4) Sensors no. 10 and no. 9 slowly equilibrated to ~ 20 kPa, while the further sensors (no. 8 to no. 6) were already in equilibrium with the imposed boundary conditions at the beginning of the phase and rested unaffected by the water table height variation throughout the whole investigated time period (Fig. 12C).

Experimental evidences and remarks
The influence of the two processes (evaporation and transpiration), acting together or separately on the evolution of suction with vertical and radial distance from the bulb, is displayed clearly in Fig. 14. In particular, Fig. 14A shows the changes in suction with depth at different instants of the test. The lower suction value recorded on day 35 at a depth of 220 mm from the root bulb, with respect to the higher depths, is due to the position of the water table and the relative increase in water content for upward seepage. Considering the suction profiles, it is clear that the closer the sensor is to the bulb, the higher the suction recorded due to RWU effects, while suction below the bulb incrementally increases as time passes. Because of the large amount of soil cover, evaporation has very limited effects at higher depths. Figure 14 B shows the evolution of suction with radial distance from the bulb at different times of the experiment. Day 35 has been considered as the initial reference condition because the lowest suction value is recorded in correspondence of that day in proximity of the bulb due to plant irrigation, while a suction increment is registered away from the plant. The downward shift of the suction profiles is due to evaporation through the soil surface, while the RWU contribution causes an increment of suction in correspondence of the first sensor The increase in suction for horizontal distances from the bulb side greater than 80 mm is essentially due to evaporation which suggests that the influence of RWU in the horizontal direction is limited to distances lower than that, as measured from sensors no. 5 and no. 4. The RWU influence zone in the vertical direction results between 80 and 130 mm from the bulb side, as shown by sensors no. 10 and no. 9.
At the end of the test (after phase 4, at day 309), the soil has been carefully removed around the plant, revealing no changes in the root bulb shape, as a result of the relatively short duration of the experiment, compared to the root growth rate of the plant. This observation follows the analysis of the suction data that showed, for the entire experiment, a reduced volume of soil interested by RWU, whose dimensions are limited to 80 and 130 mm from the bulb sides, in the horizontal and vertical directions, respectively. A greater density of roots in correspondence of the plant underside, where ponding of water occurs during irrigation events and where soil reaches higher water content for its proximity to the water table, is the most likely reason for the major root water uptake observed in the vertical direction.

Numerical modelling
The commercial code Hydrus 2D by PC-Progress (Šimůnek and van Genuchten 2008) has been used to simulate numerically the water flow within the soil mass in the different phases of the experimental program. The code has implemented a macroscopic empirical model in which a distributed sink term, representing roots extraction, is added to the governing water flow (more details in "Water movement and RWU ").

Potential evaporation and potential transpiration
Potential evapotranspiration from the reference surface (ETo) is computed using the Penman-Monteith method (Penman 1948;Monteith 1965) in the last revision by Allen et al. (1994), which requires atmospheric data as temperature, RH, solar radiation and wind speed. Two different ETo contributions have been used: -ETo _int , which is computed according to the atmospheric conditions inside the enclosed space on the top of the physical model. In this case, the RH imposed is that monitored by the sensor positioned in proximity of the soil surface, the wind speed and the solar radiation are null, the temperature is in equilibrium with the external environment and monitored by a thermometer. -ETo _ext , which is computed according to the atmospheric conditions of the external environment (the laboratory where the apparatus is stored). Also in this case, the wind speed is null, the RH and temperature are monitored by the sensor fixed to the plant stem, the solar radiation coincident with the one supplied by the grow lamp.
The terms of potential evaporation ( E p ) and potential transpiration ( T p ) are computed using the Beer's law (Beer 1852) with the partitions of the solar radiation component of the energy budget according to its interception by the canopy (Ritchie 1972), as reported in Eqs. (1) and (2): The LAI is the Leaf Area Index (dimensionless) and is defined as the average total area of leaves (one side) per unit area of ground surface (L 2 /L 2 ), and it ranges from 0 (bare ground) to 10 or more (dense forest); SCF is the soil cover fraction (dimensionless), and k is a constant governing the radiation extinction operated by the canopy (dimensionless) which is set equal to 0.463 according to Ritchie (1972).
For the present application, the values of ETo _int are used to compute E p in Eq. (1), treating the black plastic cover sealing the soil surface as a dense plant canopy with a LAI of 10, while the values of ETo _ext are used to compute T p using the LAI of the plant (1.30) in Eq. (2). In specific, the LAI of the plant has been measured by direct calculation of the leaf area (Gower et al. 1999;Jonckheere et al. 2004) which were harvested from a Melaleuca styphelioides plant of the same age and grown in the same nursery of the one in use in the apparatus. The potential transpiration ( T p ) is distributed over the root zone unevenly, according to the root architecture and, generally, to the root length density (RLD) of a specific plant (Vrugt et al. 2001;Steudle et al. 1987). These root characteristics are synthetically taken into consideration in the numerical model by means of a RWU spatial distribution function (more details in "Water movement and RWU ").

Actual evaporation and actual transpiration
The potential evaporation ( E p ), obtained as described in "Potential evaporation and potential transpiration", is used as input data for the transient top boundary condition of the FE model (i.e. the atmospheric condition) to compute the actual evaporation flux, according to the availability of water at the soil surface and to the imposed limit value of water pressure head (h crit ). As long as the (negative) pressure head at the soil-atmosphere interface is higher than h crit , the actual evaporation is equal to the potential evaporation; as this limit value is exceeded, the evaporation rate decreases abruptly because the soil permeability becomes too low to deliver the same flux rate (Simunek et al. 2012). The h crit is independent on the soil type and can be estimated from climatic variables such as temperature and RH according to the following relation (Eq. 3): where R is the universe gas constant (= 8.314 J/(mol K)), T is the absolute temperature (K) and M the molecular weight of water (= 0.018015 kg/mol). The value of h crit has been estimated, for each instant of the simulation, starting from measured atmospheric data, and it ranges between − 1.46 and − 1.58E −4 m.
The actual transpiration is computed from the potential transpiration T p according to a stress reduction function α(h), which reflects the assumption that the water absorption by roots decreases with increasing soil suction (h). (h) is a dimensionless function that ranges between 0 and 1.
The stress reduction function implemented in the code is referred to Feddes et al. (1978) and has been adopted for the presented case study. It requires 7 empirical parameters (h 1 , h 2 , h 3low , h 3high , h 4 , T p1 , T p2 ). As depicted in Fig. 15, h 1 is considered an anaerobiosis point (soil is close to saturation and for h > h 1 , the transpiration is null due to oxygen lack), and h 4 is the wilting point below which (i.e. h < h 4 ) transpiration is null. Between h 1 -h 2 and h 3 -h 4 , the water uptake increases/decreases linearly with h while in the interval h 2 -h 3 , α is equal to 1. For values of h below h 3 (defined as h 3low and h 3high ), the potential transpiration is reduced to the actual one, and the reduction start is dependent on the rate of transpiration (T p1 and T p2 , where T p1 > T p2 ). This is because the stomatal closure (and the simultaneous reduction from potential to actual transpiration) is triggered at higher value of water pressure head (h 3high ) if the transpiration rate is higher (T p1 ), while at lower h (h 3low ) if the transpiration rate is lower (T p2 ). Table 5 reports the parameters of the Feddes et al. (1978) stress reduction function adopted in the FE simulation based on previous studies of Wesseling (1991) and Taylor and Ashcroft (1972), which have provided reference values for many different plants and result the most widely used database currently available, implemented also in hydrological software as HYDRUS and SWAP.

Water movement and RWU
In order to simulate the pwp variations in the apparatus due to the upward capillary flux from the imposed water level and the outward evaporation flux from the soil surface, Richards's (1931) equation has been adopted. For the studied case, the incorporation of a sink term S(x, z, t), defined as the volume of water absorbed by roots from a unit volume of soil [L 3 L −3 T −1 ], allows to account for the RWU by the high-stem plant (Eq. 4): where h is the pressure head (L), x and z are the coordinates in the 2D space (L), the latter taken positive downward, K(h) is the unsaturated hydraulic permeability (LT −1 ), and is the volumetric water content (L 3 L −3 ). The first two terms in Eq. (4) capture the effects of capillarity associated with the diffusion transport, while the third term represents the effects of the gravity-driven flux. The soil water constitutive relations (i.e. the hydraulic conductivity function K(h) and the soil water retention function (h) ) are the subjects of the next "Soil hydraulic behaviour", while the formulation of the sink term is discussed below.
The sink term of the actual RWU , S(x, z, t) , assumes the following formulation under a non-uniform root distribution: where S p (x, z, t) represents the potential RWU rate (maximum rate occurring when water is not limiting plant transpiration) (L 3 L −3 T −1 ), and (h) is the stress reduction function ( −) introduced in "Actual evaporation and actual transpiration".
Under non-uniform root distribution, the spatial variation of S p (x, z, t) is expressed by the following formulation: where L t (L) is the width of the soil surface associated with the plant transpiration (set to 0.5 m in this case), T p (LT −1 ) is the potential transpiration, and b(x, z, t) is a normalized distribution of maximum water uptake over a soil volume of arbitrary shape (L −2 ). The normalization of b(x, z, t) ensures that it integrates to unity over the root zone Ω R , that is, The actual water uptake distribution is defined in each point of the domain introducing Eq. (6) in Eq. (5), as follows: It has to be noticed that Hydrus version 5.01, employed to simulate the test programme, does not have integrated a root growth module able to model the root development in time (i.e.b(t) changes). Despite that, given the relatively short duration of the experiment which did not lead to an evolution in root bulb configuration (as discussed in "Experimental evidences and remarks"), the b(x, z, t) function has been simplified as time independent (i. e. b(x, z)).
Under the hypothesis of axial symmetry of the plant RWU distribution, Hydrus 2D has implemented the Vrugt et al. (2001) spatial distribution function , b(x, z), defined by the x and z variables (radial and vertical distances from the tree, Fig. 15 Graphical representation of the stress response function α(h) of Feddes et al. (1978) with indication of its parameters (h 1 , h 2 , h 3low , h 3high , h 4 , T p1 , T p2 ) Table 5 The table reports the parameters of the stress response function α(h) of Feddes et al. (1978) adopted in the numerical model where Z m is the maximum rooting depth (L), X m is the maximum rooting length (L) in the radial direction; p z and p r are empirical parameters ( −);z * and x * are other empirical parameters (L). Parameters p z and p x are set to unity for values of z > z * and x > x * , respectively. Globally, the Vrugt et al. (2001) model requires the definition of six parameters ( Z m , X m , p z , p x , z * , x * ) , of which only Z m and X m have physical meaning. Integrating Eq. (6) over the whole root zone Ω R and considering Eq. (7), T p can be expressed as follows: Similarly, the cumulated actual RWU (L T −1 ), which can be assumed corresponding to the actual plant transpiration T a , can be defined by integrating Eq. (8) over the root zone Ω R :

Soil hydraulic behaviour
The Richard equation requires as mandatory input the definition of the constitutive relations between water content, pressure head and hydraulic conductivity (i.e., K(h) and (h) ) of the soils involved in the flow simulation.
The open-source computer program RETC (RETention Curve), developed by van Genuchten et al. (1991), has been used for analysing the unsaturated hydraulic and retention properties of the soils starting from a pool of observed water retention and hydraulic conductivity data and an initial set of parameters (initial guess). The code iteratively changes the parameters of the analytical function adopted to fit the observed data, trying to minimize the residuals between observed-simulated points. Durner (1994) parametric model has been chosen to represent the SWRC and fit the laboratory retention data (gravimetric water content vs suction) of the dual porosity-improved soil. The air-entry values obtained from laboratory tests have been considered fixed over the iterative process. The model fits van Genuchten's (1980) (VG) retention curves to each pore system (micropores and macropores) and superimposes them. The inflection point of the SWRC is taken as a separation between microscopic and macroscopic porosity regions. The model hence possesses two independent sets of retention parameters. The liquid phase is subdivided between macropores (inter-aggregate) and micropores (intra-aggregates) by their respective volume fraction of bulk soil ω as follows: where the subscript i = M stands for the macropores and i = m for the micropores. m (h) and M (h) have the exact formulation of the van Genuchten's SWRC function, which is reported below (van Genuchten 1980): where S e is the effective saturation degree; θ, θ s , θ r are the current, maximum and residual water content (L 3 L −3 ), respectively; α is a scale parameter equal to the reciprocal of the air entry pressure (L −1 ) and n and m are dimensionless fitting parameters ( −). The parameter n is related to the pore-size distribution of the soil while m to the overall symmetry of the SWRC.
Moreover, RETC program allows estimating the hydraulic conductivity function (HCF) using the theoretical poresize distribution of Mualem (1976), reported in Eq. (14), starting from a fixed conductivity value (in this case, the saturated hydraulic conductivity obtained from laboratory tests): where K s is the saturated hydraulic conductivity (LT −1 ) and l ( −) is a pore-connectivity parameter of the HCF that is set to 0.5 for a wide range of soils (Mualem 1976). Figure 16 A and C represent the SWRC and the HCF of the improved soil as obtained from experimental points fitting using the RETC code. The same program has been applied again to estimate the parameters of the SWRC and the HCF of the potting mix using the bimodal Durner's (1994) model (see Fig. 16B and D). In Table 6, a summary of all the hydraulic and retention parameters used in the numerical analyses is presented. The clean gravel located at the bottom of the physical model remains saturated through the experiment. Since it has no influence on water movement around the root bulb, an arbitrary saturated permeability value of 10 −4 m/s, consistent with clean gravels, has been assumed.

Geometry of the model
A 2D section of the experimental setup has been reproduced in Hydrus 2D. It is constituted by three main zones: the gravelly layer at the bottom of the container, the soil volume with the plant in the centre and the control zone on the right hand-side separated by an impervious wall (Fig. 17). The dimensions of the root bulb are equal to the ones of the pot in which the plant has been grown in the nursery. Observation points are positioned where sensors have been installed in the physical model.  Durner's (1994) model retention curve superimposed on the experimental points (black dots) of the improved soil (A) and of the potting mix (B). The hydraulic conductivity functions of the improved soil and of the potting mix are shown in C and D, respectively. "HCF" and "SWRC" in figure stand for "hydraulic conductivity function" and "soil water retention curve", respectively At each time step of the simulation, evaporation fluxes are specified in correspondence of the atmospheric boundary condition, in the upper side of the apparatus. The remaining sides are characterized by no-flux boundary (no exchange of vapour or liquid phase) apart from the bottom boundary where a variable pressure head is applied to simulate the imposed water table (Fig. 17).

Initial conditions and simulation program
The simulation starts at the beginning of phase 1A, after the third irrigation along the plant stem and one instant before the rising of the water table from 135 to 200 mm from the base of the apparatus (Table 4). Phase 0 (0-30 days) has not been simulated because RH at the atmospheric boundary was not controlled during this initial period.
The calibration of the potential RWU distribution b(x,z) uses suction data recorded in the period between day 31 (phase 1A, Table 4) and day 91 (phase 1B). Over this time interval, the Feddes et al. (1978) stress response function parameters (presented in Table 5 and defined in "Actual evaporation and actual transpiration") have a minor influence on model output, as discussed in the following "Sensitivity analysis of stress response function parameters and identification of the calibration period of the RWU distribution"; thus, suction variations in proximity of the root bulb are governed mainly by the RWU distribution.  Durner's (1994) model for the improved soil and for the potting mix. The symbol used for each parameter and its meaning are reported in the text The subscript 1 in the model parameters highlights their belonging to the retention curve in the macropore region while the subscript 2 in the micropore region. w 1 and w 2 are dimensionless weighting factors of the two retention curves and their summation is equal to unity  The validation phase uses observation data from day 92 (phase 1B) to day 231 (phase 3), when the rising of the water table causes a rapid fall of the recorded suction values. Thereafter, as seen in "Phases 3 and 4", the high water table hides the transpiration contribution of the root bulb, and the process recorded by sensors is merely a slow equilibration of the hydraulic conditions in the apparatus.
Initial conditions have been accordingly defined, in terms of pwp distribution, starting from the available suction data in the observation points which have been interpolated linearly in order to extrapolate values in the whole mesh of calculation. This allows having a good match between observed and simulated pore pressure data at the initial instant of the simulation.

Sensitivity analysis of stress response function parameters and identification of the calibration period of the RWU distribution
A sensitivity analysis has been performed on Feddes et al.'s (1978) stress response function parameters to assess their incidence on the model output over the simulation period. The simple and widely used "one factor at a time" OFAT technique (Simunek et al. 1998;Abbasi 2008) has been applied: one model parameter at a time (p j ) is varied of ± 1% its initial value, and the sensitivity coefficient, S(p j ) (Eq. 15), is determined at each instant of the simulation, to measure the local response of the model to this change: where Y p j + Δp j e j represents the model variable Y (e.g. pressure head or water content) corresponding to a 1% change in the parameter p j (Δp j = 0.01 p j ), and e j is the jth unit vector; Y(p j ) represents the value of the variable Y using the unvaried initial value of the parameter p j .
As expected, the major sensitivity coefficients are observed for the observation points located on the side of the bulb (compare no. 5 and no. 10 in Fig. 18) and closer to the bulb centre (compare no. 5 and no. 4 in Fig. 18). In the interval 31-91 days, for the same observation points, the incidence of the variations of Feddes et al.'s (1978) parameters is negligible while in correspondence of the remaining points (no. 3 to no. 1 and no. 8 to no. 6), the incidence is null over the whole simulation period.
The sensitivity analysis shows that the calibration of the b(x,z) distribution should be performed between days 31 and 91 because, in this time interval, observation data are mainly influenced by b(x,z) and only marginally by Feddes et al.'s (1978) parameters. In fact, when water content availability in proximity of the root bulb is high, plant transpiration does not experience reduction in RWU due to suction (values in between h 2 and h 3 in Fig. 15) so small variation in the parameters defining the boundaries of the response function has no effect on model results.

Calibration of the RWU spatial distribution
To correctly represent the 2D pwp distribution in the proximity of a stand-alone plant, the root distribution parameter, b(x,z), should be defined in each point of the domain. This has been achieved by: -Adopting the flexible Vrugt et al.'s (2001) root distribution function (Eq. 8) in two variables and six parameters -Adopting a user-defined root distribution function, a 2nd degree quadratic function in two variables and four parameters In both cases, the calibration procedure involves the iterative variation of the parameters of the RWU distribution functions within reasonable intervals of values and, at each iteration and in each observation point, the graphical comparison of the predicted distribution of suction to the experimental values, through a "trial and error" procedure.
The calibration of both the RWU spatial distribution functions relies on the following hypotheses and observations: 1. The compensatory uptake mechanism (Jarvis 1989), i.e. the ability of the plant to cope with soil water content distribution using a "compensatory" uptake from the wetter soils around the bulb, has been neglected in the numerical analyses. Many studies (Simunek and Hopmans 2009;Albasha et al. 2015 among others) have demonstrated that compensatory RWU is not needed when an adequate description of spatial root distribution is provided. This observation also combines with the fact that the effects of spatial root distribution may reduce if a compensatory mechanism is considered (Albasha et al. 2015), which is clearly to be avoided in a procedure of calibration of a RWU distribution. 2. As explained in "Initial conditions and simulation program", in order to isolate the effects of the RWU spatial distribution from the influence of the water stress on the uptake fluxes in the proximity of the root bulb, the calibration of b(x,z) has been carried out using observed data recorded between days 31 and 90. 3. In agreement with the observations made in "Experimental evidences and remarks", the region of influence of the RWU stretches in the horizontal and vertical direction to a distance between 50 and 80 mm and between 80 and 130 mm from the root bulb sides, respectively. From this, observation points nos. 4, 5, 9 and 10 are the only ones affected by RWU, with a b(x,z) different from zero.
The calibration is concluded when the optimal local root distribution (i.e. the b(x,z) distribution that minimizes the discrepancies between observed and simulated data) is identified. It is relevant to stress that, in both cases, the calibrated spatial distributions of b(x,z) cover a wider soil region with respect to the initial one (the pot in which the plant grew in the nursery). The use of b(x,z) is a strategy adopted to simulate the influence of the root system on soil moisture content keeping only partly a physical meaning (correlated generally to root density distribution). Many authors (Bruckler et al. 2004;Hodge 2004;Faria et al. 2010, among others) have demonstrated that the spatial distribution of RWU may differ greatly from root density, and discrepancies are expected to rise with time and in heterogeneous soil (Kuhlmann et al. 2012;Schneider et al. 2010). Vrugt et al.'s (2001) model, discussed in "Water movement and RWU " (Eq. 9), is defined by 6 parameters. Parameters X m and Z m control the extension of the area influenced by RWU. Their calibration has been performed by iteratively changing their values so that the influencing area extends between 51 and 80 mm and 81 and 130 mm from the bulb , on the side of the root bulb, and no. 9 (C) and no. 10 (D), below the root bulb, as obtained applying the one factor a time analysis (OFAT) to the Feddes et al.'s (1978) stress response function parameters (h 1 , h 2 , h 3low , h 3high , h 4 ). In each graph, the sensitivity coefficients vs time obtained varying of ± 1% each parameter value are shown sides in the horizontal and vertical direction, respectively (as discussed in point (3), "Calibration of the RWU spatial distribution"). Parameters z * and x * somewhat control the vertical and radial distance of the maximum RWU intensity with respect to the plant (origin) while p z and p x are responsible for the rate of decay of the RWU from the zone of maximum intensity outward. p x has been varied iteratively between 0 and 0.7, p z between 0 and 3 while z * and x * between 0 and 0.5. Observation point nos. 9 and 10 (below root bulb) have been used for the calibration of p z and z * while no. 5 and 4 (on the side of the bulb) for the calibration of p x and x * .

Calibration of Vrugt et al.'s (2001) root spatial distribution function
The calibrated parameters are presented in Table 7 while the obtained RWU distribution in Fig. 19. Vrugt et al.'s (2001) RWU function is potentially applicable to all existing root bulb configurations, and from this comes its complex formulation and the high number of parameters required for its definition. Other authors preferred, to a flexible but complex RWU spatial formulation, a simpler one Table 7 The table reports the calibrated parameters of the spatial distribution function of Vrugt et al. (2001) and of the user-defined one, for the presented case study Parameters for vertical distribution (Vrugt et al. 2001 with a lower number of parameters but able to describe only a specific root bulb configuration. This is the case of Coelho and Or (1996) who proposed three different 2D uptake patterns associated to three root bulb distributions resulting from different drip irrigation configurations. Most of the spatial functions found in literature have been elaborated in the agricultural fields for an efficient management of crops, and the bulb development considered results very different from that of high-stem plants, employed in earthen slope stabilization. In the laboratory apparatus, the plant is fed from the imposed water table as occurs in nature, thus resulting in a different RWU configuration with respect to one of irrigated cultures. All this spawned the idea to find the best RWU spatial distribution for this peculiar bulb configuration. The user-defined RWU distribution function relies on two hypothesis:

Calibration of the user-defined RWU spatial distribution function
1. Axial symmetry of the root distribution 2. Consistently with physical observations, the maximum value of b(x,z) is located in the centre of the root bulb, and it decreases rapidly moving away from it (Knight 1999;Doussan et al. 2006) The user-defined RWU distribution function is quadratic in two variables (x and z) and defined by three parameters ( , , ): where x and z are the radial distance and the depth from the origin of the tree, respectively.
According to hypothesis (2), the maximum RWU, b max , is imposed in the centre of the bulb: In the vertical and horizontal plane containing the origin, a hypothetical depth z 1 and radial distance x 1 in which the RWU is zero, are assumed, respectively. These assumptions provide the definition of and : From which, it follows: Since Hydrus program requires a dimensionless distribution of RWU , both b(x, z) and its variables (x and z) have to be normalized: where z m is the maximum rooting depth (L), x m is the maximum rooting length in the radial direction (L), x * and z * are the dimensionless variables ( −) and b * (x * , z * ) the normalized RWU distribution ( −).
By replacing expressions (22), (23) and (24) in (16), it follows: where Equation (26) is the dimensionless RWU spatial distribution function in 4 parameters ( x m , z m , , ) that has been used for the numerical analysis. In the calibration process, the parameters have been iteratively changed and, each time, the b(x,z) value is calculated according to Eq. (26) in all the nodes of the mesh and subsequently imported in the program in the form of triplets (x*,z*,b(x*,z*)) where x* and z* are the coordinate of a mesh node. The imported b(x,z) distribution is not normalized, i.e. it does not integrate to unity over the root zone (Eq. (7)), since this operation is performed by Hydrus internally at the start of the computation. The calibration of x m and z m has been performed similarly to X m and Z m in Vrugt et al.'s (2001) function since they have the same physical meaning (i.e. bounding of the area affected by root uptake). Parameters and are responsible for the rate of decay of b(x,z) from the bulb centre outward in the radial and vertical direction, respectively, which results extremely suited for our problem since observation nodes are positioned horizontally and vertically from the bulb. Therefore, the calibration of can rely on suction data of horizontal sensors (no. 5 and 4) while of the vertical ones (nos. 10 and 9). Both the parameters have been iteratively varied during the calibration phase in the range 0-1, and their optimized values are reported in Table 7. Figure 20 A shows the calibrated user-defined RWU spatial distribution, Fig. 20B the trend of b(x,z) in the horizontal and vertical direction from the bulb centre.

Comparison of calibrations and validation of results
In Fig. 21, the simulation output in correspondence of the observation nodes obtained by applying Vrugt et al.'s (2001) root spatial distribution is compared to the modeller-defined one in the period between 31 and 91 days (calibration phase) and 92 and 231 days (validation phase). Observed data (from installed sensors) are also shown as the reference. Both the simulations show good agreement with the observed data in the calibration and validation phases, even if the analysis characterized by the user-defined distribution locally shows a better performance (see sensors no. 5 and no. 4 in Fig. 21). Only observation points influenced by the RWU distribution (i.e. sensors no. 5, 4, 10 and 9) show differences in model results. Moreover, the greater the potential transpiration (i.e. the b(x,z) value), the higher the potential cumulated difference in suction to 231 days. Sensor no. 5, 30 mm from the bulb, shows ~ 40 kPa of difference while sensor no. 4 only ~ 5 kPa at day 231 between the two models. Sensors no. 10 and no. 9 do not show relevant differences in suction values between the two models, symptoms of an equally good representation of the RWU distribution in the vertical direction. Hence, both RWU distributions have proved to be suited for the simulation of the large-scale experiment even if the use of a less complex RWU spatial function identified by a lower number of parameters considerably facilitates the procedure, especially for a case study as the presented one, in which the optimization procedure has been performed without the use of an optimization algorithm.

Lessons learned from the large-scale apparatus experiment and the calibration procedure
The proposed calibration of the RWU spatial distribution in proximity of a high-stem plant has highlighted three important elements unmissable to guarantee good results: 1. A pressure head distribution in proximity of the root bulb able to guarantee no reduction of the potential RWU. This means that the stress response function α(h) needs to be equal to 1 (i.e. potential and actual transpiration coincides) during the calibration period. To do so, the water table height imposed to the model needs to be chosen carefully by means of a preventive numerical modelling of the experiment. 2. Installation of local suction monitoring sensors in proximity of the root bulb. As clearly seen from the experimental data, root bulbs of young plants have reduced spatial effects on pwp distribution. For this reason, it is useful to place the majority of the sensors in the zone influenced by the water uptake. The majority, but not all of them: sensors outside the "RWU influence zone" help define the boundary of the soil volume interested by the root bulb presence. It is suggested to place 3 sensors as close as possible to the root bulb in the horizontal, radial and vertical direction: 2 sensors that fall within the RWU influence area and 1 outside. In the presented Fig. 20 A Calibrated RWU spatial distribution in the investigated domain according to the user-defined function. B Trend of the root distribution parameter b(x,z) with vertical and radial distance from the root bulb centre with indication of its value in the observation nodes study, a distance of 30 mm was left between the root bulb perimeter and the first sensor to avoid the influence of the interface between two different materials (potting mix and improved soil) on suction records. The remaining sensors should be located at a distance that does not result an obstacle to the water movement around the plant (in the presented case study, a spacing of 50 mm was applied). 3. A calibration period of 60 days has been used in the presented study, and it has proven to be sufficient to obtain good results in terms of simulation of the suction distribution in proximity of the plant. Nonetheless, it is important to remember that a longer period of data acquisition is required to validate the results. In the presented case, 139 days of recorded suction data have been used to validate the RWU spatial distribution. During validation period, it is recommended to change boundary conditions to the large-scale apparatus (i.e. water table height, RH) in order to test the numerical model in conditions different from the ones used in the calibration phase.
To perform the calibration, 7 sensors have been used: nos. 5, 4, 3 (on the side of the bulb), no. 11 (at a radial distance) and nos. 10, 9, 8 (below the bulb). Sensors nos. 3, 11, 8 were used to mark the boundaries in which the root uptake is null (i.e. b(x,z) equal to 0) while sensors nos. 5, 4, 10, 9, interested by root uptake, were used to define the distribution of the b(x,z).
In light of point (2), to improve the proposed RWU spatial distribution for the Melaleuca styphelioides, two more sensors can be positioned at a radial distance but closer to the root bulb to increment the information on the radial root uptake.
Moreover, benefits would also result from the pairing of dielectric sensors (monitoring suction distribution) to soil moisture sensors: the two instruments together can give precious information on the branch of the soil water retention curve in which that point of the domain is moving.

Conclusions
The spatial pattern of RWU is function of the root architecture, and its knowledge is essential in the prediction of the spatial distribution of water contents and water fluxes in soils.
The experiment conducted in the large-scale apparatus on a 2-year-old Melaleuca styphelioides plant allowed to collect suction data in proximity of the bulb. These data Fig. 21 Graphical comparison of the model output in the observation nodes no. 3, 4, 5 (on the side of the root bulb) represented in graph (A) and nos. 9, 10 (below the root bulb) and no. 11 (radial distance from the bulb) in graph (B) using the modeller-defined RWU spatial distribution function (marked with coloured crosses) and the Vrugt et al.'s (2001) function (marked with coloured triangles) vs the observed data (marked with coloured dots) between 31 and 91 days (calibration period) and 92 and 231 days (validation period) guided the calibration of two RWU distributions, one using the spatial function of Vrugt et al. (2001), and the other is a modeller-defined function. The latter is a quadratic function in four parameters and two variables which resulted particularly suited to represent a RWU distribution resulting by a high-stem plant fed by a shallow water table. Both the RWU functions ensured a good agreement between experimental and simulated suction data although the user-defined one, given the lower number of parameters and the simpler formulation, allowed an easier implementation.
The main benefits of the experiment are the provision of controlled boundary conditions which facilitates the interpretation of the experimental data and allows a more accurate numerical simulation of the problem.
The use of dielectric sensors in the large-scale apparatus for the calibration of the RWU spatial distribution proved to be a valid alternative to destructive techniques applied to obtain information with depth on the root length density. Moreover, since root architecture changes in time for natural plant growth and external stimuli, the RWU spatial distribution needs to be updated periodically performing each time new root samplings. The use of these sensors both in the field and in laboratory tests can potentially provide the data required to update roots architecture in time with no additional costs.
The experiment can be potentially replicated with the same equipment using other plant species in order to individuate the best long-stem plant solution for the revegetation of the investigated unstable slope. Moreover, the equipment can be easily replicated in other laboratories because only commercial sensors (as TEROS 21 sensors by Meter group) and commercial devices (as the plant growing lamp) have been used, and the remaining components of the apparatus are equally easy to purchase and to put into use. Additionally, the chosen water potential sensors do not need a user calibration; they are relatively cheap and easy to install. All this, potentially, allows to compare the results obtained with different plant species in different laboratories and to create a database.
The experiment conducted on the Melaleuca styphelioides over a relatively short period (< 1 year) could be extended in time in order to investigate the development of the underground architecture. For long-stem plants that can reach considerable height, as the one in use, the laboratory-scale dimension is a viable option only up to the plant maturity, when the testing has to move necessarily to the field scale.
All things considered, the large-scale apparatus setup and the calibration strategy of the RWU spatial distribution have confirmed its applicability in those situations where a correct simulation of the transpiration contribution of a plant is required. Among them, slope stability assessment needs a high level of precision in the definition of pwp distribution in time, and this implies not only a correct simulation of the transient seepage in unsaturated conditions, but also of evaporation and transpiration fluxes. The complexity of the analysis can be overcome by breaking the problem down into parts that can be tackled more easily. Therefore, the presented work sets to a former level in which the focus is the plants, its root distribution and the suction regime in the soil volume interested by RWU. The collected pieces of information will contribute to a better understanding of the role of vegetation on the stability of shallow soil under unsaturated transient regime.