1

Simulating Nonequilibrium Movement of Water, Solutes and2Particles Using HYDRUS: A Review of Recent Applications

3
4Jiří Šimůnek1, J. Maximilian Köhne2, Radka Kodešová3, Miroslav Šejna4

5

1

6

Avenue, Riverside, CA 92521, USA, jiri.simunek@ucr.edu

7

2

8

Straße 4, D-06120 Halle (Saale), Germany

93

Department of Environmental Sciences, University of California, Riverside, 900 University

Department of Soil Physics, Helmholtz Centre for Environmental Research, Theodor-Lieser-

Department of Soil Science and Geology, Czech University of Life Sciences, Prague,

10Kamýcká 129, 16521 Prague, Czech Republic

11

4

PC-Progress, s.r.o., Anglická 28, 120 00 Prague, Czech Republic

1213

Submitted to

14
15

Soil and Water Research

16
17

Abstract

18

Water and contaminants moving through the vadose zone are often subject to a large number

19

of simultaneous physical and chemical nonequilibrium processes. Traditional modeling tools20

for describing flow and transport in soils either do not consider nonequilibrium processes at

21

all, or consider them only separately. By contrast, a wide range of nonequilibrium flow and22transport modeling approaches are currently available in the latest versions of the HYDRUS

23

software packages. The formulations range from classical models simulating uniform flow24

and transport, to relatively traditional mobile-immobile water physical and two-site chemical

25

nonequilibrium models, to more complex dual-permeability models that consider both

26

physical and chemical nonequilibrium. In this paper we briefly review recent applications of27

the HYDRUS models that used these nonequilibrium features to simulate nonequilibrium

28

water flow (water storage in immobile domains and/or preferential water flow in structured

29

soils with macropores and other preferential flow pathways), and transport of solutes

30(pesticides and other organic compounds) and particles (colloids, bacteria and viruses) in the31

vadose zone.

32

Key words: nonequilibrium flow and transport, physical nonequilibrium, chemical

33

nonequilibrium, numerical models, preferential flow, reactive transport, HYDRUS, review of34

recent applications

35
36

Introduction37

It is a pleasure to contribute to the special issue honoring Prof. Kutílek on his 80th

38birthday. Prof. Kutílek was teaching Hydropedology (a branch of science only recently

39

discovered in US, e.g., LIN et al. 2005) at the Czech Technical University in Prague, Czech

40

Republic, when I (the senior author) was a student there in the early 1980’s. His entertaining

41

and well-thought out lectures, and especially his earlier books (KUTÍLEK 1980, 1982),

42

attracted me to Soil Physics, to which I later devoted my professional career.

43

Much has changed since then. Not only our society, but also the tools that we use

44

today to address the various topics in Soil Physics, including water flow and solute transport,

45

are very different. The knowledge in Soil Physics has vastly increased. There are excellent

46

textbooks available providing an overview of the field to students and professionals alike47(e.g., KUTÍLEK & NIELSEN 1994). Moreover, while Prof. Kutílek’s lectures and books were48

primarily devoted to descriptions of various soil physical processes and their governing

49

equations, which were usually solved using relatively complex analytical and semi-analytical

50

models, today, when it is hard to imagine that IBM PC was introduced only about 25 years

2

51

ago, water flow and solute transport equations are typically solved numerically. Numerical52

methods dramatically expanded our ability to solve complex mathematical systems involving

53

a large number of simultaneous nonlinear processes that were previously unsolvable. One

54

example is nonequilibrium preferential flow and transport. Computer models based on

55

numerical solutions are now increasingly used for a wide range of applications in the research

56

and management of natural subsurface systems. The wide usage of numerical models was also

57

significantly promoted by their availability in both the public and commercial domains, and58

by the development of sophisticated graphics-based interfaces that tremendously simplify

59

their use (ŠIMŮNEK 2005).

60

The HYDRUS software packages (ŠIMŮNEK et al. 1998, 2005, 2007, 2008) are among61

the most widely used models simulating water flow and solute transport in soils. For example,

62

in March of 2007 HYDRUS-1D was downloaded more than 200 times by users from 30

63

different countries, over one thousand times in 2006, and the HYDRUS web site receives on

64

average some 700 individual visitors each day (ŠIMŮNEK et al. 2008). Traditionally, the65

HYDRUS software, as well as many other models simulating variably saturated water flow

66

and solute transport, either did not consider nonequilibrium flow and transport processes at

67

all, or only considered them separately. For example, previous versions of HYDRUS codes68

considered physical and chemical nonequilibrium separately. Physical nonequilibrium solute69

transport was accounted for by assuming a two-region, dual-porosity type formulation that70partitions the liquid phase into mobile and immobile regions (e.g., VAN GENUCHTEN &71

WIERENGA 1976). Chemical nonequilibrium solute transport was accounted for by assuming a

72

two-site sorption model, which assumes that sorption sites can be divided into two fractions

73

with sorption to different fractions of sorption sites being either instantaneous or kinetic (e.g.,

74

VAN GENUCHTEN & WAGENET, 1989).

75

Over the years, several publicly available numerical codes have been developed that3

76

consider a number of options for simulating nonequilibrium water flow and/or solute transport

77(e.g., PRUESS 1991; JARVIS 1994; VAN DAM et al. 1997; KÖHNE et al. 2008ab). Unique to the78

most recent version of the HYDRUS-1D software package (ŠIMŮNEK et al. 2005, 2008) is the79

wide range of approaches that can be selected for simulating nonequilibrium processes

80(ŠIMŮNEK & VAN GENUCHTEN 2008). The models range from classical models simulating81

uniform flow and transport, to traditional dual-porosity physical and two-site chemical

82

nonequilibrium models, to complex dual-permeability models that consider both physical and83

chemical nonequilibrium.

84

While the large number of physical and chemical nonequilibrium approaches available

85

in the latest version of HYDRUS-1D was reviewed previously by ŠIMŮNEK and VAN

86

GENUCHTEN (2008), in this paper we briefly review some recent applications of the HYDRUS

87models that used these nonequilibrium features to simulate nonequilibrium water flow and88

solute transport in the vadose zone.

8990

HYDRUS Software Packages

91

There are three different versions of the HYDRUS software packages currently in use.

92

HYDRUS-1D (ŠIMŮNEK et al. 2005), HYDRUS-2D (ŠIMŮNEK et al. 1998), and HYDRUS

93(2D/3D) (ŠIMŮNEK et al. 2006; ŠEJNA & ŠIMŮNEK 2007). While each model considers similar

94

basic processes, their main difference lies in the dimensionality of the problems they can95

address. While HYDRUS-1D considers one-dimensional problems associated with, for

96

example, soil columns, lysimeters, soil profiles and plots, HYDRUS-2D solves two-97

dimensional or axisymmetrical three-dimensional problems as encountered on the lab or field98

scale, and HYDRUS (2D/3D) calculates both two- and three-dimensional problems.

99

HYDRUS (2D/3D) is an upgrade and extension of HYDRUS-2D. In addition to basic water100

flow and solute transport processes, individual HYDRUS packages can also consider some

4

101

additional processes. For example, HYDRUS-1D can also consider the transport and102

production of carbon dioxide (ŠIMŮNEK & SUAREZ 1993) and the transport of major ions and103

major ion chemistry (ŠIMŮNEK & SUAREZ 1994; GONÇALVES et al. 2007), and HYDRUS

104(2D/3D) can simulate processes in constructed wetlands (LANGERGRABER & ŠIMŮNEK 2005,

105

2006).

106
107
108

Nonequilibrium Flow and Transport Models
Physical Nonequilibrium – Mobile-Immobile Water Model109

As already discussed, all HYDRUS codes traditionally consider physical and chemical

110

nonequilibrium separately. Physical nonequilibrium solute transport is accounted for by111

assuming a two-region, dual-porosity type formulation that partitions the liquid phase (i.e., the112

total water content θ) into mobile (flowing, inter-aggregate), θmo, and immobile (stagnant,113

intra-aggregate), θim, regions (e.g., VAN GENUCHTEN & WIERENGA 1976):114

θ = θ mo +θim

115

While water flow is assumed to be uniform in this model and the immobile water content to be

116

constant with time (i.e., there is no water transfer between two regions), the mobile-immobile

117

water concept, and thus physical nonequilibrium, is applied only to solute transport.

(1)

118
119

Chemical Nonequilibrium – Two-Site Sorption Model

120

Much like the mobile-immobile water concept (MIM), the concept of two-site sorption

121(TSS) (e.g., VAN GENUCHTEN & WAGENET 1989; model (b) in Figure 1) was already122

implemented in earlier versions of HYDRUS to permit consideration of nonequilibrium

123

adsorption-desorption reactions. The two-site sorption concept assumes that the sorption sites,124s, can be divided into two fractions:125s = se + sk
5

(2)

126

Sorption, se, on one fraction of the sites (type-1 sites) is assumed to be instantaneous, while127

sorption, sk, on the remaining (type-2) sites is assumed to be a first-order rate kinetic process.

128

Parameter f (a fraction of sorption sites at equilibrium with the solution) was used to129

discriminate between the extent of instantaneous and kinetic sorption.

130
131

Additional Chemical Nonequilibrium Models

132

While the two-site sorption model could always be simplified into a model with only133

instantaneous (f=1) or kinetic (f=0) sorption (one kinetic site model; model (a) in Figure 1),

134

there was an urgent need for a more complex sorption model for simulating transport of135

particles in porous media. The more complex, two kinetic sites model (model (c) in Figure 1)

136

was implemented into HYDRUS by SCHIJVEN and ŠIMŮNEK (2002) and ŠIMŮNEK et al.

137(2005). They also reformulated the traditional first-order sorption model∂s k= α ⎡⎣ sek - s k ⎤⎦
∂t

138
139

(3)

using the concept of attachment/detachment

ρ

140

∂s k
= kaθ c − kd ρ s k∂t141

In equations (3) and (4), sek is the sorbed concentration that would be reached at equilibrium

142

with the liquid phase concentration, sk is the sorbed concentration of the kinetic sorption sites,143

α is a first-order rate constant describing the kinetics of the sorption process, ka is the144

attachment coefficient, kd is the detachment coefficients, ρ is the bulk density, and c is the145

liquid phase concentration. ŠIMŮNEK and VAN GENUCHTEN (2008) showed that equations (3)146

and (4) are mathematically equivalent. However, the attachment/detachment model not only147

allowed for different attachment and detachment coefficients for the two sorption sites, but148

also allowed for different interpretation of those sites. That means that different processes

149

could occur at the two sorption sites.

6

(4)

150

The HYDRUS model that took into account two kinetic sites was first used by

151

SCHIJVEN and ŠIMŮNEK (2002), who simulated the removal of bacteriophages MS2 and PRD1

152

by dune recharge and removal of MS2 by deep well injection. They argued that the solid

153

phase can be divided into two fractions with different surface properties and various

154

attachment and detachment rate coefficients, reflecting the different surface characteristics of155

the soil organic carbon and ferric oxyhydroxide.

156

A different interpretation for the two sorption processes was presented by BRADFORD

157

et al. (2002, 2003, 2006ab). They used the first sorption site to represent the158

attachment/detachment process and the second sorption site represented the straining process159

of colloids. BRADFORD et al. (2003, 2006a) used a Langmuir-type blocking (a coefficient ψ160

that multiplies ka) for the attachment/detachment sites and depth-dependent blocking for161

straining sites. Both of these blocking mechanisms render equation (4) nonlinear. While162

BRADFORD et al. (2002, 2003) simulated transport of colloidal microspheres in homogeneous

163

soil columns using HYDRUS-1D, BRADFORD et al. (2004) used HYDRUS-2D to evaluate the164

effects of interfacial areas between various heterogeneities on colloid transport. This series of165

studies ultimately led to the development of a new straining model with a depth dependent

166

blocking function (BRADFORD et al. 2006a).

167

BRADFORD et al. (2006b) further used HYDRUS-1D in order to simulate transport and168

straining of E. coli O157:H7 (bacteria significantly larger than colloids) in saturated porous

169

media. For this study, HYDRUS-1D was modified to accommodate a newly developed

170

conceptual model that assumed that E. coli will aggregate when large numbers of171

monodispersed E. coli are deposited at pore constrictions or straining sites. When the172

deposited E. coli reach a critical concentration at the straining site, the aggregated E. coli

173

O157:H7 will be released into the aqueous solution as a result of hydrodynamic shearing

174

forces.

7

175

In a series of papers, GARGIULO et al. (2007ab, 2008) used the two kinetic sites model176

in HYDRUS-1D to simulate the transport and deposition of two bacteria strains (Deinococcus

177

Radiodurans and Rhodococcus rhodochrous) under both saturated and unsaturated conditions.

178

While the role of bacteria surface hydrophobicity was studied by GARGIULO et al. (2008), the179role of the matrix grain size and the bacteria surface protein was evaluated by GARGIULO et al.

180(2007b), and the effect of metabolic activity (metabolically active and stationary phase

181

Deinococcus Radiodurans) was analyzed by GARGIULO et al. (2007a). In all three studies,182

HYDRUS-1D provided an excellent tool for evaluating various factors involved in the183transport and deposition of bacteria in soils.

184
185

Additional Physical Nonequilibrium Models

186ŠIMŮNEK et al. (2003) added two additional options to HYDRUS-1D for simulating

187

nonequilibrium preferential flow. The simpler model extended the mobile-immobile water

188

concept by assuming that the immobile water content can, like the mobile water content, be

189

transient, and that there can be a water transfer between these two domains (model (2) in

190

Figure 1). They also implemented GERKE and VAN GENUCHTEN‘s (1993) full dual-

191

permeability model. The dual-permeability model assumes that the porous medium consists of192

two overlapping pore domains, with water flowing relatively fast in one domain (often called193

the macropore, fracture, or inter-porosity domain, subscripts f and F) when close to full194

saturation, and slow in the other domain (often referred to as the micropore, matrix, or intra-

195

porosity domain, subscript m and M) (model (3) in Figure 1):196

θ = θ F + θ M = wθ f + (1 − w)θ m

197

Lower case subscripts in the dual-permeability model refer to the local (pore-region) scale,

198

while upper case subscripts refer to the global (total soil medium) scale. The dual-

199

permeability model was extended further by assuming that the liquid phase of the matrix can
8

(5)

200

be further partitioned into mobile (flowing), θmo,m, and immobile (stagnant), θim,m, regions as

201

follows:202

θ m = θ mo ,m + θim ,m203

where θm is the volumetric water content of the matrix pore system (POT et al. 2005; ŠIMŮNEK

204& VAN GENUCHTEN 2008).

205

Examples of applications of dual-porosity models to a range of laboratory and field

206

data involving transient flow and solute transport are provided by ŠIMŮNEK et al. (2001),

207

ABBASI et al. (2003ab), CASTIGLIONE et al. (2003), KÖHNE et al. (2004a, 2006a) and HAWS et208

al. (2005).

209

ABBASI et al. (2003ab) applied both equilibrium and nonequilibrium (the dual-porosity

210

model) models in HYDRUS-2D to simulate water flow and tracer movement at a sandy loam

211

field plot with furrows (3 m × 3 m) in Phoenix, Arizona, USA. Soil hydraulic and solute212

transport parameters were optimized by minimizing the objective function defined using

213

water contents, infiltration rates, and solute concentrations. The similarity of the results

214

obtained assuming equilibrium and nonequilibrium flow and transport led to the conclusion

215

that equilibrium transport prevailed at this field site (ABBASI et al. 2003ab).

216

The MIM approach in HYDRUS-1D was used by KÖHNE et al. (2004a) for inverse

217

simulation of pressure heads, water contents, water outflow and Br breakthrough in the218

effluent of six aggregated soil columns with different initial water contents subject to219intermittent irrigations. The physical nonequilibrium was more pronounced for wet and dry

220

soils than for intermediate initial water contents. Inverse identification of the required MIM

221parameters was fairly successful, except for the saturated water contents in mobile and222

immobile regions, which were highly correlated (KÖHNE et al., 2004a).

223

Water and solute fluxes from the subsurface drains of two macroporous silty clay loam

224

plots (48.5 m × 60 m) in agricultural fields in West Lafayette, Indiana, were simulated using
9

(6)

225

HYDRUS-2D by HAWS et al. (2005). The drainage fluxes (hydrographs) were matched226

reasonably well by both uniform and dual-porosity (MIM) models. However, a model

227

calibrated on drainage fluxes could not reproduce solute breakthrough. HAWS et al. (2005)

228

concluded that a hydrograph fit does not guarantee a proper description of flow patterns at the229

field scale, and that consideration of solute breakthrough is needed to derive physically230

meaningful model parameters.

231

KÖHNE et al. (2006a) made similar conclusions based on their field study at the Infeld

232site in North-West Germany. The observed rapid Br effluent breakthrough at low

233

concentrations could only be simulated using the MIM approach in HYDRUS-2D. Simulation

234

results suggested that over 60% of the surface applied Br was immobilized by transfer into the235

stagnant soil water region, and that the two-dimensional flow field induced by tile drains

236

enhanced Br dispersion (KÖHNE et al. 2006a).

237

A version of HYDRUS-1D that considers GERKE and VAN GENUCHTEN‘s (1993) dual-238

permeability flow and transport model was used in the studies of ŠIMŮNEK et al. (2001), ZHANG

239

et al. (2004), KÖHNE et al. (2004b, 2006bc), and POT et al. (2005), KODEŠOVÁ et al. (2006ab,

2402008), among many others.

241

ZHANG et al. (2004) used the dual-permeability model to analyze solute and colloidal

242

tracer tests in laboratory columns that examined the hydraulic properties of a foamed

243

zeolite/iron pellet material that was developed for in situ remediation of contaminated244

groundwater. The colloidal microspheres (1 μm diameter) moved through the columns at a245much faster rate than the nonreactive solute tracer tritiated water, reflecting the inter-pellet

246

preferential flow paths in the packed material. Inverse modeling of the microsphere data using

247

a physical nonequilibrium transport model yielded the immobile water content (θim)

248

equivalent to the intra-pellet porosity (0.40), suggesting that the microspheres were excluded249

from the small intra-pellet pores and could only move through the large inter-pellet pore

10

250

spaces. The dual-permeability dual-porosity model also indicated that 6 - 11% of the total

251

porosity was preferential flow porosity, consistent with the observation of enhanced

252

microsphere transport with respect to tritiated water. Forward modeling with the dual-

253

permeability dual-porosity model suggested that the preferential flow porosity will drastically254

lower contaminant removal efficiency.

255

POT et al. (2005) used all physical nonequilibrium models from HYDRUS-1D to

256

evaluate laboratory column studies under unsaturated steady-state flow conditions generated

257

using several rainfall intensities. POT et al. (2005) needed increasingly complex models to258

describe tracer displacement column experiments for increasing fluxes. Numerical analyses

259

showed that contrasting physical nonequilibrium transport processes occurred for different

260

fluxes. Multiple (three) porosity domains (dual-permeability model with immobile zone in the261

matrix) contributed to flow at the highest rainfall intensities, including preferential flow

262

through macropore pathways. Macropores were not active any longer at the intermediate and263

lowest velocities, and dual-porosity-type models were able to describe the observed

264

preferential transport well.

265

KÖHNE et al. (2006a) studied the feasibility of the inverse (Levenberg-Marquardt)

266

identification of dual-permeability model parameters from a drainage hydrograph. The dual-

267

permeability model implemented in HYDRUS-1D was used to fit hydraulic and transport

268

parameters, either sequentially or simultaneously, using observed tile-drainage hydrographs

269

and Br concentrations. Only the simultaneous fitting procedure was successful in describing

270

Br breakthrough. From these and lab-scale results, it was inferred that a hydrograph alone is271

insufficient for the inverse identification of soil hydraulic dual-permeability model parameters

272(KÖHNE et al. 2006a).

273

KÖHNE et al. (2006b) used the dual-permeability model to analyze experiments carried274

out on a laboratory column (80 cm height, 24 cm diam.) designed for hydraulic tomography

11

275

that provided outflows separately from the matrix and the central cylindrical preferential flow276

region, as well as pressure heads and water contents in the preferential flow and matrix

277

regions at various positions. KÖHNE et al. (2006b) applied one inverse approach that relied on

278

standard (lumped) observations of infiltration and outflow, while another approach considered279

separate outflows for the matrix and the PF region. Both inverse approaches provided

280

accurate matches of bulk infiltration and outflow. However, the outflows from either the281

matrix or the preferential flow region could only be described when the dual-permeability

282

model was fit to region-specific outflow data. KÖHNE et al. (2006b) concluded that for natural283

soils where experimental data do not come in a ‘separated form’ for fast and slow flow284

regions, the domain-related hydraulic parameters of the dual-permeability model may be285

difficult to identify from water flow observations alone.

286

KODEŠOVÁ et al. (2006a) used the single-porosity and dual-permeability models in

287

HYDRUS-1D to simulate variably-saturated water movement in clay soils with and without288

macropores. Numerical simulations of water flow for several scenarios of probable macropore

289

compositions showed a considerable impact of preferential flow on water infiltration in such

290

soils. Results of numerical simulations showed that cumulative infiltration into the soil with291

and without macropores may differ by two or more orders of magnitude. Thus, the292

appropriate models must be used to describe non-equilibrium flow in such soils.

293

KODEŠOVÁ et al. (2006b) carried out a soil micromorphological study to demonstrate294

the impact of soil organisms on soil pore structure. They showed the influence of earthworms,

295

enchytraeids and moles on the pore structure of a Greyic Phaeozem by comparing two soil

296

samples either affected or not affected by these organisms. They also studied macropores

297

created by roots and soil microorganisms in a Haplic Luvisol, and subsequently affected by298

clay coatings. The dual permeability models implemented in HYDRUS-1D were applied to299

improve the numerical inversion of the multi-step outflow experiment, and to obtain

12

300

parameters characterizing multimodal soil hydraulic properties.

301

Finally, KODEŠOVÁ et al. (2008) used micromorphological images for a) characterization

302

of flow domains in three soil types, b) selection of a proper model for estimation of multimodal

303

soil hydraulic properties from multi-step outflow experiment and ponding infiltration, and c)304

numerical simulation of chlorotoluron transport within the soil profile that was experimentally305

studied in the field. They showed that chlorotoluron was regularly distributed in the highly306

connected domain of larger pores of Haplic Luvisol, from which it penetrated into the soil

307

aggregates, i.e., zones of immobile water. The highest mobility of chlorotoluron in Greyic

308

Phaeozem was caused by larger capillary pore pathways and sufficient infiltration fluxes that309

occasionally filled up these pores. The presence of clay coatings in Greyic Phaeozem that310

restricted water flow and contaminant transport between the macropore and matrix domains311

was an additional cause for this preferential transport that produced chlorotoluron penetration

312

into deeper depths. Chlorotoluron was less regularly distributed in Haplic Cambisol. Despite313

the highest infiltration rate, preferential flow only slightly affected the herbicide transport.

314

Large gravitational pores that may dominate water flow and solute transport under saturated

315

conditions were inactive during the monitored period. As a result of complex interactions

316

between meteorological conditions and the soil pore structure, the single- and dual-porosity

317

models described the herbicide behavior in Haplic Luvisol well, while the dual-permeability

318

model performed better in simulating the herbicide transport in Greyic Phaeozem and Haplic

319

Cambisol.

320
321

Physical and Chemical Nonequilibrium Models

322

However, many transport situations involve not only physical or chemical323

nonequilibrium, but both nonequilibrium processes occur simultaneously. One obvious

324

example (e.g., ŠIMŮNEK & VAN GENUCHTEN 2008) occurs during transport through an

13

325

aggregated laboratory soil column involving steady-state water flow when both a conservative326

tracer (no sorption) and a reactive solute are used. The collected tracer breakthrough curve

327

may then display a relatively rapid initial breakthrough followed by extensive tailing, which328

are both features of nonequilibrium transport. Since the tracer is non-reactive, this329

nonequilibrium must be caused by physical factors. When the reactive solute is also sorbed

330

kinetically to the solid phase (an indication of a chemical nonequilibrium), the use of a model

331

that simultaneously considers both physical and chemical nonequilibrium is required.

332

HYDRUS offers two options for simulating simultaneous physical and chemical

333

nonequilibrium. The first model, i.e., the dual-porosity model with one kinetic site (model (d)334in Figure 1), considers water flow and solute transport in a dual-porosity system, while335

assuming that sorption in the immobile zone is instantaneous. However, the sorption sites in

336

contact with the mobile zone are now divided into two fractions, subject to either

337

instantaneous or kinetic sorption, similar to the two-site kinetic sorption concept. Since the338

residence time of solutes in the immobile domain is relatively large, equilibrium likely exists

339

between the solution and the sorption complex here, in which case there is no need to consider

340

kinetic sorption in the immobile domain. On the other hand, the model, assumes the presence341

of kinetic sorption sites in contact with the mobile zone, since water can move relatively fast342

in the macropore domain and thus prevent chemical equilibrium (ŠIMŮNEK & VAN

343

GENUCHTEN 2008).

344

The last nonequilibrium option implemented into HYDRUS-1D combines chemical345

nonequilibrium with the dual-permeability model (model (e) in Figure 1). This model assumes

346

that equilibrium and kinetic sites exist in both the macropore (fracture) and micropore

347(matrix) domains. A complete list and more detailed descriptions of the different models

348

summarized here, including the specific equations used for the water flow and solute transport349

models, are given in (ŠIMŮNEK & VAN GENUCHTEN 2008).

14

350

PANG et al. (2008) used the HYDRUS-1D mobile-immobile two-region model with351

one kinetic site (model (d) in Figure 1) to evaluate the transport of fecal coliforms, Salmonella352

bacteriophage and bromide (Br), in 30 undisturbed lysimeters constructed from undisturbed

353

New Zealand soils. The dual-porosity model considered first-order attachment/detachment354

and inactivation of the microbes, first-order mass transfer between two regions for Br, and355

convection and dispersion. Model predictions matched observations reasonably well.

356GÄRDENÄS et al. (2006) compared four two-dimensional transport models (a uniform

357

flow model, MIM, a dual-porosity model, and a dual-permeability model) in HYDRUS-2D in358

order to predict preferential water flow and the leaching of the herbicide MCPA in a 50-m

359

long transect through a sloping, heterogeneous, tile-drained field soil in South Sweden. The

360

simulated time covered six weeks following the spray application. Only the dual-permeability

361

and dual-porosity models reproduced the pesticide concentration patterns in drain outflow.

362

POT et al. (2006) used various physical nonequilibrium approaches with kinetic

363

sorption to analyze the impact of different constant rainfall rates on Br-, isoproturon and364

metribuzin leaching in undisturbed soil cores collected from the grassed filter strip.

365

Observations showed a strong impact of rainfall intensity on Br- (see also the text in the366previous section) and herbicide leaching. Herbicide transport was affected by kinetic sorption367

at all flow velocities. Significantly higher estimated values for degradation rate parameters, as

368

compared to batch data, were correlated with the extent of non-equilibrium sorption (POT et369al., 2005).

370

Similarly, KÖHNE et al. (2006c) applied HYDRUS-1D to simulate the transport of371

isoproturon, terbuthylazine and Br- observed in an aggregated loamy and a macroporous372

loamy sand soil column subject to several irrigation-redistribution cycles. The early373

isoproturon breakthrough in the aggregated loamy soil could be qualitatively predicted when374

using the dual-permeability model with two-site kinetic sorption. The simulated herbicide

15

375

breakthrough curves obtained were similar when either the degradation or sorption rate376

parameters were optimized. Due to internal parameter correlations and nonuniqueness,377

inverse model applications could not accurately distinguish between degradation and sorption

378

processes for conditions involving preferential flow.

379
380

Summary381

In this paper we have summarized a wide range of nonequilibrium water flow and382

solute transport models that are available in the latest versions of the HYDRUS software

383

packages and reviewed their recent applications. The models range from classical models

384

simulating uniform water flow and solute transport, to traditional mobile-immobile water385

physical and two-site chemical nonequilibrium models, to more complex dual-permeability

386

models that consider both physical and chemical causes of nonequilibrium.

387

In general, the most common applications can be divided into two large groups. The

388

first group uses HYDRUS models to simulate transport of colloids, viruses and bacteria, i.e.,389

abiotic and biotic particles in the 0.1 μm to 10 μm size range. These applications usually use390

the chemical nonequilibrium transport capabilities of the model to simulate processes such as391

attachment/detachment to the solid phase or air-water interface, or straining. The second392

group addresses preferential water flow and transport of dissolved chemicals in structured

393

soils with macropores or other preferential flow pathways. For this problem, physical394

nonequilibrium options are used, sometimes combined with chemical nonequilibrium options,395

for simulating the transport of pesticides and other organic compounds.

396

There remains a need for more thorough studies to evaluate how much and what type397

of information is required to fully parameterize selected nonequilibrium models. Ideally,398

parameterization would be based on straightforward observations of basic soil structural and399

physicochemical properties, much in the spirit of Hydropedology, a discipline which lately

16

400

received renewed attention (LIN et al. 2005), and which was always fostered by Prof. Kutílek.

401

Until measurement techniques are fully developed, inverse parameter identification is an402

indispensable means for the application of complex nonequilibrium models. To facilitate such

403

applications, the objective function for the inverse problem in HYDRUS can be formulated in404

terms of a large number of variables involving not only boundary concentration fluxes, but405

also water and solute distributions within the soil profile and in different phases (ŠIMŮNEK &406

VAN GENUCHTEN 2008). For these reasons we believe that HYDRUS is a very attractive tool

407

for analyzing both forward and inverse flow and transport problems.

408
409

References:

410

ABBASI F., ŠIMŮNEK J., FEYEN J., VAN GENUCHTEN M. TH., SHOUSE P. J. (2003a):

411

Simultaneous inverse estimation of soil hydraulic and solute transport parameters from412

transient field experiments: Homogeneous soil, Transactions of ASAE, 46(4): 1085-1095.

413

ABBASI F., JACQUES D., ŠIMŮNEK J., FEYEN, J., VAN GENUCHTEN M. TH. (2003b): Inverse

414

estimation of the soil hydraulic and solute transport parameters from transient field

415

experiments: Heterogeneous soil, Transactions of ASAE, 46(4): 1097-1111.

416

BRADFORD S. A., YATES S. R., BETTEHAR M., ŠIMŮNEK J. (2002): Physical factors affecting the417transport and fate of colloids in saturated porous media, Water Resour. Res., 38(12): 1327,

418

doi:10.1029/2002WR001340, 63.1-63.12.

419

BRADFORD S. A., ŠIMŮNEK J., BETTEHAR M., VAN GENUCHTEN M. TH., YATES, S. R. (2003):

420

Modeling colloid attachment, straining, and exclusion in saturated porous media, Environ.

421

Sci. Technol., 37(10): 2242-2250.

422

BRADFORD S. A., BETTEHAR M., ŠIMŮNEK J., VAN GENUCHTEN M. TH. (2004): Straining and423

attachment of colloids in physically heterogeneous porous media, Vadose Zone Journal,424

3(2): 384-394.

17

425

BRADFORD S. A., ŠIMŮNEK J., BETTAHAR M., VAN GENUCHTEN M. TH., YATES S. R. (2006a):

426

Significance of straining in colloid deposition: evidence and implications, Water Resour.

427

Res., 42: W12S15, doi:10.1029/2005WR004791, 16 pp.

428

BRADFORD S. A., ŠIMŮNEK J., WALKER S. L. (2006b): Transport and straining of E. coli429

O157:H7430

doi:10.1029/2005WR004805, 12 pp.

insaturatedporousmedia.

Water

Resour.

Res.,

42:

W12S12,

431

CASTIGLIONE P., MOHANTY B. P., SHOUSE P. J., ŠIMŮNEK J., VAN GENUCHTEN M. TH., SANTINI A.

432(2003): Lateral water diffusion in an artificial macroporous system: Theory and experimental

433

evidence, Vadose Zone Journal, 2: 212-221.

434

DOUSSET S., THEVENOT M., POT V., ŠIMŮNEK J., ANDREUX F. (2007): Evaluating equilibrium435

and non-equilibrium transport of bromide and isoproturon in disturbed and undisturbed

436

soil columns, J. Contam. Hydrol., (in press).

437GÄRDENÄS A., ŠIMUNEK J., JARVIS N., VAN GENUCHTEN M. TH. (2006): Two-dimensional

438

modelling of preferential water flow and pesticide transport from a tile-drained field, J.

439

Hydrology, 329: 647-660.

440

GARGIULO G., BRADFORD S. A., ŠIMŮNEK J., USTOHAL P., VEREECKEN H., KLUMPP E.

441(2007a): Transport and deposition of metabolically active and stationary phase442

Deinococcus Radiodurans in unsaturated porous media, Environ. Sci. and Technol., 41(4):443

1265-1271.

444

GARGIULO G., BRADFORD S. A., ŠIMŮNEK J., USTOHAL P., VEREECKEN H., KLUMPP E.

445(2007b): Bacteria transport and deposition under unsaturated conditions: the role of the446

matrix grain size and the bacteria surface protein, J. Contam. Hydrology, 92: doi:447

10.1016/j.jconhyd.2007.01.009, 255-273.

448

GARGIULO G., BRADFORD S. A., ŠIMŮNEK J., USTOHAL P., VEREECKEN H., KLUMPP E. (2008):

449

Bacteria transport and deposition under unsaturated conditions: the role of bacteria surface

18

450

hydrophobicity, Vadose Zone Journal, (in press).

451

GERKE H. H., VAN GENUCHTEN M. TH. (1993): A dual-porosity model for simulating the452

preferential movement of water and solutes in structured porous media, Water Resour.

453

Res., 29: 305-319.

454

GONÇALVES M. C., ŠIMŮNEK J., RAMOS T. B., MARTINS J. C., NEVES M. J., PIRES F. P. (2006):

455

Multicomponent solute transport in soil lysimeters irrigated with waters of different

456

quality, Water Resour. Res., 42: W08401, doi:10.1029/2006WR004802, 17 pp.

457

HAWS N. W., RAO P. S. C., ŠIMŮNEK J., POYER I. C. (2005): Single-porosity and dual-porosity

458

modeling of water flow and solute transport in subsurface-drained fields using effective

459

field-scale parameters, J. of Hydrology, 313(3-4): 257-273.

460

JARVIS N. J. (1994): The MACRO model (Version 3.1), Technical description and sample

461

simulations, Reports and Dissertations 19, Dept. Soil Sci., Swedish Univ. Agric. Sci.,

462

Uppsala, Sweden, 51 pp.

463
464

KODEŠOVÁ R., KOZÁK J., ŠIMŮNEK J. (2006a): Numerical study of macropore impact onponded infiltration in clay soils, Soil and Water Research, 1(1): 16-22.

465

KODEŠOVÁ R., KODEŠ V., ŽIGOVÁ A., ŠIMŮNEK J. (2006b): Impact of plant roots and soil

466

organisms on soil micromorphology and hydraulic properties, Biologia, 61/Suppl. 19:

467

S339-S343.

468

KODEŠOVÁ R., KOČÁREK M., KODEŠ V., ŠIMŮNEK J., KOZÁK J. (2008): Impact of soil

469

micromorphological features on water flow and herbicide transport in soils, Vadose Zone

470

Journal, Special Issue “Vadose Zone Modeling”, (in press).

471

KÖHNE J. M., KÖHNE S., MOHANTY B. P., ŠIMŮNEK J. (2004a): Inverse mobile-immobile472

modeling of transport during transient flow: Effect of between-domain transfer and initial473

water content, Vadose Zone Journal, 3(4): 1309-1321.

19

474

KÖHNE J. M., MOHANTY B., ŠIMŮNEK J., GERKE H. H. (2004b): Numerical evaluation of a475

second-order water transfer term for variably saturated dual-permeability models, Water

476

Resour. Res., 40: doi:10.1029/2004WR00385.

477

KÖHNE S., LENNARTZ B., KÖHNE J. M., ŠIMŮNEK J. (2006a): Bromide transport at a tile-

478

drained field site: experiment, one- and two-dimensional equilibrium and non-equilibrium

479

numerical modeling, J. Hydrology, 321(1-4): 390-408.

480

KÖHNE J. M., MOHANTY B. P., ŠIMŮNEK J. (2006b): Inverse dual-permeability modeling of

481

preferential water flow in a soil column and implications on field-scale solute transport

482

prediction, Vadose Zone Journal, 5: 59-76.

483

KÖHNE J. M., KÖHNE S., ŠIMŮNEK J. (2006c): Multi-Process Herbicide Transport in Structured

484

Soil Columns: Experiment and Model Analysis, Journal of Contaminant Hydrology, 85:

485

1-32.

486

KÖHNE J. M., KÖHNE S., ŠIMŮNEK J. (2008a): A Review of Model Applications for Structured

487

Soils: a) Water Flow and Tracer Transport, J. Contam. Hydrology, Special Issue “Flow488

Domains”, (submitted).

489

KÖHNE J. M., KÖHNE S., ŠIMŮNEK J. (2008b): A Review of Model Applications for Structured

490

Soils: b) Pesticide Transport, J. Contam. Hydrology, Special Issue “Flow Domains”,

491(submitted).

492KUTÍLEK M. (1978): Water Management Soil Science, SNTL, Praha, (in Czech).

493KUTÍLEK M. (1984): Moisture of Porous Materials, SNTL, Praha, (in Czech).

494

KUTILEK M., NIELSEN D. R. (1994): Soil Hydrology, Catena Verlag, Cremlingen, 370 pp.

495

LANGERGRABER G., ŠIMŮNEK J. (2005): Modeling variably-saturated water flow and multi-

496

component reactive transport in constructed wetlands, Vadose Zone Journal, 4, 924-938.

20

497

LANGERGRABER G., ŠIMŮNEK J. (2006): The Multi-component Reactive Transport Module

498

CW2D for Constructed Wetlands for the HYDRUS Software Package, Manual – Version

499

1.0, HYDRUS Software Series 2, Department of Environmental Sciences, University of500

California Riverside, Riverside, CA, 72 pp.

501
502

LIN H., BOUMA J., WILDING L. P., RICHARDSON J. L., KUTÍLEK M., NIELSEN D. R. (2005):
Advances in hydropedology, Advances in Agronomy, 85: 1-97.

503

PANG L., MCLEOD M., AISLABIE J., ŠIMŮNEK J., CLOSE M., HECTOR R. (2008): Modeling504

transport of fecal coliforms and Salmonella bacteriophage in ten undisturbed New Zealand505

soils under dairy shed effluent irrigation, Vadose Zone Journal (in press).

506

POT V., ŠIMŮNEK J., BENOIT P., COQUET Y., YRA A., MARTÍNEZ-CORDÓN M.-J. (2005): Impact507

of rainfall intensity on the transport of two herbicides in undisturbed grassed filter strip

508

soil cores, J. of Contaminant Hydrology, 81: 63-88.

509
510
511
512

PRUESS K. (1991): Tough2 – A general-purpose numerical simulator for multiphase fluid andheat flow, Rep. LBL-29400, Lawrence Berkeley Laboratory, Berkeley, CA.
SCHIJVEN J., ŠIMŮNEK J. (2002): Kinetic modeling of virus transport at field scale, J. of Contam.
Hydrology, 55(1-2): 113-135.

513

ŠEJNA M., ŠIMŮNEK J. (2007): HYDRUS (2D/3D): Graphical User Interface for the HYDRUS

514

Software Package Simulating Two- and Three-Dimensional Movement of Water, Heat,

515

and Multiple Solutes in Variably-Saturated Media, published online at www.-pc-

516

progress.cz, PC-Progress, Prague, Czech Republic.

517
518
519
520
521

ŠIMŮNEK J., SUAREZ D. L. (1993): Modeling of carbon dioxide transport and production in soil:
1. Model development, Water Resour. Res., 29(2): 487-497.
ŠIMŮNEK J., SUAREZ D. L. (1994): Two-dimensional transport model for variably saturated
porous media with major ion chemistry, Water Resour. Res., 30(4): 1115-1133.
ŠIMŮNEK J., ŠEJNA M., VAN GENUCHTEN M. TH. (1998): The HYDRUS-1D software package for

21

522

simulating the one-dimensional movement of water, heat, and multiple solutes in variably-

523

saturated media, Version 2.0, IGWMC-TPS-70, International Ground Water Modeling

524

Center, Colorado School of Mines, Golden, Colorado, 202 pp.

525

ŠIMŮNEK J., WENDROTH O., WYPLER N., VAN GENUCHTEN M. TH. (2001): Nonequilibrium water526

flow characterized from an upward infiltration experiment, European J. of Soil Sci., 52(1):527

13-24.

528

ŠIMŮNEK J., JARVIS N. J., VAN GENUCHTEN M. TH., GÄRDENÄS A. (2003): Review and529

comparison of models for describing non-equilibrium and preferential flow and transport

530

in the vadose zone, Journal of Hydrology, 272: 14-35.

531

ŠIMŮNEK J. (2005): Chapter 78: Models of Water Flow and Solute Transport in the532

Unsaturated Zone, In: M. G. Anderson and J. J. McDonnell (Editors), Encyclopedia of533

Hydrological Sciences, John Wiley & Sons, Ltd., Chichester, England, 1171-1180.

534

ŠIMŮNEK J., VAN GENUCHTEN M. TH., ŠEJNA M. (2005): The HYDRUS-1D software package for535

simulating the one-dimensional movement of water, heat, and multiple solutes in variably-

536

saturated media, Version 3.0, HYDRUS Software Series 1, Department of Environmental

537

Sciences, University of California Riverside, Riverside, CA, 270 pp.

538ŠIMŮNEK J., VAN GENUCHTEN M. TH., ŠEJNA M. (2007): The HYDRUS Software Package for

539

Simulating Two- and Three-Dimensional Movement of Water, Heat, and Multiple Solutes in

540

Variably-Saturated Media, Technical Manual, Version 1.0, PC Progress, Prague, Czech

541

Republic, pp. 241.

542

ŠIMŮNEK J., VAN GENUCHTEN M. TH., ŠEJNA M. (2008): Development and applications of the543

HYDRUS and STANMOD software packages, and related codes, Vadose Zone Journal,

544

Special Issue “Vadose Zone Modeling”, (in press).

545ŠIMŮNEK J., VAN GENUCHTEN M. TH. (2008): Modeling nonequilibrium flow and transport546

with HYDRUS, Vadose Zone Journal, Special Issue “Vadose Zone Modeling”, (in press).

22

547

SUAREZ D. L., ŠIMŮNEK J. (1993): Modeling of carbon dioxide transport and production in soil:

548

2. Parameter selection, sensitivity analysis and comparison of model predictions to field data,

549

Water Resour. Res., 29(2): 499-513.

550

VAN DAM J. C., HUYGENJ., WESSELING J. G., FEDDES R. A., KABAT P., VAN VALSUM P. E. V.,

551

GROENENDIJK P., VAN DIEPEN C. A. (1997): Theory of SWAP, version 2.0. Simulation of552

water flow, solute transport and plant growth in the Soil-Water-Atmosphere-Plant

553

environment, Dept. Water Resources, WAU, Report 71, DLO Winand Staring Centre,

554

Wageningen, Technical Document 45.

555
556

VAN GENUCHTEN M. TH., WIERENGA P. J. (1976): Mass transfer studies in sorbing porous media,

I. Analytical solutions, Soil Sci. Soc. Am. J., 40: 473-481.

557

VAN GENUCHTEN M. TH., WAGENET R. J. (1989): Two-site/two-region models for pesticide

558

transport and degradation: Theoretical development and analytical solutions, Soil Sci. Soc.

559

Am. J., 53: 1303-1310.

560

ZHANG P., ŠIMŮNEK J., BOWMAN R. S. (2004): Nonideal transport of solute and colloidal

561

tracers

throughreactive562

doi:10.1029/2003WR002445.

zeolite/ironpellets,

23

Water

Resour.

Res.,

40:

Physical Equilibrium and Nonequilibrium Models
1. Uniform Flow Model

2. Dual-Porosity Model

3. Dual-Permeability Model

Water

Water

Immob.

Water

MobileSlow

SoluteSoluteSoluteImmob.

θ

Fast

Slow

Mobileθ = θim + θmoFastθ = θM + θF

Chemical Nonequilibrium Models
a. One Kinetic Site
Modelb. Two-Site Modelc. Two Kinetic Sites
Model

d. Dual-Porosity Model withOne Kinetic Sitee. Dual-permeability Modelwith Two-Site Sorption

Slow Fast

Immob. Mob.

θ

se

θ

s 1k

θ

θim

θmo

cim

cmo

smoe

s me

smok

s mk

θM

θF

cm

cf

sfe

sime

sk
c

sk

c

c

s 2k

sfk

Physical and Chemical Nonequilibrium Models
Figure 1. Physical and chemical nonequilibrium

models considered by HYDRUS-1D.
24