R

Reply to “Comment on ‘Evalua ng Interac ons between
Groundwater and Vadose Zone Using the HYDRUS-based Flow
Package for MODFLOW’” by Navin Kumar C. Twarakavi, JirkaŠimůnek, and Sophia Seo
Navin Kumar C. Twarakavi,* Jirka Šimůnek, and Sophia SeoNP(2009) commented on thecomparison of the UZF1 (Niswonger et al., 2006) andHYDRUS (Seo et al., 2007) packages and their ability to incorporate vadose zone flow processes into the MODFLOW model,presented in our manuscript (Twarakavi et al., 2008). Niswongerand Prudic’s (2009) comments give us an opportunity to better clarify the objectives of our paper (Twarakavi et al., 2008).
Our main objective was to introduce the HYDRUS package for
MODFLOW as another option for modeling coupled vadosezone and groundwater fluxes with the existing approaches, suchas the REC-ET, VSF (Thoms et al., 2006), and UZF1 packages.
Twarakavi et al. (2008) concluded that because the HYDRUS
package can account for depth-dependent root water uptake, water
fluxes in both downward and upward directions, and multiple soil
horizons with varying hydraulic properties, while the UZF1 package considers only homogeneous soil profiles, the “HYDRUS
package thus offers a good alternative to the UZF1 package whenthese factors, e.g., vegetation or multiple horizons, are significant
for a particular application.”First, Niswonger and Prudic (2009) note that the UZF1 package allows for groundwater discharge to the land surface (such as astream or a lake) whenever the water table in a particular cell exceeds
the land-surface elevation. It is true that the current version of theHYDRUS package does not address this process, since its main
objective is to account for flow in the unsaturated zone. On theother hand, the HYDRUS package implements an “atmospheric”
boundary condition for the upper boundary, which allows rainfall
or irrigation water in excess of the soil infiltration capacity to either
accumulate at the soil surface and create a dynamic surface waterlayer or generate surface runoff. We believe that the HYDRUS
package also handles evapotranspiration, an important process ofN.K.C. Twarakavi, Dep. of Agronomy and Soils, Auburn Univ., Auburn, AL36849; J. Šimůnek, Dep. of Environmental Sciences, Univ. of California Riverside, Riverside, CA 92521; S. Seo, Colorado School of Mines, 1500 Illinois St.,
Golden, CO 80401. *Corresponding author (twarakavi@auburn.edu).
Vadose Zone J. 8:820–821
doi:10.2136/vzj2008.0004L© Soil Science Society of America
677 S. Segoe Rd. Madison, WI 53711 USA.
All rights reserved. No part of this periodical may be reproduced or transmi edin any form or by any means, electronic or mechanical, including photocopying,recording, or any informa on storage and retrieval system, without permission inwri ng from the publisher.

www.vadosezonejournal.org · Vol. 8, No. 3, August 2009820

water discharge through the land surface that can be particularly
important under shallow water tables conditions, better. While theUZF1 package estimates evapotranspiration losses as a function of
the extinction depth using a simple empirical function that may
not always represent this complex process correctly, the HYDRUS
package can better approximate it by solving the Richards equationand can even simulate capillary rise. Depth-dependent calculationsof water loss by evapotranspiration used in the UZF1 package are avery rough approximation of this highly nonlinear process, which
the HYDRUS package can handle more rigorously.
Next, Niswonger and Prudic (2009) point out that the computational speed of the UZF1 package for the “Case Study 3:
Hypothetical Regional-Scale Groundwater Problem” in Twarakavi
et al. (2008) can be greatly improved and, as a result, UZF1 cangenerate results faster than the HYDRUS package. Niswonger
and Prudic (2009) note that when discretization of drainage waves
is reduced from 25 (used in Niswonger et al., 2006) to 10 increments, the UZF1 package runs two times faster than HYDRUS.
As in Niswonger et al. (2006), the HYDRUS package simulationof this case study was simply intended to illustrate the applicability of this package to larger scales without either greatly compromising the computational speed of the overall system or optimizing spatial and temporal discretization. We would also be able tosignificantly reduce the required computational time by (i) usingcoarser discretization of the HYDRUS soil profiles, (ii) decreasing
the number of zones, (iii) using coarser temporal discretization,
and/or (iv) reducing internal iteration criteria. However, it is ourbelief that it is preferable to use finer discretizations and stricter
convergence criteria, thus guaranteeing a proper description offlow processes, to optimize computation speed. We agree with
Niswonger and Prudic (2009) regarding the need to demonstrate
the applicability of the HYDRUS package to a realistic watershedscale problem. Since that was beyond the scope of the paper underdiscussion, we intend to present a few computationally demanding
and realistic case studies in the public library of HYDRUS examples on the HYDRUS Web site in the near future (http://www.
pc-progress.com/en/Default.aspx?h1d-library).
Niswonger and Prudic (2009) also noted that they do notagree with the “assumption that HYDRUS’s method for simulating one-dimensional unsaturated flow through layered sedimentis more applicable than assuming an effectively homogeneous
unsaturated zone.” Even though the current state of practice ingroundwater flow modeling still considers uniform overlaying

soils in the vadose zone, we believe that this should not remainthe status quo and that this simplification hampers the advancement of research in ground water modeling. While we agree thatdetailed information on layering may not always be available andthat considering numerous layers in such cases may not be feasible,
in situations where such information is available, it may be beneficial to consider. Most soil profiles conform to a pattern of horizons
that should be considered to improve the reliability of calculated
vadose zone fluxes. Further, Niswonger and Prudic (2009) indicatethat the HYDRUS package’s inability to simulate lateral spreadingin the unsaturated zone renders consideration of layering ineffective. We agree that no model that considers only one-dimensionalvertical flow can be effective for conditions where lateral spreading(e.g., in perched water layers) in the unsaturated zone is a process
significantly affecting groundwater recharge. For conditions wherethese processes need to be considered, fully three-dimensionalmodels, such as VSF (Thoms et al., 2006) or HYDRUS (2D/3D)(Šimůnek et al., 2006), should be used. The HYDRUS package,
which obviously can consider homogeneous conditions, is meantto be a tool for use when an assumption of one-dimensional flowin the vadose zone is reasonable. It is up to users of the HYDRUS
package to decide whether to consider layering or opt out for ahomogeneous description of the vadose zone. We note here that
almost all widely used one-dimensional models simulating waterflow and solute transport in the vadose zone currently do consider
soil layering, and that a model without such an option would likely
be rejected by both researchers and the practicing community.
Further, Niswonger and Prudic (2009) commented on waterbalance with regard to the HYDRUS package. The HYDRUS
package for MODFLOW maintains the water budget during simulations. The HYDRUS package considers the exact distribution ofwater contents above the groundwater as a solution of the Richards
equation (reflecting different textures and retention curves),
which is obtained using a mass-conservative numerical scheme,
thereby maintaining a sound mass balance. The mass balance for
the HYDRUS package simulations is outputted when requiredby the user for verification purposes. Similarly, MODFLOW
itself reports the mass balance when required by the user. Thereis obviously a slight inconsistency between MODFLOW thatrepresents the vadose zone by a constant, that is, a “specific yield,”while HYDRUS uses a complete description using full definition of retention curves. However, since both models report theirmass balances, the user has complete information about it. In cases
where the coupled mass balance is not achieved satisfactorily, wewould like to note that a version of the HYDRUS package is available that can attain a complete mass balance between the saturated
and unsaturated zones by suitably altering the unsaturated profiles
for each zone at the end of each time step. However, it is recommended that this version of the HYDRUS package should only be
used on a case-by-case basis.
We agree with Niswonger and Prudic (2009) that the numerical stability, or the lack of such stability, in the HYDRUS package
will play an important role in its successful adoption by research or
the practicing community. The numerical stability of models basedon the numerical solution of the Richards equation still depends
partly on a proper selection of spatial and temporal discretization
and related iteration criteria, and thus on the experience of theirusers. For this reason, we exert a great effort, in addition to model
development, in training HYDRUS users (short courses, online
www.vadosezonejournal.org · Vol. 8, No. 3, August 2009tutorials) or providing them with information (see, e.g., the public
library of HYDRUS projects at http://www.pc-progress.com/en/Default.aspx?h1d-library) on how to successfully use our models. The hundreds of successful applications of HYDRUS models
appearing in peer-reviewed literature (Šimůnek et al., 2008) (andlikely many more unreported) offer proof that we are on the righttrack and that the issues associated with numerical instability canbe overcome.
Although the UZF1 package has been proven to work insome case studies (Niswonger et al., 2006), since the package isbased on the Kinematic wave approach, which assumes that variably saturated water flow is driven exclusively by gravitational
potential, and ignores the effect of capillary and other forces, it
provides a poorer characterization of water flow in the vadose zone
than approaches based on the Richards equation. In addition, more
research and characterization of water flow in the vadose zone invarious soils have been performed using the Richards equation thanthe Kinematic wave approach. We believe that coupling the onedimensional Richards equation–based package to MODFLOW
balances computational speed and model efficiencies at different
scales and represents a step in the right direction.
AThe activities on which this report is based were financed in part by theDepartment of the Interior, U.S. Geological Survey, through The Regentsof the University of California and the UC Center for Water Resources.

References
Niswonger, R.G., and D.E. Prudic. 2009. Comment on “Evaluating interactions
between groundwater and vadose zone using the HYDRUS-based flow
package for MODFLOW.” Vadose Zone J. 8:818–819 (this issue).
Niswonger, R.G., D.E. Prudic, and R.S. Regan. 2006. Documentation of theunsaturated-zone flow (UZF1) package for modeling unsaturated flow
between the land surface and the water table with MODFLOW-2005.
Techniques and Methods 6-A19. Available at http://pubs.usgs.gov/tm/2006/tm6a19/ (verified 13 May 2009). USGS, Denver, CO.
Seo, H.S., J. Šimůnek, and E.P. Poeter. 2007. Documentation of the HYDRUS
package for MODFLOW-2000, the U.S. Geological Survey modular
ground-water model. GWMI 2007-01. Int. Ground Water Modeling Ctr.,
Colorado School of Mines, Golden, CO.
Šimůnek, J., M.Th. van Genuchten, and M. Šejna. 2006. The HYDRUS software
package for simulating two- and three-dimensional movement of water, heat,and multiple solutes in variably saturated media. p. 241. In Technical manual,version 1.0. PC Progress, Prague, Czech Republic.
Šimůnek, J., M.Th. van Genuchten, and M. Šejna. 2008. Development andapplications of the HYDRUS and STANMOD software packages, andrelated codes. Vadose Zone J. 7:587–600.
Thoms, R.B., R.L. Johnson, and R.W. Healy. 2006. User’s guide to the variably
saturated flow (VSF) process for MODFLOW. Techniques and Methods
6-A18. Available at http://pubs.usgs.gov/tm/2006/tm6a18/ (verified 13May 2009). USGS, Reston, VA.
Twarakavi, N.K.C., J. Šimůnek, and S. Seo. 2008. Evaluating interactions between
groundwater and vadose zone using the HYDRUS-based flow package for
MODFLOW. Vadose Zone J. 7:757–768.

821