Infiltration-soil moisture redistribution under natural conditions: experimental evidence as a guideline for realizing simulation models

The evolution in time,t , of the experimental soil moisture vertical profile under natural conditions is investigated in order to address the corresponding simulation modelling. The measurements were conducted in a plot with a bare silty loam soil. The soil water content, θ , was continuously monitored at different depths, z, using a Time Domain Reflectometry (TDR) system. Four buriable three-rod waveguides were inserted horizontally at different depths (5, 15, 25 and 35 cm). In addition, we used sensors of air temperature and relative humidity, wind speed, solar radiation, evaporation and rain as supports for the application of selected simulation models, as well as for the detection of elements leading to their improvement. The results indicate that, under natural conditions, very different trends of the θ (z, t) function can be observed in the given fine-textured soil, where the formation of a sealing layer over the parent soil requires an adjustment of the simulation modelling commonly used for hydrological applications. In particular, because of the considerable variations in the shape of the moisture content vertical profile as a function of time, a generalization of the existing models should incorporate a first approximation of the variability in time of the saturated hydraulic conductivity,K1s, of the uppermost soil. This conclusion is supported by the fact that the observed shape of θ (z, t) can be appropriately reproduced by adopting the proposed approach with K1s kept constant during each rainfall event but considered variable from event to event, however the observed rainfall rate and the occurrence of freeze-thaw cycles with high soil moisture contents have to be explicitly incorporated in a functional form for K1s(t). Correspondence to: R. Morbidelli (renato@unipg.it)


Introduction
Simulations of successive cycles of rainfall infiltration, redistribution of soil water in the no-rainfall periods and reinfiltration are needed for many hydrological applications.Usually, a simple and sufficiently accurate approach that represents on a continuous basis the evolution in time, t, of the soil moisture, θ, at different depths, z, is required.Many experiments were set up in order to show the shapes of θ (z, t) in bare soils, that were then used as guidelines for realizing the required theoretical approach.Biswas et al. (1966) examined the function θ (z, t) within redistribution periods without evaporation by laboratory experiments for different soil types as sand, silt loam and clay loam and different cumulative infiltration depth at the beginning of the process.Similar laboratory experiments were described by Young and Poulovassilis (1976) for a sandy soil.Gardner et al. (1970a,b) realized the same type of investigations on fine sandy loam soils and determined θ(z,t) during the stage of redistribution both with and without evaporation.For each of these investigations, an approximated simulation model for the reproduction of the profiles observed during redistribution was proposed (Gardner et al., 1970b).
Theoretical formulations for a continuous description of both soil moisture profiles and infiltration were proposed by several authors (see, for example, Milly, 1986;Charboneau and Asgian, 1991;Corradini et al., 1997Corradini et al., , 2000)).Milly (1986) and Charboneau and Asgian (1991) include explicitly the evapotranspiration process and represent the infiltration through the classical techniques, specifically the extended Philip equation (Chow et al., 1988) combined with the "time compression approximation" (Sivapalan and Milly, 1989) and the Mein and Larson equation (Mein and Larson, 1973), respectively.Corradini et al. (1997Corradini et al. ( , 2000)), in principle, provide a more appropriate representation of both θ (z, t) and infiltration rate through the development of composite θ profiles, even though evapotranspiration is neglected and the formulations were tested in the limits of hydrological applications involving redistribution periods up to 20 h.In any case this duration is typically acceptable when floods due to frontal rainfalls have to be considered.
An overall analysis of our discussion indicates that sufficiently simple models of θ(z, t) are available, but their effectiveness was typically examined for a limited application field and mainly for laboratory experiments where vertical homogeneous soils were generally involved.Therefore, the existing models need to be validated and generalized on the basis of experimental data carried out in natural fields for long periods and considering that sealed soils can be found (Bullock et al., 1988;Emmerich, 2003).
The main objective of this paper is to address the above issues through continuous measurements of the basic quantities determining the evolution of θ(z, t) and direct observations of the same function performed by a Time Domain Reflectometry (TDR) system.In this context, in order to quantify the possible errors generated by an inappropriate choice of the modelling, the effectiveness of a model based on the assumption of a vertically homogeneous soil (Corradini et al., 1997) with time invariant saturated hydraulic conductivity is compared to that of a model for a two-layered soil (Corradini et al., 2000).The latter is applied considering a sealed soil with the upper layer characterized by a saturated hydraulic conductivity, K 1s , time-varying as a stepwise function and the underlying layer that keeps the properties of the parent soil.The quantity K 1s is determined by calibration as a constant within each specific study period (event), while the associated variability from event to event is investigated by examining the link with the experimental values of the hydrometeorological variables observed in a selected study plot of Central Italy.This can be considered a first approximation since K 1s is really time-dependent also during a specific event (see also Assouline and Mualem, 2000) when we assume, for the sake of simplicity, a constant value.On this basis the main lines to follow for a further development of the pre-existing models are also given.Corradini et al. (1997Corradini et al. ( , 2000) ) models extended for evapotranspiration

A short description of the
The soil moisture vertical profile evolves under natural conditions because of the processes of rainfall infiltration and redistribution of soil water with evapotranspiration.These processes interact in determining the profile that can be investigated, as a first approximation, using one of the models earlier proposed for vertically homogeneous or layered soils.We have selected the models developed by Corradini et al. (1997Corradini et al. ( , 2000) ) adapted by the additional effect of evapotranspiration.
The model by Corradini et al. (2000) considers the onedimensional water flow into any two-layered soil profile under complex rainfall patterns, which produce successive infiltration-redistribution cycles.The problem formulation is simplified by assuming the initial condition for capillary head, ψ, invariant with depth and approximating in each layer the dynamic wetting profile at a given time by a distorted rectangle represented through a shape factor β (≤1) which depends on the water content at the top of the layer.
In each layer the model combines the continuity equation with the depth-integrated form of the Darcy law and uses the appropriate boundary conditions including the continuity of capillary head and water flux at the interface.Since the time when the wetting front reaches the interface, t c , we have: (1) with P L (ψ c ,t) defined as: I 2 being expressed through the difference: where we use subscripts 1 for variables in the upper layer, 2 in the underlying soil, c at the interface, 0 at the surface, i for initial conditions, and s for the saturated conditions; Z c is the upper layer thickness; C(ψ) = dθ/dψ; q water flux; K hydraulic conductivity; G net capillary drive; I the cumulative infiltration depth; the second term on the right hand side of Eq. (4) the cumulative dynamic infiltration depth in the upper layer; β, p and α quantities depending on the wetting profile shape, with β explicit function of θ (Corradini et al., 1997).Furthermore, for t < t c the wetting front moves through the homogeneous upper layer where the use of a similar procedure (Corradini et al., 1997) leads to: The quantity E in Eqs.
(1) and ( 5) is the evapotranspiration given by (Ragab, 1995): where SMI is the soil moisture index expressed by: with θ fc the field capacity defined as the bulk water content in soil at -0.33 bar pressure head and θ wp the wilting point defined as the bulk water content at −15 bar pressure head; E p is the potential evapotranspiration, given by (Maidment, 1993): or, alternatively, by the Penman-Monteith equation: with k a empirical pan coefficient; E v pan evaporation; υ latent heat of vaporization of water; gradient of e s as a function of temperature; R n net incoming radiation flux; F soil outgoing heat flux; ρ and c p density and specific heat at constant pressure of air, respectively; e s and e saturation and actual vapour pressure, respectively; r a and r st aerodynamic and stomatal resistance, respectively; γ psychrometric constant.
By adopting for θ(ψ) and K(ψ) the general functional forms used, for example, by Smith et al. (1993) and Smith (2002), numerical solutions of Eqs. ( 1)-( 2), for sealed soils, or Eq. ( 5), for homogeneous soils, provide ψ 10 (t) and ψ c (t) or ψ 10 (t), respectively.In addition, through the upper boundary conditions the cumulative dynamic infiltration depth as a function of time can be determined.
The vertical profile θ(z) at each time is represented by using the function θ (ψ) to derive the values of soil moisture at specific depths.The shape of each profile is determined by an extended form of that utilized by Corradini et al. (1997) for homogeneous soils.
An example of the results generated by Eqs.(1)-( 2) and ( 5) is given in Fig. 1a-c where, under a steady rainfall rate followed by a redistribution period, the time-varying infiltration rate and the corresponding soil moisture profiles at a few specific times are shown for both homogeneous and sealed soils.

Experimental system
We selected a natural plot characterized by slight slope with a silty loam soil (see Table 1).The application of the explanatory models described in Sect.2.1 requires measurements of a few classical hydro-meteorological quantities to be used specifically for the simulation of the actual soil moisture profiles.
Sensors of air temperature, wind speed, relative humidity and solar radiation, all operative at height of 2 m, provided the basic data for the estimate of E p ; in addition a Class A evaporation pan was used for E v of Eq. (8).A raingauge was set up in order to observe the rainfall rate that is equal to q 10 of Eqs. ( 1) and ( 5) during periods with unsaturated soil surface.In addition a TDR system was installed for the direct measurements of θ (z, t).
The water content at specific depths was monitored using the TRASE-BE device designed by the Soil Moisture Equipment Corp., Goleta, CA.Four buriable three-rod waveguides with length 20 cm were inserted horizontally at different depths along the same vertical, taking care that the metal rods were in tight contact with the soil.Specifically, the probes were installed at 5, 15, 25, 35 cm depth.Each probe provided at each time a measurement of θ at the corresponding depth, but in any case the zone of influence is a small soil volume around the same probe.The reliability of the system is supported by Zegelin et al. (1989).They showed that three-wire probes imbedded horizontally in a field profile during a rainfall event can be used to follow the wetting front during both infiltration and redistribution; in addition,  Smith et al., 1993) these probes can be used to estimate the amount of infiltrated water with error less than 10 %.The four probes were interrogated at intervals of 30 min and the universal calibration curve of Topp et al. (1980) (relationship between volumetric water content and apparent dielectric constant) was used to calculate θ from the TDR signal.

Experimental results and their analysis
Considering that small simulation errors associated to short periods could increase with time and mask the mechanisms to be clarified in this study and that the models selected were not tested for very long periods, we have chosen to investigate periods of limited duration.Therefore, many isolated events involving an infiltration period with successive redistribution/evapotranspiration of soil moisture have been investigated.Each event started from a condition of initial soil moisture that, using the measurements available, in the soil part between 5 and 35 cm has been approximated as invariant with depth.In the hypothesis of vertically homogeneous soil with constant ψ i the shape of the curve θ i (z) in the top soil, that is between the soil surface and 5 cm depth, has been represented by the same value of θ i used at larger depths, while under the condition of sealed soil by the constant value giving continuity of ψ i at the interface.
The basic quantity to be assessed in order to perform the conceptual simulations for a sealed vertical profile is the depth of the crust layer, Z c , that along the lines discussed by Mualem and Assouline (1989) and Mualem et al. (1993) has been set equal to 5 cm.In any case, a sensitivity analysis of results to the sealing layer thickness has been also carried out.The application of the models described in Sect.2.1 requires to estimate the quantities β, p and α, linked with the shape of the wetting profiles (see Corradini et al., 1997) and those involved in the functional forms θ (ψ) and K(ψ), as the Brooks-Corey pore size distribution index, λ, the air entry head, ψ b , the empirical coefficients c and d, the saturated and residual soil water content.The last two values were associated with the maximum (0.34) and the minimum (0.092) values of θ obtained during the entire period of observation.On this basis and considering the soil structure, following Corradini et al. (1997Corradini et al. ( , 2000) ) the remaining quantities were assessed (see Table 1).All the above quantities were assumed to be invariant from homogeneous to sealed soils.This implies that the continuity of ψ i leads to a constant value of θ i through the entire soil profile.Furthermore, in the application of the models the potential evapotranspiration was estimated by Eq. ( 8) with k a depending on air relative humidity, wind speed and soil cover as given by Maidment (1993).
All the significant events in the period April 2010-December 2010 have been investigated (see Fig. 2).For each specific event the observed profiles of θ have been compared with those simulated by both the homogeneous and the twolayered model.The same value of saturated hydraulic conductivity was used in the computations for the homogeneous soil and the subsoil of the sealed profile.In particular, this value was assumed time invariant and equal to that derived by calibration events occurred in the period October 2009-March 2010 (see Table 2), all associated with conditions of vertically homogeneous soil.Furthermore, the value of K 1s in the two-layered soil profile was estimated by calibration in each specific event.
The results obtained for a few events representative of the evolution in time of the processes of infiltration, redistribution/evapotranspiration during a wide study period are shown in Figs.3-9 (see also Fig. 2 and Table 2).These events were considered representative because they enable us to provide a complete interpretation of the processes involved with the time evolution of the soil moisture profile.As it can be seen  2).Saturated hydraulic conductivity 2 mm h −1 .
Table 2. Main characteristics of a few simulated events using for calibrating the saturated hydraulic conductivity of the parent soil (events 1c-4c) and for the successive representative applications of the study models (events 1-7 in Fig. 3 the curves estimated by the model for vertically homogeneous soil, with θ(z) at the end of the rainfall as well as during the redistribution/evapotranspiration period, even though in the initial infiltration stage the agreement between theoretical simulations and experimental results does not appear fairly good probably because of inaccuracies in the TDR measurements (see also Melone et al., 2006Melone et al., , 2008)).In addition, we found that the application of the two-layered model  2).Saturated hydraulic conductivity of the parent soil 2 mm h −1 and seal layer 0.07 mm h −1 .
to the same event provided representations of θ (z, t) inconsistent with the observations.On the other hand, Figs.4-6 show that the profiles simulated by the two-layered model, with the value of the saturated hydraulic conductivity in the topsoil obtained in each event by calibration and in the parent subsoil assessed in advance (K 2s = 2 mm h −1 ), fit fairly well the experimental data, even though in the limit of the lack of measurements at a depth less than 5 cm.These figures indicate also that the curves computed by the homogeneous model using K 1s = 2 mm h −1 have inadequate shapes, associated with substantial overestimates of the wetting front depth and cumulative dynamic infiltration.Furthermore, simulations carried out by the homogeneous model but with K 1s assessed by direct calibration were found in any case to provide results conflicting with the observed data.
An event occurred in the late autumn period is investigated in Fig. 7 that displays as the observed moisture contents are well simulated by the model for homogeneous soil with K 1s = 2 mm h −1 , while the existence of a crust layer, as in the event of Fig. 3, can be excluded because of the incorrect shape of the curves obtained by the sealed soil model.
Two additional events, that contribute to improve the knowledge of the mechanisms determining the time-varying vertical soil structure, are examined in Figs. 8 and 9.In spite of the results are difficult to quantify because the two events start from high soil water contents, in any case they clearly suggest a similitude with Figs. 4 and 7 in relation to the existence of a sealed soil and a homogeneous soil, respectively.
An overall analysis of our results indicates that during the high intensity summer thunderstorms there is a crust layer with a disturbed thickness of a few centimetres (see also  2).Saturated hydraulic conductivity of the parent soil 2 mm h −1 and seal layer 0.03 mm h −1 .Mualem andAssouline, 1989 andMualem et al., 1993), that keeps within no rainfall periods and is associated with high soil temperature.It becomes gradually less permeable because of successive heavy storms, until it experiences a disruption.For example, in spite of the rainfall rate associated with Fig. 6 is halved in relation to that of Fig. 4 (see Fig. 2), the ratio K 2s /K 1s has more than doubled.2).Saturated hydraulic conductivity 5 of the parent soil 2 mmh -1 and seal layer 0.12 mmh -1 .6  2).Saturated hydraulic conductivity of the parent soil 2 mm h −1 and seal layer 0.12 mm h −1 .Furthermore, in the late autumn the disruption of the crust layer with re-formation of a vertically homogeneous soil occurs.Then, during the successive months the infiltrationredistribution/evapotranspiration events can be associated to conditions of both sealed soil and homogeneous soil.
These results can be explained in the light of the coupling, proposed by Bullock et al. (1988)   experiments, between the disruption of the crust layer and occurrence of freeze-thaw cycles under the condition of high soil water content.In fact, as it can be seen in Fig. 10, the event of 24 December 2010 (Fig. 9) was preceded by many cycles of this type with air temperature, at 2 m above the ground, that dropped a few degrees below 0 • C while the soil water content at the 5 cm depth assumed values very close to soil saturation.A similar analysis can be performed for the event of 16 November 2010 (Fig. 7) when the duration of the antecedent freeze-thaw cycle was of about three days.This was a period more extended than that involving air temperature below 0 • C, because the air temperature at 2 m provides a considerable overestimate of the soil surface temperature in the nights with low wind intensity and absence of clouds.These meteorological conditions were really observed in the period aforementioned.Furthermore, in the event of Fig. 8 there was a crust with permeability larger than that of Fig. 4 because of the lower rainfall intensity involved.From the analysis of our representative events it can be deduced that the re-formation of a crust layer can occur in any month, in the absence of cycles, provided the rainfall rate is larger than a moderate threshold value.
In principle, the above results are not conflicting with the seasonal variability of the saturated hydraulic conductivity measured by Emmerich (2003) in a study area with very different climatic conditions.
Furthermore, we note that by neglecting the evapotranspiration our results experienced minor changes.The measurement errors have been also found to make an insignificant impact on our analysis.In fact the accuracy of the TDR measurements was typically within 2 % (Jones et al., 2002)  except when the zone of influence of a probe is crossed by the lower part of the wetting front.In any case the last situation interests only the results shown in Fig. 3, that have been discussed in this context.In addition the errors of the air temperature sensor, that was ±0.2 • K, cannot modify the detection of the freeze-thaw cycles we have examined in the nights with low wind intensity and absence of clouds.
Lastly, the influence of variations in sealing layer thickness on the results above discussed has been investigated considering Z c in the range 0.5-5 cm.The lower value is very close to the most common experimental results (see, for example, Sharma et al., 1981;Fohrer et al., 1999) and was frequently used in order to test a few specific models for infiltration in crusted soils (see, for example, Hillel and Gardner, 1970;Ahuja, 1983).The upper limit is linked, as pointed out by Mualem and Assouline (1989), with the extension of the region where changes in soil properties might be considered relevant, particularly under conditions of unsaturated surface.Values of the order of some centimetres were adopted by Mualem et al. (1993) and Assouline and Mualem (1997) for different soil types.In spite of the wide variability in sealing layer thickness reported in literature was explained by Assouline and Mualem (2000) through the different resolution used in measuring the soil bulk density, the uncertainty in seal thickness justifies the necessity to quantify its role in this study.For each representative event of Table 2, using the same procedure adopted with Z c = 5 cm we have simulated the soil moisture vertical profile both with Z c = 0.5 and 2 cm.Even though within the limits determined by the lack of measurements of water content above Z c = 5 cm, the observed shape of θ (z, t) appears to be reproduced fairly well independently of Z c but with the optimal value of K 1s that changes considerably with the seal thickness.In particular, our results, synthesized in Table 3, indicate that for a specific event the unit area hydraulic resistance of the sealing layer, defined as Z c /K 1s , remains unchanged while K 1s increases by a factor 10 from Z c = 5 cm to Z c = 0.5 cm.This trend of K 1s is determined by the fact that a given observed   cumulative infiltration has to be reproduced by a more compacted seal when a thinner thickness is involved.In any case, to some extent the experimental moisture content at the depth of 5 cm appears to be simulated better by adopting Z c = 5 cm.

Conclusions
Through the use of the conceptual models developed by Corradini et al. (1997Corradini et al. ( , 2000)), that were tested for "isolated" events by using the Richards equation as a benchmark, reliable simulations of the main characteristics of the natural soil moisture vertical profiles in bare soils are presented.Actually, this modelling has been adapted by incorporating the process of evapotranspiration, that in the "isolated" events examined here has a minor role but becomes significant when continuous applications for practical use have to be performed.
Our results indicate that, under summer thunderstorms and even for moderate rainfall rates occurring in other periods, a vertically homogeneous soil can change into a two-layered soil with a crust in the upper part.The latter keeps within no rainfall periods and becomes gradually less permeable because of successive storms of considerable intensity.The disruption of the sealing layer appears to occur through freezethaw cycles when high soil water contents are involved.This process was earlier suggested by Bullock et al. (1988) on the basis of laboratory experiments.Our results can be only partly compared with those provided by Emmerich (2003), obtained through direct measurements of saturated hydraulic conductivity, because of the different distribution in time of the experiments analysed.In any case, in the limits specified, the results of the two investigations are not conflicting.
The methodologies used in the development of this work and the results we have obtained suggest that: (a) the simulation of the evolution in time of the moisture content vertical profile under natural conditions, to be used in surface and sub-surface hydrological applications concerning, for example, the formation of the overland flow and the study of the atmospheric processes requiring the water content at the soil surface as boundary condition, should rely on a sealed soil approach for the vertical profile; (b) a representation of the variability in time of the sealing layer saturated hydraulic conductivity is a crucial model requirement (see also Langhans et al., 2011); (c) the model should represent appropriately the function θ (z, t) in a homogeneous soil when the hydraulic characteristics of the two layers are identical.The modeling adopted here verifies the points (a) and (c) but should be completed by incorporating point (b).
In conclusion this work, putting together measurements of some hydrometeorological variables, observed moisture vertical profiles and a conceptual modelling for their simulation in "isolated" events, addresses the problem concerning the realization of a model for a continuous prediction of θ (z, t) in bare soils by identifying its basic structure and the main elements to be used in developing a functional form for K 1s .In order to obtain a solution of the last problem it is required to extend this investigation with calibration of K 1s for event based models to a much larger number of events observed in study plots with different soil characteristics.
Figure 1.(a) Time evolution of infiltration rate for homogeneous and sealed soils.(b) and (c) 3 Vertical profiles of soil moisture at the different times indicated in the upper diagram, for 4 homogeneous and sealed soils, respectively.5

Fig. 1 .
Fig. 1.(a) Time evolution of infiltration rate for homogeneous and sealed soils.(b) and (c) Vertical profiles of soil moisture at the different times indicated in the upper diagram, for homogeneous and sealed soils, respectively.

Figure 2 .Fig. 2 .
Figure 2. Rainfall rate as a function of time observed in the study plot.The events used as 13 representative cases in this study are indicated by the numbered dashed lines.14 15 16 17 18 19 20 21 22 23 24 25 26

Fig. 8 .
Fig. 8. Water content profiles simulated by the adopted models for homogeneous soil and sealed soil and specific values observed at different depths during (a) the infiltration stage, (b) the redistribution stage.Event 6 of Fig. 2 (see also Table2).Saturated hydraulic conductivity of the parent soil 2 mm h −1 and seal layer 0.12 mm h −1 .

Figure 9 .Fig. 9 .
Figure 9. Water content profiles simulated by the adopted model for homogeneous soil and specific values observed at different depths during (a) the infiltration stage, (b) the redistribution stage.Event 7 of Fig. 2 (see also Table2).Saturated hydraulic conductivity 2 mmh -1 .

Figure 10 .Fig. 10 .
Figure 10.Observed soil water content at 5 cm depth and air temperature at 2 m above the 3 ground as functions of time in a limited interval of the study period.4 5 6 7 8 9 10 11

Table 1 .
Main characteristics of the soil used for the experiments and model parameters.The quantities θ r and θ s denote the volumetric water content residual and at saturation, respectively.For other symbols see text.All the quantities are considered invariant from the sealing layer to the parent soil.
* Involved in the characteristic equations (see ).

Table 3 .
Variability of the calibrated saturated hydraulic conductivity, K 1s , and hydraulic resistance (Z c /K 1s ) for different hypothetical seal thicknesses, Z c .All the events of Table2influenced by a sealing layer are considered.