WATER RESOURCES RESEARCH, VOL. 35, NO. 10, PAGES 2965–2979, OCTOBER 1999

Estimating unsaturated soil hydraulic properties fromlaboratory tension disc infiltrometer experiments
Jiřı́ Šimůnek
U. S. Salinity Laboratory, USDA-ARS, Riverside, California

Ole WendrothInstitute for Soil Landscape Research, Zentrum für Agrarlandschafts und LandnutzungsforschungMüncheberg, Germany

Martinus T. van Genuchten
U. S. Salinity Laboratory, USDA-ARS, Riverside, California

Abstract. Four tension disc infiltration experiments were carried out on a loamy soil inthe laboratory for the purpose of estimating the unsaturated soil hydraulic properties.
Sixteen tensiometers were installed in pairs at the following coordinate (r, z) positions:(10, 2.5), (10, 5), (10, 10), (15, 5), (15, 10), (15, 15), (15, 20), and (15, 30), where r
represents the distance from the axis of symmetry and z is the location below the soil
surface. A time domain reflectometry (TDR) probe was used to measure water contents at
a depth of 2 cm directly below the tension disc. The first three experiments involved
supply pressure heads at the disc of ⫺20, ⫺10, ⫺5, and ⫺1 cm, with the experimentlasting for ⬃5 hours. The same supply pressure heads were also used for the fourthexperiment, which lasted 6.25 days so as to reach steady state at each applied tension. The
measured data were analyzed using Wooding’s [1968] analytical solution and by numericalinversion. The parameter estimation method combined a quasi three-dimensional
numerical solution of the Richards equation with the Marquardt-Levenberg optimization
scheme. The objective function for the parameter estimation analysis was defined using
different combinations of the cumulative infiltrated volume, TDR readings, andtensiometer measurements. The estimated hydraulic properties were compared against
results obtained with an evaporation experiment as analyzed with Wind’s [1968] method.
Water contents in the retention curves were underestimated when both transient andquasi steady state experiments were analyzed by parameter estimation. Unsaturated
hydraulic conductivities obtained by parameter estimation and using Wooding’s [1968]
analysis corresponded well. Drying branches of the hydraulic conductivity function
determined by parameter estimation also corresponded well with those obtained with theevaporation method.

1.

IntroductionMore attention is increasingly being directed to accurate
measurement of the unsaturated soil hydraulic properties close
to saturation [van Genuchten et al., 1999], i.e., to moisture
conditions that are strongly affected by soil structure andmacropores. Traditional transient laboratory methods, such as
outflow or evaporation experiments, show relatively little sensitivity to the hydraulic conductivity at near-saturated conditions and hence are more suitable for estimating the hydraulicconductivity at medium saturation levels. The evaporation
method usually fails in the near-saturation range where thehydraulic conductivity is highest, leading to very small hydraulic gradients that cannot be determined with sufficient accuracy[Wendroth and Šimůnek, 1999]. Hence there is a trend toward
determining the hydraulic conductivity in the wet range withsteady state experiments, such as the tension disc infiltrometer
method [Perroux and White, 1988] or the crust method [Bouma
This paper is not subject to U.S. copyright. Published in 1999 by theAmerican Geophysical Union.
Paper number 1999WR900179.

et al., 1971], whereas transient conditions are used for the drier
range [Wendroth and Šimůnek, 1999].
Tension disc infiltrometers have recently become very popular devices for in situ measurement of the near-saturated soil
hydraulic properties [Perroux and White, 1988; Ankeny et al.,
1991; Reynolds and Elrick, 1991; Logsdon et al., 1993]. Thus far,tension infiltration data have been used primarily for evaluating saturated and unsaturated hydraulic conductivities and forquantifying the effects of macropores and preferential flow
paths on infiltration. Tension infiltration data are generally
used to evaluate the saturated hydraulic conductivity K s andthe sorptive number ␣* given by Gardner’s [1958] exponential
model (equation (10) below) of the unsaturated hydraulic conductivity using Wooding’s [1968] analytical solution. Analysesof this type require either two infiltration measurements using
two different disc diameters [Smettem and Clothier, 1989] ormeasurements using a single disc diameter but with multiple
tensions [e.g., Ankeny et al., 1991].
Detailed studies of the water flow field below the tension
disc infiltrometer have been limited. In most studies the final
water content and/or concentration profiles were obtained by

2965

2966

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

sampling after terminating the infiltration experiment [e.g.,
Clothier et al., 1992; Quadri et al., 1994]. Clothier et al. [1992]
used knowledge of the steady state water and solute (Br⫺)
distributions to estimate the mobile and immobile water fractions. Quadri et al. [1994] used a similar experimental setup to
validate their finite difference numerical model for axisymmetrical movement of both water and solute underneath aquarter-sector disc permeameter. They obtained good predictions for the observed moisture and Br⫺ profiles in the sandbox
as well as for the measured infiltration rate from the disc
[Quadri et al., 1994].
Recent developments in TDR technology have provided vadose zone hydrologists and soil scientists with a tool for quick
and accurate measurement of the water content and the bulk
electric conductivity. Kachanoski et al. [1990] used bothstraight and curved TDR transmission lines to measure soil
water distributions as a function of time at different radial
distances from the surface water source during threedimensional infiltration in a laboratory setup. A similar setup
was used by Ward et al. [1994] to characterize both water flow
and solute transport underneath a disc infiltrometer. Recently,
Vogeler et al. [1996] used horizontally and vertically installed
TDR probes to characterize simultaneously water and solute
movement below a disc permeameter in a one-dimensional
setup and successfully modeled measured data with a numerical model. The above studies were all carried out in the laboratory. In contrast, Wang et al. [1998] measured water contents with TDRs and pressure heads with tensiometers in afield experiment to develop and test different approximate
analytical infiltration models in their inverse analysis.
Although early time infiltration data can be used to estimate
the sorptivity [White and Sully, 1987] and, consequently, thematrix flux potential, only steady state infiltration rates areusually used for Wooding-type analyses. We recently suggested
using the entire cumulative infiltration curve in combination
with parameter estimation to estimate additional soil hydraulicparameters [Šimůnek and van Genuchten, 1996, 1997]. From ananalysis of numerically generated data for one supply tension
experiment we concluded that the cumulative infiltration curveby itself does not contain enough information to provide a
unique inverse solution [Šimůnek and van Genuchten, 1996].
An infinite number of combinations of the saturated hydraulicconductivity, K s and the shape factor n (see (6) and (7)) canresult in almost identical infiltration curves. A similar conclusion was reached by Russo et al. [1991] for a one-dimensional
ponded infiltration experiment. Hence additional information
about the flow process, such as the water content and/or pressure head, measured at one or more locations in the soil profile
is needed to successfully obtain unique inverse solutions forthe soil hydraulic functions.
Šimůnek and van Genuchten [1997] studied, again numerically, infiltration at several consecutive supply tensions. Šimůnek and van Genuchten considered several different scenarioswith different levels of information and concluded that the bestpractical scenario is to estimate the hydraulic parameters fromthe cumulative infiltration curve measured at several consecutive tensions applied to the soil surface, in conjunction withknowledge of the initial and final water content. Our results
suggested that one should be able to use information typically
being collected with a tension disc infiltrometer to estimate not
only unsaturated hydraulic conductivities but also, without further experiments, the soil water retention properties.
The above methodology was tested on data collected as partof the soil hydrology program of the Hydrologic Atmospheric
Pilot Experiment (HAPEX)-Sahel regional-scale experiment[Cuenca et al., 1997]. Šimůnek et al. [1998a] showed that it isindeed possible to obtain relatively reliable estimates of soil
hydraulic conductivities for near-saturated moisture conditions
from tension disc infiltration data by numerical inversion, consistent with results of the traditional Wooding analysis. However, the question of whether it is possible also to estimate
simultaneously the retention curve was not resolved satisfactorily.
Infiltration rates effectively integrate properties of the porous media underneath the disc infiltrometer, including theinfluence of local-scale heterogeneity, different soil structure,and texture irregularities, preferential pathways, layering, andanisotropy; hence infiltration rates provide a good way for
estimating the effective near-saturated soil hydraulic properties. The question of whether a more detailed description of
the flow field, for example, in terms of point measurements ofthe pressure head and water content at several locations of thesoil profile, can improve the parameter estimation process andlead to more precise estimates of the effective soil hydraulicproperties arises. Our previous studies also did not resolve
several other questions concerning the best definition of theobjective function, such as (1) how the objective function
should be defined and which measurements should be included; (2) what the importance is of different types of data
and what these can do to improve the parameterization of soil
hydraulic properties; (3) whether or not one should include allavailable measurements in the definition of the objective function, such as pressure head measurements, water content measurements, or cumulative infiltration volumes, or only some(and which) optimal combination of these; (4) if and how a
combination of different sets of measurements can augmenteach other and what weight should be given to particular measurement sets so that one measurement set does not dominate
the others; and finally, (5) whether the weighting should beclosely related to expected data errors or should reflect the
expected importance of a particular measurement set to theparameterization of the soil hydraulic properties.
In this study we will apply the parameter estimation technique to the analysis of tension disc infiltrometer laboratory
experiments in which the soil profile is instrumented with bothtime domain reflectometry (TDR) (at one location) and tensiometers (at several locations) to measure water contents andpressure heads with time. The soil hydraulic characteristics
estimated using different sets of information will be comparedagainst results of Wooding’s [1968] analysis and against independently measured soil hydraulic properties obtained with theevaporation method.

2.

Experiment

The experiments were carried out using a loamy sand soilthat was sieved through a 4 mm mesh sieve. The soil was
carefully repacked in a 1 ⫻ 1 ⫻ 1 m3 soil container in 10 cm
layers at an initial soil bulk density of 1.5 g cm⫺3. After packing
in the laboratory the soil container was taken outside andexposed to natural rainfall and evaporation conditions for ⬃2.5
years. During this time all germinated weeds were removed byhand. In December 1996 the soil container was brought backinto the laboratory, the upper 10 cm soil layer was removed,
and the soil was carefully leveled. Pressure transducer tensiometers (measurement errors ⫾0.5 cm) were installed horizontally at different depths and at different distances relativeŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

2967

Figure 1. Schematic of tension disc infiltrometer setup, including time domain reflectometry (TDR) probe
and tensiometers T1 through T8.

to the center of the permeameter (Figure 1). A small hole forthis purpose was hand dug to the intended measurement depth
and then used for inserting tensiometer from the side. The
hole was far enough from the central vertical axis below thedisc to avoid soil disturbance within the flow domain. Afterinstallation the holes were filled with soil. The tensiometercups of 6 cm length and 0.6 cm outside diameter were isolated,with the exception of their 1 cm long tips. Sixteen tensiometers
were installed in pairs at the following coordinate (r, z) positions: (10, 2.5), (10, 5), (10, 10), (15, 5), (15, 10), (15, 15), (15,20), and (15, 30), where r represents the distance from the axisof symmetry and z is the location below the soil surface. A
two-wire TDR probe (measurement error ⫾0.01 m3 m⫺3) was
horizontally installed at a depth of ⬃2 cm directly below thetension disc. The TDR rods were 10 cm long and ⬃1.5 cm
apart. Prior to the disc infiltration experiment a plexiglass ringof 2 cm height was inserted in the soil to a depth of ⬃1.5 cm.
In order to avoid smearing of the soil surface and closing poresby surface smoothing a fine sand layer was placed on the soil
surface within the ring to ensure good contact between thetension disc and the soil. In our device the tension disc wasconnected via a PVC pipe to a water reservoir placed at theintended height relative to the soil surface in the container.
The tension disc was subsequently put on the soil surface, andthe water supply was opened. The infiltrometer disc had aradius of 10 cm. A schematic of the locations of the tensiometers, the TDR probe, and the infiltrometer disc is shown in
Figure 1.
We conducted four infiltration experiments, each with consecutively increasing supply pressure heads of ⫺20, ⫺10, ⫺5,and ⫺1 cm. The first three experiments lasted between 4.5 and5.5 hours and are referred to here as the “short” or “unsteadystate” experiments. The supply tension of ⫺1 cm was appliedfor ⬃25 min, and the pressure heads of ⫺5 and ⫺10 cm wereapplied for ⬃65 min. The last experiment lasted 6.25 days andwill be referred to it as the “long” or “steady state” experiment.
The first tension in this fourth experiment was applied fornearly 4 days, the second tension was applied for ⬃2 days, andthe last two tensions were applied for 6 and 1 hour, respectively. The intent of the fourth experiment was to reach steadystate at each applied tension so that Wooding’s [1968] analytical solution could be used. Tensiometer readings were taken
automatically every 30 s. We monitored the water level in theinfiltrometer, as well as the TDR, about every 5 minutes, except for the ⫺1 cm supply head where readings were obtained
every minute. Since the infiltration tubes did hold only ⬃2 L ofwater, the experiments had to be interrupted several times
during refilling. In order to capture also the redistribution partof the experiment, tensiometers were read for an additional 2hours after infiltration had ceased. A summary of the experimental conditions is given in Table 1. Between experiments thewater in the soil container was allowed freely to redistributeand evaporate through the unprotected soil surface to recreate
approximately similar initial conditions.
Ten undisturbed soil core samples with a height of 6 cm andan inside diameter of 8 cm were used independently to measure the soil hydraulic properties. From the 10 samples, fivewere taken at 2– 8 cm depth, and the five at 10 –16 cm depth
were taken at a distance of at least 50 cm away from theinfiltration disc. The average bulk density of the samples was
1.41 g cm⫺3. Unsaturated soil hydraulic conductivities close tosaturation (at pressure heads of ⫺1, ⫺5, and ⫺10 cm) were
determined using steady state downward infiltration experiments with two tension disc permeameters [Wendroth andŠimůnek, 1999]. One infiltrometer disc was placed at the top
and another one at the bottom of the soil sample. The same

2968

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

Table 1. Summary of Experimental Conditions for the Four Tension Disc Infiltration Experiments

ExperimentI

SupplyTensions
20 cm
10 cm
5 cm
1 cm
Total
20 cm
10 cm
5 cm
1 cm
Total
20 cm
10 cm
5 cm
1 cm
Total
20 cm
10 cm
5 cm
1 cm
Total

II

III

IV

Infiltrated
Volume,L

Experiment Duration
Days

Hours

Minutes

2
1
1

39
10
5
24
18
7
6
6
20
39
3
5
6
25
39

5
2
1
1
4
2
1
1
4
23

3
2

15
48
2
5

5
1
6

6

0.084
0.522
1.213
4.304
6.12
0.155
0.647
1.502
4.224
6.53
0.213
0.470
1.549
4.465
6.70
5.92
10.12
4.15
14.33
34.5

pressure head was applied to both discs in order to establish
and maintain steady state flow conditions under a unit gradientand from this to obtain unsaturated hydraulic conductivities.
Following the infiltrometer experiments, the samples weremoved onto impermeable plates to determine water retention
and unsaturated hydraulic conductivity data for pressure heads⬍⫺10 cm using two-rate evaporation experiments [Wendroth etal., 1993; Šimůnek et al., 1998b]. The resulting data were analyzed using a modified Wind [1968] method [Wendroth et al.,1993]. After conclusion of the evaporation experiments thewater contents at pressure heads of ⫺10 and ⫺150 m were
measured in a pressure chamber on selected 0.7 cm long subsamples from the disturbed cores.

3.

Theory3.1.

Numerical ModelIn our analysis of the tension disc infiltrometer data we will
use a numerical solution of the Richards’ equation coupledwith the Levenberg-Marquardt nonlinear minimization
method [Marquardt, 1963]. The governing flow equation forradially symmetric isothermal Darcian flow in a variably saturated isotropic rigid porous medium is given by the following
modified form of the Richards’ equation:冉⭸␪ 1 ⭸
⭸h
rK
⫽
⭸t
r ⭸r
⭸r

冊 冉 冊
⫹

⭸
⭸h
K
⭸r
⭸z

⫺

⭸K
,
⭸z

(1)

where ␪ is the volumetric water content [L 3 L ⫺3 ], h is thepressure head [L], K is the hydraulic conductivity [LT ⫺1 ], r isa radial coordinate [L], z is the vertical coordinate [L] positive downward, and t is time [T]. Equation (1) was solved
numerically for the following initial and boundary conditions
applicable to a disc tension infiltrometer experiment:␪ 共r, z, t兲 ⫽ ␪ i共 z兲

t ⫽ 0,

h共r, z, t兲 ⫽ h i共 z兲

t ⫽ 0,

h共r, z, t兲 ⫽ h 0共t兲

0 ⬍ r ⬍ r0(2)
z ⫽ 0,(3)

Initial Pressure Heads, cm
2.5 cm Depth

Final Infiltration
Rates, cm3 s⫺1

30 cm Depth⫺493

0.00330
0.0767
0.379
3.81

⫺200

⫺498

0.0154
0.122
0.436
2.67

⫺210

⫺336

0.0291
0.0984
0.442
3.86

⫺182

⫺641

0.01430.0586
0.224
4.04

⫺222

⭸h共r, z, t兲
⫽1⭸zz ⫽ 0,r ⬎ r0h共r, z, t兲 ⫽ h i

r 2 ⫹ z 2 3 ⬁,
3

(4)
(5)

⫺3

where ␪ i is the initial water content [L L ], h i is the initialpressure head [L], h 0 is the time variable supply pressure head
imposed by the tension disc infiltrometer [L], and r 0 is the discradius [L]. Equation (2) specifies the initial condition in termsof either the water content or the pressure head. Boundary
condition (3) prescribes the time variable pressure head under
the tension disc permeameter, while (4) assumes a zero flux atthe remainder of the soil surface (evaporation is neglectedduring the short-duration infiltration experiments). Equation(5) states that the other boundaries are sufficiently distantfrom the infiltration source so that they do not influence theflow process. The boundary condition at the axis of symmetry(r ⫽ 0) is a no flow condition. Equation (1), subject to theabove initial and boundary conditions, was solved using a quasi
three-dimensional (axisymmetric) finite element code, HYDRUS-2D, as documented by Šimůnek et al. [1996].
3.2.

Soil Hydraulic Properties

A model of the unsaturated soil hydraulic properties must be
selected prior to application of the numerical solution of the
Richards’ equation. In this study we will limit ourselves to theunsaturated soil hydraulic functions [van Genuchten, 1980]:
S e共h兲 ⫽␪ 共h兲 ⫺ ␪ r1⫽
␪s ⫺ ␪r
共1 ⫹ 兩 ␣ h兩 n兲 m

m 2
K共 ␪ 兲 ⫽ K sS le关1 ⫺ 共1 ⫺ S 1/me 兲 兴 ,

(6)
(7)

where S e is the effective fluid saturation [ ], K s is the saturated hydraulic conductivity [LT ⫺1 ], ␪ r and ␪ s denote theresidual and saturated water contents [L 3 L ⫺3 ], respectively; lis the pore connectivity parameter [ ], and ␣ [L ⫺1 ], n [ ],and m (⫽ 1 ⫺ 1/n) [ ] are empirical shape parameters. The
pore connectivity parameter l in K( ␪ ) was estimated by
Mualem [1976] to be 0.5 as an average for many soils. Taking

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTSl ⫽ 0.5, the above hydraulic functions contain five unknown
parameters: ␪ r , ␪ s , ␣ , n, and K s . We will refer to (6) and (7)as the van Genuchten-Mualem (VGM) model.
3.3.

Table 2a. Unsaturated Hydraulic Conductivities Calculated
Using Wooding’s [1968] Analysis
Unsaturated Hydraulic Conductivity, cm s⫺1
Experimenth ⫽ ⫺15 cm

h ⫽ ⫺7.5 cm

h ⫽ ⫺3 cm

I
II
III
IV
Arithmetic mean

0.0000361
0.0000855
0.0000833
0.0000485
0.0000634

0.000388
0.000489
0.000466
0.000247
0.000398

0.003130.00268
0.00337
0.00257
0.00294

Formulation of the Inverse Problem

The objective function ⌽ to be minimized during the parameter estimation process can be formulated in terms of an arbitrary combination of cumulative infiltration data, TDRmeasured water contents, and/or tensiometer (pressure head)readings. The objective function is defined as [Šimůnek andvan Genuchten, 1996]

冘再 冘m⌽共␤, q m兲 ⫽

j⫽1

冎

nj

vj

w ij关q *j共t i兲 ⫺ q j共t i, ␤兲兴 2 ,

i⫽1

(8)

where m represents different sets of measurements (infiltration data, pressure heads, and/or water contents), n j is thenumber of measurements in a particular set, q *j (t i ) is the specific measurement at time t i for the jth measurement set, ␤ isthe vector of optimized parameters (e.g., ␪ r , ␪ s , ␣ , n, K s , andl ), q j (t i , ␤) represents the corresponding model predictions
for parameter vector ␤, and v j and w ij are weights associatedwith a particular measurement set j or a measurement i within
set j, respectively. We assume that the weighting coefficientsw ij in (8) are equal to 1; that is, the variances of the errors
inside a particular measurement set are all the same. The
weighting coefficients v j are given by
vj ⫽

1
.
n j␴ 2j(9)

The above approach views the objective function as the average weighted squared deviation normalized by measurement
variances ␴ 2j . Minimization of the objective function ⌽ is accomplished by using the Levenberg-Marquardt nonlinear minimization method [Marquardt, 1963].
3.4. Wooding’s [1968] Analysis
The traditional analysis of tension disc infiltration data
based on Wooding’s [1968] analytical solution requires twosteady state fluxes at different tensions [Ankeny et al., 1991] toyield estimates of the saturated hydraulic conductivity K s[LT ⫺1 ] and the sorptive number ␣* [L ⫺1 ] in Gardner’s [1958]
exponential model of the unsaturated hydraulic conductivity:
K共h兲 ⫽ K s exp 共 ␣ *h兲,

(10)Wooding’s solution for infiltration from a circular source witha constant pressure head at the soil surface, and with theunsaturated hydraulic conductivity described by (10), is givenby
Q共h 0兲 ⫽ ␲ r 02K共h 0兲 ⫹4r 0
K共h 0兲.
␣*(11)

where Q is the steady state infiltration rate [L 3 T ⫺1 ], r 0 is theradius of the disc [L], h 0 is the wetting pressure head [L], andK(h 0 ) is the unsaturated hydraulic conductivity [LT ⫺1 ] atpressure head h 0 . Methods for obtaining the unsaturated hydraulic conductivity in the middle of an interval between two
successively applied pressure heads are given by Ankeny et al.
[1991], Reynolds and Elrick [1991], and Jarvis and Messing[1995], among others.

2969

4.

Results and Discussion

4.1. Wooding’s [1968] Analysis
Experimental data were first analyzed using Wooding’s[1968] analytical solution (see (11)). Wooding’s analysis requires steady state infiltration rates at different supply pressure
heads. Depending upon soil texture it can take hours or evendays to reach steady state in a field experiment. Previous studies have shown that Wooding’s approach will overestimate thesoil hydraulic conductivity if steady state infiltration is notreached. Nevertheless, infiltration rates reached within 1 hour
are assumed to be the steady state rates in a majority of studies
and are used in Wooding’s analysis. The possible error is usually dismissed as being negligible as compared to the effects ofsoil heterogeneity. Another reason for dismissing errors
caused by not reaching steady state infiltration is the lack of
reproducibility of infiltration experiments. For our soil, differences between the final infiltration rates for the same type ofexperiments (i.e., the short experiments I–III) were larger than
the increased precision resulting from the use of steady state
infiltration rates as obtained with the long experiment (experiment IV) at particular supply tensions (Table 1). Actually, thefinal infiltration rates for the long experiment were sometimes
higher than the final infiltration rates for the short experiments. The problems of reproducibility of unsaturated soil flow
experiments was also raised by Cı́slerová et al. [1988], whofound that steady state ponded infiltration rates in a coarse
acid brown soil depended strongly on the initial moisture content. They contributed this effect to entrapped air that alters
the volume of available pores for flow. Such an effect of theinitial water content on infiltration rate was not evident fromour data. Hollenbeck and Jensen [1998] studied transient outflow laboratory experiments and also found that experimentswith small changes in the imposed boundary conditions could
not be reproduced well even when great care was taken to
reproduce the same initial conditions. Although experimentswith large step changes were reproduced better, the overallsystem response was virtually the same for all step levels [Hollenbeck and Jensen, 1998].
Estimated unsaturated hydraulic conductivities at pressure
heads in the middle of the imposed disc supply intervals arepresented in Table 2a. As expected, the calculated conductivities, with one exception at h ⫽ ⫺15 cm, were lowest for the
fourth experiment, where steady state infiltration rates were
presumably reached at each supply tension. The smallest deviations (of ⬃15%) between the calculated conductivities andtheir arithmetic means were for h ⫽ ⫺3 cm, while the largest
deviations, ⬃75%, occurred for h ⫽ ⫺15 cm. Overall, thedifferences between unsaturated hydraulic conductivities calculated using Wooding’s analysis at any particular supply pressure head were relatively small. Also, differences between K s

2970

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

Table 2b. Saturated Hydraulic Conductivities and Sorptive
Numbers ␣* Calculated Using Wooding’s [1968] Analysis
Saturated Hydraulic Conductivity,
cm s⫺1/Sorptive Number, ␣* cm⫺1
Experimenth(⫺20, ⫺10)h(⫺10, ⫺5)h(⫺5, ⫺1)

I
II
III
IV
Arithmetic mean

0.00404/0.315
0.00191/0.207
0.000517/0.122
0.000401/0.141
0.00172/0.196

0.00425/0.319
0.00330/0.254
0.00443/0.300
0.00185/0.268
0.00346/0.285

0.0177/0.578
0.0104/0.453
0.0171/0.542
0.0225/0.723
0.0169/0.574

values extrapolated from the same supply tension intervals fordifferent experiments were relatively small (Table 2b). The
smallest differences occurred for those calculated using thelargest pressure heads, h(⫺5, ⫺1), while the largest differences (almost 1 order of magnitude) arose for those using thelowest pressure heads, h(⫺20, ⫺10). Saturated hydraulicconductivities calculated using the largest supply tension interval were almost 1 order of magnitude lower than those obtained for the smallest tension interval (Table 2b). Similarly,
the sorptive number ␣* was much lower for the largest interval
than for the smallest one. Because the predicted values of ␣*for different tension intervals differed substantially, applicability of Gardner’s [1958] model for extrapolation of the unsaturated hydraulic conductivity function beyond a particular supply
tension interval is questionable, at least for the loamy sandused in our study.
4.2.

Numerical Inversions

The measured data were next analyzed using parameter
estimation by combining the HYDRUS-2D numerical code
[Šimůnek et al., 1996] with the Levenberg-Marquardt nonlinear
minimization method [Marquardt, 1963]. The computational
domains (r, z) for the numerical solution of either 50 ⫻ 50 cm2
or 100 ⫻ 100 cm2 were discretized into structured rectangular
finite element mesh grids containing 15 ⫻ 19 and 24 ⫻ 33 meshnodes for the first three experiments and the last experiment,
respectively. A finer discretization was used near the soil surface than at depth. A much larger transport area for the fourthexperiment was used because of the longer duration of theexperiment and the larger cumulative infiltrated volume.
As the initial condition, we used the arithmetic mean of thepressure head for all functioning tensiometers at a particulardepth (maximum 4) and linear interpolation between the averaged pressure heads at those depths. In one special case we
assumed that the initial condition was constant for the entire
soil profile and equal to the measured initial water content at
a depth of 2 cm. The assumption of homogeneity in the initialcondition likely would introduce some error in the calculations
since the initial condition was not constant (see Table 1). The
initial water content in the immediate vicinity of the tension
disc, together with the supply tension, defines the soil sorptivity
and the curvature of the cumulative infiltration curve. This
curvature is generally highest at the beginning of infiltration
experiment when the capillary forces are dominating the infiltration process. The initial condition farther away from theinfiltrometer disc has far less effect on the shape of the cumulative infiltration curve; it affects primarily the speed of theinfiltration front through the soil profile. Our ultimate objective here was to estimate the soil hydraulic parameters fromthe least amount of information possible. Although a some-what incorrect or incomplete description of the initial condition may not yield accurate descriptions of the flow field in thesoil, errors in terms of the calculated cumulative infiltration
curve may only be minimal.
The measured data were analyzed in several different ways.
We first used four different formulations for the objectivefunction. Cumulative infiltration volumes I(t), which represent an integral characteristic of the soil profile, were includedin the definition of the objective function for all optimizations.
We additionally used either all other measured data simultaneously or combined the cumulative infiltration volumes witheither only the pressure head or water content data. In one
case we defined the objective function in terms of the measured cumulative infiltration and the final water content asmeasured with the TDR at 2 cm depth. This measured final
water content ␪ f was assumed to define one point of the retention curve together with the final supply pressure head (⫺1cm). This optimization option, in combination with the initialcondition being defined in terms of the water content, represents the optimization scenario that was suggested by Šimůnek
and van Genuchten [1997] as being the best practical setup for
predicting the soil hydraulic parameters from a tension infiltration experiment.
Next we optimized the soil hydraulic parameters in the
VGM model with and without assuming hysteresis according
to Kool and Parker [1987] (referred to here as the VGM-H
model). In a final scenario we also optimized the pore connectivity number l (VGM-Hl model). Below we will concentrate
mainly on experiment III, representing the short experiments,and steady state experiment IV. The soil hydraulic parameters
obtained with different optimization options for experiments
III and IV are summarized in Tables 3 and 4, respectively. In
our analysis we used 44 water content, 22,752 pressure head,
and 62 cumulative infiltration data points for experiment III
and 91 water content, 12,564 pressure head, and 100 cumulative infiltration data points for experiment IV.
4.2.1. Experiment III. Selected results for different parameter optimization options for experiment III (the shortduration infiltration experiment that lasted ⬃4.5 hours) aresummarized in Table 3. In this table, ⌽ *l , ⌽ *h , and ⌽*␪ represent contributions to the total value of the objective function ⌽of the weighted residuals between measured and optimized
values of the cumulative infiltration, pressure head, and watercontent measurement sets, respectively. Symbols in front of asemicolon within the argument of ⌽ (Table 3) represent measurement sets used in the optimization, while parameters listed
after the semicolon define the optimized unknowns. The best
fit for any particular measurement set (i.e., cumulative infiltration, water content, or pressure heads data) and the best overall fit as characterized by the highest R 2 and the lowest valueof the objective function ⌽ are shown. For comparison the lastcolumn of the Table 3 gives results obtained with a direct(forward) numerical solution using soil hydraulic parameters
fitted to independently measured hydraulic data obtained withthe evaporation method, the steady state unsaturated verticalflow measurements, and the pressure chamber.
Notice that only ␪ s was optimized when the objective function was defined in terms of cumulative infiltration data andthe tensiometer readings. Optimizations can lead to unrealistic
values of the optimized parameters ␪ r and/or ␪ s when no
information about the water content is included in the objective function and the initial condition is given in terms of thepressure head. Only the water content range can then be op-

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

2971

Table 3. Summary of Parameter Estimation Results for Experiment III
Objective
Function⌽(I, ␪ f ; ␪ r , ␪ s ,␣ w , n, K s )⌽(I, ␪ ; ␪ r , ␪ s ,␣ w , n, K s )⌽(I, h; ␪ s ,␣ w , n, K s )⌽(I, h, ␪ ; ␪ r , ␪ s ,␣ w , n, K s )⌽(I, h, ␪ ; ␪ r , ␪ s ,␣ w , n, K s , ␣ d )⌽(I, h, ␪ ; ␪ r , ␪ s ,␣ w , n, K s , ␣ d , l )DirectSolution

⌽
⌽ *I
⌽ *h
⌽ *␪
␪r
␪s
␣ w , cm⫺1n
K s , cm min⫺1
␣ d , cm⫺1l
R2

0.00526
0.00502䡠䡠䡠
0.000240
0.000
0.325
0.231
1.46
1.79
䡠䡠䡠
0.5b
0.998

0.185
0.0175
䡠䡠䡠
0.167
0.001
0.324
0.206
1.43
1.42
䡠䡠䡠
0.5b
0.992

0.152
0.00580.146a䡠䡠䡠
0.000b
0.246
0.196
2.07
0.623
䡠䡠䡠
0.5b
0.953

0.390
0.0046
0.195
0.191
0.015
0.334
0.257
1.44
1.98
䡠䡠䡠
0.5b
0.956

0.379
0.0048
0.184
0.191
0.016
0.332
0.263
1.42
2.11
0.147
0.5b
0.959

0.303a0.0046a0.152
0.146a
0.002
0.310
0.171
1.43
1.59
0.111
1.97
0.966a

3.858
0.557
0.604
2.697
0.0761b0.348b䡠䡠䡠
1.46b0.105b
0.0347b
0.5b
0.775

a
b

Best fit for any particular measurement set and the best fit as characterized by the highest R2 and the lowest value of the objective function ⌽.
Parameters fixed during a particular optimization; ␣ w and ␣ d are shape parameters for the wetting and drying cycle, respectively.

timized since the residual and saturated water contents aremutually correlated. For this reason we fixed ␪ r equal to zeroand optimized only ␪ s . Although this approach provided thebest overall fit of the tensiometer readings, the optimization
resulted in unrealistically low values of ␪ s (or the water content
interval ␪ s ⫺ ␪ r ). Also, the parameter n differed significantly
from those obtained with the other optimizations.
Figure 2 shows measured and calculated (fitted) water contents and cumulative infiltration volumes for the third experiment. The VGM soil hydraulic parameters were optimized
using the objective function defined in terms of all three measured variables, i.e., cumulative infiltration volumes, tensiometer measurements, and TDR measurements (Table 3, fifth
column). Notice the excellent agreement between measured
and calculated values for cumulative infiltration. While more
scatter appears in the measured water content values (Figure2b), the final fit correctly reflects the water content pattern.
The predicted water contents reach the final value imposed byeach supply tension faster than the measured values, which
tend to increase more gradually with time and do not quite
reach constant values during the short-duration experiment.
This may be partly due to space averaging by the TDR probe
as compared to the use of nodal (point) values in the numericalsolution. We will discuss this issue later in more detail. Figure
2b shows a small underprediction of the water content for
supply tensions of 20, 10, and 5 cm. Notice from Table 3 thatthe results presented in Figure 2 do not represent the best fitof the measured water contents and that significantly better fitswere obtained with other combinations of optimized parameters.
Figure 3 shows measured and calculated pressure heads ateight locations below the disc infiltrometer. Since a few tensiometers were not working correctly during the entire experiment, we used only one tensiometer reading for locations (10,10), (15, 15), (15, 20), and (15, 30). The differences in pressure
head readings at two tensiometers located at the same radial
coordinates, for example, for (10, 2.5) and (15, 5), reflect local
soil heterogeneity. Almost identical readings were obtained forthe duplicate readings at coordinates of (10, 5) and (15, 10).
The calculation underpredicted the arrival of the water front at
both (10, 2.5) and (10, 5) in the immediate vicinity of thetension disc. Agreement was better at location (10, 10) deeper
in the profile, while the agreement was excellent at depths of
10 and 20 cm when r was equal to 15 cm.
We continued the measurements for about 2 more hours
after we terminated the infiltration of water. The ensuing redistribution part of the experiment could not be described wellwith the VGM model that ignored hysteresis (Figure 3, time ⬎340 min). Figure 4 shows measured and calculated pressure
heads at the same eight tensiometer locations when an additional soil hydraulic parameter ␣ d was included into the optimization (Table 3, sixth column) to account for hysteresis. For
the hysteretic case we used the model of Kool and Parker[1987], who coupled the VGM model with the simplified scal-

Table 4. Summary of Parameter Estimation Results for Experiment IV
Objective
Function⌽(I, ␪ f ; ␪ r , ␪ s ,␣ w , n, K s )⌽(I, ␪ ; ␪ r , ␪ s ,␣ w , n, K s )⌽(I, h; ␪ s ,␣ w , n, K s )⌽(I, h, ␪ ; ␪ r ,␪ s , ␣ w , n, K s )⌽(I, h, ␪ ; ␪ r , ␪ s ,␣ w , n, K s , ␣ d )⌽(I, h, ␪ ; ␪ r , ␪ s ,␣ w , n, K s , ␣ d , l )DirectSolution

⌽
⌽ *I
⌽ *h
⌽ *␪
␪r
␪s
␣ w , cm⫺1n
K s , cm min⫺1
␣ d , cm⫺1l
R2

0.01310.0130䡠䡠䡠
0.0001
0.105
0.302
0.201
1.56
1.01
䡠䡠䡠
0.5b
0.995

0.107
0.0241䡠䡠䡠
0.0830
0.127
0.307
0.219
1.44
1.08
䡠䡠䡠
0.5b
0.991

0.244
0.00034a0.243a䡠䡠䡠
0.000b
0.470
0.854
1.22
34.4
䡠䡠䡠
0.5b
0.991

0.424
0.09080.272
0.0609
0.001
0.306
0.134
1.30
0.463
䡠䡠䡠
0.5b
0.961

0.348
0.04660.252
0.0495
0.017
0.332
0.197
1.30
0.764
0.129
0.5b
0.975

0.339a0.03720.252
0.0490a
0.025
0.332
0.228
1.28
1.13
0.150
0.106
0.979a

30.18
25.42
0.70
4.06
0.0761b0.348b䡠䡠䡠
1.46b0.105b
0.0347b
0.5b
0.901

a
b

Best fit for any particular measurement set.
Parameters fixed during a particular optimization; ␣ w and ␣ d are shape parameters for the wetting and drying cycle, respectively.

2972

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

Figure 2. Measured and fitted (a) cumulative infiltration volumes and (b) water contents for the short-duration experimentIII. Fitted curves were obtained using optimized van Genuchten-Mualem (VGM) parameters (see Table 3, fifth column).
ing approached of Scott et al. [1983] to describe the scanning
curves. Scott et al. [1983] assumed that the shape parameters
for all drying scanning curves are the same as those for themain drying curve and, similarly, the shape parameters for all
wetting scanning curves are the same as those for the main
wetting curve. Scanning curves are then calculated by varying
the residual and saturated water contents for the wetting anddrying scanning curves, respectively. Kool and Parker [1987]
also assumed that the shape parameter n is the same for bothwetting and drying, thus decreasing the number of required
parameters. Using the additional restrictions that ␪ r and ␪ s are
the same for both drying and wetting, the shape parameter ␣ dfor the drying retention curve ( ␣ w is used for wetting) is the
only additional parameter describing hysteresis. An advantageof the model by Scott et al. [1983] is that one can obtain the
shape parameters of the drying and wetting curves from knowledge of any two main or primary, secondary, or higher-orderscanning drying and wetting curves. While the infiltration partof the experiment in Figure 4 is similar to that in Figure 3, asignificantly better description of the redistribution part of theexperiment is obtained with the simple VGM-H hysteresis
approach. The optimized soil hydraulic parameters are very
similar to those for Figure 3 (see Table 3). Interestingly, theoptimized value of ␣ d (0.147 cm⫺1) is close to ␣ w / 2 (⫽ 0.132cm⫺1), an approximation often thought to be reasonable for simlating hysteresis [Kool and Parker, 1987; Luckner et al., 1989].
Overall, all optimizations that included in the objective function either water content data or only the initial and final
values of the water content at a depth of 2 cm provided similar
sets of optimized parameters. The use of pressure head data
did not significantly improve the results. The main benefit ofusing pressure head measurements is the obtaining of information about both infiltration and redistribution and, consequently, about hysteresis in the soil hydraulic properties.
4.2.2. Experiment IV. Parameter estimation results forexperiment IV, the 6 day slow infiltration experiment for which
steady state was reached for each supply pressure head, aresummarized in Table 4. The best fit of cumulative infiltrationand pressure head readings was obtained when the objectivefunction did not include water content measurements. The
predicted and measured cumulative infiltration volumes were
essentially identical (figure not shown). Nevertheless, the estimated soil hydraulic parameters for this optimization were, as
for experiment III, unrealistic, producing a high saturated water content and very high values of ␣ and K s (Table 4). These
results suggest that the optimization requires some additional
information about the water content. Additional information isusually available for transient laboratory experiments (such as
outflow or evaporation experiments) in terms of the total water
volume in the sample but generally not for field experiments.
For example, Inoue et al. [1998] analyzed a multistep soil water
extraction field experiment and obtained unrealistically highestimates of the residual water content ( ␪ s was fixed) whenusing only pressure head data and cumulative extraction rates
despite having an excellent fit of these data. Estimation of theretention curve dramatically improved only when water content data were included in the optimization [see Inoue et al.,1998, Figure 10]. Similarly, estimates of ␪ s ( ␪ r was fixed) werevery unreliable when only pressure head and cumulative infiltration data were analyzed in a modified cone penetrometer
experiment [Kodešová et al., 1998; Šimůnek et al., 1999]. Fixing
␪ s or including water content data in the latter example improved the estimation. Consistent with these results, Šimůnek
and van Genuchten [1997] found out that the identifiability ofsoil hydraulic parameters from numerically generated disc infiltrometer data dramatically improved when the initial condition was expressed in terms of the water content instead of thepressure head.
Figure 5 shows measured and calculated water contents andcumulative infiltration rates for the overall best fit, i.e., when
the objective function was defined in terms of all measured
variables (Table 4, seventh column). Some interruptions occurred in the collection of cumulative infiltration and watercontent data during the nights; however, tensiometer readings
continued since they were logged in automatically. The watercontent again gradually increased with time for each supply
pressure head during the entire experiment and did not reacha constant value. The water content increased long after thenumerical model reached a constant water content and in spiteof the infiltration rates being already constant. Again, as forthe third example, this could be caused by water content averaging over a relatively large volume when using the TDR
probe, while the numerical solution only assumed nodal values.
An equally or more probable explanation is that water moves
only gradually into relatively small intra-aggregate or dead-end
pores that do not significantly contribute to water flow. Water
initially may move primarily through the larger pores (depend-

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

2973

Figure 3. Measured and fitted pressure heads at eight locations below the disc infiltrometer for the shortduration experiment III. Calculated curves were obtained using optimized VGM parameters (see Table 3, fifth
column).
ing upon the supply tension) that provide the main paths forwater flow. These larger pores are the most important in termsof defining the overall flow process and determining the tensiometer pressure heads; water hence may bypass the smaller
pores, which are vulnerable to air entrapment. This situation
would cause the pressure head and the infiltration rates toreach apparent steady state conditions much faster than thewater content. Notice from Figure 6 that the measured pressure heads indeed reached equilibrium conditions faster than
the numerically calculated pressure heads. A similar water flow
nonequilibrium situation was observed by Wang et al. [1998],

who measured increases in the water content under the tension
infiltrometer long after tensiometers had reached equilibrium
at the same spatial coordinates. While the agreement between
measured and fitted cumulative volumes in Figure 5 was worse
than the best fit mentioned above, the match was still relatively
good. For the optimized parameters the numerical code overpredicted and underpredicted the infiltration rates for supply
tensions of 10 and 1 cm, respectively.
Figure 6 shows again the measured and calculated pressure
heads at eight tensiometer locations. Inserts give details about
the pressure head at early times. Notice again an underpredic-

2974

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

Figure 4. Measured and fitted pressure heads at eight locations below the disc infiltrometer for the shortduration experiment III. Calculated curves were obtained using optimized VGM parameters with hysteresis(VGM-H) (see Table 3, sixth column).
tion of the moisture front arrival at a depth of 2.5 cm and a
relatively good agreement at other locations having a radial
coordinate of 10 cm. There is an especially good correspondence between the measured and calculated pressure heads forlocations with a radial coordinate of 15 cm at depths of 5 and10 cm, whereas the model predicts later arrivals of the moisture front to depths of 20 and 30 cm. Once again, all optimizations that included in the objective function either watercontent measurements or only the initial and final values ofwater content at a depth of 2 cm provided similar sets of
optimized parameters.
4.3.

Measured and Optimized Soil Hydraulic Properties

Soil water retention data determined with the evaporation
method for h ⬎ ⫺700 cm and in a pressure chamber for h ⫽⫺10 and ⫺150 m for 10 samples are presented in Figure 7a.
Corresponding hydraulic conductivities estimated using theevaporation method for h ⬍ ⫺10 cm and those obtained using

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

the tension disk infiltrometer for h ⫽ ⫺1, ⫺5, and ⫺10 cm
are shown in Figure 7b. Overall, the results obtained for all 10
samples using the evaporation method show a high level ofconsistency. Also, the K(h) relations in the wet range estimated from the infiltration experiments correspond reasonably
well with those in the drier range determined using the evaporation method. The soil hydraulic parameters were obtainedby simultaneously fitting (6) and (7) to all ␪ (h) and K(h) datausing the retention curve (RETC) code [van Genuchten et al.,1991]. We obtained an excellent fit when the parameter l was
fixed at 0.5 for the evaporation data (Figure 7, fit 1). The fit,
however, overpredicted water contents at the higher tensions
and underpredicted the hydraulic conductivity at a pressure
head of ⫺1 cm. A summary of the obtained soil hydraulicparameters, together with R 2 , is given in Table 5. A better
description of the retention curve at higher tensions was obtained when the parameter l was also considered to be an
unknown. Since the first fit provided a better description of thesoil hydraulic data within the experimental range for our laboratory tension disc experiments, we used these parameters forthe direct solutions given in Tables 3 and 4.
Figure 8 shows retention curves obtained with the parameter
estimation method for experiments III and IV and, for comparison purposes, results produced with the evaporation
method. Even though different optimization options wereused, the resulting retention curves for each experiment arevery similar; they define a relatively narrow band and converge
to approximately the same saturated water content. Notice thatwe plotted both the wetting and drying curves when hysteresis
was considered. The retention curves estimated from experiment IV gave higher water contents over nearly the entire hrange than those estimated from experiment III. The greatest
difference between the two tests is ⬃0.04 and occurs at apressure head of ⬃⫺100 cm. Although all retention curves
have similar shapes, the curves obtained by parameter estimation do underestimate water contents at particular pressure
heads relative to curves based on the evaporation method. This
underestimation is more important for retention curves determined from the short-duration experiment III than for the
long-duration experiment IV. The average difference between
retention curves determined from experiment IV and the evaporation method is again ⬃0.04 cm3 cm⫺3. This difference canbe explained partly by hysteresis, partly by entrapped air, andpartly by observed nonequilibrium behavior.
Since the soil samples used for the evaporation method were
initially saturated with deionized water, the resulting retention
curve represents the main drying curve with ␪ s approximately
equal to porosity. It is generally accepted that the fieldsaturated (or satiated) water content is much smaller thanporosity because of entrapped air. Part of the difference canalso be explained by the presence of flow irregularities anddeviations from equilibrium flow theory (such as a gradually
increasing water content even when the infiltration rate andthe pressure head reached steady state). For both experiments
the water content was 0.185 after 2 hours of infiltration. For
experiment IV the water content increased another 0.044 units
before the end of the tension infiltration step (⫺20 cm). The
water content increased ⬃0.025 even after tensiometers located at deeper depth than TDR reached a constant value.
Since the water content measurements are included in thedefinition of the objective function, but somehow lag in timebehind the pressure heads, the optimization routines must
underestimate the water content at particular pressure heads.

2975

Figure 5. Measured and fitted (a) cumulative infiltration volumes and (b) water contents for the long-duration experiment IV.
Calculated curves were obtained using optimized VGM parameters with hysteresis (VGM-Hl) (see Table 4, seventh column).

Figure 9 shows hydraulic conductivities determined directly
using the evaporation method, those measured during thesteady state infiltration, and estimated conductivity functions
assuming different optimization options for experiments III
and IV. The optimized hydraulic conductivities fall in the samerange as those measured during steady state infiltration andthose calculated using Wooding’s [1968] analysis. Good correspondence between conductivities obtained using Wooding’s
approach and the inverse results were expected since bothmethods analyze the same data. While hydraulic conductivities
measured during steady state infiltration correspond closelywith the K(h) functions estimated from experiment IV, Wooding’s conductivities are closer to the results from experimentIII. The hydraulic conductivities estimated from experiment III
were higher than those estimated from experiment IV at thelower tensions and tended to converge toward higher tensions.
The average K s estimated from experiments III and IV were⬃1.75 and 0.9 cm min⫺1, respectively. The evaporation
method provided slightly higher values of the conductivities atparticular tensions than those estimated by numerical inversion. Nevertheless, the drying branches of the hydraulic conductivity functions are very close to point values obtained withthe evaporation method.

5.

Concluding Remarks

The parameter estimation technique proved to be a powerful method for analyzing transient water flow experiments.

2976

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

Figure 6. Measured and fitted pressure heads at eight locations below the disc infiltrometer for the longduration experiment IV. Calculated curves were obtained using optimized VGM parameters with hysteresis(VGM-Hl) (see Table 4, seventh column). Inserts give details of pressure heads at the early time.

Unsaturated hydraulic conductivities obtained with parameter
estimation corresponded well with results of Wooding’s [1968]
analysis as well as with those determined with steady state
infiltration. Drying branches of the K(h) function determinedby parameter estimation corresponded relatively well with hydraulic conductivities obtained by the evaporation method.
Retention curves estimated from the infiltration experiments
were quite different from those determined using the evapo-

ration method. They, nevertheless, may prove to be more useful in describing infiltration and the transport of contaminants
through the vadose zone than retention curves determined bysteady state methods or from transient processes of a totallydifferent nature.
An advantage of axisymmetric experiments, such as those
discussed in this paper, is that one can have several measurements at locations with the same radial coordinates. However,ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

2977

Figure 8. Soil water retention curves determined using different optimization options and data collected for experiments
III and IV. Retention data were measured using the evaporation method and the pressure chamber.

Figure 7. (a) Soil water retention curves determined with theevaporation method for h ⬎ ⫺700 cm and as determined witha pressure chamber for h ⫽ ⫺10 and ⫺150 m and (b) hydraulic conductivity functions determined with the evaporation
method for h ⬍ ⫺10 cm and using a tension infiltration
method for h ⫽ ⫺1, ⫺5, and ⫺10 cm.

fronts for both cases would probably result in fairly similar
optimized shape parameters ( ␣ , n), while the optimized valueof K s could be quite different between the two optimizations.
The arithmetic or geometric mean of these two values could
then be taken as the effective hydraulic conductivity of the soil.
One could also use mean value of the measured variable, or
one could define the objective function in terms of both measurements, in which case the optimization probably will resultin the same mean hydraulic conductivity but with widely different optimized shape parameters ␣ and n. The optimized
parameters will then predict a much more gradual moisture
front as shown by the dotted line in Figure 10. Such a wetting
front would be characteristic of a more fine-textured soil (especially a lower n value). This tendency toward predicting soil
hydraulic parameters of finer soils is a general trend caused bythis type of averaging of measured pressure heads. A more
appropriate way of averaging measured variables at the samethese measurements do not necessarily provide the same results because of local soil heterogeneity. Figure 10 shows idealized but typical tensiometer responses of such duplicate measurements in the soil profile with each tensiometer showing a
relatively sharp curve. An interesting question is how best to
analyze these curves when they differ significantly. If the response curves were analyzed independently, the sharp wetting
Table 5. Soil Hydraulic Parameters Fitted to the
Combined Data Measured With the Evaporation Method,
the Pressure Chamber, and the Steady State Tension
Infiltration
Parameter

Fit 1

Fit 2

␪r
␪s
␣ d , cm⫺1n
K s , cm d⫺1
lR2

0.0761
0.348
0.0347
1.46
151.4
0.5a
0.967

0
0.337
0.0303
1.31
255.1
2.23
0.968

a

Parameter fixed during optimization.

Figure 9. Hydraulic conductivity functions determined using
different optimization options and data collected for experiments III and IV. Hydraulic conductivities were measured
using the evaporation method and making use of Wooding’s[1968] analysis of tension disc experiments III and IV as well as
steady state infiltration experiments.

2978

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS

Figure 10. Hypothetical pressure head profiles at a nodal
point with coordinates (r, z) and possible averaging at given
values of time or the pressure head.
radial coordinates seems to be to use the mean time for thesame value of the measured variable (dash-dotted line in Fig.
10) rather than the mean for the same time. This alternative
type of averaging would lead to the correct shape parameters( ␣ , n) and to approximately the correct (mean) K s .
Acknowledgment. Thanks are due to N. Romano for his many
useful review suggestions and comments.

References
Ankeny, M. D., M. Ahmed, T. C. Kaspar, and R. Horton, Simple field
method for determining unsaturated hydraulic conductivity, Soil Sci.
Soc. Am. J., 55, 467– 470, 1991.
Bouma, J., D. I. Hillel, F. D. Hole, and C. R. Amerman, Field measurement of unsaturated hydraulic conductivity of unsaturated soil
near saturation, Soil Sci. Soc. Am. Proc., 35, 362–364, 1971.
Cı́slerová, M., J. Šimůnek, and T. Vogel, Changes of steady-state
infiltration rates in recurrent ponding infiltration experiments, J.
Hydrol., 104, 1–16, 1988.
Clothier, B. E., M. B. Kirkham, and J. E. McLean, In situ measurement of the effective transport volume for solute moving throughsoil, Soil Sci. Soc. Am. J., 56, 733–736, 1992.
Cuenca, R. H., et al., Soil measurements during HAPEX-Sahel intensive observation period, J. Hydrol., 188-189, 224 –266, 1997.
Gardner, W. R., Some steady-state solutions of the unsaturated moisture flow equation with application to evaporation from a watertable, Soil Sci., 85, 228 –232, 1958.
Hollenbeck, K. J., and K. H. Jensen, Experimental evidence of randomness and nonuniqueness in unsaturated outflow experiments
designed for hydraulic parameter estimation, Water Resour. Res., 34,
595– 602, 1998.
Inoue, M., J. Šimůnek, J. W. Hopmans, and V. Clausnitzer, In-situ
estimation of soil hydraulic functions using a multi-step soil water
extraction technique, Water Resour. Res., 34, 1035–1050, 1998.
Jarvis, N. J., and I. Messing, Near-saturated hydraulic conductivity insoils of contrasting texture measured by tension infiltrometers, Soil
Sci. Soc. Am. J., 59, 27–34, 1995.
Kachanoski, R. G., I. J. Van Wesenbeeck, P. Von Bertoldi, A. Ward,and C. Hamplen, Measurement of soil water content during threedimensional axial-symmetric water flow, Soil Sci. Soc. Am. J., 54,
645– 649, 1990.
Kodešová, R., M. M. Gribb, and J. Šimůnek, Estimating soil hydraulicproperties from transient cone penetrometer data, Soil Sci., 163,
436 – 453, 1998.

Kool, J. B., and J. C. Parker, Development and evaluation of closedform expressions for hysteretic soil hydraulic properties, Water Resour. Res., 23, 105–114, 1987.
Logsdon, S. D., E. L. McCoy, R. R. Allmaras, and D. R. Linden,
Macropore characterization by indirect methods, Soil Sci., 155, 316 –
324, 1993.
Luckner, L., M. T. van Genuchten, and D. R. Nielsen, A consistent setof parametric models for the two-phase flow of immiscible fluids inthe subsurface, Water Resour. Res., 25, 2187–2193, 1989.
Marquardt, D. W., An algorithm for least-squares estimation of nonlinear parameters, SIAM J. Appl. Math., 11, 431– 441, 1963.
Mualem, Y., A new model for predicting the hydraulic conductivity ofunsaturated porous media, Water Resour. Res., 12, 513–522, 1976.
Perroux, K. M., and I. White, Design for disc permeameters, Soil Sci.
Soc. Am. J., 52, 1205–1215, 1988.
Quadri, M. B., B. E. Clothier, R. Angulo-Jaramillo, M. Vauclin, and
S. R. Green, Axisymmetric transport of water and solute underneatha disk permeameter: Experiments and numerical model, Soil Sci.
Soc. Am. J., 58, 696 –703, 1994.
Reynolds, W. D., and D. E. Elrick, Determination of hydraulic conductivity using a tension infiltrometer, Soil Sci. Soc. Am. J., 55,
633– 639, 1991.
Russo, D., E. Bresler, U. Shani, and J. C. Parker, Analysis of infiltration events in relation to determining soil hydraulic properties by
inverse problem methodology, Water Resour. Res., 27, 1361–1373,
1991.
Scott, P. S., G. J. Farquhar, and N. Kouwen, Hysteresis effects on net
infiltration, in Advances in Infiltration, Publ. 11-83, pp. 163–170, Am.
Soc. of Agric. Eng., St. Joseph, Mich., 1983.
Šimůnek, J., and M. T. van Genuchten, Estimating unsaturated soil
hydraulic properties from tension disc infiltrometer data by numerical inversion, Water Resour. Res., 32, 2683–2696, 1996.
Šimůnek, J., and M. T. van Genuchten, Estimating unsaturated soil
hydraulic properties from multiple tension disc infiltrometer data,
Soil Sci., 162, 383–398, 1997.
Šimůnek, J., M. Šejna, and M. T. van Genuchten, The HYDRUS-2D
software package for simulating water flow and solute transport in
two-dimensional variably saturated media: Version 1.0, IGWMCTPS-53, 167 pp., Int. Ground Water Model. Cent., Colo. Sch. ofMines, Golden, 1996.
Šimůnek, J., R. Angulo-Jaramillo, M. G. Schaap, J.-P. Vandervaere,
and M. T. van Genuchten, Using an inverse method to estimate thehydraulic properties of crusted soils from tension disc infiltrometer
data, Geoderma, 86, 61– 81, 1998a.
Šimůnek, J., O. Wendroth, and M. T. van Genuchten, A parameter
estimation analysis of the evaporation method for determining soil
hydraulic properties, Soil Sci. Soc. Am. J., 62, 894 –905, 1998b.
Šimůnek, J., R. Kodešová, M. M. Gribb, and M. T. van Genuchten,
Estimating hysteresis in the soil water retention function from modified cone penetrometer test, Water Resour. Res., 35, 1329 –1345,
1999.
Smettem, K. R. J., and B. E. Clothier, Measuring unsaturated sorptivity and hydraulic conductivity using multiple disk permeameters, J.
Soil Sci., 40, 563–568, 1989.
van Genuchten, M. T., A closed-form equation for predicting thehydraulic conductivity of unsaturated soils, Soil Sci. Soc. Am. J., 44,892– 898, 1980.
van Genuchten, M. T., F. J. Leij, and S. R. Yates, The RETC code forquantifying the hydraulic functions of unsaturated soils, EPA/600/291-065, Off. of Res. and Dev., U.S. Environ. Prot. Agency, Washington, D. C., 1991.
van Genuchten, M. T., F. J. Leij, and L. Wu (Eds.), Characterizationand Measurement of the Hydraulic Properties of Unsaturated PorousMedia, Univ. of Calif., Riverside, Calif., 1999.
Vogeler, I., B. E. Clothier, S. R. Green, D. R. Scotter, and R. W.
Tillman, Characterizing water and solute movement by time domainreflectometry and disk permeametry, Soil Sci. Soc. Am. J., 60, 5–12,
1996.
Wang, D., S. R. Yates, and F. F. Ernst, Determining soil hydraulicproperties using tension infiltrometers, time domain reflectometry,
and tensiometers, Soil Sci. Soc. Am. J., 62, 318 –325, 1998.
Ward, A. L., R. G. Kachanoski, and D. L. Elrick, Laboratory measurements of solute transport using time domain reflectometry, Soil Sci.
Soc. Am. J., 58, 1031–1039, 1994.
Wendroth, O., and J. Šimůnek, Soil hydraulic properties determinedfrom evaporation and tension infiltration experiments and their use

ŠIMŮNEK ET AL.: LABORATORY TENSION DISC INFILTROMETER EXPERIMENTS
for modeling field moisture status, in Characterization and Measurement of the Hydraulic Properties of Unsaturated Porous Media, editedby M. T. van Genuchten, F. J. Leij, and L. Wu, pp. 737–748, Univ.
of Calif., Riverside, Calif., 1999.
Wendroth, O., W. Ehlers, J. W. Hopmans, H. Kage, J. Halbertsma,
and J. H. M. Wösten, Reevaluation of the evaporation method for
determining hydraulic functions in unsaturated soils, Soil Sci. Soc.
Am. J., 57, 1436 –1443, 1993.
White, I., and M. J. Sully, Macroscopic and microscopic capillary
length and time scales from field infiltration, Water Resour. Res., 23,1514 –1522, 1987.
Wind, G. P., Capillary conductivity data estimated by a simple method,in Water in the Unsaturated Zone, edited by P. E. Rijtema and H.

2979

Wassink, vol. 1, pp. 181–191, Int. Assoc. of Sci. Hydrol., Gentbrugge,
Germany, 1968.
Wooding, R. A., Steady infiltration from large shallow circular pond,
Water Resour. Res., 4, 1259 –1273, 1968.
J. Šimůnek and M. T. van Genuchten, U. S. Salinity Laboratory,
USDA-ARS, 450 West Big Springs Road, Riverside, CA 92507.
(JSimunek@ussl.ars.usda.gov)
O. Wendroth, Institute for Soil Landscape Research, ZALF,
Eberswalder Straße 84, D-15374 Müncheberg, Germany.
(Received October 9, 1998; revised May 10, 1999;
accepted June 1, 1999.)2980