HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS		QUICK SEARCH:   [advanced] Author: Keyword(s): Year:  Vol:  Page:

This Article Abstract Figures Only Full Text (PDF) Free Alert me when this article is cited Alert me if a correction is posted Services Similar articles in this journal Similar articles in ISI Web of Science Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via ISI Web of Science (2) Citing Articles via Google Scholar Google Scholar Articles by Ramos, T. B. Articles by Pires, F. P. Search for Related Content PubMed Articles by Ramos, T. B. Articles by Pires, F. P. GeoRef GeoRef Citation Agricola Articles by Ramos, T. B. Articles by Pires, F. P. Related Collections Soil Physics Inverse Procedures/Parameter Estimation

ORIGINAL RESEARCH

Estimation of Soil Hydraulic Properties from Numerical Inversion of Tension Disk Infiltrometer Data

^a Estação Agronómica Nacional, Quinta do Marquês, 2784-505 Oeiras, Portugal
^b USDA-ARS, George E. Brown, Jr., Salinity Lab., 450 W. Big Springs Road, Riverside, CA 92507-4617
^* Corresponding author (Tiago_Ramos{at}netcabo.pt)

Received 27 June 2005.

	ABSTRACT

TOP EXECUTIVE SUMMARY ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

 
Many applications involving variably saturated flow and transport require estimates of the unsaturated soil hydraulic properties. Numerical inversion of cumulative infiltration data during transient flow, complemented with initial or final soil water content data, is an increasingly popular approach for estimating the hydraulic curves. In this study, we compared Mualem–van Genuchten (MVG) soil hydraulic parameters obtained from direct laboratory and in situ unsaturated hydraulic conductivity measurements with estimates using numerical inversion of tension infiltration data of four coarse- to medium-textured soils in Alentejo (Portugal). The laboratory methods used were suction tables, pressure plates, and the evaporation method as applied to undisturbed soil samples collected from the surface horizons of four different soil profiles. Field measurements were taken with a tension disk infiltrometer using consecutive supply pressure heads of –15, –6, –3, and 0 cm. Six MVG parameters (residual soil water content [_r], saturated soil water content [_s], empirical shape factors, , and , and saturated hydraulic conductivity [K_s]) were estimated from the field data by numerical inversion using the HYDRUS-2D software package, and compared with values estimated from the laboratory data. Macroporosity was also determined. The laboratory- and field-measured water retention curves were found to agree closely for most experiments as reflected by relatively high values of the coefficient of determination, the modified coefficient of efficiency, and the modified index of agreement (always >0.9949, 0.8412, and 0.8931, respectively). The unsaturated hydraulic conductivity curves were predicted less accurately, although good estimates of K_s were obtained.

Abbreviations: MVG, Mualem–van Genuchten

	INTRODUCTION

TOP EXECUTIVE SUMMARY ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

 
MATHEMATICAL MODELS are increasingly used to address a broad range of variably saturated flow and contaminant transport problems. Such simulations are generally based on numerical solutions of the Richards equation, which in turn require knowledge of the unsaturated soil hydraulic properties. These properties consist of the water retention curve, which relates the volumetric water content () to the soil water pressure head (h), and the hydraulic conductivity curve, which relates the conductivity (K) to h or . While a large number of laboratory and field methods are available for direct measurement of the hydraulic properties (e.g., Dirksen, 1991; Dane and Topp, 2002), most techniques remain time consuming and costly, especially for hydraulic conductivity (van Genuchten and Nielsen, 1985) of fine-textured soils. Moreover, the hydraulic properties frequently show significant variations in space and time due to subsurface heterogeneity, agricultural activities, shrink–swell phenomena of fine-textured soils, the effects of particle dispersion and soil crusting, and changes in the concentration and ionic composition of the soil solution (van Genuchten and imnek, 1996). This implies that many samples are required to quantify those properties for most large-scale applications.

Among the laboratory methods available for measuring the soil water retention curve are those based on porous media principles such as the use of classical sand and sand plus kaolin boxes (Stakman, 1974; Romano et al., 2002) and the pressure plate extractor (Dane and Hopmans, 2002a), as well as the hanging water column method (Dane and Hopmans, 2002b). Laboratory methods for the unsaturated hydraulic conductivity include steady-state procedures based on direct inversion of Darcy's law such as the long-column method (Corey, 2002) and the crust method (Bouma et al., 1983), as well as transient procedures that involve some type of approximation or simplification of the Richards equation, such as the horizontal infiltration method (Bruce and Klute, 1956), the hot-air method (Arya et al., 1975), and the evaporation method (Wind, 1968). Evaporation methods also allow simultaneous measurement of both the water retention function and the hydraulic conductivity.

Field methods are usually considered more realistic than laboratory methods because of the larger volume of soil involved and because of continuity in the soil profile vs. depth. Popular field methods include the instantaneous profile (Vachaud and Dane, 2002), and the internal drainage and zero-flux plane methods (Vachaud et al., 1978; Arya, 2002). Methods employing tension disk infiltrometers have recently also become very popular, especially for in situ measurements of the near-saturated (h > –35 cm) soil hydraulic properties (Perroux and White, 1988; Ankeny et al., 1991; imnek et al., 1999a). Tension infiltrometers are especially useful for quantifying the effects of macropores and preferential flow paths on infiltration in the field. The method requires only minimal disturbance of the soil, is relatively rapid, and functions most effectively for pressure heads close to saturation where soil macropores are hydraulically the most active (Ankeny et al., 1991). Also, tension disk infiltration rates integrate various properties of the porous medium underneath the infiltrometer, such as local-scale heterogeneity, soil structure, textural irregularities and soil layering, preferential pathways, and possibly anisotropy (Mohanty et al., 1997; imnek et al., 1999a; National Research Council, 2001; Young et al., 2004). Tension infiltrometry additionally is useful for characterizing the water flux of macroporous soils in terms of two-domain (dual-porosity or dual-permeability) models in which one domain pertains to the soil matrix where Darcy's equation can be applied, and the second domain is dominated by preferential flow through the macropores (e.g., imnek et al., 2003).

Tension disk infiltration data are traditionally analyzed using Wooding's (1968) analytical solution for unconfined steady-state infiltration from a disk. As an alternative to using Wooding's analysis, imnek and van Genuchten (1996) proposed a methodology that allows indirect determination of the hydraulic parameters from transient tension infiltrometer data using inverse modeling. The unknown unsaturated soil hydraulic properties are then estimated from observed cumulative infiltration data by numerical inversion of the Richards equation. Inverse methods are based on the minimization of a suitable objective or likelihood function, which expresses the discrepancy between the observed values and the predicted system response. Initial estimates of the parameters are iteratively improved during the minimization process until a desired precision is obtained (imnek and van Genuchten, 1996; Hopmans et al., 2002).

Interest in the use of inverse methods has increased dramatically during the last few years, even though no standard procedures have yet been adopted. Still, the procedure has many advantages in that most or all available information from an experiment (such as soil water contents, pressure heads, cumulative infiltration data, and independent water retention and hydraulic conductivity data), and potentially also soft data, can be included in the objective function (Durner et al., 1999; Abbaspour et al., 1997; Hopmans et al., 2002; Yeh and imnek, 2002; Wang et al., 2003).

imnek et al. (1998a, 1998b) were first to apply the inverse methodology to field data, while imnek et al. (1999a) applied the approach to laboratory infiltration data in conjunction with simultaneously measured in situ tensiometer data and soil water contents measured with a time domain reflectometry probe. In the latter study, relatively close agreement was obtained between near-saturated hydraulic conductivities estimated using inverse modeling and Wooding's analytical method; however, the simultaneously estimated soil water retention curve using inverse modeling deviated from independent steady-state soil water retention data obtained with a pressure chamber. imnek et al. (1999a) noted that water retention data determined from a transient tension disk infiltrometer should be more useful for field conditions than those obtained from steady-state laboratory methods. Still, few studies exist were the various direct laboratory and field inverse modeling approaches have been compared.

The main objective of this study was to further test the inverse modeling approach of imnek and van Genuchten (1996) by using the method to characterize the hydraulic properties of four field sites in Portugal. We compared the resulting hydraulic properties with independent estimates using Wooding's analysis, as well as with direct laboratory measurements. Results were further interpreted in terms of separate soil macroporosity measurements.

	MATERIALS AND METHODS

TOP EXECUTIVE SUMMARY ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

 
Wooding's Analytical Approach
The traditional approach for analyzing tension disk infiltrometer data is based on Wooding's (1968) analytical solution for unconfined steady-state infiltration from a disk, given by

[1]

where Q is the steady-state infiltration rate (L³ T^–1), R is the radius of the disk (L), K is the hydraulic conductivity (L T^–1), h₀ is the supply wetting pressure head (L), and M is the matrix flux potential (L² T^–1). The first term on the right represents the effect of gravitational forces and the second term the effect of capillary forces. Several approaches are possible for analyzing steady-state tension infiltrometer data using Eq. [1]. In this study, we used the relatively simple approach of Ankeny et al. (1991), which requires knowledge of two steady-state fluxes, Q(h₁) and Q(h₂), at two tensions, h₁ and h₂, obtained with the same disk infiltrometer. Their method leads to two equations containing four unknowns:

[2]

[3]

A third equation can be obtained by assuming a constant K(h)/M(h) ratio throughout the pressure range between h₁ and h₂. Support for using a constant relationship between K and M can be found in Philips (1985). Alternatively, one could also assume an exponential relationship between K and h (Gardner, 1958):

[4]

in which (L^–1) is the sorptivity parameter.

Wooding's analysis requires steady-state infiltration rates at different supply pressure heads. Depending on soil texture, it can take hours or even days to reach steady state in a field experiment. Previous studies (e.g., Bagarello et al., 2000) have shown that Wooding's approach tends to overestimate the hydraulic conductivity if steady-state infiltration is not reached. Nevertheless, a majority of studies using Wooding's analysis assume that steady-state conditions are obtained within 1 h (e.g., imnek et al., 1999a). The possible error is usually dismissed as being negligible relative to errors related to soil heterogeneity or lack of reproducibility of the infiltration experiments.

Inverse Solution Approach
Inverse analyses of tension infiltrometer data require numerical solutions of the following modified Richards equation for radially symmetric Darcian flow:

[5]

subject to initial and boundary equations of the form (Warrick, 1992)

[6a]

[6b]

[6c]

[6d]

where t is time (T), r is the radial coordinate (L), z is the vertical coordinate (L), being positive upward with z = 0 corresponding to the soil surface, h_i is the initial pressure head (L), _i is the initial water content (L³ L^–3), and h₀(t) is the imposed supply pressure head (L). imnek et al. (1999b) developed a quasi-three-dimensional finite element code, HYDRUS-2D, to solve the above set of equations.

The inverse analysis further requires a parameterization for the unsaturated soil hydraulic properties. If the MVG equations (van Genuchten, 1980) are used, the soil water retention (h) and hydraulic conductivity K(h) are given by

[7]

[8]

where S_e is effective saturation, _r and _s denote the residual and saturated water contents, respectively (L³ L^–3), K_s is the saturated hydraulic conductivity (L T^–1), and (L^–1), , and are empirical shape factors (van Genuchten, 1980; van Genuchten and Nielsen, 1985). imnek and van Genuchten (1996), among others, pointed out that the selected soil hydraulic model, and the number of parameters being optimized, generally influences the identification, uniqueness, and stability of the inverse solution. Also, tension disk infiltration in general is a wetting process, which suggests that the hydraulic parameters in Eq. [7] and [8] should represent wetting branches of the unsaturated hydraulic properties.

The inverse modeling approach by imnek and van Genuchten (1996) is based on minimization of the following objective function, :

[9]

where m represents the number of different sets of measurements (e.g., cumulative infiltration data, pressure heads, or additional information) used in the analysis, n_j is the number of measurements in a particular set, q_j*(t_i) is the specific measurement at time t_i for the jth measurement set, ß is the vector of optimized parameters (i.e., _r, _s, , , , and K_s), q_j(t_i,ß) represents the corresponding model predictions for parameter vector ß, and v_j and w_ij are weights associated with a particular measurement set j or measurement i within set j, respectively. imnek and van Genuchten (1996) used values of 1 for the weighting coefficients w_i,j in Eq. [9], thus assuming that variances of the errors inside a particular measurement set are all the same. The weighting coefficients v_j are used to minimize differences in weighting between different data types because of different absolute values and numbers of data involved, and are given by

[10]

This approach represents the objective function as the average squared deviation normalized by measurement variances _j². The different measurement sets could consist of cumulative infiltration data, unsaturated hydraulic conductivities obtained by Wooding's analysis, in situ determined pressure heads, or the final water content. Minimization of the objective function in HYDRUS-2D was accomplished using the Levenberg–Marquardt nonlinear minimization method (Marquardt, 1963).

In a test of the above inverse methodology, imnek and van Genuchten (1996) found that cumulative infiltration rates measured with a tension disk infiltrometer at one particular tension did not provide enough information to estimate more than two MVG soil hydraulic parameters. To obtain at least three parameters (i.e., , , and K_s), additional information was needed, such as the soil water content or pressure head corresponding to the last measured tension under the disk infiltrometer. imnek and van Genuchten (1997) concluded that a combination of multiple-tension cumulative infiltration data with measured values of the initial and final water contents yielded unique solutions for the unknown parameters.

Following imnek et al. (1998a, 1998b), our objective functions (Q, _i, _f) were defined in terms of the measured cumulative infiltration data (Q) at multiple pressure heads, and the initial and final soil water contents (_i and _f, respectively). The weighting coefficients w_ij in Eq. [9] for the different infiltration data points, as well as for the initial water content, were all assumed to be 1 since the observation errors of the measurements were unknown; however, the final water content was given a weight of 10 to guarantee a reasonable effect on the final results relative to the cumulative infiltration data. The six Mualem–van Genuchten parameters (_r, _s, , , , and K_s) were estimated simultaneously by numerical inversion from the data using HYDRUS-2D, and compared with those derived from the laboratory data using the RETC fitting program of van Genuchten et al. (1991).

Field Tension Infiltration Experiments
The field tension infiltration experiments were performed in Aljustrel and Alvalade (Alentejo), Portugal, in two experimental areas cropped with maize (Zea mays L.) and irrigated with a center-pivot irrigation system.

In Aljustel, the infiltration experiments were performed on three different Gleyic Luvisols (LVgl), and in Alvalade on a Haplic Fluvisol (FLha) (soil classification according to ISSS-ISRIC-FAO, 1998). The field tension infiltrometer measurements were performed twice for each of the four soils, designated as Run A and Run B. The two runs in each case were performed at a distance of 1 m from each other using tension infiltrometers with the disks detached from the supply and tension control tubes. A nylon mesh was attached to the disks (all having a radius of 10 cm) to improve hydraulic contact with the soil surface. The infiltrometers used in this study were capable of completing the infiltration tests at multiple tensions without interruptions, with the supply tube hence containing enough water to complete each set of experiments. A level was used to ensure that the disk and the infiltrometer base were always at the same level, as was the case during calibration of the tensiometer in the laboratory. This was to make sure that the pressure heads in the bubbling outlet at the bottom of the water supply tube and in the disk membrane were always the same. A fine layer of silica sand, having a particle diameter between 0.2 and 0.3 mm and a much higher K_s than the soil, was used to obtain good contact between the disk membrane and the soil. A porous mesh was additionally used under the sand to avoid possible blockage of soil macropores. The sand was moistened immediately before placing the disk membrane on the soil to further improve contact between the disk membrane and the sand, and to prevent air entry into the disk (Cameira et al., 2002).

All four infiltration experiments were conducted with consecutive supply pressure heads of –15, –6, –3, and 0 cm. Readings of the water supply tube were done visually. Figure 1 shows measured cumulative infiltration rates vs. time. The readings were taken as soon as the tension infiltrometer was installed in the field and the closure clamps at the air bubble entry tube had been opened.

View larger version (21K): [in this window] [in a new window]   Fig. 1. Measured cumulative infiltration rates and corresponding fitted values from the final estimates of the inverse solution, at consecutive supply pressure heads of –15, –6, –3, and 0 cm for four surface horizons of soils in the Alentejo area of Portugal.

Undisturbed samples (100 and 630 cm³) were collected from the top 2 to 10 cm of the A horizon of each profile at a distance of at least 1 m from the infiltration disk to determine the soil water retention and hydraulic conductivity properties, as well as the dry bulk density. Water retention data were determined in the laboratory on three samples of 100 cm³ each per horizon using suction tables with sand for suctions of 2.5, 10.0, 31.6, 63.1, and 100.0 cm, and with sand and kaolin for suctions of 199.5, 316.2, and 501.2 cm (Stakman, 1974), while a pressure plate apparatus was used for suctions between 1000 and 15 850 cm. The evaporation method (Wind, 1968; Halbertsma and Veerman, 1994) was further used to simultaneously estimate water retention and hydraulic conductivity data between pressure heads of approximately –50 and –800 cm. Two samples of 630 cm³ (10 cm diameter by 8 cm high) each were used for the evaporation method, with tensiometers placed at depths of 1, 3, 5, and 7 cm. The evaporation data were analyzed using procedures documented by Halbertsma and Veerman (1994). The same samples had been previously used to determine the K_s using a constant-head method (Stolte, 1997).

The water retention and hydraulic conductivity data were next analyzed in terms of Eq. [7] and [8] using the RETC computer program (van Genuchten et al., 1991). The dry bulk density was obtained by drying volumetric soil samples (100 cm³) at 105°C. Total porosity was determined from the maximum gravimetric water contents of the soil samples and the bulk density. The particle size distribution was obtained using the pipette method for particles having diameters <20 µm (clay and silt fractions), and by sieving for particles between 200 and 2000 µm (coarse sand) and between 20 and 200 µm (fine sand). These textural classes follow the Portuguese classification system (Gomes and Silva, 1962) and are based on international soil particle limits (Atterberg scale). The OM (organic matter) content, which quantifies the organic fraction of the soil on a mass basis, was estimated from the OC (organic carbon) content determined by the Walkley–Black method, using the relation OM = 1.724 x OC (Nelson and Sommers, 1982). The main physical and hydraulic property results are given in Table 1.

Statistical Analysis
Results obtained with the different approaches were compared using several statistical parameters, including the mean absolute error, the root mean square error, and the coefficient of determination. The mean absolute error (MAE), given by

[11]

describes the difference between the observations (O_i) and the model predictions (P_i) in the units of the variable, with N being the number of observations. The root mean square error (RMSE), given by

[12]

is the square root of the mean square error in the units of the variable. In general, RMSE MAE. The degree to which a RMSE value exceeds the MAE is usually a good indicator of the presence and extent of outliers, or the variance of the differences between the modeled and observed values (Legates and McCabe, 1999).

The coefficient of determination (R²) is the square of Pearson's product moment correlation coefficient (r), describing the proportion of the total variance in the observed data that can be explained by the model, given by

[13]

in which SSR is the regression sum of squares and SSQ is the total sum of squares. Values of R² range from 0.0 to 1.0, with higher values indicating better agreement.

According to Legates and McCabe (1999), correlation and correlation-based measures are very sensitive to extreme values (outliers), but relatively insensitive to additive and proportional differences between model predictions and observations. They suggested the use of alternative goodness-of-fit tests that overcome many of the limitations of correlation-based measures. Alternative goodness-of-fit tests used in our study are the coefficient of efficiency, the modified coefficient of efficiency, the index of agreement, and the modified index of agreement. The coefficient of efficiency (E) is given by (Nash and Sutcliffe, 1970)

[14]

in which

[15]

Values of E ranges from – to 1.0, with higher values indicating better agreement. The modified coefficient of efficiency (E₁) is given by

[16]

which is a modification of E to reduce the effect of squared terms. We also used the index of agreement (d) given by (Wilmott, 1981)

[17]

which ranges from 0.0 to 1.0, with higher values indicating better agreement between the model and the data, similar to the interpretation of the coefficient of determination. Finally, the modified index of agreement (d₁) given by (Wilmott et al., 1985)

[18]

also reduces the effect of squared terms, similarly to E₁.

In this study, we considered as observed data the MVG soil hydraulic (water retention and conductivity) functions fitted with RETC to the laboratory data. These functions always showed excellent agreement with the observations. Soil hydraulic curves generated by numerical inversion of the field tension infiltrometer measurements were used as the predicted values.

Macroporosity Estimates
Tension infiltrometry permits rapid measurements of the hydraulic properties near saturation where water flow is determined primarily by the macroporosity of a soil (Mohanty et al., 1997; Haws and Rao, 2004). While several methods have been used to characterize macroporosity and soil structure, one approach suggested by Watson and Luxmoore (1986) is to use the maximum number of effective pores per unit area (N), which can be calculated from the minimum pore radius, R(L), in a particular class, and application of the capillary equation in conjunction with Poiseuille's law to give

[19]

where µ is the viscosity of water (M L^–1 T^–1), _w the density of water (M L^–3), and K_d the difference in K between two consecutive tensions (L T^–1). Consistent with Eq. [19], the effective porosity (L³ L^–3) is given by

[20]

Use of the capillary rise equation indicates that a pressure head > –3 cm corresponds to the macropore class (r > 0.5 mm) following the classification suggested by Wilson and Luxmoore (1988). Supply pressure heads of –6 and –300 cm similarly hold for pore radii of 0.25 and 0.005 mm, respectively, which define the boundaries of the mesopore class using the same classification. Unfortunately, the mesopore class cannot be studied across the complete range of pore sizes (0.005 < r < 0.5 mm) using tension disk infiltrometry, since tension infiltrometers can be used only for pressure heads close to saturation, even with special equipment (e.g., Castiglione et al., 2005). To still quantify the effects of the larger pores, we established within the mesopore class two pore subclasses: Mesoporosity 1, defined by 0.25 < r < 0.5 mm, and Mesoporosity 2, given by 0.1 < r < 0.25 mm.

We note there that Eq. [19] holds for laminar flow and assumes that macropores are completely water filled and not interconnected, and that the effects of tortuosity and the presence of pore necks on flow can be neglected. Because of these assumptions, N represents merely an equivalent number of macropores, not the true value. While not completely accurate, this equivalent value can still provide a relative estimate (e.g., Logsdon et al., 1993) of the number of hydraulically active macropores within relatively small depth intervals (Watson and Luxmoore, 1986; Wilson and Luxmoore, 1988; Cameira et al., 2002).

	RESULTS AND DISCUSSION

TOP EXECUTIVE SUMMARY ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

 
Figure 2 shows the water retention and hydraulic conductivity curves obtained by parameter estimation from the tension disk infiltration data at consecutive supply pressure heads of –15, –6, –3, and 0 cm. Each plot contains results for the replicated infiltration tests (Run A and Run B). Also included in Fig. 2 are the separately measured laboratory water retention and hydraulic conductivity data, the RETC curves fitted to the laboratory data, and independently analyzed field infiltrometer data using the approach of Ankeny et al. (1991).

View larger version (28K): [in this window] [in a new window]   Fig. 2. Water retention, (h), and hydraulic conductivity, K(h), curves obtained through numerical inversion of the field-measured tension disk infiltrometer data (Runs A and B), by means of separate laboratory measurements (suction tables and pressure plates) followed by analysis with RETC, using the evaporation method on laboratory samples, and using Wooding's analysis. Results are for the surface layers of four soil profiles (LVgl1, LVgl2, LVgl3, and FLha) used in the experiments.

The MVG soil hydraulic parameters (_r, _s, , , , and K_s) for the various curves in Fig. 2, as well as confidence intervals for the parameters estimated by numerical inversion, are listed in Table 2. The laboratory data were found to correspond closely with the fitted RETC water retention and hydraulic conductivity curves as reflected by the coefficients of determination aways being >0.9940. This excellent agreement between observed and fitted curves made it possible to generate water retention and hydraulic conductivity data points for direct comparison with the field data, thus permitting a better comparison of the different methods involved. Table 3 presents values of the various statistical parameters used to compare the different methods.

View this table: [in this window] [in a new window]   Table 3. Mean absolute error (MAE), RMSE, r, R2, coefficient of efficiency (E), modified coefficient of efficiency (E1), index of agreement (d), and the modified index of agreement (d1) comparing the observed laboratory data and the numerical inversion results shown in Fig. 2 for water retention, (h), and hydraulic conductivity, K(h) for four soils.

The residual water content, _r, showed far less consistency among the laboratory and field measurements. Nevertheless, this parameter was determined satisfactorily for half of the situations (i.e., all of Run A); however, some problems are apparent for the LVgl2 and FLha B runs, as reflected by negative values of the lower confidence limit (which has no physical meaning). The upper limits of the confidence intervals for LVgl2 and FLha (0.2488 and 0.2957, respectively) are also unrealistically high for medium-textured soils. The _r parameter for LVgl3 was impossible to estimate and had to be fixed to zero. This poor definition of _r for the LVgl2 and LVgl3 runs is not surprising because of the relatively small water content ranges involved, with the second LVgl3 showing a difference between the initial and final measured water content of only 0.02 cm³ cm^–3. The narrow range of pressure heads (–15 h 0) and associated water contents made it very difficult to accurately determine the slope of the retention curve in this case. While additional measurements at the dry end would have been helpful (e.g., the wilting point) to better define _r, the majority of retention curves in Fig. 2 all seem to converge to the same general curve, which suggests that the range of pressure heads used in our measurements will be sufficient for most applications.

The parameters and defining the shape of the retention curves showed good agreement between the numerical inversions and the laboratory data, reflected in part by their relatively small confidence intervals. One complication in comparing the laboratory and field data is the hysteretic nature of the retention curve. Since our laboratory methods represented drying processes and the infiltration experiments wetting processes, hysteresis should be present in the retention curves. When Kool and Parker (1987) coupled the MVG model with the simplified scaling approach of Scott et al. (1983), only was used to describe hysteresis. They used _w for the wetting branch and _d for the drying curve of the retention, while keeping all other parameters constant. Their approach could not be tested with our data, in part because the initial water content, _i, had been measured during desorption, after irrigation, following high evaporative demand conditions. The initial water content, as we will show below, was found to be essential for defining the shape of the curves.

The statistical analysis resulted in high values of R², E, and d, with the values always exceeding 0.9923, 0.9887, and 0.9984, respectively (Table 3); however, visual analysis still showed some differences (Fig. 2) between the RETC-derived curves and the curves obtained for Run B of the LVgl2 and FLha infiltration experiments. Legates and McCabe (1999) considered values of E₁ and d₁ of particular interest. The advantage in using of E₁ and d₁ is that errors and differences are given more appropriate weighting, not inflated by their squared values. While squaring in statistics is useful since squares are easier to manipulate mathematically than absolute values, the use of squares places relatively more weight on the larger values. An analysis of our retention curves showed also high values of E₁ and d₁, although lower than the squared coefficients, with E₁ between 0.8081 and 0.9263, and d₁ between 0.8918 and 0.9618, except for the LVgl2 and FLha B runs. For these two experiments, the agreement was not as good as for the other examples, with the LVgl2 curves starting to diverge from the laboratory results for log₁₀ h > 2.7, leading to E₁ and d₁ values of 0.6116 and 0.7724, respectively, while FLha started to differ from the laboratory curve for log₁₀ h > 4.2, giving E₁ = 0.6667 and d₁ = 0.7930. The values of RMSE were generally very low and close to MAE for the retention curves, thus indicating good agreement between the two field and laboratory data sets.

Hydraulic Conductivity Curves
Estimated values of K_s also corresponded satisfactorily among the various methods, i.e., the constant-head method, the method of Ankeny et al. (1991), and numerical inversion. The laboratory K_s closely reproduced the field-measured values, with the exception of LVgl1, for which we found a value of 53.4 cm d^–1 using the constant-head method, and values of 17.3 and 15.1 cm d^–1 using the infiltration tests. The close agreement was somewhat surprising since K_s is generally the most difficult parameter to quantify in view of field-scale soil spatial and temporal variability and the use of small (10-cm-diam.) undisturbed soil cores in our laboratory tests. While tension infiltration experiments never reached complete saturation (e.g., imnek and van Genuchten, 1996), at least for homogeneous soil profiles, because of the imposition of zero or negative boundary supply pressure heads, we considered the _s and K_s values obtained with the tension infiltrometer with a pressure head of 0 cm to be excellent estimates of the saturated water content and the saturated hydraulic conductivity, considering all other errors that typically occur during field experimentation and equipment calibration.

Estimates of the shape parameter using numerical inversion did not agreed well with those obtained from the laboratory evaporation data for most cases, with the confidence limits showing considerable uncertainty resulting in large differences between the upper and lower limits.

Similarly as for the retention curve, the coefficients obtained in the statistical analysis of the hydraulic conductivity curves also produced relatively high values. The E₁ and d₁ produced more realistic values. For LVgl1, the hydraulic conductivity curves showed relatively large differences between the numerical inversion results and the laboratory evaporation data, leading to E₁ values of 0.5129 and 0.5058 for Runs A and B, respectively, and d₁ values of 0.5649 and 0.5649 for these same replicates; LVgl3 actually produced a negative value (–0.0628) for E₁. According to Wilcox et al. (1990), negative values indicate that the observed mean is a better predictor than the model. Again, this poor result was probably due to the small difference between the initial and final water content data in the field. By comparison, the hydraulic conductivity for FLha resulted in good agreement between the numerical inversion estimates and the laboratory evaporation results, leading to relatively high values for both E₁ (0.8783 and 0.8814 for Runs A and B, respectively) and d₁ (0.7516 and 0.8004 for Runs A and B, respectively).

Figure 2 shows also results obtained with Wooding's analysis using the methodology of Ankeny et al. (1991). The K(h) values determined at pressure heads of –15, –6, –3, and 0 cm, by numerical inversion and using the traditional approach, are displayed in Table 4. Note the good agreement between the two methods of analysis.

View larger version (15K): [in this window] [in a new window]   Fig. 3. Effect of a lower initial water content on the inverse estimated soil hydraulic functions [Run (i)], results obtained when water contents at –100 and –15000 kPa are included in the objective function [Run (pF)], and analysis with RETC.

View this table: [in this window] [in a new window]   Table 5. Soil hydraulic parameters, with respective lower and upper confidence limits, estimated from numerical inversion of the tension disk infiltrometer data using a low value for the initial soil water content (I) in the objective function [Run (i)], and with independently measured water contents at –100 and –15 850 cm were included in the objective function [Run (pF)].

View this table: [in this window] [in a new window]   Table 6. Mean absolute error (MAE), RMSE, r, R2, coefficient of efficiency (E), modified coefficient of efficiency (E1), index of agreement (d), and the modified index of agreement (d1) comparing the observed laboratory data and the estimated numerical inversion results shown in Fig. 3 for water retention, (h), and hydraulic conductivity, K(h).

imnek et al. (1998b)

imnek et al. (1999a)

Macroporosity Estimation
Table 7 shows estimates of the macroporosity as calculated with the approach of Watson and Luxmoore (1986) from the measured K(h) curves at supply pressures of –15, –6, –3, and 0 cm. Very similar results were obtained using the inversely estimated K(h) curves and the curves based on Wooding's analysis. Results for the LVgl2 and LVgl 3 A runs, which had less water infiltrated during the experiments (Fig. 1), showed some differences in the value of Mesoporosity 2. By contrast, FLha's macroporosity was considerably larger than the macroporosities of the Luvisols, as was expected since the cumulative infiltration rates for the Luvisols were much higher.

View this table: [in this window] [in a new window]   Table 7. Number of pores (N) and effective porosity () values associated with macropore size class, where r is pore size. Results are for hydraulic conductivity, K(h), measurements at different soil water pressure heads, h, estimated by numerical inversion (NI) and using the Ankeny et al. (1991) method (AN).

	CONCLUSIONS

TOP EXECUTIVE SUMMARY ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

imnek and van Genuchten, 2000

The water retention curves obtained by numerical inversion closely matched the laboratory-measured curves for the four surface horizons where the infiltration experiments were performed, presenting high modified coefficients of efficiency and modified indices of agreement. The hydraulic conductivity curves were predicted less accurately, although good estimates of K_s were obtained. Hydraulic conductivities obtained with the inversely estimated MVG parameters also corresponded well with results obtained using Wooding's traditional approach of disk infiltrometer data following the methodology of Ankeny et al. (1991). This correspondence was further reflected by the very similar estimates of the macroporosity, Mesoporosity 1, and Mesoporosity 2 pore classes we calculated from the estimated K(h) curves using the approach of Watson and Luxmoore (1986).

One major limitation of the numerical inversion method is its extreme dependence on the field-measured water content values. Due to soil spatial variability, some problems may arise, especially with the initial soil water content, which (unlike the final water content) cannot be determined at exactly the same location where the tension infiltrometer measurements are performed. This was shown in this study for the LVgl2 and LVgl3 field measurements, which produced curves that deviated considerable from the laboratory-derived curves. By comparison, the measured final water contents corresponded very well with the saturated water contents measured in the laboratory. Contributing to this good match was the fact that the final water contents were determined within the top 2 cm of the soil immediately after wetting of the sand layer. The HYDRUS-2D simulations indicated that the wetting front had already reached depths of >10 cm at all times when samples were collected.

We also obtained more reliable results when independently measured water contents at –100 and –15 850 cm were added to the objective function. While this will require more time and effort (and additional equipment), and as such negate some of the advantages of numerical inversion methods (i.e., speed and ease of use), the results will be more reliable in cases where more complete laboratory data are not available for comparison purposes. An alternative approach for especially fine-textured soils would be to simply fix the residual water content using pedotransfer functions.

	ACKNOWLEDGMENTS

 
This work was financed by PEDIZA II 1462.1 and AGRO 14, and also by the SAHRA Science Technology Center as part of NSF grant EAR-9876800. We acknowledge Dr. Gonçalo Jacinto of the Mathematics Department of Évora's University for his help with the numerical inversion methods.

	REFERENCES

TOP EXECUTIVE SUMMARY ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

 

• Abbaspour, K.C., M.Th. van Genuchten, R. Schulin, and E. Schläppi. 1997. A sequential uncertainty domain inverse procedure for estimating subsurface flow and transport parameters. Water Resour. Res. 33:1879–1892.[CrossRef]
• Ankeny, M.D., A. Mushtaque, T. Kaspar, and R. Horton. 1991. Simple field method for determining unsaturated hydraulic conductivity. Soil Sci. Soc. Am. J. 55:467–470.[Abstract/Free Full Text]
• Arya, L.M. 2002. Wind and hot-air methods. p. 916–926. In J.H. Dane and G.C. Topp (ed.) Methods of soil analysis. Part 4. Physical methods. SSSA Book Ser. 5. SSSA, Madison, WI.
• Arya, L.M., D.A. Farrel, and G.R. Blake. 1975. A field study of soil water depletion patterns in presence of growing soybean roots. I. Determination of hydraulic properties of the soil. Soil Sci. Soc. Am. J. 45:1023–1030.
• Bagarello, V., M. Iovino, and G. Tusa. 2000. Factors affecting measurement of the near-saturated soil hydraulic conductivity. Soil Sci. Soc. Am. J. 64:1203–1210.[Abstract/Free Full Text]
• Bouma, J., C. Belmans, L.W. Dekker, and W.J.M. Jeurissen. 1983. Assessing the suitability of soils with macropores for subsurface liquid waste disposal. J. Environ. Qual. 12:305–311.
• Bruce, R.R., and A. Klute. 1956. The measurement of soil moisture diffusivity. Soil Sci. Soc. Am. Proc. 20:458–462.
• Cameira, M.R., R.M. Fernando, and L.S. Pereira. 2002. Soil macropore dynamics affected by tillage and irrigation for a silty loam alluvial soil in southern Portugal. Soil Tillage Res. 70:131–140.[CrossRef]
• Castiglione, P., P.J. Shouse, B.P. Mohanty, D. Hudson, and M.Th. van Genuchten. 2005. Improved tension infiltrometer for measuring low fluid flow rates in unsaturated fractured rock. Vadose Zone J. 4:885–890.[Abstract/Free Full Text]
• Corey, A.T. 2002. Long column. p. 899–903. In J.H. Dane and G.C. Topp (ed.) Methods of soil analysis. Part 4. Physical methods. SSSA Book Ser. 5. SSSA, Madison, WI.
• Dane, J.H., and J.W. Hopmans. 2002a. Pressure plate extractor. p. 688–690. In J.H. Dane and G.C. Topp (ed.) Methods of soil analysis. Part 4. Physical methods. SSSA Book Ser. 5. SSSA, Madison, WI.
• Dane, J.H., and J.W. Hopmans. 2002b. Hanging water column. p. 680–683. In J.H. Dane and G.C. Topp (ed.) Methods of soil analysis. Part 4. Physical methods. SSSA Book Ser. 5. SSSA, Madison, WI.
• Dane, J.H., and G.C. Topp. 2002. Methods of soil analysis. Part 4. Physical methods. SSSA Book Ser. 5. SSSA, Madison, WI.
• de Vos, J.A., J. imnek, P.A.C. Raats, and R.A. Feddes. 1999. Identification of the hydraulic characteristics of a layered silty loam. p. 783–798. In M.Th. van Genuchten et al. (ed.) Characterization and measurement of the hydraulic properties of unsaturated porous media. Univ. of California, Riverside.
• Dirksen, C. 1991. Unsaturated hydraulic conductivity. p. 209–269. In K.A. Smith and C.E. Mullins (ed.) Soil analysis; physical methods. Marcel Dekker, New York.
• Durner, W., B. Schultze, and T. Zurmühl. 1999. State-of-the-art in inverse modelling of inflow/outflow experiments. p. 661–681. In M.Th. van Genuchten et al. (ed.) Characterization and measurement of the hydraulic properties of unsaturated porous media. Univ. of California, Riverside.
• Gardner, W.R. 1958. Some steady state solutions of the unsaturated moisture flow equation with the application to evaporation from water table. Soil Sci. 85:228–232.
• Gomes, M.P., and A.A. Silva. 1962. Um novo diagrama triangular para a classificação básica da textura do solo. (In Portuguese). Garcia da Orta 10:171–179.
• Halbertsma, J.M., and G.J. Veerman. 1994. A new calculation procedure and simple set-up for the evaporation method to determine soil hydraulic functions. Rep. 88. DLO Winand Staring Centre, Wageningen, the Netherlands.
• Haws, N.W., and P.S.C. Rao. 2004. The effect of vertically decreasing macropore fractions on simulations of non-equilibrium solute transport. Vadose Zone J. 3:1300–1308.[Abstract/Free Full Text]
• Hopmans, J.W., J. imnek, N. Romano, and W. Durner. 2002. Inverse methods. p. 963–1008. In J.H. Dane and G.C. Topp (ed.) Methods of soil analysis. Part 4. Physical methods. SSSA Book Ser. 5. SSSA, Madison, WI.
• ISSS-ISRIC-FAO. 1998. World reference base for soil resources. World Soil Resour. Rep. 84. FAO, Rome.
• Kool, J.B., and J.C. Parker. 1987. Development and evaluation of closed-form expressions for hysteretic soil hydraulic properties. Water Resour. Res. 23:105–114.
• Legates, D.R., and G.J. McCabe, Jr. 1999. Evaluating the use of "goodness-of-fit" measures in hydrologic and hydroclimatic model validation. Water Resour. Res. 35:233–241.
• Logsdon, S.D., E.L. McCoy, R.R. Allmaras, and D.R. Linden. 1993. Macropore characterization by indirect methods. Soil Sci. 155:316–324.
• Marquardt, D.W. 1963. An algorithm for least-squares estimation of non-linear parameters. SIAM J. Appl. Math. 11:431–441.[CrossRef]
• Martins, J.C., T.B. Ramos, F.P. Pires, A.V. Oliveira, M.A. Bica, N.L. Castanheira, J.L. Reis, and F.L. Santos. 2005. Runoff and soil erosion in a soil with maize irrigated by center-pivot under different tillage technique and polyacrylamide application. (In Portuguese, with English abstract). Rev. Ciênc. Agrárias 28:49–62.
• Mohanty, B.P., R.S. Bowman, J.M.H. Hendrickx, and M.Th. van Genuchten. 1997. New piecewise-continuous hydraulic functions for modeling preferential flow in an intermittent flood-irrigated field. Water Resour. Res. 33:2049–2063.[CrossRef]
• Nash, J.E., and J.V. Sutcliffe. 1970. River flow forecasting through conceptual models. I. A discussion of principles. J. Hydrol. 10:282–290.
• National Research Council. 2001. Conceptual models of flow and transport in the fractured vadose zone. National Acad. Press, Washington, DC.
• Nelson, D.W., and L.E. Sommers. 1982. Total carbon, organic carbon, and organic matter. p. 539–579. In A.L. Page et al (ed.) Methods of soil analysis. Part 2. Chemical and microbiological properties. Agron. Monogr. 9. ASA and SSSA, Madison, WI.
• Perroux, K.M., and I. White. 1988. Design for disc permeameters. Soil Sci. Soc. Am. J. 52:1205–1215.[Abstract/Free Full Text]
• Philips, J.R. 1985. Reply to "Comments on the steady infiltration from spherical cavities." Soil Sci. Soc. Am. J. 49:788–789.[Free Full Text]
• Romano, N., J.W. Hopmans, and J.H. Dane. 2002. p. 692–698. In J.H. Dane and G.C. Topp (ed.) Methods of soil analysis. Part 4. Physical methods. SSSA Book Ser. 5. SSSA, Madison, WI.
• Scott, P.S., G.J. Farquhar, and N. Kouwen. 1983. Hysteresis effects on net infiltration. Adv. Infiltration 11:163–179.
• Schwartz, R.C., and S.R. Evett. 2003. Conjunctive use of tension infiltrometry and time domain reflectometry for inverse estimation of soil hydraulic properties. Vadose Zone J. 2:530–538.[Abstract/Free Full Text]
• imnek, J., and M.Th. van Genuchten. 1996. Estimating unsaturated soil hydraulic properties from tension disc infiltrometer data by numerical inversion. Water Resour. Res. 32:2683–2696.[CrossRef]
• imnek, J., and M.Th. van Genuchten. 1997. Estimating unsaturated soil hydraulic properties from multiple tension disk infiltrometer data. Soil Sci. 162:383–398.