Journal of Hydrology, 105 (1989) 225-251
Elsevier Science Publishers B.V., Amsterdam -- Printed in The Netherlands

225

[2]
EXPERIMENTAL
AND THEORETICAL ANALYSIS OF SOLUTE
TRANSPORT FROM A DIFFUSE SOURCE OF POLLUTION

H.C. VAN OMMEN 1., M. TH. VAN GENUCHTEN2, W.H. VAN DER MOLEN 1, R. DIJKSMA !
and J. HULSHOF1
1Department of Land and Water Use, Agricultural University, Nieuwe Kanaal 11, 6709 PA
Wageningen (The Netherlands)
U.S. Salinity Laboratory, U.S.D.A./A.R.S., 4500 Glenwood Drive, Riverside, CA 92501 (U.S.A.)(Received June 16, 1988; accepted for publication June 26, 1988)

ABSTRACTVan Ommen, H.C., Van Genuchten, M. Th., Van der Molen, W.H., Dijksma, R. and Hulshof, J.,
1989. Experimental and theoretical analysis of solute transport from a diffuse source ofpollution. J. Hydrol., 105: 225--251.
Transport processes governing soil and groundwater contamination by a diffuse (non-point)source were studied by means of a field ~racer experiment in the Hupselse Beek catchment area in
The Netherlands. A bromide tracer was applied to the intake area (11 m x 260m) of a drain, andsubsequently monitored in the unsaturated zone, as well as in the drain efl]uent. Data from theunsaturated zone were analyzed using both the classical convection-dispersion solute transport
equation (CDE) and a regional stochastic model (RSM) which assumes lognormal distributions of
the pore water velocity and the dispersion coefficient across the field. Both models described thedata equally well, with the CDE model giving a somewhat better mass balance of the recovered
bromide. Calculated ~3olute fluxes leaving the unsaturated were used as input into a mixing cell
type model for transport in the saturated zone. Drain efliuent concentrations predicted with thecoupled unsaturated--saturated transport model were compared with the observed breakthroughcurve; results indicate the presence of preferential flow in the unsaturated zone.

INTRODUCTION
P u b l i c a w a r e n e s s a n d c o n c e r n s a b o u t t h e l e a c h i n g of n u t r i e n t s a n d
p e s t i c i d e s f r o m a g r i c u l t u r a l l a n d s h a s n e v e r b e e n h i g h e r t h a n it is t o d a y . T h i s
is n o t s u r p r i s i n g as m o r e a n d m o r e e v i d e n c e s h o w e d t h a t diffuse ( n o n - p o i n t )s o u r c e s of p o l l u t i o n p l a y a m a j o r r o l e i n t h e d e g r a d a t i o n of t h e q u a l i t y of o u rsoil a n d g r o u n d w a t e r s y s t e m s ( B r u y n , 1984; V a n d e r Veen, 1986; F u n k , 1987).
G r o u n d w a t e r p o l l u t i o n s t e m m i n g f r o m diffuse a g r i c u l t u r a l s o u r c e s t y p i c a l l y
i n v o l v e s first t r a n s p o r t t h r o u g h t h e u n s a t u r a t e d zone, f o l l o w e d by t r a n s p o r t in
t h e u n d e r l y i n g s a t u r a t e d g r o u n d w a t e r system. T h e t i m e scale for t r a n s p o r t in
t h e u n s a t u r a t e d z o n e is u s u a l l y m u c h l a r g e r t h a n t h a t for t r a n s p o r t i n t h e
u n s a t u r a t e d zone, t h u s f u r t h e r c o n t r i b u t i n g to t h e t i m e - c o n s u m i n g a n d c o s t l yn a t u r e of field t r a c e r e x p e r i m e n t s i n v o l v i n g c o u p l e d u n s a t u r a t e d - s a t u r a t e d* Present address: Grontmij Ltd., Engineers and Consultants, P.O. Box 203, 3730 AE, De Bilt (TheNetherlands).
0022-1694/89/$06.65© 1989 Elsevier Science Publishers B.V.

226

systems. Because of these considerations, field experimental studies of the
migration of solutes from diffuse sources are best carried out in shallow
aquifers with relatively short residence times. An aquifer of this type is locatedin the Hupselse Beek catchment area of The Netherlands. Soils in the region
consist of an approximately 2m thick surface layer of thin fluviatile andaeolian deposits,underlain by a thick Miocene clay layer of low permeability.
ity.
In this paper we present the results of a field transport experiment in theabove Hupselse Beek groundwater basin using bromide as a tracer. The
experiments were carried out in a I ha fieldwith recently installed subsurface
drains. A diffuse source of pollution was simulated by uniformly applying asmall pulse of KBr solution to a strip of land, 260 m long and having a width
twice the drain spacing (11 m). Main objective of the experiment was to study
and evaluate field-scalesolute transport processes operative i~ the unsaturated
zone, and in a coupled unsaturated-saturated soil/groundwater system.
Transport in the unsaturated zone will be analyzed by means of the classicalconvection-dispersion solute transport equation, and a regional stochastic
model accounting for field-scale heterogeneity in the transport process.
Transport in the saturated zone is described using a recently developed mixing
cell type model (Van Ommen, 1986).
THEORETICAL

Transport in the unsaturated zone

Numerous theoretical studies of solute transport in the vadose zone havebeen carried out over the last several decades. Initially, solute transport wasoften simulated by applying the classical Fickian-type convection-dispersion
equation (CDE) to field conditions. The CDE for nonreactive transport is givenby:~c~--~ =

~2c
~x2

D---

V ~c
~x

(1)

where t is time [T],x is soil depth [L],c is the solution concentration [ML-3],
D is the dispersion coefficient[L2T -'],and V is the average pore water velocity[LT-'].
The C D E has been shown to describe solute transport under laboratory
conditions quite satisfactorily(Danckwerts, 1953; Nielsen and Biggar, 1962;
Lindstrom et al.,1967). Its extension to field conditions affects primarily thevalue of the dispersion coefficient,which presumably relatesto the scale of theheterogeneities present in the flow field (Anderson, 1979; Dagan, 1986). This
flow field is often considered to be at steady state with respect to water flow,
and homogeneous with respect to water content. This was shown to be justified
by experimental and theoretical studies by Cassel et al. (1975) and Wierenga(1977). The assumption of steacLy-statewater flow associated with an average

227water content for the soil gave similar solute distributionsas those obtained bymodeling solute transport during transient wster flow.
The dispersion coefficientD in eqn. (1) is often expressed mathematically by:
D = Do + ~V [L2T -1](2)

where Do is the porous medium ionic or molecular diffusion coefficient [L2T - 1],and • is the dispersivity [L]. The applicability of this expression to transport inunsaturated soils has often been questioned (e.g., Beese and Wierenga, 1983).
Still there is reasonably strong evidence that eqn. (2) holds if the scalardispersivity is equal to a few centimeters (Sposito et al., 1986). In most cases theeffect of Do on the value of the dispersion coefficient is negligible as comparedto the effect of the dispersivity term in eqn. (2) (e.g., Biggar and Nieisen, 1976;
Schulin et al., 1987a). In practice, D is often used as an empirical parameter
representing the effects of all solute spreading mechanisms, such as transientwater flow, or various physical and chemical nonequilibrium conditions(Nielsen et al., 1986). In our case we will assume that the dispersion coefficient
is linearly related to the local pore water velocity, i.e.:
D = ~V(3)

There is now increasing evidence that solute transport in field soils cannot beadequately described by the CDE (Biggar and Nielsen, 1976; Van de Pol et al.,
1977; Bowman and Rice, 1986; Richter and Jury, 1986). This appears speciallytrue during surface ponding, when local variations in the pore water velocity
can become so large that the dispersion process can no longer be described witha physically realistic value of the dispersivity. Dagan and Bresler (1979)
introduced a method to determine variations in the pore water velocity interms of the statisticaldistribution of the hydraulic conductivity. From thesevariations in the velocity,first-and second-order spatial and temporal momentsof the concentration distributioncan be derived. An application of this theory
was given by Bresler and Dagan (1979) for constant input concentration. Theseauthors compared solute profilesbased on spatiallyvariable data with solute
profiles computed with the classical CDE. For both cases they derived
expressions for the length of the transition zone, defined as the difference
between the depths at which the concentration equalled 0.1 and 0.9 of the input.
In case of the CDE, this zone expands as the square root of time, whereas in the
Bresler-Dagan model the length of this zone increases linearly with time. The
profile calculated with the proposed model could also be obtained with the
CDE, provided the dispersivity used for predicting the concentration at I mdepth ranged from 16.1cm for low recharge cases to 334.0 cm for high recharge
rates. From this the authors concluded that pore-scale dispersion is not significant for transport in heterogeneous fields.
T w o years later Bresler and Dagan (1981)incorporated pore-scaledispersion
into their earlier model. This had very littleeffect on the calculated solute
distributions, as already had been foreseen in the 1979 paper. Finally, Bresler
and Dagan (1983b) developed expressions for the first- and second-order

228
moments of the concentration distribution as a function of depth and time,using a simplified solution (Dagan and Bresler, 1983; Bresler and Dagan, 1983a)of the nonsteady water flow equation. The transient water flow description
gave approximately the same results as the steady-state approach.
Jury (1982) introduced the concept of the transfer function model as a meansof describing solute transport in heterogeneous fields. In this approach, a soil
system is characterized uniquely by the way in which it transforms an input
function into an output function, and assuming a log-normal distribution of
travel times and no local dispersion. Jury et al. (1986) later gave a more generalderivation of transfer function models for transport of solutes in soils. Solutes
in the model may be conservative, or may undergo physico.chemical or
biological transformations in a soil with spatially or temporally variable flow
properties. The authors conclude that any mechanistic solute transport :~,aodel
which maintains mass balance can be made consistent with the transfer
function concept.
A slightly different approach was followed by Amoozegar-Fard et al. (1982)
who used Monte Carlo simulation techniques to obtain various solute concentration distributions as affected by variabilities in the pore water velocity
and the dispersion coefficient. This approach assumes that the field can be
represented by a large number of independent noninteracting soil columns inwhich vertical transport is described by the CDE. Each column has its own
pore water velocity V, which is considered to be lognormally distributed among
the individual columns. The lognormal distribution, p(V), with mean ~ andstandard deviation a is given by:
1p(V) = Va~/(2n) exp[-(ln(V)-~)
2a 2

2]
(4)

The first moment of this distribution is the field average (global) velocity, ( V )given by:(V)

= exp(l~ + a ~/2)

(5)

The probability density function of the dispersion coefficient D is considered tobe completely defined in terms of p(V) and deterministic ~.
The spreading of a solute pulse in a "field-averaged" soil profile can thus be
attributed to two mechanisms: (1) local dispersion in the columns themselves;and (2) differences in the pore-water velocities among the columns.
Results of the Monte Carlo simulations by Amoozegar-Fard et al. (1982)
show that the variability of Vhas a dominant effect on the field-averaged solute
profile. In the limiting case of zero dispersivity, the approach is analogous to
the transfer function model of Jury (!982).
In our study we will consider two methods of describing solute transport in
the unsaturated soil profile: the classical CDE, and the method introduced byAmoozegar-Fard et al. (1982) commonly referred to as regional-stochastic
model (RSM) (Parker and Van Genuchten, 1984a; Schulin et al., 1987b). To
partially account for the transient flow regime in the system, we will use the229amount of drainage water, rather than time, as independent variable. Thisapproach is somewhat similar to the linearization method used by Parker and
Van Genuchten (1984).
T r a n s p o r t in t h e s a t u r a t e d z o n eOnce the bromide tracer leaves the unsaturated zone to reach the groundwater table, transport occurs in the saturated zone along flow lines in the
direction of the drains. Transport in this part of the system is simulated byconsidering the saturated zone to be a perfectly mixed reservoir. Previous
examples of such an analogy are given by Eldor and Dagan (1972), Ernst (1973),
Gelhar and Wilson (1974), Raats (1983) and Van Ommen (1986), among others.
The approach assumes that: (1) groundwater flow is at steady state in whichcase the drain discharge from the system becomes equal to the recharge rate
multiplied with the total recharge area; (2) ~he thickness of the aquifer and its
porosity are constant; and (3) groundwater flow is horizontal, thus making the
Dupuit-Forchheimer assumptions valid. The concentration in the drainage
water for a unit step input of a diffuse pollution source is then given by:
C(t)

=

1 -

(6)

e -N¢~H~

where C is the concentration of the drainage water, t is the time, N is the steadygroundwater recharge rate [LT -1 ] and e [L3L -3] is the porosity of the aquiferof thickness H [L]. As before, the total amount of outflow, Nt, can be used asa substitute variable for time if the flow regime is transient.
Whereas flow in the direct vicinity of the drain is essentially radial the effectof this on the solute breakthrough curve is only of importance for values of
H I 5 > 0.1, where L is the drain distance (Ernst, 1973). Groenend~jk (1986)~h~wed that the concentration of the drainage water can be corrected for radial
fl~w within a distance H from the drain by using a modified travel time
distribution as first derived by Ernst (1973). Groenendijk's analysis leads to:
c(t)

=

1 -

(7)

A

where:A=

1 -2 x/(2Nt)/(~eL) -for t ___ (~H2)/(2NL)

and:
A = (1 - 2 H ] L ) exp [ - ( N t ) / ( e H ) + ( ~ H ) / ( 2 L ) ] ,for t > ( ~ e H 2 ) / ( 2 N L )Although H / L is generally smaller than 0.1 in our study, we still included theabove correction for radial flow near the drain into our solution.
Drain water concentrations as a function of the cumulative amount ofdrainage is now calculated by the convolution integral of the impulse responseof the saturated transport system with the solute loading rate leaving theunsaturated zone. We refer to the study of Van Ommen (1985) for a moredetailed description of this technique.

230

MATERIALS A N D M E T H O D S

Geographical setting
The experimental field was located in the Hupsel catchment area in the
Eastern part of The Netherlands. The field met the requirements with respect
to the thickness of the saturated zone (see Introduction). Drains were first
installed in the area in 1982. Soil type is a loamy sand Haplaquod (WSsten etal., 1985) with a top soil having a finer texture than the deeper soil layers(Dekker, pers. commun., 1986). Beneath these layers at a depth ranging from 1.5
to 2.5 m, a 20-30 m thick impervious clay layer of Miocene age was found. This
clay layer forms a natural boundary for the unconfined groundwater flow
system (Studygroup Hupselse Beek, 1974). One half of the drained field wasunder permanent grass and the remainder under corn during the growing
season.

Tracking the drains
The bromide tracer was applied to the catchment area of one drain after thisarea had been determined in the field.Drains in the fieldwere located by using
the T R A C K A system (manufacturer Hormann) which consisted of a radio wave
emitter attached to a special jethead normally used for cleaning drainage
tubes. At places where the emitted signals were most clearly received, a marker
was put in the soil. Having tracked a whole drainage tube in this way, a cord
was stretched along the markers. Three drains were jetted and tracked in thisway: the drain to be sampled and the two adjacent ones. The average drain
distance between the center drain and its two neighbors was taken to be half
the distance between the outer drains. This resulted in an average drain
spacing of about 10.85m. Because the drains were approximately 260m long,the catchment area comprised about 2850 m s. The markers of the outer drains
also served as boundaries for the trace~-application. After application of the
tracer the center drain was dug up at several places to determine its depth in
the field.Figure la gives a schematic plane view of the experimental field;avertical cross-section is shown in Fig. lb.

Application of the tracer
The tracer used in our study was the negatively charged bromide ion, often
recommended as a soil water tracer (Jury et al.,1982; Bowman, 1984). Application took place on December 12, 1985 in the form of a 4 . 5 m m K B r solution(2.55g I-I) using a sprinkler irrigation system (Fig. 2) with a high application
uniformity (the coefficient of variation was less than 2 % when not driving,
Institute for Soil Fertility,pers. commun., 1985). T w o strips of 5 m width eachon both sides of the drain were sprayed four times, resulting in a total application rate of 60kgBr- to a surface area of 2 0 m × 260m. Thus, we treated an

231maize field(a)grassland

~L,.

drain

E]

cGtchmentof center
drain
booth

~i

e

m

m

drain

t~

In

row of piezometersand microfi|ters

2pplication stripof Dec.12,1985
IO m

Q

row of piezometersand microfilters

0

50 m

==so~lsampling plot~/~baseflow

maize
130 m

(b)

mm
T 1400-

drain

/

~l

grassland
130 m

I

1000
i

•

I

:
~

600
drain

200
BOOTH,
200,""l[Irrnnrnrnrr/B;nf/rr/Tr//"-50 cm600

NfF1~'
1I'I'I'F
fI 0

1000

1400

~.TrrlrlTzrf?
TrrlTf/77TflTf/77r/77nz
n'

50m I0_..,nffll fllf'

1800
Fig. 1. (a) Schematic plane view of the experimental field. The small squares indicate the sampling
spots; Co) vertical cross-section of the experimental field.

232

Fig. 2. Application of the tracer solution.

area nearly twice as large as the catchment area of the center drain. This wasdone to account for possible irregularities i,l the water divide between thedrains. No ponding or surface runoff was observed during the application.

Soil sampling
To quantify solute transport processes in the unsaturated zone, soil samples
after 69.5, 167.0, 225.0 and 304.5 mm of drain discharge were collected. At eachsampling time, the soil was sampled at twenty plots regularly distributed along
the drain: ten plots in the corn field and an equal number in the grassland. Eachplot wa~ located in the center of the first "application strip" next to the central
drain. The first plot was located in the corn field at the right-hand side of thecentral drain, the second plot at the left-hand side, etc. (Fig. la). Sampling
within a plot consisted of taking soil samples every 5 cm with a 6 cm diameter
auger up to a depth of 70 cm. Deeper sampling would sometimes lead to samplesof the saturated zone (see Fig. lb). Distances within each sampling plot neverexceeded 30cm. This procedure thus resulted in 20 × 14 = 280 soil samples at
each sampling time. Analysis of these samples in the laboratory involved the
preparation of 1:2 soil-water extracts, followed by bromide concentration determination.

Piezometers and sampling of the saturated zone
Rows of piezometers were installed perpendicular to the drain in the middleof the corn and grass fields (Fig. la). Bromide concentrations of the ground-233

water at these same locations were monitored by means of microfilters (3/4"
inner diameter and 10-cm length) which had been installed at different depths
in the saturated zone. Groundwater samples were withdrawn after 109.5, 193.5,
278.0 and 353.0mm of drain discharge. Piezometer readings were taken at
irregular intervals.

Drain discharge measurements and sampling
Drain water was collected using a tank which was partially sunken into thebottom of the ditch at the drain outlet. Water was pumped automatically fromthis tank in increments of 0.089 m 3 whenever the water level reached a floating
contact. Once every 16th time of pumping, an automatic water sampler(Manning $4400) received a signal to operate. In this way (16 x 0.089m3/
2850m 2 soil surface area) a water sample was taken after every 0.5 mm of draindischarge. Samples were thus collected on a discharge-proportional basis. The
entire system was controlled with a Hewlett Packard 75B programmablecalculator. The equipment was installed in a booth with an electrical heating
system to protect it against severe colds during the winter season.

Bromide analysis
Soil water extracts, drainage water samples and groundwater samples wereanalyzed for bromide using a high performance ion-chromatographic methodwith UV-detection at a wavelength of 193 nm (Neele, 1987).

Bromide uptake by plants
In order to account for possible bromide uptake during the winter season,
grass samples were taken at the end oi" March, 1986. Grass was sampled at six
plots, each randomly located in the application strip. The corn field remained
uncultivated until the end of May, 1986.
RESULTS AND DISCUSSION

Local water balance
During the summer and early autumn preceding the experiment, the drain
did not discharge and the ditches were dry. The bromide tracer was applied afew weeks after the drain started to flow. Precipitation data were obtained fromthe meteorological station Assink, located some 900 m from the experimental
field (Warmerdam, 1981). Figures 3 and 4 show the precipitation and ~hecumulative amounts of precipitation and drain discharge from day I (December
12, 1985) to day 191 (June 20, 1986) at the end of the percolation period whendrainflow again ceased. The very dry period from day 47 to day 84 correspondedwith a long period of frost (February 1986 was the dryest month of February on

234

2,I

precipifofion during the experiment

_ 2/+
E
E 20
t-

._o 164--•~. 12-

8

I

II
r!l! !,,,r,I,,,I',t, ,I
,.,,,

I

20

Ill

40

. . .60. . . .,...........
. .80.!!.................
. 100II, ,!f...........................
120

140

160

180

day number
Fig. 3. Precipitation during the experiment.
...- 500gf"

cumutnfive rain
cumutot'ive drain discharge

b 4oo
e-U

/
.J
i.e~ .... •°./
f
t

..,

e-

_i, - ~ ~ ~

300

\
(,.,

~/,°*
........

200

','t~

/

E

~ loo

tI #J'/
20

40

60

80

100
120
day number

140

160

180

Fig. 4. Cumulstlve precipitation and drain discharge during the experiment.

record in The Netherlands). The differences between cumulative 8mounts ofdischarge and precipitation increased during the experiment because ofincreased evapotranspiration in the spring.
The experiment was continued until the summer of 1987. At the beginningof the experiment we assumed that drain discharge would result exclusively
from water infiltrating the application strip between the water divides at eachside of the center drain. Unfortunately, piezometer readings in the corn field
later indicated that a regional baseflow was superimposed upon the local flows
between the drains, at least during part of the experL~nental period (Fig. la). A
simple hydrological calculation showed that this base~ow had a flux of at mosta few cm per day, at least in the corn field. No head gradient was observed
perpendicular to the drain in the grass area. Because the tracer was applied to235a strip that extended 5m beyond the water divide between the drains,
horizontal inflow of bromide by water from the catchment of the central drain
would only be guaranteed for some 250 days. Therefore, only results for the first
percolation season ending in June, 1986, will be considered here. Drain water
samples in autumn 1986 showed bromide cc~ncentrations of 1.4 mg 1-1 for the
upstream drain, 14.0mg1-1 for the centre drain, and ll.8mg1-1 for thedownstream drain, thus confirming the presence of the baseflow. Note that the
average of the upstream and downstream drainage water concentrations isabout half the concentration in the center drain, as would have been expected.

Soil samples and unsaturated transport
Soil bromide concentrations were grouped according to land use: grass orcorn. They will be treated separately, even though about half of the concentration distributions from the grass field were statistically not significantly
different from those of the corn land, as predicated by the Kolmogorov-Smirnov statistic (based on the maximum difference between two empiricalcumulative distribution functions), and the Wilcoxon test (based on the rankof the data in both groups). Because the grass and corn fields were treatedseparately, the number of soil sample replications at each time and depth wasten. Distributions of the concentration~ were tested for normality and lognormality by the tests of Kolmogorov-Smirnov, Cramer-Von Mises (Rao et al.,1979) and Shapiro-Wilk (Shapiro and Wilk, 1965). No clear distinction could bemade from these analyses: sometimes the distributions were better described bya normal probability density function (PDF), sometimes by a log-normal PDF,
and sometimes neither of the two couL'! adequately describe the concentration
data. Hence, spatially-averaged (global) bromide concentrations at given
depths and times were determined by simply taking arithmetic mean of theconcentrations. Confidence intervals (90%) were calculated using the bootstrapping procedures (Dane et al., 1986) which do not require a-priori
assumptions regarding the probability distribution. The number of bootstrap
replications was chosen at 5000. Concentration profiles for the four samplingtimes are presented in Figs. 5 (corn) and 6 (grass), together with the optimized
curves for the CDE and RSM transport models.
The nonlinear least-squares optimization program developed by Parker and
Van Genuchten (1984a) was used for the analysis of the field data. This program(called CXTFIT,version 1985) can be used to estimate the unknown parameters
in several one-dimensional solute transport models from spatial distributionsof the resident concentration (Parker and Van Genuchten, 1984b). In our study
where the flow rate varied greatly in time, the controlling independent variable
was taken to be the amount of drain discharge, expressed in mm. By doing so,
the pore water velocity MVis given ~vith respect to the flux density of the draindischarge. Thus, a computed velocity M V of, say, 4.0 means that the pore watervelocity is four times larger than the water flux density. Its reciprocal(MV-1 _ 0.25), corresponds then with the mobile water content of the system.

236
concentration (rag t-1)
6O 80 120 160 200

conrc,,rrotion (rag I -I)
I.O 80 120 160 200

01 ~%--c-7

o.'1

'

0.2

,

, <

E
E

0-51 ~
• _I

~

0.6.

COE_modet.... R S

sompIlng 2
CDE. mode[o.61 i-~_ mode[02
c0ncentroh0n

0

~0

80

120

(rag t-1)
160 200

0

40
•

•

'

.modet.... Rs][Z]O'7

concentration (rag [-1)
80 120 160 200
•

,

•

'

•

•

r

,

-

-

•

~

~

,

,

i

0.1
0.2

0.2

g a3

~,0.3.
..,= •
~Qz,.
1;3

(15
i ~
Q61 ~ .
0,TJ

~COE.modeL.... RS .modet

"II

Q6

sampling z,
C D E .modet.... RS . modet

0,7

Fig. 5. Concentration profiles in the corn field at different sampling times.

CDE results
Parameter estimation results for the convection-dinpersion equation (CDE)
are given in Table ]. The sampling number refers to the accumulated draindischarge when the soil was sampled (i.e.,sampling 1 corresponds with adischarge of 69.5turn, 2 with 167.0mm, 3 with 225.0mm, and 4 with 304.5mm).
Because the application of the tracer solution was uniform, the applied tracer
pulse length towas 4.5 mm. This value was not corrected for bromide uptake bygrass as only an estimated 3 % of the amount of bromide applied was recoveredfrom the grass.
Initially,we optimized the pore water velocity M V and the dispersion coefficient D, of the C D E model. The mass balance was verified by carrying out asecond optimization with unknown tobut with M V and D fixed at their values
from the firstrun. The optimized tovalue was divided by the k n o w n tovalue toyield the recovered fractions of the tracer as viewed by the C D E model. These237
concentration {rag 1-1}

i0.,. ~, ~, ~toncentrotion {rng t 1}

I00 . 120

0

20

~

60

80

100

120

0.2

1

~03.

,i

~0.~.
I '/ I
i / I

05.
(16.

I'/I

snmp[ing 1- - [ D E_ model.... RS _ model

0.7[ DE _ mode[..--~_. RS _model
0.7.

concentration (mg t "1 }

20

~

i

0.1,
0.2
_ 0.3

0.4,
0.5
(16.
O]

i

60
,

80
,

|

100
•

'

120
'

0

, 20

concentration ling t -I }, ~, ~, , SO ,IO.O

.I~0

|

0.1
-I

02.

03.
I

1~. O.t+,_0,5sampling 3[DE. model...... RS _ modet

0,6.

sampling/-,
ODE.model.....RS. modet

--

0.'~

Fig. 6. Concentration profiles in grassfield at the different sampling times.

mass balance ratios are also indicated in Table 1. Values range between 0.911
and 1.008 times the amount applied. Optimizations were carried out for all four
samplings, and for various combinations of these as shown in Table 1. Optimizations were also performed using data from all plots separately, instead ofusing areally-averaged values. Results for the individual plots were generally
very similar to the global estimates except for substantially lower correlation
coefficients between the fitted and observed concentrations because of fewer
sampling replicates and the inclusion of spatial variability.
Because the relative pore water velocity in this study is the inverse of thevolumetric water content, one should be able to estimate the average mobile
volumetric water contents of the corn and grass plots from the optimized
parameter values for MV in Table 1. Results compare well with the measured
gravimetric water content distributions for the grass and corn fields shown in
Figs. 7 and 8, respectively. These globally averaged distributions, with their

238
TABLE 1
Parameter estimation results for the convection dispersion equation (CDE)
Sampling

Mass balance

Number ofobservations

MV

(-)

(ram)

14
14
14
14
56
42
42

3.27
5.40
3.77
3.68
3.19
3.33
4.31

37
196
95
95
62
42
178

0.981
0.988
0.992
1.00
0.911
0.957
0.953

3.04
7.16
3.28
3.78
3.03
3.09
4.15

145
322
211
173
198
154
353

0.966
0.984
1.01
0.997
0.931
0.951
0.934

L a n d use: corn

1
2
3
4
I-4
1,3,4
2,3,4

Land use: p e r m a n e n t grass

1
2
3
4
1-4
1,3,4
2,3,4

14
14
14
14
56
42
42

90% confidence intervals, were based on measured gravimet:~ic water contents
at the four sampling times. Dry soil bulk densities derived from open pit
measurements ranged from 1.3 gcm -3 for the topsoil to about 1.7 for the deeper
layers. At a few locations we found some thin peat layers with dry bulk
densities of about 0.gg cm -'~. Unfortunately we were unable to measure thevolumetric moisture content and the tracer concentration of the same soil
sample. Confidence intervals of the water contents in Figs. 8 and 9 were basedon the normal distribution of the gravimetric water contenlt. Validity of thisassumption was verified with the same three statistical tests used previously.
The inverses of the M V values in Table 1 correspond well with the measured
water content profiles provided the b u ~ density factor of 1..3-1.7 is taken into
account. An exception to this is the second sampling when the optimization
procedure resulted in an unreasonably high pore water velocity. Note that porewater velocities for the corn field seem to be somewhat hi~;her than those for~he grass because of lower water contents.
The dispersivitie~ ~ show an analogous behavior: for corn they are in
reasonable agreement with values reported in the literature, except for thesecond sampling when a relatively large value resulted. The dispersivitiesfor
the grassfield were always larger than those for the corn fieldwith the second
sampling again showing a much larger value. The larger values of ~ for grass
may be attributed to more prominent preferential flow paths in the soilbecause

239

gro~ water content (g.woter/g.soi[)0.10 020 0.30 0.~ 0.50 0.60

0

,,,,

0.1

'

'

,~..__:1
I I,

' '1

'

O2

0.2

~0.3

Q3

sampling I

0.5.

gray. vwoter content (g.wnter/g.soi[)0.10 0.20 030 OH 0,50 0.60

05 t

L~

sampting 2

0.6,
0.7.

0.7

gray. water content (g. water/g.soit )00.10 0.20

0.30 0.~0 050

Of&) 0

0.1.

01

0.2

0.2
I

"~ 0.3

g a3

II

¢,.

"~ o.~..
0.5,

grav. water content (g.~ter/g.soit )o~, 0.20 0.3o, o.~ o,so 0.~

I

~ or,
sampling 3

0.5

0,6,

0.6

0.'/,

o7

'

I
I

I
'

~

I
'

I
J

sompting~,

Fig. 7. Gravimetric moisture content profiles in the corn field at different sampling times (90%confidence intervals are all indicated).

of continuous plant root activity in the winter, which in turn could haveresulted in a larger variation of pore flow velocities.
Taking all four samplings together, (1-4 in Table 1) yields values for MV
which are in good agreement with the measured water content profiles, both forthe corn field and the grassland. Note, however, that the first sampling has adisproportionally large influence on the results. This is because sampling ! hasby far the largest contribution to the sum of squares than the other samplings
with much lower concentrations.
The parameter optimizations above could or also have been executed withall three transport parameters as unknown (i.e., MV, D and the amount ofapplied solute to). This would have given a more direct check on the massbalance. Such a procedure is only successful if the solute concentration profiles
exhibit clear maxima. Otherwise, a nonunique set of transport parametersgrov. wofer content (g.woter/g.soit )
0.10, (120, 0.30, 0.~ 0.50 0.60

gray woter content (g.woter/g, soil)o.lo 0.20 0.30 o.~,o o.s,o 0.60
|

|

i

0.1.

sompling I

~,

0.2

[,

.=
~0.4

I
I

0.5.

I
I

'
'

'

I

Q6.

[
I

Lr ,-q~ 0.,.,.

t

,

I

'

I

~ 0.3

Q2.
e 0.3,

sompting 2

_I

0.6.

_

"

't

0.5.

I

D I

I

II

'
'

'

I,

'

I

I

0.7,

grov. woter content (g.woter/g. soil )0.I0 0.20 030 0.~ 0.50 (160

t,I

i

i

0.1.

i

~~ I

0.2.
-E

03

z:

~. 0.~,

i

i

I
L

'

I

A
U

grov. woter con~enf (g.woter / g.soi[ )0.10 0.20 0.30 0.~ 0.50 0.60
Im

sampling3 0,1
0.2

I

-

I,

,

1

I

I

I

ttl

(16

0.7.

~]

0.7

..I.

I

I

I

[l

0.6.

|

sompLing/+

g a3.
=~ 0.4,

0.5

I

,

I

I
J

Fig. 8. G r a v i m e t r i c m o i s t u r e c o n t e n t profiles in g r a s s f i e l d at t h e d i f f e r e n t s a m p l i n g times.

soil surfoce
unsoturofed

zone

sa,~urofedzone

groundwoter
toble

Fig. 9. Vertical cross-sectionof the coupled transport systcm.
(notably for to) can describe the concentration profile. For example during the
optimizations M V and t0 will be positively correlated in case of increasing
concentrations with depth, and negatively correlated for decreasing concentrations. Hence, the three-parameter optimization can only be applied to241
TABLE 2Mass balance results as based on data from the first sampling time

Convection dispersion equation
M V (-)• (mm)

to (mm)Mass balance

36

4.398

0.977

101

3.812

0.847

• (ram)aln{MV}

Mass balance

92

0.210

0.812

L a n d use: corn

331

~ a n d use: permanent g r a s s l a n d741?egional stochastic model:MV~ (-)
L a n d use: grassland?.61

the first sampling time when a distinct maximum is present. Results for thiscase are given in Table 2; they show relatively good mass recoveries. Velocities
and dispersivities agree quite well with those found earlier, especially for thecorn field.
R S M results
Parameter estimation results for the regional stochastic model (XSM) arepresented in Table 3. Optimizations in this case were carried out for threetransport parameters P~ncMV),aln{MV)and ~, instead of two as was done for the
CDE model. A mass balance check was also made, similarly as for the CDE.
Some results have been omitted since gin{My)at times converged towards zero
during the optimization procedure, thus giving the same values of the transportparameters as the CDE optimization. Most values of the average pore watervelocity- ( M V ) are somewhat higher than those calculated with the convection-dispersion approach. The velocities calculated for the second sampling
are again remarkably higher than the others. The mass balance checks
generally showed a substantial decrease in the sum of squares when theamount of applied bromide (to) was also optimized. This amount was always lessthan the actually applied amount (4.5ram x 2.55g1-1 = l l . 4 8 g B r - m-2), thus
resulting in mass balances of less than 100%.
A direct check on the mass balance was only carried out for the firstsampling for the grass, again resulting in the relatively low recovery of 81.2%(Table 2). No explanation for this low recovery could be found.

242
TABLE 3
Parameter estimation results for the regional stochastic model
Sampiing
numberNumber ofobservations~MV) (-)

14
14
14
14
56
42
42

.
8.26
4.71
.
4.~
3.53
6.42

~ (mm)#ln~Mv~

Massbalance

0.651
0.457

0.65
0.801

0.411
0.186
0.607

0.748
0.917
0.665

9.303
1.010

0.779
0.352

0.647
0.521
0.300
0.840

0.628
0.651
0.825
0.463

L a n d use: corn

1
2
3
4
1-4
1,3,4
2,3,4

.

.

.
64
33
.
15
33
46

.

.

Land use: grassland
1
2
3
4
1-4
1,3,4
2,3,4

14
14
14
14
56
42
42

3.85
20.03
.
5.97
4.71
3.70
8.85

.

85
0.2
.
70
61
103
33

.

Statistical properties of the transport parameters
The transport parameters were also determined for every vertical concentration profile separately, using a two.parameter fit for M V and D for the C D E(Table 4).The resulting distributions of velocitiesand dispersivitiesof the C D E
optimization process were again tested for normality and log-normality by the
Kolmogorov-Smirnov, Cramer-Von Mises and Shapiro-Wilk tests. In somecases we could not discriminate between these distributions as often noither
the normal, the lognormal, or both, were found to be acceptable. This isindicated with a question mark in Table 4. The tab!e also shows correlations
between M V and D, and between ln(MV) and ln(D). The values for (MV)assume a log-normal distribution of M V as based on eqn. (5). Correlations areoften close to one, with some exceptions.
Linear relations between V and D assuming zero intercepts are presented inTable 5. If the intercept was not set equal to zero, the regression line would
have crossed the D-axis often at values below the origin, which is physicallynot correct. Table 5 shows that very high values of the dispersivity will occur
when the power of the velocity in eqn. (5) is set equal to one, as is often assumed
in the literature.
The calculated values for ~MV) and ~ln~Mv~(Table 4) at samplings 2 and 3(corn) and 2 and 4 (grass) corresponded rather well with the values as
determined from the field-averaged profiles (Table 3) and the RSM parameter

243
TABLE 4
Statistical properties of the fitted transport parameters
Sampling

Number ofobservations

Distribution of M V

Distribution cf D

I0
I0
10
10

normallog-normal
log-normal?

log-normal?
log-normal
log-normal

40

9

9

Correlation
between M V - D

Correlation
between ln(MV)-ln(D)-0.194
0.993
0.933
0.997
0.915

0.045
0.844
0.900
0.969
0.861

( M V ) (-)

O'in(MV)

3.04
7.55
4.51
8.67
5.49

0.358
0.906
0.506
1.096
0.778

Land use: corn

I
2
3
4
1-4

1
2
3
4
1-4

1
2
3
4
1-4

10
10
10
10
40

10
10
10
10
40

Land use: permanent grass
Distribution of M V

Distribution of D

1

I0

?

2

10

log.normal

3
4

10
10

log.normal
?

1-4

40

?

log.normal
log-normal
log-normal
log-normal
log.normal

Correlation
between M V - D

Correlation
between ln(MV)-ln(D)

1
2
3
4
1-4

1
2
3
4
1-4

10
10
10
10
40

10
10
10
10
40

0.111
0.806
0.921
0.925
0.831

0.490
0.957
0.720
0.865
0.024

( M ~ (-)

~71n{MV )

3.06
14.89
4.38
5.16
5.90

0.467
1.202
0.614
0.560
0.819

244
TABLE 5R e l a t i o n s b e t w e e n t h e fitted d i s p e r s i o n coefficient a n d t h e l o c a l r e l a t i v e p o r e w a t e r v e l o c i t y
Sampling

N u m b e r ofobservations

D

--

• MV

D

=

• (MV) ~

10
10
10
10
40

D
D
D
D
D

=
=
=
=
=

46.9MV758 M V309 M V
3770 M V
2870 M V

D
D
D
D
D

=
=
=
=
=

70.1 M V °'14
10.2 M V T M
3 . 3 7 M V 2"~
7 . 4 6 M V 249
8 . 3 6 M V 223

D
D
D
D
D

=
=
=
-=

182 M V
1520 M V
452 M V
~27 M V
1380MV

D
D
D
D
D

=
=
=
=
=

82
M V L23
18.9 M V z15
22.3 M V 2"1°
16.1 M V 2"~523
M V °'°e

L a n d use: corn

1
2
3
4
1-4

L a n d use: p e r m a n e n t g r a s s l a n d

1
2
3
4
1-4

10
10
10
10
4O

estimation process. Note that in the latter case, the effect of dispersion has
been omitted, so that the results are not quite comparable. However, if we
assume that the effect of dispersion is negligible with respect to the variability
in the velocity as characterized by a~n~MV),then the comparison is justified.
We concluded that both the convection-dispersion equation and the regional.stochastic model describe the field data reasonably well. Because of thesurprisingly large local variations across a seemingly quite homogeneous andsmall area, a definite preference for one of the two models cannot be establishedfrom the available data. On the other hand, the CDE, performs somewhat better
than the RSM in terms of the mass balance accuracy. Therefore, in thefollowing we will only use the CDE description of the processes in the unsaturated zone.

Samples of the saturated zone
When the field tracer experiment was first initiated we planned to construct
two-dimensional contour maps of the measured groundwater bromide concentrations perpendicular to the drain, and to hopefully use the plots fordiscerning information about the transport of bromide along streamlines
towards the drains. Unfortunately, groundwater samples withdz.awn from themicrofilters installed at different depths in the aquifer yielded such an erratic
concentration pattern that any attempt to analyze them in terms of the above
objective was judged to be useless. W e were also unable to clearly deduce thepresence of a base flow in the saturated zone from the measured solute concentrations, especially under the corn field.Still,as expected, concentrations
at greater depths tended to increase somewhat with time during the245experiment. The large variability of the measured groundwater concentrations
may be attributed to a large variability in the permeability of the aquifer at thescale of the sampling program, a feature which was recently described by
Sudicky (1986), and/or to irregular bromide solute fluxes across the catchmentfrom the unsaturated zone into the groundwater system.

Drainage water sampling and coupled breakthrough curves
Figure 10 shows the measured concentration distribution of the drainwater.
The plot reflects the behavior of the coupled unsaturated-saturated transportsystem upon an instantaneous input of bromide across the catchment area.
Careful analysis of this showed that about 53% of the total amount of applied
bromide was recovered during the first 6 months of the experiment.
Having determined the parameters of the two transport models for solutetransport in the unsaturated zone, and assuming that transport in thesaturated zone obeys eqn. (7), it should be possible to predict the bromide
concentration of the drainage water as a function of the drain discharge. The
calculation procedure is schematically shown in Fig. 9, which gives a verticalcross section of the coupled unsaturated-saturated transport system perpendicular to the drain. The input concentration as a function of the amount (mm)of drain discharge for the saturated transport system is obtained by calculating
the breakthrough curve from the unsaturated zone at the groundwater table.
In this case, flux concentrations as defined by Parker and Van Genuchten

o.

I"M

J
0

Z
0
,.-,,
I.---

o
lid

PREDICTED/MEASURED BTCs
MFASURED BTC
oPREDICTED FROM T2
0PREDICTED FROM T4~ P R E D I C T E D FROM T 2 , T 3 , T 4+PREDICTED FROM T ; - 4

n," o
i.i.

Z

.

~

W~
0
Z

od

o

C.3

•

60.0

90.0

120.0

150.0

DRAIN

210.0

240.0

DISCHARGE

180.0

(MM)

270.0

300.0

330

0

300.0

Fig. 1O. Measured and predicted drain effluent BTC's using the fitted unsaturated solute transportparameters.

246(1984b) were used, and not the resident concentrations employed in thep~rameter estimation process of the unsaturated transport models. The concentration of the drainage water was obtained by numerically evaluating theconvolution integral of the time-dependent input concentration with theproposed impulse response for the saturated transport system. Because thethicknesses of the unsaturated and saturated zones were not constant along
the drain, the concentration of the drain water was calculated as the averagesolute concentration of twenty equally-spaced locations, each having its own
coupled transport system. The breakthrough curve for a certain location along
the drain then depends on land use, the thickness of the unsaturated zone, thethickness of the saturated zone, and its water content. The thickness of theunsaturated zone at a given location could be equate~, to the drain depth,
except that piezometric observations during winter indicated that the unsaturated zone was about 10 to 20cm less. This correction i~creased thethickness of the saturated zone as well.
The liquid-filled porosity of the aquifer was determined from soil samples
taken from the deepest layers (about 110 cm below the soil surface) exposed inten profile pits; from measured dry bulk densities, the mean porosity wasestimated to be 0.33; this value will be used in all subsequent calculations.
Figure 10 compares the measured breakthrough curve (I~TC) with four
predicted curves which were chosen such that the full range of transport
velocities as determined by CXTFXTwas considered, including the highest
velocities for sampling 2, and the lowest velocity resulting from the simul.
taneous fit of all four samr~ing times. The simultaneous fit considerably underestimated the drain wa~er concentrations, especially at early times, e v e nthough the M V values were in good agreement with the measured soil watercontents, while at the same time the dispersivity ~ agreed well with literature
values. The predicted BTC based on the combined fit of samplings 2, 3 and 4shows a slightly better fit although the breakthrough still starts too early, andthe concentrations are still too low after about 120 mm of drain discharge. The
fastest BTC obtained for sampling 2 considerably overestimates the drainwater concentrations up to about 130ram of discharge, after which value theconcentrations are being somewhat underestimated.
From Fig. 10 we conclude that solute transport in this soil system cannot be
described adequately by assuming one average water flow velocity for theentire unsaturated zone. We believe that preferential flow may have caused
some of the observed discrepancies.
Although solute transport in the topsoil was described qv.ite well with the
CDE, we speculate that at a certain depth water flow co~ld have becof~eunstable, causing water to flow downwards along preferential flow paths (Glasset el., 1987). This type of preferential flow would result i~ an accelerated
transport through the subsoil. Thus, "slow" transport in the topsoil may havebeen followed by relatively "fast" transport in the subsoil, thus resulting in a
BTC that fits the measured data better.
Given the above considerations, we tried to predict the drain efIluent curve

247

by assuming a layered soilprofileconsisting of a uniform topsoil characterized
by relatively slow bror.lidetransport, and a more heterogeneous subsoil withrelatively fast transport caused by some type of shortcircuiting.Transport inthis layered system wa'.:calculated with a numerical finiteelement solution(Van Genuchten, 19~;) of the steady-state transport equation for layered
profiles. For transport in the topsoil we used the fitted parameter values aslisted in Table I for sampling 1. However, transport velocitiesfor the subsoil
were calculated from the fittedM V value of the second sampling. Since M V for
this second sampling also included the residence time of the tracer in the
topsoil, its value was corrected accordingly to yield a relative pore watervelocity M V of 7.21 for the subsoil of the corn field,and 11.6for the grass field.
Dispersivities were assumed to be the same as for the firstsampling time.
The B T C obtained in this way is shown in Fig. 11. The curve was found to
be rather insensitive to changes of +/- 5 cm in the thickness of the topsoil;
hence only the B T C for a topsoil of 25 cm is presented. Close examination ofFig. 11 shows that the peak concentration after about 1 5 9 m m discharge is still
somewhat underestimated but in general the measured B T C is predicted surprisingly well. Thus we conclude that the assumption of accelerated transport
in the subsoil is justified. It still remains an open question, whether thetransport parameters from the second, third and fourth samplings have aphysical basis, or should be considered as a set of parameters describing fast
breakthrough only.

o;
oO
IM

~'.

__1

t",d

Z~
i,-.,
I,--

/

--

PREDICTED/MEASURED BTC
MEASURED BTC• PREDICTED BTC

rf o
I--•
Zc~
C.~
Z o

°ii

CJ

.ILj

" 30.0

J

J

l

I

60.0

I

I

i

I

90.0

i

l

i

I

120,0

I

I

m

I

150.0

t~--.~

180.0

I ,~

I

I

210-0

,

i

I

I

240,0

I

I

'

!

270.0

l

J

I

I

300,0

z

~

,330.0,360,0DRAIN DISCHARGE (MM)Fig. 11. Measured and the predicted drain effluent bTC using a numerical model for a layered soil
system.

248

Additional evidence of accelerated (preferential) transport in the subsoil
can be found in the literature, especially in cases where a fine-textured layer
overlies a coarse-textured horizon. Starr et al. (1978, 1986) found similar results
for field measurements, while Glass et al. (]987) showed this fingering
phenomenon to be present during elaborate laboratory experiments. The
occurrence of preferential flow in our study is indirectly supported by: (1) the
presence of a textural change from a low permeable topsoil to a more permeable
subsoil at approximately 25 cm depth (see Materials and Methods), which may
cause unstable wetting fronts (Raats, 1973), and thus contribute to preferential
flow even during wet conditions (Glass et al., 1987); and (2) the presence 9f
water repellent layers throughout the profile; as shown by Hendrickx et al.
(1988), these layers can also trigger the formation of preferential flow paths.
Experimental evidence of the presence of preferred flow zones in this field
was further presented by Van Ommen et al. (1989, this volume), who were able
to visualize these zones by means of a iodide-starch color technique describedin Van Ommen et ai. (1988).
S U M M A R Y A N D CONCLUSIONS

Transport of a nonreactive tracer (bromide) was studied in an unsaturatedsaturated transport system. Two models were used to describe unsaturated
transport: the convection--dispersion equation (CDE) and a regional stochastic
model (RSM). Both transport models could describe the observed datareasonably well, although the CDE model performed slightly better inpredicting the amount of recove~,ed bromide.
Solute migration velocities and dispersivities for transport in the topsoil
generally agreed well with measured water contents and dispersivity values
found in the literature. In the subsoil, however, transport velocities were notin agreement with observed water contents, especially for the grassfield. The
accelerated breakthrough of the subsoil was probably caused by preferential
flow (shortcircuiting) originating at the transition between topsoil and subsoil.
Using measured tracer concentration distributions for the unsaturated zone,
the drain water quality could be predicted accurately only when the effects offast transport in the subsoil was accounted for in the unsaturated-saturated
transport system.
ACKNOWLEDGEMENTS

The experimental and theoretical objections of this study could only berealized thanks to the much-appreciated inputs by Mr. Stricker, Mr. Post, Mr.
van den Abeele and Mr. Beem of the Agricultural University Wageningen; Mr.
Verstrate of the Institute fur Inland Water Management and Waste Water
Treatment; Mr. van der Plicht and Mr. De Klein of the Water Quality Control
Service East-Gelderland; and Mr. Reijnders and Mr. Neele of the National
Institute of Public Health and Environmental Hygiene. We would like to249

acknowledge also the Province of Gelderland and the Institute for Inland
Water Management and Waste Water Treatment for financially supporting the
project, and Mr. Veehof and Mr. Beusink for their permission to use their fields
during the experiment. Finally, Mr. Hopmans is acknowledged for hiscomments on the manuscript, and Mrs. Arissen for typing and editing the draft
versions of the paper.
REFERENCES

Amoozegar-Fard, A., Nielsen, D.R. and Warrick, A.W., 1982. Soil solute concentration distributions for spatially varying pore water velocities and apparent diffusion coefficients. Soil Sci.
Soc. Am., J., 46: 3-9.
Anderson, M.P., 1979. Using models to simulate the movement of contaminants through groundwater flow systems. Dep. Geol. Geophys., University of Wisconsin, Madison, Wisc., CRC
Critical Reviews in Environmental Control, pp. 97-156.
Beese, F. and Wierenga, P.J., 1983. The variability of the apparent diffusion coefficient in undisturbed soil columns. Z. Pflantzenernaehr. Bodenk., 146: 302-3!5.
Biggar, J.W. and Nielsen, D.R., 1976. Spatial variability of the leaching characteristics of a fieldsoil. Water Resour. Res., 12(1): 78-84.
Bowman, R.S., 1984. Evahiation of some new tracers for soil water studies. Soil Sci. Soc. Am., J.,
48: 987-993.
Bowman, R.S. and Rice, R.C., 1986. Transport of conservative tracers in the field under intermittent
flood irrigation. Water Resour. Res., 22(11): 1531-1536.
Bresler, E. and Dagan, G., 1979. Solute dispersion in unsaturated heterogeneous soil at field scale:
II. Applications. Soil Sci. Soc. Am., J., 43: 467-472.
Bresler, E. and Dagan, G., 1981. Convective and pore scale dispersive s~lute transport in unsaturated heterogeneous fields. Water Resour. Re 3. 17(6): 1688-1993.
Bresler, E. and Dagan, G., 1983a. Unsaturated flow in spatially variable fields. 2. Application ofwater flow models to various fields. Water Resour. Res., 19(2): 421-428.
Bresler, E. and Dagan, G., 1983b. Unsaturated flow in spatially variable fields. 3. Transport models
and their application to two fields. Water Resour. Res., 19(2): 429-435.
Bruyn, J., 1984. Groundwater quality - manuring: problems with nitrate in East-Gelderland. H20,
17:502 (in Dutch).
Cassel, D.K., Van Genuchten, M.Th. and Wierenga~ P.J., 1975. Predicting anion movement indisturbed and undisturbed soils. Soil Sci. Soc. Am., Proc., 39: 1015-1019.
Dagan, G., 1986. Statistical theory of groundwater flow and transport, pore to laboratory,laboratory to formation, and formation to regional scale. Water Resour. Res., 22(9): 120S-134S,
Dagan, G. and Bresler, E., 1979. Solute dispersion in unsaturated heterogeneous soil at fie!el scale:
I. Theory. Soil Sci. Soc. Am., J, 43: 461-467.
Dagan, G. and Bresler, E., 1983. Unsaturated flow in spatially variable fields 1. Derivation of~.aodels of in£1tration and redistribution. Water Resour. Res., 19(2): 413-420.
Danckwerts, P.V., 1953. Continuous flow systems: distribution of residence times. Chem. Eng. Sci.,
2: 1-17.
Dane, J.H., Reed, R.B. and Hopmans, J.W., 1986. Estimating soil parameters and sample size bybootstrapping. Soil Sci. Soc. Am. J., 50: 283-287.
Eldor, M. and Dagan, G., 1972. Solutions of hydrodynamic dispersion in porous media. Water
Resour. Res, 8(5): 1316-1331.
Ernst, L.F., 1973. The determination of residence times in case of saturated groundwater flow. Nota
755 !.C.W., "rst. Land Water Manage. Res., Wageningen, 42 pp. (in Dutch).
Funk, W.H, 1987. Water quantity and quality: the growing crisis. Washington's Land and People,
Vok 1 (1): 5-7.
Gelhar, L.W. and Wilson, J.L., 1974. Groundwater quality modeling. Ground Water, 12: 399-408.

250
Glass, R.J,,Parlange, J.-Y. and Steenhuis, T.S., 1987. In: Yu-Si Fok (Editor),InfiltrationDevelop.
ment and Application. Proc. Int. Conf. InfiltrationDev. Appl., pp. 66-81.
Groenendijk, P.: 1986. Onderzoek naar de effecten van wateraanvoer en peilveranderingen inagrarische gebieden op de waterkwaliteit in natuurgebieden. Deel 2: Modellering van twee.
dimensionaal transport van een conservatieve stofin de verzadigde zone. Nota I C W 1736, 38 pp.
Hendrickx, J.M.H., Dekker, L.W., Bannink, M.H. and Van Ommen, H.C., 1988. Significance of soil
survey for agrohydrological studies. Agric. Water Manage., 14: 195-208.
Jury, W.A., 1982. Simulation of solute transport using a transfer function model. Water ~¢sour.
Res., 18(2):363--368.
Jury, W.A., Sposito,G. and White, R.E., 1986.A transferfunction model of solute transport throughsoil 1. Fundamental concepts. Water Resour. Res., 22(2):243-247.
Lindstrom, F.T., Hague, R., Freed, V.H. and Boersma, L., 1967. Theory on the movement of someherbicides in soils:linear diffusion and convection of chemicals in soil.Environ. Sci. Tech., 1:562-565.
Neele, J., 1988. Ion-chromatographic determination of bromide in soil water extracts and groundwater. I~.tern. Rep., Natl. Inst. Public Health Environ. Hygiene, Bilthoven (in prep.).
Nielsen, D.R. and Biggar, J.W., 1962. Miscible displacement: III. Theoretical considerations. Soil
Sci. So::. Am., Proc., 26: 216-221.
Nielsen, D.R., Van Genuchten, M.Th. and Biggar, J.W., 1986. Water flow and solute transport
processes in th~ unsaturated zone. Water Resour. Res., 22(9): 89S-108S.
Par~er, J.C. and V~.~ Genuchten, M.Th., 1984a. Determining transport parameters from laboratory
and field tracer experiments. Va. Agric. Exp. Stn. Bull., 84-3, 96pp.
Parker, J.C. and Van Genuchten, M.Th., 1984b. Flux-averaged and volume.averaged concentrations in continuum approaches to solute transport. Water Resour. Res., 20(7): 866-872.
Raats, P.A.C., 1973. Unstable wetting fronts ir uniform and non-uniform soils. Soil Sci. Soc. Am.,
Proc., 37: 681-685.
Raats, P.A.C., 1983. Responses of hydrological systems to ~hanges in water quantity and quality.
Relatiov of groundwater quantity and quality. Proc. Hamburg Syrup., IAHS Publ. No. 146.
Rao, P.V., Rao, P.S.C., Davidson, J.M. and Hammond, L.C., 1979. Use of goodness-of-fit tests forcharacterizing the spatial variability of soil properties. Soil Sci. Soc. Am., J., 43: 274-278.
Richter, G. and Jury, W.A., 1986. A microlysimeter field study of solute transport through astructured sandy loam soil. Soil Sci. Soc. Am., J., 50: 863-868.
Schulin, R., Wierenga, P.J., Fl~ihler,H. and Leuenberger, J., 1987a. Solute transport through astony soil.Soil Sci. Soc. Am., J., 51'.36-42,
Schulin, R., Van Genuchten, M.Th., Fliihler,H. and Ferlin,P., 1987b. Experimental and theoreticalanalysis of solute transport in a stony fieldsoil.Water Resour. Res., 23(4):1785-1794.
Shapiro, S,S, and Wilk, M.B., 1965, An analysis of variance test for normality (complete samples).
Biometrika, 52: 591-611.
Sposito,G., Jury, W.A. and Gupta, V.K., 1986. Fundamental problems in the stochasticconvectiondispersion model of solute transport in aquifers and fieldsoils.Water Resour. Res., 22(1):77-88.
Starr, J.L., DeRoo, H.C., Frink, C.R. and Parlange, J.-Y., 1978. Leaching characteristics of alayered fieldsoil.Soil Sci. Soc. Am., J., 42: 386-391.
Starr, J.L.,Parlange, J.-Y.and Frink, C.R., 1986. Water and chloride movement through a layered
fieldsoil.Soil Sci. Soc. Am., J., 50: 1384-1390.
Studygroup Hupselse Beek, 1974. Third report of the hydrologic research in the catchment of the
Hupselse Beek. Dep. Hydraul. Catchment Hydrol., Agric. Univ., Wageningen.
Sudicky, E.A., 1986. A natural gradient experiment on solute transport in a sand aquifer: Spatial
variabilityof hydraulic conductivity and itsrole in the dispersion process.Water Resour. Res.,
22(13): 2069-2082.
Van Genuchten, M.Th., 1986. A numerical model for water and solute movement in and below theroot zone. Model description and user manual. Unpubl. Res. Rep., U.S. Salinity Lab., USDA/
ARS, Riverside, Calif.
Van Ommen, H.C., 1985. Systems approach to a coupled unsaturated/saturated groundwater
quality model, including adsorption, decomposition and bypass. Agric. Water Manage., 10:251
193-203.
Van Ommen, ~.C., 1986. Influence of diffuse sources of contamination on the quality of outflowing
groundwater including non-equilibrium adsorption and decomposition. J. Hydrol., 88" 79-95.
Van Ommen, H.C., Dekker, L.W., Dijksma, R., Hulshof, J. and Van der Molen, W.H., 1988. A newtechnique for evaluating the presence of preferential flow paths in non.structured soils. Soil Sci.
Soc. Am., J., 52: 1192-1193.
Van Ommen, H.C., Dijksma, R., Hendrickx, J.M.H., Dekker, L.W., Hulshoff, J. and Van den Heuvel,
M., 1989. Experimental assessment of preferential flow paths in a field soil. J. Hydrol., 105:253-262 (this volume).
Van de Pol, R.M., Wierenga, P.J. and Nielsen, D.R., 1977. Solute movement in a field soil. Soil Sci.
Soc. Am., J., 41: 10-13.
Van der Veen, C., 1986. Tsjernobyl and drinking water supply. H20, 19: p. 266. (in Dutch).
Warmerdam, P.M.M., 1981. Effect of wind on rainfall measurement. Comparison of raingauges.
H20, 14: p. 16 (in Dutch).
Wierenga, P.J., 1977. Solute distribution ]~rofiles computed with steady state and transient water
movement models. Soil Sci. Soc. Am. J., 41: 1050-1055.
W~sten, J.H.M., Bouma, J. and Stoffelsen, G., 1985. Use of Soil Survey data for regional soil water
simulation models. Soil Sci. Soc. Am., J., 49: 1238-1244.