6.6 Solute Transport During Variably Saturated Flow—
Inverse MethodsJIÌÍ ŠIMæNEK AND
MARTINUS TH. VAN GENUCHTEN, USDA-ARS, George E. Brown, Jr. Salinity
Laboratory, Riverside, California
DIEDERIK JACQUES, Belgian Nuclear Research Centre, Boeretang, Belgium
JAN W. HOPMANS, University of California, Davis, California
MITSUHIRO INOUE, Arid Land Research Centre, Tottori University, Tottori, JapanMARKUS FLURY, Washington State University, Pullman, Washington

6.6.1 Introduction
The use of parameter estimation techniques for determining soil hydraulic properties is well established (Kool et al., 1987; Hopmans & Šimçnek, 1999; Sections 1.7
and 3.6.2). The approach has been widely used for various laboratory and field experiments. Among others, laboratory experiments include one-step and multistep
outflow experiments, upward flux or head-controlled infiltration, the evaporation
method, and others (Section 3.6.2). In separate lines of research, solute transportparameters are often obtained from column experiments assuming steady-state waterflow (e.g., Nkedi-Kizza et al., 1984; Sections 6.3–6.5), and using parameter estimation codes such as CFITIM (van Genuchten, 1981), CXTFIT (Toride et al., 1995),
or STANMOD (Šimçnek et al., 1999) for fitting analytical solutions of the transport equation to experimental breakthrough curves. Obtaining solute transport parameters for conditions for which no analytical solutions exist, such as for nonlinear adsorption, can be accomplished using numerical solutions (Kool et al., 1989;Šimçnek & van Genuchten, 1999).
The above parameter estimation efforts for water flow and solute transport
have thus far remained relatively disjointed. Although numerous studies exist thatcombine estimation of flow and transport parameters for groundwater flow problems (e.g., Sun & Yeh, 1990; Medina & Carrera, 1996; Weiss & Smith, 1998), only
a few studies have used combined transient variably saturated water flow andsolute transport experiments for simultaneous estimation of soil hydraulic andsolute transport parameters (Mishra & Parker, 1989; Abbaspour et al., 1997; Inoueet al., 2000; Jacques et al., 2001).
Different strategies in combined estimation of water flow and solute transport parameters can be followed (Inoue et al., 2000). First, one could use water flow
information only (e.g., pressure heads and/or fluxes) to estimate the soil hydraulicparameters, followed by estimation of the transport parameters using information
1435

1436

CHAPTER 6

from the transport part of the experiment (e.g., solute concentrations). Alternatively,
combined water flow and transport information can be used to estimate soil hydraulicand solute transport parameters in a sequential manner. Finally, combined water flow
and transport information can be used to simultaneously estimate both the soil hydraulic and solute transport parameters. This last approach is the most beneficial
since it uses crossover effects between state variables and parameters (Sun & Yeh,1990), and it takes advantage of all available information since concentrations area function of water flow (Medina & Carrera, 1996). Mishra and Parker (1989)
showed that simultaneous estimation of hydraulic and transport properties yields
smaller estimation errors for model parameters than sequential inversion of hydraulicproperties from water content and matric pressure head data, followed by inversionof transport properties from concentration data.
In this section we focus on the application of parameter estimation methods
to indirect estimation of the soil hydraulic and solute transport parameters from variably saturated transient flow experiments. Since the approach requires the use ofnumerical methods, we also include a relatively simple example involving transport of solute with nonlinear chemical reaction (during steady water flow) for
which no analytical solution is available. Parameter estimation methods coupled withcomprehensive numerical models can provide extremely useful tools for analyzing
experimental data that may not be evaluated optimally using conventional tools suchas analytical solutions (Sections 6.3 through 6.5). Because of the generality of numerical models, these tools are very attractive for analyzing conveniently and accurately a broad range of steady-state and transient laboratory and field flow andtransport experiments.
General reviews of the application of inverse modeling was presented by Hopmans and Šimçnek (1999) and in Section 1.7. Specific applications to variably saturated flow problems only are discussed in Section 3.6.2. We will briefly present
the governing equations for simultaneous variably saturated water flow and nonlinear nonequilibrium solute transport, define the inverse problem, and then givefour examples demonstrating the application of parameter estimation techniques tosolute transport problems of increasing complexity. The first example deals withnonlinear solute transport in a laboratory column. The second example applies theparameter estimation method to a laboratory column experiment that used a flowinterruption technique to evaluate rate-limited mass transfer. The third example deals
with variably saturated water flow and solute transport during an infiltration experiment in a laboratory column. Finally, the fourth example describes a relatively
complex field experiment involving time-dependent boundary conditions, a layeredsoil profile, and nonequilibrium solute transport. We will devote little attention toquestions of uniqueness, stability, parameter confidence intervals, and strengths andweaknesses of the approach, since these topics are extensively discussed in Sections
1.7 and 3.6.2.
6.6.2 Theory of Flow, Transport, and Optimization
Isothermal one-dimensional Darcian water flow in a variably saturated rigid
isotropic porous medium is given by the following modified form of the Richards
equation:

MISCIBLE SOLUTE TRANSPORT

1437

∂θ
∂q
∂
∂h
__ = − __− S = __ K __ + 1 − S∂t∂z∂z — ‰ ∂z

[6.6–1]

where θ is the volumetric water content, h is the pressure head, q is the volumetricflux density, z is the vertical coordinate positive upwards, t is time, S is a sink term,and K is the unsaturated hydraulic conductivity function. Various analytical expressions can be used to describe the retention and hydraulic conductivity functions(Section 3.3.4).
Solute transport during transient water flow in a variably saturated rigid
porous medium is given by:∂θc ___∂ρs __
∂
∂c
∂qc
___
+
=
θD __ − ___ − Scr − µwθc − µsρs + γwθ + γsρ∂t∂t∂z ‰ ∂z ∂z[6.6–2]where c and s are solute concentrations in the liquid and solid phases, respectively;
ρ is the soil bulk density, µw and µs are first-order rate constants for solutes in theliquid and solid phases, respectively; γw and γs are zero-order rate constants for theliquid and solid phases, respectively; S is the sink term in the flow equation (Eq.
[6.6–1]), cr is the concentration of the sink term, and D is the dispersion coefficientfor the liquid phase:
D = Dwτ + λ|v|

[6.6–3]where Dw is the molecular diffusion coefficient, τ is a tortuosity factor, λ is the longitudinal dispersivity, and v is the average pore water velocity (= q/θ).
Interactions between the solution (c) and adsorbed (s) concentrations can be
described in terms of equilibrium and nonequilibrium models. For equilibrium adsorption, the isotherm relating s and c may be described by a generalized nonlinear equation of the forms = kdcβ/(1 + ηcβ)[6.6–4]

where kd, β, and η are empirical coefficients. The Freundlich, Langmuir, and linear isotherms are special cases of Eq. [6.6–4].
Nonequilibrium sorption is usually described by a mass transfer reaction
where the transfer between solutes in the solution and the sorbed phases is governedby the concentration differences and a mass transfer coefficient αmt:∂s/∂t = αmt[kdcβ/(1 + ηcβ) − S]

[6.6–5]It has often been observed that sorption sites can be separated into two types. Type1 sorption sites (fraction f) undergo instantaneous sorption, described by Eq.
[6.6–4], whereas the remaining type-2 sorption sites (fraction 1 − f) undergo nonequilibrium sorption, described by Eq. [6.6–5]. This concept is known as the twosite sorption model (e.g., Selim et al., 1977; van Genuchten & Wagenet, 1989). The
general mass balance equation for the type-2 sorption sites is then given by

1438

CHAPTER 6∂sk
kdcβ − s___= αmt (1 − f) ______k∂t—
1 + ηcβ[6.6–6]

where sk is the sorbed phase concentration at the kinetically controlled sorption sites.
Physical nonequilibrium transport may also exist (Section 6.3). Physicalnonequilibrium may be implemented using the concept of two-region dual-porosity type transport (van Genuchten & Wierenga, 1976) by assuming that the liquid
phase can be partitioned into mobile (flowing), θm, and immobile (stagnant), θim,regions, and that solute exchange between the two liquid regions can be modeled
as a first-order process; that is,
∂θ
∂f′ρs
∂
∂cm
∂qmcm
mcm
_____
+ _____m = __ θmDm ___− _____− Scr − α′mt (cm − cim)
∂t
∂t
∂z ‰
∂z 
∂z
β−1
kdβcim∂cim____
θim + ρ(1 − f′) _________= α′mt (cm − cim)
β 2—
(1 + ηcim)  ∂t

[6.6–7][6.6–8]where cim is the concentration of the immobile region, α′mt the mass transfer coefficient for exchange between the two liquid regions, f′ is the fraction of exchangesites in contact with the immobile solution phase, and the subscript “m” refers to
the mobile region. Notice that the parameters f′ and α′mt in Eq. [6.6–8] have a slightly
different meaning than f and αmt in Eq. [6.6–6].
Equations [6.6–1] through [6.6–8] must be solved numerically for variably
saturated flow and/or nonlinear solute transport problems for a given set of initialand boundary equations. When the van Genuchten (1980) analytical expression forthe soil hydraulic properties is used, the optimized soil hydraulic parameters may
include the residual, θr, and saturated, θs, water contents, the saturated hydraulicconductivity, Ks, and the shape parameters n, α, and l (Section 3.3.4). The optimized
solute transport parameters for solute transport as described with Eq. [6.6–2]
through [6.6–8] could include the following parameters: Dw, f, f′, kd, β, θim, η, αmt,
α′mt, τ, and/or λ.
The objective function, Φ, to be minimized during the parameter estimation
process may be defined as
m

nj

j=1

i=1

Φ(b,y) = Σvj Σ wi,j[yj*(z,ti) − yj(z,ti,b)]2[6.6–9]where the right-hand side represents the residuals between the measured (yj*) andcorresponding model predicted (yj) space–time variables, using the vector of optimized soil hydraulic and solute transport parameters, b. The first summation sign
sums residuals for different measurement types, m, whereas the variable nj in thesecond summation denotes the number of measurements for a certain measurement
type. The factors w and v are weighting factors for the individual measurements andthe measurement types, respectively. For instance, w can represent the standard deviation or analytical errors for measurements. We assume weighting factors w tobe equal to unity, and factors v are defined as follows:
vj = 1/(njσj2)[6.6–10]MISCIBLE SOLUTE TRANSPORT

1439

where σj2 is the variance of the measurements yj*. Measured values can include waterflow variables—such as water contents, pressure heads, and/or water fluxes; solutetransport variables—such as concentrations and/or solute fluxes; or integral variables—such as electrical conductivity. Minimization of the objective function (Eq.
[6.6–9]) can be carried out using optimization techniques discussed in Section 1.7.
6.6.3 Examples
We demonstrate the use of the parameter estimation method with four examples of increasing complexity, using the variably saturated water flow and solutetransport numerical model HYDRUS-1D (Šimçnek et al., 1999) for each example.
However, any suitable numerical flow and transport code could be used for this purpose.
The first example deals with transport of solute with nonlinear chemical reaction in a fully saturated laboratory column. Parameters of the nonlinear Freundlich
adsorption isotherm are obtained by analyzing a measured breakthrough curve. In
the second example, we applied parameter estimation to a laboratory column experiment that used the flow interruption technique to evaluate rate-limited masstransfer. The third example involves variably saturated water flow and equilibrium
solute transport during infiltration in a laboratory column. The soil hydraulic andsolute transport parameters are optimized simultaneously in this example. The
fourth and final example uses data from a relatively complex field experiment involving atmospheric boundary conditions, a layered soil profile, and nonequilibrium solute transport. The soils hydraulic and solute transport parameters for individual soil layers are obtained sequentially. This example considers nonequilibrium
solute transport using the mobile–immobile dual porosity concept.
6.6.3.1 Steady-State Laboratory Flow Experiment with Nonlinear Transport
The first example demonstrates the use of the inverse procedure for the estimation of nonlinear sorption parameters from a column breakthrough curve reported
by Selim et al. (1987). A 10.75-cm-long soil column was first saturated with a 0.005
M CaCl2 solution. The transport experiment consisted of applying a 14.26 pore volume pulse (t = 358.05 h) of a 0.005 M MgCl2 solution, followed by the original CaCl2
solution. The flow rate was equal to 6.495 cm d−1. The adsorption isotherm was determined independently with the help of batch experiments (Selim et al., 1987), andfitted with the Freundlich equation, to yield kd = 1.69 cm3 g−1 and β = 1.62.
Only the dispersion coefficient (D) and the Freundlich coefficients (Eq.
[6.6–4] with η = 0) kd and β were optimized in this example. Since the governing
solute transport equation is nonlinear, no analytical solution is available, so that one
must resort to a numerical model. The soil profile was discretized using 44 nodes,
and a third-type boundary condition was used at the top of the column. The observed
and fitted Mg breakthrough curves are shown in Fig. 6.6–1. The results indicate areasonable prediction of the measured breakthrough curve for the final estimatesof the optimized solute transport parameters (λ = 3.76 ± 1.40 cm, kd = 1.44 ± 1.14
cm3 g−1, and β = 1.62 ± 0.35). A very high negative correlation of r 2 = −0.99 was
found between the sorption parameters kd and β.

1440

CHAPTER 6

Fig. 6.6–1. Measured and fitted Mg breakthrough curves for Abist loam (data taken from Selim et al.,1987).

6.6.3.2 Laboratory Transport Subject to Flow Interruption
In the second example, we analyzed a column breakthrough curve (Fortin etal., 1997) obtained using the flow interruption technique. Flow interruption techniques are often used to elucidate rate-limited sorption processes (Murali & Aylmore, 1980; Brusseau & Rao, 1989; Fortin et al., 1997). This particular experiment
was carried out in a 15-cm-long repacked soil column (loamy sand). The bulk density was 1.41 g cm−3 and the saturated water content 0.47 cm3 cm−3. The experiment was performed under steady-state, fully saturated flow conditions in a vertical direction at a flow rate of 0.674 cm h−1, corresponding to a pore water velocity
v of 1.43 cm h−1. Bromide and the herbicide simazine (6-chloro-N,N′-diethyl1,3,5-triazine-2,4-diamine) were dissolved in a 0.01 M CaSO4 solution at concentrations of 50 mg L−1 and 0.025 µg L−1, respectively, and introduced into the column as step inputs. Flow was stopped for 185 h after five pore volumes of input.
Fortin et al. (1997) analyzed the bromide and simazine breakthrough curves using
equilibrium and nonequilibrium solute transport models, respectively. They used
three nonequilibrium models (i.e, a two-stage, two-rate; a two-stage, one-rate; andrate-limited sorption models) and concluded that the two-stage, two-rate model, in
which the first rate proceeds considerably faster than the second one, could describesimazine breakthrough data the best.
Here we illustrate the parameter optimization procedure using the two-site
sorption model (Eq. [6.6–2] and [6.6–6]). First, the bromide breakthrough curves
were used to estimate the hydrodynamic dispersion of the column system. Flow interruption experiments indicated that bromide was not subject to chemical or physical nonequilibrium processes. The fitting of the steady-state solute transport equation to the bromide breakthrough data yielded a dispersivity of λ = 0.341 ± 0.011
cm (bromide breakthrough curve not shown). This dispersivity was then used as afixed parameter in the subsequent estimation of the sorption parameters. It is as-

MISCIBLE SOLUTE TRANSPORT

1441

Fig. 6.6–2. Measured and optimized simazine breakthrough curves (vs. pore volume and real time) obtained using the flow interruption technique. Optimization I uses dispersivity estimated from the bromide effluent data; Optimization II fits dispersivity independently (data taken from Fortin et al., 1997).

sumed that, under the experimental conditions, the sorption isotherm for simazine
is linear (Fortin et al., 1997). Figure 6.6–2 shows relative concentration values (vs.
pore volume and in real time) of the two-site sorption model fitted to the measured
simazine breakthrough curve. The fitted simulation (Optimization I) shows a reasonable representation of the simazine breakthrough for the final solute transportparameters (kd = 0.814 ± 0.062 cm3 g−1, f = 0.533 ± 0.046, and αmt = 0.0121 ± 0.0044h−1). Correlation coefficients, r2, among optimized parameters were below 0.7. Note,however, that the spreading of the breakthrough curve is not very well captured by

1442

CHAPTER 6

the model simulation. Figure 6.6–2 also shows results (Optimization II) when four
solute transport parameters, including dispersivity, were optimized (λ = 1.08 ± 0.246
cm, kd = 1.49 ± 2.57 cm3 g− 1, f = 0.328 ± 0.564, and αmt = 0.00128 ± 0.00406 h−1).
The data are much better represented by Optimization II; however, the increased
number of fitting parameters resulted in very high correlation coefficients amongkd, f, and αmt (r 2 = ± 0.99), indicating low confidence of these parameters (large
confidence intervals). More details on the parameter optimization using this twosite two-rate model can be found in Fortin et al. (1997), but the general proceduresof the parameter estimation are identical as employed in the example presented here.
6.6.3.3 Transient Laboratory Experiment with Equilibrium
Solute Transport
The third example demonstrates the simultaneous estimation of soil hydraulic and solute transport parameters from a column infiltration experiment in a
30-cm-long and 5-cm-i.d. laboratory soil column (Inoue et al., 2000), containing arepacked coarse-textured soil (Tottori sand). The infiltration experiment was carried out by simultaneously increasing the solute concentration (from 0.02 to 0.1 mol
L−1 NaCl) and the infiltration rate (from 0.00032 to 0.0026 cm s−1). Matric pressure heads, h, and bulk soil electrical conductivities, σa, were measured using automated minitensiometer and four-electrode sensors, respectively, at the 23-cm
depth. The measured bulk soil electrical conductivity is a variable that integrates
information on both water flow and solute transport, and can thus be beneficiallyused to simultaneously estimate soil hydraulic and solute transport parameters.
The soil hydraulic (van Genuchten, 1980) and solute transport parameters
were estimated both sequentially and simultaneously. The objective function for thesequential optimization was defined first in terms of pressure heads to estimate soil
hydraulic parameters (θr, θs, Ks, n, α, l; Section 3.3.4) and then in terms of electrical conductivities to estimate solute transport parameters (λ). For the simultaneous
optimization of both soil hydraulic and solute transport parameters, the objectivefunction was defined in terms of both measured pressure heads and electrical conductivities. The results for the sequential optimization are shown in Fig. 6.6–3. It
is specifically noted that the bulk electrical conductivity data show two separate
fronts, demonstrating that water and solute fronts travel at different velocities depending on the ratio of initial and final water content values. Figure 6.6–3 also shows
excellent correspondence between both optimized and measured pressure head andelectrical conductivity values.
The sequentially estimated soil hydraulic and solute transport parameters(Table 6.6–1) were compared with independently measured solute dispersion (Fig.
6.6–4), soil water retention, and unsaturated hydraulic conductivity (Fig. 6.6–5) data.
Notice from Table 6.6–1, that although the final parameter estimates for both sequential and simultaneous optimizations are very similar, the confidence intervals
for the simultaneous optimization of soil hydraulic and solute transport parameters,
as well as the final value of the objective function, are smaller. Dispersion coefficients (data points in Fig. 6.6–4) were obtained independently by fitting breakthrough curves derived from steady-state water flow column experiments (Inoueet al., 2000). The retention data of Fig. 6.6–5 were obtained from an additional in-

MISCIBLE SOLUTE TRANSPORT

1443

Fig. 6.6–3. Measured and optimized pressure heads and bulk electrical conductivities for a column infiltration experiment (data taken from Inoue et al., 2000).

Table 6.6–1. Soil hydraulic and solute transport parameters and their confidence intervals obtained using
sequential and simultaneous optimizations.
Parameterθr (cm3 cm−3)
θs (cm3 cm−3)
Ks (cm s−1)n (-)α (cm−1)l (-)λ (cm)
Φ

Sequential optimization

Simultaneous optimization

0.0265 ± 0.0286
0.310 ± 0.0386
0.138 ± 0.038
2.014 ± 0.318
0.0446 ± 0.0388−1.15 ± 0.178
0.221 ± 0.024
0.003470.0206 ± 0.0198
0.310 ± 0.0248
0.137 ± 0.024
1.969 ± 0.186
0.0570 ± 0.0322−0.816 ± 0.220
0.207 ± 0.020
0.00170

Fig. 6.6–4. Dispersion coefficient D as function of the mean pore water velocity v obtained by inverse
optimization and from analysis of steady-state data (data taken from Inoue et al., 2000).

1444

CHAPTER 6

Fig. 6.6–5. Comparison between soil hydraulic properties estimated inversely (lines) and obtained fromanalysis of independent infiltration data (data points; data taken from Inoue et al., 2000).

filtration experiment in which the solute concentrations were kept constant, so thatthe measured electrical conductivities could be converted directly to water contents
and related to the simultaneously measured pressure head values. The independent unsaturated hydraulic conductivity data were measured during steady-state experiments at different water flow rates. Correspondence between both optimized
and independently measured soil hydraulic and dispersion properties is very good
within the experimental range (Fig. 6.6–4 and 6.6–5, respectively). Further details
about the above-described experiment, including calibration of the four-electrodeprobe, are given in Inoue et al. (2000).
6.6.3.4 Field Experiment with Nonequilibrium Solute Transport
The fourth example uses data collected during an elaborate field experimentthat involved atmospheric boundary conditions, a layered soil profile, and nonequilibrium solute transport (Jacques et al., 2001). The soil hydraulic and solutetransport parameters for individual soil layers were obtained sequentially. This example considers nonequilibrium solute transport using the mobile–immobile dual
porosity concept.
Measured data were obtained from a 5.5-m-long and 1.2-m-deep transect located at Bekkevoort, Belgium (Jacques et al., 2001). Water contents and resident
concentrations were measured laterally each 50 cm at five depths (15, 35, 55, 75,
and 95 cm) using two-rod time domain reflectometry (TDR) probes (25 cm long,
0.5-cm rod diameter, 2.5-cm rod spacing). Matric pressure head values were measured using porous cups (6-mm diameter, 25 mm long) located at a horizontal distance of 10 cm from each TDR probe. TDRs and pressure transducers calibration
procedures are described in Jacques et al. (2000).
The soil profile was described as having four horizons: Ap (0–25 cm), C1(25–55 cm), and C2 (55–100 cm) overlying a Bt horizon. A thin layer of gravel was

MISCIBLE SOLUTE TRANSPORT

1445

Fig. 6.6–6. Hourly precipitation rates during the Bekkevoort field experiment.

spread over the soil surface to minimize evaporation, prevent erosive effects by raindrops, and decrease weed growth. The soil surface was exposed to the atmosphere(Fig. 6.6–6). The tracer (CaCl2) was applied on 28 Aug. 1998 in three pulses at aconcentration of 75 g L−1.
Field-averaged water contents, pressure heads, and resident concentrations
were calculated for each depth by averaging measurements at particular depths. Resident concentrations were normalized as follows:θ(zi,t)c r(zi,t)c r*(zi,t) = ____________t∫0 θ(zi,t)c r(zi,t)dt[6.6–11]

where θ(zi,t) is the average water content at depth zi and time t, c r(zi,t) is the corresponding average resident solute concentration, and c r*(zi,t) is the time-integral-normalized resident concentration (Vanderborght et al., 1996). Only four daily observations for each variable were used in the optimization. Data for 384 d were usedin numerical simulations of water flow (11 Mar. 1998–30 Mar. 1999). Observations
between Day 160 and 384, with solute applied on Day 168.5, were used for combined water flow and solute transport simulations. The layered soil profile was divided into four different horizons (0–15, 15–35, 35–55, and 55–150 cm, corresponding with the measurement depths), with layer-specific hydraulic and transport
properties.
The following procedure was followed in order to minimize problems of
uniqueness in the optimized parameters. The soil hydraulic parameters were optimized sequentially for each layer starting with the top layer. First, we assumed thatthe soil profile was homogeneous and optimized the soil hydraulic parameters (vanGenuchten, 1980) against depth-averaged water contents and pressure heads measured at the first depth. Next, we fixed the parameters of the first layer, defined a two-

1446

CHAPTER 6

Fig. 6.6–7. Example of optimized pressure head and water content time series at a depth of 75 cm for
the Bekkevoort field experiment (data taken from Jacques et al., 2001).

layer soil profile (0–15 and 15–150), and optimized parameters of the second layer
against data from the second depth. Then a three-layer soil profile was defined, andparameters of the third layer were optimized and so on. A typical example ofmeasured and optimized pressure head and water content time series is presented
in Fig. 6.6–7 for the 75-cm depth.
Once a satisfying correspondence was obtained between measured and calculated water content and pressure head profiles, normalized resident concentrations (Eq. [6.6–11]) were used to optimize parameters for the physical nonequilibrium transport model, that is, the parameters θim (the immobile water content),
λ (longitudinal dispersivity), and α′mt (mass transfer coefficient). A similar proce-MISCIBLE SOLUTE TRANSPORT

1447

Fig. 6.6–8. Measured and optimized normalized resident concentrations at depths of (a) 15 cm, (b) 35
cm, (c) 55 cm, (d) 75 cm, and (e) 95 cm for the Bekkevoort field experiment (data taken from Jacqueset al., 2001).

dure as for water flow was subsequently followed also for solute transport. Measured and optimized normalized resident concentrations at five different depths areshown in Fig. 6.6–8, while the optimized parameters are given in Table 6.6–2. Notice that the nonequilibrium transport model was able to adequately describe concentrations at all five depths, including the early solute arrival and long tailing at
deeper depths. Our calculations showed that a large part of the porous space did notcontribute to solute transport, and that only between 5 and 12% of the pore space

1448

CHAPTER 6

Table 6.6–2. Optimized solute transport parameters for the Bekkevoort field experiments.
Layer

θs

θim

λ (cm)

α′mt(d−1)

1
2
3
4

0.378
0.371
0.370
0.391

0.32
0.32
0.26
0.26

0.1
6.8
6.0
14.6

0.035
0.080
0.014
0.043

was mobile. Also, notice from Table 6.6–2 that the dispersivity increased with depth
and that some variations existed in the transfer coefficients of the various layers.
Jacques et al. (2000, 2001) gives additional details about this field experiment andthe invoked analyses.
6.6.4 Conclusions
As was also discussed in Sections 1.7 and 3.6.2, the parameter estimation
methods coupled with comprehensive numerical models can provide extremely useful tools for analyzing experimental data that may not be evaluated optimally using
conventional tools such as analytical solutions. Various numerical codes (e.g., HYDRUS-1D) have been developed for identifying, either simultaneously or sequentially, soil hydraulic and solute transport parameters from unsaturated flow and transport data in one-, two-, and quasi-three-dimensional porous media. We demonstrated
the utility of such tools using four examples of increasing complexity. Because of
their generality (in terms of the definition of the objective function, the possible combination of different boundary and initial conditions, and options for considering
multilayered systems), these tools are very attractive for analyzing conveniently abroad range of steady-state and transient laboratory and field flow and transport experiments.
6.6.5 References
Abbaspour, K.C., M.Th. van Genuchten, R. Schulin, and E. Schläppi. 1997. A sequential uncertainty
domain inverse procedure for estimating subsurface flow and transport parameters. Water Resour. Res. 33:1879–1892.
Brusseau, M.L., and P.S.C. Rao. 1989. Sorption nonideality during organic contaminant transport inporous media. Crit. Rev. Environ. Control 19:33–99.
Fortin, J., M. Flury, W.A. Jury, and T. Streck. 1997. Rate-limited sorption of simazine in saturated soilcolumns. J. Contam. Hydrol. 25:219–234.
Hopmans, J.W., and J. Šimçnek. 1999. Review of inverse estimation of soil hydraulic properties. p.
643–659. In M.Th. van Genuchten et al. (ed.) Characterization and measurement of the hydraulicproperties of unsaturated porous media. University of California, Riverside, CA.
Inoue, M., J. Šimçnek, S. Shiozawa, and J.W. Hopmans. 2000. Estimation of soil hydraulic and solutetransport parameters from transient infiltration experiments. Adv. Water Resour. 23:677–688.
Jacques, D. 2000. Analysis of water flow and solute transport at the field-scale. Ph.D. thesis, no. 454.
Faculteit Landbouwkundige en Toegepaste Biologische Wetenschappen, K.U. Leuven, Belgium.
Jacques, D., J. Šimçnek, A. Timmerman, and J. Feyen. 2001. Calibration of the Richards’ and convection–dispersion equations to field-scale water flow and solute transport under rainfall conditions.
J. Hydrol. 259:15–31.
Kool, J.B., J.C. Parker, and M.Th. van Genuchten. 1987. Parameter estimation for unsaturated flow andtransport models—A review. J. Hydrol. 91:255–293.
Kool, J.B., J.C. Parker, and L.W. Zelazny. 1989. On the estimation of cation exchange parameters fromcolumn displacement experiments. Soil Sci. Soc. Am. J. 53:1347–1355.

MISCIBLE SOLUTE TRANSPORT

1449

Medina, A., and J. Carrera. 1996. Coupled estimation of flow and solute transport parameters. Water
Resour. Res. 32:3063–3076.
Mishra, S., and J.C. Parker. 1989. Parameter estimation for coupled unsaturated flow and transport. Water
Resour. Res. 25:385–396.
Murali, V., and L.A.G. Aylmore. 1980. No-flow equilibration and adsorption dynamics during ionic transport in soils. Nature 283:467–469.
Nkedi-Kizza, P., J.W. Biggar, H.M. Selim, M.Th. van Genuchten, P.J. Wierenga, J.M. Davidson, and
D.R. Nielsen. 1984. On the equivalence of two conceptual models for describing ion exchange
during transport through an aggregated oxisol. Water Resour. Res. 20:1123–1130.
Selim, H.M., R. Schulin, and H. Flühler. 1987. Transport and ion exchange of calcium and magnesiumin an aggregated soil. Soil Sci. Soc. Am. J. 51:876–884.
Šimçnek, J., M. Šejna, and M.Th. van Genuchten. 1998. The HYDRUS-1D software package for simulating water flow and solute transport in two-dimensional variably saturated media. Version 2.0.
IGWMC-TPS-70. International Ground Water Modeling Center, Colorado School of Mines,
Golden, CO.
Šimçnek, J., and M.Th. van Genuchten. 1999. Using the HYDRUS-1D and HYDRUS-2D codes for estimating unsaturated soil hydraulic and solute transport parameters. p. 1523–1536. In M.Th. van
Genuchten et al. (ed.) Characterization and measurement of the hydraulic properties of unsaturated porous media. University of California, Riverside, CA.
Šimçnek, J., M.Th. van Genuchten, M. Šejna, N. Toride, and F.J. Leij. 1999. The STANMOD computer
software for evaluating solute transport in porous media using analytical solutions of convection–dispersion equation. Versions 1.0 and 2.0. IGWMC-TPS-71. International Ground Water
Modeling Center, Colorado School of Mines, Golden, CO.
Sun, N.-Z., and W.W.-G. Yeh. 1990. Coupled inverse problems in groundwater modeling. 1. Sensitivity analysis and parameter identification. Water Resour. Res. 26:2507–2525.
Toride, N., F.J. Leij, and M.Th. van Genuchten. 1995. The CXTFIT code for estimating transport parameters from laboratory or field tracer experiments. Version 2.0. Research Report No. 137. U.S.
Salinity Laboratory, USDA, ARS, Riverside, CA.
Vanderborght, J., M. Vanclooster, D. Mallants, J. Diels, and J. Feyen. 1996. Determining convective lognormal solute transport parameters from resident concentrations. Soil Sci. Soc. Am. J.
60:1306–1317.
van Genuchten, M.Th. 1980. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Soc. Am. J. 44:892–898.
van Genuchten, M.Th. 1981. Non-equilibrium transport parameters from miscible displacement experiments. Research Report No. 119. U. S. Salinity Laboratory, USDA, ARS, Riverside, CA.
van Genuchten, M.Th., and R.J. Wagenet. 1989. Two-site/two-region models for pesticide transport anddegradation: Theoretical development and analytical solutions. Soil Sci. Soc. Am. J.
53:1303–1310.
van Genuchten, M.Th., and P.J. Wierenga. 1976. Mass transfer studies in sorbing porous media. I. Analytical solutions. Soil Sci. Soc. Am. J. 40:473–481.
Weiss, R., and L. Smith. 1998. Parameter space methods in joint parameter estimation for groundwater flow models. Water Resour. Res. 34:647–661.