Isabelle Braud
email: isabelle.braud@cemagref.fr

Isotopic composition of bare soil evaporated water vapor. Part II: Modeling of RUBIC IV experimental results

Keywords: Soil evaporation, Soil water, Water stable isotope, Kinetic fractionation factor, RUBIC IV experiment, Modeling

Stable water isotopes such as oxygen 18, are tracers of water movement within the soilvegetation -atmosphere system. They provide useful information for a better understanding of evaporation and water vapor transport within soils. In part I of this paper, we presented a novel control experimental set up under non steady conditions, dedicated to the measurement of the evaporation flux and corresponding isotopic composition from six bare soil columns.

Data analysis raised several questions about the soil depth controlling the isotopic composition of the evaporated water vapor, suggesting different behavior before and after the appearance of back diffusion. Experimental data also suggested a time variable value of the kinetic fractionation factor. The present paper presents the modeling of the experimental results using the coupled heat, water and stable isotope transfer model SiSPAT_Isotope.

Model results were used for investigating the above questions more in details. For this purpose, model parameters were calibrated for each soil column in order to reproduce the data. Then model results were inverted to estimate the kinetic fractionation factor. The results show that the hypothesis that the kinetic fractionation factor varies in time is plausible but the uncertainty is too large to derive firm conclusions. The largest uncertainty is found when the soil relative humidity is lower than one but water vapor is still negligible. When back diffusion has occurred, model results are the most robust and confirm that the isotopic composition of the evaporated water vapor is controlled by the soil isotopic composition of the liquid water at the peak. In this case, the retrieved kinetic fractionation factor is close to 18.9°/ °°, corresponding to laminar flow.

Introduction

Measurements of stable isotopic composition of soil water can provide useful information for better quantifying water transfer within soils (mainly vapor diffusion coefficients and hydraulic properties). We developed a physically based model, called SiSPAT_Isotope (Braud et al., 2005a) for bare soil, representing the full interactions between the atmosphere, the soil and stable isotope species. The model was evaluated against two sets of laboratory data (Braud et al., 2005b). The data were composed of soil columns which were let evaporating freely in the atmosphere. This first study evidenced lots of uncertainty in the modeling and interpretation of stable isotopic composition in terms of evaporation. The uncertainty was related to a lack of control of the experimental conditions, especially of the atmospheric relative humidity and to a lack of knowledge of the kinetic fractionation factor, which controls the isotopic flux close to the surface, for unsaturated soils. In Part I of this paper, we presented a novel control experimental set up, dedicated to the measurement of the evaporation flux and corresponding isotopic composition of evaporated water vapor from bare soil columns. The experiment allowed a more accurate determination of atmospheric conditions and evaporation flux than in previous studies. For the first time, the isotopic composition of soil evaporated water under quasi-saturated and unsaturated conditions was measured. In the present paper, the SiSPAT_Isotope model was used for interpreting the experimental results, focusing on the determination of the kinetic fractionation factor, as deduced from model results. Section 2 presents the model, the methodology retained to run it and calibrate its parameters, as well as the calculation of the kinetic fractionation factor from model results. Section 3 discusses the results and conclusions and perspectives are drawn in section 4.

Material and methods

The experimental set up

The RUBIC IV reactor was developed in order to determine the isotopic composition of the water vapor released by an evaporating soil and to monitor its time evolution as long as the soil was drying. The leak tight experimental set up allowed us to send a dry air gas flow simultaneously over six evaporating soil columns and to trap continuously all the water vapor released by evaporation. The water samples were then analyzed using a mass spectrometer to determine their content in oxygen 18. Therefore, the only water vapor source was that of the evaporating soil columns, without any contamination from external sources. Furthermore, a cryoscopic trapping ensured that all the water vapor coming from the soil evaporation was completely trapped. In the design of the RUBIC IV experimental set up (see Part I, Braud et al., submitted), the following constraints were taken into account:

-A continuous non destructive sampling for the water vapor phase and a destructive sampling for the liquid phase within the soil columns. It meant to dismantle the soil columns at successive dates, cut them into slices and analyze the isotopic composition of all samples, as well as their gravimetric water content -A tight experimental set up with regards to external contaminations and loss of vapor or water from the system through a leak -The continuous measurement and the regulation of climatic variables -The continuous measurement of the most important variables describing the soil moisture status in the soil (water pressure and water content)

The experiment was launched on April 11 2005 (Day of the Year, DoY 101) and lasted 338 days. At successive dates during the drying process corresponding to respectively 31,53,91,163,235 and 338 days of drying (DoD), one of the column was dismantled for analysis (Table 1). All the details, including the description of the reactor and the measurements which were carried out, are given in Part I of the paper [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF]

The SiSPAT_Isotope model

The coupled heat, water and stable isotope transfer model (SISPAT_Isotope) was based on a Soil -Vegetation -Atmosphere Transfer (SVAT) model (SiSPAT, Braud et al., 1995;[START_REF] Braud | SiSPAT users's manual[END_REF][START_REF] Braud | SiSPAT User's Manual Update[END_REF]. The bare soil version of the model was used in this study and the isotope transport module was extensively described in Braud et al. (2005a).

The isotope transport module is based on the mass conservation equation for the isotope species. Both the liquid and vapor concentrations as well as the liquid and vapor isotope fluxes are considered. The number of unknowns is reduced to the liquid isotope concentration only, by using the hypothesis of instantaneous equilibrium between the liquid and vapor phases. Details on the final set of equations, which can be written as a convection-diffusion type transport equation, as well as its numerical discretization are given in Braud et al. (2005a).

To close the system, upper and lower boundary conditions must be specified. For the bottom boundary condition, the user can chose between a constant value of the liquid water isotopic concentration or a bottom flux calculated as the product of the bottom liquid flux times the bottom liquid water isotopic concentration. For the upper boundary condition, the soil-atmosphere isotope flux was written as:

( ) a v ia s l is e w i a K a i q R q R M M r E - = * α α ρ (1)
where v is R and v ia R are the liquid isotopic ratio at the soil surface and the water vapor isotopic ratio in the atmosphere, respectively; M i (kg) and M w (kg) are the heavy and ordinary water molecular mass, respectively; q s and q a (kg kg -1 ) are the specific humidity at the soil surface and in the atmosphere, respectively; ρ a (kg m -3 ) is the air volumetric mass; r a (s m -1 ) is the aerodynamic resistance to heat and water vapor transfer; * e α (-) is the liquid-vapor isotopic fractionation factor at equilibrium, and α K (-) is the so-called kinetic isotopic fractionation factor. Several formulae of this factor, taken from the literature were previously reviewed and tested by Braud et al. (2005a) for synthetic data sets, showing a high sensitivity, which will be investigated in this paper using experimental data and model results. The model is able to perform calculation for oxygen 18 and deuterium. Only oxygen 18 was considered in this paper.

Model set up

Each 1D soil column (12 cm in diameter, 35 cm of soil in height) was divided into 132 cells, with a logarithmic increase of the thickness from about 0.1 mm close to the surface and to the bottom and 5 mm in the center of the column. The soil was assumed homogeneous. The retention curve, relating the volumetric water content θ (m 3 m -3 ) to the soil water pressure h (m) was described by the Van Genuchten (1980) expression and the hydraulic conductivity K (m s -1 ) by the [START_REF] Brooks | Hydraulic properties of porous media[END_REF] equation, respectively:

n m with h h m n 2 1 1 VG sat - = ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + = - θ θ (2) η θ θ θ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ = sat sat ) ( K K (3)
where θ sat (m 3 m -3 ) is the saturated water content, h VG (m) is a scaling parameter for soil water pressure, K sat (m s -1 ) is the saturated hydraulic conductivity and m, n, and η are shape parameters related mainly to soil texture [START_REF] Haverkamp | Soil water retention: I. Introduction of a shape index[END_REF]. In addition, soil water vapor diffusivity depends on the soil tortuosity τ (-) which is another parameter of the model.

The residual water content was set at zero, which is consistent with the measurements (minimum measured value was 0.01 m 3 m -3 ). The saturated water content was also known from the dry bulk density and initial water volume used to saturate each column from their base (see Table 2 in Braud et al., 2009 for details). For the retention curve of Eq. ( 2), the simultaneous measures of soil water pressure h and volumetric water content θ at -2.5 cm depth were used to fit the value of the h VG and n parameters (see Fig. 5 in [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF].

The other parameters were considered for optimization (see below).

To run the model, initial and boundary conditions must also be specified. They are given in Table 2. Let's concentrate on the upper boundary condition. In the original model (Braud et al, 2005a), surface evaporation and isotope fluxes were calculated from the climatic conditions at a reference level and the surface conditions by using a resistance type approach.

In the context of the RUBIC IV experiment, these fluxes were measured. The first idea was to force the model using these measured values of the fluxes. However, we found that, if the soil hydraulic properties were not correctly prescribed, numerical instability could be generated in the model, especially when it was forced by evaporation fluxes the soil was not able to provide. This solution could therefore not be used in the optimization process, to avoid failure of the algorithms. For the consistency of the model in the optimization procedure, the model was let calculate the evaporation flux using the following formula:

av a s a r q q E - = ρ (4)
ρ a and q a being deduced from the measured air temperature and humidity at the outlet of the column, which can be assumed to be representative of the conditions prevailing above the soil surface [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF]; q s was calculated by the model. The unknown was therefore the aerodynamic resistance r av . Unfortunately, the wind speed above the column was not measured continuously. Assumptions were therefore required for the estimation of the resistance. They are detailed hereinafter in the description of the optimization process. Once the parameters were optimized, the model was also run by forcing it with the measured values of the evaporation flux, to verify the consistency of the results.

When running the isotope transport module, the measured isotopic composition v iout R of the outlet water vapor was used and the surface isotope flux was forced as:

E R M M E v iout w i i =
(5)

Another method could have been used to provide the isotope transport module surface boundary condition, at least before the appearance of back diffusion. Equation (15) of Part I (Braud et al., 2008) could have been used to estimate the soil surface liquid isotopic composition from the measured isotopic composition of the evaporated water vapor. However this solution was not retained as (i) it required the choice of a value for the kinetic fractionation factor, we wanted to estimate independently and (ii) it was not applicable for the whole experiment.

Optimization of parameters of the water and heat transfer modules

The first application of the SiSPAT_Isotope model presented in Braud et al. (2005b) showed that a proper simulation of water transport was a prerequisite to a correct simulation of stable isotope concentrations. As the major objective of the present study was to understand the processes involved in evaporation, we chose first to optimize the missing soil and aerodynamic parameters in order to reproduce the observations as well as possible. In a second step, this allowed us to use model results for testing various hypotheses about active processes and for interpreting and deriving the kinetic fractionation factor. For the optimization, the HAUS59 software, based on the iterative Gauss-Marquardt method and the maximum likelihood criterion [START_REF] Bard | Non-linear parameter estimation[END_REF] was used. The optimization function F(P aram ) was defined as follows: y are the corresponding observations. The k index refers to the various series which were used. For column 1, the series of volumetric water content θ at depths -2.5, -8, -16 and -24 cm and the evaporation flux E (kg m -2 s -1 ) were used. For the other five soil columns, only the evaporation flux was available. As the observed variables were of different nature, the differences between modeled and observed values were normalized using k obs y , the average of each observation series. Indices i=1, n refer to the various time steps where observations were available. For the efficiency of the optimization, we used daily average values in Eq. ( 6), although the original observations were hourly values.

2 5 1 mod 1 ) ( ) ( ∑ ∑ = = ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - = k k
The optimization process was conducted in several steps, using only the coupled heat and water transport modules (without isotopes). Several strategies were tested on column 1 (which was provided with the largest number of observations and the longest time series, see Table 1) before the final strategy, presented below was retained.

In a first step, we tried to optimize three soil parameters, namely K sat , η and τ as well as a constant value of the aerodynamic resistance r av as used by [START_REF] Schneider | Novel evaporation experiment to determine soil hydraulic properties[END_REF]. We also tested an optimization including h VG and n, but the algorithm results showed that these parameters were correlated with K sat and η, without significant gain in the optimization function, but with much less robustness in the optimized parameters. Therefore the values of h VG and n, fitted on the measurements, were used in the remainder of the study. Figure 1.

shows the results of this first optimization for column 1. When using a constant aerodynamic resistance, Fig. 1 (top) shows that the calculated time series of the evaporation flux was not consistent with the measured one. Furthermore, Figure 1 (bottom) shows that the aerodynamic resistance, calculated by inverted Eq. ( 4) using the measured E and q s calculated by the model was not constant, as assumed in the optimization. Contrarily to the experimental conditions of [START_REF] Schneider | Novel evaporation experiment to determine soil hydraulic properties[END_REF], it was thus not possible to assume that the aerodynamic resistance was constant during the whole experiment. Figure 1 (bottom) suggests that the aerodynamic resistance can be approximated by straight lines segments. Except for the first one, the dates at which the slope of the lines was changing were corresponding with the dates where the constraints imposed on the air relative humidity (h a ) above the soil columns were modified as well (Table 3). As explained in Part I of the paper [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF], it was necessary to modify the constraint on the relative air humidity above the soil column to fulfill experimental constraints on the maximum value of the dry air flux flowing above the columns. The change in relative air humidity had of course an impact on the evaporation flux and soil isotope composition profiles. However, all these changes were monitored and taken into account directly into the modeling though the use of the measured evaporation flux and isotopic composition of the evaporated water vapor as surface boundary conditions.

A second set of optimization was therefore realized using the function given by Eq. ( 7) for the aerodynamic resistances. To avoid a large number of optimized values, only three segments were considered, with t 1 =DoY 207 16h 45 (first change in air relative humidity), t 2 = DoY 306 12 h (last change in air relative humidity), t 3 =DoY 439 12 h (end of the simulation). The following three resistances r av0 , r av1 , r av2 were considered for optimization.

( ) ( ) Several trials and errors were conducted to determine the sensitive parameters and get a robust optimization. They showed that only the Brooks and Corey soil parameter, K sat , (Eq.

( ) ( ) ( ) ( ) ( ) ⎪ ⎪ ⎪ ⎩ ⎪ ⎪ ⎪ ⎨ ⎧ ≤ ≤ = ≤ ≤ - - - + = ≤ ≤ - + = 3 2 2 2 1 1 1 2 1 2 1 1 1 0 1 0 0 t t t
(3)) was independent and the other one, η, was fixed at η = 6.67. In the same way, the soil tortuosity was fixed atτ = 0.4. Finally, only the three aerodynamic resistances and the saturated hydraulic conductivity were optimized. The optimization algorithm was run several times with different initial values of the corresponding parameters in order to control that the model had reached the absolute optimum. For the columns with the shortest time series (columns 5 and 6), the evaporation flux was found to be quite insensitive to the saturated hydraulic conductivity (evaporation was more atmosphere controlled than soil controlled during these periods). Therefore, the K sat value of column 1 was used and only the first two aerodynamic resistances (r av0 and r av1 ) were optimized (see Table 4). After the optimization, the model was run using the measured evaporation flux as the upper boundary condition with the optimized soil parameters in order to verify the consistency of the results.

Calibration of the isotope transport module

Results of Braud et al. (2005a, b) have shown that the value of the soil tortuosity τ was tightly related to the depth of the peak within the isotope concentration profiles. The value τ=0.4 which was retained in section 2.4 was fulfilling this requirement.

A second parameter to be specified in the isotope transport module is the dispersivity length Λ (m) which appears in the total liquid diffusivity * l i D (m 2 s -1 ) for isotope i. It is given by Eq. ( 8)

where lo i D (m 2 s -1 ) is the molecular diffusivity of isotope i in liquid water, θ l (m 3 m -3 ) is the volumetric water content and q l (m s -1 ) is the liquid water flux.

l l l lo i l i q D D θ τθ / * Λ + = (8)
The sensitivity of model results to this parameter was first investigated. Then, the value was calibrated by hand for each column, by matching observed and simulated maximum or surface isotopic composition within the soil profiles (see section 3.2).

Calculation of the kinetic fractionation factor using model results

In Part I of the paper, [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF] showed that, given the experimental conditions and the expression of the surface isotopic flux (Eq. ( 1)), the kinetic fractionation factor can be retrieved as:

( ) ( ) v iout a s a v iout s l is e K R q q q R q R - - = * α α (9)
To apply Eq. ( 9), we used the measured data for the specific humidity of the air above the soil column q a , and its isotopic composition v iout R , as well as the values calculated by the model for the soil surface specific humidity q s and the isotopic composition of surface liquid water content l is R .

As the depth controlling the kinetic fractionation within the soil profile remains unknown for unsaturated soils, additional calculations were performed using the specific humidity q s_front and the isotopic composition of the soil liquid content l front is R _ at the evaporation front. The latter can easily be identified within the soil profile as it corresponds to the depth at which liquid and vapor fluxes are equal [START_REF] Barnes | The distribution of deuterium and 18O in dry soils: 1[END_REF]Braud et al., 2005a). This depth also corresponds to the depth of the peak within the soil isotopic composition profile. The appearance of this peak is referred to as the appearance of back diffusion in the remainder of the paper. The concept of evaporation front appears in modeling approaches of evaporation proposed by [START_REF] Bastiaanssen | A new simulation model of bare soil evaporation in arid regions (EVADES)[END_REF] and revisited by [START_REF] Gowing | Evaporative flux from a shallow water table: the influence of a vapour-liquid phase transition[END_REF]. The soil humidity controlling evaporation is taken at the evaporation front with a soil resistance r soil added to the aerodynamic one to take into account the resistance to water vapor transfer within the dry soil layer. Therefore, Eq. ( 4) can be rewritten as:

soil av a front s a r r q q E + - = _ ρ (10)
where the subscript _front refers to the value at the depth of the evaporation front. The inversion of this equation allows the determination of r av +r soil .

For the isotopic surface flux, the analogy is more complicated. Between the evaporation front and the surface, the isotope flux occurs mostly in the vapor phase. Braud et al. (2005a) defined the vapor isotopic flux within the soil, v i q (kg m -2 s -1 ), as:

z C D q C q v i * v i v v i v i ∂ ∂ - = (11)
where q v (m s -1 ) is the vapor water flux, * v i D is the total vapor diffusivity (m 2 s -1 ) for isotope i.

v i C (kg m -3 ) is the vapor isotopic concentration, which is related to the liquid concentration l i C by:

l i i l i w v e v l i e w i v v i w i v i C C R M M R M M C * * * β ρ ρ α ρ α ρ = = = = (12)
with

w v e i ρ ρ α β * * = (13)
and * e α (-) is the liquid-vapor isotopic fractionation factor at equilibrium given by [START_REF] Majoube | Fractionnement en oxygène-18 et en deutérium entre l'eau et sa vapeur[END_REF] as a function of temperature T (K), ρ (kg m -3 ) is the volumetric mass of water either in the liquid, ρ w , or in the water vapor phase, ρ v .

If the diffusive term of the isotopic vapor flux is negligible as compared to the convective term in Eq. ( 11), and if the vapor flux is uniform between the evaporation front and the surface, thus, the isotopic composition of the evaporated water vapor would not be modified very much between the evaporation front and the surface. Eq. ( 9) can thus be applied with the values at the evaporation front and the estimation of the kinetic fractionation factor α iK becomes:

( ) ( ) v iout a front s a v iout front s l front is e K R q q q R q R - - = _ _ _ * α α (14)
Eqs. ( 9) and ( 14) should lead to the same results before the appearance of back diffusion.

Hypotheses leading to Eq. ( 14) are consistent with data analysis performed by [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF] suggesting that the isotopic composition of the evaporated water could be controlled by the peak liquid concentration within the soil profile. The validity of this hypothesis will be evaluated in the next section.

Results

Optimization of parameters and simulation of water transfer

Values of the optimized parameters are given in Table 4, together with their estimated confidence interval, provided by the optimization algorithm. The figures show that most of the parameters are estimated with a narrow uncertainty, except the third resistance r av2 . The uncertainty on the third resistance r av2 is larger for columns 2 and 4 because only the evaporation flux time series was used in the optimization, whereas both evaporation flux and soil moisture time series were used for column 1. Table 5 shows the correlation between the estimated values (Eq. ( 4) with the optimized resistances) and observed hourly values of the evaporation flux for all the columns. Figure 2 provides an example of the time series for columns 1 and 5. Results are much more satisfactory for the longest drying periods. The measured evaporation flux exhibits high frequency fluctuations, certainly related to time variations in the air temperature, which contrarily to the air relative humidity was not regulated. These fluctuations made the modeling and data interpretation more difficult. Furthermore, simulations were performed assuming vertically homogeneous soil profiles, whereas the real ones certainly exhibit some vertical heterogeneity due to heterogeneities in the soil filling and saturation.

Subsequent results analysis shows that it was relevant to optimize a parameter set for each soil column in order to get the best estimates of the evaporation flux. Indeed, model results interpretation in terms of isotope transfer appear very sensitive to the values of this flux.

Figure 3 and Table 5 provide the comparison between measured and simulated soil water content profiles both using the evaporation flux calculated by Eq. ( 4) and the measurement flux as surface boundary condition. The agreement is satisfactory, especially for the last profiles. When using the measured evaporation flux, results were very similar to those obtained by Eq. ( 4), except for column 3 (91 DoD, DoY 192) where the drying is much less when using the measured evaporation flux (77.2 mm versus 85.4 mm). The comparison of cumulated evaporation derived from three methods exposed in [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF] showed cumulative values varying between 77.2 and 80.4 mm. However, the simulated profile is closer to the observed one using Eq. ( 4). For column 4 (235 DoD, DoY 336), the model diverged when using the measured flux, due to an evaporation flux larger than the value the soil was able to provide. In the following, results of column 4 were discarded as they were inconsistent, putting some doubts on the quality of the measurements. The comparison of the simulations using the evaporation flux calculated by Eq. ( 4) and the measurement flux as surface boundary condition, shows that the hypotheses on the aerodynamic resistances and the parameter optimization were reasonable. These columns are representative of the results obtained before (column 5) and after (columns 2 and 1) the establishment of a back diffusion profile. It appears that model results are very sensitive to the value of the dispersivity length. An increase in Λ leads to a decrease of the peak value and conversely. Results also show that a zero value of Λ leads to a large overestimation of the peak, except for column 4 (see Fig. 5) where it was still underestimated with Λ=0. Before back diffusion occurs, model results are very sensitive to the dispersivity length with a range of about 15 °/°° for the surface values (column 5). The sensitivity of the peak value is smaller after the establishment of back diffusion and decreases as long as the drying period length increases: the range is about 10 °/°° for column 2 and 5 °/°° for column 1. The range of the simulated surface values is larger than for the peak (about 10 °/°° ) for these two columns.

Simulation of isotope concentrations profiles

After the establishment of back diffusion, it was not clear if the optimization of Λ should favor a better representation of the surface or of the peak values (as it was not possible to get a value leading to a correct simulation of both, except for column 2). The retained optimal value for the dispersivity presented in Table 4 for column 1 was chosen so that the surface value was correctly simulated, although it was leading to an overestimation of the peak value (Fig. 4,right).

Table 4 provides the optimized values of Λ for the six columns. Because dispersivity is wellknown to be sensitive to heterogeneities at small scales, this may explain the more or less erratic values obtained for the various columns and the absence of a logical evolution amongst them. For column 4, even a zero value provides an underestimation of the peak and very low values close to the surface. As mentioned before, the results of this column are thought to be non reliable. For column 1, none of the values of the dispersity was leading to a correct simulation of both the surface and peak value, whereas it was possible for column 2, but with a very large value of the dispersity. It must be mentioned that the sampling of column 1 was very difficult because the soil was very dry and very friable. Furthermore, as the water content close to the surface was very low, the measure of the oxygen 18 isotopic composition is probably less reliable than for column 2, as can be seen from the shape of the profile which is much more irregular.

Figure 5 shows the simulated oxygen 18 composition soil profiles obtained using the "optimal" value of the dispersivity provided in Table 4 for each soil column and an evaporation flux either calculated using Eq. (4) or measured. They are compared with the observations. Table 5 provides the statistics between the simulated and observed profiles in both cases. The agreement is very good for columns 6, 5, 2 and 1 (R 2 larger than 0.95) and results are similar for the two choices of the boundary conditions, with a slight decrease in performance when the measured evaporation flux was used. For column 3, model performance is poorer (R 2 =0.90 with a 0.73 slope) with a large decreases in statistical criteria when the measured evaporation was used, corresponding to an underestimation of the surface values. It must be related to the slower drying mentioned in the previous section. As for water content, we noticed a divergence of the model for column 4 forced by the measured evaporation flux.

Simulation of water and isotope fluxes

In order to get insight into the model results, we investigated possible processes which were involved and verified if hypotheses stated in section 2.6 hold. Figure 6 shows the liquid, vapor and total fluxes of oxygen 18 for columns 5, 2 and 1. It can be seen that below the peak of oxygen 18 isotopic content (located at the depth where liquid and vapor fluxes are equal, [START_REF] Barnes | The distribution of deuterium and 18O in dry soils: 1[END_REF]Braud et al., 2005a,b), liquid isotope (and water, not shown) flux is dominant whereas above the peak, vapor flux is dominant. Figure 6 also shows that the columns are not in a stationary state, as the total isotope flux increases from zero at the bottom to a maximum value at the surface (column 5). When water vapor flux is dominant (dashed lines in Fig. 6), the isotope flux is almost constant between the soil surface and the depth of the peak of oxygen 18 isotopic composition (columns 2 and 1). Thus, one of the assumptions made in section 2.6, is shown to hold. This non stationary characteristic of the experimental set up is important to keep in mind, as the traditional results in isotopic geochemistry, based on isotopic stationary states hypotheses, do not hold. The results obtained after the establishment of back diffusion can corroborate the findings of [START_REF] Saravanapavan | Analysis of rate-limiting processes in soil evaporation with implications for soil resistance models[END_REF]. They showed that, when the evaporation flux was limited by the capacity of capillary rises to provide the water, the depth of the evaporation front was adapting to allow the soil to transfer the water flux from deeper layers towards the surface. Our experimental and modeling results confirm this finding as the vapor flux is almost constant between the evaporation front and the soil surface.

Figure 7 (left) provides the ratio between the diffusive and total liquid oxygen 18 flux. This figure shows that, when the liquid flux is dominant (i.e. below the peak), the diffusive component of the flux is only a few percents of the total flux, which is dominated in this case by convective transport. This fraction increases from the bottom towards the top of the profile, up to the depth of the maximum oxygen 18 content. The diffusive flux is directed downwards whereas convective flux is directed upwards (not shown). Due to this opposite sign, the enrichment of the isotopic profile close to the surface is limited. The relative contribution of the diffusive part of the flux increases as long as the soil dries (Fig. 7 shows an increasing contribution from column 5 to 2 and 1).

The analysis of the vapor isotope flux (Fig. 7, right) only makes sense after the establishment of back diffusion, i.e. for columns 2 and 1 (and above the peak). The contribution of the diffusive part to the total vapor isotope flux is still restricted to a few percents. This confirms the hypothesis made in section 2.6, that convective transport of water vapor should be dominant between the evaporation front and the surface. Thus the isotopic composition of the evaporated water should be related to the isotopic composition of the soil liquid water at the peak. From this analysis of the model results, it can be concluded that Eq. ( 14) should be suited to derive an estimation of the kinetic fractionation factor, as will be shown in the following section.

Kinetic fractionation factor

Figure 8 shows the calculated kinetic fractionation factor before the establishment of back diffusion (columns 6, 5, 3) (top) and the corresponding soil surface relative humidity (bottom). The retrieved values are very different from one column to the other. Column 5 was exhibiting the most consistent results when the simulation was performed either using Eq. ( 4)

or the measured values and can be considered as the most reliable column. For all the columns, the calculated range for the kinetic fractionation factor is larger than the range reported in the literature and the estimation can sometimes become negative (column 3).

Apart from the differences in absolute values, Figure 8 shows that the three columns exhibit similar behaviors: larger values are calculated at the very beginning of the experiment, but results are certainly less reliable due to time required by the system to reach its equilibrium.

Then, the kinetic fractionation factor seems to stabilize (up to about DoY 150) around a value which remains very uncertain, and then begins to decrease to reach smaller values than those given in the literature. [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF] show that an isotopic stationary state was noticeable before DoY 150 (Fig. 9 in Part I) and it corresponds to a relative humidity very close to 1. The calculated kinetic fractionation factor decreases when the soil surface relative humidity begins to decrease, although the latter only reaches a value of 0.9975. Figure 8 seem to confirm one of the hypotheses raised in [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF] from the data analysis that the kinetic fractionation factor could vary in time. Stable values are obtained when the relative humidity is close to one, and then the kinetic fractionation factor seems to decrease, even below values reported in the literature.

Figure 9 provides the calculated kinetic fractionation factor for columns 2 and 1 for the optimum dispersivity length (0.046m for column 2 and 0.055m for column 1 -see Table 4) and a value of 0.01 m. For both cases, Figure 9 gives estimates derived from Eq. ( 9), i.e. using the calculated soil surface value (full line) and Eq. ( 14), i.e. using the calculated values at the evaporation front (dashed line). For column 2, Figure 9 (top) shows that the calculated kinetic fractionation factors are very low and part of the time negative when using the optimum dispersivity of Table 4. More realistic values were obtained with a dispersivity of 0.01 m, which was leading to an overestimation of the isotopic peak value (see Fig. 4). For column 1 and the optimized dispersivity (0.055m), the calculated kinetic fractionation factor falls in the range of literature values. For both columns, the shape of the curve is similar to that of columns 6, 5, 3 (Fig. 8, top) at the beginning of the drying process. When Eqs. ( 9) and ( 14) provide different values, it means that back diffusion process has begun to occur. It is remarkable that, for column 1, the calculated kinetic fractionation factor is very stable after back diffusion, whatever the value of the dispersivity is, and converges towards a kinetic fractionation factor of 18.9 (per mil), when using Eq. ( 14) (Fig. 9 bottom, dashed lines). This result confirms the hypothesis of an isotopic composition of the water vapor controlled by the value of the isotopic composition of the soil liquid water at the evaporation front after back diffusion has taken place. On the other hand, the use of Eq. ( 9) and of the surface value to retrieve the kinetic fractionation factor leads to lower values than in the literature (full lines Fig. 9 bottom) and the result depends on the dispersivity length. This result shows that the use of Eq. ( 26) given in Braud et al. (2005a):

( ) a v ia s v is w i a K a i q R q R M M r E - = α ρ (15)
to derive the surface isotopic flux, might either not be correct when back diffusion has occurred, or the formulation of the isotope transport in the vapor phase might be deficient.

Discussion

The results presented in the paper show that experimental data analysis and modeling are complementary. Data analysis allow to raise hypothesis about processes. These hypotheses can be tested using the model which also provides a way to "interpolate" between measurements. Model results also give access to variables which cannot be measured such as water and isotope fluxes and which are important to understand the processes. Indeed, in [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF], data analysis allowed raising some hypotheses about the processes involved in soil evaporation and the isotopic composition of the evaporated water vapor.

These hypotheses were that (i) the isotopic signature of soil evaporation (under non stationary conditions) might be controlled by the value of the isotopic content at the depth of the evaporation front and (ii) the kinetic fractionation factor was not constant during the whole drying process. The modeling results presented in this paper allowed partly confirming these hypotheses. The results rely on a careful calibration of the water transport module, using measured evaporation flux and hypotheses about aerodynamic resistances. Then, sensitivity tests were performed to assess the impact of the dispersivity length on the simulated oxygen 18 profiles. Although this parameter was optimized, it remains largely uncertain due to uncertainties on the measured soil oxygen 18 profiles, evaporation fluxes and isotopic composition of the evaporated water vapor. Therefore, the retrieved kinetic fractionation factors are not reliable in absolute values, except when the back diffusion is well established and the evaporation flux is located several cm below the soil surface (column 1). It is remarkable that all the columns show the same characteristic time evolution of this factor.

Larger values are retrieved at the beginning of the drying process, when the soil surface relative humidity is close to one. When the soil surface dries and back diffusion is not established, the kinetic fractionation factor decreases, even towards values lower than those reported in the literature. It is certainly the period, where our conclusions are the less certain.

Then, when back diffusion has established, the retrieved kinetic fractionation factor is quite stable and the calculated value is close to 1.0189 for oxygen 18.

As far as modeling is concerned, the surface isotopic flux was calculated using the measured isotopic signature of the water vapor (see Eq. ( 5)). If the latter was not known, Eq. ( 15) should be used instead. But the difference between Eqs. ( 9) and ( 14) in terms of retrieved kinetic fractionation factor, shows that this formula should be used with lower values than those of the literature and our results do not tell us how it could vary and how it could be formulated. Eq. ( 15) with values taken at the evaporation front could be used, but it is not clear if it would be consistent with the model formulation of vapor isotopic transfer.

As far as field studies are concerned, non stationary isotopic situations are likely to be encountered. The experiment results reported in Part I [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF] and in this paper show that the traditional method which derives the isotopic composition of evaporation from the [START_REF] Craig | Deuterium and oxygen 18 variations in the ocean and the marine atmosphere[END_REF] formula and the sampling of an arbitrary thickness (i.e. [START_REF] Yakir | The use of stable isotopes to study ecosystem gas exchange[END_REF] for the soil surface isotopic composition is not correct when back diffusion occurs. In the absence of back diffusion, we have shown that the isotopic composition of the evaporated water vapor is controlled by the value within a very thin soil surface layer, which is almost impossible to sample. The traditional method is therefore likely to underestimate the isotopic composition of evaporated water vapor. Furthermore, our results

show that the picture is much more complex and further efforts in modeling and results

analysis are required to improve the modeling approach and/or better assess the dependency of the kinetic fractionation factor with the soil surface relative humidity, water content and isotope fluxes profiles.

Conclusions and perspectives

This paper and Part I (companion paper by [START_REF] Braud | Isotopic composition of bare soil evaporated water vapour. Part I: Experimental setting and results[END_REF] From the RUBIC IV experiment and the subsequent modeling exercise, it can be concluded that for evaporating unsaturated soils, the isotopic enrichment of the liquid water within the soil is dependent on the kinetic fractionation factor. The value of this factor is suspected not to be constant in time but the accuracy of the experiment was not sufficient to draw firm conclusions about this hypothesis. Smaller values below those reported in the literature seem to be possible. The experiment allowed drawing reliable conclusions for very dry soils with a well established back diffusion. In this case, the composition of evaporated water vapor is controlled by the isotopic concentration at the depth of the peak and the kinetic fractionation factor is close to 1.0189, i.e. the [START_REF] Merlivat | Molecular diffusivity of H 2 16 O, HD 16 O and H 2 18 O in gases[END_REF] value for laminar flux. When back diffusion has not occurred and the soil has dried enough to have a relative humidity significantly smaller than 1, the experiment set-up did not allow complete understanding of all the processes involved in the isotopic composition of evaporated water vapor. Nevertheless, the experimental design appears to be promising but should be improved as stated above to allow deriving less uncertain conclusions in this case.

The study we have performed was restricted to controlled conditions, because it was the only way to isolate the isotopic composition of soil water from external contaminations. Under field conditions, the Keeling plot [START_REF] Keeling | The concentration and isotopic abundances of carbon dioxide in rural and marine air[END_REF] showing the isotopic composition of the water vapor as a function of the water vapor isotopic composition (values are measured at various altitudes from the soil surface up to about 10-12m height) is able to give elements showing the possible impact of soil evaporation on the isotopic composition of the air water vapor. If the soil surface is bare, the ordinate at the origin provides the isotopic composition of the evaporated water. This value can be compared with the use of the Craig and Gordon model which must be applied with caution. Our results show that the application of this formula using the isotopic composition of the soil water averaged over several cm thickness is likely to be underestimated, leading to biased estimation of the partition between bare soil evaporation and regional contribution.

In the presence of vegetation, the picture is more complex. If there are only two sources of water vapor (regional and transpiration) the ordinate at the origin should be the isotopic composition of transpired water (itself being equal to that of the tiller of the plants under stationary isotopic conditions). It this hypothesis does not hold, thus a contribution from bare soil evaporation must be invoked. In general the isotopic composition of the evaporation water vapor is obtained from the Craig and Gordon model. As before, our results show that the use of this formula using the isotopic composition of the soil water averaged over several cm thickness is likely to be underestimated, leading to biased estimation of the partition between bare soil evaporation and transpiration.

List of figures

Figure 1 Top: Measured evaporation fluxes (10 -6 kg m -2 s -1 ) (light line) versus calculated ones with a constant value of the aerodynamic resistance (black line) for column 1. Bottom:

Optimized constant aerodynamic resistance (black horizontal line) versus aerodynamic resistance recalculated using the measured evaporation fluxes and the calculated soil surface specific humidity (light line) for column 1. The prescribed relative humidity h a is also provided (see Table 3). Figure 8 Top: Time evolution of the calculated kinetic fractionation factor (using Eqs. ( 9) and ( 14) which provide the same value) for the columns 6 (dotted line), 5 (dashed light line), 3

(full line). The horizontal lines are the literature values of [START_REF] Merlivat | Molecular diffusivity of H 2 16 O, HD 16 O and H 2 18 O in gases[END_REF] corresponding to α K =14.0, 18.9 and 28.5 (per mil). Bottom: Soil surface relative humidity for the same columns.

Figure 9 Time evolution of the calculated kinetic fractionation factor using Eqs. (9) (full line) and ( 14) (dotted line) for the columns 2 (top) and 1 (bottom). Black lines correspond to the optimum value of the dispersivity length (0.046m for column 2 and 0.0055m for column 1see Table 4) and light lines are for Λ=0.01 m. The horizontal lines are the literature values of [START_REF] Merlivat | Molecular diffusivity of H 2 16 O, HD 16 O and H 2 18 O in gases[END_REF] corresponding to α K =14.0, 18.9 and 28.5 (per mil).

Figure 4

 4 Figure4shows the comparison of simulated and observed oxygen 18 soil water isotopic

  investigated the isotopic composition of bare soil evaporated water vapor. A novel controlled experiment was set up, which allowed measuring both the evaporation flux and its isotopic composition for six columns. Model results analysis showed that, although the air relative humidity was regulated, the variations in the air temperature were leading to fluctuations in the measured fluxes, which were complicating data and model analysis. It would be advisable to perform again the experiments by regulating the air temperature, to avoid diurnal fluctuations in the evaporation flux. As much as possible, the set up should be improved to get a more constant value of the air relative humidity. The experimental design allowed measuring the isotopic composition of evaporated water vapor with a good accuracy. However, the measurements of the soil isotopic composition of the liquid water profiles remain quite uncertain, due to difficulties in sampling a very fine thickness of soil. A better accuracy would have possibly been obtained by freezing the soil columns.

Figure 2

 2 Figure 2 Comparison of the calculated evaporation fluxes after optimization of the parameters

Figure 3

 3 Figure 3 Comparison of the observed and calculated soil water content profiles when the six

Figure 4

 4 Figure 4 Sensitivity of the simulated soil oxygen 18 ratio to the dispersivity length for

Figure 5

 5 Figure 5 Comparison of the observed and calculated oxygen 18 concentration profiles when

Figure 6

 6 Figure 6 Profiles of simulated oxygen 18 fluxes at the end of the drying period of columns 5

Figure 7

 7 Figure 7 Left: Ratio of diffusive over total liquid oxygen 18 fluxes for columns 5 (53 DoD)

Acknowledgements

The "Programme National de Recherche en Hydrologie" of the French "ECosphère COntinentale" program is acknowledged for providing the financial support of the experiment. A. Chanzy is thanked for providing the HAUS59 optimization code.


Ibrahim Mubarak
email: ibrahim.mubarak@cemagref.fr

Jean Claude Mailhol

Rafael Angulo-Jaramillo

Pierre Ruelle

Pascal Boivin

Mohammadreza Khaledian

Temporal variability in soil hydraulic properties under drip irrigation

Keywords: Soil hydraulic properties, Beerkan infiltration method, Structural pores, Drip irrigation, Wetting bulb

Predicting soil hydraulic properties and understanding their temporal variability during the irrigated cropping season are required to mitigate agro-environmental risks. This paper reports field measurements of soil hydraulic properties under two drip irrigation treatments, full (FT) and limited (LT). The objective was to identify the temporal variability of the hydraulic properties of field soil under high-frequency water application during a maize cropping season. Soil hydraulics were characterized using the Beerkan infiltration method. Seven sets of infiltration measurements were taken for each irrigation treatment during the cropping season between June and September 2007. The first set was measured two weeks before the first irrigation event. The results demonstrated that both soil porosity and hydraulic properties changed over time. These temporal changes occurred in two distinct stages. The first stage lasted from the first irrigation event until the root system was well established. During this stage, soil porosity was significantly affected by the first irrigation event, resulting in a decrease in both the saturated hydraulic conductivity K s and the mean pore effective radius ξ m and in an increase in capillary length α h . These hydraulic parameters reached their extreme values at the end of this stage. This behavior was explained by the "hydraulic" compaction of the surface soil following irrigation. During the second stage, there was a gradual increase in both K s and ξ m and a gradual decrease in α h when the effect of irrigation was overtaken by other phenomena. The latter was put down to the effects of wetting and drying cycles, soil biological activity and the effects of the root system, which could be asymmetric as a result of irrigation with only one drip line installed for every two plant rows.

The processes that affected soil hydraulic properties in the two irrigation treatments were similar. No significant change in ξ m and α h was observed between FT and LT. However, as a result of daily wetting and drying cycles, which were strongest in LT, the soil in this treatment was found to be more conductive than that of FT. This showed that most of the changes in pore-size distribution occurred in the larger fraction of pores.

The impact of these temporal changes on the dimensions of the wetting bulb was studied using a simplified modeling approach. Our results showed that there were marked differences in the computed width and depth of wetting bulb when model input parameters measured before and after irrigation were used. A temporal increase in capillary length led to a more horizontally elongated wetting bulb. This could improve both watering and fertilization of the root zone and reduce losses due to deep percolation. As a practical result of this study, in order to mitigate agro-environmental risks we recommend applying fertilizers after the restructuration of tilled soil. Further studies using improved models accounting for temporal changes in soil hydraulic properties are needed.

INTRODUCTION

Drip irrigation has become quite common thanks to its great potential to use less water and to localize chemical applications, thereby enhancing the efficiency of irrigation and fertilization and reducing the risk of pollution. However, these objectives can only be achieved if the irrigation system is correctly designed (e.g. emitter discharge rate, emitter spacing, tape lateral spacing, diameter and length of the lateral system) and well managed (e.g. irrigation scheduling and fertilization strategy) for any given set of soil, crop and climatic conditions.

In contrast to surface or sprinkler systems, the frequency of the water application under drip irrigation is high. This means the infiltration period is a very important stage of the irrigation cycle [START_REF] Rawlins | Principals of managing high frequency irrigation[END_REF]. A good knowledge of the soil hydraulic properties involved in the multidirectional infiltration process during the course of this cycle is required to optimize water applications. The ability to estimate the dimensions of the wetting bulb i.e., water extending laterally and vertically away from an emitter is an important criterion for the design of drip systems to ensure efficient irrigation and to avoid the movement of water beyond the root zone [START_REF] Bresler | Analysis of trickle irrigation with application to design problems[END_REF][START_REF] Zur | Control of irrigation amounts using velocity and position of wetting front[END_REF][START_REF] Zur | Wetted soil volume as a design objective in trickle irrigation[END_REF][START_REF] Revol | Infiltration from a Surface Point Source and Drip Irrigation 2. An Approximate Time-Dependent Solution for Wet-Front Position[END_REF]. Because analytical models provide a rapid way of determining the position of the wetting front [START_REF] Revol | Infiltration from a Surface Point Source and Drip Irrigation 2. An Approximate Time-Dependent Solution for Wet-Front Position[END_REF][START_REF] Cook | Infiltration from surface and buried point sources: the average wetting water content[END_REF][START_REF] Thorburn | Soil-dependent wetting from trickle emitters: implications for trickle design and management[END_REF], researchers have tried to develop a simple model to describe the soil wetting pattern with micro-irrigation systems. [START_REF] Schwartzman | Emitter spacing and geometry of wetted soil volume[END_REF] developed a simplified semi-empirical model of wetted soil geometry with surface trickle irrigation, which depends on specific parameters i.e., soil type (saturated hydraulic conductivity), emitter discharge per unit length of laterals, and total amount of water in the soil. Al-Qinna and Abu-Awwad ( 2001) estimated an exponential function with a water application rate to describe the horizontal width and the vertical depth of the advance of the wetting front. In their field study, [START_REF] Revol | Infiltration from a Surface Point Source and Drip Irrigation 2. An Approximate Time-Dependent Solution for Wet-Front Position[END_REF] found that the infiltration solutions of [START_REF] Philip | Travel times from buried and surface infiltration point sources[END_REF] provided good estimations of the radial (r) and vertical (z) distance of the wetted zone from the water sources. [START_REF] Warrick | Soil water dynamics[END_REF] reviewed many analytical solutions describing water infiltration from point and line sources.

In a particular soil-water-plant system and climatic conditions, the transport properties of the soil surface layer can change during the growing season. This temporal variation is likely due to modifications in surface soil conditions resulting from tillage practices [START_REF] Mohanty | Spatial analysis of hydraulic conductivity measured using disk infiltrometers[END_REF][START_REF] Cameira | Soil macropore dynamics affected by tillage and irrigation for a silty loam alluvial soil in southern Portugal[END_REF], and to the effects of the rooting system [START_REF] Shirmohammadi | Effect of surface conditions on infiltration for shallow water table soils[END_REF][START_REF] Rasse | Alfalfa root and shoot mulching effects on soil hydraulic properties and aggregation[END_REF][START_REF] Iqbal | Spatial Variability Analysis of Soil Physical Properties of Alluvial Soils[END_REF]. Wetting and drying cycles and the irrigation system can also alter the soil structure. Hydrodynamic behavior is consequently affected primarily by the current state of the soil structure, as well as its texture [START_REF] Mapa | Temporal variability of soil hydraulicproperties with wetting and drying subsequent to tillage[END_REF][START_REF] Messing | Temporal variation in the hydraulic conductivity of a tilled clay soil as measured by tension infiltrometers[END_REF][START_REF] Angulo-Jaramillo | Seasonal variation of hydraulic properties of soils measured using a tension disk infiltrometer[END_REF][START_REF] Cameira | Soil macropore dynamics affected by tillage and irrigation for a silty loam alluvial soil in southern Portugal[END_REF][START_REF] Mailhol | Simulation of furrow irrigation practices (SOFIP): a field-scale modeling of water management and crop yield for furrow irrigation[END_REF]. However, the above-mentioned studies dealt with traditional irrigation systems that supply large amounts of water to the soil system at low frequencies. Only [START_REF] Mapa | Temporal variability of soil hydraulicproperties with wetting and drying subsequent to tillage[END_REF] addressed the effects on soil hydraulics of the wetting and drying cycles caused by drip irrigation following tillage. These authors found that soil hydraulic properties changed significantly after only one wetting/drying cycle in a silty clay loam and in a clay loam. However, in their study, each wetting/drying cycle included wetting (18 h of irrigation) followed by drying (7-10 days). In our opinion, this irrigation schedule thus more resembled that of traditional irrigation systems. To our knowledge, no study has tried to identify temporal variations in soil hydraulic properties during a cropping season under highfrequency water application.

The objectives of the present study were (i) to characterize temporal variability of soil hydraulic properties due to changes in soil structure under high-frequency drip irrigation, using the Beerkan infiltration method, and (ii) to illustrate the effects of temporal variability on the geometry of the wetting pattern generated by emitters, i.e., radius and depth, using the 1 infiltration solution of [START_REF] Philip | Travel times from buried and surface infiltration point sources[END_REF]. 2

MATERIALS AND METHODS

Experimental Site, Soil and Agricultural Practices

Field experiments were conducted on a loamy soil containing an average of 43% sand, 40% silt and 17% clay in the plowed layer with a relatively small coefficient of variation. The experimental field is located at the Cemagref Experimental Station in Montpellier, France (43°40' N, 3°50' E) where there is a fully equipped meteorological station. In the 2006-2007 cultivation years, the field was plowed to a depth of 35 cm on November 15 with a moldboard plow. The seed bed (top 8 cm) was prepared in April with a rotary harrow to fragment the soil.

On April 24, maize (Pionner PR35Y65) was sown with a row spacing of 0.75 m and at a plant density of 100000 ha -1 . Water was delivered through a line drip irrigation system with an emitter spacing of 30 cm and a flow rate of 3.67 L h -1 m -1 . One drip line was installed for every two plant rows leading to a lateral spacing of 1.5 m. Two different irrigation treatments, i.e. a full treatment (FT) and a limited treatment (LT), were investigated. On the basis of estimated crop water demand, irrigation was applied daily for three hours in FT. As it was not possible to irrigate during weekends, supplementary irrigation of three hours were applied on Fridays and Mondays. This irrigation schedule was maintained until the last week of July, after which irrigation was reduced to two hours per day with no supplementary irrigation at weekends. In LT, irrigation water was applied for three hours every other day (only on three days per week) until the end of the irrigation season. In both treatments, the crop was first irrigated for approximately 7.5 hr on June 19 to maintain the wetting bulb throughout the cropping season.

Potential crop evapotranspiration ETc was calculated from the daily reference evapotranspiration [START_REF] Allen | Crop evapotranspiration: Guidelines for computing crop requirements[END_REF] and the crop coefficient Kc (Doorenbos and Pruit, 1977). The maximum value of Kc was 1.2. Figure 1 shows crop water demand and water inputs from rainfall and irrigation over the cropping season. The cumulative amount of irrigation water applied in FT was approximately equal to estimated ETc except at the end of irrigation season due to a breakdown in the pumping station. The cumulative irrigation amount applied in FT was 430 mm. In LT, the total water amount was 300 mm, i.e. approximately 60% of estimated crop evapotranspiration.

Infiltration Measurements and Data Collection

The Beerkan infiltration method, (a 3D axi-symmetric infiltration method) was used with a simple annular ring [START_REF] Haverkamp | Hydrological and thermal behavior of the vadose zone in the area of Barrax and Tomelloso (Spain): experimental study, analysis and modeling[END_REF][START_REF] Lassabatère | Beerkan Estimation of Soil Transfer parameters through infiltration experiments-BEST[END_REF]. Seven sets of infiltration measurements were performed during the cropping season. Major variations in soil structure correspond to the reorganization of soil particles following the first irrigation event [START_REF] Mapa | Temporal variability of soil hydraulicproperties with wetting and drying subsequent to tillage[END_REF][START_REF] Mailhol | Simulation of furrow irrigation practices (SOFIP): a field-scale modeling of water management and crop yield for furrow irrigation[END_REF]. For this reason, in the present study we tried to characterize the soil following the first irrigation event in detail. Measurements were taken two weeks before irrigation started, five days after the 1 st irrigation event, one week later, one week later around the time the root system reached maximum, one month later, two weeks later, and just prior to harvest.

At each measurement date, six to eight infiltration tests were performed on a preselected plot in the center of each treatment area. Each plot was 4.5×4.5 m 2 square and included three adjacent drip lines. Infiltration tests were performed under the irrigation pipelines midway between irrigated plant rows, and spaced 1.5 m apart.

Infiltration tests were performed using a 65 mm-radius cylinder inserted to a depth of about 1cm to avoid lateral water losses. A fixed volume of water (100 ml) was poured into the cylinder at time zero, and the time needed for the infiltration of the known volume of water was measured. When the first volume had completely infiltrated, a second known volume of water was added to the cylinder, and the time needed for it to infiltrate was added to the previous time. This procedure was repeated until apparent steady state was reached, i.e. until two consecutive infiltration times were identical, and the cumulative infiltration was then recorded [START_REF] Lassabatère | Beerkan Estimation of Soil Transfer parameters through infiltration experiments-BEST[END_REF]. In addition, before the infiltration test, a soil sample was extracted in the vicinity of the infiltration ring to determine the initial soil gravimetric water content and particle-size distribution (< 2 mm), the latter being analyzed using the sedimentation method. Another 200 cm 3 sample was taken to measure dry bulk density (ρ d ). Saturated volumetric water content (θ s ) was calculated as total soil porosity considering the density of the solid particles to be 2.65 g cm -3 .

Soil Hydraulic Characterization and Data Analysis

The BEST algorithm (Beerkan Estimation of Soil Transfer parameters) [START_REF] Lassabatère | Beerkan Estimation of Soil Transfer parameters through infiltration experiments-BEST[END_REF] was used to determine soil hydraulic properties. This specific algorithm determines characteristic hydraulic curves that take into account the van Genuchten equation for the water retention curve, h(θ), (Eq. 1a) with the Burdine condition (Eq. 1b) and the Brooks and Corey relation (Eq. 2) for hydraulic conductivity curve, K(θ), [START_REF] Burdine | Relative permeability calculations from pore size distribution data[END_REF][START_REF] Brooks | Hydraulic properties of porous media[END_REF][START_REF] Van Genuchten | A closed form equation for predicting the hydraulic conductivity of unsaturated soils[END_REF] :

[1a] [1b] [2]
where θ r and θ s are the residual and saturated volumetric water content [L 3 L -3 ], respectively; n and m are shape parameters, and h g is the scale parameter [L] of the water retention curve h(θ), K s is saturated hydraulic conductivity [L T -1 ] and η is the shape parameter of the K(θ) relationship. θ r is usually very low and was thus considered to be zero. θ s was calculated as total soil porosity considering the density of the solid particles to be 2.65 g cm -3 . The hydraulic properties were thus represented using five parameters: θ s , n, h g , K s , η. Following [START_REF] Haverkamp | Hydrological and thermal behavior of the vadose zone in the area of Barrax and Tomelloso (Spain): experimental study, analysis and modeling[END_REF], the shape parameters n and η were assumed to mainly depend on soil texture, while θ s , h g and K s were assumed to mainly depend on soil structure.

The BEST algorithm estimates the shape parameters from particle-size distribution by classical pedotransfer functions. BEST derives saturated hydraulic conductivity and sorptivity S (L T -0.5 ) by modeling the 3D infiltrations performed at zero water pressure head (i.e. Beerkan infiltration method). The 3D cumulative infiltration I(t) and the infiltration rate q(t) can be approached by the explicit transient two-term and steady-state equations given by [START_REF] Haverkamp | Three-dimensional analysis of infiltration from the disc infiltrometer. 2. Physically based infiltration equation[END_REF]. The reader is referred to the study of [START_REF] Lassabatère | Beerkan Estimation of Soil Transfer parameters through infiltration experiments-BEST[END_REF] for more details on fitting experimental data on infiltration on analytical expressions that can provide estimations of scale parameters. [START_REF] Lassabatère | Beerkan Estimation of Soil Transfer parameters through infiltration experiments-BEST[END_REF] also described the main characteristics of the BEST algorithm which is coded with MathCAD 11 (Mathsoft Engineering and Education, 2002).

The capillary length (α h ) was then estimated from sorptivity (S) and the other hydraulic parameters through [START_REF] Haverkamp | Soil Properties and Moisture Movement in the Unsaturated Zone[END_REF]):

( ) s s s p h K c S ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - - = η θ θ θ θ α 0 0 2 1 [3] n m 2 1- = ( ) η θ θ θ θ θ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - - = r s r s K K m n g r s r h h - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + = - - 1 θ θ θ θ
where θ 0 is initial water content and c p is a function of the shape parameters for the van Genuchten (1980) water retention equation (see [START_REF] Haverkamp | Soil Properties and Moisture Movement in the Unsaturated Zone[END_REF][START_REF] Lassabatère | Beerkan Estimation of Soil Transfer parameters through infiltration experiments-BEST[END_REF]. The relationship between the capillary length and the characteristic microscopic pore "radius" i.e., the mean characteristic dimension of hydraulically functional pore, ξ m , is known to be: [4] where σ is the surface tension [MT -2 ], ρ w is the density of water [ML -3 ] and g is the gravitational acceleration [LT -2 ]. Taking the properties of pure water at 20 °C as appropriate, Equation (4) reduces to [5] where α h and ξ m are expressed in mm.

Standard statistical analysis was used in this study. In both treatments, the Kolgomorov-Smirnov test was used to check the assumption of normality of the data sets. The effects of time on each soil structure parameter, K s , α h , ξ m and the comparison between the two treatments, FT and LT, were evaluated by analysis of significance using the t-test at 0.05 probability level.

Determining the Horizontal and Vertical Components of the Wetting Front

The analytical simulation model developed by [START_REF] Philip | Travel times from buried and surface infiltration point sources[END_REF] was used to estimate the surface radius (r) and the depth (z) of the wetted soil volume. Because this model relies on assumptions such as the homogeneity of soil hydraulic properties, we used it as indicative rather than perspective. [START_REF] Philip | Travel times from buried and surface infiltration point sources[END_REF] found that the travel time of the wetting front away from a point source radially and vertically is given for dimensionless time (T), vertical distance (Z) and radial distance (R) as:

( )

Z Z Z T + + - = 1 ln 2 2 [6] and 2 2 exp 2 3 2 - ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + - = R R R T [7]
where

3 16 h qt T α θ π Δ = [8a] Z z h α 2 = [8b] R r h α 2 = [8c]
where q is the emitter flow rate, t is time, α h is the capillary length and Δθ is the difference between the average volumetric water content in the wetted soil and the initial volumetric water content.

RESULTS AND DISCUSSION

Soil Hydraulic Properties

Combining analysis of particle size distribution with modeling of the 3D infiltration experiments enabled us to fully determine the hydraulic parameters of the water retention and unsaturated hydraulic conductivity curves. Table 1 summarizes the statistical parameters of data sets of physical and hydraulic parameters.

The shape parameters of h(θ) and K(θ) varied little over time. This low variability is consistent with the assumption that the shape parameters mainly depend on soil texture [START_REF] Haverkamp | Hydrological and thermal behavior of the vadose zone in the area of Barrax and Tomelloso (Spain): experimental study, analysis and modeling[END_REF].

The K s , α h , and ξ m data sets were consistent with a log normal distribution (Kolgomorov-Smirnov test). Similar results were reported in studies of soil hydraulic properties (Mulla and McBratney 2002). Statistical analysis was thus performed on logtransformed values. At almost every measurement date, K s showed the highest variability with coefficients of variation (CV), ranging from 11% to 43%. The data sets of both ξ m and α h showed low CV values ranging from about 14% to 24% over the whole measurement period. This confirms the homogeneity of soil preparation, the precision of the Beerkan infiltration method, and low spatial variability at small scales.

Figures 2, 3 and4 show the values of K s , α h , and , ξ m obtained during the cropping season in the two treatments, FT and LT, as a function of time and as a function of cumulative water application depth. Each graph can be divided into two separate stages. The first stage lasted from the first irrigation event until July 10 around the time the root system was well established. The second stage lasted until the end of the cropping season. During the 1 st stage, hydraulic properties were significantly affected by the first irrigation event (according to the ttest) (Table 2). The hydraulic parameters reached their extreme values at the end of this stage. In FT, the mean value of K s decreased sharply over time and with cumulative water application depth. Its value dropped from 8.5×10 -3 mm s -1 to 1.6×10 -3 mm s -1 at the end of this stage. It subsequently increased gradually during the second stage to reach a mean value of 5.3×10 -3 mm s -1 just prior to harvest (Fig. 2 and Table 3). In LT, K s reached its minimal value (3.3×10 -3 mm s -1 ) only five days after the beginning of the irrigation season. Subsequently, the mean value of K s did not significantly vary, according to the t-test, for a period of some weeks before rising significantly again to reach 7.7×10 -3 mm s -1 at the end of the second stage (Fig. 2 and Table 3).

In both FT and LT, ξ m decreased from 0.14 mm just before the first irrigation event to about 0.06 mm at the end of the first stage. Subsequently, this parameter gradually increased to reach approximately 0.10 mm at the end of the second stage (Fig. 4 and Table 3). As the mean effective pore radius is inversely proportional to capillary length, the latter displayed inverse behavior with respect to ξ m (Fig. 3 andTable 3). It increased significantly from about 54 mm before the start of the irrigation season to a maximal value of 120 mm at the end of the first stage. It then decreased gradually to reach a value of about 72 mm just before harvest. The temporal variations in soil hydraulic properties were in agreement with the temporal change in dry bulk density during the cropping season. Dry bulk density showed a significant increase five days after irrigation started, and then stabilized somewhat before starting a downward trend at the end of the second stage (Fig. 5).

Temporal variations in soil hydraulic properties were, in general, similar in FT and LT but with different apparent response intensity (Figs 2, 3 and 4). In the second stage, K s values were significantly higher in LT than in FT, while in the first stage the values did not differ (according to the t-test) (Table 3). No significant differences were found in α h and ξ m between the two treatments at each measurement date (t-test) (Table 3). As a function of cumulative water application depth, FT showed that during the irrigation season, the three parameters reached their extreme values at the same cumulative water application depth of 144 mm. This value was much higher than that in LT (81 mm) (Figs 2, 3 and 4).

In the two treatments, the temporal changes in K s were generally associated with changes in ξ m at all measurement dates and inversely to changes in both ρ d and α h . The mean characteristic pore dimension and the capillary length are actual estimates of the capillary component of water transfer into the soil. Thus, the higher ξ m , the greater the effect of gravity compared to capillarity, as the infiltration driving force.

Temporal changes in the soil hydraulic properties can be commented as follows. At the beginning of the irrigation season, temporal variations in soil hydraulic properties could be due to the "hydraulic" compaction of soil following water application. This is in agreement with the concomitant increase in the dry bulk density (Fig. 5 and Table 1). The volume of water delivered during the first irrigation event and the following frequent irrigations caused restructuration of the fragile structural porosity created by soil preparation operations [START_REF] Mapa | Temporal variability of soil hydraulicproperties with wetting and drying subsequent to tillage[END_REF][START_REF] Angulo-Jaramillo | Seasonal variation of hydraulic properties of soils measured using a tension disk infiltrometer[END_REF][START_REF] Cameira | Soil macropore dynamics affected by tillage and irrigation for a silty loam alluvial soil in southern Portugal[END_REF]. During the second part of the irrigation season, when the irrigation rate decreased due to lower water requirements, the effect of irrigation was overtaken by the effects of wetting and drying cycles associated with the effects of the rooting system. The rooting system could be asymmetric as a result of irrigation with one drip line installed for every two plant rows. Some roots could grow and develop preferentially in the vicinity of emitters, creating new channels or continuity between existing pores. [START_REF] Shirmohammadi | Effect of surface conditions on infiltration for shallow water table soils[END_REF] observed that in cropped soils, hydraulic conductivity was much higher than in bare soils due to the effect of the rooting system. Similar results were also reported by [START_REF] Cameira | Soil macropore dynamics affected by tillage and irrigation for a silty loam alluvial soil in southern Portugal[END_REF], [START_REF] Iqbal | Spatial Variability Analysis of Soil Physical Properties of Alluvial Soils[END_REF] and [START_REF] Rasse | Alfalfa root and shoot mulching effects on soil hydraulic properties and aggregation[END_REF]. The latter authors found that the crop rooting system tended to increase water flow resulting in higher K s values and a larger macroporosity. In addition, as both α h and ξ m are characteristic microscopic lengths estimating the capillary component of water transfer, their temporal variations indicated a change in structure of the fine soil fraction. An increase in α h and a decrease in ξ m during the first stage indicates that the structural pattern of the fine fraction of soil changed from a poorly connected porous network before irrigation started to a more interconnected network in the middle of irrigation season. The presence of connected pores in that period may be the result of soil biological activity that is locally stimulated by the permanent high humidity of the wetted zone of the soil. [START_REF] White | Macroscopic and microscopic capillary length and time scales from field infiltration[END_REF] linked the higher value of the mean characteristic pore dimension to soil biological activity which is expected to be maximum in the surface soil layer. In addition, during the first stage, the change in these two parameters indicates a decrease in soil porosity leading to an increase in dry bulk density (Fig. 5). The reverse occurred during the second stage, confirmed by an increase in soil porosity at the end of irrigation season, and by a change in structural porosity from an interconnected porous network to a less connected one [START_REF] Angulo-Jaramillo | Seasonal variation of hydraulic properties of soils measured using a tension disk infiltrometer[END_REF].

The processes affecting soil hydraulic properties were similar in the two irrigation treatments. Variations in α h and ξ m over time were similar in FT and LT. Although we used two different irrigation schedules i.e., two different water amounts with two different application frequencies, the two microscopic lengths reached their extreme values at the same date i.e., at the end of the first stage, but not at the same cumulated water application depth. In our soil and climatic context, an added cumulated water application depth of from 81 mm to 144 mm did not lead to a significant change in the mean values of these two parameters. This indicates that both irrigation management systems play similar roles in the soil matrix. In other words, a change in irrigation management does not appear to adversely influence the fine fraction of the soil.

In contrast, K s values were significantly higher in LT than in FT during the second stage (Fig. 2 and Table 3). In LT, the soil was thus more conductive than in FT. This can be attributed to irrigation management. The frequency of water application used in LT (every other day) increased the importance of the alternating daily effects of wetting and drying cycles leading to the formation of the moderate cracks on the soil surface observed with this treatment during measurement campaigns. Thus, we can hypothesize that structural porosity is responsible for the largest water transport fraction [START_REF] Cameira | Soil macropore dynamics affected by tillage and irrigation for a silty loam alluvial soil in southern Portugal[END_REF]. This shows that different irrigation management systems mainly affect the larger fraction of pores. This confirms the effects of irrigation, the alternative cycles of wetting and drying, biological activity, root development on soil porosity.

The concomitance of extreme values of soil hydraulic properties in the two treatments may be due to the fact that these hydraulic properties are related to the same internal pore geometry. The difference between the two treatments in the cumulated water application depth at which our three hydraulic parameters reached their extreme values can be explained by the different volumes of water applied and by the irrigation schedule that stimulated or not the development of structural porosity in the surface soil layer.

Effect of Temporal Changes on Simulated Components of the Wetting Front

Based on the temporal changes in measured soil hydraulic parameters identified above, three sets of hydraulic properties were used to examine the effects of this phenomenon on the dimensions of the wetted soil volume. The first set represents the fragmented soil before irrigation started, so the initial values of soil hydraulic properties were used. The second set represents the restructuration of the soil some days after irrigation started and the average situation over the irrigation season. The third set represents the extreme situation of "hydraulic" compaction after three weeks of irrigation, that is to say at the 4 th measurement date (see Figs. 2, 3 and 4 from June 1 until July 8). The input parameters of the mathematical equations of Philip's model are (i) two given irrigation periods: t = 12 hr and 24 hr (ii) q = 2 l h -1 as emitter flow rate, (iii) the difference (Δθ = 0.2) between initial water content and mean value of water content of the wetted soil volume which are assumed to be uniform, and (iv) the mean value of capillary length for each hydraulic parameter set. Table 4 shows the width (r) and depth (z) of the wetting pattern computed from Eq. 6, 7 and 8 for the three sets and for the two irrigation periods. For the three parameter sets, it can be seen that the volume of wetted soil and of the water applied to the soil are similar. Results showed that there were differences in the dimensions of the wetting pattern. Initial soil hydraulic properties tended to form slightly vertically elongated bulb. The increase over time in the capillary length, i.e., characteristic of the capillary component of water transfer, led to an increase in the horizontal component of the wetting bulb (about 30%) and to a decrease in its vertical component (about 25%) for the two periods of irrigation concerned. In other words, it led to a horizontally elongated bulb. Differences in the dimensions were due only to temporal variability associated with water applications for the two given irrigation periods.

CONCLUSIONS

In this study, the behavior of a loamy soil under drip irrigation was analyzed using the Beerkan infiltration method to identify the temporal variability of its hydraulic properties caused by high-frequency irrigation during a maize cropping season. Two different irrigation treatments, a full (FT) and a limited (LT), were investigated. Our results demonstrated that both soil porosity and hydraulic properties varied over time. Soil behavior can be divided into two separate stages. The first stage lasted from the first irrigation event (56 days after sowing of the crop) until the root system was well established (75 days after sowing). During this stage, soil porosity was significantly affected by the first irrigation event, resulting in a decrease in both saturated hydraulic conductivity K s and in the mean pore effective radius ξ m , and in an increase in capillary length α h . These hydraulic parameters reached their extreme values at the end of this stage. This behavior was explained by the "hydraulic" compaction of the surface soil following irrigation.

Later in the season, during the second stage, a gradual increase in both K s and ξ m and a gradual decrease in α h were observed in both FT and LT, when the effect of irrigation was overtaken by other phenomena. The latter was put down to the effects of wetting and drying cycles, soil biological activity, and the effect of the root system, which could be asymmetric as a result of irrigation with only one drip line installed for every two plant rows.

The processes that affect soil hydraulic properties in the two irrigation treatments were found to be similar. No significant change in mean effective pore radius and in capillary length was observed between FT and LT. However, as a result of daily wetting and drying cycles, which were very marked in LT, the soil was found to be more conductive than in the FT. This shows that different irrigation management systems mainly affect the larger fraction of pores i.e., structural pores.

Our work raises questions about the usefulness of taking this phenomenon into account to improve the efficiency of both water and fertilizer under drip irrigation. To answer this question, we studied the impact of temporal changes in soil hydraulic properties on the dimension of the wetting bulb using a simplified modeling approach. Our results demonstrated that there were major differences in the computed width and depth of the wetted bulb when model input parameters measured before and after irrigation were used. An increase in capillary length over time led to the more horizontally elongated bulb. This change can improve the efficiency of both irrigation and fertilization of the root zone and reduce losses of both water and solute due to deep percolation. As a practical result of this study, we recommend fertigation after the restructuration of tilled soil in order to mitigate agroenvironmental risks and to improve the efficiency of the use of solutes. The temporal changes in soil hydraulic properties identified in our study should be taken into account in future studies when simulating soil water transfer under drip irrigation in order to improve irrigation scheduling practices. NS indicates no significant difference.

1 st June: before irrigation started.

24 th June: some days after the 1 st irrigation event.

8 th July: rooting system approximately reaches its maximal value 10 th September: just prior to harvest. 24 th June: some days after the 1 st irrigation event.

8 th July: rooting system approximately reaches its maximal value

Fig. 1 Fig. 4

 14 Fig. 1 Cumulative crop water use (ETc) and water inputs from rainfall and irrigation over crop season.

Evolution of dry bulk density, ρ d (arithmetic mean with standard deviation) for both treatments with time as well as with cumulative water application depth.

	Initial
	value
	1.6						1.6
	1.2 1.4 30-M ay value Fig. 5 Initial ρd (g/cm3)	19-Jun 1 st irrigation	9-Jul	29-Jul	18-A ug	7-Sep FT LT	1.2 1.4	0	50 Cum ulated w ater application depth (m m ) 100 150 200 250 300 350 400 F T 450 LT

Table 1 . Statistical parameters of dry bulk density and hydraulic properties for both treatments FT and LT

 1

	Date	Treatment FT	ρd gcm -3	θs m 3 m -3	n	m	η	Ks mm s -1	αh mm	ξm mm
	1 st June	Mean(8) +	1.383	0.478	2.201	0.0915	12.93	8.49E-03*	53.9*	0.276*
		CV %	6.7	7.3				20.6	18.0	16.3
	24 th June	Mean(8) +	1.550	0.415	2.210	0.0951	12.52	2.79E-03*	81.2*	0.183*
		CV %	1.9	2.4				33.3	15.0	16.4
	1 st July	Mean(6) +	1.550	0.415	2.210	0.0951	12.52	1.79E-03*	111.1*	0.134*
		CV %	1.3	1.9				31.4	21.2	24.2
	8 th July	Mean(8) +	1.540	0.419	2.210	0.0948	12.54	1.60E-03*	121.0*	0.123*
		CV %	1.4	1.9				31.7	18.1	15.6
	6 th August	Mean(6) +	1.532	0.422	2.210	0.0946	12.57	2.76E-03*	85.6*	0.174*
		CV %	1.6	2.0				31.8	14.1	13.8
	22 th August	Mean(6) +	1.500	0.434	2.210	0.0939	12.65	4.17E-03*	73.5*	0.202*
		CV %	2.9	4.0				9.1	15.6	14.8
	10 th September	Mean(8) +	1.445	0.455	2.204	0.0926	12.80	5.26E-03*	72.2*	0.206*
		CV %	2.7	3.3				26.1	19.8	19.5
	Date	Treatment LT	ρd gcm -3	θs m 3 m -3	n	m	η	Ks mm s -1	αh mm	ξm mm
	1 st June	Mean(8) +	1.383	0.478	2.201	0.0915	12.93	8.49E-03*	53.9*	0.276*
		CV %	6.7	7.3				20.6	18.0	16.3
	24 th June	Mean(8) +	1.446	0.454	2.206	0.0933	12.72	3.25E-03*	70.5*	0.211*
		CV %	2.3	2.7				21.1	8.7	8.4
	1 st July	Mean(6) +	1.517	0.428	2.210	0.0950	12.53	3.35E-03*	88.4*	0.168*
		CV %	2.2	3.0				11.6	7.7	7.6
	8 th July	Mean(8) +	1.525	0.425	2.210	0.0952	12.51	4.07E-03*	107.1*	0.139*
		CV %	2.6	3.3				23.5	6.0	5.8
	6 th August	Mean(6) +	1.490	0.438	2.210	0.0943	12.60	6.15E-03*	79.7*	0.187*
		CV %	0.9	1.0				16.4	11.1	11.0
	22 th August	Mean(6) +	1.452	0.452	2.206	0.0934	12.70	6.25E-03*	74.3*	0.200*
		CV %	3.2	3.7				18.6	22.3	22.2
	10 th September	Mean(6) +	1.451	0.452	2.206	0.0934	12.70	7.72E-03*	69.7*	0.214*
		CV %	3.4	4.2				19.1	11.4	11.6

+ (Test size)

*Geometric mean

Table 2 . Comparison of soil hydraulic properties in each treatment i.e., Full Treatment and Limited Treatment Full treatment Limited treatment Date Ks mm s -1 αh mm ξm mm Ks mm s -1 αh mm ξm mm 1 st June

 2

		8.49E-03* 53.9	0.276	8.49E-03	53.9	0.276
	24 th June 2.79E-03	81.2	0.183	3.25E-03	70.4	0.211
	Difference	5.70E-03 + -27.3 + 0.093 +	5.24E-03 +	-16.5 + 0.065 +
	24 th June 2.79E-03	81.2	0.183	3.25E-03	70.4	0.211
	8 th July 1.60E-03 121.0 0.123	4.07E-03	107.1 0.139
	Difference	1.19E-03 + -39.8 + 0.060 + -8.21E-04NS -36.6 + 0.072 +
	8 th July 1.60E-03 121.0 0.123	4.07E-03	107.1 0.139
	10 th September 5.26E-03	72.2	0.206	7.72E-03	69.7	0.214
	Difference	-3.66E-03 + 48.7 + -0.083 + -3.65E-03 +	37.4 + -0.075 +
	1 st June 8.49E-03	53.9	0.276	8.49E-03	53.9	0.276
	10 th September 5.26E-03	72.2	0.206	7.72E-03	69.7	0.214
	Difference	3.24E-03 + -18.3 + 0.070 +	7.73E-04 +	-15.7 + 0.062 +

*Geometric mean value. Difference = Date 1 -Date 2 + Significance at 0.05 probability level.

Table 3 . Comparison of soil hydraulic properties between Full Treatment and Limited Treatment.

 3

	Date	Treatment	Ks mm s -1	αh mm	ξm mm
	24 th June	FT	2.79E-03*	81.2	0.183
		LT	3.25E-03	70.4	0.211
		Difference -4.61E-04NS 10.7NS -0.028NS
	8 th July	FT	1.60E-03	121.0	0.123
		LT	4.07E-03	107.1	0.139
		Difference -2.47E-03 + 13.9NS -0.016NS
	10 th September	FT	5.03E-03	72.2	0.206
		LT	7.72E-03	69.7	0.214
		Difference -2.69E-03 +	2.6NS -0.008NS
	*Geometric mean value. Difference = Date 1 -Date 2

+ Significance at 0.

05 probability level. NS indicates no significant difference.

Table 4 . Comparison of wetting pattern dimensions between different dates for input parameters: q = 2 Lh -1 as emitter flow rate, and the difference (Δθ=0.2) between initial water content and mean value of water content of the wetted soil volume supposed uniform t=12hr t=24hr 1 st June 24 th June 8 th July 1 st June 24 th June 8 th July The radius of the wetted soil volume, r (cm

 4

	)=	22.4	25.6	28.5	25.85	30.02	33.85
	The depth of the wetted soil volume, z (cm) =	67.3	58.8	52.7	92.5	79.73	70.51
	Applied water volume (liter) = qt	24	24	24	48	48	48
	Volume of the wetted soil (litre )= qt / Δθ	120	120	120	240	240	240
	1

st June: before irrigation started.

st irrigation

ACKNOWLEDGEMENT

The AEC of Syria is greatly acknowledged for the Ph.D. scholarship granted to Ibrahim Mubarak.