Special Section:
Preferential Flow

Radka Kodešová*Jií Šimnek
Antonín Nikodem
Veronika Jirk

Estimation of the DualPermeability Model Parameters
using Tension Disk Inltrometer
and Guelph Permeameter
A determinaon of the parameters describing the soil hydraulic properes of matrix andmacropore domains and mass exchange between these domains is crucial when preferenal water ow in structured soils is simulated using the dual-permeability model. This study
focused on esmang the parameters of the radially symmetric dual-permeability modelfrom cumulave inltraon measured in the surface horizon of a Haplic Luvisol. Whileparameters obtained from the numerical inversion of the tension disk inltraon, using thesingle-porosity ow model in HYDRUS 2D/3D, were used to describe the matrix domain, theparameters characterizing the macropore domain and mass exchange between domains
were esmated using the Guelph permeameter inltraon and the dual-permeability owmodel in HYDRUS 2D/3D. The mass transfer coecient between the two pore domainsaected the simulated water regime considerably, and subsequently, the calibrated valueof the saturated hydraulic conducvity in the macropore domain, Ksf. A less signicant
impact of the aggregate shape factor was observed due to a low range of possible values
compared with the other two parameters, which either varied within orders of magnitude(the eecve saturated hydraulic conducvity of the interface between the two poredomains, Ksa) or were squared (the characterisc length of an aggregate, a). The Ksf values
increased when mass exchange decreased (when a increased and Ksa decreased). Sinceboth parameters are mutually correlated and therefore have a similar impact on simulateddata, we suggest that a be determined independently, and Ksa and Ksf should be simultaneously opmized when the parameters of the dual-permeability model are evaluated using
the presented experimental procedure.
Abbreviaons: MSO, mulstep oulow; PI, ponding inltraon.

Hydraulic and mass exchange term
parameters for the HYDRUS 2D/3D
dual-permeabilit y model were
estimated using data from micromorphological images and dye tracer
experiments, and numerical inversionof transient infiltration data measured with tension disk inltrometer
and Guelph permeameter methods inthe surface horizon of Haplic Luvisol.
R. Kodešová, A. Nikodem, V. Jirk, Czech Univ.
of Life Sciences Prague, Faculty of Agrobiology, Food and Natural Resources, Dep. of
Soil Science and Soil Protecon, Kamýcká
129, 16521 Prague 6, Czech Republic; and J.
Šimnek, Dep. of Environmental Sciences,
900 University Av., Univ. of California, Riverside, CA 92521. *Corresponding author(kodesova@af.czu.cz).
Vadose Zone J. 9:213–225
doi:10.2136/vzj2009.0069
Received 3 June 2009.
Published online 3 May 2010.
© Soil Science Society of America
5585 Guilford Rd., Madison, WI 53711 USA.
All rights reserved. No part of this periodical may
be reproduced or transmi ed in any form or by anymeans, electronic or mechanical, including photocopying, recording, or any informaon storage andretrieval system, without permission in wring fromthe publisher.

2010, Vol. 9

Soil and groundwater contaminaon by various chemicals used in agriculture and industry present a serious environmental problem that has been widely studied.
Review studies by Šimůnek et al. (2003), Gerke (2006), Jarvis (2007), Clothier et al. (2008),
Šimůnek and van Genuchten (2008), and Köhne et al. (2009a,b) have summarized anddocumented that water flow and contaminant transport in soils is frequently influenced
by either water and solute temporal immobilization or preferential flow. Many numerical models have been developed to describe such nonequilibrium water flow and solutetransport in soils. Commonly used physically based models simulating water flow andsolute transport in soils use either single continuum or bi- and multicontinuum approaches.
Bi- and multimodal concepts assume that the soil porous system is divided into two (ormore) domains, while each domain is characterized by its own set of transport properties
and equations for describing local flow and transport processes. Only two-domain models
have been widely tested, however.
The dual-porosity (mobile–immobile) approach (Philip, 1968; Šimůnek et al., 2003)
describes water flow and solute transport in systems consisting of domains with both
mobile and immobile water. The dual-porosity formulation is based on a set of equations
describing water flow and solute transport in the mobile domain (Richards and advection–
dispersion equations), with mass balance equations describing soil water and solute contents
in the immobile domain. The dual-permeability approach assumes that water flow andsolute transport occur in both domains. The dual-permeability formulation is based on a setof equations that describe water flow and solute transport separately in each domain (thematrix and macropore domains). While the Richards equation is typically used to describewater flow in the matrix domain, various equations are used to describe water flow in themacropore domain (Šimůnek et al., 2003; Gerke, 2006). The gravity-driven kinematic

www.VadoseZoneJournal.org | 213

wave approach was applied by Germann (1985), Germann andBeven (1985), and Jarvis (1994). The Richards equation was usedby Othmer et al. (1991) and Gerke and van Genuchten (1993,
1996). There are also other approaches describing preferential flowin discrete macropores, which are based on Poiseuille’s equation(Ahuja and Hebson, 1992), boundary layer theory for viscous flow
(Germann, 1990), and the Chezy–Manning equation for turbulent
flow in macropores (Chen and Wagenet, 1992a). The transport in
the macropore domain can be described using either the advection–dispersion equation (Gerke and van Genuchten, 1993) or apiston-type transport (Jarvis, 1994).
There are also many approaches for describing water and solutetransfer between the mobile and immobile domains and between
the macropore and matrix domains (Šimůnek et al., 2003; Gerke,2006). Pressure-head-based (Gerke and van Genuchten, 1993;
Zimmerman et al., 1996) or saturation-based (Philip, 1968;Šimůnek et al., 2003; Jarvis, 1994) mass transfer terms may be
used to describe water transfer. The Green–Ampt infiltration equation (Ahuja and Hebson, 1992), the Philip infi ltration equation(Chen and Wagenet, 1992a), or the Richards equation (Chen andWagenet, 1992b) are applied to simulate water infi ltration fromthe discreet macropores into the matrix. Depending on the masstransfer model used to describe water flow, diff usion or advection
mass transfer terms are used to describe solute transfer.
The HYDRUS-1D and HYDRUS 2D/3D codes (Šimůnek et al.,2008), which were used in this study, include single-porosity, dualporosity, and dual-permeability models all based on the numericalsolution of the Richards and advection–dispersion equations fordifferent flow domains. Pressure-head-based or saturation-basedmass transfer terms are used to characterize water transfer between
the mobile and immobile domains. Pressure-based terms areapplied for simulating water transfer between the macropore andmatrix domains of the dual-permeability model.
To apply two- (or multi-) domain models, the domain fractions
and the properties characterizing water flow and solute transport
within and between the domains must be specified. Numericalinversion of the observed water flow and tracer transport data is
usually used to obtain the desired information (see review articles
referenced above); however, many parameters must be measured
independently or determined from the literature to obtain reliable results. A summary of the various approaches for independent
parameter determination was given by Köhne et al. (2009a). The
hydraulic parameters of multimodal soil porous systems are usually determined based on the transient flow and transport dataobtained under various laboratory conditions (Pot et al., 2005;
Köhne et al., 2006a,b; Kodešová et al., 2006a, 2008, 2009).
Tension disk infiltrometers and Guelph permeameters are frequently
used to measure the unsaturated, K(h), and saturated, Ks, hydraulicconductivities of the three-dimensional soil–pore system in the field.

The analytical expressions presented by Wooding (1968) and Zhang
et al. (1998) are usually used to calculate K values from the tension
disk infiltrometer and the Guelph permeameter, respectively. An
inversion of models based on the numerical solution of the Richards
equation and the van Genuchten–Mualem expressions describing
soil hydraulic properties has also been applied to evaluate the othersoil hydraulic parameters from the tension disk infiltrometer transient flow data (Šimůnek and van Genuchten, 1996, 1997; Schwartz
and Evett, 2002; Ramos et al., 2006).
On the other hand, data obtained using the Guelph permeameter
have not yet been numerically analyzed to estimate soil hydraulicparameters. Numerical simulations have been used only to define
the shape factor, which was utilized in the analytical evaluationof the saturated hydraulic conductivity (Hayashi and Quinton,2004); however, parameter estimation was applied to a new tool
called a Restricted Interval Guelph permeameter (Freifeld andOldenburg, 2000). Numerical inversion was also utilized toevaluate the soil hydraulic parameters from a cone permeameter
experiment (Kodešová et al., 1998, 1999; Šimůnek et al., 1999),
which had a geometry of the flow domain similar to the Guelph
permeameter experiment. None of these studies assumed nonequilibrium water flow.
The main objective of this study was to take advantage of different
boundary flow conditions applied during tension disk infi ltrometer (performed under unsaturated conditions) and Guelph
permeameter (performed under ponding infi ltration conditions)experiments to assess preferential flow in macropores, i.e., poreswith diameters larger than the diameter corresponding to theapplied pressure head of the tension disk infi ltrometer. The single-porosity and dual-porosity models in HYDRUS 2/3D wereapplied to estimate the parameters of the soil hydraulic functionsof both models using the tension disk infiltrometer and the Guelph
permeameter transient flow data. The goals of this study were: (i)assessment of the impact of various retention curve parameters(obtained independently using two different laboratory methods)
on optimized Ks using the single-porosity model for both field
experiments; (ii) evaluation of the possible influence of macropores on water flow and, consequently, on the estimated valuesof the Ks by comparing the optimized values obtained using thesingle-porosity model and different field experiments (disk orpermeameter); (iii) Estimation of the Ksf value for the macropore
domain of the dual-permeability model using the Guelph permeameter transient flow data. In the case of this third goal, the soil
hydraulic parameters obtained using the tension disk infiltrometerdata were assumed to characterize the matrix domain. The parameters characterizing dual-domain geometry, water retention in themacropore domain, and the hydraulic conductivity of the interface
between the macropore and matrix domains were also unknown
and had to be estimated. Therefore, additional goals of this study
were: (iv) evaluation of the impact of various parameters describing aggregate structure on the estimated Ksf values; (v) evaluation

www.VadoseZoneJournal.org | 214

of the influence of various soil water retention parameters in themacropore domain on the estimated Ksf values; and (vi) evaluationof the impact of various effective saturated hydraulic conductivitiesof the interface between the macropore and matrix domains on
the estimated Ksf values.

6 Materials and Methods
Field Experiments
The study was performed at the experimental station of the Crop
Research Institute in Hněvčeves, the Czech Republic. The soil
studied was Haplic Luvisol (loess parent material). Four soil
diagnostic horizons were identified in the soil profi le: Ap, 0 to
35 cm; Bt1, 35 to 57 cm; Bt2, 57 to 93 cm; and Ck, 93 to 129 cm.
Conventional tillage was used at this location. Crops were planted
within small areas that allowed the soil properties to be studied
every year under similar conditions. Despite a 5-yr rotation system,
the properties were studied using the same crop (winter barley[Hordeum vulgare L.] or wheat [Tritcum aestivum L.]) and with areasonable distance between the sampling and experimental plots.
Before the study began, the hydraulic properties of the diagnostic
horizons were studied in the laboratory using multistep outflow
experiments performed on undisturbed 100-cm3 soil samples(soil core height of 5.1 cm and cross-sectional area of 19.60 cm2)
(Kodešová et al., 2008). Initially, fully saturated soil samples placedin the Tempe cells were slowly drained using nine pressure headsteps (a minimum pressure head of −1000 cm) during a 3-wk
period. Ponding infi ltration tests on undisturbed 1125-cm3 soilcolumns (soil core height of 13 cm and cross-sectional area of 86.54
cm2) were then performed to study the hydraulic properties of different diagnostic horizons (Kodešová et al., 2009). The ponding
depth was 1.5 cm, and the experiments lasted 150 min. In bothstudies, micromorphological images were used to evaluate themacroporosity, which was defined as pores with diameters >40 μm.
The single- and dual-porosity models in HYDRUS-1D (Šimůnek
et al., 2008) were successfully applied to simulate water flow in the
Ap1, Ap2, and Bt1 horizons for the multistep outflow experiments(experiments performed mostly under unsaturated conditions) andto optimize the soil hydraulic parameters in all mobile domains.
The single-porosity and dual-permeability models in HYDRUS-1D
were used to simulate water flow in the Ap, Bt1, and Bt2 horizons
under ponding infiltration conditions (experiments involving both
unsaturated and saturated conditions) and to estimate the soil
hydraulic parameters and solute transport parameters in all domains.
While the dual-permeability model, in combination with information obtained from the micromorphological images, successfully
described water flow and solute transport in the subsurface horizons,
it was less successful for describing flow in the Ap horizon. An analysis of the soil pores showed bimodality in the detected macropore
system, with a high fraction of gravitational pores. The influenceof the gravitational pores on the observed water flow and solutetransport could not be simulated successfully using the assumptionof one macropore domain representing both large capillary and gravitational pores. To remedy this, new experiments were performed
directly in the field to assess this phenomenon.
Field experiments were performed immediately after the harvest of
winter barley in 2008 (at approximately the time of the previous
years’ sampling), at a close distance from the previous years’ sampling plots. Measurements were performed randomly within anarea of 10 m2 , avoiding the locations impacted by either the wheelsof heavy machinery or by the experimental work performed in thearea. A tension disk infi ltrometer (with a disk radius of 10 cm) wasused to measure the cumulative water infi ltration under unsaturated conditions created using a pressure head of −2 cm (Watsonand Luxmoore, 1986). Measurements were performed at a depthof 5 cm. A Guelph permeameter was used to measure the cumulative water flux under surface ponding conditions. The depth of thedrilled well was 10 cm, the well radius was 3 cm, and the well ponding depth was 5 cm. The depths of the experiments were selectedbased on the following observations: (i) it was observed in the fieldthat after excluding the top 1 cm of the soil crust, the humic Ap
horizon was relatively homogeneous to a depth of 25 cm; (ii) waterflow under the tension disk infi ltrometer was mostly downward
but also horizontal and therefore the wetted domain was expected
to be between depths of 5 and approximately 15 cm; and (iii) using
the Guelph permeameter, water flow is downward, horizontal, andalso upward and therefore the wetted domain was expected to be
between depths of 1 cm (i.e., above a depth of 5 cm, which was the
position of a water table in the well) and 20 cm. The maximum
wetted horizontal cross-section was expected to be similar to thearea of the disk infi ltrometer.
The soil surface was carefully leveled using a sharp knife before
the tension disk infi ltrometer tests. Then, a 1-mm layer from thesame soil, sieved through a 2-mm sieve, was formed to ensure close
contact between the soil and the disk. The standard procedure
recommended in the Guelph permeameter manual (Soilmoisture
Equipment Corp., 2008) was used to prepare the infi ltration well.
A soil auger was used to drill the well to the desired depth. A sizing
auger was used for the final well shaping, and the well brush wasapplied to remove the smear layer on the well sides. The tension
disk infi ltrometer tests lasted at least 60 min and the Guelph
permeameter tests at least 30 min. Data were collected manually. Cumulative infi ltration was measured using a Mariotte tube
(shaded from sunlight if necessary), and times were measured witha stopwatch. The results of three replicates with each permeameterare presented here.
A dye infi ltration experiment was performed in 2009, using a procedure similar to that described by Sander and Gerke (2007). One
hundred liters of solution with the food color Brilliant Blue FCF(5 kg m−3) was infi ltrated into a 1- by 1-m plot (applying an initial
ponding depth of 10 cm) immediately after wheat harvest. On the

www.VadoseZoneJournal.org | 215

next day, one half of the plot was sliced horizontally and the otherhalf vertically to study the dye distribution in the soil profi le toa depth of 100 cm. An example of a staining pattern within thevertical and horizontal (to a depth of 20 cm) sections is shown in
Fig. 1 and was used to estimate the dual-domain geometry. Figure
1 shows that a less permeable crust layer 1 cm thick reduced waterand dye infi ltration at the top of the surface Ap1 horizon (0–25cm). Below the crust layer, the dye tracer was transported mostly
regularly, and only some isolated uncolored domains were visible
in this horizon. The dye transport in the Ap2 horizon (25–35 cm,
the plow pan) slowed down due to the very compact matrix structure. The dye transport through the plow pan was caused mainly
by gravitational biopores.

Water Flow Models
The single-porosity and dual-permeability models implemented
in HYDRUS 2D/3D (Šimůnek et al., 2008) were used to simulate the observed cumulative water fluxes and to estimate thesoil hydraulic parameters via numerical inversion. The Richards
equation, which describes isothermal Darcian flow in a variably
saturated rigid porous medium, is used in both models. The flow
domain was assumed to be isotropic, as previous micromorphological studies (Kodešová et al., 2008, 2009) showed that the soilmaterial consisted of regular angular blocky soil aggregates withoutlayering or other features that could cause anisotropic conditions.
The Richards equation for the radially symmetric, single-porosity
system may be applied to simulate both transient experiments: 1  h   h= rK ( h )  + K (h )+ K ( h )
 t r  r 
 r   z z[1]

where  is the volumetric soil water content [L3 L−3], h is the pressure head [L], K is the hydraulic conductivity [L T−1], t is time [T],r is the radial coordinate [L], and z is the vertical axis [L]. Equation[1] is solved for the entire flow domain using one set of soil water
retention and hydraulic conductivity functions.
In the dual-permeability model, the Richards equation is applied
separately to each of the two pore regions, i.e., the macropore (fractures, a domain of larger pores) and matrix domains (Gerke andvan Genuchten, 1993). The Richards equation for the radially symmetric, dual-permeability system may be written as
 f
h 
1  
 rK f ( hf ) f =r  r 
t
r 

h
 
w
 K f ( hf ) f + K f ( hf )
+ fw
 z 
z m 1  h  rK m ( hm ) m =r r 
t
r 

h
 
 K m ( hm ) m + K m ( hm ) +
+

 1
z 
z

[2]

w

fw

where the subscripts f and m refer to the macropore and matrix
domains, respectively, w is the mass transfer term for waterexchange between the macropore and matrix domains [T−1], andfw is the macropore domain fraction (dimensionless), i.e., the ratioof the volumes of the macropore and total flow domains. The matrix
domain fraction is specified as 1 − fw. Equations [2] are solved using
two sets of soil water retention and hydraulic conductivity functions,which are defined for each domain. The total soil water content andthe hydraulic conductivity are defined as the sum of the soil watercontents or hydraulic conductivities of each domain multiplied bycorresponding domain fractions.
The mass exchange between the matrix and macropore regions,is calculated as
w =

a2

wKa

( ha )( hf

and

w,

hm )
[3]

K a ( ha ) = K sa K rm ( ha )Fig. 1. Distribution of Brilliant Blue FCF dye within the Ap1 and Ap2 horizons after
infiltration in a 1- by 1-m plot.

where Ka is the effective hydraulic conductivity of theinterface between the two pore domains [L T−1], Ksa
is the saturated hydraulic conductivity of the interface,
K rm(h a) is the matrix relative unsaturated hydraulic conductivity function calculated using the van
Genuchten model (see below), and ha is the pressure
head at the interface [taken as (hf + hm)/2]. The parameters describing the aggregate shapes (Gerke and van
Genuchten, 1996) are the dimensionless shape factor(15 for spherical aggregates, 3 for cubic aggregates),
the characteristic length of an aggregate, a [L] (a sphereradius or half size of the cube edge), and the dimensionless scaling factor w (= 0.4).

www.VadoseZoneJournal.org | 216

The analytical expressions proposed by van Genuchten (1980) for
the soil water retention curve, (h), and the hydraulic conductivity
function, K(), are used in both models:e = (h ) r
=
s r

1

(

1 + |  h| n

m

)

e =1

K (  ) = K s  le 1
K () = K s

h<0
[4]
h 0

m 2

(1 1/e m ) 

h<0

[5]

h 0

where  e is the effective soil water content (dimensionless), Ks is thesaturated hydraulic conductivity [L T−1],  r and  s are the residual
and saturated soil water contents [L3 L−3], respectively, l is thepore-connectivity parameter (dimensionless),  is the reciprocalof the air-entry pressure head [L−1], n (dimensionless) is related
to the slope of the retention curve at the inflection point, and m= 1 − 1/n.

Soil Hydraulic Parameter EsmaonFirst, the single-porosity model was applied to estimate the soil
hydraulic parameters of the Ap1 horizon from the cumulative
water infi ltration measured using both the tension disk infi ltrometer and the Guelph permeameter. The finite element mesh for an
axisymmetrical quasi-three-dimensional flow domain is shown in
Fig. 2 for the Guelph permeameter. The flow domain either contained an axisymmetrical well of a radius of 3 cm and a depth of
10 cm for the Guelph permeameter or was flat at the top for thetension disk infi ltrometer (not shown). The right side and bottom
boundaries of both flow domains were set far enough away so as

not to influence the simulated transient flow data. An entire simulated domain was assumed to consist of only a single material. The
initial conditions were specified using pressure heads that were
calculated for the average measured initial soil water content (0.33
cm3 cm−3) using the applied soil water retention curve (which is
specified below). A constant pressure head of −2 cm was set at theboundary where the tension disk infi ltrometer came into contact
with the soil. A constant hydrostatic pressure head (5 cm at thebottom) was set at the boundary where ponding infi ltration wasapplied for the Guelph permeameter. A free-drainage boundary
condition was specified at the bottom of the transport domain.
Because it was assumed that boundary conditions at the remaining
domain boundaries do not affect the simulated flow, a zero-flux
boundary condition was set at these boundaries (see, for instance,Šimůnek and van Genuchten, 1996, 1997; Šimůnek et al., 1999;
Kodešová et al., 1998, 1999).
Because it has been shown that cumulative infi ltration for the onestep disk infi ltrometer experiment without additional transient
flow data does not provide enough information for the numerical
optimization of all soil hydraulic parameters (Šimůnek and van
Genuchten, 1996), the following simplifications were used. The
parameters of the soil water retention curve ( r,  s, , and n) were
set to values obtained in previous studies (Table 1) performed on
the same soil. In these studies, multistep outflow (MSO) experiments (Kodešová et al., 2008) and ponding infi ltration (PI) tests(Kodešová et al., 2009) were analyzed using the HYDRUS-1D
single-porosity model. It should be mentioned that while theparameters  r, , n, and Ks were optimized via numerical inversionof both the MSO and PI transient flow data, the saturated watercontents  s were set to measured values. The  s value for the MSO
data analysis was obtained from a saturation pan. The  s value for
the PI data analysis was equal to the sample porosity. A lower  svalue for the MSO test compared with that obtained for the PI test
resulted because samples are usually taken in such a way as to avoid
large (gravitational) pores or cracks, and the probability that large
pores or cracks are present in the soil core increases with its size.
Thus, a smaller fraction of macropores could be included in smaller
samples, resulting in a lower value of the sample porosity comparedwith that in larger soil samples (100- vs. 1125-cm3 soil samples). In
addition, the  s value obtained from the saturation pan might havebeen slightly lower than the actual sample porosity. The values ofthe  and n parameters ( MSO <  PI and nMSO < nPI) indicate
Table 1. Soil hydraulic parameters (van Genuchten, 1980) obtained inprevious studies.
Experiment

r

s



n

—— cm3 cm−3 —— cm−1

Fig. 2. A finite element mesh generated for the Guelph permeameter experiments. Coordinates are in centimeters. Blue nodes at theboundary represent a constant-head boundary condition. Coordinates of the A, B, C, D, E, and F points in the radially symmetricsystem are in centimeters.

Ks
cm min−1

Ponding infi ltration† 0

0.428

0.01242.001

0.0377

Multistep outflow‡

0.363

0.0072

1.758

3.12 × 10−4

0.217

† Kodešová et al. (2009).
‡ Kodešová et al. (2008).

www.VadoseZoneJournal.org | 217

that the shape of the soil water retention curve obtained from the
MSO test is more gradual than that obtained from the PI test. Notethat  from the wetting experiment (PI) is larger than from thedrying experiment (MSO), while the n values are relatively similar,
as would be expected considering the hysteresis model of Kool andParker (1987). Parameters resulting from different experiments
were chosen to assess the impact of various phenomena, such as
different experimental history (drying vs. wetting experiments)and flow conditions (mostly unsaturated vs. mostly saturated), on
optimized K s values. The value of l was set equal to 0.5 (Mualem,1976). Saturated hydraulic conductivities, Ks, were optimized via
numerical inversion using the single-porosity model in HYDRUS
2D/3D. The same procedure was applied to analyze the cumulative
infi ltration measured using the Guelph permeameter.

deeper progress of the wetting front in the macropore domaincompared with that in the matrix domain). These parameters (atleast in this case) influenced cumulative infi ltration only marginally, however. Kodešová et al. (2009) used values of = 5 and a =0.2 cm to characterize small blocky aggregates. The ponding dye
infi ltration experiment performed in 2009 (Fig. 1) showed thatwater and the dye tracer moved mostly uniformly through thesoil; however, some pore domains remained isolated. Therefore,
several scenarios were simulated (see Table 2). The shape factor,, was either 5 or 15 to allow study of the impact of the aggregate
shape (blocky or spherical). The characteristic length of aggregates,a, was either 0.2 or 1.5 cm to simulate water flow in systems witheither small (for a region with regular dye distribution and for avalue assessed using the micromorphological images) or larger (forregions less impacted by dye tracer) aggregates.

Next, the dual-permeability model in HYDRUS 2D/3D wasused to analyze the cumulative infi ltration measured using the
Guelph permeameter. The same finite element mesh as was used
for the single-porosity model was utilized, and the same boundary
conditions were applied for both domains. The initial conditions
were set similar to those set for the single-porosity model. To haveinitially no mass transfer between the two domains, the same pressure head was specified in both domains. The initial pressure heads
were calculated using the applied soil water retention curves (which
are specified below for each scenario; Tables 1 and 2), assumingthat the sum of the initial soil water contents in both domains
multiplied by their domain fractions equals the average measured
initial soil water content (0.33 cm3 cm−3). The fraction of themacropore domain ( fw = 0.0277) was determined, based on the
micromorphological study of Kodešová et al. (2009), as the ratioof the image porosity for pores >1000 μm (0.0119) and the sample
porosity (0.428).

The retention curve parameters ( r,  s,  and n) optimized forthe laboratory ponding infi ltration (PI inputs in Table 1) and Ks
optimized for the tension disk infi ltrometer experiment (Table 2,tension disk infi ltrometer Replicate 1) were used to characterize
the parameters ( r,m,  s,m,  m, nm, and Ks,m) of the matrix domain.
The data from the laboratory ponding infi ltration were chosen due
to their similar experimental history (a wetting experiment). In
addition, Gribb et al. (2004) showed that for sandy soils the shapesof hydraulic properties obtained using the cone permeameter (alsoan infi ltration experiment and a flow domain geometry similar
to that of the Guelph permeameter) were less gradual than thosemeasured using an MSO experiment.
The parameters of the soil water retention curve for the macropore
domain were set to ensure that macropores were fi lled with water
only for pressure heads close to zero (a step-like shape curve). This
assumption was clearly demonstrated in the analysis of bimodal
soil water retention curves by Othmer et al. (1991). Previous studies also showed that the saturated water content of the macropore
domain,  s,f, may be the same as the  s,m of the matrix domain(Gerke and van Genuchten, 1996; Šimůnek et al., 2003; Kodešová
et al., 2008, 2009) or larger (Gerke and Köhne, 2004; Gerke et al.,2007; Kodešová et al., 2006b, 2008), depending on the available
data and parameter analysis. The residual soil water content of themacropore domain,  r,f, was commonly assumed to be close to zeroin these studies. The parameters  f and nf varied from 0.05 to 0.15
cm−1 and from 2 to 4, respectively. On the other hand, Vogel et al.

Parameters characterizing the aggregate geometry were set basedon the study of Kodešová et al. (2006b), the micromorphologicalstudies of Kodešová et al. (2009), and the ponding dye infi ltration
experiment performed in 2009. Kodešová et al. (2006b) discussed
the impact of the variable shape factor ( = 3 and 15) and thecharacteristic length of the aggregates (a = 1 and 2.5 cm) on the
one-dimensional numerical simulation of 1-h ponding infi ltration
into clay soils. They showed that a lower value and a greater avalue increased the final wetting front depth. In addition, a greatera value increased the irregularity of the wetting front shape (aTable 2. Saturated hydraulic conductivities obtained from tension disk infi ltrometer and Guelph permeameter tests.

Infi ltration test

Soil-water retention curveparameters obtained from

Replicate 1
Ks

Replicate 2
R2

cm min−1
Tension disk infi ltrometer
Guelph permeameter

Ks

Replicate 3
R2

cm min−1

Ks

R2

cm min−1

Avg. KsSD

cm min−1

ponding infi ltration

1.25 × 10−4

0.995

1.02 × 10−4

0.999

0.92 × 10−4

0.990

1.06 × 10−4

0.57 × 10−9

multistep outflow

1.82 × 10−4

0.998

1.53 × 10−4

0.998

1.35 × 10−4

0.994

1.57 × 10−4

1.12 × 10−9

ponding infi ltration

1.61 × 10−3

0.999

1.07 × 10−3

0.994

0.16 × 10−3

0.980

9.47 × 10−4

1.07 × 10−6

multistep outflow

1.66 × 10−3

0.998

1.15 × 10−3

0.990

0.21 × 10−3

0.974

10.07 × 10−3

1.08 × 10−6

www.VadoseZoneJournal.org | 218

(2000), who did not assumed a step-like shaped soil water retention
curve for the macropore domain, used the same soil water retention parameters in both domains. A one-dimensional numericalsimulation of ponding infi ltration (Kodešová et al., 2006b) indicated that  s,f influenced the wetting front depth but did not havea significant impact on cumulative infi ltration. On the other hand,an increased  f value (= 0.07 or 0.15) slightly decreased cumulative
infiltration. As was expected, fw and Ks,f had the greatest impact on
cumulative fluxes. Increasing values of these parameters resulted
in increased cumulative infi ltration and the wetting front moving
to deeper depths. Because the fraction of the macropore domain
was defi ned here as a ratio between the image porosity and thetotal porosity,  s,f was equal to  s,m (both equal to the measured
porosity). Th is means that the total saturated soil water content
was equal to the sum of the saturated soil water content of thematrix domain multiplied by its domain fraction (0.972), and thesaturated soil water content of the macropore domain multipliedby its domain fraction (0.028). The residual water content of themacropore domain,  r,f, was set to zero.
Since no independent measurements or predictions of  f and nf
parameters were available, several scenarios were simulated toassess the impact of retention curve shape parameters (see Table2). Parameter  f was set to 0.05, 0.1, and 0.15 cm−1. Parameternf was set to 2, 3, and 4. Simulations were performed only for= 15 (spherical aggregates, which are usually expected in the
Ap horizon) and a = 0.2 cm. Four scenarios (for all the scenarios
described above) with various values of the effective saturated
hydraulic conductivity, Ksa (Table 3), were studied. The first Ksa
value was set to the value obtained by numerical inversion using
the dual-permeability model in HYDRUS-1D (Kodešová et al.,2009). The other three values of Ksa were sequentially decreased
by one order of magnitude.
The value of Ks,f was calibrated manually because HYDRUS 2D/3D
does not provide inverse capabilities for the dual-permeability

model. The objective function , to be minimized during theparameter estimation process, was formulated using the cumulative infi ltration, q (cm3):N ( L, q ) =  vi  q* (t i ) q (t i , L )

2

[6]

i =1

where κ is the vector of estimated parameters, N is the numberof measurements of cumulative infiltration, q*(ti) are specificmeasurements at time ti, qi(ti,κ) are the corresponding model predictions for parameter vector κ, and vi is the weight associated witha particular measurement point (vi = 1 for all data).
In addition, the possible impact of the less permeable plow pan (theAp2 horizon) between depths of 25 and 35 cm was investigated.
The same geometry structure parameters and soil hydraulic parameters as for the Ap1 horizon, except for the Ks,f value, were assumed
for this layer. The Ks,f value was set lower (0.01 cm min−1) compared with that in the Ap1 horizon. No impact on the estimatedKs,f of the Ap1 horizon, simulated cumulative infi ltration, or the
shape of the wetting front was found. Therefore, these simulations
are not discussed further.

6 Results and Discussion
Single-Porosity Model
The optimized Ks values, obtained using the single-porosity modelby numerically inverting the data for three replications of both thetension disk infi ltrometer and the Guelph permeameter experiments, are shown in Table 2. Different applied soil-water retention
curve parameters ( r,  s, , and n), which resulted from the previous MSO and PI tests (Table 1), noticeably affected the Ks values
obtained from both experiments. In both cases, greater Ks values
were obtained when the MSO parameters were used than when
PI parameters were applied. The first reason for this is the lower
MSO  s value and the subsequently larger initial relative saturation

Table 3. Saturated hydraulic conductivities of the macropore domain, Ks,f, correlation coefficients R 2 , and objective functions  obtained for Test 1
with the Guelph permeameter for various scenarios with different effective saturated hydraulic conductivities of the interface between the two poredomains, Ksa , including aggregate shape factors, , characteristic lengths of an aggregate, a, and the soil water retention parameters  and n of themacropore domain.
Ksa = 5.3 × 10−5 cm min−1
Scenario

a

f

cm

cm−1

nf

R2

Ks,f



cm min−1

Ksa = 5.3 × 10−6 cm min−1
Ks,f

R2



cm min−1

Ksa = 5.3 × 10−7 cm min−1
Ks,f

R2



cm min−1

Ksa = 5.3 × 10−8 cm min−1
Ks,f

R2



cm min−1

A

5

0.2

0.1

3

0.158

0.999

691.7

0.171

0.999

101.9

0.218

0.996

623.6

0.240

0.997

424.8

B

5

1.5

0.1

3

0.160

0.999

627.8

0.239

0.997

420.7

0.245

0.997

374.1

0.244

0.997

373.9

C

15

0.2

0.1

3

0.156

0.999

845.8

0.165

0.999

293.3

0.195

0.997

614.4

0.235

0.997

471.3

D

15

1.5

0.1

3

0.182

0.998

366.8

0.232

0.998

572.4

0.245

0.997

567.0

0.240

0.997

548.9

E

15

0.2

0.05

3

0.100

0.999

733.0

0.105

0.999

218.0

0.125

0.997

602.5

0.150

0.997

437.7

F

15

0.2

0.15

3

0.197

0.999

843.6

0.205

0.999

350.4

0.245

0.996

695.7

0.291

0.996

552.2

G

15

0.2

0.1

2

0.200

0.999

857.1

0.208

0.999

315.6

0.249

0.996

695.7

0.299

0.996

564.4

H

15

0.2

0.1

4

0.140

0.999

893.7

0.148

0.999

304.1

0.175

0.997

592.0

0.215

0.997

495.9

www.VadoseZoneJournal.org | 219

(given by the effective water content always calculated using thesame initial soil water content and either MSO or PI values of r and  s) for the MSO parameters than the PI parameters. The
second reason is the more gradual shape of the retention curvewhen the MSO parameters were used, which resulted in faster flow
and higher Ks values than when the PI parameters were used tosimulate the same cumulative infi ltration. Calculated differences
between the Ks values obtained from the tension disk infi ltrometer using either MSO or PI parameters (5.7 × 10−5 [Replicate 1],
5.1 × 10−5 [Replicate 2], and 4.3 × 10−5 cm min−1 [Replicate 3])
were similar to those calculated for the Guelph permeameter (5.0
× 10−5 [Replicate 1], 8.0 × 10−5 [Replicate 2], and 5.0 × 10−5 cm
min−1 [Replicate 3]). While the K s values obtained from the tension disk infiltrometer using MSO parameters were 1.45 (Replicate1), 1.50 (Replicate 2), and 1.46 (Replicate 3) times higher thanthose using the PI parameters, however, the Ks values obtained
from the Guelph permeameter using MSO parameters were only
1.03 (Replicate 1), 1.07 (Replicate 2), and 1.31 (Replicate 3) times
higher than those using the PI parameters. The impact of different
soil hydraulic properties on the predicted soil water content distributions for a single Guelph permeameter experiment (Replicate 1)
is illustrated in Fig. 3. A more gradual shape of the soil water retention curve, with a lower value of the saturated soil water content
obtained from the multistep outflow experiment (compared to thelaboratory ponding infi ltration), resulted in a deeper wetting front
and more gradual soil water distributions.
The Ks values obtained from the Guelph permeameter were oneorder of magnitude higher than those obtained from the tension
disk infi ltrometer (Table 2). In addition, the analysis of three
replicates for each scenario shows that the variability of the Ks
values obtained using the Guelph permeameter is greater than thatobtained using the tension disk infiltrometer. The greater Ks values
obtained using the Guelph permeameter and their larger variability

may be explained by gravitational macropores. It has been shownin many studies (e.g., Jarvis, 2007), and it is demonstrated by thedye distribution in Fig. 1, that the presence of gravitational macropores considerably increases water and solute fluxes in soils undersaturated conditions. Their occurrence depends on factors suchas root system development, worm activity, and soil shrinking,and may thus vary significantly within the studied area. On theother hand, the Ks values obtained using the tension disk infi ltrometer represent properties of the well-consolidated higher orderaggregates (the matrix soil structure) at the end of the vegetative
period, which reflect similar aggregation conditions and processes
characteristic for a particular soil type and climate. As a result,experiments performed under unsaturated conditions could not
aid in predicting Ks values characterizing saturated conditions ina soil with macropores.
The K s values obtained from the tension disk infiltrometer(Table 2) were on the same order of magnitude as the Ks value
obtained from the MSO experiment (Table 1). Th is similarity
was probably caused by the fact that both experiments were performed primarily under unsaturated conditions, that is, whenlarger gravitational pores were not involved in water flow. The
Ks values obtained from the Guelph permeameter experiments(Table 2) were closer to those obtained using laboratory ponding infi ltration (Table 1). Th is was probably due to the fact thatboth experiments involved saturated conditions. The Ks values
obtained from the Guelph permeameter, however, were one orderof magnitude lower than those for the ponding infi ltration. Th isdifference may be caused by the spatial and temporal variabilityof gravitational pores in the Ap horizon or the possible closing
or disruption of macropores during the drilling of the Guelph
permeameter well (which could occur despite the precautions
taken during the well preparation).

Dual-Permeability Model

Fig. 3. Final soil water contents (cm3 cm−3) for the Guelph permeameter test (Replicate 1),
assuming the single-porosity model and soil water retention parameters obtained from thelaboratory ponding infiltration (left) and the multistep outflow (right) experiments. Coordinates of the A, B, C, G, H, and I points in the radially symmetric system are in centimeters.

The optimized Ks values obtained for all scenarios using the dual-permeability model by
numerical inversion of data from the Guelph
permeameter (Replicate 1) are shown in Table 3.
Comparisons of Scenarios A with B and C withD show that larger calibrated Ks,f values wereobtained in all cases for the larger characteristic length of aggregates, a = 1.5 cm (Scenarios
B and D), than for the lower a = 0.2 cm value(Scenarios A and C). In the case of the larger
a value, a contact area between the matrix andmacropore domains was smaller than for the
lower a value. Therefore, for a = 1.5 cm, lesswater penetrated into the matrix domain fromthe macropore domain, and the wetting front
in the macropore domain moved faster (higher
Ks,f values) and farther from the source compared with a = 0.2 cm, while at the same time

www.VadoseZoneJournal.org | 220both scenarios produced the same cumulative infi ltration as measured (Fig. 3). In addition, the lower mass exchange between both
domains was compensated by the large direct cumulative infi ltration into the matrix domain (from the infi ltration well), as shown
for Ksa = 5.3 × 10−6 cm min−1 (scenarios with the lowest valuesof the objective function) in Fig. 4.
Different saturations of the matrix and macropore domains at theend of the simulated infi ltration (e.g., final wetting front
positions) for the same Ksa value (= 5.3 × 10−6 cm min−1)
and Scenarios C and D are clearly documented in Fig. 5
and 6. The fi nal soil water contents in the matrix andmacropore domains are shown in Fig. 5 and 6, respectively. The fi nal mass transfer between the macropore
and matrix domains is presented in Fig. 7. It should bepointed out that the total soil water content is equal
to the sum of the soil water contents of the matrix andmacropore domains multiplied by their corresponding
domain fractions (0.972 and 0.028). The figures showthat pressure heads, and consequently the soil water contents, were more or less in equilibrium in both domains
for scenarios with smaller aggregates and correspondingly
faster mass exchange. Therefore, mass exchange occurred
mainly at the wetting front of the macropore domain. On
the other hand, for scenarios with larger aggregates, nonequilibrium in pressure heads between the two domains
caused mass exchange throughout almost the entire flow
domain impacted by macropore flow. Comparisons of
Scenarios A with C and B with D show a very low impactof the aggregate shape parameter (due to its small rangeof possible values) on the calibrated Ks,f values. The Ks,fvalues for = 5 (except for scenarios with Ksa = 5.3 ×
10−5 cm min−1) were only slightly higher than those for= 15. Figure 4 shows, however, that evidently impacted
more mass exchange between the macropore and matrix
domains and the soil water content distribution for a =1.5 cm than for a = 0.2 cm.
Comparisons of Scenarios E, C, and F demonstrate theimpact of the  f parameter (0.05, 0.1, and 0.15 cm−1),
which is frequently interpreted as an inverse value of the
air-entry pressure head, on the calibrated Ks,f values. Thus,
a lower value of the  f parameter produces an increased
air-entry value and a correspondingly larger range of pressure heads for which the soil water contents and hydraulicconductivities are close to saturated values. Larger Ks,fvalues were thus calibrated for the macropore domaincharacterized by a higher  f to obtain the same measured cumulative infi ltration. The limited impact of the f value on simulated cumulative fluxes for Scenarios C,
E, and F with Ksa = 5.3 × 10−6 cm min−1 can be seen in
Fig. 4. Wetting fronts for the same scenarios (for Scenario
C shown in Fig. 5 [left] and 6 [left] and for Scenarios

E and F, which are not shown) moved in both regions to similar
distances from the infi ltration well. Water content distributions
at the wetting front were slightly more gradual for a lower valueof the  f parameter.
Comparisons of Scenarios G, C, and H show the impact of thenf parameter (= 2, 3, and 4) and the correspondingly increasing
degree of the curve S-shape (from gradual to a step-like shape) on

Fig. 4. Measured and simulated (using the dual-permeability model) cumulative
infiltrations for the Guelph permeameter test (Replicate 1) for various scenarios(Table 3) and an effective saturated hydraulic conductivity of the interface between
the two pore domains Ksa = 5.3 × 10−6 cm min−1. Cumulative macropore flux is the
direct (from the well) cumulative infiltration into the macropore domain; cumulative matrix flux is the direct cumulative infiltration into the matrix domain; sum of
cumulative fluxes is the sum of the direct cumulative fluxes; cumulative macropore/
matrix exchange is the cumulative infiltration from the macropore domain into thematrix domain.

www.VadoseZoneJournal.org | 221

Fig. 5. Final soil water contents (cm3 cm−3) in the matrix domain (simulated using
the dual-permeability model) for the Guelph permeameter test (Replicate 1) for
an effective saturated hydraulic conductivity of the interface between the twopore domains Ksa = 5.3 × 10−6 cm min−1 and for Scenarios B (left) and D (right).
Coordinates of the A, B, C, G, H, and I points in the radially symmetric system
are in centimeters.

Fig. 6. Final soil water contents (cm3 cm−3) in the macropore domain (simulatedusing the dual-permeability model) for the Guelph permeameter test (Replicate 1) for
an effective saturated hydraulic conductivity of the interface between the two poredomains Ksa = 5.3 × 10−6 cm min−1 and Scenarios B (left) and D (right). Coordinatesof the A, B, C, G, H, and I points in the radially symmetric system are in centimeters.

Fig. 7. Final mass transfer (min−1) between the macropore and matrix domains(simulated using the dual-permeability model) for the Guelph permeameter test
(Replicate 1) for an effective saturated hydraulic conductivity of the interface
between the two pore domains Ksa = 5.3 × 10−6 cm min−1 and Scenarios B (left)and D (right). Coordinates of the A, B, C, G, H, and I points in the radically symmetric system are in centimeters.

calibrated Ks,f values. As for increased air-entry pressure
heads, the more step-like shape of the retention curve (i.e.,higher nf values) increased the range of pressure heads withsoil water contents and hydraulic conductivities close tosaturated values. Larger calibrated Ks,f values were thus
calibrated for the macropore domain characterized withlower nf parameters. Again, Fig. 4 shows the limited impactof nf on simulated cumulative fluxes for Ksa = 5.3 × 10−6
cm min−1. Wetting fronts in both domains were again very
similar in all cases (for Scenario C shown in Fig. 5 [left]
and 6 [left] and Scenarios G and H, which are not shown).
The lowest Ks,f values were obtained for Scenario F and thelargest Ks,f values for Scenario G (for all Ksa values).
A comparison of results for different values of Ksa showsthat Ks,f increased with decreasing Ksa except for Scenarios
B and D (for a = 1.5 cm). Almost the same Ks,f values wereobtained for Ksa values of 5.3 × 10−7 and 5.3 × 10−8 cm
min−1. Evidently, the reduction in the mass exchange
between the macropore and matrix domains due to theaggregate size (discussed above) and low Ksa values was sohigh that cumulative infi ltration was no longer influenced
by the mass exchange. The impact of Ksa on particular
cumulative infi ltration fluxes is shown for Scenario A in
Fig. 8. Figure 8 shows that relatively high exchange rates(Ksa = 5.3 × 10 −5 and 5.3 × 10 −6 cm min−1) between
the macropore and matrix domains reduced the simulated cumulative infi ltration into the matrix domain andincreased the cumulative infi ltration into the macropore
domain. On the other hand, lower exchange rates (Ksa =5.3 × 10−7 and 5.3 × 10−8 cm min−1) increased the cumulative infi ltration into the matrix domain and decreased
the cumulative infi ltration into the macropore domain. As
with the impact of aggregate size on mass transfer, the wetting front in the macropore domain moved faster (due tohigher Ks,f values) and farther from the source for lower
values of Ksa than for higher values while producing thesame cumulative infi ltration as was measured.
The flow regimes for Scenario A and various Ksa values
are documented in Fig. 9, 10, and 11. While the final soil
water contents in the matrix and macropore domains areshown in Fig. 9 and 10, respectively, the final mass transfer
between the macropore and matrix domains is presented
in Fig. 11. Pressure heads, and consequently soil water contents, were more or less in equilibrium in both domains for
Ksa values of 5.3 × 10−5 and 5.3 × 10−6 cm min−1 and correspondingly faster mass exchange (Fig. 9 and 10, top). Mass
exchange for these scenarios took place close to the wetting
front in the macropore domain. Figures 9 and 10 (bottom)
for Ksa values of 5.3 × 10−7 and 5.3 × 10−8 cm min−1 show
similar water content distributions in the matrix domain
close to the Guelph permeameter well, which was largely

www.VadoseZoneJournal.org | 222

impacted by direct matrix infi ltration. The faster flow
and deeper percolation of the wetting front in themacropore domain, compared with the wetting front
in the matrix domain, are clearly visible for scenarioswith high Ksa values. While for the scenario with a
Ksa value of 5.3 × 10−7 cm min−1, the mass exchange
between the macropore and matrix domains caused asmall increase in water contents in the matrix domain
at larger distances from the well, no visible (at least on
the scale applied) increase in water contents was noticeable for the scenario with a Ksa value of 5.3 × 10−8 cm
min−1. Correspondingly, the final mass exchange wassimulated for almost the entire flow region impacted
by macropore flow in the former case, and almost no
visible mass exchange was simulated for the latter case.
It is evident from the values of the calculated objectivefunction that the best fit was obtained for Scenario
A and K sa = 5.3 × 10 −6 cm min−1 ( = 101.9).
Additionally, a considerably better fit was obtainedusing the dual-permeability model ( = 101.9) than
the single-porosity model with the PI parameters (= 835.0). It should be noted, however, that a better fit
could be obtained for the single-porosity model using
different  r,  s, , and n parameters. For instance, the value for the single-porosity model with the MSO
parameters was only 192.4, but such gradual soil water
retention curves are not characteristic for wetting
experiments (Gribb et al., 2004). After comparing thewater content distributions around the infiltration
well simulated using the dual-permeability model (Fig.
9 and 10) with the results of the single-porosity model(Fig. 3), it is obvious that macropore flow and mass
exchange between the two domains produced more
uniform water content distributions.
Close correlation among parameters , a, and Ksa
reflects their direct effect on mass transfer (Eq. [3]).
Independent evaluation of aggregate shape parameters
is recommended because new imaging techniques, suchas x-ray computed tomography, magnetic resonance,
dye tracer, and micromorphological imaging (Köhne etal., 2009a; Kodešová, 2009), provide relatively reliable
information about the soil structure. The Ksa value may
also be independently evaluated. For instance, Gerke
and Köhne (2002) measured the hydraulic propertiesof the aggregate skin from different adsorption rates
between the aggregates with and without the skin. Thismethod requires a spatial apparatus, however, and may
be applied only for aggregates with well-developed
skins. In addition, optimized Ksa values may compensate for the possible inaccuracies of independently
determined a and values.

Fig. 8. Measured and simulated (using the dual-permeability model) cumulative
infiltrations for the Guelph permeameter test (Replicate 1) for an effective saturated
hydraulic conductivity of the interface between the two pore domains Ksa = 5.3 ×
10−5 (top left), Ksa = 5.3 × 10−6 (top right), Ksa = 5.3 × 10−7 (bottom left), and Ksa =5.3 × 10−8 cm min−1 (bottom right). Cumulative macropore flux is the direct (fromthe well) cumulative infiltration into the macropore domain; cumulative matrix flux isthe direct cumulative infiltration into the matrix domain; sum of cumulative fluxes is
the sum of the direct cumulative fluxes; cumulative macropore/matrix exchange is the
cumulative infiltration from the macropore domain into the matrix domain.

Fig. 9. Final soil water contents (cm3 cm−3) in the matrix domain (simulated using thedual-permeability model) for the Guelph permeameter test (Replicate 1) for an effective saturated hydraulic conductivity of the interface between the two pore domains
Ksa = 5.3 × 10−5 (top left), Ksa = 5.3 × 10−6 (top right), Ksa = 5.3 × 10−7 (bottom
left), and Ksa = 5.3 × 10−8 cm min−1 (bottom right), and for Scenario A. Coordinatesof the A, B, C, G, H, and I points in the radially symmetric system are in centimeters.

www.VadoseZoneJournal.org | 223

6 Conclusions
Tension disk infiltrometer and Guelph permeameter tests
were used to estimate the saturated hydraulic conductivity, Ks,
using the single-porosity model. It was shown that different
soil water retention curves, obtained by either drying or wetting laboratory experiments, impacted the saturated hydraulicconductivities optimized from the field experiments. The Ks
values obtained using parameters from the ponding infi ltration experiment were suggested to better characterize the field
experiments because this experiment represented a wetting process, and thus higher consistency was expected between them.
In addition, greater Ks values were obtained from the Guelph
permeameter tests than from the tension disk infi ltrometer
tests due to an occurrence of preferential flow caused by thegravitational pores. Therefore, the Guelph permeameter tests
were used to estimate the saturated hydraulic conductivity ofthe macropore domain, Ksf, using the dual-permeability model.
Fig. 10. Final soil water contents (cm3 cm−3) in the macropore domain (simulated using the dual-permeability model) for the Guelph permeameter test
(Replicate 1) for an effective saturated hydraulic conductivity of the interface
between the two pore domains Ksa = 5.3 × 10−5 (top left), Ksa = 5.3 × 10−6
(top right), Ksa = 5.3 × 10−7 (bottom left), and Ksa = 5.3 × 10−8 cm min−1(bottom right), and for Scenario A. Coordinates of the A, B, C, G, H, and I
points in the radially symmetric system are in centimeters.

Fig. 11. Final mass transfer (min−1) between the macropore and matrix
domains (simulated using the dual-permeability model) for the Guelph permeameter test (Replicate 1) for an effective saturated hydraulic conductivityof the interface between the two pore domains Ksa = 5.3 × 10−5 (top left),
Ksa = 5.3 × 10−6 (top right), Ksa = 5.3 × 10−7 (bottom left), and Ksa = 5.3
× 10−8 cm min−1 (bottom right) and for Scenario A. Coordinates of the A,
B, C, G, H, and I points in the radially symmetric system are in centimeters.

The soil hydraulic parameters determined from the tension disk
infi ltrometer, obtained using the single-porosity model, wereused as parameters for characterizing the matrix domain in the
dual-permeability model. Since many parameters characterizing
the dual-domain geometry, water retention in the macropore
domain, and the effective saturated hydraulic conductivity ofthe domain interface were unknown and could not be estimatedsimultaneously, we have evaluated multiple scenarios with different fixed parameter sets. The use of these different parameter sets
produced different estimates of Ksf. While the effective saturated
hydraulic conductivities, Ksa, characterizing the exchange term
hydraulic conductivity differed by orders of magnitude, however,
the estimated saturated hydraulic conductivities of the macropore domain, Ksf, varied only within one order of magnitude(that is, within the natural range of spatial heterogeneity for thesaturated hydraulic conductivity [Nielsen and Wendroth, 2003]).
In addition, although simple infiltration experiments without
additional measurements of transient flow data such as watercontents or pressure heads in the soil are not expected to provide
enough information for estimating other parameters (Šimůnek
and van Genuchten, 1996), the objective function showed highsensitivity to at least two estimated parameters.
It is evident that the shape parameters of the soil water retention curve of the macropore domain influenced the estimatedK sf value. When no information is available about theseparameters, however, the average  f and nf values may serve as
reference values. It is possible that additional measurements ofpressure heads in the soil may help improve parameter estimation. The impact of the aggregate shape parameter, , was, incomparison with the effect of other parameters, very small. The
characteristic length of aggregates, a, noticeably influenced theestimated Ksf values. Decreasing Ksa and increasing a similarly
decreased mass transfer between the macropore and matrix

www.VadoseZoneJournal.org | 224domains and subsequently caused a large nonequilibrium between
the two domains and a large irregularity at the wetting front. Our
results document the applicability of the approach presented forthe assessment of some of these parameters (preferably Ksf and Ksa)of the dual-permeability model.

Acknowledgments
We acknowledge the financial support of the Grant Agency of the Czech Republic (Grant no. GA CR 526/08/0434) and the Ministry of Education, Youth andSports of the Czech Republic (Grant no. MSM 6046070901). We also thank KarelNěmeček and Marcela Mühlhanselová for their help with the field work. We would
like to thank to Zdeněk Tomíček and his colleagues from the experimental stationof the Crop Research Institute in Hněvčeves for their support.

References
Ahuja, L.R., and C. Hebson. 1992. Root zone water quality model. GPSR Tech.
Rep. 2. USDA-ARS Great Plains Syst. Res., Fort Collins, CO.
Chen, C., and R.J. Wagenet. 1992a. Simulaon of water and chemicals inmacropore soils: 1. Representaon of the equivalent macropore inuence and its eect on soil water ow. J. Hydrol. 130:105–126.
Chen, C., and R.J. Wagenet. 1992b. Simulaon of water and chemicals
in macropore soils: 2. Applicaon of linear lter theory. J. Hydrol.
130:127–149.
Clothier, B.E., S.R. Green, and M. Deurer. 2008. Preferenal ow and transport in soil: Progress and prognosis. Eur. J. Soil Sci. 59:2–13.
Freifeld, B.M., and C.M. Oldenburg. 2000. Restricted interval Guelph permeameter: Theory and applicaon. Water Resour. Res. 36:1373–1380.
Gerke, H.H. 2006. Preferenal ow descripons for structured soils. J. Plant
Nutr. Soil Sci. 169:382–400.
Gerke, H.H., J. Dušek, T. Vogel, and J.M. Köhne. 2007. Two-dimensional dualpermeability analyses of a bromide tracer experiment on a le-drainedeld. Vadose Zone J. 6:651–667.
Gerke, H.H., and J.M. Köhne. 2002. Esmang hydraulic properes of soil
aggregate skins from sorpvity and water retenon. Soil Sci. Soc. Am. J.
66:26–36.
Gerke, H.H., and J.M. Köhne. 2004. Dual-permeability modeling of preferenal bromide leaching from a le-drained glacial ll agricultural eld. J.
Hydrol. 289:239–257.
Gerke, H.H., and M.Th. van Genuchten. 1993. A dual-porosity model for
simulang the preferenal movement of water and solutes in structured
porous media. Water Resour. Res. 29:305–319.
Gerke, H.H., and M.Th. van Genuchten. 1996. Macroscopic representaon of
structural geometry for simulang water and solute movement in dualporosity media. Adv. Water Resour. 19:343–357.
Germann, P.F. 1985. Kinemac wave approach to inltraon and drainage
into and from soil macropores. Trans. ASAE 28:745–749.
Germann, P. 1990. Preferenal ow and generaon of runo: 1. Boundary
layer ow theory. Water Resour. Res. 26:3055–3063.
Germann, P.F., and K. Beven. 1985. Kinemac wave approximaon to inltraon into soils with sorbing macropores. Water Resour. Res. 21:990–996.
Gribb, M.M., R. Kodešová, and S.E. Ordway. 2004. Comparison of soil hydraulic property measurement methods. J. Geotech. Geoenviron. Eng.
130:1084–1095.
Hayashi, M., and W.L. Quinton. 2004. A constant-head well permeameter
method for measuring eld-saturated hydraulic conducvity above animpermeable layer. Can. J. Soil Sci. 3:255–264.
Jarvis, N.J. 1994. The MACRO model: Technical descripon and sample simulaon. Rep. Diss. 19. Dep. of Soil Sci., Swedish Univ. of Agric. Sci., Uppsala.
Jarvis, N.J. 2007. Review of non-equilibrium water and solute transport insoil macropores: Principles, controlling factors and consequences for water quality. Eur. J. Soil Sci. 58:523–546.
Kodešová, R. 2009. Soil micromorphology use for modeling of a non-equilibrium water and solute movement. Plant Soil Environ. 55:424–428.
Kodešová, R., M.M. Gribb, and J. Šimnek. 1998. Esmang soil hydraulicproperes from transient cone permeameter data. Soil Sci. 163:436–453.
Kodešová, R., M. Ko árek, V. Kodeš, J. Šimnek, and J. Kozák. 2008. Impact ofsoil micromorphological features on water ow and herbicide transport
in soils. Vadose Zone J. 7:798–809.
Kodešová, R., V. Kodeš, A. Žigová, and J. Šimnek. 2006a. Impact of plantroots and soil organisms on soil micromorphology and hydraulic properes. Biologia 61:S339–S343.
Kodešová, R., J. Kozák, and J. Šimnek. 2006b. Numerical study of macropore
impact on ponded inltraon in clay soils. Soil Water Res. 1:16–22.

Kodešová, R., S.E. Ordway, M.M. Gribb, and J. Šimnek. 1999. Esmaon ofsoil hydraulic properes with the cone permeameter: Field Studies. Soil
Sci. 164:527–541.
Kodešová, R., N. Vignozzi, M. Rohošková, T. Hájková, M. Ko árek, M. Pagliai,
J. Kozák, and J. Šimnek. 2009. Impact of varying soil structure on transport processes in dierent diagnosc horizons of three soil types. J. Contam. Hydrol. 104:107–125.
Köhne, J.M., S. Köhne, and J. Šimnek. 2006a. Mul-process herbicide transport in structured soil columns: Experiments and model analysis. J. Contam. Hydrol. 85:1–32.
Köhne, J.M., S. Köhne, and J. Šimnek. 2009a. A review of model applicaons for structured soils: a) Water ow and tracer transport. J. Contam.
Hydrol. 104:4–35.
Köhne, J.M., S. Köhne, and J. Šimnek. 2009b. A review of model applicaons
for structured soils: b) Pescide transport. J. Contam. Hydrol. 104:36–60.
Köhne, J.M., B.P. Mohanty, and J. Šimnek. 2006b. Inverse dual-permeability
modeling of preferenal water ow in a soil column and implicaons foreld-scale solute transport. Vadose Zone J. 5:59–76.
Kool, J.B., and J.C. Parker. 1987. Development and evaluaon of closed-form
expressions for hysterec soil hydraulic properes. Water Resour. Res.
23:105–114.
Mualem, Y. 1976. A new model for predicng the hydraulic conducvity ofunsaturated porous media. Water Resour. Res. 12:513–522.
Nielsen, D.R., and O. Wendroth. 2003. Spaal and temporal stascs. Catena-Verlag, Reiskirchen, Germany.
Othmer, H., B. Diekkrüger, and M. Ku lek. 1991. Bimodal porosity and unsaturated hydraulic conducvity. Soil Sci. 152:139–149.
Philip, J.R. 1968. The theory of adsorpon in aggregated media. Aust. J. Soil
Res. 6:1–19.
Pot, V., J. Šimnek, P. Benoit, Y. Coquet, A. Yra, and M.-J. Marnez Cordun.
2005. Impact of rainfall intensity on the transport of two herbicides inundisturbed grassed lter strip soil cores. J. Contam. Hydrol. 81:63–88.
Ramos, T.B., M.C. Gonçalves, J.C. Marns, M.Th. van Genuchten, and F.P.
Pires. 2006. Esmaon of soil hydraulic properes from numerical inversion of tension disk inltrometer data. Vadose Zone J. 5:684–696.
Sander, T., and H.H. Gerke. 2007. Preferenal ow pa erns in paddy elds
using a dye tracer. Vadose Zone J. 6:105–115.
Schwartz, R.C., and S.R. Eve . 2002. Esmang hydraulic properes of a netextured soil using a disc inltrometer. Soil Sci. Soc. Am. J. 66:1409–1423.
Šimnek, J., N.J. Jarvis, M.Th. van Genuchten, and A. Gärdenäs. 2003. Review
and comparison of models for describing non-equilibrium and preferenal ow and transport in the vadose zone. J. Hydrol. 272:14–35.
Šimnek, J., R. Kodešová, M.M. Gribb, and M.Th. van Genuchten. 1999. Esmang hysteresis in the soil water retenon funcon from cone permeameter test data. Water Resour. Res. 35:1329–1345.
Šimnek, J., and M.Th. van Genuchten. 1996. Esmang unsaturated soil
hydraulic properes from tension disc inltrometer data by numericalinversion. Water Resour. Res. 32:2683–2696.
Šimnek, J., and M.Th. van Genuchten. 1997. Esmang unsaturated soil
hydraulic properes from mulple tension disc inltrometer data. Soil
Sci. 162:383–398.
Šimnek, J., and M.Th. van Genuchten. 2008. Modeling nonequilibrium ow
and transport with HYDRUS. Vadose Zone J. 7:782–797.
Šimnek, J., M.Th. van Genuchten, and M. Šejna. 2008. Development and applicaons of the HYDRUS and STANMOD so ware packages and related
codes. Vadose Zone J. 7:587–600.
Soilmoisture Equipment Corp. 2008. Model 2800K1 Guelph Permeameter
Operang Instrucons. Soilmoisture Equipment Corp., Santa Barbara, CA.
van Genuchten, M.Th. 1980. A closed-form equaon for predicng the hydraulic conducvity of unsaturated soils. Soil Sci. Soc. Am. J. 44:892–898.
Vogel, T., H.H. Gerke, R. Zhang, and M.Th. van Genuchten. 2000. Modelingow and transport in a two-dimensional dual-permeability system withspaally variable hydraulic properes. J. Hydrol. 238:78–89.
Watson, K.W., and R.J. Luxmoore. 1986. Esmang macroporosity in a forest
watershed by use of a tension inltrometer. Soil Sci. Soc. Am. J. 50:578–582.
Wooding, R.A. 1968. Steady inltraon from a shallow circular pond. Water
Resour. Res. 4:1259–1273.
Zhang, Z.F., P.H. Groenevelt, and G.W. Parkin. 1998. The well-shape factor forthe measurement of soil hydraulic properes using the Guelph permeameter. Soil Tillage Res. 49:219–221.
Zimmerman, R.W., T. Hadgu, and G.S. Bodvarsson. 1996. A new lumped-parameter model for ow in unsaturated dual-porosity media. Adv. Water
Resour. 19:317–327.

www.VadoseZoneJournal.org | 225