Click
Here

WATER RESOURCES RESEARCH, VOL. 42, W08401, doi:10.1029/2005WR004802, 2006

for

Full
ArticleMulticomponent solute transport in soil lysimeters irrigatedwith waters of different quality
Maria C. Gonçalves,1 Jirka Šimůnek,2 Tiago B. Ramos,1 José C. Martins,1
Maria J. Neves,1 and Fernando P. Pires1
Received 12 December 2005; revised 3 April 2006; accepted 21 April 2006; published 2 August 2006.

[1] A variety of analytical and numerical models have been developed during thepast several decades to predict water and solute transfer processes between the soil surface
and the groundwater table. While many models quantifying solute transport in soils
usually consider only one solute and severely simplify various chemical interactions,others such as the geochemical module of HYDRUS-1D consider multiple solutes andtheir mutual interactions. In this study we use HYDRUS-1D to analyze water flow andsolute transport in three soil lysimeters (1.2 m2  1 m) irrigated during the summer
months with waters of different quality that were used to evaluate salinization andalkalization hazards. The soil monoliths were constructed in a Eutric Fluvisol in Alentejo,
Portugal. The electrical conductivity (EC) of irrigation water varied between 0.4 and3.2 dS m1, and the sodium adsorption ratio (SAR) varied between 1 and 6 (mmol(c)L1)0.5, while maintaining a ratio of Ca:Mg equal to 1:2. The soil monoliths were
subjected to regular rainfall and leaching during the rest of the year. Water contents andfluxes, concentrations of individual ions (Na+, Ca2+, and Mg2+), electrical conductivity ofthe soil solution, SAR, and exchangeable sodium percentage (ESP) indices weremonitored from May 2001 to September 2004 at four depths (10, 30, 50, and 70 cm) in all
three soil monoliths. Irrigation water with EC up to 1.6 dS m1 did not cause salinizationor alkalization hazards. The rainfall water leached the salts accumulated during theirrigation period down to a depth of 100 cm. Rainfall, however, did not restore the salinity
and sodicity of the soil to its original values below a depth of 60 cm for the lysimeterirrigated with water having an EC equal to 3.2 dS m1. HYDRUS-1D successfully
described field measurements of the water content (R2 = 0.60), overall salinity (R2 = 0.65),
and the concentration of individual soluble cations (R2 ranged between 0.62 and 0.78) as
well as the sodium adsorption ratio (R2 = 0.87) and the exchangeable sodium percentage(R2 = 0.76).
Citation: Gonçalves, M. C., J. Šimůnek, T. B. Ramos, J. C. Martins, M. J. Neves, and F. P. Pires (2006), Multicomponent solutetransport in soil lysimeters irrigated with waters of different quality, Water Resour. Res., 42, W08401, doi:10.1029/2005WR004802.

1. Introduction[2] Salinity is increasingly a concern in many irrigatedareas of the Alentejo region of Southern Portugal where
Mediterranean conditions prevail with high summer temperatures and scarce rainfall. Salinization (by irrigation water)
is a process whereby soluble salts from the irrigation water
accumulate in soil due to inadequate leaching, high watertables and/or high evaporation rates [Keren, 2000]. The
water quality in several reservoirs of the region at present
time varies between about 0.4 (Campilhas, Alto Sado,
Odivelas) and 1.6 dS m1 (Roxo), although higher values(to about 2.6 dS m1) have been frequently observed in
Roxo over the years. Many irrigated areas in the region will1
Department of Soil Science, Estação Agronómica Nacional, Oeiras,
Portugal.
2
Department of Environmental Sciences, University of California,
Riverside, California, USA.

Copyright 2006 by the American Geophysical Union.
0043-1397/06/2005WR004802$09.00

be integrated in the near future within the Alqueva project,such that additional 110,000 ha of agricultural land will be
irrigated from a newly constructed dam on the Guadiana
river. This river is one of the largest in the southern Iberian
Peninsula, with flow rates that exhibit large seasonal andinterannual fluctuations. Large changes in river water quality are due to changes in the flow rate, the local climate withfrequent dry years, the existence of 300,000 ha of irrigatedfields in Spain’s upstream catchment area, and the existenceof several urban agglomerations. As water is being used andreused, risks of soil salinization for irrigated soils undoubtedly will worsen. A sound irrigation policy must be established to mitigate these risks. Such policy must be based ona quantitative understanding of the subsurface movement ofwater and dissolved chemicals. Modeling of subsurface
water flow and the transport of major soluble ions inand below the root zone is essential for predicting groundwater quality, implementing better irrigation and fertilization practices, and quantifying salinization and alkalization
hazards.

W08401

1 of 17GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401

W08401

Table 1. Chemical Species Included in the Major Ion Module of HYDRUS-1D
Number of Species

Ions

Aqueous components
Complexed species

610Precipitated species

6

Sorbed (exchangeable) species
CO2-H2O species

47Ca2+, Mg2+, Na+, K+, SO42, Cl
CaCO3, CaHCO3+, CaSO4, MgCO3,MgHCO+3 , MgSO4, NaCO
3,
NaHCO3, NaSO
4 , KSO4
CaCO3, CaSO4 2H2O, MgCO3 3H2O,
Mg5(CO3)4(OH)2 4H2O, Mg2Si3O7.5(OH) 
3H2O, CaMg(CO3)2
Ca, Mg, Na, K+
PCO2, H2CO3*, CO2
3 , HCO3 , H , OH , H2O[3] The past several decades have seen considerable
progress in quantifying water flow and solute transport
processes in the unsaturated zone [van Genuchten andŠimůnek, 2004]. A large number of analytical and numericalmodels are now available to predict flow and transport
processes between the soil surface and the groundwatertable. These models are typically based on the numericalsolutions of the Richards equation for variably saturated
flow, and analytical or numerical solutions of the Fickianbased convection-dispersion equation for solute transport.
Deterministic solutions of these classical equations likely
will remain important tools for predicting water and solute
movement in the vadose zone, and for analyzing specificlaboratory or field experiments involving unsaturated waterflow and solute transport. Models of this type are also
helpful tools for extrapolating information from a limited
number of field experiments to different soil types, crop andclimatic conditions, as well as to different tillage and water
management schemes. Unfortunately, evaluation of thesemodels under field conditions is still very limited, in largepart because of their need for a variety of input data,
including soil hydraulic properties, solute transport parameters, parameters characterizing the partitioning between thesolid phase and the soil solution, and meteorological andcrop related information.
[4] While many models have been developed over theyears, most consider the transport of only one solute andseverely simplify various chemical interactions [Šimůnek
and Valocchi, 2002]. The relatively complex processes ofadsorption and cation exchange are often accounted for bymeans of empirical linear or nonlinear adsorption isotherms.
Other processes such as precipitation/dissolution and biodegradation are frequently simulated by invoking simplified
first- or zero-order rate equations. During the last two
decades several models have been developed that canconsider multiple solutes and their various interactions,such as precipitation/dissolution and competition for sorption sites. Šimůnek and Valocchi [2002] divided these
multicomponent solute transport models into two large
groups: general models and models with specific chemistry.
Models with generalized chemistry provide users considerable freedom in terms of designing their particular chemicalsystem, and thus permit a relatively broad range of applications. On the other hand, models with specific chemistry
are much more restrictive in terms of the prescribed (and
often simplified) chemical system, and thus are generally
constrained to very specific applications. However, they areoften much easier to use and computationally much more
efficient than general models.

[5] Models simulating the transport of major ions are
typical examples of models with specified chemistry. For
example, Robbins et al. [1980a, 1980b] developed chemical
precipitation-dissolution and cation exchange subroutines
using equilibrium chemistry and coupled them with a onedimensional unsaturated water flow model. Their model
was tested against experimental data from a lysimeter study.
Dudley et al. [1981] further evaluated this same model forfield conditions under cropped and uncropped conditions.
They obtained adequate predictions of total salinity but notof individual ion concentrations. The salinity model ofRobbins et al. [1980a] was later used also by Russo[1986] to theoretically investigate the leaching of a gypsiferous-sodic soil using different water qualities. Robbins’
equilibrium chemistry model additionally formed a basis for
the LEACHM numerical code of Wagenet and Hutson[1987].
[6] Šimůnek et al. [1996] and Šimůnek and Suarez [1994]
later coupled a major ion chemistry module (Table 1 lists
considered species) with one- and two-dimensional variably
saturated water flow models, respectively, while also considering solute transport, carbon dioxide transport, and heat
flow. The resulting UNSATCHEM models considered theeffects of CO2 producing microbiological activity and CO2transport in the soil on geochemical transport. Contrary tomodels based on the geochemical module of Robbins et al.
[1980a], the UNSATCHEM models can consider kinetic
chemical reactions (precipitation/dissolution of calcite anddissolution of dolomite). They can also be used to evaluate
chemical reactions for solutions having very high ionicstrengths since they evaluate activity coefficients either
using an extended version of the Debye-Hückel equation[Truesdell and Jones, 1974] for dilute to moderately saline
solutions, or by means of the Pitzer expressions [Pitzer,
1979]. The major ion chemistry and carbon dioxide modules of UNSATCHEM were recently incorporated in the
HYDRUS-1D software package [Šimůnek et al., 2005].
This model thus represents a relatively powerful tool forevaluating major ion chemistry in the subsurface, forassessing the effects of irrigation water quality on groundwater recharge, and for quantitative predictions of theamount of water and amendment required in order toreclaim sodic and saline soils to desired levels of salinity
and ESP (exchangeable sodium percentage) [Šimůnek andSuarez, 1997; Šimůnek and Valocchi, 2002]. While thewater flow and solute transport components of HYDRUS1D have been used widely in both research and management, application of the major ion chemistry module has
been much more limited. A majority of applications dealt

2 of 17GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401

W08401

Table 2. Ionic Composition of Irrigation Waters Applied to the Soil Monoliths and Their Classification According to U.S. Salinity
Laboratorya
Monolith
Waters IA
B
C
Waters II
A0
B0
C0

EC,
dS m1SAR,(mmol(c) L1)0.5

Ca2+,mmol(c) L1

Mg2+,
mmol(c) L1

Na+,
mmol(c) L1

Cl,
mmol(c) L1

USSL
Classificationb

0.3
0.8
1.6

1.0
3.0
6.0

1.00
1.28
1.93

1.00
2.56
3.86

1.00
4.16
10.21

3.00
8.00
16.00

C2-S1
C3-S1
C3-S2

0.8
1.6
3.2

1.5
3.0
6.0

1.85
3.16
5.10

3.65
6.32
10.20

2.50
6.52
16.70

8.00
16.00
32.00

C3-S1
C3-S1
C4-S2

a

Richards [1954].
C2, medium-salinity water (EC 0.25 – 0.75 dS m1); C3, high-salinity water (EC 0.75 – 2.25 dS m1); C4, very high salinity water(EC >2.25 dS m1); S1, low-sodium water [SAR 0 – 10 (mmol(c) L1)0.5]; S2, medium-sodium water [SAR 10 – 18 (mmol(c) L1)0.5].
bmainly with hypothetical examples [e.g., Šimůnek andSuarez, 1994, 1997] because of a lack of appropriate
experimental data.
[7] The objective of this study was to carry out field
experiments to quantify salinization and alkalization risksof a medium-textured Eutric Fluvisol, irrigated with waters ofdifferent qualities, in Alvalade-Sado (Alentejo) region of
Portugal. We also evaluated the effectiveness of the
HYDRUS-1D software package [Šimůnek et al., 2005] topredict water contents and fluxes, concentrations of individual ions, the overall salinity given by the electrical conductivity (EC), as well as SAR (sodium adsorption ratio) and ESP(exchangeable sodium percentage) under field conditions.

2. Materials and Methods
2.1. Field Experimental Conditions[8] Salinization and alkalization processes were studied
on soils irrigated with three waters of different quality in
Alvalade do Sado, Alentejo, Portugal. Three soil monoliths(1.2 m2  1.0 m deep) were for this purpose constructed ina Eutric Fluvisol [Food and Agricultural Organization,
1998]. The monoliths were laterally isolated with plastic
to prevent lateral water and solute fluxes and subjected toatmospheric conditions at the top and free drainage conditions at the bottom. The soil monoliths were covered byannual spontaneous vegetation (gramineous plants, leguminosae and compositae). TDR probes using waveguides fromthe Trase System (Soil Moisture Equipment Corp., Goleta,
CA) and ceramic cups were installed in each soil monolithin two replicates at depths of 10, 30, 50 and 70 cm.
[9] The monoliths were manually irrigated from 2001 to
2004 during the normal irrigation period between May and
September. Application amounts were 10 mm each, for a
total of 500 mm per year. The monoliths were exposed tonatural atmospheric conditions (rainfall and evapotranspiration) during the rest of the year. Because of variable
Mediterranean conditions, large annual differences existed
in the cumulative rainfall between September and May,
which corresponds with the soil’s main leaching cycle.
Rainfall totals in this time period ranged from 445 mm in2001 – 2002, 587 mm in 2002 – 2003, to 481 mm in 2003 –2004. In order to easily collect soil solution samples, the soil
water content was kept between 25 and 35% (pressure
heads of 500 and 40 cm, respectively) during theirrigation period. This water content interval corresponds
approximately with the soil field capacity. Irrigation watercompositions, named A, B, and C according to the monolith
to which they were applied, are described in Table 2. Watersthat were used in the first 2 years (waters I: B and C) wereobtained by adding increasing concentrations of NaCl,
CaCl2 and MgCl2 to water A available in the region, usinga ratio of 1:2 for Ca:Mg. The electrical conductivity (EC) of
the three waters was increased in 2003 and 2004, while
maintaining the same balance between the cations (waters
II: A0, B0, and C0).
[10] Soil solutions were collected from ceramic cups,
while soil water contents were monitored using TDRs atfour depths twice a week during the irrigation periods andonce a week during the remaining months. At the beginningof the experiment, at the end of each irrigation period, andafter the following winters, soil samples were collected at 5
depths (0 – 20, 20– 40, 40– 60, 60– 80, and 80– 100 cm) tomeasure the exchangeable cations, and the cation exchange
capacity (CEC).
2.2. Input Data
2.2.1. Evapotranspiration[11] The experimental field was located 10 m from a
meteorological station. Collected meteorological data wereused to determine the potential evapotranspiration rate(ET0) using the Penman-Monteith method [Allen et al.,1998]. A crop coefficient (Kc), which accounts for both
crop transpiration and soil evaporation, was used to calculate crop evapotranspiration. The crop evapotranspiration
rate (ETc), being the product of ET0 and Kc, was calculatedusing Kc values for pasture (0.4 for early growth stages,
0.85 for periods without irrigation [Allen et al., 1998] and1.0 when irrigated, since during each irrigation period thesoil water content was kept at the field capacity). The leaf
area index (LAI) (a mean value of 3.5 m2 m2) and acorresponding soil cover factor (SCF) were used to accountfor different stages of the soil cover and to divide the ETcdaily values into crop transpiration (T) and soil evaporation(E) rates.
2.2.2. Root Distribution[12] Evaporation and plant transpiration exert a major
influence on water and solute distributions in near-surface
environments. These processes concentrate salts by decreasing the amount of water in the soil, which, when combined

3 of 17GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401

Table 3. van Genuchten – Mualem Parameters for the Soil
Hydraulic Functionsa
Depth, cm3

3

qr, cm cm
qs, cm3 cm3a, cm1N‘
Ks, cm d1
R2

0 – 4848 – 85

85 – 170

0.050
0.427
0.029
1.21
4.39
18.2
0.995

0.108
0.428
0.108
1.16
5.91
99.3
0.994

0.000
0.373
0.040
1.15
6.91
21.4
0.998

a 2

R is the regression coefficient between measured and fitted soil
hydraulic property data.

with irrigation in arid regions, can lead to highly salineconditions. For this reason it is important to consider theroot distribution in the soil profile. Roots of the existing
vegetation in the top of the soil monolith were not very
deep, except for a few plants. On the basis of fieldobservations that showed a relatively uniform rooting depthof 30 cm, with only a small fraction of roots extending
below this depth, we assumed that 80% of all roots were
distributed in the first 30 cm of the soil profile and that theremaining 20% extended to a depth of 70 cm.
2.2.3. Soil Hydraulic Properties[13] Disturbed and undisturbed samples (100, 200, and4700 cm3) were collected at the beginning of the experimentfrom different soil layers to measure the unsaturated soil
hydraulic properties (i.e., the soil water retention andhydraulic conductivity functions). The soil water retention
curve was determined using suction tables with sand or
kaolin for suctions below 500 cm and a pressure plateapparatus for suctions above 1000 cm. The hydraulicconductivity curve was obtained by combining measurements on the 4700 and 200 cm3 samples. A constant head
method was used to measure the saturated hydraulic conductivity. The crust method [Bouma et al., 1983] was used
to measure the unsaturated hydraulic conductivity on larger
samples, while the hot air method [Arya et al., 1975] wasused at greater suctions to measure unsaturated hydraulicconductivities on 200 cm3 samples. The soil hydraulicfunctions were described using the van Genuchten-Mualem
equations [van Genuchten, 1980], with their parametersoptimized using the RETC code [van Genuchten et al.,1991]. Table 3 lists the van Genuchten-Mualem parametersthat describe the soil hydraulic functions for the analyzed
soils. The saturated hydraulic conductivities in Table 3 wereassumed to represent optimum hydraulic conditions; they
were reduced during the calculations using the reduction
functions of McNeal [1968] to account for the adverse
effects of the solution composition.
2.2.4. Solute Transport Parameters[14] The solute transport parameters (dispersivities) wereobtained from solute displacement experiments carried out
on undisturbed 9040 cm3 cylindrical samples with a crosssectional area of 452 cm2. A 0.05 M KCl pulse was applied
during steady state saturated flow. The sampling and preparation of the soil columns were carried out as described byMallants et al. [1994]. All experimental procedures are
explained by Gonçalves et al. [2001]. The chloride breakthrough curves were expressed using the dimensionless

W08401concentration as a function of the number of pore volumes
leached through the soil column. Dispersivity values(Table 4) were obtained with the nonlinear parameter
estimation code CXTFIT 2.1 of Toride et al. [1995] byfitting analytical solutions of the CDE to observed breakthrough data.
2.2.5. Physical and Chemical Analysis[15] The particle size distribution was obtained using thepipette method for particles with diameters less than 20 mm(clay and silt fractions) and by sieving for particles between
200 and 2000 mm (coarse sand) and between 20 and 200 mm(fine sand). The dry bulk density (rb) was obtained bydrying volumetric soil samples at 105C. Concentrations ofsoluble cations Na+, Ca2+, and Mg2+ in the soil solutions
collected with the ceramic cups and from saturation extracts
were measured using atomic absorption spectrophotometry.
The electrical conductivity and Cl concentrations were
determined potentiometrically. Exchangeable cations andcation exchange capacity (CEC) were determined with amodified Melich method [Melich, 1948] using a solution of
BaCl2+triethanolamine at pH 8.1. The pH was measured ona suspension of soil and distilled water (1:2.5 w/v soil/
solution) using a potentiometric method [Van Lierop, 1990].
Gapon selectivity coefficients [White and Zelazny, 1986]
describing the partitioning between the solid phase and thesolution were calculated from the initial soil conditions.
Table 4 presents the main physical and chemical characteristics of the soil.

3. Data Analyses[16] The experimental data were analyzed using version
3.0 of the HYDRUS-1D software package [Šimůnek et al.,Table 4. Physical and Chemical Soil Characteristics (Initial
Conditions)a
Depth, (c)0 – 48
1

Coarse sand, g kg
62
532
Fine sand, g kg1
1
296
Silt, g kg
110
Clay, g kg1
TextureSilty-loam
3
1.49
Bulk density, g cm
0.20
EC, dS m1
5.94
pH (H2O)
Soluble cations, mmol(c) L1
2+
Ca
0.753
0.613
Mg2+
0.625
Na+
1.991
Cl,b mmol(c) L1
1
Exchangeable cations, mmol(c) kg
2+
Ca
46.63
21.88
Mg2+
1.45
Na+
69.96
CEC, mmol(c) kg1
1 0.5
0.756
SAR, (mmol(c) L )ESP, %
2.00
1.93
K (Na/Ca), (mol L1)1/2
K, Mg/Ca
0.52
Dispersivity, cm
5.36

48 – 85

85 – 170

51
468
292
189
Silty-loam
1.51
0.23
6.58

61
428
282
229
Loam
1.61
0.23
6.74

0.673
0.560
0.990
2.223

0.670
0.520
1.000
2.190

58.50
27.00
2.16
87.66
1.260
2.46
1.37
0.50
0.68

62.00
28.13
2.17
92.30
1.296
2.35
1.28
0.51
12.18

a
EC, electrical conductivity; CEC, cation exchange capacity; SAR,
sodium adsorption ratio; ESP, exchangeable sodium percentage; K, Gapon
selectivity coefficient.
b
Calculated to maintain the charge balance.

4 of 17GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W084012005]. In addition to one-dimensional variably saturated
water flow, and solute and heat transport, this program also
numerically simulates the transport of carbon dioxide andmajor ions, and their mutual interactions. The program
numerically solves the Richards equation describing variably saturated water flow and various forms of the convection-dispersion equation describing the transport of heat,solutes, and carbon dioxide. The flow equation incorporates
a sink term to account for water uptake by plant roots. Waterand salinity stress response functions that reduce the potential root water uptake can be defined according to functions
proposed by Feddes et al. [1978] or van Genuchten [1987].
Although the program considers multiple analytical models
to represent the soil hydraulic properties, in this study weused only van Genuchten’s [1980] analytical model. Thismodel requires six parameters, i.e., the residual, qr, andsaturated, qs, water contents, the saturated hydraulic conductivity Ks, the tortuosity factor l, and two shape parameters a and n. Atmospheric and free drainage boundary
conditions were used at the surface and bottom of the soil
profile, respectively. Atmospheric boundary conditions require specification of daily values of precipitation, irrigation, evaporation, and transpiration. The governing
equations for all simulated processes, possible initial andboundary conditions, and other details can be obtained fromŠimůnek et al. [2005].
[17] HYDRUS-1D was recently updated with a major ion
chemistry module extracted from the UNSATCHEM model[Šimůnek and Suarez, 1994; Šimůnek et al., 1996]. The
geochemical module considers the transport of seven major
ions (six aqueous components in Table 1, and alkalinity)and various temperature-dependent equilibrium chemical
reactions between these ions, such as aqueous complexation, precipitation and dissolution of several solid phases,
and cation exchange. Although the geochemical module has
been extended to other species (e.g., boron and silica
species by Suarez and Šimůnek [1997] or Fe and As by
Decker et al. [2006]), in this study we will consider only
major ions. The aqueous complexes and mineral species
included in the geochemical module of HYDRUS-1D arecompiled in Table 1, while thermodynamic coefficients,
kinetic reaction constants, and other details are provided
in Šimůnek et al. [1996, 2005].
[18] Since during our experiment and modeling analyses
we did not observe any precipitation or dissolution ofmineral phases, we will not discuss these further. The twomost important chemical reactions in our application wereaqueous complexation and cation exchange. Equations forthe aqueous complexation reactions can be obtained using
the law of mass action. Partitioning between the solid
exchange phase and the solution phase (cation exchange
reaction) is described in HYDRUS by the Gapon exchange
equation [White and Zelazny, 1986]. As an example, for thecase of exchange of the cations Ca2+ and Na+, i.e.,
2þ
þ
2þ
2 cþ
Na þ cCa , cCa þ 2 cNa

ð1Þ

we used the following Gapon equation
KCa=Na ¼

gNa cþ
Na1=2
þ 
cNa gCa c2þ
Ca

c2þ
Ca

ð2Þ

W08401

It is further assumed that the cation exchange capacity cT (orCEC) was constant and independent of pH. Whenexchanger phase concentrations are expressed in equivalents per mass of solid, then
cT ¼

Mz
X

ci

ð3Þ

i¼1

where Mz is the number of ion exchange species.
[19] The sodium adsorption ratio (SAR) was calculated asfollows:
ð Naþ ÞSAR ¼ rﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
ðCa2þ þ Mg 2þ Þ
2

ð4Þ

where Na+, Ca2+, and Mg2+ are the cation concentrations inmmol(c) L1 (millimol of charge per liter is an SI unit
corresponding to milliequivalent per liter). We further usedthe exchangeable sodium percentage (ESP), which isdefined as the ratio of exchangeable Na to the sum of
exchangeable cations, expressed in equivalents per mass ofsolid:
þ

ESP ¼

Na
100
ðCECÞ

ð5Þ

The electrical conductivity in HYDRUS was determinedfrom individual anions and cations following the method of
McNeal et al. [1970], while the osmotic coefficient f wascalculated using the semiempirical equation of Pitzer [1973].
[20] The accumulation of monovalent cations, such as
sodium and potassium, or the use of high-quality water,
often leads to clay dispersion or swelling. These processes
can have an adverse effect on the soil hydraulic properties
including hydraulic conductivity, infiltration rates and soil
water retention as a result of swelling and clay dispersion.
These negative effects are usually explained based on the
diffuse double layer theory [Russo and Bresler, 1977;
Russo, 1988]. They become more pronounced with decreasing salt concentration and the valence of the adsorbed ions
[Shainberg and Levy, 1992]. HYDRUS accounts for theseadverse effects on the hydraulic conductivity by using
reduction functions developed by McNeal [1968] (seedescription by Šimůnek et al. [1996, 2005]):
K ðh; SAR; C0 Þ ¼ rðSAR; C0 ÞK ðhÞ

ð6Þ

where K(h) is the hydraulic conductivity function, C0 is thetotal salt concentration of the ambient solution inmmolcL1, and r is a scaling factor which represents theeffect of the solution composition (SAR and dilution) on the
final hydraulic conductivity (dimensionless), and which isrelated to SAR and salinity. The hydraulic conductivity
without the scaling factor r can be assumed to be theoptimal value under favorable chemical conditions in termsof optimal SAR and salinity.
[21] In addition to a visual check, field-measured values
were compared with results of the HYDRUS-1D predictions
using simple regression analysis (R2). Regression analyses
between measured and simulated data were carried out for

5 of 17

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

water contents, electrical conductivities of the soil solution(representing the overall salinity), for the individual ion
concentrations (Na+, Ca2+, Mg2+), and for the SAR and ESP
relationships. The analyses were carried out for a particularvariable within the overall data set (all depths combined), as
well as using data of the individual depths.

4. Results and Discussion
4.1. Volumetric Water Contents[22] The experiments started on 23 May 2001 (day 1). As
initial conditions for the HYDRUS-1D simulations we used
TDR water content readings of 0.1775, 0.2140, 0.2605,
0.2828 cm3 cm3 at 10, 30, 50, and 70 cm depths,
respectively, immediately before the start of the irrigation
period. Other input for HYDRUS-1D consisted of the soil
hydraulic parameters for the three soil horizons, water
applied to the top of three monoliths either from irrigation
or rainfall, and daily ETc values divided into separate
evaporation (E) and transpiration (T) components. Simulations were carried for a time period of four years. Although
the initial water contents varied slightly between the threemonoliths, all other parameters and properties were verymuch the same, which allowed us to assume in our initial
calculations that one single HYDRUS-1D simulation could
well characterize all three monoliths. We thus neglected any
possible effects of different irrigation waters on root water
uptake and plant growth. Figure 1 shows for the threemonoliths the mean values and standard deviations of thewater contents measured with TDRs at depths of 10, 30,
50, and 70 cm, and compares these values with results ofthe HYDRUS-1D simulation between 23 May 2001 and
30 September 2004 (i.e., 1227 days).
[23] During the four irrigation periods (i.e., from day 1 to76; 388 to 479; 737 to 820; and 1128 to 1186), TDR
measurements showed water content values near the field
capacity (0.25 to 0.35 cm3 cm3) and then gradually
decreasing values until the beginning of the rainy seasons(days 134 – 323, 480 – 713, and 861 – 1084) when they
started to increase again. Overall, simulated water contents
closely mirrored measured values at all four depths (Figure 1).
When comparing measured and simulated values, one must
realize that water contents were measured only sporadicallyand never during or immediately after rainfall events, thus
avoiding the highest water content values. On the other
hand, HYDRUS-1D calculates the entire process and predicts continuous water contents during all events, including
possible extremes.
[24] The relationship between measured and simulated
water contents at four depths in all three monoliths wasanalyzed by fitting a linear model to the data. Regression
coefficients were calculated for each layer and for all layers
combined. An overall regression coefficient (R2) of 0.60
was obtained for all n = 1180 observations combined, whilevalues for the individual layers decreased with depth from0.66 to 0.42 (Table 5). Better correspondence between
measured and simulated values probably could have beenobtained using the inverse (parameter estimation) option of
HYDRUS-1D [Šimůnek et al., 2005]. However, since ourgoal was to use as much as possible independently measured parameters, the parameter estimation option of
HYDRUS-1D was not used for our water flow simulations.

W08401[25] In subsequent simulations involving the transport ofmajor ions, we considered the possible effects of different
irrigation waters (and soil solution compositions) on hydraulic conductivities and root water uptake.
4.2. Overall Salinity[26] Since the soil used in our experiments (a medium
textured Fluvisol) had favorable hydraulic characteristics(saturated hydraulic conductivities of about 18, 100, and22 cm d1 in the 0– 48, 48– 75, and 75– 120 cm depths,
respectively; see Table 3) and because of leaching during
the winter rain seasons, the use of irrigation waters up to
1.6 dS m1 (monoliths A and B) did not lead to soil
salinization (Figure 2). Only at the end of each irrigation
cycle did the soil salinity increase to 5 – 6 dS m1 in theupper soil layers (0 – 40/50 cm). The salinity of the deeper
horizons increased during and after the rainfall periods,
often reaching 3 dS m1, due to leaching of salts from theupper layers.
[27] The application of irrigation water C0 (EC = 3.2 dS m1)during the last two irrigation periods had the biggest effecton the soil solution salinity, leading to a value of about
12 dS m1 in the 0 – 20 cm layer. During the subsequent
rainfall period (between 864 to 1084 days), salinity values
increased to about 10 dS m1 in layers below the 50-cm
depth due to leaching of salts from upper layers. SuchEC values can lead to significant yield reductions for a large
number of crops, especially vegetables [e.g., Ayers andWestcot, 1985; Mass, 1990; Steppuhn et al., 2005].
[28] HYDRUS-1D simulations of the total soil salinity
resulted in a generally good agreement with the observed
patterns in all three soil monoliths. Figure 2 shows observed
and simulated total soil salinities for four layers of monoliths A, B and C. We often lacked measured data between
the end of the irrigation period and the beginning of therainfall period. During this time period, which lasted about
2 months, soil water contents decreased significantly below
field capacity due to relatively high air temperatures and soil
evaporation rates, which prevented us from extracting soil
solution samples with ceramic cups. Although the measured
data hence could not validate the simulated EC peaks, thecalculated values are reasonable and to be expected during
periods of decreasing or low water contents caused byevapotranspiration. A value of 0.65 for R2 (for n = 1039
observations) was obtained for the relation between measured and simulated EC values (Table 5). Table 5 also
provides regression coefficients calculated for each depth
from the observed and simulated data.
4.3. Individual Cations[29] Measured and simulated concentrations of soluble
Na+, Ca2+ and Mg2+ for the three monoliths at four depths(10, 30, 50, and 70 cm) are presented in Figures 3, 4, and 5,
respectively. The general behavior of sodium is the same asthat of the electrical conductivity in Figure 2. Na+ concentrations increased gradually at all depths during the irrigation cycles, and peaked during the following dry periods
before the winter rainfalls. The highest concentrations werereached in monolith C with the application of the lowestquality water (C0). The peak Na+ concentration was about
40 mmol(c) L1, which is more than double the concentration of the irrigation water. Sodium concentrations subsequently decreased during the rainfall periods in the two

6 of 17

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401Figure 1. Measured (TDR) and simulated (Hydrus) volumetric water contents at 10, 30, 50, and 70 cm
depths. I and R denote the irrigation and rainfall periods, respectively.
surface layers (10 and 30 cm), reaching values similar to theinitial conditions. The deeper depths (50 and 70 cm) showed
an increase in sodium content due to the leaching from theupper layers. In addition to the sudden increases in Na+
concentrations during the winter periods, a gradual increasein Na occurred in deeper layers during the entire experiment.
[30] HYDRUS-1D simulations of the sodium concentrations resulted in a generally very good agreement withmeasured values in all three monoliths. A value of 0.78
for R2 was obtained between measured and simulated
sodium concentrations for the 1080 observation data points(Table 5).
[31] Figures 4 and 5 present results for calcium andmagnesium concentrations, respectively. Once again, thegeneral dynamics of the Ca and Mg concentrations during
the experiments was similar to those for sodium and EC.

7 of 17GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401

W08401

Table 5. Results of Regression Analyses Between Measured and Simulated Soil Water Contents, Soluble Na+, Mg2+, and Ca2+
Concentrations, Electrical Conductivities of the Soil Solution, Sodium Adsorption Ratios, and Exchangeable Sodium Percentagesa
Lysimeter Layers
10
Parameters

R

2

Water content
Na+
Ca2+
Mg2+
EC
SAR
ESP

0.66
0.84
0.73
0.79
0.87
0.88
0.81

30
n

2

R

259
277
271
293
274
264
24

0.69
0.71
0.51
0.53
0.53
0.84
0.84

50
n

R

2

310
269
249
256
261
292
24

0.53
0.82
0.54
0.46
0.60
0.87
0.77

70
n

R

2

286
272
255
256
254
312
24

0.42
0.78
0.76
0.68
0.70
0.74
0.87

90
n

R

2

325
262
240
250
250
311
24

0.49

Overall
n

2

R

n

24

0.60
0.78
0.63
0.62
0.65
0.87
0.76

1180
1080
1015
1055
1039
1180
120

a
Regression coefficients R2 were calculated both for the different layers and for all layers combined; n is the number of data points used in the regression
analyses.

However, the correspondence between measured and simulated values for these two cations was somewhat less thanfor sodium. The calculated R2 value was equal to 0.63 forcalcium (n = 1015) and 0.62 for magnesium (n = 1055)(Table 5). The largest differences were found at a depth of
70 cm between days 952 and 1073 for both calcium andmagnesium in all three monoliths. Smaller but noticeable
differences were also observed for sodium at the same times
and depth. This disagreement could be caused by experimental errors in the soil solution extraction process. However, a more likely explanation is that roots during this timeinterval reached this deeper soil horizon and extracted here
more water than was considered by the model. Notice thatthe calculated water contents were overpredicted during thistime period at this depth, resulting in underpredictions of
the cation concentrations. Nevertheless, visual inspection
and the obtained R2 values suggest relatively good overallcorrespondence between measured and simulated data.
4.4. Sodium Adsorption Ratio[32] Figure 6 presents results for the measured andsimulated SAR values. For this variable we obtained thelargest regression coefficient (R2 = 0.87 for n = 1180)
between measured and simulated values (Table 5). Relatively large SAR values were observed during time periods
when lysimeters were irrigated with C and C0 waters at
depths of 10 and 30 cm. The maximum SAR value, about 9(mmol(c) L1)0.5, was obtained for the first layer of monolith
C when irrigated with water C’. During the following winterseason, the SAR of this first layer decreased to about 5(mmol(c) L1)0.5 due to the leaching by rainwater. The valueof 5 (mmol(c) L1)0.5 is about 6 times larger than the initial
SAR values in the soil profile (i.e., 0.8 (mmol(c) L1)0.5).
Winter leaching, however, led to an increase in SAR in the
deeper horizons. SAR gradually increased with time at
depths of 50 and 70 cm of monolith C without any visible
fluctuations.
[33] Similar trends were observed for the other twomonoliths (A and B), although the SAR values were mostly
much lower. A gradual increase in SAR with time at the 30,
50 and 70 cm depths was observed in monolith B, especially after day 800 when the quality of the irrigation water
decreased. SAR fluctuations at 10 cm depth were similar as
in the monolith C, exhibiting larger values during theirrigation periods and lower values during the winter rainfall
seasons.

[34] SAR is an integral variable that characterizes saltaffected soils and provides information on the comparative
concentrations of Na+, Ca2+, and Mg2+ in soil solutions.
This variable takes into consideration that the adverse
effects of sodium are moderated by the presence of calcium
and magnesium ions. SAR reflects relatively complex
processes in the soil profile that, at least in our study, were
dominated by cation exchange. SAR gradually increased
with time at different depths after 800 days when the qualityof the irrigation water changed. However, the increase wasdelayed since the SAR front moves slower than the concentration fronts of nonreactive solutes, such as Cl, due tothe effects of the cation exchange.
4.5. Exchangeable Sodium Percentage[35] Soil samples were collected at 5 depths at the end ofeach irrigation period and after each winter season tomeasure the exchangeable cations and the cation exchange
capacity (CEC), and to calculate the exchangeable sodium
percentage (ESP). ESP was used to evaluate sodification
hazards.
[36] In monoliths A and B, the highest ESP values (about7%) occurred after the irrigation periods (Figure 7). Thisvalue does not imply sodification risks since ESP values
must be larger than 15% to lead to soil sodicity. However,lower values of ESP (about 10%) in fine-textured soils maylead to clay dispersion and thus affect soil structural
stability, reflected by a decrease in hydraulic conductivities
and overall deterioration of soil properties. The highest ESP
value in the surface layer of monolith C, irrigated with thepoorest quality water (C0), was 17% after the last irrigation
cycle (September 2003). ESP decreased to about 7% during
the following winter (March 2004), still much higher value
than the initial condition (Figure 7). Although ESP values
decreased in all soil monoliths after the winter rainfall
periods, there was a gradual increase in ESP when considered on an annual basis.
[37] Figure 7 shows a comparison between ESPs simulated with HYDRUS-1D and calculated from measured exchangeable sodium concentrations and the CEC at selected
times after the irrigation and winter rainfall periods. A good
agreement between measured and simulated values wasfound, particularly for monolith C. Exceptions were ESP
values from monolith B which showed relatively large
differences for the last three irrigation and rainfall periods(March 2003, September 2003, and March 2004). Still, a

8 of 17

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401Figure 2. Measured and simulated soil solution electrical conductivities in monoliths A, B, and C. I and
R correspond to the irrigation and rainfall periods, respectively.
relatively large R2 value of 0.76 (n = 120 observations) was
obtained for all observations.
4.6. Precipitation-Dissolution of Solid Phases[38] The highest concentrations of individual cations andanions were recorded in monolith C after the fourth irriga-

tion period and before the start of the rainy season. The
HYDRUS model produced during this time period thelowest values of pIAP (negative logarithm of the ion
activity product) for both calcite and gypsum. Since theminimum observed pIAP values for calcite and gypsum
were 12.6 and 30.3, respectively, while precipitation occurs

9 of 17

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

Figure 3. Measured and simulated soluble sodium concentrations in monoliths A, B, and C. I and R
correspond to the irrigation and rainfall periods, respectively.

10 of 17

W08401

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401Figure 4. Measured and simulated soluble calcium concentrations in monoliths A, B, and C. I and R
correspond to the irrigation and rainfall periods, respectively.

11 of 17

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

Figure 5. Measured and simulated soluble magnesium concentrations in monoliths A, B, and C. I and R
correspond to the irrigation and rainfall periods, respectively.

12 of 17

W08401

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401Figure 6. Measured and simulated sodium adsorption ratios (SAR) for monoliths A, B, and C. I and R
correspond to the irrigation and rainfall periods, respectively.

13 of 17

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401

optimal between pressure heads h2 and h3, whereas forpressure head between h3 and h4 (or h1 and h2), water
uptake decreases (or increases) linearly with h. In this study,
we used parameters suggested by Wesseling et al. [1991] forpasture, i.e., h1 = 10 cm, h2 = 25 cm, h3 = 200 to800 cm, h4 = 8000 cm.
[40] To evaluate the effects of the water and salinity stresson root water uptake we reran the simulation for monolith C
without considering osmotic head reductions. Cumulative
potential transpiration for the entire 4-year period was equal
to 331 cm. This total potential value was reduced due to thepressure head (water) stress to 240 cm and due to combinedpressure head and osmotic head stress to 228 cm. The
osmotic stress effects on the cumulative root water uptake(transpiration) were thus relatively minor in our lysimeters.
The effects of salinity stress on predicted Na concentrations
are shown in Figure 8. Since osmotic stress reduces root
water uptake, the model predicts larger water contents, andsmaller Na concentrations (Figure 8). Predicted Na concentrations are especially lower during the dry periods whenlower water contents lead to larger solution concentrations,and thus to larger osmotic heads and more reductions in root
water uptake.

Figure 7. Measured and simulated exchangeable sodium
percentages (ESP) for monoliths A, B, and C.

at values of 8.37 and 4.848 [Truesdell and Jones, 1974], no
precipitation occurred in our lysimeters.

5. Sensitivity Analysis
5.1. Effects of the Osmotic Stress on Root
Water Uptake[39] For the simulations reported above we assumed thatthe potential root water uptake was reduced due to bothwater stress (not sufficient supply of water) and salinity
stress (high levels of salinity). We further assumed that theeffects of the water and salinity (osmotic) stress are additive[van Genuchten, 1987]. To describe the root water uptake
reductions, we used the model of Feddes et al. [1978],which assumes that root water uptake is zero close tosaturation (for pressure heads h > h1) and below the wilting
point pressure head (h < h4). Water uptake was considered

5.2. Effects of Solution Composition on the
Hydraulic Conductivity[41] Our simulations thus far assumed that hydraulicconductivities in the soil profile are reduced due to the effectsof the solution composition as described by equation (6).
To evaluate the effects of these reductions in the hydraulicconductivity on predicted concentrations of the majorions, we reran the simulation for monolith C without
considering the hydraulic conductivity reductions. Whilein our original simulations the hydraulic conductivities were
often reduced by up to two orders of magnitude as a resultof solution composition changes, the reductions did not leadto dramatically different water content (not shown) or ion
concentration (Figure 8) profiles. While large reductions in
the hydraulic conductivity certainly will have a dramaticeffect on many short-term processes, such as especiallyponded infiltration, the effects on long-term (seasonal)
processes appear relatively minor when flux-type atmospheric boundary conditions are used (as was the case in
our study).
5.3. Effects of Model Complexity[42] The HYDRUS-1D model as used in this study withthe major ion chemistry module is a relatively complex code
which requires a large number of input parameters that are
not always readily available. This poses the obvious question whether it is necessary to use such a complex model todescribe the transport of major ions, or if one could use
simpler models and obtain similar results. Our results
indicate that neglecting osmotic stress on root water uptake
or the effects of salinity on the hydraulic conductivity didlead to very similar results. Since no solid phases precipitated during our experiments, or were added as amendments, we did not need to consider precipitation anddissolution. Could we further simplify the problem by notconsidering cation exchange and using a much simpler
linear Kd (the distribution coefficient between the liquidand solid phase) approach that assumes that particular

14 of 17

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401Figure 8. Effects of osmotic stress and reductions in the hydraulic conductivity due to changes in thesolution composition on calculated Na+ concentrations in lysimeter C.

solutes are independent of each other? To answer thisquestion we reran again the simulation for monolith C,
but now assumed that cations are transported independentlyof each other and used distribution coefficients (Kd) equal toeither zero (a nonadsorbing tracer) or values obtained byrelating the liquid and sorbed concentrations from collected

samples. The Kd coefficients obtained in this way for
different soil layers varied between 4.8 and 5.55 L kg1
for Na, between 86.5 and 147 L kg1 for Mg, and between
149 and 267 L kg1 for Ca. Sodium concentrations at a
depth of 10 cm obtained with three simulations that considered either cation exchange, no sorption (tracer-like

15 of 17

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

W08401

W08401Figure 9. Calculated Na+ concentrations in lysimeter C using the Gapon exchange equation(equation (2)), using the distribution coefficient Kd approach, and assuming no sorption (tracer).
behavior), or linear sorption are shown in Figure 9. Noticethat the three runs produced dramatically different results.
Although the solute concentrations obtained with simulations that considered either cation exchange or no sorption
showed similar patterns, absolute values of the concentrations were significantly larger for the tracer. Also, thesimulations assuming linear sorption completely failed todescribe the observed concentration patterns. Hence, while
not all options of a complex model always need to be used
for a particular data set or application, any simplification of
the model by neglecting or simplifying relevant processes
must be done with extreme caution.

6. Conclusion[43] The use of irrigation water with EC values up to
1.6 dS m1 between May 2001 to September 2004 did notlead to salinization/alkalization of the medium-textured
Fluvisol because of favorable hydraulic characteristics andwinter precipitation. Although the soil salinity at the end ofeach irrigation cycle gradually increased in the surface
layers (0 – 40/50 cm), reaching a maximum value of about
5 dS m1, EC values in the upper layers after each rainfall
period were reduced down to their initial values. Soilsalinity steadily increased below a depth of 40 cm due toleaching of salts from the upper layers, but did not exceed
EC values of 3 dS m1. SAR and ESP values also increased
after the application of irrigation water and decreased during
the rainfall leaching, but showed a gradual increase fromyear to year.
[44] Irrigation with water having an EC of 3.2 dS m1
can cause significant soil salinization and alkalization. After
two irrigation cycles (in 2003 and 2004), the EC of the soil
solution of the surface soil horizon (0 – 20 cm) increased toabout 12 dS m1, while SAR increased to about 8 (mmol(c)L1)0.5, and ESP to 17%. Although the winter rainfall
leached most salts from the surface layers, this was not
sufficient to restore the soil to its initial conditions bellow adepth of 40 cm. EC values reached 3 dS m1, SAR
7 (mmol(c) L1)0.5 and ESP 9% at a depth of 50 cm inMarch 2004, immediately after the third rainfall season. The

trend of soil salinization/alkalization below 40/50 cm must
be taken into account when using continuously poor quality
irrigation waters, especially in dry years.
[45] HYDRUS-1D successfully simulated the water regime, as well as the effects of different irrigation waters on
the geochemistry of the studied Fluvisol. The correspondence between observed and simulated variables was remarkable, considering that simulations were carried out topredict field measurements over a considerable time period(4 years) without any calibration, and with all input variables (i.e., the soil hydraulic properties, solute transportparameters, atmospheric demand, Gapon constants, physicaland chemical characteristics of the soil, LAI and root depth)
measured independently. The agreement between measured
and simulated values was best for soluble sodium concentrations (R2 = 0.78), SAR (R2 = 0.87) and ESP (R2 = 0.76),
while water contents (R2 = 0.60), ECs (R2 = 0.65), andsoluble magnesium (R2 = 0.62) and calcium concentrations(R2 = 0.63) were predicted slightly less well.
[46] Carried out sensitivity analysis showed that model
predictions for our experimental conditions were not significantly affected by the reduction of hydraulic conductivities in the soil profile due to the effects of the solution
composition or by the reduction of the root water uptake
due to osmotic stress. Our results thus indicate that it would
be possible for our experimental conditions to simplify themodel and to neglect the effect of the osmotic stress on root
water uptake or the effect of salinity on the hydraulicconductivity. On the other side, the solute concentrations
could not be well predicted without considering cationexchange. Simulations that considered either no or linear
sorption completely failed to describe the observed concentration patterns. Although it may not be always necessary toconsider all options of a complex model for a particular data
set or application, any simplification of the model byneglecting or simplifying relevant processes must be donewith extreme caution.
[47] In spite of the considerable input data demands,
HYDRUS-1D proved to be an effective tool that may
become very useful for irrigation management, and forpredicting the effects of irrigation water quality on soil

16 of 17

W08401

GONÇALVES ET AL.: SOLUTE TRANSPORT IN SOIL LYSIMETERS

and groundwater quality. This may be particularly importantwhen designing new irrigation areas. Models such as
HYDRUS-1D, after proper calibration and validation, hence
should be considered for establishing sound irrigation
policies and for mitigating associated environmental risks.

References
Allen, R. G., L. S. Pereira, D. Raes, and M. Smith (1998), Crop evapotranspiration: Guidelines for computing crop water requirements, Irrig.
Drain. Pap., 56.
Arya, L. M., D. A. Farrel, and G. R. Blake (1975), A field study of soil
water depletion patterns in presence of growing soybean roots, I. Determination of hydraulic properties of the soil, Soil Sci. Soc. Am. J., 45,1023 – 1030.
Ayers, R., and D. Westcot (1985), Water quality for agriculture, Irrig.
Drain. Pap., 29.
Bouma, J., C. Belmans, L. W. Dekker, and W. J. M. Jeurissen (1983),
Assessing the suitability of soils with macropores for subsurface liquid
waste disposal, J. Environ. Qual., 12, 305 – 311.
Decker, D. L., J. Šimůnek, S. W. Tyler, C. Papelis, and M. Logsdon (2006),
Variably saturated reactive transport of arsenic in heap leach facilities,
Vadose Zone J., 5, 430 – 444.
Dudley, L. M., R. J. Wagenet, and J. J. Jurinak (1981), Description of soil
chemistry during transient solute transport, Water Resour. Res., 17,1498 – 1504.
Feddes, R. A., P. J. Kowalik, and H. Zaradny (1978), Simulation of Field
Water Use and Crop Yield, John Wiley, Hoboken, N. J.
Food and Agricultural Organization (1998), World reference base for soil
resources, World Soil Resour. Rep., 84.
Gonçalves, M. C., F. J. Leij, and M. G. Schaap (2001), Pedotransfer functions for solute transport parameters of Portuguese soils, Eur. J. Soil Sci.,
52, 563 – 574.
Keren, R. (2000), Salinity, in Handbook of Soil Science, edited by M. E.
Summer, pp. G1 – G26, CRC Press, Boca Raton, Fla.
Mallants, D., M. Vanclooster, M. Meddahi, and J. Feyen (1994), Estimating
solute transport in undisturbed soil columns using time domain reflectrometry, J. Contam. Hydrol., 17, 91 – 109.
Mass, E. V. (1990), Crop salt tolerance, in Agricultural Salinity Assessment
and Management, Manual Eng. Pract., vol. 71, edited by K. K. Tanji,chap.13, pp. 262 – 304, Am. Soc. of Civ. Eng., Reston, Va.
McNeal, B. L. (1968), Prediction of the effect of mixed-salt solutions on
soil hydraulic conductivity, Soil Sci. Soc. Am. Proc., 32, 190 – 193.
McNeal, B. L., J. D. Oster, and J. T. Hatcher (1970), Calculation of electrical conductivity from solution composition data as an aid to in-situ
estimation of soil salinity, Soil Sci., 110, 405 – 414.
Melich, A. (1948), Determination of cations and anions exchange properties
of soils, Soil Sci., 66, 429 – 445.
Pitzer, K. S. (1973), Thermodynamics of electrolytes I: Theoretical basis
and general equations, J. Phys. Chem., 77, 268 – 277.
Pitzer, K. S. (1979), Activity coefficients in electrolyte solutions, in Activity
Coefficients in Electrolyte Solutions, edited by R. M. Pitkowicz, pp.157 –
208, CRC Press, Boca Raton, Fla.
Richards, L. A. (1954), Diagnosis and Improvement of Saline and Alkali
Soils, Agric. Handb., vol. 60, U.S. Dep. of Agric., Washington, D. C.
Robbins, C. W., R. J. Wagenet, and J. J. Jurinak (1980a), A combined salt
transport-chemical equilibrium model for calcareous and gypsiferous
soils, Soil Sci. Soc. Am. J., 44, 1191 – 1194.
Robbins, C. W., R. J. Wagenet, and J. J. Jurinak (1980b), Calculating cationexchange in a salt transport model, Soil Sci. Soc. Am. J., 44, 1195 – 1200.
Russo, D. (1986), Simulation of leaching of a gypsiferous-sodic desert soil,
Water Resour. Res., 22, 1341 – 1349.
Russo, D. (1988), Numerical analysis of the nonsteady transport of interacting solutes through unsaturated soil: I. Homogeneous systems, Water
Resour. Res., 24, 271 – 284.
Russo, D., and E. Bresler (1977), Analysis of the saturated – unsaturated ahydraulic conductivity in a mixed Na – Ca soil system, Soil Sci. Soc. Am.
J., 41, 706 – 710.
Shainberg, I., and G. J. Levy (1992), Physico-chemical effects of salts upon
infiltration and water movement in soils, in Interacting Processes in Soil

W08401

Science, edited by R. J. Wagenet, P. Baveye, and B. A. Stewart, CRC
Press, Boca Raton, Fla.
Šimůnek, J., and D. L. Suarez (1994), Two-dimensional transport model forvariably saturated porous media with major ion chemistry, Water Resour.
Res., 30, 1115 – 1133.
Šimůnek, J., and D. L. Suarez (1997), Sodic soil reclamation using multicomponent transport modeling, J. Irrig. Drain. Eng., 123(5), 367 – 376.
Šimůnek, J., and A. J. Valocchi (2002), Geochemical transport, in Methodsof Soil Analysis, part 1, Physical Methods, 3rd ed., edited by J. H. Daneand G. C. Topp, chap. 6.9, pp. 1511 – 1536, Soil Sci. Soc. of Am.,
Madison, Wis.
Šimůnek, J., D. L. Suarez, and M. Šejna (1996), The UNSATCHEM software package for simulating one-dimensional variably saturated water
flow, heat transport, carbon dioxide production and transport, and multicomponent solute transport with major ion equilibrium and kinetic chemistry, version 2.0, Res. Rep. 141, 186 pp., U. S. Salinity Lab., Agric. Res.
Serv., Riverside, Calif.
Šimůnek, J., M. T. van Genuchten, and M. Šejna (2005), The HYDRUS-1D
software package for simulating the one-dimensional movement of water,
heat, and multiple solutes in variably-saturated media, version 3.0, HYDRUS Software Ser. 1, 270 pp., Dep. of Environ. Sci., Univ. of Calif.,
Riverside.
Steppuhn, H., M. T. van Genuchten, and C. M. Grieve (2005), Root-zone
salinity: II. Indices for tolerance in agricultural crops, Crop Sci., 45,
221 – 232.
Suarez, D. L., and J. Šimůnek (1997), UNSATCHEM: Unsaturated waterand solute transport model with equilibrium and kinetic chemistry, Soil
Sci. Soc. Am. J., 61, 1633 – 1646.
Toride, N., F. J. Leij, and M. T. van Genuchten (1995), The CXTFIT codefor estimating transport parameters from laboratory or field tracer experiments, version 2.0, Res. Rep. 137, U.S. Salinity Lab., Riverside, Calif.
Truesdell, A. H., and B. F. Jones (1974), Wateq, a computer program forcalculating chemical equilibria of natural waters, J. Res. U.S. Geol. Surv.,
2(2), 233 – 248.
van Genuchten, M. T. (1980), A closed form equation for predicting thehydraulic conductivity of unsaturated soils, Soil Sci. Soc. Am. J., 44,892 – 898.
van Genuchten, M. T. (1987), A numerical model for water and solute
movement in and below the root zone, Res. Rep. 121, U.S. Salinity
Lab., Agric. Res. Serv., Riverside, Calif.
van Genuchten, M. T., and J. Šimůnek (2004), Integrated modeling ofvadose zone flow and transport processes, in Unsaturated Zone Modelling: Progress, Challenges and Applications, Frontis Ser., vol. 6, editedby R. A. Feddes, G. H. de Rooij, and J. C. van Dam, pp. 37 – 69,
Springer, New York.
van Genuchten, M. T., F. J. Leij, and S. R. Yates (1991), The RETC codefor quantifying the hydraulic functions of unsaturated soils, Rep. EPA/
600/2-91-065, U.S. Environ. Prot. Agency, Washington, D. C.
Van Lierop, W. (1990), Soil pH and lime requirement determination, in
Methods of Soil Analysis, part 2, 2nd ed., edited by A. L. Page, R. H.
Miller, and D. R. Keeney, pp. 76 – 120, Soil Sci. Soc. of Am., Madison,
Wis.
Wagenet, R. J., and J. L. Hutson (1987), LEACHM: Leaching Estimation
and Chemistry Model, a process-based model of water and solute movement, transformations, plant uptake and chemical reactions in the unsaturated zone, continuum 2, report, Water Resour. Inst., Cornell Univ.,
Ithaca, N. Y.
Wesseling, J. G., J. A. Elbers, P. Kabat, and B. J. van den Broek (1991),
SWATRE: Instructions for input, report, Winand Staring Cent., Wageningen, Netherlands.
White, N. L., and L. M. Zelazny (1986), Charge properties in soil colloids,
in Soil Physical Chemistry, edited by D. L. Sparks, pp. 39 – 81, CRC
Press, Boca Raton, Fla.



M. C. Gonçalves, J. C. Martins, M. J. Neves, F. P. Pires, and T. B.
Ramos, Department of Soil Science, Estação Agronómica Nacional,
Avenida República, P-2784-505 Oeiras, Portugal. (mc.goncalves@netc.pt)
J. Šimůnek, Department of Environmental Sciences, University ofCalifornia, Riverside, CA 92521, USA.

17 of 17