Published 1999

9

Principles of Solute Transport
FEIKE J. LEU AND MARTINUS TH. VAN GENUCHTEN

u.s. Salinity Laboratory
USDA -ARS, Riverside, California

I. INTRODUCTION
Solute transport in soils has traditionally been described with the convection-dispersion equation (CDE). This equation incorporates two constitutive
transport processes: (i) movement as a result of liquid flow, and (ii) spreading asa result of known and unknown processes such as diffusion and small-scale variations in the water flow velocity. In Part II we will review several modes of solutetransport; convection, diffusion in free liquids and subsequently in soils, mechanical dispersion, and (hydrodynamic) dispersion. The governing transport equation is derived from mass balance principles in Part III. Many other processes
may affect the movement and fate of solutes in soils. Part IV is devoted to theimportant process of linear and nonlinear solute adsorption.

II. THE CONVECTION-DISPERSION EQUATION
Consider transport of a chemical species in a three-phase, soil-air-water
system and assume that the chemical species (the solute) is completely misciblewith the water (the solvent). At the macroscopic level, the mass balance equation
for a solute species in a control volume of the soil can be given as

aec =-"'J',,+ R," +R C

at

V

'

II '

(I

[IJ

where e is the volumetric water content (L 3L-3), C is the concentration of thesolute expressed in solute mass per solvent volume (ML -3), t is time (T), V' . js isthe divergence of the solute flux (ML- 3r l), H, denotes arbitrary solute sinks «0)or sources (>0) (ML- 3r\ and the last term denotes injection (>0) or pumping«0) of water with concentration Co into or from the control volume at a rate RII

Copyright © 1999 American Society of Agronomy, Crop Science Society of America, Soil Science
Society of America. 677 S. Segoe Rd .. Madison. WI 53711. USA. Agricultural Drainage. Agronomy
Monograph no. 38.

331

LEU & VAN GENUCHTEN

332(L 3L-3r i ) . Classical descriptions of solute transport usually distinguish the following two components of the solute flux,js

js = qC + jD

[2]

where q is a vector quantifying the Darcy water flux (Lr\ and jD is anautonomous solute flux (e.g., Bolt, 1982b) to quantify transport caused by a gradient in the solute concentration (ML- 2r i ). Both solute fluxes in Eq. [2] areexpressed per unit area of porous medium.
The passive movement of a solute with flowing water, as quantified by thesolute flux qC, is referred to as convection or also as advection. We assume herethat the Darcy water flux, q, is known . This flux is generally a function of timeand position; many (analytical) studies pertaining to solute transport simplify the
problem by assuming steady-state one-dimensional flow or by changing independent variables (e.g., cumulative infiltration instead of time). We will often treat
the velocity as a scalar rather than a vector. Because dissolved substances movein a passive fashion, convective transport can be readily quantified when the solvent flux, q, is known. Hence, solute transport cannot be studied without adequately quantifying the flow regime. Although macroscopic values for the water
flux may be available, unknown microscopic variations in water velocity willlead to unequal solute movement in the flow direction. This movement is quantified with the dispersive flux. Macroscopic variations in flow, which manifest
themselves as preferential flow, require modifications of the transport model. Leijet al. (1998b, Chapter II) discuss a relatively simple model for transport in a soil
where a mobile and immobile flow domain can be distinguished.
Miscible displacement experiments involving steady-state flow in watersaturated soil columns have shown that if the solute concentration of the solution at the inlet of the column is changed abruptly at time t = 0, the observed
breakthrough of solute at the column outlet at times t > 0 will not exhibit anequivalent abrupt change (Nielsen & Biggar, 1961). Instead, the solute concentration will change more gradually with time due to hydrodynamic dispersion.
Hydrodynamic dispersion is caused by mechanical dispersion, resulting fromflow in a porous medium, and molecular or ionic diffusion. Whereas flow variations at the pore scale can ordinarily not be observed directly, they can be quantified with the use of a tracer. Hydrodynamic dispersion has been discussed at
great length in the literature, and we will only discuss some of the more important underlying concepts of the dispersion process. We will use the diffusion
process as a starting point, because of the analogy between molecular diffusion
and mechanical dispersion.
A. Diffusion
Molecular or ionic diffusion is an important mechanism for solute transport
in the liquid phase in directions of low mean pore-water velocities with relatively
little or no water flow. Solute molecules in a free solution possess random thermal motion which causes an exchange of molecules between adjacent volume
elements. A net transfer of molecules of a solute species usually occurs when theconcentration of the species differs in adjacent volume elements, i.e., more parti-

PRINCIPLES OF SOLUTE TRANSPORT333

Table 9-1. Diffusion coefficients of strong electrolytes in water for three electrolyte
levels (after Weast et aI., 1984).
Do (cm 2 d~l)Solute

O.OIM

O.IM

1.0M

NaCl
NaBr
NaI
KCl
KEr
KI
KN0 3
KClO 4CaCl 2
BaCl 2
Na2S04
MgS0 4

1.335

1.281
1.311
1.313
1.593
1.619
1.611

1.282
1.379
1.436
1.635
1.706
1.784

0.959
1.001

1.039
1.019

1.656

1.595
1.547
1.026
1.093
0.970
0.613

cles move from the elements with higher concentrations to those with lower concentrations than vice versa. The resulting process is referred to as diffusion. The
mass flux, hi (ML 2r'), due to molecular diffusion is given by Fick's first law
which, for one-dimensional diffusion, yields

de~Iit = -Do-

dX

[31

where x is distance (L), and Do is the coefficient of molecular diffusion for a freeor bulk solution (L 2!'), whose value depends on properties of both the soluteand solvent. The experimentally observed proportionality between ii/if and theconcentration gradient can be described at the molecular or ionic level with a balance of forces. The driving force for particle movement from higher to lower
concentrations is a gradient in the chemical potential. For mixtures with ideal
behavior the chemical potential can be expressed in terms of its mole fraction.
For solutions with nonideal mixing behavior the activity coefficient of the soluteneeds to be determined. Ionic activity coefficients can be estimated with theextended Debye-Hiickel equation or the Davies extension for solutions up to 0.5
M (cf. Novozamsky et aI., 1978). Activity coefficients for a greater concentration
range-up to 16 M-are estimated with the Pitzer virial equations (Pitzer, 1979;
Harvie & Weare, 1980). Table 9-1 contains values for Do as published by Weast
et aI. (1984) for some common salts in water at three different molarities. Manyother sources exist (cf. Table 43-1 in Kemper, 1986; Table I in Ellis et al., 1983).
For large spherical particles moving through a viscous solvent, the frictional force can be expressed with Stokes' law. The diffusivity may then be givenby the Stokes~Einstein approximation
[41

334

LEU & VAN GENUCHTEN

where Rg is the gas constant (ML 2r 2K- l mol- l ; 8.314 JK-Imol- I), T is theabsolute temperature (K), Ai is Avogadro's number (mor l ; 6.022 1023 mol-I), rp
is the particle radius (L), and 11 is the kinematic viscosity (ML-lr l). Several
other theories exist for diffusion in gases and liquids (e.g., Bird et aI., 1960).
Equation [4] was derived for colloidal particles in a dilute suspension.
Codiffusion or counterdiffusion occurs in systems with multiple ion species.
Diffusion rates for individual species as predicted solely by Fick's model would
violate the electroneutrality principle. The ionic diffusion flux consists of terms
for ordinary Fickian diffusion and electric transference of ions. The corresponding diffusion coefficient is related to the ionic mobility using the Nernst-Planck
equation (cf. Helfferich, 1962). Additional approaches have been used for diffusion in multicomponent systems, such as the Maxwell-Stefan equation (Standart
et aI., 1979).
To characterize diffusion in a porous medium such as soil, the diffusivity ina free solution is typically adjusted with terms accounting for a reduced solution
phase (a smaller cross-sectional area available for diffusion), and an increased
path length. A general treatment of diffusion in soils can be found in Olsen and
Kemper (1968), Nye and Tinker (1977), and Nye (1979). The macroscopic diffusive flux per unit area of soil can be written as[5]
where D. is the coefficient of molecular or ionic diffusion for the liquid phase of
the soil. The diffusion coefficients for the soil liquid and a free liquid are relatedby (e.g., Epstein, 1989)

D =
*

(

D()

D

)2 =--'2..=Dr- r 2 - () {/
L*/L

[6]

where L* is the actual path length for diffusion in the soil (which depends on 8),
L is the (straight) length of the soil, 't = LIL is known as the tortuosity, 't 2 is the
tortuosity factor, while 'to = (UL)2 is designated as an apparent tortuosity factor.
Note that the tortuosity LIL appears twice in Eq. [6] so as to account forchanges in the concentration gradient along the streamline with length L" andthe variability in travel distance, compared to the fixed concentration difference
and travel distance along a straight path with length L in free water (Olsen &Kemper, 1968). Values for the tortuosity are typically greater than unity. The
tortuosity factor 't 2 is the ratio of the diffusivities in a free solution and theporous medium.
There is no uniformity in the literature, especially among different disciplines, regarding the use of the terms tortuosity and tortuosity factor. The use of
tortuosity for 't and tortuosity factor for 't2 is relatively standard in the geochemical and chemical engineering literature (e.g., Dykhuizen & Casey, 1989; Epstein,1989). In the soil science and hydrology literature the terms tortuosity (Bear &Verruijt, 1987; Luckner & Schestakow, 1991) or tortuosity factor (Hillel, 1980;
Bresler et a!., 1982; van Rees et aI., 1991) are used for 'to in Eq. [6]. The corre-

335

PRINCIPLES OF SOLUTE TRANSPORT

sponding values are less than unity. Adding to the confusion is the fact that someauthors include the water content in their definition of tortuosity (e.g., Dykhuizen& Casey, 1989), while sometimes the effects of solute adsorption (retardation)
are implicitly included in either 1 or the porous medium diffusion, D .]. (Nye,1979; Robin et aI., 1987).
Electrical conductance measurements have frequently been used to characterize the pore structure of porous media (cf. Sauer et aI., 1955). In geophysics,
I I( 1,,9) is often referred to as the formation factor and because of the saturated
conditions the water content is typically replaced by the porosity. The formation
factor is analogous to the cell constant in conventional conductometry. Sauer etal. (1955) applied a model to distinguish between the conductance by the solution
and the solid surface of a porous medium. The model was used by Gast (1966) toestimate diffusion in dilute bentonite pastes while a similar model was employed
by Rhoades et al. (1976) for soils. Diffusion coefficients for cationic species
undergoing adsorption have been estimated from soil electrical conductivity by
Palmer and Blanchar (1980) and Rhue (1992).
For unsaturated conditions it is convenient to quantify the dependency ofthe diffusion coefficient on water content. This may be accomplished by modifying the tortuosity term, for instance according to Millington (1959) and Millington and Quirk (1961). Their expression has been used frequently to describe
gaseous diffusion in soils but it also can be applied to diffusion in the liquid
phase of variably saturated porous media[7]
where £ is the soil porosity (L 3L-\
Diffusion coefficients in soil systems are usually determined by mathematically analyzing solute concentration profiles in the soil as a function of time orposition. Such concentration profiles will develop after systems spiked with different solute concentrations have been brought in contact with each other. Van
Rees et al. (1991) measured diffusivities by allowing diffusion from: (i) a spiked
solution into the soil, (ii) a spiked soil into a solution, and (iii) a spiked soil into
the soil. The first method proved to be quite convenient fo r in situ laboratory
measurements.
Kemper and van Schaik (1966) measured diffusivities of NaCI and CaClc
solutions in a bentonite-clay using different electrolyte levels and water contents.
They found that the influence of the electrolyte level was relatively small. Theirfindings pertaining to the influence of water content were described by Olsen and
Kemper (1968) using the relationshipDaD.'." == -'-'
8 exp(b8)[8]in which a and b are empirical coefficients- a is a function of the surface area, avalue that typically ranges between 0.00 I and 0.005, and b '" I O. Robin et al.
(1987) observed an increase in diffusivity with an increasing electrolyte level andLEU & VAN GENUCHTEN

336

Table 9-2. Selected diffusion coefficients for aqueous solutions in clay and soil materials.
D,

Pb
g/cm3

Comments

1.46
0.64
0.70
1.47
0.79
1.66
1.46
0.69
0.97
0.46
Robin et al. (1987)+

1.42
1.42
1.35
1.35
1.42

Spiked water on top of sediment (e = 0.42)
Spiked water on top of sediment (e = 0.42)
Spiked water on top of sediment (8 = 0.45)
Lake water on top of spiked sediment (e = 0.45)
Sediment on top of spiked sediment (e = 0.42)

0.19
0.11
0.20
0.11
0.36
0.11
0.40
0.12
0.54
0.11
0.11
0.56
Oscarson et al. (l992)§

1.63
1.61
1.62
1.61
1.63
1.64

25°C
25°C
60°C
60°C
90°C
90°C

To

cm /d
van Rees et al. (1991)t
2

0.19
0.15
0.10
0.05

0.33
0.27
0.17
0.08

0.90
1.12
1.31
1.50

t Water eH 20) tracer in litoral sediment.
t Chlorine (36CI) tracer in bentonite-sand mixture using a spiked (Co'" 0.27 M) andunspiked soil plug.
25 l) tracer in compacted bentonite using spiked and un spiked clay plugs.

§ Iodineetemperature. Table 9-2 provides a few values for diffusion coefficients in clays
and soils with estimates for 'ta o More soil diffusivity data can be found by consulting the reviews of Hamaker (1972), Nye (1979), and Ellis et al. (1983),
among others.

B. Mechanical Dispersion
Variations in water flow in a porous medium will lead to mechanical dispersion. Several factors contributing to mechanical dispersion are illustrated withhypothetical tracer particles in Fig. 9-1. Dispersion may occur because of: (a) the
development of a velocity profile in an individual pore such that the highest
velocity occurs in the center of the pore, and presumably little or no flow at thepore walls; (b) different mean flow velocities in pores of different sizes; (c) themean water flow direction in the porous medium is different from the actual
streamlines within individual pores, which differ in shape, size, and orientation;and (d) particles originating from different pores end up in the same pore, andvice versa. All of these processes contribute to increased spreading, initially steep
concentration fronts will become smoother during movement along the main
flow direction .
The effects of dispersion may be illustrated by considering a hypotheticallaboratory experiment in which water and a dissolved tracer are applied to an ini-

PRINCIPLES OF SOLUTE TRANSPORT

337

(a)

(b)

( c)

(d)

Direction of average flow
Fig. 9-1. Schematic concepts contrihuting to mechanical dispersion.

tially tracer-free, uniformly packed soil column of length L (Fig. 9-2). The column is subjected to steady-state water flow with a uniform water content. As
more of the tracer is added, the initially very sharp concentration front near thesoil surface becomes spread out because of dispersion. Eventually a smooth andsigmoidally shaped effluent curve can be monitored at the column exit as shownin Fig. 9-2d. In the absence of dispersion, the front of a perfect inert tracer willtravel as a square wave through the column (a process often called piston flow) toreach the bottom of the column at time t = Uv, where v is the average pore-wateror interstitial velocity. This velocity is the ratio of the Darcian water flux density,q, and the volumetric water content, 8. Notice that v is expressed per area of fluid
whereas q is defined per area of (bulk) soil. We emphasize this difference since
the literature is not always clear whether a reported velocity includes the volumetric water content or porosity-the issue also is important for diffusion or dispersion coefficients. For piston flow the tracer reaches the column exit after one
pore volume of tracer solution has been injected (or collected at the column exit).
Pore volume is defined as the amount of water stored in that column.
The degree of spreading is usually related to the solute travel time,
although some constraints do exist on the amount of spreading. Dispersion, asFig. 9-2. Hypothetical laboratory tracer experiment: (a) column of soil, (b) influent curve, (e) concentration distributions inside the column, and (d) breakthrough curves with and without dispersion(modified after van Genuchten, 1988).

sketched in Fig. 9-1 a, is limited because of transverse molecular diffusion whichcauses solutes to move from the center of a tube to areas near the pore walls, andvice versa. This transverse diffusion counteracts spreading caused by the variations in the longitudinal flow velocity. Dispersion also is limited since capillaries
in a porous medium generally are not independent cylindrical tubes, but branch
and join or rejoin each other at distances characteristic of the pore or particle-size
distribution of the medium. This promotes lateral mixing of solutes from different pores at certain mixing lengths (Bolt, 1982b) as sketched in Fig. 9-1d.
Although molecular diffusion and mechanical dispersion are different
processes, the macroscopic solute flux due to mechanical dispersion is often conveniently expressed with Fick's first law of diffusion (Scheidegger, 1974). For
one-dimensional dispersion in a uniform isotropic porous medium this leads toPRINCIPLES OF SOLUTE TRANSPORT

339

[9]
where J"i' is the dispersive solute flux (ML- 2r l) and D"i' is the coefficient of
mechanical dispersion (L2rl). For three-dimensional dispersion, D"i, should formally be written in tensor form. The components of the symmetric dispersion
tensor for an isotropic porous medium can be given as (Bear & Verruijt, 1987)[10]
where Ivl denotes the magnitude of the pore-water velocity with Vi as the ith component (Lr\ Oij is the Kronecker delta (oij = I if i =) and oij = 0 if it:- i), and aL
and aT are respectively the longitudinal and transverse dispersivity (L).
Mechanical dispersion, as sketched in an idealized fashion in Fig. 9-1 a,can be reversed by changing the flow direction to make a smooth front steep
again. Some reversibility has indeed been observed (Heller, 1960; Scotter &Raats, 1968). In actual porous media, however, dispersion is not reversible since
mixing erases antecedent concentration distributions, as illustrated in Fig. 9- 1d.
Absolute values are therefore used for V in Eq. [10].
In a Cartesian coordinate system where the pore-water velocity v is givenby (vx , V \, vJ , the components of the dispersion tensor can be written as (Bear &Verruijt, 1987)

Du =aTlvl+(aL-aT)v~/lvl=[aT(v~ +vn+aLv~lIlvl
D,r = Dr, = (a L- aT)V, vl'/Ivl

D,~ = D~., =(al-aT)v, v~/lvl
D"r == arlvl+ (a L-aT)v; /Ivl = [aT(v~ +vn+ aLV~ lIlvl[11]

D,~ == D~\ == (a L- aT) V, ,,:)Ivl

D~~ = aTlvl +(a L - aT)v; /Ivl == [aT(v; +vn+ aLv; lIlvl
For uniform water flow parallel to the x axis of the coordinate system only threecomponents remain, i.e.,[12lwhere D" is the coefficient of longitudinal (mechanical) dispersion, and Dl'r andD~~ are the coefficients of transverse dispersion.

Relatively little work has been done to quantify the processes of diffusionand dispersion in unsaturated soils. No simple relationship exists between watercontent and the coefficient of mechanical dispersion (Krupp & Elrick. 1968). Overall , dispersion willlike\y increase with decreasing water content because of microscopic variations in the sequences of fully and partially filled pores, and because of
the increased importance of transport in surface films at low water contents.

LEU & VAN GENUCHTEN

340

A foundation for the Fickian description of mechanical dispersion is provided in the classical paper by Taylor (1953), who treated fluid flow in a cylindrical tube. Taylor's analysis provides an explanation of porous medium dispersion.
We will review his analysis since it has often been used to investigate transport in
repacked soil columns in the laboratory or in field soils (e.g., Nielsen & Biggar,
1961; Gupta & Bhattacharya, 1983). Taylor (1953) considered laminar flow in acircular tube, which was filled with water, using the parabolic velocity profile[ 13]
where q" is the maximum velocity of a microscopic volume element of water at
the axis (Lr'), r is the radial distance from the axis (L), and ro is the tube radius(L). It can be shown that the mean velocity over a cross section of the tube, <q>,is equal to qJ2. The continuity equation, with C(x,r, t) as dependent variable, forcylindrical coordinates and for a constant diffusion coefficient Do yields (Bird etaI., 1960)[14]
In order to quantify the effect of the variable fluid velocity and molecular
diffusion on solute spreading, Taylor assumed that diffusion is only significant
for radial transport, but negligible for longitudinal transport. Deviations from themean velocity and the mean concentration in the transverse direction are givenby (Fischer et aI., 1979)

q'=q(r)-<q>[15]

C' = C( r)- < C >[16]

where the averaged values are defined asf",<q>= ~ rq(r)dr=!lsL
r,) 0
2

f
r"

< C> = ~ rC(r)dr
I'v 0

with

ekC> =0
dr

[17]

[18]

In order to quantify dispersion, a solution needs to be obtained for C', i.e.,the deviation in concentration from the cross-sectional mean <C>. With the helpof Eqs. [15] and [16], Eq. [14] without the longitudinal diffusion terms can be
written as

PRINCIPLES OF SOLUTE TRANSPORT

341

where variations of <C> in the radial direction can obviously be ignored. Equation [19] may be rewritten using the following coordinate transformation~ = x- < q > t

r =t

[20]

Z= r/t;)to give

This equation provides a fairly detailed description of solute transport in
the tube, further simplifications are needed to obtain an analytical solution. Aftersome time a steady state is reached. Taylor (1953) assessed the magnitude of thedifferent terms in Eq. [21] and neglected the left-hand side as well as the last termon the right-hand side (cf. Fischer et aI., 1979). In this case longitudinal convective and transverse diffusive transport will balance each other:a2 e' I ae' q'r} a < C>--+--=----az 2 Z az
D(I
a~

[221

To solve the steady-state profile, q' in Eq. [22] is written in terms of qo with
Eqs. [13] and [15]. After integrating twice, the following local deviation from themean concentration is obtained

e'= t;,2 q() a<C>
DII

a~
."

[£_ Z4J+a
8

16

InZ+aoI-[23]

where (XI and (X2 are integration constants. Because no solute transfer occurs
across pore walls, there is a zero concentration gradient at Z = I and hence (XI =O. The parameter (X2 is determined by the initial condition.
The dispersive flux, i.e., mass transport relative to the coordinate axis moving with velocity <q>, for an initially solute-free tube can be given byj.
'dIS=~J';' r 'e'dr=- q(~I(~ a<C>
r2
(I

0

q

192D

0

a~

[24]

."

This flux may be expressed in terms of Fick's law for macroscopic dispersive
transport as[25]
and the following expression for the coefficient of mechanical dispersion can beobtained:

LEU & VAN GENUCHTEN

342
2 2

2

= r" qo = r" < q>

D.
alS

192Do

2

48D()

[26]

Note that D dis is inversely proportional to the coefficient of molecular diffusion,

DO" Also note that Eq. [26] is similar to Eq. [12] since <q> is inversely proportional to r~.
Taylor's approach is restricted to flow and transport in a cylindrical tube
after a certain travel time has elapsed whereas longitudinal molecular diffusion isneglected. Aris (1956) dropped these assumptions and investigated solute dispersion for flow in systems of arbitrary geometries using spatial moments, mp(p=0,1,2,oo.)[27]

Solute spreading increases linearly with time for Fickian dispersion

dm2 -2- ma Ddt[28]

where D is an apparent dispersion coefficient, ma denotes the total amount ofsolute mass per cross-sectional area, and m2 is equivalent to the variance of theconcentration curve in the longitudinal direction. The apparent dispersion coefficient accounts for different mechanisms contributing to spreading-in the case of
mechanical dispersion and molecular diffusion D is referred to as the coefficientof hydrodynamic dispersion. Aris's analysis yields expressions for m2 for arbitrary geometries and times, including the early stage for which Taylor's analysis
is not applicable. Aris (1956) included molecular diffusion and obtained the following expression for the apparent dispersion coefficient2

D=D + /0;, <q>

o

Do

2

[29]

where K is a constant depending on the geometry of the cross section and the relative variation of velocity and diffusion coefficient.
Research toward a better understanding of mechanical dispersion in porous
media is still evolving. The development of a comprehensive theory for such dispersion is considerably more difficult than for solute movement in, for instance, awater-filled tube. One major reason is that microscopic flow pathways and thecomplex geometry of the solid phase cannot be quantified precisely. Some of theearlier models for dispersion in porous media simplify the pore structure as a network of capillary tubes (e.g., de Josselin de Jong, 1958; Saffman, 1960). Most ofthese models assume a random-walk process as the solute particle moves throughthe capillary network. Reviews of the dispersion process have been given byFried and Combarnous (1971), Scheidegger (1974), and Dullien (1992). A rigorous theory for dispersion in spatially periodic porous media was developed fromthe concepts of Brownian motion theory by Brenner and coworkers (Brenner &Gaydos, 1977; Brenner, 1980; Brenner & Adler, 1982; Edwards et aI., 1991).

PRINCIPLES OF SOLUTE TRANSPORT

343

C. Hydrodynamic Dispersion
The macroscopic similarity between diffusion and mechanical dispersion
has led to the practice of describing both processes with one coefficient of hydrodynamic dispersion (i.e. , D = D",,+ D *).This practice is consistent with resultsfrom typical laboratory and field experiments, which do not permit a distinction
between mechanical dispersion and molecular diffusion. The following contributions from molecular diffusion (Eq. [5]) and mechanical dispersion (Eq . [9]) to ahydrodynamic dispersive flux , 1D (Eq. [2]) , can be distinguished

10 = ~lit +~Ii.'[30]

Diffusion is insignificant compared to mechanical dispersion for many transport
problems. The terms mechanical and hydrodynamic dispersion are therefore
often used interchangeably and solute spreading is described simply in terms of
one dispersion coefficient. We will often use the term dispersion with the understanding that hydrodynamic dispersion is meant. Sometimes, nonequilibrium
phenomena caused by a macroscopically nonuniform porous medium also areincluded in such a dispersion coefficient (Bolt, 1982b).
Saffman (1960) and Pfannkuch (1963) investigated the relationshipbetween the dimensionless dispersion coefficient, DJD", and the Peclet number
for molecular diffusion, P = vd/D", where DL is the longitudinal coefficient ofhydrodynamic dispersion and d is the mean particle diameter. Figure 9- 3 shows
this widely publicized relationship between DJD" and P. Notice that molecular
diffusion is relatively important at low pore-water velocities , whereas mechanical dispersion is mostly responsible for solute mixing when P becomes greaterthan five . A similar relationship was shown for the transverse dispersion coefficient, Dr, by Crane and Gardner (1961). Values for the transverse dispersivity,._ - - -10-- ------- -/-

I:;;::-.-~~_T--I --- ~~- - - - - - - 1-----

10

P=vd/D.

Fig. 9- 3. The dimensionless longitudinal di spersion coeffici ent, OdD", as a function of the Peelet
number. P = vdlDo (modified after Pfannkuch, 1963).

LEU & VAN GENUCHTEN

344

Table 9-3. Selected information on field-observed dispersivities (after Gelhar et aI., 1992).
Ev

Glaciofluvial
sandMedium to coarse
sand overlying
silty sand and till
Brecciated basalt

0.33

crnJd
9

0.39

43

Layered medium
sand
Sand and gravelwith clay lenses
Sand andsandstone
Sand

0.35

Referencesite name

Aquifer material

Freyberg (1986)
Borden, Canada
Garabedian et ai. (1988)
Cape Cod, MA
Gelhar (1982)
Hanford, WA
Huyakorn et ai. (1986)
Mobile,AL
Lau et ai. (1957)
Berkeley, CA
Mercado (1966)
Yavne, rsrael
Molinari and Peaudecerf(1977), Sauty (1977)
Bonnaud, France
Rajaram and Gelhar(1991) Borden, Canada
Valocchi et ai. (1981)
Palo Alto, CA

SoluteScale

BrCIBr-

aL

aT

9000

43

3.9

25000

96

1.8

131r

1710

60

Br-

3830

400

em

0.30

700

CI-

1900

200-300

0.233

84

60 Co

<1150050-150

1300
3320

79
194

11

9000

50

5

1600

100

10

270
200

Glaciofluvial sand

0.33

9

Sand, gravel, silt

0.25

2700

CIr131r
BrCICI-aL, are often 10 to 100 times greater than aT (Crane & Gardner, 1961; Harleman& Rumer, 1963), albeit this difference can become much smaller for certain field

experiments.
Gelhar et al. (1992) reviewed published field-scale dispersivities determined in aquifer materials. Their analysis was based on publications documenting a variety of experimental conditions, including different types of solute,application methods, monitoring systems, observation scales, and flow configurations. Table 9-3 is a summarized version of Table 1 of Gelhar et al. (1992) showing values for the macroscopic longitudinal (aL) and transverse (a r ) dispersivity.
Table 9-3 also lists the reference and site name, aquifer material, effective porosity (f), the (horizontal) pore-water velocity (v), type of solute, and the scale of theexperiments. Further information on the table can be found in Gelhar et al.
(1992) and in the original publications. Information on laboratory-scale dispersion is relatively abundant and can be readily found in the soil science literature(cf. Kutflek & Nielsen, 1994).

III. MASS BALANCE FORMULATION
The mass balance equation, Eq. [I], was formulated in a general manner byconsidering the accumulation or depletion of solutes in a control volume overtime as the net result of a solute inflow and outflow, solute reactions, and injection or extraction of solute along with the solvent. If the last two terms of Eq. [1]are ignored and we assume solute transport by convection and dispersion, Eq. [I]PRINCIPLES OF SOLUTE TRANSPORT

345

can be written as the classical convection-dispersion equation (CDE), also sometimes referred to as the advection-dispersion equation, i.e.,

aec = -V ·(qC-eD· Vc)

-

at

[31 ]

If we apply the water conservation equation, Eq. [31] can be expressed asaC

e - = -qVC+ V ·(eD· Vc)

at

[32]

where D is the hydrodynamic dispersion tensor (L 2r'). Other forms of the CDE
also are possible. The CDE may be classified as a second-order parabolic differential equation, it is a special form of the Fokker-Planck equation (Risken, 1984).
In case of a time and space invariant water content, Eq. [31] transforms into

ac = -V ·(vC-D· Vc)
at

-

[33]

with v = q/8.
In Cartesian coordinates, the governing equation for three-dimensionaltransport during one-dimensional steady water flow in an isotropic soil profile is[34]
where x, y, and z are positions along the coordinates (L), and v denotes the velocity in the longitudinal x direction. Ellsworth and Butters (1993) showed how thetransport problem can be solved for arbitrary Cartesian coordinates (i.e. , thevelocity is not parallel to the flow direction) . The coefficients of hydrodynamic
dispersion may be expressed in a similar manner as the mechanical dispersion
coefficients (cf. Eq. [12])
[35]

[36]
where the coefficient n is usually close to unity (Scheidegger, 1974). Dispersivities are ordinarily considered to be scalars, although evidence suggests that for
many porous media the dispersivity is actually a second-rank tensor (Brenner,1980).
The previous formulation of the CDE is based on a Eulerian approach, i.e.,
solute transport is observed from a fixed coordinate position. In the Lagrangian
approach solute transport is quantified with respect to a moving coordinate system, which is attached to a designated particle or concentration moving throughthe medium. A special example is the use of the material or substantial derivative
based on the average pore-water velocity

LEU & VAN GENUCHTEN

346

DC dC
dC
dC
-=-+v·V'C=-+vDt
dt
dt
dx

[37]

Using this total derivative in a CDE such as Eq. [34] amounts to formulating the mass balance equation with respect to a coordinate system moving at the
speed of the (mean) pore-water velocity.
This approach is attractive for analytical studies since the convectiondispersion problem can now be transformed into an equivalent diffusion problem(see Leij et aI., 1998a, Chapter 10, Part II.A). Problems of this kind, including
those for heat conduction, have been studied extensively (Carslaw & Jaeger,
1959; Luikow, 1968; Crank, 1975). Mixed Eulerian-Lagrangian formulations
have recently become quite popular for numerically solving transport problems(e.g., Galeati et aI., 1992; Huang et aI., 1994). Such formulations offer several
advantages in numerical studies, for instance to avoid numerical smearing andundesired oscillations when convective transport dominates dispersion.

IV. REACTIVE SOLUTES
A variety of solute source or sink terms may be implemented in the CDE.
The most common source/sink term is due to adsorption/desorption and ion
exchange stemming from chemical and physical interactions between the soluteand the soil solid. Many other processes such as radioactive decay, aerobic andanaerobic transformations, volatilization, photolysis, precipitation/dissolution,
reduction/oxidation, and complexation also may affect the solute concentration.
Such processes may all be included somehow in the term R, of Eq. [I]. In practice, Rs is usually approximated by a combination of first- and zero-order rate
expressions. A further refinement of the transport model is necessary in the caseof nonuniform interactions between the solute and the soil, or if there is adsorption on moving particles and colloids. The discussion here will be devoted
mainly to interactions at the solid-liquid interface. These interactions have a profound effect on the movement of reactive solutes and should therefore be considered in solute transport models.
Dissolved substances in the liquid phase can interact with several soil constituents such as primary minerals, oxides, and inorganic and organic colloids.
Reviews of the mechanisms and characteristics of reactions in solid-solutionsolute systems are given by Sposito (1981,1984), and Bolt (l982a), among others. Dissolved ions in the soil solution counterbalance the surface charge of soil
particles caused by isomorphous substitution of one element for another in the
crystal lattice of clay minerals, by the presence of hydronium or hydroxyl ions at
the solid surface, and other mechanisms. The net surface charge of an assemblageof soil particles produces an electric field that affects the distribution of cations
and anions within water films surrounding the soil particles.
Several models have been proposed to describe the ionic distribution in thesoil solution near such a charged solid. One of them is the diffuse double-layer(DDL) model, which relates the electric strength of a planar solid of uniform
charge density with the solute concentration using the Poisson-Boltzmann equa-

347

PRINCIPLES OF SOLUTE TRANSPORT

tion. Figure 9-4 shows distributions of monovalent ions in the soil solution as afunction of distance from a particle surface according to the DOL theory. For
concentrated solutions with CII equal to 0.1 N or mol L- 1, the impact of the electrical field is not evident at distances greater than about 5 nm , while for dilute
solutions (0.00 I N) the electrical field extends further than 20 nm from the surface. The thickness of the electrical double layer that neutralizes the excess
charge of soil particles is not only affected by the total electrolyte concentration,
but also by the mineralogical composition of soil particles, and the valency andhydration of ions in the soi l solution. For example, the extent of the double layer
as shown in Fig. 9-4 tends to be less for divalent than monovalent ions, and lessfor larger than for smaller, generally more hydrated, ions. Moreover, changes inwater content may alter the thickness of the water films on the soil particle surface and affect the extent of the double layer.
The term adsorption is used here to refer to the net accumulation of solute,the adsorbate, at or near the fluid-solid interface, the adsorbent. Adsorption
resulting from a charge distribution may be viewed as a physical redistribution
process rather than a conventional chemical reaction, although other types ofreactions (e.g., chemisorption) also contribute to adsorption. Some ion species
are adsorbed preferentially from the liquid phase by the solid, whereas othersorbed species are more readily exchanged for different dissolved species. Since
most natural soil-water systems contain many solute species, a large number ofdifferent solute interactions may need to be considered. Precipitation of the

~

3

2
1.0
.~

>-

I-I

.2

e!

0.1

cr

.O~

::!:
0
Z

I

Z
0

.02

C."O.OI

0 .01

l-

e!

cr

I-

Z

W

U

Z
0

0 .001

U

20

40

60

80

100

DISTANCE (AI FROM PARTICLE

120

180

SURFACE

Fig. 9-4. Distribution of monovalent cations and anions near Ihe surface of a montmorillonite parti cle (modified after Nielsen et al.. 1972).

LEU & VAN GENUCHTEN

348

solute at the fluid-solid interface is not considered adsorption, although a clear
distinction between the two processes can often not be made.
The CDE for the one-dimensional transport of an adsorbed solute during
steady flow may be written as[38]

where S is the adsorbed concentration, expressed as mass of solute per mass ofdry soil (MM- '). The above equation can be expressed in terms of one dependent
variable by assuming a suitable relationship between the adsorbed and liquid
concentrations. This is typically done with a simple adsorption isotherm toexpress the sorbed concentration as a function of the liquid concentration at aconstant temperature. Besides temperature, the adsorption isotherm also generally is affected by the solution composition, total concentration, the pH of thebulk solution, and sometimes the method used for measuring the isotherm. Alternative methods to quantify solute adsorption are mass action models or surface
complexation models (e.g., Weber et aI., 1991).
A. Linear Adsorption

From a mathematical point of view it is important to distinguish between
linear and nonlinear isotherms since nonlinear isotherms generally do not permit
the derivation of relatively simple closed-form analytical solutions of the CDE.
We will consider linear adsorption here, while Subpart B briefly deals with nonlinear exchange.
A further distinction can be made between equilibrium isotherms, where achange in the solution phase induces an instantaneous change in the sorbed phase(or vice versa), and nonequilibrium isotherms, where a change in C is accompanied by a "delayed" change in S. The adsorption rate can be described assuming
first-order kinetics~~ = kf(C,S)[39]

where k is a rate parameter (r') and f is a function to quantify how far theadsorption or desorption process is removed from equilibrium. Amacher (1991)
reviewed methods of obtaining and analyzing kinetic data for soil systems,
including batch and flow methods for adsorption. A single-valued isotherm forequilibrium adsorption, C), is used to define f( C,S) according to

rc

fCC,S) = reC) - S[40]

For equilibrium adsorption k ~and hence fiC,S) ~ 0, which implies
that S = c). For a linear adsorption isotherm, the relation between rand C cansimply be given as

rc

00

[41]

349

PRINCIPLES OF SOLUTE TRANSPORT

where Kd is the partition coefficient, often referred to as the distribution coefficient, expressed in volume of solvent per mass of soil (L 3M-I). For r = S, substitution of Eq. [41] into Eq. [38] leads to the common expression for the CDE todescribe transport of a solute that undergoes linear equilibrium exchange

l42]
in which the retardation factor R is given by[43]

As a result of adsorption the convective and dispersive fluxes are reduced by afactor R, and the solute is said to be retarded with respect to the average porewater velocity. If there is no interaction between the solute and the soil (Kd = 0)
the value for R is equal to unity.
The effects of linear adsorption on solute transport in a homogeneous soil
profile are shown in Fig. 9-5 where solution and adsorbed concentrations areplotted 4 d after the start of a J-d application of influent with a unit solute concentration to an initially solute-free soil profile subject to steady saturated water3
3
7
flow. Other parameters are q = 10 cmld, 8 = 0.40 Cill' Ieill' , and D = 62.5 cm-/d.
The pore-water velocity, v =q18, is hence 25 cmld and (XL = 2.5 cm. Distributions
are given for three values of the retardation factor, R. Notice from Fig. 9-5a thatwhen R is increased from 1.0 to 2.0, the apparent solute velocity (vIR) is reduced
by one-half, causing a shallower penetration of the solute pulse into the profile.

R=2

E

~
N

.c

a.
Q)

R=l

Q

'0 100

100

120

120

C/)

,.R=l

(0)

140
160

.1

.2

.3

.4

.5

Solution Concentration, C

( b)

140

.6

160

.05

.10

.15

.20

.25

.30

Adsorbed Concentration, S

Fig. 9-5. Effect of adsorption. as accounted for by the retardation factor R. on solution (C) andadsorbed (5) concentration distributions in a homogeneous soil profile.

350

LEIJ & VAN GENUCHTEN

At the same time the area under the curve in Fig. 9-5a also is reduced by onehalf. When R = 4, the apparent solute velocity and the area under the curve areagain reduced by half. Figure 9-5b shows similar distributions for the adsorbed
concentration, S, which increases from zero (no adsorption) when R = 1, to a
maximum when R = 4. Assuming a soil bulk density, Ph, of 1.25 gJcm 3 and thesame water content as before (8 = 0.40 cm 3Jcm\ one may calculate, using Eq.
[43] that the distribution coefficient Kd equals 0, 0.32, and 0.96 Llkg, respectively, for R = 1, 2, and 4, respectively.
In some cases part of the liquid phase, especially near the solid, does notparticipate in the transport process. This occurs when the chemical is subject toanion exclusion, or when relatively immobile liquid regions are present in thesoil, for example inside aggregates (see Leij et ai., 1998, Chapter I1). Anion
exclusion results from the repulsion of anions from negatively charged surfacesof clays and ionizable organic matter; anions are repelled from such surfaces andaccumulate in the center of pores (cf. Fig. 9-4). Because water flow velocities are
zero at pore walls and maximum in the center of pores (Fig. 9-1 a), the averageanion movement may be faster than the average water movement. Many displacement experiments also suggest faster anion than water movement simply
because the displacement volume pertaining to anions is smaller for anions thanwater. The quantity (I - R) is the relative anion exclusion volume. The exclusion
volume per unit mass of soil also can be estimated according to (cf. Sposito,
1989)
[44]
where Vex is the exclusion volume (L3M- I ), Ci(x) is the concentration of an anion
i at position x (ML- 3), Co its bulk concentration (ML-3), and V is the entire vol-ume encompassing the liquid phase. Instead of using R < 1, anion transport may
be modeled with a transport model, with R = 1, which restricts the accessible liquid volume (cf. Krupp et ai., 1972). Further discussions on anion exclusion canbe found in Bresler et al. (1982), Bolt and de Haan (1982), and Corapcioglu andLingam (1994).
Just as cations, anions also may be adsorbed by the soil. This occurs
through surface complexation and adsorption onto positively charged areas of thesolid matrix. If the effect of adsorption exceeds exclusion, the anion will moveslower than water. Frequently, one value for the retardation factor (R > 0) is used
as an effective parameter in equations like Eq. [42] to model transport of a solute
species, which may undergo a variety of adsorption and exclusion processes.
Figure 9-6 shows breakthrough curves typical of the transport of an
excluded cr, a nonreactive solute (tritiated water, 3Hz O), and an adsorbed Ca 2+.
The first two tracers pertain to transport through 30-cm columns containing disturbed Glendale clay loam soil (fine silty mixed thermic Typic Torrifluvent) (P. J.
Wierenga, personal communication; see also van Genuchten & Cleary, 1982),
while the Ca 2+ data are for transport through a 30-cm long column containing a
Troup loam (loamy siliceous thermic Groassarenic Paleudult) and a Savannahfine loam (fine-loamy, siliceous, thermic Typic Frafiudult) (Leij & Dane, 1989).
Analysis of the three breakthrough curves in terms of the CDE, using inverse

PRINCIPLES OF SOLUTE TRANSPORT

351

o.o~~+---~~~--~--~--~--~--~--~--~

0.0

2.0

4.0

6.0

8.0

10.0

Pore Volumes, T
Fig. 9-6. Observed and fitted CDE breakthrough curves for three tracers typifying the transport ofanions (Cl-), a (nearly) nonreactive solute (.1 H20 ), and an adsorbing solute (Ca 2+).

procedures (Parker & van Genuchten, 1984), yielded R values of 0.681, 1.027,
and 4.120 for Cl-, .lH20, and Ca 2+. Hence, the cr curve was strongly affected byanion exclusion, while 3H20 transport was subject to relatively minor adsorption/exchange.

B. Nonlinear Adsorption

In many cases the adsorption and the retardation factor cannot be describedusing a simple "K" approach." For nonlinear equilibrium adsorption R is given asR( C) = 1+ Ph dr
e dC

[45]

Two common nonlinear adsorption isotherms are the Langmuir and Freundlich
equation[46]

r = k 1C" Freundlich[47]

where kJ, k2' k3, and n are empirical constants. Many other equations for adsorption exist (e.g., Travis & Etnier, 1981; van Genuchten & Cleary, 1982), including

352

LEU & VAN GENUCHTEN

cases where apparent hysteresis has been observed between the adsorption anddesorption isotherms of mostly organic chemicals (Hornsby & Davidson, 1973;
van Genuchten et ai., 1974).
We use the Freundlich isotherm here to illustrate the effects of nonlinear
equilibrium adsorption on solute transport. In order to keep the calculations simple, the value of k3 in Eq. [47] is taken to be 0.64, while three different values ofthe exponent n are used: 0.5, 1.0, and 1.5 to demonstrate "favorable," linear, and"unfavorable" adsorption (Fig. 9-7). Figure 9-8 showed calculated distributionsof the solution (C) and adsorbed (5) concentrations vs. soil depth 8 d after application of a 4-d long solute pulse to the soil surface. The same pore-water velocity
is used as in the example illustrated in Fig. 9-5, but a smaller dispersion coefficient is used, i.e., D = 25 cm 2/d (aL = 1 cm). Notice that, as before in Fig. 9-5, thesolution concentration distribution for n = 1 (linear adsorption) has a nearly symmetrical shape vs. depth. This is not the case for the other two n values.
When n = 0.5, a very sharp concentration front develops, while the curvenear the soil surface becomes more dispersed. The sharp front can be explainedby considering the retardation factor R for nonlinear adsorption (Eq. [45]), which
for n = 0.5, Ph = 1.25 g/cm 3, e = 0.40, and k3 = 0.64 yields R = 1 + 11K Thisshows that R increases rapidly when C decreases, becoming infinite when C
approaches zero. Consequently, the apparent solute velocity Va = viR is very small
at the lower liquid concentrations, but increases at higher values. Of course,
higher concentrations cannot move faster than lower concentrations; front sharpening will lead to a steep solute front. This front never quite becomes a step function because the large concentration gradient across the front will create a large
diffusion/dispersion flux. When n is much smaller than one, an estimate of thefront can be obtained from the average slope of the isotherm between the initialconcentration, C;, and the maximum concentration, Co. Because in the present

IC""'0>"'
EC/).6

k = 0.64.5°O~--.~2---47---~B--~.8~~1.0
Solution Concentration C (g mOl

Fig. 9-7. Freundlich equilibrium plots for kJ = 0.64 and three values of the exponent 11.

PRINCIPLES OF SOLUTE TRANSPORT

353

0
20

\

\

\

,,

,
...........

,,

(b)

,,

,

_----..:: , , ,

,'"
, .

~.

40

/'

60

~ 80

o

a-

100

/

,I

120
140

/

/

I

!

,/

/

. ,,

/'

\

./

----/-- -- --/

/

/

,

J

80

120

k = 0.64

1 = 8 days

n=0.5

10=4 days

140

n = 1.0n = 1.5

160O~-.-'-2-----L:4----J.6'---..l-8--'I.O 160 O!::--.~2--':.4-~.6:------!.8
S (mg kg')
Fig. 9-8. Effect of nonlinear adsorption on solution (C) and adsorbed (S) concentrations in a deep
homogeneous soil profile during steady-state flow. The distributions were obtained for the threeisotherms shown in Fig. 9-7.

°

example Ci = and Co "" I, the average slope of the isotherm is exactly the sameas the linear distribution coefficient, i.e., dr/dC = 0.64. Substituting this value
into Eq. [45] yields R = 3. Hence, the apparent solute velocity Va equals 25/3 or8.33 cm/d, and the solute front after 8 d is located at a depth of about 67 cm (cf.
Fig. 9-8).
A reverse scenario occurs when n > I (i.e., unfavorable exchange). Adsorption at the lower concentrations is now relatively small and the toe of the front
moves through the profile at a velocity nearly equal to that of an inert solute.
Adsorption at the higher concentrations, on the other hand, is now much more
extensive, resulting in a lower apparent solute velocity in the higher range of concentrations. As a result, the concentration front now becomes increasingly dispersed in time (Fig. 9-8). Ignoring dispersion, the velocity of the solute front, V"=viR, at any given value of C is given by
V

(/

V

=---=

1+3..JC[481

while the depth of the solute front can be approximated by[491

where z( C,G) is the location of a solute concentration with value Cat t = 0.

LEU & VAN GENUCHTEN

354

The above discussion mostly pertained to adsorption of a single ion
species. Of course, the typical cation exchange process in transport studies
involves at least two species. The simplest case arises when two cations of thesame valency and total concentration are considered, e.g., Ca2+ and Mg2+. The
resulting exchange process is then approximately linear for relatively smallchanges in the composition of the soil solution. Exchange between Na+ and Ca2+,
on the other hand, is considerably more nonlinear.
The exchange of two cationic species A and B with valency x and y may be
represented as

Al/xx+ + ~ BY+ ~ Ii Y+ + ~ Bx+
~ l/v
Y

.

x

[50]

where the bar denotes the adsorbed species. The Gapon equation (Gapon, 1933)
can be written as (cf. White & Zelazny, 1986)[51 ]where adsorbed concentrations are expressed in equivalents (or mole) whilemolar values are used for the solution activities. The Gapon equation is often successfully employed to describe the Na and Ca exchange (cf. Amrhein & Suarez,
1991), although its use has been questioned when the complex becomes saturated
with Na (Bolt, 1967).
Many other exchange equations have been reported such as the Kerr,
Vanselow, and Eriksson equations. Exchange according to the Gapon equation
also can be qualitatively described with a Freundlich or Langmuir isotherm. For
example, a Freundlich isotherm with n > 1 can describe Na+ adsorption while
Ca2+ adsorption is consistent with n < 1 (see the example in Leij et aI., 1998,
Chapter 11, Section III.B). Addition of sodium-rich irrigation water to a soil thatis initially high in Ca leads therefore to a highly skewed Na distribution like thecurve for n = 1.5 shown in Fig. 9-8, while the addition of Ca to an initially sodic
soil leads to very sharp Ca concentration fronts. These processes of increasing
and decreasing soil sodicity provide an example of nonfavorable and favorable
exchange, respectively (Bolt, 1978).
A more formal treatment of ion exchange requires knowledge of theexchange isotherm which provides the charge fraction of adsorbed ions as a function of the charge fraction of ions in the liquid phase for a particular species.
Cation exchange for Species A and B with valence ZA and ZB, respectively, can be
represented as (Bolt, 1982a)
1-

I

1

1-

ZA

ZB

ZA

ZB

-A+-B~-A+-B[52]

The thermodynamic equilibrium constant for this reaction can be given by[53]

PRINCIPLES OF SOLUTE TRANSPORT

355

with a as the activity of the adsorbed and dissolved species relative to a chosen
standard state (cf. Gaines & Thomas, 1953). The solute concentration in theadsorbed phase is often expressed in equivalent or charge fractions[54]
where n signifies the number of moles of charge in the adsorbed phase per massof dry soil, and N is the equivalent or charge fraction. The sum of the adsorbed
concentrations is sometimes referred to as the cation exchange capacity, CEC,
expressed in centimoles of charge per kilogram, which corresponds to the traditional milliequivalent per 100 g. At equilibrium the equivalent sorbed fraction is
written as[55J

Since activity coefficients for the sorbed phase are impractical, Eq. [531 may be
written as/lc

8Kex = KN (YCA

wherefi =ax /~, h =afj/NB

[56)

·A

where the equilibrium constant consists of a part with activity coefficients for thesorbed phase, fA and iii, and a part with an observable selectivity coefficient
expressed as[57J
For pragmatic reasons most of the soil chemistry research involving cationexchange has been directed toward the selectivity coefficient, K N , rather than thethermodynamic equilibrium constant, KeX" Bruggenwert and Kamphorst (1982)
provided selectivity coefficients and thermodynamic constants for binary
exchange in many clay and soil systems.
Selectivity coefficients may be convenient to formulate the solute retardation factor, R. However, both selectivity coefficients and retardation factors may
only be constant for a limited range of concentrations. If exchange data are available, the value for R can be determined by differentiating the experimental
isotherm (cf. Eq. [45]) as for the example illustrated in Fig. 9-8. Typically, an
exchange isotherm provides the dimensionless solute fraction in the sorbed
phase, NA , vs. the equivalent fraction in the liquid phase, XA- For a concave
isotherm (cf. n = 0.5 in Fig. 9-7) the exchanger prefers the "incoming" ion B over
the resident ion A, whereas for a convex isotherm the exchange complex prefers
the competing "resident" ion. These cases were already referred to as favorable(;iNAldX/ < 0) and unfavorable exchange (d 2NAldX/ > 0). In most soils theexchange does not behave according to Fig. 9-7, some exchange sites may havea high preference for a particular species whereas the remaining sites favor thecompeting species. A somewhat broader classification is based on the ratio NAIXA

LEU & VAN GENUCHTEN

356

instead of the second derivative, for favorable exchange NAIXA > 1 and for unfavorable exchange NAIXA < 1.
The effect of nonlinear exchange on the concentration profile vs. distanceor time was discussed by Bolt (1982b), among others (cf. Fig. 9-8). Front sharpening will occur for a favorably exchanged incoming solute; ion exchange counteracts spreading due to hydrodynamic dispersion. The solute displacement canbe in the form of traveling waves, i.e., the movement and the shape of the front
do not change with time. On the other hand, front spreading will take place for an
unfavorably exchanged ion. In both cases the nonlinearity prevents an exact analytical description. Analytical modeling of solute movement is possible with avariety of approximate techniques (e.g., van der Zee, 1990). An effective K"
approach may sometimes be useful (Valocchi, 1984) but the approach may beoverly simplistic (Reardon, 1981). Usually numerical methods need to be
employed to model transport of reactive solutes in case of nonlinear exchange.
The exchange and adsorption models described or referred to in this chapter may
be used for this purpose.

REFERENCES
Amacher, M.C. 1991. Methods of obtaining and analyzing kinetic data. p. 19-59. In D. L. Sparks and
D. L. Suarez (ed.) Rates of soil chemical processes. SSSA Spec. Publ. 27. SSSA, Madison, WI.
Amrhein, C , and D.L. Suarez. 1991. Sodium-calcium exchange with anion exclusion and weatheringcorrections. Soil Sci. Soc. Am. J. 55:698-706.
Aris, R. 1956. On the dispersion of a solute in a fluid flowing through a tube. Proc. R. Soc. London A
235:67-77.
Bear, J., and A. Verruijt. 1987. Modeling groundwater flow and pollution. Kluwer Acad. Pub!., Norwell , MA.
Bird, R.B., W.E. Stewart, and E.N. Lightfoot. 1960. Transport phenomena. John Wiley & Sons, Inc.,New York.
Bolt, G.H. 1967. Cation-exchange equations used in soil science-a review. Neth. J. Agric. Sci.
15:81-103.
Bolt, G.H. 1978. Transport and accumulation of soluble soil components. p. 126-140. In G.H. Bolt
and M. G. M. Bruggenwert (Ed.). Soil Chemistry. A. Basic Elements. Elsevier, Amsterdam.
Bolt, G. H. 1982a. Thermodynamics of cation exchange. p. 27-46. In G.H. Bolt (Ed.). Soil Chemistry. B. Physico-chemical models. Elsevier, Amsterdam.
Bolt, G. H. 1982b. Movement of solutes in soils: Principles of adsorption/exchange chromatography. p.
285-348. In G.H. Bolt (Ed.). Soil Chemistry. B. Physico-chemical models. Elsevier, Amsterdam.
Bolt, G. H., and F. A. M. de Haan. 1982. Anion exclusion in soil. p. 233-257. In G. H. Bolt (Ed.). Soil
Chemistry. B. Physico-chemical models. Elsevier, Amsterdam.
Brenner, H. 1980. Dispersion resulting from flow through spatially periodic porous media. Phil.
Trans. R. Soc. London. Ser. A297:81-133.
Brenner, H., and P. M. Adler. 1982. Dispersion resulting from flow through spatially periodic porous
media. II. Surface and intraparticle transport. Phil. Trans. R. Soc. London. A 307: 149-200.
Brenner, H. J., and L. J. Gaydos. 1977. The constrained Brownian movement of spherical particles incylindrical pores of comparable radius. J. Colloid Interface Sci. 58:312-356.
Bresler, E., B. L. McNeal, and D. L. Carter. 1982. Saline and sodic soils. Springer, Berlin.
Bruggenwert, M. G. M., and A. Kamphorst. 1982. Survey of experimental information on cationexchange in soil systems. p. 141-203. In G. H. Bolt (Ed.). Soil Chemistry. B. Physico-chemical
models. Elsevier. Amsterdam.
Carslaw, H. S., and J. C Jaeger. 1959. Conduction of heat in solids. Clarendon Press, Oxford.
Corapcioglu M. Y. and R. Lingam. 1994. The anion exclusion phenomenon in the porous media flow:
A review. p. 151-168. In M. Y. Corapcioglu (Ed.) Advances in porous media. Vol. 2. Elsevier,
Amsterdam.
Crane, F. E., and G. H. F. Gardner. 1961. Measurements of transverse dispersion in granular media. J.
Chem. Eng. 6:283-287.

PRINCIPLES OF SOLUTE TRANSPORT

357

Crank, J. 1975. The mathematics of diffusion. Clarendon Press, Oxford.
de Josselin de Jong, G. 1958. Longitudinal and transverse diffusion in granular deposits. Trans. Amer.
Geophys. Union. 39:67-74.
Dullien, F. A. L. 1992. Porous Media: Fluid transport and pore structure. 2nd ed. Academic Press, San
Diego.
Dykhuizen, R. c., and W. H. Casey. 1989. An analysis of solute diffusion in rocks. Geochimica et
Cosmochimica Acta 53:2797-2805.
Edwards , D. A., M. Shapiro, H. Brenner, and M. Shapira. 1991. Dispersion of inert solutes in spatiallyperiodic, two-dimensional model porous media. Transp. Porous Media 6:337-358.
Ellis, B. G., B. D. Knezek, and L. W. Jacobs. 1983. The movement of micro nutrients in soils. In Nelson et aI. (ed.) Chemical mobility and reactivity in soil systems. SSSA Special Publication II.
SSSA-ASA, Madison, WI.
Ellsworth, T. R., and G. L. Butters. 1993. Three-dimensional analytical solutions to the advection-dispersion equation in arbitrary Cartesian coordinates. Water Resour. Res. 29:3215-3225.
Epstein, N. 1989. On tortuosity and the tortuosity factor in flow and diffusion through porous media.
Chern. Eng. Sci. 44(3):777-780.
Fischer, H. B., E. List, R. C. Y. Koh, J. 1mberger, and N. H. Brooks. 1979. Mixing in Inland andCoastal Waters. Academic Press. Orlando, FL.
Freyberg, D. L. 1986. A natural gradient experiment on solute transport in a sand aquifer. 2. Spatial
moments and the advection and dispersion of nonreactive tracers. Water Resour. Res.
22:2031-2046.
Fried, 1. J., and M. A. Combarnous. 1971. Dispersion in porous media. In Adv. Hydrosci. 9: 169-282.
Academic, New York.
Gaines, G. L., and H. C. Thomas. 1953. Adsorption studies on clay minerals. II. A formulation of thethermodynamics of exchange adsorption. J. Chern. Phys. 21 :714-718.
Galeati, G., G. Gambolati, and S. P. Neuman. 1992. Coupled and partially coupled EulerianLagrangian model of freshwater-seawater mixing. Water Resour. Res. 28: 149-165.
Gapon, E. N. 1933. On the theory of exchange adsorption in soil. 1. Gen. Chern. U.S.S.R. 3: 144-152.
Garabedian, S. P., L. W. Gelhar, and M. A. Celia. 1988. Large-scale dispersive transport in aquifers:
Field experiments and reactive transport theory. Rep. 315. Ralph M. Parsons Lab .. MIT. Cambridge, MA.
Gast, R. G. 1966. Applicability of models to predict rates of cation movement in clays. Soil Sci. Soc.
Am. Proc. 30:48-52.
Gelhar, L. W. 1982. Analysis of two-well tracer tests with a pulse input. Rep. RHO-BW-CR-131 P,
Rockwell IntI., Richland, WA.
Gelhar, L. w., C. Wei thy, and K. R. Rehfeldt. 1992. A critical review of data on field-scale dispersion
in aquifers. Water Resouf. Res. 28: 1955-1974.
Gupta, Y. K., and R. N. Bhattacharya. 1983. A New Derivation in the Taylor-Aris Theory of Solute
Dispersion in a Capillary. Water Resour. Res. 19:945-951.
Hamaker, J. W. 1972. Diffusion and volatilization. p. 341-397. In C.A.!. Goring and 1.w. Hamaker(ed.) Organic chemicals in the soil environment. Marcel Dekker, Inc. , New York.
Harleman, D. R. F., and R. R. Rumer. 1963. Longitudinal and Lateral Dispersion in an Isotropic
Porous Medium. J. Fluid Mech. 16:385-394.
Harvie, C. E., and J. H. Weare. 1980. The prediction of mineral solubilities in natural waters: The NaK-Mg-Ca-CI-S04 -H 20 system from zero to high concentration at 25°C. Geochim. Cosmochim.
Acta 44:981-997.
Helfferich, F. 1962. Ion Exchange. McGraw-Hill. New York.
Heller, J. P. 1960. An "unmixing" demonstration. Am. J. Phys. 28:348-353.
Hillel, D. 1980. Fundamentals of Soil Physics. Avademic Press. New York, NY. 413 p.
Hornsby, A. G., and J. M. Davidson. 1973. Solution and adsorbed fluometuron concentrations distrihution in a water-saturated soil: experimental and predicted evaluation. Soil Sci. Soc. Am. Proc.
37:823-828.
Huang, K., R. Zhang, and M. Th. van Genuchten. 1994. An Eulerian-Lagrangian approach with anadaptively corrected method of characteristics to simulate variably saturated water flow. Water
Resour. Res. 30:497-507.
Huyakorn. P. S., P. F. Anderson. F. J. Molz. O. Giiven, and J. G. Melville. 1986. Simulations of twowell tracer tests in stratified aquifers at the Chalk River and the Mobile sites. Water Resour. Res.
22: 1016-1030.
Kemper. W. D. 1986. Solute DifTusivity. In A. Klute (ed.) Methods of Soil Analysis. I. Physical andMineralogical Methods. ASA-SSSA. Madison. WI.
Kemper. W. D., and J. C. van Schaik. 1966. Diffusion of salts in clay-water systems. Soil Sci. Soc.
Amer. Proc. 30:534-540.

358

LEU & VAN GENUCHTEN

Krupp, H. K., and D. E. Elrick. 1968. Miscible Displacement in an Unsaturated Glass Bead Medium.
Water Resour. Res. 4:809- 815.
Krupp, H. K. , J. W. Biggar, and D. R. Nielsen. 1972. Relative Flow Rates of Salt and Water in Soil.
Soil Sci. Soc. Amer. Proc. 36:412-417.
Kut'lek, M. , and D. R. Nielsen. 1994. Soil Hydrology. Catena, Cremlingen-Destedt, Germany.
Lau, L. K., W. 1. Kaufman, and D. K. Todd. 1957. Studies of dispersion in a radial flow system. Canal
seepage research: Dispersion phenomena in flow through porous media. Progress Rep. 3, I.E.R.
Ser. 93 , 3. Sanit. Eng. Res. Lab., Dep. of Engineering and School of Public Health, Univ. of California, Berkeley.
Leij, F. J. , and J. H. Dane. 1989. Determination of transport parameters from step and pulse displacement of cation and anions. Ala. Agric. Exp. Sta. Dept. of Agronomy and Soils Series 135.
Luckner, L., and W. Schestakow. 1991. Migration Processes in the Soil and Ground-Water Zone.
Lewis Publ., Chelsea, MI. 485 p.
Luikow, A. Y. 1968. Analytical heat diffusion theory. Academic Press, New York.
Mercado, A. 1966. Recharge and mixing tests at Yavne 20 well field . Underground Water Storage
Study Tech. Rep. 12, Publ. 611. Tahal Water Plann. for Israel, Tel Aviv.
Millington, R. 1. 1959. Gas diffusion of porous solids. Science 130: 100-102.
Millington, R. J. , and 1. P. Quirk. 1961. Permeability of porous solids. Trans. Faraday Soc.
57:1200- 1206.
Molinari , J., and P. Peaudecerf. 1977. Essais conjoints en laboratoire et sur Ie terrain en vue d ' uneapproach simplefiee de la prevision des propagations de substances miscibles dans les aquiferes
reels. Symposium on Hydrodynamic Diffusion and Dispersion in Porous Media . Int. Assoc.
Hydraul. Res., Pavis, Italy.
Nielsen, D. R. , and J. W. Biggar. 1961. Miscible Displacement in Soils. I. Experimental Information.
Soil Sci. Soc. Amer. Proc. 25: 1-5.
Nielsen, D. R. , R. D. Jackson , J. W. Carey, and D. D. Evans (eds.). 1972. Soil Water. American Society of Agronomy, Madison, WI. 175 p.
Novozamsky, I., J. Beek, and G. H. Bolt. 1978. Chemical equilibria. p. 13-42. In G.H. Bolt and
M.G.M. Bruggenwert (Ed.). Soil Chemistry. A. Basic Elements. Elsevier, Amsterdam.
Nye, P. H . 1979. Diffusion of Ions and Uncharged Solutes in Soils and Soil Clays. Adv. Agron. 31:
225- 272.
Nye, P. H., and P. B. Tinker. 1977. Solute movement in the soil-root system. University of Cali fomi a
Press, Berkeley, CA.
Olsen , S. R. , and W. D. Kemper. 1968. Movement of nutrients to plant roots. Advances in Agronomy
20:91-151.
Oscarson, D. W., H. B. Hume, N. G. Sawatsky, and S. C. H. Cheung. 1992. Diffusion of iodide incompacted bentonite. Soil Sci. Soc. Am. J. 56: 1400-1406.
Palmer, C. J. , and R. W. Blanchar. 1980. Prediction of diffusion coefficients from the e lectrical conductance of soil. Soil Sci. Soc. Am. J. 44:925- 929.
Parker, J. c., and M. Th. van Genuchten. 1984. Determining Transport Parameters from Laboratory
and Field Tracer Experiments. Virgo Agric. Exp. Stat. Bull. 84--5.
Pfannkuch, H. O. 1963. Contribution a l'Etude des Deplacements de Fluides dans un Milieu poreux.
Rev. Inst. Fr. Petr. 18:215- 270.
Pitzer, K. S. 1979. Activity coefficients in electrolyte solutions. CRC Press, Boca Raton, FL.
Rajarem, H. , and L. W. Gelhar. 1991. Three-dimensional spatial moments analysi s of the Borden
tracer test. Water Resour. Res. 27: 1239- 1251.
Reardon, E. J. 1981. Kd's-can they be used to describe ion sorption reactions in contaminant migration? Ground Water 19:279- 286.
Rhoades, J. D., P. A. C. Raats, and R. J. Prather. 1976. Effects of liquid-phase electrical conductivity,
water content, and surface conductivity on bulk soil electrical conductivity. Soil Sci. Soc. Am. J.
40:651 - 655.
Rhue, D. 1992. Coupled diffusion of exchangeable cations in soil. Soil Sci. Soc. Am. J. 56:683-689.
Risken, H. 1984. The Fokker-Planck Equation. Methods of Solution and Application. Springer- Verlag. New York.
Robin, M. J. L., R. W. Gillham, and D. W. Oscarson. 1987. Diffusion of strontium and chloride incompacted clay-based materials. Soil Sci. Soc. Am. J. 51 (5): 1102-1108.
Saffman, P. G. 1960. Dispersion due to molecular diffusion and macroscopic mixing in flow througha network of capi lIaries. 1. Fluid. Mech. 2: 194- 208.
Sauer, M. C. , P. F. Southwick, K. S. Spiegler, and M. R. J. Wyllie. 1955. Electrical conductance ofporous plugs: Ion exchange resin solution systems. Ind. Eng. Chern. 47:2187-2193.

PRINCIPLES OF SOLUTE TRANSPORT

359

Sauty, J. P. 1977. Contribution a l'identification des parameters de dispersion dans les aquifieres parinterpretation des experiences de tracage. Diss. Univ. Sci. et Med. et 1nst. Natl. Poly tech. de
Grenoble, Grenoble, France.
Scheidegger, A. E. 1974. The Physics of Flow Through Porous Media. University of Toronto Press.
Scotter, D. R., and P. A. C. Raats. 1968. Dispersion in Porous Mediums Due to Oscillating Flow.
Water Resour. Res. 6:1201-1206.
Sposito, G. 1981. The thermodynamics of soil solutions. Oxford Univ. Press, New York.
Sposito, G. 1984. The surface chemistry of soils. Oxford Univ. Press, New York.
Sposito, G. 1989. The chemistry of soils. Oxford Univ. Press, New York.
Standart, G. L., R. Taylor, and R. Krishna. 1979. The Maxwell-Stefan formulation of irreversible
thermodynamics for simultaneous heat and mass transfer. Chern. Eng. Comm. 3:277-289.
Taylor, G. I. 1953. Dispersion of Soluble Matter in Solvent Flowing through a Tube. Proc. R. Soc.
London A 219: 186-203.
Travis, C. G., and E. C. Etnier. 1981. A survey of sorption relationships for reactive solutes in soils. J.
Environ. Qual. 10:8-17.
Va10cchi , A. J. 1984. Describing the transport of ion-exchanging contaminants using an effective KJ
approach. Water Resour. Res. 20:499-503.
Valocchi, A. J., P. V. Roberts, G. A. Parks, and R. L. Street. 1981. Simulation of the transport of ionexchanging solutes using laboratory-determined chemical parameter values. Ground Water
19:600-607.
van der Zee, S. E. A. T. M. 1990. Analytical traveling wave solutions for transport with nonlinear andnonequilibrium adsorption. Water Resour. Res. 26:2563-2578.
van Genuchten, M. Th. 1988. Solute transport. In S. P. Parker (ed.), 1989 McGraw-Hili Yearbook of
Science and Technology, pp. 360-362. McGraw-Hill Book Co., New York, NY.
van Genuchten, M. Th., and R. W. Cleary. 1982. Movement of solutes in soils: Computer-simulated
and laboratory results. p. 349-386. In G.H. Bolt (Ed.). Soil Chemistry. B. Physico-chemical
models. Elsevier. Amsterdam.
van Genuchten, M. Th., J. M. Davidson, P. J. Wierenga. 1974. An evaluation of kinetic and equilibrium equations for the prediction of pesticide movement through porous media. Soil Sci. Soc.
Am. Proc. 38:29-35.
van Rees, K. C. J., E. A. Sudicky, P. S. C. Rao, and K. R. Reddy. 1991. Evaluation of laboratory techniques for measuring diffusion coefficients in sediments. Env. Sci. Tech. 25(9): 1605-1611.
Weast, R. c., M. J. Astle, and W. H. Beyer. 1984. Handbook of chemistry and physics. CRC Press.
Boca Raton. FL.
Weber, W. J., P. M. McGinley, and L. E. Katz. 1991. Sorption phenomena in subsurface systems:
Concepts, models and effects on contaminant fate and transport. Water Res. 25:499-528.
White, G. N., and L. W. Zelazny. 1986. Charge properties of soil colloids. p. 39-81. In D. L. Sparks(ed.) Soil Physical Chemistry. CRe Press, Boca Raton, FL.