Documentation of the Unsaturated-Zone Flow (UZF1) 

Package for Modeling Unsaturated Flow Between the 

Land Surface and the Water Table with MODFLOW-2005

Techniques and Methods 6-A19

A Product of the Ground-Water Resources Program

U.S. Department of the Interior

U.S. Geological Survey

Chapter 19 of

Section A, Ground Water, of

Book 6, Modeling Techniques

Evapotranspiration (ET)

ET extinction depth

Deep water table

One-dimensional unsaturated flow

Three-dimensional ground-water flow

Land surface

Overland flow 

to streams and lakes

B

Evapotranspiration (ET)

ET extinction depth

Shallow water table

One-dimensional unsaturated flow

Three-dimensional ground-water flow

Land surface

Overland flow 

to streams and lakes

A

Evapotranspiration (ET)

ET extinction depth

Deep water table

One-dimensional unsaturated flow

Three-dimensional ground-water flow

Land surface

Overland flow 

to streams and lakes

B

Evapotranspiration (ET)

ET extinction depth

Shallow water table

One-dimensional unsaturated flow

Three-dimensional ground-water flow

Land surface

Overland flow 

to streams and lakes

A

A. B.
Documentation of the Unsaturated-Zone 

Flow (UZF1) Package for Modeling 

Unsaturated Flow Between the Land Surface 

and the Water Table with MODFLOW-2005

By Richard G. Niswonger, David E. Prudic, and R. Steven Regan

Chapter 19 of

Section A, Ground Water, of

Book 6, Modeling Techniques

Techniques and Methods 6-A19

U.S. Department of the Interior

U.S. Geological Survey
U.S. Department of the Interior

Gale A. Norton, Secretary

U.S. Geological Survey

P. Patrick Leahy, Acting Director

U.S. Geological Survey, Reston, Virginia: 2006

For sale by U.S. Geological Survey, Information Services

Box 25286, Denver Federal Center

Denver, CO 80225

For more information about the USGS and its products:

Telephone: 1-888-ASK-USGS

World Wide Web: http://www.usgs.gov/

Any use of trade, product, or firm names in this publication is for descriptive purposes only and does not imply 

endorsement by the U.S. Government.

Although this report is in the public domain, permission must be secured from the individual copyright owners to 

reproduce any copyrighted materials contained within this report.

Suggested citation:

Niswonger, R.G., Prudic, D.E., and Regan, R.S., 2006, Documentation of the Unsaturated-Zone Flow (UZF1) Package 

for modeling unsaturated flow between the land surface and the water table with MODFLOW-2005: U.S. Geological 

Survey Techniques and Methods 6-A19, 62 p. 
iii

PREFACE

A new Unsaturated-Zone Flow (UZF1) Package was written for use with the U.S. Geological 

Survey (USGS) MODFLOW-2005 ground-water model. The UZF1 Package is designed to simulate 

percolation through an unsaturated zone between land surface and the water table. The 

performance of this computer program has been tested in models of hypothetical ground-water 

flow systems; however, future applications of the programs may reveal errors that were not 

detected in the test simulations. Users are requested to notify the USGS if errors are found in 

the documentation report or in the computer program. Correspondence regarding the report or 

program should be sent to:

U.S. Geological Survey, WRD

2730 North Deer Run Road

Carson City, Nevada 89701 

Attention: Richard G. Niswonger or David E. Prudic

Although the computer program has been written and used by the USGS, no warranty, 

expressed or implied, is made by the USGS or the United States Government as to the accuracy 

and functionality of the program and related program material. Nor shall the fact of distribution 

constitute any such warranty, and no responsibility is assumed by the USGS in connection 

therewith.

MODFLOW-2005, the UZF1 Package, and other ground-water programs are available from the 

USGS at the following World Wide Web (WWW) address:

http://water.usgs.gov/software/ground_water.html
iv

Figures

Figure 1. Diagram showing one-dimensional unsaturated-zone flow coupled to

three-dimensional ground-water flow ................................................... 3

Figure 2. Diagram showing a wetting front moving through a uniform column of

unsaturated material affected by a decrease in surface flux, and results from 

a kinematic-wave solution of the wetting front represented by leading and 

trailing waves .............................................................................. 6

Figure 3. Graph showing kinematic-wave approximation of a wetting front moving 

through a uniform column of unsaturated material affected by a constant

evapotranspiration rate of 4.35 millimeters per day .................................... 7

Figure 4. Graphs showing comparison of results from the kinematic-wave 

approximation used to simulate unsaturated flow (UZF1) with a 

two-dimensional finite-difference solution of Richards’ equation (VS2DT; 

Healy, 1990) assuming evapotranspiration demand rates of: (A) 5 × 10-8 meter 

per second, (B) 5 × 10-7 meter per second, and (C) 1 × 10-6 meter per second...... 11

Contents

PREFACE ................................................................................................... iii

Abstract ................................................................................................... 1

Introduction................................................................................................ 1

Purpose and Scope ................................................................................. 2

Acknowledgments ................................................................................. 2

Theory and Conceptualization of Unsaturated-Zone Flow .......................................... 4

Solution Method .................................................................................... 4

Representation of Wetting and Drying Fronts ............................................. 5

ET in the Unsaturated Zone .................................................................. 7

Infiltration Boundary ........................................................................ 8

Coupling UZF1 to MODFLOW-2005 ............................................................... 8

Ground-Water Discharge to Land Surface ................................................ 9

Steady-State Simulations .................................................................. 10

Assumptions and Limitations ..................................................................... 10

Comparison of UZF1 and VS2DT Simulations ......................................................... 10

Description .......................................................................................... 10

Results................................................................................................ 11

Test Simulations .......................................................................................... 12

Test Simulation 1: Sagehen Watershed ......................................................... 12

Test Simulation 2: Simple Hypothetical Problem ................................................ 20

Implementation of UZF1 With MODFLOW-2005 ...................................................... 26

Summary ................................................................................................... 26

References Cited.......................................................................................... 27

Appendix 1: Data Input Instructions for the Unsaturated-Zone Flow (UZF1) Package ......... 29

Appendix 2: Selected Input Datasets and Printed Results for Test Simulation 2 ............... 33
v

Figure 5. Map showing distribution of hydraulic conductivity from ground-water

model used for test simulation 1 ......................................................... 13

Figure 6. Map showing distribution of steady-state ground-water depths used for test

simulation 1 ................................................................................. 15

Figure 7. Graph showing infiltration rate and evapotranspiration demand used for test

simulation 1 ................................................................................. 16

Figure 8. Map showing distribution of infiltration factors used for test simulation 1 ......... 17

Figure 9. Graph showing volumetric rates of infiltration, recharge, and

evapotranspiration summed over the model domain for test simulation 1 ......... 18

Figure 10. Graph showing streamflow at the basin outlet, ground-water discharge to

land surface plus saturation excess, and ground-water discharge to

streams for test simulation 1 ............................................................... 18

Figure 11. Maps showing distribution of ground-water seepage to land surface for (A)

steady state, and (B) July 1. ............................................................... 19

Figure 12. Graph showing water-content profiles in the unsaturated zone of a model

cell located in an upland area (column 9, row 36) at selected dates for test

simulation 1 ................................................................................. 20

Figure 13. Diagram showing well locations, elevation of top of consolidated rocks

beneath basin-fill aquifer, and distribution of hydraulic conductivity and

specific yield of the basin-fill aquifer used in test simulation 2 ........................ 21

Figure 14. Diagram showing hypothetical basin-fill aquifer with model grid,

land-surface contours, active cells, and stream segment and reach

numbering for test simulation 2 ............................................................ 22

Tables

Table 1. Unsaturated zone variables used in comparison of results from the

Unsaturated-Zone Flow (UZF1) Package in MODFLOW-2005 with those from

the Variably-Saturated Two-Dimensional Flow and Transport (VS2DT)

Model ....................................................................................... 10

Table 2. MODFLOW-2005 Packages and files used for test simulation 1 ..................... 14

Table 3. Hydraulic properties and other selected variables used in the

Layer-Property Flow (LPF) and the Unsaturated-Zone Flow (UZF1) Packages

for test simulation 1 ........................................................................ 14

Table 4. Method for calculating stream depth and width, inflow rates, streambed

properties, stream dimensions, and Manning’s roughness coefficients for

stream segments used in test simulation 2 ............................................. 23

Table 6. Hydraulic properties and other selected variables used in the

Layer-Property Flow (LPF) and Unsaturated-Zone Flow (UZF1) Packages for

test simulation 2.............................................................................. 24

Table 5. Specified infiltration and pumping rates for the 12 stress periods used in test

simulation 2 ................................................................................. 24

Table 7. Computed ground-water budget for test simulation 2 ................................. 25

Table 8. Computed unsaturated-zone budget for test simulation 2.............................. 25

Figures—Continued
vi

Conversion Factors, Datums, and Acronyms

Conversion Factors

Multiply By To obtain

cubic centimeter (cm3

) 0.06102 cubic inch (in3)

cubic meter (m3

) 264.2 gallon (gal)

cubic meter (m3

) 35.31 cubic foot (ft3)

cubic meter per day (m3

/d) 35.31 cubic foot per day (ft3/d)

cubic meter per second (m3

/s) 35.31 cubic foot per second (ft3/s)

cubic meter per year (m3

/yr) .000811 acre-foot per year (acre-ft/yr)

kilometer (km) .6214 mile (mi)

meter (m) 3.281 foot (ft)

meter per day (m/d) 3.281 foot per day (ft/d)

meter per second (m/s) 3.281 foot per second (ft/s)

millimeter (mm) .03937 inch (in.)

square kilometer (km2

) 247.1 acre

square kilometer (km2

) .3861 square mile (mi2)

square meter (m2

) .0002471 acre 

square meter (m2

) 10.76 square foot (ft2)

Datums

Vertical coordinate information is referenced to the North American Vertical Datum of 1988 

(NAVD 88).

Horizontal coordinate information is referenced to the North American Datum of 1983 (NAD 83).

Altitude, as used in this report, refers to distance above the vertical datum.

Acronyms

Acronyms Meaning

BCF Block-Centered Flow

ET evapotranspiration

GWF Ground-Water Flow

LPF Layer-Property Flow

MODFLOW modular finite-difference ground-water flow model

SFR2 Streamflow-Routing

SIP Strongly Implicit Procedure

USGS U.S. Geological Survey

UZF1 Unsaturated-Zone Flow

VSF Variably-Saturated Flow

VS2DT Variably-Saturated Two-Dimensional Flow and Transport Model
Abstract 

Percolation of precipitation through unsaturated zones is 

important for recharge of ground water. Rain and snowmelt at 

land surface are partitioned into different pathways including 

runoff, infiltration, evapotranspiration, unsaturated-zone 

storage, and recharge. A new package for MODFLOW2005 called the Unsaturated-Zone Flow (UZF1) Package 

was developed to simulate water flow and storage in the 

unsaturated zone and to partition flow into evapotranspiration 

and recharge. The package also accounts for land surface 

runoff to streams and lakes.

A kinematic wave approximation to Richards’ equation 

is solved by the method of characteristics to simulate vertical 

unsaturated flow. The approach assumes that unsaturated 

flow occurs in response to gravity potential gradients only 

and ignores negative potential gradients; the approach further 

assumes uniform hydraulic properties in the unsaturated zone 

for each vertical column of model cells. The Brooks-Corey 

function is used to define the relation between unsaturated 

hydraulic conductivity and water content. Variables used by 

the UZF1 Package include initial and saturated water contents, 

saturated vertical hydraulic conductivity, and an exponent 

in the Brooks-Corey function. Residual water content is 

calculated internally by the UZF1 Package on the basis of the 

difference between saturated water content and specific yield. 

The UZF1 Package is a substitution for the Recharge 

and Evapotranspiration Packages of MODFLOW-2005. 

The UZF1 Package differs from the Recharge Package in 

that an infiltration rate is applied at land surface instead of a 

specified recharge rate directly to ground water. The applied 

infiltration rate is further limited by the saturated vertical 

hydraulic conductivity. The UZF1 Package differs from 

the Evapotranspiration Package in that evapotranspiration 

losses are first removed from the unsaturated zone above 

the evapotranspiration extinction depth, and if the demand 

is not met, water can be removed directly from ground 

water whenever the depth to ground water is less than the 

extinction depth. The UZF1 Package also differs from the 

Evapotranspiration Package in that water is discharged directly 

to land surface whenever the altitude of the water table 

exceeds land surface. Water that is discharged to land surface, 

as well as applied infiltration in excess of the saturated vertical 

hydraulic conductivity, may be routed directly as inflow 

to specified streams or lakes if these packages are active; 

otherwise, this water is removed from the model. 

The UZF1 Package was tested against the U.S. 

Geological Survey’s Variably-Saturated Two-Dimensional 

Flow and Transport Model for a vertical unsaturated flow 

problem that includes evapotranspiration losses. This report 

also includes an example in which MODFLOW-2005 with 

the UZF1 Package was used to simulate a realistic surfacewater/ground-water flow problem that includes time and space 

variable infiltration, evapotranspiration, runoff, and groundwater discharge to land surface and to streams. Another 

simpler problem is presented so that the user may use the input 

files as templates for new problems and to verify proper code 

installation.

Introduction

Historically, hydrologic processes that occur near the land 

surface, such as infiltration and evapotranspiration (ET) have 

been modeled separately from regional ground-water flow 

processes. Watershed-runoff models simulate near-surface 

hydrology and non-physically-based shallow unconfined 

ground-water flow, but rarely simulate physically-based threedimensional regional ground-water flow. Similarly, regional 

ground-water modeling studies often simplify near-surface 

hydrologic processes such that recharge is specified to the 

model on the basis of external calculations independent of 

ground-water levels. Recently, however, regional ground-water 

flow models have been coupled to watershed-runoff models 

(Sophocleous and Perkins, 2000). Coupled models are used 

to simulate the effects of near-surface hydrology on regional 

ground-water flow by partitioning precipitation into various 

pathways among the atmosphere, surface, and subsurface. In 

this manner, recharge in regional ground-water models can 

be simulated in concert with changing climate. Moreover, 

ground-water flow and watershed-runoff modeling is improved 

by more realistically including the effects of ground-water 

seepage to land surface and streams, and by simulating the 

rejection of infiltration caused by high ground-water levels. 

Documentation of the Unsaturated-Zone Flow (UZF1) 

Package for Modeling Unsaturated Flow Between the 

Land Surface and the Water Table with MODFLOW-2005

By Richard G. Niswonger, David E. Prudic, and R. Steven Regan
Modeling flow through the unsaturated zone has been 

a major obstacle to coupling watershed-runoff and regional 

ground-water flow models. One approach is to use Richards’ 

equation to model both saturated and unsaturated flow 

(Freeze, 1971). Richards’ equation is highly non-linear and 

more difficult to solve than the conventional ground-water 

flow equation. Models that incorporate Richards’ equation 

require much finer spatial grids and time steps, which limits 

their applicability for simulating regional ground-water 

flow problems. Additionally, modeling both saturated and 

unsaturated flow with the numerical solution of Richards’ 

equation requires that unsaturated flow be simulated in 

three dimensions, which is often unnecessary for basinscale applications because the direction of unsaturated flow 

averaged over large grid cells is usually vertical (Mantoglou, 

1992; Chen and others, 1994). 

Another approach to coupling regional ground-water flow 

models to watershed-runoff models is to simulate unsaturated 

flow using Richards’ equation in only one dimension and 

three-dimensional ground-water flow. This approach has been 

used to model unsaturated flow over large areas (Pikul and 

others, 1974) and to develop the integrated surface-water and 

ground-water flow model MIKESHE (Refsgaard and Storm, 

1995). The approach used by Pikul and others (1974) and 

Refsgaard and Storm (1995) is similar to the approach used 

in developing the Unsaturated-Zone Flow (UZF1) Package 

for MODFLOW-2005, which solves the three-dimensional 

ground-water flow equation using finite-difference techniques 

(Harbaugh, 2005; fig. 1). The one-dimensional form of 

Richards’ equation is approximated by a kinematic-wave 

equation that is solved by the method of characteristics (Smith, 

1983). The method of characteristics solution for unsaturated 

flow precludes the need to develop a structured grid of 

the unsaturated zone for numerical stability and simplifies 

handling of the moving boundary defined by the interface 

between the water table and unsaturated zone.

Purpose and Scope

This report presents a new package for simulating 

vertical unsaturated flow in MODFLOW-2005 (Harbaugh, 

2005). The package is intended to provide an efficient means 

of simulating recharge in MODFLOW-2005 that considers 

the effects of flow, ET, and storage in the unsaturated zone; 

it is also intended to be used for coupling MODFLOW to 

precipitation-runoff models. The theoretical development 

of the kinematic-wave approximation to Richards’ equation 

and the method of characteristic solution was presented by 

Smith (1983) and Smith and Hebbert (1983). The method of 

characteristics solution for unsaturated flow is extended in 

this report to consider ET losses. The method was initially 

incorporated with MODFLOW-2000 through the simulation 

of unsaturated flow beneath intermittent streams without 

considering ET (Niswonger and Prudic, 2004). 

This report also describes how the UZF1 Package is 

coupled to MODFLOW-2005 and how the UZF1 Package 

calculates ground-water seepage to land surface. A comparison 

is presented between MODFLOW-2005 using the UZF1 

Package and the U.S. Geological Survey (USGS) VariablySaturated Two-Dimensional Flow and Transport Model 

(VS2DT; Healy, 1990) for a simple problem that involves 

vertical unsaturated flow with ET. Two test simulations also 

are described: the first test simulation shows results from a 

transient model of a 27 km2 watershed in which infiltration, 

ET, and ground-water discharge vary through time. The 

second test simulation is a simple problem that can be used to 

verify proper code installation and as a template for setting up 

new problems. 

A more complex method for simulating flow through 

the unsaturated zone is available for MODFLOW-2000. The 

Variably-Saturated Flow (VSF) Process (Thoms and others, 

2006) uses a finite-difference approximation to solve the threedimensional form of Richards’ equation. VSF replaces the 

standard Ground-Water Flow (GWF) Process in MODFLOW2000 (Harbaugh and others, 2000). The process is designed for 

problems where horizontal flow within the unsaturated zone is 

important. VSF includes packages that can simulate seepage 

faces, infiltration with ponding, bare soil evaporation, and 

plant transpiration. The VSF Process offers a more rigorous 

and more computationally demanding treatment of flow 

through the unsaturated zone than that provided by the UZF1 

Package.

Acknowledgments

The U.S. Geological Survey Ground-Water Resources 

Program and the Office of Ground Water funded this 

work. The original idea to develop the UZF1 Package for 

MODFLOW began in the summer of 1995 and was inspired 

by interest in coupling MODFLOW to a precipitation-runoff 

model. The development of this work began with many 

conversations with Jon P. Fenske and Arlen D. Feldman 

(retired) of the Army Corps of Engineers Hydrologic 

Engineering Center, Davis, California, and with Stanley A. 

Leake and Richard W. Healy of the U.S. Geological Survey. 

The Hydrologic Engineering Center provided the initial 

support to develop a plan for incorporating precipitation 

runoff models with MODFLOW (Fenske and Prudic, 1998). 

The authors are grateful to Steven L. Markstrom of the U.S. 

Geological Survey National Research Program for his help 

in the final design of UZF1, and the technical reviewers of 

this report, including Paul M. Barlow, Randall T. Hanson, 

Arlen W. Harbaugh, Devin L. Galloway, and Richard W. 

Healy of the U.S. Geological Survey.

2 Documentation of the Unsaturated-Zone Flow Package
Figure 1. One-dimensional unsaturated-zone flow coupled to three-dimensional 

ground-water flow.

NV_UZF1_fig01

Evapotranspiration (ET)

ET extinction depth

Deep water table

One-dimensional unsaturated flow

Three-dimensional ground-water flow

Land surface

Overland flow 

to streams and lakes

B

Evapotranspiration (ET)

ET extinction depth

Shallow water table

One-dimensional unsaturated flow

Three-dimensional ground-water flow

Land surface

Overland flow 

to streams and lakes

A

Introduction 3
Theory and Conceptualization of 

Unsaturated-Zone Flow

Vertical flow through a homogeneous unsaturated 

zone can be approximated with kinematic waves (Colbeck, 

1972; Smith, 1983; Smith and Hebbert, 1983; Charbeneau, 

1984). This approximation is made by simplifying Richards’ 

equation, which can be written in the vertical dimension as:

¶

¶ = ¶

¶ - = ¶

¶

¶

¶ - é

ë

ê

ê

ù

û

ú

ú

- q q) q q) t

q

z

i

z

D

z

( ( K i , (1)

where

q is the volumetric water content (volume of water per 

volume of rock);

q is the water flux (volume of water per time per unit 

area);

z is the elevation in the vertical direction (length);

D( )q is the hydraulic diffusivity (length squared per time);

K( )q is the unsaturated hydraulic conductivity as a function 

of water content (length per time);

i is the ET rate per unit depth (length per time per 

length); and

t is time.

If equation 1 is simplified to remove the diffusive term, 

assuming that the vertical flux is only driven by gravitational 

forces, the vertical flux owing to gravity is:

q K = - ( )q , (2)

where

q is shown as negative downward.

Substituting equation 2 into equation 1, neglectingD

z

( )q ¶q

¶

(the diffusive term), and assuming evaporation is removed 

from the soil profile instantaneously and all flow is downward 

vertical yields:

¶

¶ +

¶

¶ + =

q q

t

K

z

i ( ) 0 . (3)

Solution Method

The method of characteristics solution to equation 3 is 

obtained by rewriting equation 3 so that q is taken partially in 

terms of t and z:

¶

¶ +

¶

¶

¶

¶ = - q qq

q

t

K

z

i ( ) , (4)

and by introducing the equation of variation (Abbott, 1966):

¶

¶ +

¶

¶ = q q qt

dt

z

dz d . (5)

Equations 4 and 5 written in matrix form yields:

1 ¶

¶

é

ë

ê

ê

ê

ê

ù

û

ú

ú

ú

ú

¶

¶

¶

¶

é

ë

ê

ê

ê

ê

ê

ù

û

ú

ú

ú

ú

ú

=

éë

ê

ê

ê

ù

û

K

dt dz

t

z

i

d

( )q

q

q

q q

ú

ú

ú

. (6)

¶ ¶ q / t and ¶ ¶ q / z must be indeterminate along characteristic 

lines (Abbott, 1966). This occurs only when the determinant 

of the coefficient matrix in equation 6 and the determinants 

of the matrices obtained by substituting the right hand side 

vector for each column in the coefficient matrix are all zero, as 

described in equation 7:

1

0

¶

¶ =

K

dt dz

( )q

q , (7a)

1

0

-

=

i

dt dq , (7b)

- ¶

¶ = i K

d dz

( )q

q

q

0 . (7c)

Expanding the determinants yields:

dz

dt

K = v¶

¶ = ( ) ( ) qq

q , (8a)

d

dt

i q = - , (8b)

d

dz

i

v

q

q = -

( ) , (8c)

where

v( )q is the characteristic velocity restricted to the downward 

(positive z) direction (length per time) and 

q decreases according to the ET rate along a 

characteristic path.

4 Documentation of the Unsaturated-Zone Flow Package
The characteristic equations 8a, b, and c define the velocity 

of a wave, the change in water content of the wave with time, 

and the change in water content with depth behind the wave, 

respectively. This wave represents a wetting front in the 

unsaturated zone that was generated by a sudden increase in 

infiltration at land surface and that is being affected by ET. 

Without ET, the water content of the wave is constant over the 

wave profile. However, the inclusion of ET creates a 

non-linear slope ( ) d

dz

q in water content over the wave profile

if a non-linear function for K( )q is used, such as the 

Brooks-Corey function (Brooks and Corey, 1966). 

The derivative ¶

¶

K( )q

q is discontinuous over a sharp 

wetting front that results from neglecting the diffusive term 

in equation 1, such that the spatial derivative in equation 3 is 

not defined in the absence of hydraulic diffusion. An analytic 

solution for the velocity may still be derived by considering 

the effects of diffusion and substituting an equivalent sharp 

wetting front of equal mass (Smith, 1983; Charbeneau, 1984). 

The solution to equation 1 for a wetting front that considers 

hydraulic diffusion can be found by integrating over a control 

volume containing a single wetting front (Charbeneau, 1984):

d

dt

dz K D

z z

z

z

z

q q q q+ -

¶

¶ ò ( ( ) ( ) ) =1

2

1

2

0 , (9)

where

z1 and z2 are points above and below the wetting front at

distances far enough that ¶

¶ »q

z

0, respectively.

D( )q can be neglected because ¶

¶ »q

z

0 at z1 and z2 such that:

d

dt

dz K K

z

z

ò q q + - ( )2 1 ( ) q = 0

1

2

, (10)

where

K( ) q1 and K( ) q2 are the values of K( )q at depths z1 and z2.

Integrating over a profile containing a sharp front with 

equivalent mass gives:

q q q

z

z

dz f f z z z z

1

2

ò = - 1 1 + - 2 2 ( ) ( ), (11)

where

zf is the depth of the sharp front.

Combining equations 10 and 11 gives:

dz

dt

u f K K = = s

( )- ( )

- ( , q q ) q q

q q 1 2

1 2

1 2

, (12)

where

us is the velocity of a sharp wetting front (length per time);

q1 is the volumetric water content above a depth zf ; and

q2 is the volumetric water content below a depth zf .

Representation of Wetting and Drying Fronts

An increase in the infiltration rate will cause a wetting 

front to form, which is represented by a lead wave. A decrease 

in the infiltration rate will cause a drying front to occur, which 

is represented by a trailing wave. Thus, waves are used to 

represent both wetting and drying fronts. Attenuation of a lead 

wave occurs as a trailing wave of higher velocity overtakes 

it. When a trailing wave overtakes a lead wave, the water 

content of the lead wave becomes equal to the water content 

of the trailing wave. Consequently, this reduces the velocity 

of the lead wave. Conversely, when a lead wave overtakes 

a trailing wave or another lead wave of lower velocity, the 

overtaken wave is removed, and the water content and flux 

of the uppermost lead wave are maintained, resulting in 

rewetting. The process of waves intercepting each other is 

termed a shock, and this process is important to the method of 

characteristics because wave properties, such as water content 

and flux, become discontinuous at the instant a shock occurs 

(Abbott, 1966).

In contrast to a wetting front that stays sharp, trailing 

waves elongate with time owing to gravity. Consequently, 

trailing waves must be divided into a series of incremental 

waves or be represented by a function that describes internal 

drainage with time. In either case, internal drainage is 

determined based on the water content and relative location 

of points along a trailing wave. For example, points A and B

shown in figure 2 represent two points along the trailing wave 

where the velocities may be determined. The velocity at points 

A and B is calculated from the characteristic solution (equation 

8a) and the ratio of their velocities is equal to the ratio of their 

depths:

z

z

u

u

dK d

dK d

B

A

B

A

B

A

= = =

=

( )

( )

( )/

( )/

q

q

q q

q q

q q

q q

, (13)

where

Theory and Conceptualization of Unsaturated-Zone Flow 5
zA is the depth of the trailing wave increment below land 

surface represented by point A in figure 2, and 

zB is the depth of trailing wave increment below land 

surface represented by point B in figure 2.

In contrast to a lead wave that forms a discontinuity 

in the water content profile, the derivative ¶ ¶ K( )q q / is 

continuous over a trailing front. Accordingly, the BrooksCorey unsaturated hydraulic conductivity function can be 

used to evaluate ¶ ¶ K( )q q / . The Brooks-Corey function for 

unsaturated hydraulic conductivity can be expressed as:

K Ks r

s r

( )q q q

q q = --

é

ë

ê

ê

ù

û

ú

ú

e

, (14)

where

Ks is the saturated hydraulic conductivity;

qr is the residual water content;

qs is the saturated water content; and

e is the Brooks-Corey exponent.

The drainable water from an unconfined or water-table 

aquifer, referred to as the specific yield or Sy, is equal to the 

porosity minus the specific retention ( ) Sr , where Sr is the 

volume of water per unit volume of rock that is retained when 

the rock is drained by gravity (Meinzer, 1923, p. 28). Specific 

yield is used by MODFLOW-2005 to calculate the amount 

of storage within an unconfined aquifer. Thus, continuity 

between the unsaturated zone and unconfined aquifers in 

MODFLOW-2005 is maintained through the specific yield by 

approximating S S r r y s r with q q , and with - q . Residual 

water content generally is smaller than Sr, such that caution 

should be used when representing unsaturated flow parameters 

in the UZF1 Package with measured values.

The velocity of the deepest point along a trailing wave 

is determined by taking the derivative of equation 14 with 

respect to q:

v

K

S S

s

y

r

y

( )q q q = - é

ë

ê

ê

ê

ù

û

ú

ú

ú

- e e 1

, (15)

where

v and q are the velocity and water content, respectively, at 

the deepest point along a trailing wave.

A relation between the deepest point along a trailing wave 

and all other points along a trailing wave is determined by 

substituting the velocity defined by equation 15 into equation 

13 and results in:

z z r

r

( )q q q

q q = --

é

ë

ê

ê

ù

û

ú

ú

-

0

0

e 1

, (16)

where

z( )q is the depth of a point on a trailing wave;

q is the water content of a point on a trailing wave; and

q0 and z0 are the water content and depth of the deepest 

point along a trailing wave, respectively.

For example, the depths of points along a trailing wave may 

be generated by choosing water content values between q2

and q3 shown in figure 2 and by using equation 16 to calculate 

the depth at a particular time following a decrease in the 

infiltration rate. Thus, equation 15 may be used to calculate 

the velocity and the corresponding depth of the front of a 

trailing wave and then equation 16 can be used to determine 

the drainage profile above the trailing front.

During a pulse of infiltration, a lead wave’s velocity 

and water content will decay during a subsequent period of 

less infiltration as trailing waves overcome the lead wave. 

The spreading of a trailing wave caused by gravity is defined 

by equations 15 and 16; however, these equations are not 

Figure 2. A wetting front moving through a uniform column 

of unsaturated material affected by a decrease in surface flux, 

and results from a kinematic-wave solution of the wetting front 

represented by leading and trailing waves.

(Modified from Smith and Hebbert, 1983.)

6 Documentation of the Unsaturated-Zone Flow Package

Water table

Qi Q2 Q1

Trailing wave 

represented as a 

step function at 

time t2

Actual wetting front 

at time equals t1

Actual wetting front 

at time equals t2

Leading wave with

water content Q

 

at time t2

 Leading wave with

water content Q

at time t1

WATER CONTENT (Q)

DEPTH (z)

A

B

Q
applicable after a trailing wave has intercepted a lead wave. 

Analytic relations can be derived for routing trailing waves 

during interception with a lead wave but these equations are 

difficult to solve. A simpler approach is to divide the trailing 

wave into steps or increments and calculate velocities of each 

increment based on a finite-difference approximation (Smith, 

1983):

v

K K = - - D

D

( )q q( )q

q

, (17)

where

Dq is the change in water content between two adjacent 

locations along a trailing wave.

The UZF1 Package relies on equations 15 and 16 to route 

trailing waves until the trailing waves intercept another lead 

wave, at which point the trailing wave is discretized into 

multiple waves that are routed according to equation 17.

The finite-difference approximation of the true trailing 

wave as depicted in figure 2 does not conserve mass 

perfectly when a trailing wave intercepts a lead wave and 

vice versa. However, for problems involving rapidly varying 

infiltration rates, mass is conserved within an error of less 

than 0.05 percent when each trailing wave is represented 

by 15 increments, although 10 increments are sufficient for 

most problems to achieve a mass balance error of less than 

0.05 percent. The user specifies the number of increments 

used to represent a trailing wave with the UZF1 Package input 

parameter NSTRAIL.

As described, waves are generated by changes in the 

infiltration rate and will have differing velocities leading to the 

possibility of shocks. Determination of the timing of shocks is 

required for modeling because the wave properties 

(v, q) become discontinuous in time and have to be 

recalculated in order to continue the solution. Thus, the UZF1 

Package routinely calculates the shortest time for any two 

waves to intersect. Waves are routed continuously between 

wave intersections, while wave properties remain constant. 

Consequently, the UZF1 Package requires a unique time 

discretization that is specific to the problem being solved; 

the UZF1 Package time steps are different than those used in 

MODFLOW-2005.

ET in the Unsaturated Zone

Evaporation and uptake by roots in the unsaturated zone 

can cause negative potential gradients, such that diffusive 

forces can be important. Evaporation can cause water to move 

upward in the unsaturated zone by drying out the soil at land 

surface. These conditions cannot be modeled with kinematic 

waves because diffusive forces are neglected. However, 

infiltrating water affected by evaporation and uptake by roots 

can be approximated with kinematic waves during relatively 

wet conditions because negative potential gradients are less 

important. This approach assumes that evaporation and uptake 

by roots can be grouped together as ET, and that they occur as 

an instantaneous loss of water over a depth interval equal to 

the depth of root uptake. 

Assuming that i is independent of q, equations 8b and 8c 

can be integrated to find:

q q t t +t = -it , (18a)

z

K K iz

i 2 = 2 1 - - 1-

[ (q q ) ( )] , (18b)

where

qt is the water content at the head of a lead wave above 

the ET extinction depth after time t, and

qt+t is the water content at the head of a lead wave after 

time t +t subject to the ET rate per unit depth i.

ET is specified as a rate (length per time) within the UZF1 

Package (input variable PET) and is converted internally to 

a rate per unit depth by dividing the specified ET rate by the 

extinction depth (EXTDP). 

Equation 18a is used to calculate the water content at 

the head of a lead wave through time and equation 18b is 

used to calculate the depth of points along the leading-wave 

profile between q q = q and q q = t+t, where qq is the water 

content at land surface corresponding to the infiltration rate q, 

and q1 and q2 are the water contents at arbitrary points along 

the lead-wave profile with q2 being deeper than q1 (fig. 3). 

The water content of a trailing wave is decreased according 

to the ET rate and equation 18a. If the water table elevation 

is above the extinction depth and the ET demand is not met 

by the unsaturated zone, then ET is removed directly from 

ground water. ET is removed from ground water using the 

same method implemented in the MODFLOW ET Package 

(McDonald and Harbaugh, 1988, p. 10-1). 

Figure 3. Kinematic-wave approximation of a wetting front moving 

through a uniform column of unsaturated material affected by a constant 

evapotranspiration rate of 4.35 millimeters per day.

8

6

4

2

0

0.20 0.21 0.22 0.23 0.24 0.25

WATER CONTENT, IN CUBIC METER PER CUBIC METER

dQ

v dz

d i

dz

Q - = Q	

DEPTH, IN METERS BELOW LAND SURFACE

Theory and Conceptualization of Unsaturated-Zone Flow 7
Infiltration Boundary

The specified infiltration rate is converted to water 

content in order to apply the characteristic solutions for 

unsaturated flow. A relation between water content and 

infiltration rate can be determined assuming the flux is equal 

to the unsaturated hydraulic conductivity using equations 2 

and 14 (Charbeneau, 1984):

q q = æ

è

ç

ç

çç

ö

ø

÷

÷

÷

÷ + < £ q

K

in S q K

s y r in s

1

0

/

( )

e

 , (19a)

q q = < s s K qin , (19b)

where

q is the water content of a wave generated from 

infiltration, and

qin is the infiltration rate (length per time).

The water content is set to the saturated water content when 

the infiltration rate specified in the UZF1 Package input file 

exceeds the saturated vertical hydraulic conductivity. If the 

specified infiltration rate is greater than the saturated vertical 

hydraulic conductivity, then the difference is multiplied by 

the model cell plan-view area and this volumetric rate of 

water can be added to a user-specified stream reach or lake by 

setting the UZF1 Package variable IRUNFLG to greater than 

zero. Equation 19 ignores the time to ponding and infiltration 

rates that exceed the saturated vertical hydraulic conductivity 

(Smith and Parlange, 1978). 

An algorithm was written in FORTRAN using 

equations 12, 14, 15, 16, 17, 18 and 19 to simulate vertical 

unsaturated flow and ET. This algorithm keeps track of the 

location, water content, velocity, and flux of all waves in the 

unsaturated zone through time and the interaction of waves 

overtaking one another. Additionally, changes in the thickness 

of the unsaturated zone owing to changes in the elevation of 

the water table are considered with respect to waves reaching 

the water table and contributing to ground-water storage.

Coupling UZF1 to MODFLOW-2005

Flow from the unsaturated zone as recharge to the 

water table in an unconfined aquifer is subtracted from the 

right-hand side of the system of equations that are solved by 

MODFLOW-2005: 

AX = - B QUZF1 , (20)

where

A is a matrix containing the coefficients of the 

conductance equations (HCOF array) that are 

solved by MODFLOW-2005;

X is a one-dimensional vector containing the groundwater heads that are solved by MODFLOW-2005;

B is a one-dimensional vector containing all known 

terms in the conductance equations that are not 

multiplied by unknown head values (RHS array); 

and 

QUZF1 is the volumetric rate (volume per time) of recharge 

computed in a given model cell from the UZF1 

Package.

Recharge that is simulated using the UZF1 Package is 

dependent on the water table elevation (ground-water head). 

Thus, the UZF1 Package is coupled to MODFLOW-2005 

within the iteration loop. However, because ground-water 

recharge is dependent on the amount of storage and the 

downward flux rate in the unsaturated zone, there is no 

general equation that defines the relation between recharge 

and ground-water head, such as the conductance equation in 

the Streamflow-Routing (SFR1) Package (Prudic and others, 

2004). Recharge from the UZF1 Package does not affect the 

HCOF array. 

Based on the test simulations presented in this report, 

the UZF1 Package does not cause additional instabilities or 

slow convergence of MODFLOW-2005 when ground water 

is not discharging to land surface. Ground-water levels that 

rise above land surface add an additional head dependency 

when using the UZF1 Package owing to the calculation of 

ground-water discharge to land surface, and consequently, 

convergence may be slowed significantly when ground 

water is discharging to land surface from many model cells. 

MODFLOW-2005 simulations using the UZF1 Package may 

increase simulation time compared to using the Recharge 

Package owing to the additional calculations required to 

simulate flow through the unsaturated zone.

An infiltration rate at land surface (length per time) is 

specified for each vertical column in the active model domain 

(appendix 1). The infiltration rate is multiplied by the planview area of the cell to produce a volumetric flow rate (volume 

per time). Three options are used in the UZF1 Package to 

determine where flow in the unsaturated zone recharges 

ground water (UZF1 variable NUZTOP). The three options 

are those used in the Recharge Package (McDonald and 

Harbaugh, 1988, p. 7-6; Harbaugh and others, 2000, p. 67-68). 

The first option specifies that all recharge is to the top model 

layer only; the second option allows for recharge to be 

specified for a particular layer in each vertical column; and the 

third option allows for recharge to the uppermost active cell 

8 Documentation of the Unsaturated-Zone Flow Package
in a vertical column. When the last option is used in the UZF1 

Package, the program will continue to route water through 

the unsaturated zone until it reaches the water table of the 

uppermost active cell. However, the hydraulic properties of the 

unsaturated zone will remain constant even if the unsaturated 

zone spans more than one cell in the vertical direction. If there 

are no active cells in a column beneath land surface, then 

infiltration and recharge do not occur and a warning statement 

is printed to the output file. If a cell becomes active during a 

simulation and the corresponding element in the UZF1 array 

IUZFBND is nonzero, then infiltration and recharge will start 

to occur. 

As stated previously, mass balance is achieved by 

approximating the residual water content variable with the 

specific retention computed as the saturated water content 

minus the specific yield. Thus, if the water table rises into the 

unsaturated zone, the aquifer will yield an amount of water 

equal to the water table rise multiplied by the difference 

between the water content and the specific retention. 

Alternately, if the water table falls due to ground-water 

discharge or pumping at some other location, then the water 

content in the unsaturated zone between the old and new water 

table is set equal to the specific retention. 

Recharge occurring during a rising water table over a 

time step is calculated by first routing waves through the 

unsaturated zone using the base of the unsaturated zone as 

the water-table elevation at the end of the previous time step. 

Flow across the base is accumulated for the time step and the 

volume of water in the unsaturated zone through which the 

water table rose over the time step is added to the volume that 

flowed past the base. This sum is divided by the time step to 

obtain a volumetric recharge rate. Once the flow equation has 

been solved for the time step, a new base of the unsaturated 

zone in each active cell is computed. Recharge occurring 

during a declining water table is computed based on the flux 

of water from the unsaturated zone during the time step using 

the base of the unsaturated zone predicted at the previous 

iteration or time step. 

Ground-Water Discharge to Land Surface

The UZF1 Package also allows for ground-water 

discharge to the land surface whenever the water table in a cell 

is higher than the land-surface elevation. The volumetric rate 

of ground-water discharge to land surface is calculated on the 

basis of the following equation:

Q A K h celtop celthk gw cel = - l v ( )/ . 0 5 ) , (21)

where

h is the water-table elevation;

celtop is the land-surface elevation;

Acell is the plan-view area of the model cell (equal to 

the column length times the row length of the 

model grid);

Kv is the vertical hydraulic conductivity of the model 

cell; and

celthk is the thickness of the model cell.

The value that is subtracted from the HCOF array in 

MODFLOW-2005 (corresponds to the A matrix in equation 20) 

is equal to A K celthk cell v /(0 5. ). The value subtracted from the 

RHS array in MODFLOW-2005 (corresponds to the B matrix in 

equation 20) is equal to A K celtop celthk cell v /(0 5. ). 

The UZF1 Package includes an option to specify a twodimensional array of values (UZF1 variable IRUNBND) that 

identify where ground-water discharge to land surface will be 

added to a stream or lake if the SFR2 and LAK Packages are 

active. Additionally, the UZF1 Package variable IRUNDFLG

must be greater than 1 for water to be added to streams or 

lakes. Ground-water discharge to land surface and infiltration 

in excess of the saturated vertical hydraulic conductivity are 

added instantaneously to specified streams and lakes. The 

SFR2 Package is documented by Niswonger and Prudic (2006) 

and the LAK3 Package is documented by Merritt and Konikow 

(2000). 

A water budget for the unsaturated zone is tracked 

independently of the ground-water budget within 

MODFLOW-2005. The mass balance for the unsaturated zone 

accounts for all water between the land surface and the water 

table. The mass balance for the saturated zone is accounted for 

in the ground-water budget. When the water table rises above 

the top of a model cell, then the unsaturated-zone mass balance 

ceases to be calculated for that cell, and unsaturated storage 

and the change in unsaturated-zone storage for that cell are set 

to zero. 

There are no procedural limitations regarding stress period 

and time step lengths associated with the UZF1 Package. 

However, because the UZF1 Package delays recharge to the 

water table, changes in stress do not necessarily correspond to 

the beginning of a stress period. Consequently, caution should 

be used when simulating time steps that increase during a 

stress period (specifying the variable TSMULT>1.0 in the data 

input for the Discretization file; Harbaugh, 2005) because 

water percolating through the unsaturated zone may reach the 

water table at the end of a stress period, when the time step 

length is a maximum.

Theory and Conceptualization of Unsaturated-Zone Flow 9
Steady-State Simulations

Steady state can be specified for the initial stress period 

of a simulation when using the UZF1 Package. When steady 

state is specified for the initial stress period, the ground-water 

recharge rate computed in UZF1 is set equal to the infiltration 

rate and ET is not removed from the unsaturated zone. The 

water content is constant within each model cell during steady 

state. The specified ET rate is subtracted from the infiltration 

rate to calculate a steady-state recharge rate. ET can also 

occur from ground water when the water table is above the ET 

extinction depth. A uniform water content profile is computed 

for each active cell, which is used automatically as the initial 

conditions for a following transient stress period.

Steady-state water content profiles may be printed to a 

file for specified locations. Initial water content values are not 

required when the initial stress period is specified as steady 

state. Thus, initial water-content values must not be specified 

within the UZF1 input file when the first stress period is 

steady state. Another option for developing initial water 

content profiles is to run the model for several years while 

repeating a cycle of infiltration rates. This step may be carried 

out subsequent to a steady-state simulation. 

Assumptions and Limitations

The kinematic-wave approximation does not account 

for capillary-induced infiltration at the onset of wetting 

when negative pressure gradients are large relative to gravity 

potential gradients. However, infiltration approaches the 

unsaturated hydraulic conductivity when the wetting front 

has advanced sufficiently far beneath the surface (Mein and 

Larson, 1973; Childs and Bybordi, 1969). Negative pressure 

gradients also may result in lateral and vertical redistribution 

in the unsaturated zone. Because negative pressure gradients 

are ignored, UZF1 likely under predicts the advancement of 

the wetting front at early times.

Because capillary-induced infiltration is neglected, there 

is no consideration for the period prior to ponding when the 

infiltration should be represented by a constant flux equal to 

the rainfall rate (Smith and Parlange, 1978), rather than the 

constant water content assigned at the boundary (equation 19). 

Thus for greater rainfall or snowmelt rates, the model will 

underestimate infiltration during the onset of infiltration. 

UZF1 does not simulate a capillary fringe above the water 

table and it assumes that ground water is instantaneously 

released from or taken into storage. 

The initial version of this package (UZF1) does not 

simulate unsaturated flow through multiple layers of varying 

hydraulic properties. The model can simulate unsaturated flow 

through inactive cells to deeper active model cells. However, 

the hydraulic properties of an unsaturated column between the 

land surface and the uppermost active cell are uniform during 

the simulation and are always equal to their values at the onset 

of the simulation.

Comparison of UZF1 and VS2DT 

Simulations

The kinematic-wave approximation used in the UZF1 

Package was compared to the USGS’s Variably-Saturated 

Two-Dimensional Flow and Transport (VS2DT) Model 

(Healy, 1990). VS2DT solves Richards’ equation numerically 

in terms of pressure head. The results were limited to the 

simulation of vertical flow through a column.

Description

Percolation through the unsaturated zone was simulated 

with VS2DT as a 30-m long vertical column of materials 

having hydraulic properties typical of a sandy clay loam, 

whereas percolation through the unsaturated zone was 

simulated with the UZF1 Package and MODFLOW-2005 

for a single cell (row 10, column 4) from the model in test 

simulation 2, described herein. The initial water table in row 

10, column 4 was 22 m below land surface and the water table 

declined less than 1 m during the simulation period. The same 

variables that describe the unsaturated zone were used in both 

model simulations (table 1). The wetting front did not reach 

the water table for both simulations because the comparison 

was done to test differences in the propagation of the wetting 

front with varying ET demand rates. VS2DT simulates 

evaporation separate from transpiration by root uptake; only 

transpiration was used in the simulation. Total simulation time 

was 17.7 days. A variable time step was used in the simulation 

with VS2DT, whereas a constant time step was used in the 

UZF1/MODFLOW-2005 simulation. 

Variable Value

Saturated vertical hydraulic conductivity (m/s) 4×10-6

Saturated water content (m3 of water per m3 of aquifer) 0.4

Residual water content (m3 of water per m3 of aquifer) 0.2

Brooks-Corey exponent (dimensionless) 3.5

Air entry pressure (m)1 -0.05

Infiltration rate (m/s) 2×10-6

Evapotranspiration extinction depth (m) 2

Evapotranspiration rate (m/s) 

 Simulation 1 5×10-8

 Simulation 2 5×10-7

 Simulation 3 1×10-6

1Value was not used in the simulation with the Unsaturated-Zone Flow 

(UZF1) Package in MODFLOW-2005.

Table 1. Unsaturated zone variables used in comparison of results from 

the Unsaturated-Zone Flow (UZF1) Package in MODFLOW-2005 with those 

from the Variably-Saturated Two-Dimensional Flow and Transport (VS2DT) 

Model.

[Abbreviations: m, meter; m/s, meter per second; m3, cubic meter] 

10 Documentation of the Unsaturated-Zone Flow Package
A constant infiltration rate was applied at the upper 

boundary of the column and was 2.0 ×10-6 m/s. The 

infiltration rate was half the saturated vertical hydraulic 

conductivity (table 1). The constant infiltration rate at land 

surface was converted to water content within UZF1 on the 

basis of equation 19. The initial water content was set equal 

to the residual water content throughout the unsaturated zone. 

The ET extinction depth was set at a depth of 2 m below land 

surface, and an ET demand rate of 5.0 ×10-8 m/s, 5.0 ×10-7

m/s, and 1.0 ×10-6 m/s was specified for simulations 1, 2 and 3, 

respectively (table 1). Although the two highest ET demand 

rates used in this example are too high to represent any actual 

conditions, they serve to demonstrate the model’s accuracy 

during extreme situations.

Results

Advancement and magnitude of the wetting front for 

all three simulations of differing specified ET rates were 

nearly the same for both models (fig. 4). The kinematicwave approximation in UZF1 slightly under-predicted the 

advancement of the wetting front at early time and overpredicted advancement at later time. The effect was more 

pronounced when the ET demand rate was 1.0 ×10-6 m/s 

(fig. 4C). Discrepancies between UZF1 and VS2DT for each 

simulation were related to neglecting the negative pressure 

gradients with the kinematic-wave approximation.

Figure 4. Comparison of results from the kinematic-wave approximation used to simulate unsaturated 

flow (UZF1) with a two-dimensional finite-difference solution of Richards’ equation (VS2DT; Healy, 1990) 

assuming evapotranspiration demand rates of: (A) 5 × 10-8 meter per second, (B) 5 × 10-7 meter per 

second, and (C) 1 × 10-6 meter per second.

0

5

10

15

20

0.15 0.2 0.25 0.3 0.35 0.4

WATER CONTENT, IN CUBIC METER PER CUBIC METER

7.5 days

2.5 days

12.6 days

17.7 days

2.5 days

7.5 days

12.6 days

17.7 days

2.5 days

7.5 days

12.6 days

17.7 days

B

A

C

0

5

10

15

20

0

5

10

15

20

UZF1

VS2DT

DEPTH, IN METERS BELOW LAND SURFACE

Comparison of UZF1 and VS2DT Simulations 11
Test Simulations

Two hypothetical test simulations are used to illustrate 

the capabilities of the UZF1 Package. The first is an 

application of the model to a watershed that encompasses an 

area of 27 km2. This simulation was designed to demonstrate 

all the capabilities of UZF1. However, because the first test 

simulation involved relatively large data input and model 

output files that are not easily printed, a second simpler 

test simulation is described that can be used as a guide for 

data input and model output structure and to verify proper 

installation of the UZF1 Package in MODFLOW-2005. 

Detailed results are described for the first test simulation in 

this document, but the data input and model output are not 

included as part of the printed document. However, data input 

and model output for the second test simulation are included 

in appendix 2 along with a summary of the unsaturated-zone 

and ground-water budgets in the discussion of the second test 

simulation. Both test simulations and this report are available 

at http://water.usgs.gov/software/.

Test Simulation 1: Sagehen Watershed

The first test simulation is a conceptual yet realistic 

simulation of the Sagehen watershed (fig. 5), a USGS 

Hydrologic Benchmark Network Basin located on the 

eastern slope of the northern Sierra Nevada, near Truckee, 

California (Mast and Clow, 2000). The basin drains an area 

of 27 km2 and ranges in altitude from 1,926 to 2,663 m above 

mean sea level. The areally averaged annual precipitation is 

about 970 mm and the annual hydrograph is dominated by 

snowmelt. The Sagehen watershed comprises volcanic rocks 

overlain by a veneer of alluvium. 

The test simulation consisted of a single model layer 

with 73 rows and 81 columns in which all cells had a constant 

width and length equal to 90 m (fig. 5). The top elevations of 

model cells were set equal to estimated land-surface altitudes 

for each cell; the elevation of the bottom of the aquifer 

(bottom of model cells) ranged from 1,830 to 2,480 m above 

mean sea level. The simulation period was 1-year beginning 

on December 1 and ending November 30. The simulation 

period was divided into 365 1-day stress periods. The first 

stress period was steady state, and subsequent stress periods 

were transient. One time step was simulated for each stress 

period. No-flow conditions were simulated across the bottom 

and sides of the model except for three cells beneath and 

adjacent to the stream at the basin outlet, which were specified 

as constant head cells.

Several packages and optional files were used in the 

simulation of the Sagehen watershed (table 2). The LayerProperty Flow (LPF) Package (Harbaugh and others, 2000, 

p. 59) was used and both horizontal and vertical hydraulic 

conductivity values were specified. Because one layer was 

specified in the model, the vertical hydraulic conductivity of 

the aquifer only affected ground-water seepage to land surface. 

All cells were specified as convertible (unconfined). The 

horizontal hydraulic conductivity ranged from 0.005 m/d on 

the ridges to 0.3 m/d in the valleys (table 3). A lower hydraulic 

conductivity was specified on the ridges where volcanic rocks 

are near the surface (fig. 5). The hydraulic conductivity within 

each cell was assumed isotropic. Specific storage was set to 

5 ×10-7 m-1, and the specific yield was specified as 0.05 on the 

ridges and 0.25 in the valleys near streams. 

The vertical hydraulic conductivity of the unsaturated 

zone had the same values as the horizontal and vertical 

hydraulic conductivity used in the LPF Package (table 3). A 

constant Brooks-Corey exponent of 4 and a constant saturated 

water content of 0.3 were assigned to all cells.

The streams were represented by 15 segments made up 

of 201 reaches using the Streamflow-Routing (SFR2) Package 

(Niswonger and Prudic, 2006) and included the simulation 

of unsaturated flow beneath streams. All stream reaches 

were set to a constant width of 3 m (ICALC=1) and ranged 

in length from 30 to 114 m. The streambed and unsaturated 

zone hydraulic conductivity were equal to 0.3 m/d (table 3). 

Unsaturated-zone properties were constant among all stream 

reaches. Saturated water content beneath the stream was 0.30. 

The Brooks-Corey exponent was 3.5. 

The distribution of hydraulic conductivity was created 

initially on the basis of the surface geology. The distribution 

of hydraulic conductivity was adjusted during a steady-state 

calibration that assumed that the spatially varying groundwater recharge was proportional to the distribution of mean 

annual precipitation. Precipitation was distributed according 

to 860 mm below an altitude of 2,100 m and increased linearly 

above 2,100 m to a maximum of 1,140 mm at an altitude of 

2,600 m. The range in ground-water recharge for the steadystate simulation was determined by approximating the mean 

daily discharge at the outlet of the Sagehen watershed in early 

December. The steady-state simulation resulted in a calculated 

water table that was as much as 80 m below land surface 

along the ridges and at, or slightly above, land surface in the 

lowest parts of the valleys next to streams (fig. 6). Although 

no observation wells have been drilled on the ridges in the 

watershed, the maximum depth to ground water was based 

on depths measured in wells elsewhere in the northern Sierra 

Nevada.

12 Documentation of the Unsaturated-Zone Flow Package
Figure 5. Distribution of hydraulic conductivity from ground-water model used for test simulation 1.

Elevation 1,935 meters

Elevation 2,653 meters

Sagehen watershed North

Sierra Nevada

California

Nevada

Model boundary

0.005

0.01

0.03

0.1

0.3

Hydraulic Conductivity—

meter per day

Sagehen

Creek

EXPLANATION 

Independence Lake

Shaded relief base from USGS 10-meter National Elevation Data, 

illumination from the northwest at 45°.

Albers Equal-Area Conic projection 

Standard parallels 29°30'N and 45°30'N, central meridian 120°00'W

North American Datum of 1983

Perspective Information:

Altitude is 7,300 meters above land surface

Viewing angle is 25 degrees

Vertical exaggeration times 2

Test Simulations 13
The initial water content of the unsaturated zone was 

calculated on the basis of the simulated steady-state recharge 

rate, which equaled the infiltration rate at land surface. An ET 

extinction depth of 2.5 m was assigned to all cells (table 3) 

and ET from the upper 2.5 m was simulated on the basis 

of the assigned ET demand rate. Water was first removed 

from the unsaturated zone and then from the saturated zone 

only when the water table was within 2.5 m of land surface. 

The water content for ET extinction was set to a constant of 

0.1002. Infiltration and ET demand rates were specified for 

each stress period and were estimated assuming infiltration 

and ET were at a minimum from December 1 to mid-April, 

Table 3. Hydraulic properties and other selected variables used in the Layer-Property Flow (LPF) and the Unsaturated-Zone Flow (UZF1) Packages for 

test simulation 1.

[Abbreviations: m, meter; m/d, meter per day; m3, cubic meter]

Variable Value

Variables assigned to the Layer-Property Flow (LPF) Package

Horizontal hydraulic conductivity of aquifer (m/d) 0.005 to 0.3

Vertical hydraulic conductivity of aquifer (m/d) 0.005 to 0.3

Specific storage (per m of aquifer) 5 ×10-7

Specific yield (m3 of water per m3 of aquifer) 0.05 to 0.25

Variables assigned to the Unsaturated-Zone Flow (UZF1) Package 

Vertical hydraulic conductivity of the unsaturated zone beneath land surface (m/d) 0.005 to 0.3

Brooks-Corey exponent for unsaturated zone beneath land surface (unitless) 4.0

Saturated water content of unsaturated zone beneath land surface (m3 of water per m3 of unsaturated zone) 0.3

Initial water content of unsaturated zone beneath land surface (m3 of water per m3 of unsaturated zone) 0.2

Steady-state infiltration rate (m/d) 0.0004 to 0.00104

Evapotranspiration extinction depth (m) 2.5

Evapotranspiration extinction water content (m3 of water per m3 of unsaturated zone) 0.1002

Variables assigned to the Streamflow-Routing (SFR2) Package

Hydraulic conductivity of streambed (m/d) 0.3

Streambed thickness (m) 1.0

Vertical hydraulic conductivity of unsaturated zone beneath streams (m/d) 0.3

Brooks-Corey exponent for unsaturated zone beneath streams (unitless) 3.5

Saturated water content of unsaturated zone beneath streams (m3 of water per m3 of unsaturated zone) 0.3

MODFLOW-2005 Packages and Files

Basic (BAS6) Package

Output Control (OC) Option file

Discretization (DIS) file

Layer-Property Flow (LPF) Package

Preconditioned Conjugate-Gradient (PCG) Solver Package 

Unsaturated-Zone Flow (UZF1) Package

Streamflow-Routing (SFR2) Package

Gage (GAG) Package

Table 2. MODFLOW-2005 Packages and files used for test simulation 1.

14 Documentation of the Unsaturated-Zone Flow Package
Figure 6. Distribution of steady-state ground-water depths used for test simulation 1.

N

Depth to water 

 table, in meters 

 below land 

 surface

Greater than 25

12.5

0.0

Sagehen Creek watershed

Model boundary

0 1 2 4 KILOMETERS

infiltration increased to a maximum during May, and the daily 

ET demand increased to a maximum during June (fig. 7). 

The specified infiltration rates varied in the model and were 

three times greater at the tops of the ridges as compared to the 

valley lowlands (fig. 8). The range in the applied infiltration 

rates were estimated on the basis of the average difference 

between the two precipitation gages located in the upper and 

lower portions of the watershed.

The Preconditioned Conjugate-Gradient (PCG) Package 

was used to solve for ground-water head throughout the model 

domain. The head-closure criterion used in the transient 

simulation was 0.0009 m and the flow-closure criterion 

was 0.005 m3/d. The absolute mass-balance error for the 

unsaturated zone ranged from 0 to 0.08 percent for all stress 

periods and the saturated zone ranged from 0 to 0.22 percent. 

The cumulative mass-balance error at the end of the 365-

day simulation was 0 percent for the unsaturated zone and 

0.10 percent for the saturated zone. 

Test Simulations 15
Results from the transient simulation of ET and 

percolation through the unsaturated zone indicate that recharge 

patterns are a result of the cumulative effects of ET, spatially 

varying recharge, hydraulic conductivity, topography, and 

thickness of the unsaturated zone. Model calculated groundwater recharge summed over the entire basin indicates two 

distinct patterns: a nearly immediate response to peaks in 

infiltration in the valley lowlands where the unsaturated 

zone is thin, and a diffuse response that varies slowly over 

several months in upland areas where the unsaturated zone is 

relatively thick (fig. 9). 

The simulated volumetric ET rate was correlated to the 

specified ET demand rate. The cumulative volumetric ET 

rate for the watershed was limited by water availability in 

the unsaturated zone and the depth to the water table, except 

during periods of high infiltration. Because the water table 

generally was deeper than the ET extinction depth for most of 

the watershed, a larger proportion of ET was removed from 

the unsaturated zone (fig. 9). However, the quantity of ET 

loss is sensitive to the extinction depth and extinction water 

content, such that attention should be placed on estimating 

the extinction depth and extinction water content to reduce 

uncertainty in recharge estimates. 

Most of the simulated streamflow during the first 

six months was from ground-water discharge to the stream 

(fig. 10). Total overland flow, including ground-water 

discharge to land surface plus saturation excess, exceeded 

ground-water discharge to streams beginning mid-May until 

mid-August. After mid-August, total overland flow was 

diminished and ground-water discharge to streams dominated 

the hydrograph (fig. 10). 

The outlet streamflow at the end of the simulation period 

was greater by 5,900 m3/d than that at the beginning; whereas 

overland flow caused by ground-water discharge to land 

surface (spring flow) was greater by 1,700 m3/d (fig. 10). The 

greater ground-water discharge at the end of the simulation 

is consistent with the greater mean recharge of 46,700 m3/d 

for the one-year simulation compared with a mean recharge 

of 10,400 m3/d for the steady-state period. The results of the 

simulation show how a year of well-above-mean recharge may 

have a prolonged effect on ground-water discharge to land 

surface and to streams. 

0

0.005

0.010

0.015

0.020

0.025

Dec Jan Feb Mar Apr May June July Aug Sept Oct Nov

Infiltration

Evapotranspiration demand

VOLUMETRIC RATE, IN CUBIC METER

 PER DAY PER SQUARE METER

MONTH

Figure 7. Infiltration rate and evapotranspiration demand used for test simulation 1.

The infiltration rate was multiplied by the spatial distribution of infiltration factors shown in figure 8. ET 

demand was constant over the model and was satisfied by storage in the unsaturated zone or from ground 

water.

16 Documentation of the Unsaturated-Zone Flow Package
Figure 8. Distribution of infiltration factors used for test simulation 1. 

The factors were multiplied by infiltration rates that varied daily as shown in figure 7. Infiltration factors are based on land-surface altitude 

and were used to account for variations of precipitation in relation to altitude.

Infiltration and ground-water recharge during the 

simulations resulted in a general increase of ground-water 

seepage to land surface between the beginning of the 

simulation and the end of June when ground-water recharge 

was at its peak (fig. 10). The increase in the area over which 

ground water discharged to land surface caused the increase in 

overland runoff (ground-water discharge to land surface plus 

saturation excess) shown in figure 11.

Water-content profiles beneath individual models cells 

indicate that the water content follows a seasonal envelope that 

becomes narrower with depth because of the seasonal pattern 

of infiltration and ET (fig. 12). The water-content profiles 

illustrate the utility of the kinematic-wave approximation in 

computing vertical flow through the unsaturated zone. The 

wetting front created by many different infiltration rates took 

longer than a year to reach the water table in areas where the 

saturated vertical hydraulic conductivity was low and depth to 

ground water was greater than 20 m below land surface.

The marked difference among infiltration, recharge, 

and streamflow for the example of the Sagehen watershed 

suggests that the unsaturated zone could have a strong effect 

on fluctuations in ground-water levels and streamflow during 

and after storms. Years in which infiltration rates exceed the 

mean annual rate may increase ground-water storage even 

in relatively poorly permeable rocks and result in increased 

discharge to streams that may extend well beyond the period 

of infiltration. The UZF1 Package in MODFLOW-2005 

provides a method in which daily-averaged values of overland 

runoff and saturation excess runoff can be simulated in a 

complex, yet simplified, manner. 

N

Infiltration 

rate factor 

3

2

1

0 1 2 4

KILOMETERS

Sagehen Creek watershed

Model boundary

Test Simulations 17
NV_UZF1_fig09

0

100

200

300

400

500

600

Dec Jan Feb Mar Apr May June July Aug Sept Oct Nov

Infiltration

Unsaturated-zone evapotranspiration

Ground-water evapotranspiration

Recharge

VOLUMETRIC RATE, IN THOUSANDS OF CUBIC METERS PER DAY 

MONTH

Figure 9. Volumetric rates of infiltration, recharge, and evapotranspiration summed over the model domain 

for test simulation 1.

0

20

30

40

50

60

Dec Jan Feb Mar Apr May June July Aug Sept Oct Nov

VOLUMETRIC RATE, IN THOUSANDS OF CUBIC METERS PER DAY

Streamflow at basin outlet

Ground-water discharge to land surface

plus saturation excess

Ground-water discharge to stream

10

MONTH

Figure 10. Streamflow at the basin outlet, ground-water discharge to land surface plus saturation 

excess, and ground-water discharge to streams for test simulation 1.

18 Documentation of the Unsaturated-Zone Flow Package
N

N

0 1 2 4 KILOMETER

Model boundary

Sagehen Creek watershed

Sagehen Creek watershed

Model cell with

 no discharge

Model cell with discharge—discharge

 increased from light blue to red

A. Steady state

B. July 1

Figure 11. Distribution of ground-water seepage to land surface for (A) steady state, and (B) July 1.

Test Simulations 19
NV_UZF1_fig12

0

5

10

15

20

0.08 0.10 0.12 0.14 0.16

February 1

June 1

August 1

October 1

December 1

DEPTH, IN METERS BELOW LAND SURFACE

WATER CONTENT, IN CUBIC METER PER CUBIC METER

Steady state 

Test Simulation 2: Simple Hypothetical Problem

A second smaller test problem was developed for use 

as a data-input template and testing proper code installation. 

Abbreviated listings of the data-input and output files for test 

simulation 2 are included in appendix 2. 

Test simulation 2 used the hypothetical stream-aquifer 

interaction problem described in Prudic and others (2004; 

p. 13–19), except the UZF1 Package was used to replace the 

ET and Recharge Packages. The units used in the hypothetical 

stream-aquifer interaction problem were in English units, and 

were not converted to SI units for test simulation 2. Thus, 

the data input and printed output listed in tables 4–8 and in 

appendix 2 are in English units. 

The hypothetical stream-aquifer interaction problem 

used in this test simulation was developed for an alluvial 

basin in a semiarid region. The principal aquifer comprises 

unconsolidated deposits of mostly sand and gravel. Recharge 

to the aquifer primarily is leakage from the streams that enter 

the basin from the mountains on the northwest, northeast, and 

southwest (fig. 13). The main stream in the southern part of 

the valley is perennial, whereas all other streams in the valley 

are intermittent with small drainage areas. 

Streamflow entering the model domain and diversions 

from streams were the same as that used by Prudic and others 

(2004) and remained the same for all stress periods (table 4). 

Different methods were used to represent the relations 

between width, depth and flow in the streams (see ICALC

in table 4 and ICALC option descriptions in Prudic and 

others, 2004). Infiltration rates in the UZF1 Package file and 

pumping rates in the Well Package file were varied during 12 

stress periods (table 5) but the distribution of infiltration and 

pumping among grid cells did not change. The ET demand 

was specified as constant over the entire modeled area equal 

to 1.6 ft/yr to a maximum depth of 15 ft below land surface. 

The LPF Package was used to specify hydraulic properties for 

the aquifer. The hydraulic conductivity and specific yield were 

173 ft/d and 0.2, respectively, in the vicinity of the stream 

channels and the hydraulic conductivity and specific yield 

were 35 ft/d and 0.1, respectively, elsewhere in the alluvial 

basin (fig. 14; table 6). 

Figure 12. Water-content profiles in the unsaturated zone of a model cell located in an upland area (column 9, row 

36) at selected dates for test simulation 1.

20 Documentation of the Unsaturated-Zone Flow Package
Figure 13. Well locations, elevation of top of consolidated rocks beneath basin-fill aquifer, and 

distribution of hydraulic conductivity and specific yield of the basin-fill aquifer used in test simulation 2.

(Modified from Prudic and others, 2004.)

NV_UZF1_fig13

0 10,000 20,000 FEET N

1000

900

800

700

600

700

800

900

Green Valley

COLUMNS

ROWS

1

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

2 3 4 5 6 7 8 9 10

900

EXPLANATION

Top of consolidated rock contour—Shows elevation of top of

consolidated rock. Contour interval is 100 feet

Basin-fill boundary

Well

Hydraulic conductivity = 35 feet per day; specific yield = 0.1

Hydraulic conductivity = 173 feet per day; specific yield = 0.2

Test Simulations 21
Figure 14. Hypothetical basin-fill aquifer with model grid, land-surface contours, active cells, and 

stream segment and reach numbering for test simulation 2.

(Modified from Prudic and others, 2004).

22 Documentation of the Unsaturated-Zone Flow Package

Ditch

1100

1100

1060

1080

1040

1040

1060

1080

1020

1000

0 10,000 20,000 FEET N

1,1

1,2

1,3

1,4

2,1

2,2

2,3

2,4

2,5

6,1

6,2

6,3

6,4

6,5

7,1

7,2 7,3 8,3 8,2 8,1 8,4 8,5

8,6

5,4

5,3

5,2

5,1

4,3

4,2 4,1

3,6

3,5

3,4

3,1 3,2 3,3

1,1

1,2

1,3

1,4

2,1

2,2

2,3

2,4

2,5

6,1

6,2

6,3

6,4

6,5

7,1

7,2 7,3 8,3 8,2 8,1 8,4 8,5

8,6

5,4

5,3

5,2

5,1

4,3

4,2 4,1

3,6

3,5

3,4

3,1 3,2 3,3

8,4 Segment number and reach number

EXPLANATION

Contour—Shows altitude of land surface.

 Contour interval is 20 feet 

Basin-fill boundary

1000

Active model cell

Green Valley

COLUMNS

ROWS

1

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

2 3 4 5 6 7 8 9 10

Blue River

Green Creek

Little

Creek
The model grid was divided into uniformly spaced cells 

5,000 ft on each side (fig. 14). The strongly implicit procedure 

(SIP) was used to solve the flow equation for test simulation 2. 

The head-closure criterion was 0.0002. The option to add 

overland runoff to stream segments was specified in the UZF1 

Package; however, no ground-water seepage to land surface 

was simulated in the test. The number of stress periods in the 

Discretization file was 12. Each stress period was 2.628 ×106

seconds or 30.42 days. The first stress period was steady 

state and had one time step. This allowed for the calculation 

of initial water contents for the unsaturated zone on the basis 

of the steady-state infiltration rate. The remaining 11 stress 

periods were transient and were divided into 15 time steps that 

increased sequentially by a factor of 1.1. 

Table 4. Method for calculating stream depth and width, inflow rates, streambed properties, stream dimensions, and Manning’s roughness coefficients 

for stream segments used in test simulation 2.

[Stream: Streambed properties and stream dimensions for each stream segment are specified for the first and last reach in each segment. A linear interpolation 

is used to determine each value at the midpoint of each reach. Stream segment no: Segment 3 is a diversion from segment 2 assuming all available flow to the 

rate specified is diverted. Method determines how stream depth and width are calculated for each reach in a segment. Stream depth and width are specified 

when ICALC is 0; stream depth is calculated from Manning’s equation assuming a wide rectangular channel when ICALC is 1; stream depth and width are 

calculated from Manning’s equation using an eight-point cross section when ICALC is 2; depth and width are calculated from a power function when ICALC

is 3; and depth and width are calculated from a table of values when ICALC is 4. See appendix 1 in Prudic and others (2004). Inflow rate: Inflow calculated 

from 1 or more tributary segments. Manning’s roughness coefficient is not used for calculating stream depth or width when ICALC is 0, 3, or 4, and Manning’s 

roughness coefficient for the overbank is not used when ICALC is 1. Abbreviations: ft, foot; ft/s, foot per second; ft3/s, cubic foot per second]

Stream

Method 

(ICALC)

Inflow rate 

(ft3/s)

Streambed 

hydraulic 

conductivity 

(ft/s)

Streambed 

thickness 

(ft)

Elevation of 

streambed 

(ft)

Stream 

depth (ft)

Stream 

width (ft)

Manning’s roughness 

coefficient (unitless)

Segment

No.

First and 

last reach Channel Over bank

1 1 4 25 0.00003 3 1,095 – – – –

4 .00003 3 1,075 – –

2 1 1 – .00003 3 1,075 – 12 0.030 –

5 .00003 3 1,050 – 12

3 1 0 10 .00003 2 1,075 2 10 – –

6 .00003 2 1,060 1 6

4 1 1 10 .00003 3 1,080 – 10 .030 –

3 .00003 3 1,060 – 10

5 1 3 – .00003 3 1,060 – – – –

4 .00003 3 1,045 – –

6 1 1 – .00003 3 1,045 – 12 .030 –

5 .00003 3 1,025 – 12

7 1 2 150 .00003 3 1,040 – – .025 0.045

3 .00003 3 1,025 – –

8 1 2 – .00003 3 1,025 – – .025 .045

6 .00003 3 990 – –

Water budgets at the end of selected stress periods for 

ground water and the unsaturated zone are summarized in 

tables 7 and 8, respectively. The results indicate that changes 

in flow and storage in the unsaturated zone affect groundwater recharge and, when combined with stream interaction, 

ground-water ET, and ground-water pumping, the effect 

of the unsaturated zone can be important to the overall 

assessment of ground-water flow. The absolute mass-balance 

errors were less than 0.05 percent for both the unsaturatedzone and ground-water budgets, and the cumulative massbalance error was 0 percent for the unsaturated-zone budget 

and -0.01 percent for the ground-water budget.

Test Simulations 23
24 Documentation of the Unsaturated-Zone Flow Package

Table 6. Hydraulic properties and other selected variables used in the Layer-Property Flow (LPF) and Unsaturated-Zone Flow (UZF1) Packages for test 

simulation 2.

[English units used for test simulation 2 are consistent with units of example problem as described by Prudic and others (2004) and were not changed to conform 

to the Standard International units used elsewhere in this document. Abbreviations: ft/s, foot per second; ft, foot; ft3, cubic foot]

Variable Value

Variables assigned to the Layer-Property Flow (LPF) Package

Horizontal hydraulic conductivity (ft/s) 4 ×10-4 to 2 ×10-3

Vertical hydraulic conductivity (ft/s) 1 ×10-6

Specific storage (per ft) 1.0 ×10-6

Specific yield (ft3 of water drained per ft3 of aquifer) 0.1 to 0.2

Variables assigned to the Unsaturated-Zone Flow (UZF1) Package

Saturated water content of unsaturated zone (ft3 of water per ft3 of bulk volume) 0.3

Brooks-Corey exponent (unitless) 3.5

Vertical hydraulic conductivity of the unsaturated zone (ft/s) 1 ×10-6

Steady-state infiltration rate (ft/s) 1 ×10-9 to 8 ×10-9

Evapotranspiration rate (ft/s) 5 ×10-8

Evapotranspiration extinction depth (ft) 15

Evapotranspiration extinction water content (ft3 of water per ft3 of bulk volume) 0.100 to 0.205

Table 5. Specified infiltration and pumping rates for the 12 

stress periods used in test simulation 2.

[English units used for test simulation 2 are consistent with units 

of example problem as described by Prudic and others (2004) 

and were not changed to conform to the Standard International 

units used elsewhere in this document. Infiltration rate was 

multiplied by the area of each cell (2.5 × 108 square feet) and 

by a factor that ranged from 1 to 8. Pumping rate was assigned 

to 10 model cells (appendix 2). Abbreviations: ft/s, foot per 

second; ft3/s, cubic foot per second]

Stress period Infiltration rate

(ft/s)

Pumping rate

(ft3/s)

1 1 ×10-9 -2

2 8 ×10-9 -2

3 1 ×10-8 -1

4 2 ×10-8 -1

5 5 ×10-9 -3

6 1 ×10-8 0

7 3 ×10-9 -2

8 5 ×10-9 0

9 3 ×10-9 -1

10 2 ×10-9 0

11 1 ×10-9 -3

12 1 ×10-9 0
Table 7. Computed ground-water budget for test simulation 2.

[English units used for test simulation 2 are consistent with units of example problem as described by Prudic and others 

(2004) and were not changed to conform to the Standard International units used elsewhere in this document. All units are in 

cubic foot per second, rounded to nearest tenth]

Steady 

state End of stress

Period 1 Period 3 Period 5 Period 7 Period 9 Period 12

Inflow

Leakage from streams 43.7 42.8 42.0 42.4 42.7 42.9

Recharge 8.5 16.0 8.3 9.3 11.4 11.3

Underflow beneath streams entering basin .2 .2 .2 .2 .2 .2

Total Inflow 52.4 59.0 50.5 51.9 54.3 54.4

Outflow

Wells 20.0 10.0 30.0 20.0 10.0 0.0

Evapotranspiration from ground water 19.5 13.6 17.4 17.0 17.2 18.6

Leakage to streams 12.6 14.1 14.2 13.9 13.9 13.6

Underflow beneath stream exiting basin .3 .3 .3 .3 .3 .3

Ground-water seepage to land surface .0 .0 .0 .0 .0 .0

 Total outflow 52.4 38.0 61.9 51.2 41.4 32.5

 Removed from storage .0 .2 14.9 5.5 .8 .8

 Added to storage .0 21.1 3.6 6.3 13.8 22.8

Table 8. Computed unsaturated-zone budget for test simulation 2.

[English units used for test simulation 2 are consistent with units of example problem as described by Prudic and others 

(2004) and were not changed to conform to the Standard International units used elsewhere in this document. All units are in 

cubic foot per second, rounded to nearest tenth]

Steady 

state End of stress

Period 1 Period 3 Period 5 Period 7 Period 9 Period 12

Inflow

Infiltration 8.5 84.8 42.4 25.4 25.4 8.5

Outflow

Unsaturated-zone evapotranspiration – 109.5 61.8 69.4 58.1 29.8

Ground-water recharge 8.5 16.0 8.3 9.3 11.4 11.3

 Storage change 0.0 -40.7 -27.7 -53.3 -44.1 -32.6

Test Simulations 25
Implementation of UZF1 With 

MODFLOW-2005

The UZF1 Package comprises 13 subroutines and one 

Fortran module that are used with the Ground-Water Flow 

(GWF) Process of MODFLOW-2005 (Harbaugh, 2005). 

Each subroutine or module incorporated in the main listing 

of MODFLOW-2005 is identified first by the Ground-Water 

Flow (GWF2) Process, followed by the UZF1 Package, 

and then followed by the division of procedures used in 

MODFLOW-2005. The Fortran module is used to share 

and store data. Memory allocation, reading input data and 

initialization of variables, formulating characteristic and finitedifference equations, and computing mass balances occur 

within GWF2UZF1AR, GWF2UZF1RP, GWF2UZF1FM, 

and GWF2UZF1BD subroutines, respectively. The subroutine 

GWF2UZF1VKS was added to define the saturated vertical 

hydraulic conductivity in the unsaturated zone, either as a 

variable read into the UZF1 Package in the GWF2UZF1RP 

subroutine or obtained from the LPF or Block-Centered 

Flow (BCF) Packages (Harbaugh, 2005). In addition to 

the subroutines called from the main program, several 

new subroutines were included in the UZF1 Package for 

routing kinematic waves through the unsaturated zone. 

These subroutines are UZF1FLOW2, LEADWAVE2, and 

TRAILWAVE2. The subroutine GWF2UZF1OLF was 

included to add ground-water discharge to land surface 

and overland flow to streams and lakes. The subroutine 

TRANSPIRATION was included to account for ET loss 

from the unsaturated zone. Finally, memory deallocation, 

setting pointer for all data for the package, and saving pointer 

locations occurs within GWF2UZF1DA, GWF2UZF1PNT, 

and GWF2UZF1PSV, respectively.

Additionally, the MAIN program of MODFLOW-2005 

has been modified to designate “UZF” in position 56 in 

CUNIT. Whenever “UZF” is used at the beginning of a line in 

the NAME file, the program will use the specified unit number 

and file name to obtain information for the UZF1 Package. 

Output from the UZF1 Package can be printed to any or all of 

(1) the main listing file; (2) a separate listing file; (3) separate 

files for each specified cell that provides a time series of 

selected information; and (4) separate files for budget items 

summed for all model cells.

The UZF1 Package has been programmed to apply 

infiltration to a designated underlying cell. If the water level 

in the underlying cell becomes less than the bottom elevation 

of the cell, the cell becomes inactive. When this happens, 

unsaturated flow will be routed vertically to the uppermost 

active cell. However, the unsaturated zone hydraulic properties 

remain unchanged even if the properties change across 

model layers. If all cells in the vertical column are inactive, 

no interaction with ground water is allowed and water is not 

routed through the unsaturated zone and a warning is printed 

in the main output file.

Summary

The U.S. Geological Survey modular finite-difference 

ground-water flow model (MODFLOW-2005) has been used 

to model surface-water and ground-water interaction over a 

wide range of spatial scales. The capability to simulate vertical 

unsaturated flow was added to MODFLOW-2005 for the 

purpose of simulating the delay between infiltration at land 

surface and recharge at the water table. Unsaturated flow is 

simulated using a kinematic-wave approximation of Richards’ 

equation and assuming one-dimensional gravity-driven 

vertical flow. The method of characteristics is used to solve the 

kinematic-wave approximation. The approach simulates flow 

and storage in the unsaturated zone and time-delayed recharge, 

while maintaining the applicability of MODFLOW-2005 to 

basin-scale problems.

Unsaturated flow is simulated independently of saturated 

flow within each model cell in which the water table is below 

land surface. The relation between the unsaturated hydraulic 

conductivity and water content in the unsaturated zone is 

defined on the basis of the Brooks-Corey function. Variables 

used to simulate unsaturated flow include saturated and initial 

water content, saturated vertical hydraulic conductivity, and 

the Brooks-Corey exponent. Additionally, the residual water 

content is approximated by the aquifer specific retention in 

the Brooks-Corey function. Specific retention is calculated 

as the porosity represented by the saturated water content 

minus specific yield. These variables can be different for 

each MODFLOW column of cells but cannot change between 

model layers. 

The maximum infiltration rate (flow per unit area) 

across the land surface is limited to the vertical saturated 

hydraulic conductivity. Recharge is added to ground water 

only when water in the unsaturated zone drains to the water 

table. This differs from the standard Recharge Package, in 

which the recharge rate must be determined separately from 

MODFLOW. When the water table rises, the quantity of water 

stored in the unsaturated zone above the residual water content 

is added to recharge. Thus, the water table can rise rapidly 

when the amount of water stored in the unsaturated zone 

is at or near saturation. A water budget for the unsaturated 

zone is tracked independently of the ground-water budget in 

MODFLOW. 

Evapotranspiration may be simulated over a specified 

depth in the unsaturated zone. The volumetric rate of 

evapotranspiration is dependent on the volume of water stored 

in the unsaturated zone above the specified extinction depth 

and extinction water content. Evapotranspiration also may be 

supplied by ground water when the water table is above the 

specified extinction depth. Variables that must be specified 

to simulate evapotranspiration include the evapotranspiration 

demand rate, the extinction depth, and the extinction water 

content. 

26 Documentation of the Unsaturated-Zone Flow Package
Ground water may discharge to the land surface when 

using the UZF1 Package. This water may be added to a 

specified stream segment in the Streamflow-Routing (SFR2) 

Package or to a specified lake in the Lake (LAK3) Package, 

provided the packages are being used in the simulation. Water 

discharging to land surface is added to stream segments or 

lakes without delay during each time step, regardless of the 

distance between the cell where ground water is discharged to 

the nearest stream segment or lake.

Results from the kinematic-wave approximation in 

the UZF1 Package were compared with results from U.S. 

Geological Survey’s Variably-Saturated Two-Dimensional 

Flow and Transport (VS2DT) Model that solves Richards’ 

equation. The test simulation demonstrates the capability of 

the UZF1 Package to simulate the advancement of a wetting 

front assuming three different rates of evapotranspiration. 

The simulation was in good agreement in the magnitude and 

downward progression of the wetting front with respect to the 

results from the VS2DT. An application to a 27 km2 watershed 

was used to demonstrate the capabilities of the UZF1 Package. 

A 1-year transient simulation in which infiltration and 

evapotranspiration demand changed daily demonstrated that 

simulating flow through the unsaturated zone is important for 

estimating recharge. Years in which infiltration rates exceed 

the mean annual ET rate may increase unsaturated-zone and 

ground-water storage, even in relatively poorly permeable 

rocks, and can result in increased discharge to streams that 

may extend well beyond the period of infiltration. A second 

and simpler problem is included to provide an example of the 

input-file structure and model output and to provide testing of 

the package after installation. 

References Cited

Abbott, M.B., 1966, An introduction to the method of 

characteristics: New York, Elsevier Publishing Company, 

Inc., 243 p. 

Brooks, R.H., and Corey, A.T., 1966, Properties of porous 

media affecting fluid flow: American Society of Civil 

Engineers, Journal of Irrigation and Drainage, v. 101, 

p. 85-92.

Charbeneau, R.J., 1984, Kinematic models for soil moisture 

and solute transport: Water Resources Research, v. 20, 

no. 6, p. 699-706.

Chen, ZhiQiang, Govindaraju, R.S., and Kavvas, M.L., 1994, 

Spatial averaging of unsaturated flow equations under 

infiltration conditions over areally heterogeneous fields—1. 

Development of models: Water Resources Research, v. 30, 

no. 2, p. 523-534.

Childs, E.C., and Bybordi, M., 1969, The vertical movement 

of water in stratified porous material, 1. Infiltration: Water 

Resources Research, v. 5, no. 2, p. 446-459.

Colbeck, S.C., 1972, A theory of water percolation in snow: 

Journal of Glaciology, v. 2, no. 63, p. 369-385.

Fenske, J.P., and Prudic, D.E., 1998, Development of HMS/

MODFLOW for simulation of surface and groundwater 

flow, in Poeter, Eileen, Zheng, Chunmiao, and Hill, Mary, 

eds., MODFLOW ’98: Golden, Colo., Colorado School of 

Mines, October, 4-8, 1998, Proceedings, v. 1, p. 463-470.

Freeze, R.A., 1971, Three-dimensional, transient, saturatedunsaturated flow in a groundwater basin: Water Resources 

Research, v. 7, no. 5, p. 153-171. 

Harbaugh, A.W., 2005, MODFLOW-2005, the U.S. 

Geological Survey modular ground-water model—the 

Ground-Water Flow Process: U.S. Geological Survey 

Techniques and Methods 6-A16, variously paginated.

Harbaugh, A.W., Banta, E.R., Hill, M.C., and McDonald, 

M.G., 2000, MODFLOW-2000, the U.S. Geological Survey 

modular ground-water model-modularization concepts and 

the ground-water flow process: U.S. Geological Survey 

Open-File Report 00-92, 121 p.

Healy, R.W., 1990, Simulation of solute transport in variably 

saturated porous media with supplemental information on 

modification of the U.S. Geological Survey’s computer 

program VS2D: U.S. Geological Survey Water-Resources 

Investigations Report 90-4025, 125 p.

McDonald, M.G., and Harbaugh, A.W., 1988, A modular 

three-dimensional finite-difference ground-water flow 

model: U.S. Geological Survey Techniques of WaterResources Investigations, Book 6, Chapter A1, 586 p.

Mantoglou, Aristotelis, 1992, A theoretical approach for 

modeling unsaturated flow in spatially variable soils—

Effective flow models in finite domains and nonstationarity: 

Water Resources Research, v. 28, no. 1, p. 251-267.

Mast, M.A., and Clow, D.W., 2000, Environmental 

characteristics and water quality of hydrologic benchmark 

network stations in the Western United States, 1963-95: 

U.S. Geological Survey Circular 1173-D, 114 p.

Mein, R.G., and Larson, C.L., 1973, Modeling infiltration 

during a steady rain: Water Resources Research, v. 9, no. 2, 

p. 384-394.

Meinzer, O.E., 1923, Outline of ground-water hydrology, with 

definitions: U.S. Geological Survey Water-Supply Paper 

494, 71 p.

Merritt, M.L., and Konikow, L.F., 2000, Documentation of a 

computer program to simulate lake-aquifer interaction using 

the MODFLOW ground-water flow model and the MOC3D 

solute-transport model: U.S. Geological Survey WaterResources Investigations Report 00-4167, 146 p.

Niswonger, R.G., and Prudic, D.E., 2004, Modeling variably 

saturated flow using kinematic waves in MODFLOW, 

in Hogan, J.F., Phillips, F.M., and Scanlon, B.R., eds., 

Groundwater Recharge in a Desert Environment—The 

Southwestern United States: Washington, D.C., American 

Geophysical Union, Water Science and Application Series, 

v. 9, p. 101-112. 

References Cited 27
Niswonger, R.G., and Prudic, D.E., 2006, Documentation 

of the Streamflow-Routing (SFR2) Package to include 

unsaturated flow beneath streams—A modification to 

SFR1: U.S. Geological Survey Techniques and Methods, 

Book 6, Chap. A13, 57 p.

Pikul, M.F., Street, R.L., and Remson, Irving, 1974, A 

numerical model based on coupled one-dimensional 

Richards and Boussinesq equations: Water Resources 

Research, v. 10, no. 2, p. 295-302.

Prudic, D.E., Konikow, L.F., and Banta, E.R., 2004, A new 

Streamflow-Routing (SFR1) Package to simulate streamaquifer interaction with MODFLOW-2000: U.S. Geological 

Survey Open-File Report 2004-1042, 95 p.

Refsgaard, J.C., and Storm, B., 1995, MIKE SHE, in

Computer Models of Watershed Hydrology: Singh, V.P., 

ed., Highlands Ranch, Colo., Water Resources Publications, 

p. 809-846.

Smith, R.E., 1983, Approximate sediment water movement by 

kinematic characteristics: Soil Science Society of America 

Journal, v. 47, p. 3-8.

Smith, R.E., and Hebbert, R.H.B., 1983, Mathematical 

simulation of interdependent surface and subsurface 

hydrologic processes: Water Resources Research, v. 19, 

no. 4, p. 987-1001.

Smith, R.E., and Parlange, J.Y., 1978, A parameter-efficient 

hydrologic infiltration model: Water Resources Research, 

v. 14, no. 3, p. 533-538.

Sophocleous, Marious, and Perkins, S.P., 2000, Methodology 

and application of combined watershed and groundwater 

model in Kansas: Journal of Hydrology, v. 236, p. 185-201.

Thoms, R.B., Johnson, R.L., and Healy, R.W., 2006, User’s 

Guide to the Variably Saturated Flow (VSF) Process for 

MODFLOW: U.S. Geological Survey Techniques and 

Methods Book 6, Chap. A18, 58 p.

28 Documentation of the Unsaturated-Zone Flow Package
Appendix 1: Data Input Instructions for the Unsaturated-Zone Flow (UZF1) 

Package 

MODFLOW Name File

Percolation of water through an unsaturated zone is activated by including a record in the MODFLOW name file using the 

file type (Ftype) “UZF” to indicate that relevant calculations are to be made in the model and to specify the related input data 

file. The UZF1 Package is compatible with MODFLOW-2005 (Harbaugh, 2005), but not with earlier versions of MODFLOW. 

The user can optionally specify unsaturated-zone water budgets and water content profiles for selected model cells by including 

a record in the MODFLOW name file using the file type (Ftype) “DATA” that specifies the relevant output data file name for 

each model cell. The row and column cell indices for each selected model cell are included in the UZF output file. Three types 

of information may be printed to the specified file for each model cell depending on the OUTTYPE option specified. The three 

OUTTYPE options are: option 1 prints volumes of water entering, leaving, and stored in the unsaturated zone; option 2 prints 

volumes and rates for water entering, leaving, and stored within the unsaturated zone; and option 3 prints the water content 

profile between land surface and the water table. Additionally, a time series of infiltration, unsaturated zone evapotranspiration, 

recharge, and ground-water discharge summed over the model domain may be printed to a specified file. 

Unsaturated-Zone Flow Package Input Data

The UZF1 Package input file consists of items numbered from 0 through 16, each consisting of one or more records. 

These data are used to specify information about the hydraulic properties of the unsaturated zone, the infiltration rate, 

evapotranspiration, and certain output control options. All input variables are read using free format unless indicated otherwise. 

In free format, variables are separated by one or more spaces or by a comma and optionally one or more spaces. A blank field 

cannot be used to set a variable value to zero. Units of length (L) and time (T) are used to define the dimensions of each variable.

The UZF1 Package relies on the specific yield values as specified in the Layer-Property Flow (LPF) Package or the BlockCentered Flow (BCF) Package. Thus, the option for unsaturated flow is only available when either the LPF or BCF Package is 

used in the simulation. When the option to use vertical hydraulic conductivity in the LPF Package is specified, the layer(s) that 

contain cells where unsaturated flow will be simulated must be specified as convertible. That is, the variable LAYTYP specified 

in the LPF Package must not be equal to zero, otherwise the model will print an error and stop execution. The variable VKA

in LPF may be specified as the ratio of horizontal to vertical hydraulic conductivity (LAYVKA is not zero) or as the vertical 

hydraulic conductivity (LAYVKA=0). When using the BCF Package, the right digit of LTYPE (LAYCON) must be greater 

or equal to one or the model will print an error and stop execution. Neither parameter estimation nor solute transport can be 

simulated with the UZF1 Package. 

For Each Simulation

0. Data: [#Text]

Text A character variable (up to 199 characters) that is printed when the file is read. The “#” character must be in 

 column 1, and, accordingly, the variable starts in column 2. Any characters can be included in Text.

Note 1: Item 0 can be repeated multiple times.

1. Data: NUZTOP IUZFOPT IRUNFLG IETFLG IUZFCB1 IUZFCB2 NTRAIL2 NSETS2 NUZGAG SURFDEP 

NUZTOP An integer value used to define which cell in a vertical column that recharge and discharge is simulated.

 1 Recharge to and discharge from only the top model layer. This option assumes land surface is defined as 

 top of layer 1.

 2 Recharge to and discharge from the specified layer in variable IUZFBND. This option assumes land surface 

 is defined as top of layer specified in IUZFBND.

 3 Recharge to and discharge from the highest active cell in each vertical column. Land surface is determined 

 as top of layer specified in IUZFBND. A constant head node intercepts any recharge and prevents deeper 

 percolation.

Appendix 1 29
IUZFOPT An integer value equal to 0, 1, or 2. A value of zero indicates that flow will not be routed through the 

 unsaturated zone, rather infiltration will be applied directly to the saturated zone. A value of 1 indicates that 

 the vertical hydraulic conductivity will be specified within the UZF1 Package input file using array VKS. 

 A value of 2 indicates that the vertical hydraulic conductivity will be specified within either the BCF or LPF

 Package input file.

IRUNFLG An integer value that specifies whether ground water that discharges to land surface will be routed to stream

 segments or lakes as specified in the IRUNBND array (IRUNFLG not equal to zero) or if ground-water

 discharge is removed from the model simulation and accounted for in the ground-water budget as a loss of

 water (IRUNFLG=0). The Streamflow-Routing (SFR2) and(or) the Lake (LAK3) Packages must be active if

 IRUNFLG is not zero.

IETFLG An integer value that specifies whether or not evapotranspiration (ET) will be simulated. ET will not be

 simulated if IETFLG is zero, otherwise it will be simulated.

IUZFCB1 An integer value used as a flag for writing ground-water recharge, ET, and ground-water discharge to land

 surface rates to a separate unformatted file using subroutine UBUDSV. If IUZFCB1>0, it is the unit number

 to which the cell-by-cell rates will be written when “SAVE BUDGET” or a non-zero value for ICBCFL is

 specified in Output Control. If IUZFCB1 ≤ 0, cell-by-cell rates will not be written to a file. 

IUZFCB2 An integer value used as a flag for writing ground-water recharge, ET, and ground-water discharge to land

 surface rates to a separate unformatted file using module UBDSV3. If IUZFCB2>0, it is the unit number to

 which cell-by-cell rates will be written when “SAVE BUDGET” or a non-zero value for ICBCFL is specified

 in Output Control. If IUZFCB2≤0, cell-by-cell rates will not be written to file. 

NTRAIL2 An integer value equal to the number of trailing waves used to define the water-content profile following a

 decrease in the infiltration rate. The number of trailing waves varies depending on the problem, but a range

 between 10 and 20 is usually adequate. More trailing waves may decrease mass-balance error and will increase

 computational requirements and memory usage. This variable is not specified if IUZFOPT is zero or less.

NSETS2 An integer value equal to the number of wave sets used to simulate multiple infiltration periods. The number

 of wave sets should be set to 20 for most problems involving time varying infiltration. The total number of

 waves allowed within an unsaturated zone cell is equal to NTRAIL2*NSETS2. An error will occur if the

 number of waves in a cell exceeds this value. This variable is not specified if IUZFOPT is zero or less.

NUZGAG An integer value equal to the number of cells (one per vertical column) that will be specified for printing

 detailed information on the unsaturated zone water budget and water content. A gage also may be used to print

 the budget summed over all model cells.

SURFDEP A real value equal to the average undulation depth within a finite-difference cell.

Note 2: The variables NTRAIL2 and NSETS2 are used for dimensioning arrays.

2. Data: [IUZFBND (NCOL, NROW)] -- U2DINT

IUZFBND An array of integer values used to define the aerial extent of the active model in which recharge and discharge

 will be simulated.

If IRUNFLG > 0: Read item 3

3. Data: [IRUNBND (NCOL, NROW)] -- U2DINT

IRUNBND An array of integer values used to define the stream segments within the Streamflow-Routing (SFR2) Package

 or lake numbers in the Lake (LAK3) Package that overland runoff from excess infiltration and ground-water

 discharge to land surface will be added. A positive integer value identifies the stream segment and a negative

 integer value identifies the lake number.

30 Documentation of the Unsaturated-Zone Flow Package
IF IUZOPT = 1: Read item 4.

4. Data: [VKS (NCOL, NROW)] -- U2DREL

VKS An array of positive real values used to define the saturated vertical hydraulic conductivity of the unsaturated

 zone (LT-1).

Note 3: Vertical hydraulic conductivity is optional when either the LPF or BCF Package is used to define vertical hydraulic

conductivity of cells (IUZFOPT=2).

5. Data: [EPS (NCOL, NROW)] -- U2DREL

EPS An array of positive real values for each model cell used to define the Brooks-Corey epsilon of the unsaturated

 zone. Epsilon is used in the relation of water content to hydraulic conductivity (Brooks and Corey, 1966). This 

 variable is not specified if IUZFOPT is zero or less.

6. Data: [THTS (NCOL, NROW)] -- U2DREL

THTS An array of positive real values used to define the saturated water content of the unsaturated zone in units of

 volume of water to total volume (L3

L-3). This variable is not specified if IUZFOPT is zero or less.

7. Data: [THTI (NCOL, NROW)] -- U2DREL

THTI An array of positive real values used to define the initial water content for each vertical column of cells in 

 units of volume of water at start of simulation to total volume (L3L-3). THTI should not be specified for

 steady-state simulations. This variable is not specified if IUZFOPT is zero or less.

If NUZGAG>0: Item 8 is repeated NUZGAG times

8. Data: {IUZROW} {IUZCOL} IFTUNIT {IUZOPT}

IUZROW An integer value equal to the row number of the cell that unsaturated-zone information will be printed for each

 time step. The variable is not used when IFTUNIT is negative.

IUZCOL An integer value equal to the column number of the cell that unsaturated-zone information will be printed for

 each time step. The variable is not used when IFTUNIT is negative.

IFTUNIT An integer value equal to the unit number of the output file. A positive value is for output of individual cells

 whereas a negative value is for output that is summed over all model cells. The summed output includes

 applied infiltration, runoff, actual infiltration, ground-water discharge to land surface, ET from the unsaturated

 zone, ET from ground water, recharge, and change in unsaturated-zone storage.

IUZOPT An integer value that is a flag for the type of expanded listing desired in the output file. The variable is not

 used when IFTUNIT is negative.

 1 Prints time, ground-water head, and thickness of unsaturated zone, and cumulative volumes of

 infiltration, recharge, storage, change in storage and ground-water discharge to land surface.

 2 Same as option 1 except rates of infiltration, recharge, change in storage, and ground-water discharge

 also are printed.

 3 Prints time, ground-water head, thickness of unsaturated zone, followed by a series of depths and

 water contents in the unsaturated zone. This option cannot be used when IUZFOPT is zero or less.

Note 4: Item 8 must include exactly NUZGAG lines (records) of data. If NUZGAG > 1, it is permissible to repeat the same cell if 

a summary of volumes and rates and water contents are desired for a given cell. The records within item 8 can be listed in any 

order.

Note 5: A unique unit number for IFTUNIT must be specified for each line or record in item 8 and each number must 

match that used in the MODFLOW name file as Ftype “DATA” file types and file names (see Harbaugh, 2005).

Appendix 1 31
For Each Stress Period

9. Data: NUZF1 

NUZF1 An integer value for reusing or reading infiltration rates that can change each stress period. If NUZF1 > 0, 

 then infiltration rates at land surface are specified. If NUZF1 < 0, then infiltration rates from the previous 

 stress period are used.

If NUZF1 > 0:

10. Data: [FINF (NCOL, NROW)] – U2DREL

FINF An array of positive real values used to define the infiltration rates (LT-1) at land surface for each vertical

 column of cells.

Note 6: If FINF is specified as being greater than the vertical hydraulic conductivity then FINF is set equal to the vertical 

unsaturated hydraulic conductivity. Excess water is routed to streams or lakes when IRUNFLG is not zero, and if SFR2 or LAK3 

is active.

If IETFLG > 0: Read Items 11 through 16

11. Data: NUZF2 

NUZF2 An integer value for reusing or reading ET demand rates that can change each stress period. If NUZF2>0, then

 ET demand rates are specified. If NUZF2<0, then ET demand rates from the previous stress period are used.

If NUZF2 > 0:

12. Data: [PET (NCOL, NROW)] – U2DREL

PET An array of positive real values used to define the ET demand rates (L1T-1) within the ET extinction depth

 interval for each vertical column of cells.

13. Data: NUZF3

NUZF3 An integer value for reusing or reading ET extinction depths that can change each stress period. If NUZF3>0,

 then ET extinction depths are specified. If NUZF3<0, then depths from the previous stress period are used.

If NUZF3 > 0:

14. Data: [EXTDP (NCOL, NROW)] – U2DREL

EXTDP An array of positive real values used to define the ET extinction depths.

Note 7: The quantity of ET removed from a cell is limited by the volume of water stored in the unsaturated zone above 

the extinction depth. If ground water is within the ET extinction depth, then the rate removed is based on a linear decrease 

in the maximum rate at land surface and zero at the ET extinction depth. The linear decrease is the same method used in the 

Evapotranspiration Package (McDonald and Harbaugh, 1988, chap. 10). 

15. Data: NUZF4 

NUZF4 An integer value for reusing or reading the extinction water content that can change each stress period. If

 NUZF4>0, then extinction water contents are specified. If NUZF4<0, then the extinction water contents from

 the previous stress period are used. 

If NUZF4 > 0:

16. Data: [EXTWC (NCOL, NROW)] – U2DREL

EXTWC An array of positive real values used to define the extinction water content below which ET cannot be removed

 from the unsaturated zone. This variable is not specified if IUZFOPT is zero or less.

Note 8: EXTWC must have a value between (THTS-Sy) and THTS, where Sy is the specific yield specified in either the 

LPF or BCF Package.

32 Documentation of the Unsaturated-Zone Flow Package
Appendix 2: Selected Input Datasets and Printed Results for Test Simulation 2

Selected sections of several key input and output data files are shown below for test simulation 2—sometimes with 

annotations; gaps in the listings are indicated by an ellipsis. A complete set of the files is available for distribution over the 

Internet as discussed in the Preface. Contents of some files are enclosed in a border and explanations are noted outside of the 

border; for other files, explanations are sometimes included as comments following a semicolon on the line being explained. 

Some brief annotations were added to the selected output listings to help the reader understand the purpose of various sections of 

output. These annotations are written in bold italics. Font sizes in the following listings are sometimes reduced so that lines will 

fit within page margins. Information pertaining specifically to the UZF1 Package is highlighted by shading.

Listing of Selected Input Data Sets for Test Simulation 2

Following (enclosed in a border) are the contents of the MODFLOW name file for test simulation 2. Shading indicates 

input files in which values are used in the UZF1 Package and output files in which information from the UZF1 Package are 

printed or saved:

File name: UZFtest2.nam

LIST 9 Out\UZFtest2.lst ← Output file for MODFLOW

BAS6 75 Data\UZFtest2.ba6 ← Input file for Basic Package

LPF 7 Data\UZFtest2.lpf ← Input file for Layer-Property Flow Package

DIS 8 Data\UZFtest2.dis ← Input file for Discretization file

SIP 13 Data\UZFtest2.sip ← Input file for Strongly Implicit Procedure Package

OC 14 Data\UZFtest2.oc ← Input file for Output Control option 

SFR 15 Data\UZFtest2.sfr ← Input file for Streamflow-Routing Package

WEL 16 Data\UZFtest2.wel ← Input file for Well Package

GHB 17 Data\UZFtest2.ghb ← Input file for General-Head Boundary Package

UZF 19 Data\UZFtest2.uzf ← Input file for Unsaturated-Zone Flow Package

GAGE 32 Data\UZFtest2.gag ← Input file for Gage Package

DATA 58 Out\UZFtest2hd.out ← Output file for calculated heads

DATA(BINARY) 61 Out\UZFtest2.uzfot ← Output file for recharge and ground-water discharge

DATA 65 Out\UZFtest2.uzf1 ← Output file unsaturated zone cell 1 (row 3, column 6) 

DATA 66 Out\UZFtest2.uzf2 ← Output file unsaturated zone cell 2 (row 6, column 3)

DATA 67 Out\UZFtest2.uzf3 ← Output file unsaturated zone cell 3 (row 10, column 5)

DATA 68 Out\UZFtest2.uzf4 ← Output file of times series of unsaturated-zone water budgets

DATA 81 Out\UZFtest2.flw ← Auxiliary output file for the Streamflow-Routing Package 

DATA 83 Out\UZFtest2.sg1 ← Output file for last reach in stream segment 1

DATA 84 Out\UZFtest2.sg2 ← Output file for last reach in stream segment 2

DATA 85 Out\UZFtest2.sg3 ← Output file for last reach in stream segment 3

DATA 86 Out\UZFtest2.sg4 ← Output file for last reach in stream segment 4

DATA 87 Out\UZFtest2.sg5 ← Output file for last reach in stream segment 5

DATA 88 Out\UZFtest2.sg6 ← Output file for last reach in stream segment 6

DATA 89 Out\UZFtest2.sg7 ← Output file for last reach in stream segment 7

DATA 90 Out\UZFtest2.sg8 ← Output file for last reach in stream segment 8

DATA 91 Out\UZFtest2dv.sg9 ← Output file for first reach in stream segment 3

 ↑ ↑ ↑

 1 2 3

1 Ftype (the type of file).

2 Unit number.

3 File name (name chosen to reflect contents of file).

Appendix 2 33
Following (enclosed in a border) are the contents of the Basic Package input file for test simulation 2; explanations 

are noted outside of border:

File name: UZFTest2.ba6

#Test Simulation 2-- Example problem using the UZF1 Package ← 1

#Modified from example Problem 1 of SFR1 documentation report ← 1

FREE ← 2

INTERNAL 1 (10I4) 6 IBOUND layer 1 ← 3

 1 1 0 0 0 0 0 0 0 0 ← 3

 1 1 1 0 0 1 1 1 0 0

 0 1 1 1 1 1 1 1 0 0

 1 1 1 1 1 1 1 1 1 0

 1 1 1 1 1 1 1 1 1 1

 1 1 1 1 1 1 1 1 1 1

 1 1 1 1 1 1 1 1 1 0

 0 1 1 1 1 1 1 1 1 0

 0 1 1 1 1 1 1 1 1 0

 0 1 1 1 1 1 1 1 1 1

 0 1 1 1 1 1 1 1 1 0

 1 1 1 1 1 1 1 1 1 0

 1 1 1 1 1 1 1 1 1 0

 0 1 1 1 1 1 1 1 1 0

 0 0 1 1 1 1 0 0 0 0

6999.00 HNOFLO ← 4

INTERNAL 1.000E+00 (10F8.0) 0 SHEAD layer 1 ← 5

 1092.9 1089.4 6999.0 6999.0 6999.0 6999.0 6999.0 6999.0 6999.0 6999.0 ← 5

 1089.9 1086.7 1081.4 6999.0 6999.0 1072.3 1071.5 1071.1 6999.0 6999.0

 6999.0 1079.4 1077.5 1073.6 1072.5 1070.6 1069.1 1068.3 6999.0 6999.0

 1069.1 1070.7 1070.2 1067.0 1065.5 1063.6 1063.9 1065.9 1066.1 6999.0

 1064.9 1063.9 1062.9 1061.5 1060.8 1060.8 1060.2 1062.0 1065.7 1068.1

 1059.8 1058.3 1057.2 1056.4 1056.3 1056.9 1057.2 1059.1 1061.9 1065.6

 1056.1 1053.0 1051.7 1051.2 1051.3 1051.9 1053.2 1055.5 1056.5 6999.0

 6999.0 1046.3 1045.7 1045.5 1045.8 1046.5 1047.8 1049.9 1050.4 6999.0

 6999.0 1039.1 1039.1 1039.5 1040.0 1041.0 1041.6 1042.9 1044.0 6999.0

 6999.0 1031.1 1031.8 1033.3 1034.3 1035.4 1036.5 1037.6 1039.3 1041.6

 6999.0 1021.4 1023.1 1026.8 1028.7 1030.5 1032.1 1033.6 1035.7 6999.0

 1001.2 1009.2 1014.5 1019.6 1023.1 1026.2 1028.6 1030.4 1034.3 6999.0

 994.5 1000.0 1006.7 1012.9 1018.7 1024.3 1029.1 1034.1 1036.1 6999.0

 6999.0 1004.0 1009.0 1013.7 1018.9 1024.0 1029.5 1035.1 1039.3 6999.0

 6999.0 6999.0 1011.9 1015.1 1019.1 1022.4 6999.0 6999.0 6999.0 6999.0

1 Two header lines of comments.

2 Options line.

3 IBOUND array.

4 Head value assigned to inactive cells.

5 Initial head; rounded to the nearest tenth of foot.

34 Documentation of the Unsaturated-Zone Flow Package
Following (enclosed in a border) are the contents of the Discretization input file for test simulation 2; explanations are 

noted outside of border:

File name: UZFtest2.dis

#Test Simulation 2-- Example problem using the UZF1 Package

#Modified from example test problem 1 of SFR1 documentation report

1 15 10 12 1 1 NLAY,NROW,NCOL,NPER,ITMUNI,LENUNI ← 1

 0 ← 2

CONSTANT 5.000E+03 DELR ← 3

CONSTANT 5.000E+03 DELC ← 4

INTERNAL 1.000E+00 (10F7.0) 12 TOP ← 5

1105.0 1105.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ← 5

1105.0 1095.0 1105.0 0.0 0.0 1100.0 1100.0 1110.0 0.0 0.0

 0.0 1100.0 1085.0 1080.0 1078.0 1076.0 1073.0 1090.0 0.0 0.0

1110.0 1100.0 1090.0 1085.0 1076.0 1080.0 1080.0 1070.0 1100.0 0.0

1110.0 1100.0 1085.0 1080.0 1078.0 1074.0 1077.0 1066.0 1075.0 1080.0

1105.0 1095.0 1080.0 1077.0 1073.0 1069.0 1070.0 1063.0 1085.0 1095.0

1105.0 1090.0 1075.0 1072.0 1070.0 1065.0 1065.0 1069.0 1090.0 0000.0

 0.0 1080.0 1070.0 1068.0 1063.0 1060.0 1063.0 1067.0 1090.0 0000.0

 0.0 1075.0 1065.0 1062.0 1058.0 1052.0 1055.0 1065.0 1090.0 0000.0

 0.0 1070.0 1060.0 1055.0 1046.0 1046.0 1044.0 1060.0 1088.0 1100.0

 0.0 1055.0 1050.0 1045.0 1041.0 1042.0 1040.0 1050.0 1075.0 0.0

1027.0 1028.0 1028.0 1032.0 1032.0 1033.0 1032.0 1040.0 1060.0 0.0

997.0 1004.0 1008.0 1014.0 1020.0 1027.0 1032.0 1038.0 1050.0 0.0

 0.0 1030.0 1040.0 1045.0 1050.0 1050.0 1050.0 1050.0 1044.0 0.0

 0.0 0.0 1080.0 1080.0 1080.0 1080.0 0.0 0.0 0.0 0.0

INTERNAL 1.000E+00 (10F7.0) 0 BOT layer 1 ← 6

1000.0 1000.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ← 6

1000.0 900.0 900.0 0.0 0.0 950.0 950.0 950.0 0.0 0.0

 0.0 850.0 670.0 680.0 680.0 770.0 770.0 890.0 0.0 0.0

 925.0 850.0 650.0 650.0 650.0 680.0 740.0 850.0 950.0 0.0

 925.0 800.0 650.0 600.0 600.0 600.0 650.0 720.0 850.0 950.0

 925.0 800.0 650.0 600.0 550.0 600.0 650.0 720.0 850.0 970.0

 925.0 800.0 650.0 600.0 550.0 600.0 650.0 720.0 950.0 0.0

 0.0 850.0 700.0 650.0 570.0 600.0 700.0 850.0 950.0 0.0

 0.0 850.0 700.0 650.0 540.0 580.0 700.0 850.0 950.0 0.0

 0.0 850.0 700.0 630.0 530.0 580.0 650.0 800.0 900.0 950.0

 0.0 850.0 700.0 620.0 580.0 620.0 700.0 800.0 900.0 0.0

 900.0 830.0 710.0 640.0 600.0 650.0 700.0 800.0 900.0 0.0

 900.0 810.0 750.0 700.0 650.0 700.0 750.0 800.0 900.0 0.0

 0.0 950.0 850.0 850.0 850.0 850.0 850.0 850.0 850.0 0.0

 0.0 0.0 950.0 950.0 950.0 950.0 0.0 000.0 000.0 0.0

 2.628E+06 1 1.00 SS STRESS PERIOD 1 ← 7

 2.628E+06 15 1.10 TR STRESS PERIOD 2 ← 7

... ← 7

 2.628E+06 15 1.10 TR STRESS PERIOD 12 ← 7

1 Grid characterization data. 

2 Flags for quasi-3D simulation.

3 Row discretization data.

4 Column discretization data.

5 Top elevations for layer 1 read from following data.

6 Bottom elevations for layer 1 read from following data.

7 Time step and stress period information.

Appendix 2 35
Following (enclosed in a border) are the contents of the Layer-Property Flow Package input file for test simulation 2; 

explanations are noted outside of border. Shading indicates that the specific yield specified within the Layer-Property Flow 

Package is used to calculate the residual water content in the UZF1 Package.

File name: UZFtest2.lpf

#Layer-Property Flow Package input for test simulation 2 ← 1

0 0. 0 ILPFCB, HDRY, NPLPF ← 2

1 LAYTYP ← 3

0 LAYAVG ← 3

1.0 CHANI ← 3

0 LAYVKA ← 4

0 LAYWET ← 5

INTERNAL 2.000E-04 (10F5.0) 0 HY layer 1 ← 6

10.0 10.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ← 6

 2.0 10.0 10.0 0.0 0.0 2.0 2.0 2.0 0.0 0.0

 0.0 2.0 10.0 10.0 2.0 2.0 2.0 2.0 0.0 0.0

 2.0 2.0 2.0 10.0 10.0 10.0 2.0 2.0 2.0 0.0

 2.0 2.0 2.0 10.0 10.0 10.0 10.0 10.0 10.0 10.0

 2.0 2.0 2.0 10.0 10.0 10.0 10.0 10.0 2.0 2.0

 2.0 2.0 2.0 10.0 10.0 10.0 10.0 10.0 2.0 0.0

 0.0 2.0 2.0 10.0 10.0 10.0 10.0 2.0 2.0 0.0

 0.0 2.0 2.0 10.0 10.0 10.0 10.0 2.0 2.0 0.0

 0.0 2.0 2.0 10.0 10.0 10.0 10.0 2.0 2.0 2.0

 0.0 2.0 2.0 10.0 10.0 10.0 10.0 2.0 2.0 0.0

 2.0 2.0 10.0 10.0 10.0 10.0 10.0 2.0 2.0 0.0

10.0 10.0 10.0 10.0 10.0 10.0 10.0 10.0 10.0 0.0

 0.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 10.0 0.0

 0.0 0.0 2.0 2.0 2.0 2.0 0.0 0.0 0.0 0.0

CONSTANT 1.E-6 VKA ← 7

CONSTANT 1.E-6 Ss ← 8

INTERNAL 0.2 (10F4.0) 0 Sy ← 9

 1.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ← 9

 0.5 1.0 1.0 0.0 0.0 0.5 0.5 0.5 0.0 0.0

 0.0 0.5 1.0 1.0 0.5 0.5 0.5 0.5 0.0 0.0

 0.5 0.5 0.5 1.0 1.0 1.0 0.5 0.5 0.5 0.0

 0.5 0.5 0.5 1.0 1.0 1.0 1.0 1.0 1.0 1.0

 0.5 0.5 0.5 1.0 1.0 1.0 1.0 1.0 0.5 0.5

 0.5 0.5 0.5 1.0 1.0 1.0 1.0 1.0 0.5 0.0

 0.0 0.5 0.5 1.0 1.0 1.0 1.0 0.5 0.5 0.0

 0.0 0.5 0.5 1.0 1.0 1.0 1.0 0.5 0.5 0.0

 0.0 0.5 0.5 1.0 1.0 1.0 1.0 0.5 0.5 0.5

 0.0 0.5 0.5 1.0 1.0 1.0 1.0 0.5 0.5 0.0

 0.5 0.5 1.0 1.0 1.0 1.0 1.0 0.5 0.5 0.0

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 0.0 0.5 0.5 0.5 0.5 0.5 0.5 0.5 1.0 0.0

 0.0 0.0 0.5 0.5 0.5 0.5 0.0 0.0 0.0 0.0

1 Comment lines. 

2 Flag for cell-by-cell flow terms; head in cell for wet-dry simulations; and flag related to wet-dry simulations.

3 Layer type; method for computing transmissivity; and horizontal anisotropy.

4 Flag to determine vertical hydraulic conductivity.

5 Flag to determine if rewetting is active. 

6 Hydraulic conductivity data for layer 1 read from following data. 

7 Ratio of vertical to horizontal hydraulic conductivity.

8 Specific storage for layer 1.

9 Specific yield for layer 1.

36 Documentation of the Unsaturated-Zone Flow Package
Following (enclosed in a border) are the contents of the Unsaturated-Zone Flow Package input file for test simulation 2; 

explanations are noted outside of border:

File name: UZFtest2.uzf

#UZF Package input for test simulation 2 ← 1

1 1 1 1 -1 61 25 20 4 ← 2

INTERNAL 1 (20I4) 6 UZFBND ← 3

 1 1 0 0 0 0 0 0 0 0 ← 3

 1 1 1 0 0 1 1 1 0 0

 0 1 1 1 1 1 1 1 0 0

 1 1 1 1 1 1 1 1 1 0

 1 1 1 1 1 1 1 1 1 1

 1 1 1 1 1 1 1 1 1 1

 1 1 1 1 1 1 1 1 1 0

 0 1 1 1 1 1 1 1 1 0

 0 1 1 1 1 1 1 1 1 0

 0 1 1 1 1 1 1 1 1 1

 0 1 1 1 1 1 1 1 1 0

 1 1 1 1 1 1 1 1 1 0

 1 1 1 1 1 1 1 1 1 0

 0 1 1 1 1 1 1 1 1 0

 0 0 1 1 1 1 0 0 0 0

INTERNAL 1 (20I4) 6 IRUNBND ← 4

 1 1 0 0 0 0 0 0 0 0 ← 4

 1 1 1 0 0 3 3 3 0 0

 0 1 1 1 3 3 3 3 0 0

 8 8 2 2 2 2 2 3 3 0

 8 8 8 8 2 2 2 3 4 4

 8 8 8 8 8 2 2 4 4 4

 8 8 8 8 8 2 5 5 5 0

 0 8 8 8 8 2 5 5 5 0

 0 8 8 8 8 6 6 6 6 0

 0 8 8 8 8 6 6 6 6 6

 0 8 8 8 8 8 6 6 7 0

 8 8 8 8 8 8 6 7 7 0

 8 8 8 8 8 8 7 7 7 0

 0 8 8 8 8 8 7 7 7 0

 0 0 8 8 8 8 0 0 0 0

INTERNAL 1.0E-6 (10F4.0) 0 VKS ← 5

 1.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ← 5

 1.0 1.0 1.0 0.0 0.0 1.0 1.0 1.0 0.0 0.0

 0.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0 0.0

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 0.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 0.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 0.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0

 0.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 0.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 0.0 0.0 1.0 1.0 1.0 1.0 0.0 0.0 0.0 0.0

Appendix 2 37
File name: UZFtest2.uzf (continued)

CONSTANT 3.5 BROOKS-COREY EPSILON ← 6

CONSTANT 0.30 THTS ← 7

 3 6 65 1 IUZROW IUZCOL IFTUNIT IUZOPT ← 8

 6 3 66 2 IUZROW IUZCOL IFTUNIT IUZOPT

10 5 67 3 IUZROW IUZCOL IFTUNIT IUZOPT

-68 IFTUNIT

1 NUZF1 STRESS PERIOD 1 ← 9

INTERNAL 1.000E-09 (10F4.0) 0 FINF STRESS PERIOD 1 ← 9

 1.0 2.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ← 9

 2.0 1.0 2.0 0.0 0.0 8.0 8.0 5.0 0.0 0.0

 0.0 5.0 1.0 1.0 8.0 8.0 5.0 5.0 0.0 0.0

 8.0 5.0 1.0 1.0 1.0 1.0 1.0 1.0 2.0 0.0

 8.0 5.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0

 8.0 5.0 5.0 2.0 2.0 1.0 1.0 1.0 1.0 5.0

 8.0 5.0 5.0 2.0 2.0 1.0 1.0 2.0 2.0 0.0

 0.0 5.0 5.0 2.0 2.0 1.0 1.0 2.0 5.0 0.0

 0.0 5.0 5.0 2.0 2.0 1.0 1.0 2.0 5.0 0.0

 0.0 5.0 5.0 2.0 2.0 1.0 1.0 2.0 5.0 8.0

 0.0 1.0 1.0 1.0 1.0 1.0 1.0 2.0 5.0 0.0

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 2.0 0.0

 1.0 1.0 2.0 2.0 2.0 2.0 2.0 1.0 1.0 0.0

 0.0 2.0 5.0 5.0 5.0 5.0 5.0 2.0 1.0 0.0

 0.0 0.0 8.0 8.0 8.0 8.0 0.0 0.0 0.0 0.0

1 NUZF2 STRESS PERIOD 1 ←10

CONSTANT 5.0E-08 PET STRESS PERIOD 1 ←10

1 NUZF3 STRESS PERIOD 1 ←11

CONSTANT 15.0 EXTDP STRESS PERIOD 1 ←11

1 NUZF4 STRESS PERIOD 1 ←12

INTERNAL 0.100005 (10F4.0) 0 EXTWC STRESS PERIOD 1 ←12

 1.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ←12

 2.0 1.0 1.0 0.0 0.0 2.0 2.0 2.0 0.0 0.0

 0.0 2.0 1.0 1.0 2.0 2.0 2.0 2.0 0.0 0.0

 2.0 2.0 2.0 1.0 1.0 1.0 2.0 2.0 2.0 1.0

 2.0 2.0 2.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0

 2.0 2.0 2.0 1.0 1.0 1.0 1.0 1.0 2.0 2.0

 2.0 2.0 2.0 1.0 1.0 1.0 1.0 1.0 2.0 0.0

 0.0 2.0 2.0 1.0 1.0 1.0 1.0 2.0 2.0 0.0

 0.0 2.0 2.0 1.0 1.0 1.0 1.0 2.0 2.0 0.0

 0.0 2.0 2.0 1.0 1.0 1.0 1.0 2.0 2.0 2.0

 0.0 2.0 2.0 1.0 1.0 1.0 1.0 2.0 2.0 0.0

 2.0 2.0 1.0 1.0 1.0 1.0 1.0 2.0 2.0 0.0

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0

 0.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 1.0 0.0

 0.0 0.0 2.0 2.0 2.0 2.0 0.0 0.0 0.0 0.0 

...

Input data for stress periods 2-11 are not shown.

38 Documentation of the Unsaturated-Zone Flow Package
File name: UZFtest2.uzf (continued)

1 NUZF1 STRESS PERIOD 12 ←9

INTERNAL 1.000E-9 (10F4.0) 0 FINF STRESS PERIOD 12 ←9

 1.0 2.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ←9

 2.0 1.0 2.0 0.0 0.0 8.0 8.0 5.0 0.0 0.0 

 0.0 5.0 1.0 1.0 8.0 8.0 5.0 5.0 0.0 0.0 

 8.0 5.0 1.0 1.0 1.0 1.0 1.0 1.0 2.0 0.0 

 8.0 5.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 

 8.0 5.0 5.0 2.0 2.0 1.0 1.0 1.0 1.0 5.0 

 8.0 5.0 5.0 2.0 2.0 1.0 1.0 2.0 2.0 0.0 

 0.0 5.0 5.0 2.0 2.0 1.0 1.0 2.0 5.0 0.0 

 0.0 5.0 5.0 2.0 2.0 1.0 1.0 2.0 5.0 0.0 

 0.0 5.0 5.0 2.0 2.0 1.0 1.0 2.0 5.0 8.0 

 0.0 1.0 1.0 1.0 1.0 1.0 1.0 2.0 5.0 0.0 

 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 2.0 0.0 

 1.0 1.0 2.0 2.0 2.0 2.0 2.0 1.0 1.0 0.0 

 0.0 2.0 5.0 5.0 5.0 5.0 5.0 2.0 1.0 0.0 

 0.0 0.0 8.0 8.0 8.0 8.0 0.0 0.0 0.0 0.0

-1 NUZF2 STRESS PERIOD 12 ←10

-1 NUZF3 STRESS PERIOD 12 ←11

-1 NUZF4 STRESS PERIOD 12 ←12

 1 Comment lines. 

 2 Flags for unsaturated flow options; print options; number of trailing waves; trailing sets; and the number of 

 unsaturated zone cells for printing time-series data.

 3 Boundary array for unsaturated flow.

 4 Array that specifies stream segments for adding ground-water discharge to streams.

 5 Vertical hydraulic conductivity of the unsaturated zone.

 6 Brooks-Corey exponent. 

 7 Saturated water content. Note–no THTI because the first stress period is steady state.

 8 Row and column indices; file unit number; and gage option for printing unsaturated flow results.

 9 Applied infiltration rate. 

10 Evapotranspiration demand rate.

11 Evapotranspiration extinction depth.

12 Evapotranspiration extinction water content.

Appendix 2 39
Following (enclosed in a border) are the contents of the Well Package input file for test simulation 2; explanations are noted 

outside of border:

File name: UZFtest2.wel

#Well Package input for test simulation 2 ←1

10 -1 MXACTW IWELCB ←2

10 0 ITMP NP -- Stress Period 1 ←3

 1 6 4 -2.00 Layer Row Column Q ←3

 1 6 5 -2.00

 1 7 4 -2.00

 1 7 5 -2.00

 1 8 4 -2.00

 1 8 5 -2.00

 1 9 4 -2.00

 1 9 5 -2.00

 1 10 4 -2.00

 1 10 5 -2.00

...

10 0 ITMP NP -- Stress Period 11 ←3

 1 6 4 -3.00 Layer Row Column Q ←3

 1 6 5 -3.00

 1 7 4 -3.00

 1 7 5 -3.00

 1 8 4 -3.00

 1 8 5 -3.00

 1 9 4 -3.00

 1 9 5 -3.00

 1 10 4 -3.00

 1 10 5 -3.00

 0 0 ITMP NP -- Stress Period 12 ←3

1 Comment lines. 

2 Maximum number of wells used during stress period; and flag for printing.

3 Number of non-parameter wells read for current stress period; number of parameters used in the 

 current stress period; and layer, row, and column for well and pumping rate.

Following (enclosed in a border) are the contents of the Strongly Implicit Procedure Package input file for test simulation 2; 

explanations are noted outside of border:

File name: UZFtest2.sip

 300 5 MXITER NPARM ←1

 1.0 2.0E-04 1 0.0 1 ←2

1 Maximum iterations; and number of iteration parameters.

2 Acceleration parameter; head change criterion; flag for seed; seed; and printout interval for SIP.

40 Documentation of the Unsaturated-Zone Flow Package
Following (enclosed in a border) are the contents of the General-Head Boundary Package input file for test simulation 2; 

explanations are noted outside of border:

File name: UZFtest2.ghb

 2 -1 MXACTB IGHBCB ← 1

 2 0 ITMP NP -- Stress Period 1 ← 2

 1 13 1 988.0 0.038 ← 3

 1 14 9 1045.0 0.038

 -1 0 ITMP NP -- Stress Period 2 ← 2

... 

 -1 0 ITMP NP -- Stress Period 12 ← 2

1 Maximum number of general-head boundary cells; and unit number for storing cell-by-cell flow terms.

2 Read flags for stress period.

3 General-head boundary cells for stress period 1.

Following (enclosed in a border) are the contents of the Streamflow-Routing Package input file for test simulation 2; 

explanations are noted outside of border:

File name: UZFtest2.sfr

#SFR2 Package input for Test simulation 2 ← 1

36 8 0 0 1.486 0.000001 -1 81 ← 2

 1 1 1 1 1 4500. ← 3

 1 2 2 1 2 7000. ← 3

 1 3 3 1 3 6000. ← 3

 1 3 4 1 4 5550. ← 3

 1 4 5 2 1 6500. ← 3

 1 5 6 2 2 5000. ← 3

 1 6 6 2 3 5000. ← 3

 1 7 6 2 4 5000. ← 3

 1 8 6 2 5 5000. ← 3

 1 3 5 3 1 5000. ← 3

 1 3 6 3 2 5000. ← 3

 1 3 7 3 3 4500. ← 3

 1 4 8 3 4 6000. ← 3

 1 5 8 3 5 5000. ← 3

 1 6 8 3 6 2000. ← 3

 1 5 10 4 1 2500. ← 3

 1 5 9 4 2 5000. ← 3

 1 6 8 4 3 3500. ← 3

 1 6 8 5 1 4000. ← 3 

 1 7 7 5 2 5000. ← 3

 1 8 7 5 3 3500. ← 3 

 1 8 6 5 4 2500. ← 3

 1 9 6 6 1 5000. ← 3

 1 10 7 6 2 5000. ← 3

 1 11 7 6 3 5000. ← 3

 1 12 7 6 4 5000. ← 3

 1 13 7 6 5 2000. ← 3

 1 14 9 7 1 5000. ← 3

 1 13 8 7 2 5500. ← 3

 1 13 7 7 3 5000. ← 3

 1 13 6 8 1 5000. ← 3

 1 13 5 8 2 5000. ← 3

 1 13 4 8 3 5000. ← 3

 1 13 3 8 4 5000. ← 3

 1 13 2 8 5 5000. ← 3

 1 13 1 8 6 3000. ← 3

Appendix 2 41
File name: UZFtest2.sfr (continued)

8 0 0 ← 4

 1 4 2 0 11 25.0 0.0 0.0E-08 0.0E-08 ← 5

 0.00003 3. 1095. ← 6

 0.00003 3. 1075. ← 6

 0.5 1.0 2.0 4.0 7.0 10.0 20.0 30.0 50.0 75.0 100. ← 7

 0.25 0.4 0.55 0.7 0.8 0.9 1.1 1.25 1.4 1.7 1.8 ← 7

 3.0 3.5 4.2 5.3 7.0 8.5 12. 14.0 17.0 20.0 22.0 ← 7

 2 1 6 0 0.0 0.0 0.0E-08 0.0E-08 0.030 ← 5

 0.00003 3. 1075. 12. ← 6

 0.00003 3. 1050. 12. ← 6

 3 0 5 1 0 10.0 0.0 0.0E-08 0.0E-08 ← 5

 0.00003 2. 1075. 10. 2. ← 6

 0.00003 2. 1060. 6. 1. ← 6

 4 1 5 0 10.0 0.0 0.0E-08 0.0E-08 0.030 ← 5

 0.00003 3. 1080. 10. ← 6

 0.00003 3. 1060. 10. ← 6

 5 3 6 0 0.0 0.0 0.0E-08 0.0E-08 0.3 0.35 3.8 0.6 ← 5

 0.00003 3. 1060. ← 6

 0.00003 3. 1045. ← 6

 6 1 8 0 0.0 0.0 0.0E-08 0.0E-08 0.030 ← 5

 0.00003 3. 1045. 12. ← 6

 0.00003 3. 1025. 12. ← 6

 7 2 8 0 150.0 0.0 0.0E-08 0.0E-08 0.025 0.045 ← 5

 0.00006 3. 1040. ← 6

 0.00006 3. 1025. ← 6

 0.0 10. 80. 100. 150. 170. 240. 250. ← 8

 20. 13. 10. 2. 0. 10. 13. 20. ← 8

 8 2 0 0 0.0 0.0 0.0E-08 0.0E-08 0.025 0.045 ← 5 

 0.00006 3. 1025. ← 6 

 0.00006 3. 990. ← 6

 0.0 10. 80. 100. 150. 170. 240. 250. ← 8 

 25. 17. 13. 4. 0. 10. 16. 20. ← 8

-1 0 0 Stress Period 2 ← 9

... -1 0 0 Stress Period 12 ← 9

1 Optional text line for comments. 

2 Basic stream network description.

3 Stream reach descriptions.

4 Read and print flags for stress period.

5 General stream segment data.

6 Stream segment data for upstream and downstream ends.

7 Stream segment data when relation of flow to depth and width read from table.

8 Stream segment data for 8-point cross section.

9 Read and print flags for stress periods 2–12.

42 Documentation of the Unsaturated-Zone Flow Package
Following (enclosed in a border) are the contents of the Gage package input file test simulation 2; explanations are noted as 

comments at the end of each record:

File name: UZFtest1.gag

9 NUMGAGE (number of gaging stations)

1 4 83 4 Stream Segment & Reach #; Output file for Gage 1; Output option

2 3 84 1 Stream Segment & Reach #; Output file for Gage 2; Output option

3 6 85 1 Stream Segment & Reach #; Output file for Gage 3; Output option

4 3 86 1 Stream Segment & Reach #; Output file for Gage 4; Output option

5 4 87 4 Stream Segment & Reach #; Output file for Gage 5; Output option

6 5 88 1 Stream Segment & Reach #; Output file for Gage 6; Output option

7 3 89 4 Stream Segment & Reach #; Output file for Gage 7; Output option

8 6 90 4 Stream Segment & Reach #; Output file for Gage 8; Output option

3 1 91 5 Stream Segment & Reach #; Output file for Gage 9; Output option

Following (enclosed in a border) are the contents of the Output-Control Option input file for test simulation 2. 

File name: UZFtest2.oc

HEAD PRINT FORMAT 4

DRAWDOWN PRINT FORMAT 0

HEAD SAVE FORMAT (10F10.3) LABEL 

HEAD SAVE UNIT 58 

PERIOD 1 STEP 1 

 PRINT BUDGET 

 SAVE BUDGET 

PERIOD 1 STEP 5 

 PRINT BUDGET 

 SAVE BUDGET 

PERIOD 1 STEP 10 

 PRINT BUDGET 

 SAVE BUDGET 

PERIOD 1 STEP 15 

 PRINT BUDGET 

 PRINT HEAD 

 SAVE HEAD 

 SAVE BUDGE

... 

PERIOD 12 STEP 5 

 PRINT BUDGET 

 SAVE BUDGET 

PERIOD 12 STEP 10 

 PRINT BUDGET 

 SAVE BUDGET 

PERIOD 12 STEP 15 

 PRINT BUDGET 

 PRINT HEAD 

 SAVE HEAD 

 SAVE BUDGET 

Appendix 2 43
Listing of Selected Output Files for Test Simulation 2

Following are the contents of the main MODFLOW listing file for test simulation 1:

File name: UZFtest2.lst

 MODFLOW-2005

 U.S. GEOLOGICAL SURVEY MODULAR FINITE-DIFFERENCE GROUND-WATER FLOW MODEL

 VERSION 1.00 12/15/2005

 This model run combines GLOBAL and LIST output into this single file.

 GLOBAL LISTING FILE: OUT\UZFtest2.lst

 UNIT 9

 OPENING DATA\UZFtest2.ba6

 FILE TYPE: BAS6 UNIT 75 STATUSES: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

 OPENING DATA\UZFtest2.lpf

 FILE TYPE: LPF UNIT 7 STATUSES: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

 OPENING DATA\UZFtest2.dis

 FILE TYPE: DIS UNIT 8 STATUSES: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

 OPENING DATA\UZFtest2.sip

 FILE TYPE: SIP UNIT 13 STATUS: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

 OPENING DATA\UZFtest2.oc

 FILE TYPE: OC UNIT 14 STATUSES: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

 OPENING DATA\UZFtest2.sfr

 FILE TYPE: SFR UNIT 15 STATUSES: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

 OPENING DATA\UZFtest2.wel

 FILE TYPE: WEL UNIT 16 STATUSES: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

 OPENING DATA\UZFtest2.ghb

 FILE TYPE: GHB UNIT 17 STATUSES: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

 OPENING DATA\UZFtest2.uzf

 FILE TYPE: UZF UNIT 19 STATUSES: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

 OPENING DATA\UZFtest2.gag

 FILE TYPE: GAGE UNIT 32 STATUS: OLD 

 FORMAT: FORMATTED ACCESS: SEQUENTIAL 

44 Documentation of the Unsaturated-Zone Flow Package
 OPENING OUT\UZFtest2hd.out

 FILE TYPE: DATA UNIT 58 STATUS: UNKNOWN

 FORMAT:FORMATTED ACCESS:SEQUENTIAL 

 OPENING OUT\UZFtest2.uzfot

 FILE TYPE:DATA(BINARY) UNIT 61 STATUS:UNKNOWN

 FORMAT:BINARY ACCESS:SEQUENTIAL 

 OPENING OUT\UZFtest2.uzf1

 FILE TYPE:DATA UNIT 65 STATUS:UNKNOWN

 FORMAT:FORMATTED ACCESS:SEQUENTIAL 

 OPENING OUT\UZFtest2.uzf2

 FILE TYPE:DATA UNIT 66 STATUS:UNKNOWN

 FORMAT:FORMATTED ACCESS:SEQUENTIAL 

 OPENING OUT\UZFtest2.uzf3

 FILE TYPE:DATA UNIT 67 STATUS:UNKNOWN

 FORMAT:FORMATTED ACCESS:SEQUENTIAL 

 OPENING OUT\UZFtest2.uzf4

 FILE TYPE:DATA UNIT 68 STATUS:UNKNOWN

 FORMAT:FORMATTED ACCESS:SEQUENTIAL 

...

 DISCRETIZATION INPUT DATA READ FROM UNIT 8

 #Test problem 2-- Discretization input 

 #Example problem using the UZF1 Package

 1 LAYERS 15 ROWS 10 COLUMNS

 12 STRESS PERIOD(S) IN SIMULATION

MODEL TIME UNIT IS SECONDS

MODEL LENGTH UNIT IS FEET

 Confining bed flag for each layer:

 0

 DELR = 5000.00 

 DELC = 5000.00 

...

Appendix 2 45
TOP ELEVATION OF LAYER 1

 READING ON UNIT 8 WITH FORMAT: (10F7.0) 

 1 2 3 4 5 6 7 8 9 10

 ...........................................................................

 1 1105. 1105. 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 

 2 1105. 1095. 1105. 0.000 0.000 1100. 1100. 1110. 0.000 0.000 

 3 0.000 1100. 1085. 1080. 1078. 1076. 1073. 1090. 0.000 0.000 

 4 1110. 1100. 1090. 1085. 1076. 1080. 1080. 1070. 1100. 0.000 

 5 1110. 1100. 1085. 1080. 1078. 1074. 1077. 1066. 1075. 1080. 

 6 1105. 1095. 1080. 1077. 1073. 1069. 1070. 1063. 1085. 1095. 

 7 1105. 1090. 1075. 1072. 1070. 1065. 1065. 1069. 1090. 0.000 

 8 0.000 1080. 1070. 1068. 1063. 1060. 1063. 1067. 1090. 0.000 

 9 0.000 1075. 1065. 1062. 1058. 1052. 1055. 1065. 1090. 0.000 

 10 0.000 1070. 1060. 1055. 1046. 1046. 1044. 1060. 1088. 1100. 

 11 0.000 1055. 1050. 1045. 1041. 1042. 1040. 1050. 1075. 0.000 

 12 1027. 1028. 1028. 1032. 1032. 1033. 1032. 1040. 1060. 0.000 

 13 997.0 1004. 1008. 1014. 1020. 1027. 1032. 1038. 1050. 0.000 

 14 0.000 1030. 1040. 1045. 1050. 1050. 1050. 1050. 1044. 0.000 

 15 0.000 0.000 1080. 1080. 1080. 1080. 0.000 0.000 0.000 0.000 

 MODEL LAYER BOTTOM EL. FOR LAYER 1

 READING ON UNIT 8 WITH FORMAT: (10F7.0) 

 1 2 3 4 5 6 7 8 9 10

 .........................................................................

 1 1000. 1000. 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 

 2 1000. 900.0 900.0 0.000 0.000 950.0 950.0 950.0 0.000 0.000 

 3 0.000 850.0 670.0 680.0 680.0 770.0 770.0 890.0 0.000 0.000 

 4 925.0 850.0 650.0 650.0 650.0 680.0 740.0 850.0 950.0 0.000 

 5 925.0 800.0 650.0 600.0 600.0 600.0 650.0 720.0 850.0 950.0 

 6 925.0 800.0 650.0 600.0 550.0 600.0 650.0 720.0 850.0 970.0 

 7 925.0 800.0 650.0 600.0 550.0 600.0 650.0 720.0 950.0 0.000 

 8 0.000 850.0 700.0 650.0 570.0 600.0 700.0 850.0 950.0 0.000 

 9 0.000 850.0 700.0 650.0 540.0 580.0 700.0 850.0 950.0 0.000 

 10 0.000 850.0 700.0 630.0 530.0 580.0 650.0 800.0 900.0 950.0 

 11 0.000 850.0 700.0 620.0 580.0 620.0 700.0 800.0 900.0 0.000 

 12 900.0 830.0 710.0 640.0 600.0 650.0 700.0 800.0 900.0 0.000 

 13 900.0 810.0 750.0 700.0 650.0 700.0 750.0 800.0 900.0 0.000 

 14 0.000 950.0 850.0 850.0 850.0 850.0 850.0 850.0 850.0 0.000 

 15 0.000 0.000 950.0 950.0 950.0 950.0 0.000 0.000 0.000 0.000 

46 Documentation of the Unsaturated-Zone Flow Package
 STRESS PERIOD LENGTH TIME STEPS MULTIPLIER FOR DELT SS FLAG

 ----------------------------------------------------------------------------

 1 2628000. 1 1.000 SS

 2 2628000. 15 1.100 TR

 3 2628000. 15 1.100 TR

 4 2628000. 15 1.100 TR

 5 2628000. 15 1.100 TR

 6 2628000. 15 1.100 TR

 7 2628000. 15 1.100 TR

 8 2628000. 15 1.100 TR

 9 2628000. 15 1.100 TR

 10 2628000. 15 1.100 TR

 11 2628000. 15 1.100 TR

 12 2628000. 15 1.100 TR

 COMBINED STEADY-STATE AND TRANSIENT SIMULATION

...

 LPF -- LAYER PROPERTY FLOW PACKAGE, VERSION 7, 5/2/2005

 INPUT READ FROM UNIT 7

 #Layer-Property Flow Package input for test simulation 2 

 HEAD AT CELLS THAT CONVERT TO DRY= 0.0000 

 No named parameters

 LAYER FLAGS:

 LAYER LAYTYP LAYAVG CHANI LAYVKA LAYWET

 ---------------------------------------------------------------------------

 1 1 0 1.000E+00 0 0

 INTERPRETATION OF LAYER FLAGS:

 INTERBLOCK HORIZONTAL DATA IN

 LAYER TYPE TRANSMISSIVITY ANISOTROPY ARRAY VKA WETTABILITY

 LAYER (LAYTYP) (LAYAVG) (CHANI) (LAYVKA) (LAYWET)

 ---------------------------------------------------------------------------

 1 CONVERTIBLE HARMONIC 1.000E+00 VERTICAL K NON-WETTABLE

 WETTING CAPABILITY IS NOT ACTIVE IN ANY LAYER

Appendix 2 47
 HYD. COND. ALONG ROWS FOR LAYER 1

 READING ON UNIT 7 WITH FORMAT: (10F5.0) 

 1 2 3 4 5 6 7 8 9 10

 .............................................................................

 1 2.E-03 2.E-03 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 

 2 4.E-04 2.E-03 2.E-03 0.000 0.000 4.E-04 4.E-04 4.E-04 0.000 0.000 

 3 0.000 4.E-04 2.E-03 2.E-03 4.E-04 4.E-04 4.E-04 4.E-04 0.000 0.000 

 4 4.E-04 4.E-04 4.E-04 2.E-03 2.E-03 2.E-03 4.E-04 4.E-04 4.E-04 0.000 

 5 4.E-04 4.E-04 4.E-04 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03

 6 4.E-04 4.E-04 4.E-04 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 4.E-04 4.E-04

 7 4.E-04 4.E-04 4.E-04 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 4.E-04 0.000 

 8 0.000 4.E-04 4.E-04 2.E-03 2.E-03 2.E-03 2.E-03 4.E-04 4.E-04 0.000 

 9 0.000 4.E-04 4.E-04 2.E-03 2.E-03 2.E-03 2.E-03 4.E-04 4.E-04 0.000 

 10 0.000 4.E-04 4.E-04 2.E-03 2.E-03 2.E-03 2.E-03 4.E-04 4.E-04 4.E-04

 11 0.000 4.E-04 4.E-04 2.E-03 2.E-03 2.E-03 2.E-03 4.E-04 4.E-04 0.000 

 12 4.E-04 4.E-04 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 4.E-04 4.E-04 0.000 

 13 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 2.E-03 0.000 

 14 0.000 4.E-04 4.E-04 4.E-04 4.E-04 4.E-04 4.E-04 4.E-04 2.E-03 0.000 

 15 0.000 0.000 4.E-04 4.E-04 4.E-04 4.E-04 0.000 0.000 0.000 0.000 

 VERTICAL HYD. COND. = 1.000000E-06 FOR LAYER 1

 SPECIFIC STORAGE = 1.000000E-06 FOR LAYER 1

 SPECIFIC YIELD FOR LAYER 1

 READING ON UNIT 7 WITH FORMAT: (10F5.0) 

 1 2 3 4 5 6 7 8 9 10

 ..........................................................................

 1 0.200 0.200 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 

 2 0.100 0.200 0.200 0.000 0.000 0.100 0.100 0.100 0.000 0.000 

 3 0.000 0.100 0.200 0.200 0.100 0.100 0.100 0.100 0.000 0.000 

 4 0.100 0.100 0.100 0.200 0.200 0.200 0.100 0.100 0.100 0.000 

 5 0.100 0.100 0.100 0.200 0.200 0.200 0.200 0.200 0.200 0.200 

 6 0.100 0.100 0.100 0.200 0.200 0.200 0.200 0.200 0.100 0.100 

 7 0.100 0.100 0.100 0.200 0.200 0.200 0.200 0.200 0.100 0.000 

 8 0.000 0.100 0.100 0.200 0.200 0.200 0.200 0.100 0.100 0.000 

 9 0.000 0.100 0.100 0.200 0.200 0.200 0.200 0.100 0.100 0.000 

 10 0.000 0.100 0.100 0.200 0.200 0.200 0.200 0.100 0.100 0.100 

 11 0.000 0.100 0.100 0.200 0.200 0.200 0.200 0.100 0.100 0.000 

 12 0.100 0.100 0.200 0.200 0.200 0.200 0.200 0.100 0.100 0.000 

 13 0.200 0.200 0.200 0.200 0.200 0.200 0.200 0.200 0.200 0.000 

 14 0.000 0.100 0.100 0.100 0.100 0.100 0.100 0.100 0.200 0.000 

 15 0.000 0.000 0.100 0.100 0.100 0.100 0.000 0.000 0.000 0.000 

...

48 Documentation of the Unsaturated-Zone Flow Package
UZF1 -- UNSATURATED FLOW PACKAGE, VERSION 1.2, 1/14/2006

 INPUT READ FROM UNIT 19

 #UZF Package input for test simulation 2

...

 OPTION 1 -- UZF RECHARGE/ET IN TOP LAYER ONLY 

 AREAL EXTENT OF UZ FLOW

 READING ON UNIT 19 WITH FORMAT: (20I4) 

 1 2 3 4 5 6 7 8 9 10 

 ....................................................

 1 1 1 0 0 0 0 0 0 0 0

 2 1 1 1 0 0 1 1 1 0 0

 3 0 1 1 1 1 1 1 1 0 0

 4 1 1 1 1 1 1 1 1 1 0

 5 1 1 1 1 1 1 1 1 1 1

 6 1 1 1 1 1 1 1 1 1 1

 7 1 1 1 1 1 1 1 1 1 0

 8 0 1 1 1 1 1 1 1 1 0

 9 0 1 1 1 1 1 1 1 1 0

10 0 1 1 1 1 1 1 1 1 1

11 0 1 1 1 1 1 1 1 1 0

12 1 1 1 1 1 1 1 1 1 0

13 1 1 1 1 1 1 1 1 1 0

14 0 1 1 1 1 1 1 1 1 0

15 0 0 1 1 1 1 0 0 0 0

 ROUTING OVERLAND RUNOFF

 READING ON UNIT 19 WITH FORMAT: (20I4) 

 1 2 3 4 5 6 7 8 9 10

 ....................................................

 1 1 1 0 0 0 0 0 0 0 0

 2 1 1 1 0 0 3 3 3 0 0

 3 0 1 1 1 3 3 3 3 0 0

 4 8 8 2 2 2 2 2 3 3 0

 5 8 8 8 8 2 2 2 3 4 4

 6 8 8 8 8 8 2 2 4 4 4

 7 8 8 8 8 8 2 5 5 5 0

 8 0 8 8 8 8 2 5 5 5 0

 9 0 8 8 8 8 6 6 6 6 0

10 0 8 8 8 8 6 6 6 6 6

11 0 8 8 8 8 8 6 6 7 0

12 8 8 8 8 8 8 6 7 7 0

13 8 8 8 8 8 8 7 7 7 0

14 0 8 8 8 8 8 7 7 7 0

15 0 0 8 8 8 8 0 0 0 0

Appendix 2 49
 SATURATED VERTICAL K

READING ON UNIT 19 WITH FORMAT: (10F4.0) 

 1 2 3 4 5 6 7 8 9 10

 ...........................................................................

 1 1E-06 1E-06 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 

 2 1E-06 1E-06 1E-06 0.000 0.000 1E-06 1E-06 1E-06 0.000 0.000 

 3 0.000 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 0.000 0.000 

 4 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 0.000 

 5 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06

 6 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06

 7 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 0.000 

 8 0.000 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 0.000 

 9 0.000 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 0.000 

 10 0.000 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06

 11 0.000 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 0.000 

 12 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 0.000 

 13 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 0.000 

 14 0.000 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 1E-06 0.000 

 15 0.000 0.000 1E-06 1E-06 1E-06 1E-06 0.000 0.000 0.000 0.000 

 BROOKS-COREY EPSILON = 3.50000 

 SATURATED WATER CONTENT = 0.300000 

WATER CONTENT PROFILES WILL BE PRINTED TO SEPARATE FILES FOR SELECTED MODEL CELLS, 

AND TIME STEPS AS DEFINED BY OUTPUT CONTROL PACKAGE

SELECTED MODEL CELLS ARE: 

 ROW NUMBER COLUMN NUMBER FORTRAN UNIT NUMBER OUTPUT OPTION

 3 6 65 1

 6 3 66 2

 10 5 67 3

 0 0 68 4

...

 STRESS PERIOD NO. 1, LENGTH = 2628000. 

 -----------------------------------------------

 NUMBER OF TIME STEPS = 1

 MULTIPLIER FOR DELT = 1.000

 INITIAL TIME STEP SIZE = 2628000. 

...

 SOLUTION BY THE STRONGLY IMPLICIT PROCEDURE

 -------------------------------------------

 MAXIMUM ITERATIONS ALLOWED FOR CLOSURE = 300

 ACCELERATION PARAMETER = 0.80000 

 HEAD CHANGE CRITERION FOR CLOSURE = 0.20000E-03

 SIP HEAD CHANGE PRINTOUT INTERVAL = 1

 CALCULATE ITERATION PARAMETERS FROM MODEL CALCULATED WSEED

50 Documentation of the Unsaturated-Zone Flow Package
 AREAL INFILTRATION RATE

 READING ON UNIT 19 WITH FORMAT: (10F4.0) 

 1 2 3 4 5 6 7 8 9 10

 ...........................................................................

 1 1E-09 2E-09 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 

 2 2E-09 1E-09 2E-09 0.000 0.000 8E-09 8E-09 5E-09 0.000 0.000 

 3 0.000 5E-09 1E-09 1E-09 8E-09 8E-09 5E-09 5E-09 0.000 0.000 

 4 8E-09 5E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 2E-09 0.000 

 5 8E-09 5E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09

 6 8E-09 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 1E-09 1E-09 5E-09

 7 8E-09 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 1E-09 2E-09 0.000 

 8 0.000 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 2E-09 5E-09 0.000 

 9 0.000 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 2E-09 5E-09 0.000 

 10 0.000 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 2E-09 5E-09 8E-09

 11 0.000 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 2E-09 5E-09 0.000 

 12 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 2E-09 0.000 

 13 1E-09 1E-09 2E-09 2E-09 2E-09 2E-09 2E-09 1E-09 1E-09 0.000 

 14 0.000 2E-09 5E-09 5E-09 5E-09 5E-09 5E-09 2E-09 1E-09 0.000 

 15 0.000 0.000 8E-09 8E-09 8E-09 8E-09 0.000 0.000 0.000 0.000 

 ET DEMAND RATE = 5.000000E-08

 ET EXTINCTION DEPTH = 15.0000 

 EXTINCTION WATER CONTENT

 READING ON UNIT 19 WITH FORMAT: (10F4.0) 

 1 2 3 4 5 6 7 8 9 10

 .........................................................................

 1 0.10 0.10 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 

 2 0.20 0.10 0.10 0.00 0.00 0.20 0.20 0.20 0.00 0.00 

 3 0.00 0.20 0.10 0.10 0.20 0.20 0.20 0.20 0.00 0.00 

 4 0.20 0.20 0.20 0.10 0.10 0.10 0.20 0.20 0.20 0.10 

 5 0.20 0.20 0.20 0.10 0.10 0.10 0.10 0.10 0.10 0.10 

 6 0.20 0.20 0.20 0.10 0.10 0.10 0.10 0.10 0.20 0.20 

 7 0.20 0.20 0.20 0.10 0.10 0.10 0.10 0.10 0.20 0.00 

 8 0.00 0.20 0.20 0.10 0.10 0.10 0.10 0.20 0.20 0.00 

 9 0.00 0.20 0.20 0.10 0.10 0.10 0.10 0.20 0.20 0.00 

 10 0.00 0.20 0.20 0.10 0.10 0.10 0.10 0.20 0.20 0.20 

 11 0.00 0.20 0.20 0.10 0.10 0.10 0.10 0.20 0.20 0.00 

 12 0.20 0.20 0.10 0.10 0.10 0.10 0.10 0.20 0.20 0.00 

 13 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.00 

 14 0.00 0.20 0.20 0.20 0.20 0.20 0.20 0.20 0.10 0.00 

 15 0.00 0.00 0.20 0.20 0.20 0.20 0.00 0.00 0.00 0.00 

Appendix 2 51
SOLVING FOR HEAD 

 AVERAGE SEED = 0.00926840

 MINIMUM SEED = 0.00187768

 5 ITERATION PARAMETERS CALCULATED FROM AVERAGE SEED:

 0.000000E+00 0.689722E+00 0.903727E+00 0.970129E+00 0.990732E+00

 27 ITERATIONS FOR TIME STEP 1 IN STRESS PERIOD 1

...

 OUTPUT CONTROL FOR STRESS PERIOD 1 TIME STEP 1

 PRINT BUDGET

 SAVE BUDGET

 WELLS PERIOD 1 STEP 1

 WELL 1 LAYER 1 ROW 6 COL 4 RATE -2.00000 

 WELL 2 LAYER 1 ROW 6 COL 5 RATE -2.00000 

 WELL 3 LAYER 1 ROW 7 COL 4 RATE -2.00000 

 WELL 4 LAYER 1 ROW 7 COL 5 RATE -2.00000 

 WELL 5 LAYER 1 ROW 8 COL 4 RATE -2.00000 

 WELL 6 LAYER 1 ROW 8 COL 5 RATE -2.00000 

 WELL 7 LAYER 1 ROW 9 COL 4 RATE -2.00000 

 WELL 8 LAYER 1 ROW 9 COL 5 RATE -2.00000 

 WELL 9 LAYER 1 ROW 10 COL 4 RATE -2.00000 

 WELL 10 LAYER 1 ROW 10 COL 5 RATE -2.00000 

 HEAD DEP BOUNDS PERIOD 1 STEP 1

 BOUNDARY 1 LAYER 1 ROW 13 COL 1 RATE -0.262024 

 BOUNDARY 2 LAYER 1 ROW 14 COL 9 RATE 0.209432

 

 UBDSV3 SAVING “ GW ET “ ON UNIT 61 AT TIME STEP 1, STRESS PERIOD 1

 UBDSV3 SAVING “ UZF RECHARGE “ ON UNIT 61 AT TIME STEP 1, STRESS PERIOD 1

 UBDSV3 SAVING “ SURFACE LEAKAGE “ ON UNIT 61 AT TIME STEP 1, STRESS PERIOD 1

UNSATURATED ZONE PACKAGE VOLUMETRIC BUDGET FOR TIME STEP 1 STRESS PERIOD 1

 ------------------------------------------------------------------------------

 CUMULATIVE VOLUMES L**3 RATES FOR THIS TIME STEP L**3/T

 ------------------ ------------------------

 IN: IN:

 --- ---

 

 INFILTRATION 22272299.5597 INFILTRATION = 8.4750

 OUT: OUT:

 ---- ----

 UZ ET = 0.0000E+00 UZ ET = 0.0000E+00

 UZF RECHARGE = 22272299.5597 UZF RECHARGE = 8.4750

 IN - OUT = 0.0000E+00 IN - OUT = 0.0000E+00

 STORAGE: STORAGE:

 -------- --------

 STORAGE CHANGE = 0.0000E+00 STORAGE CHANGE = 0.0000E+00

 PERCENT DISCREPANCY IS DIFFERENCE BETWEEN IN-OUT MINUS CHANGE IN STORAGE

 DIVIDED BY THE AVERAGE OF IN AND OUT TIMES 100

 PERCENT DISCREPANCY = 0.00 PERCENT DISCREPANCY = 0.00

52 Documentation of the Unsaturated-Zone Flow Package
...

 VOLUMETRIC BUDGET FOR ENTIRE MODEL AT END OF TIME STEP 1 IN STRESS PERIOD 1

 ------------------------------------------------------------------------------

 CUMULATIVE VOLUMES L**3 RATES FOR THIS TIME STEP L**3/T

 ------------------ ------------------------

 IN: IN:

 --- ---

 STORAGE = 0.0000 STORAGE = 0.0000

 CONSTANT HEAD = 0.0000 CONSTANT HEAD = 0.0000

 WELLS = 0.0000 WELLS = 0.0000

 HEAD DEP BOUNDS = 551734.6875 HEAD DEP BOUNDS = 0.2099

 UZF RECHARGE = 22272300.0000 UZF RECHARGE = 8.4750

 GW ET = 0.0000 GW ET = 0.0000

 SURFACE LEAKAGE = 0.0000 SURFACE LEAKAGE = 0.0000

 STREAM LEAKAGE = 116098504.0000 STREAM LEAKAGE = 44.1775

 TOTAL IN = 138922544.0000 TOTAL IN = 52.8625

 OUT: OUT:

 ---- ----

 STORAGE = 0.0000 STORAGE = 0.0000

 CONSTANT HEAD = 0.0000 CONSTANT HEAD = 0.0000

 WELLS = 52560000.0000 WELLS = 20.0000

 HEAD DEP BOUNDS = 686396.1875 HEAD DEP BOUNDS = 0.2612

 UZF RECHARGE = 0.0000 UZF RECHARGE = 0.0000

 GW ET = 53700992.0000 GW ET = 20.4342

 SURFACE LEAKAGE = 0.0000 SURFACE LEAKAGE = 0.0000

 STREAM LEAKAGE = 31978764.0000 STREAM LEAKAGE = 12.1685

 TOTAL OUT = 138926160.0000 TOTAL OUT = 52.8638

 IN - OUT = -3616.0000 IN - OUT = -1.3771E-03

 PERCENT DISCREPANCY = 0.00 PERCENT DISCREPANCY = 0.00

 TIME SUMMARY AT END OF TIME STEP 1 IN STRESS PERIOD 1

 SECONDS MINUTES HOURS DAYS YEARS

 -----------------------------------------------------------

 TIME STEP LENGTH 2.62800E+06 43800. 730.00 30.417 8.32763E-02

 STRESS PERIOD TIME 2.62800E+06 43800. 730.00 30.417 8.32763E-02

 TOTAL TIME 2.62800E+06 43800. 730.00 30.417 8.32763E-02

Appendix 2 53
...

 STRESS PERIOD NO. 3, LENGTH = 2628000. 

 -----------------------------------------------

 NUMBER OF TIME STEPS = 15

 MULTIPLIER FOR DELT = 1.100

 INITIAL TIME STEP SIZE = 82713.07 

...

 AREAL INFILTRATION RATE

 READING ON UNIT 19 WITH FORMAT: (10F4.0) 

 1 2 3 4 5 6 7 8 9 10

 ..........................................................................

 1 1E-08 2E-08 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 

 2 2E-08 1E-08 2E-08 0.000 0.000 8E-08 8E-08 5E-08 0.000 0.000 

 3 0.000 5E-08 1E-08 1E-08 8E-08 8E-08 5E-08 5E-08 0.000 0.000 

 4 8E-08 5E-08 1E-08 1E-08 1E-08 1E-08 1E-08 1E-08 2E-08 0.000 

 5 8E-08 5E-08 1E-08 1E-08 1E-08 1E-08 1E-08 1E-08 1E-08 1E-08

 6 8E-08 5E-08 5E-08 2E-08 2E-08 1E-08 1E-08 1E-08 1E-08 5E-08

 7 8E-08 5E-08 5E-08 2E-08 2E-08 1E-08 1E-08 2E-08 2E-08 0.000 

 8 0.000 5E-08 5E-08 2E-08 2E-08 1E-08 1E-08 2E-08 5E-08 0.000 

 9 0.000 5E-08 5E-08 2E-08 2E-08 1E-08 1E-08 2E-08 5E-08 0.000 

 10 0.000 5E-08 5E-08 2E-08 2E-08 1E-08 1E-08 2E-08 5E-08 8E-08

 11 0.000 1E-08 1E-08 1E-08 1E-08 1E-08 1E-08 2E-08 5E-08 0.000 

 12 1E-08 1E-08 1E-08 1E-08 1E-08 1E-08 1E-08 1E-08 2E-08 0.000 

 13 1E-08 1E-08 2E-08 2E-08 2E-08 2E-08 2E-08 1E-08 1E-08 0.000 

 14 0.000 2E-08 5E-08 5E-08 5E-08 5E-08 5E-08 2E-08 1E-08 0.000 

 15 0.000 0.000 8E-08 8E-08 8E-08 8E-08 0.000 0.000 0.000 0.000 

USING ET DEMAND RATE FROM PREVIOUS STRESS PERIOD. CURRENT PERIOD IS: 3

USING ET EXTINCTION DEPTH FROM PREVIOUS STRESS PERIOD. CURRENT PERIOD IS: 3

USING EXTINCTION WATER CONTENT FROM PREVIOUS STRESS PERIOD.

CURRENT PERIOD IS: 3

54 Documentation of the Unsaturated-Zone Flow Package
... 

 UBDSV3 SAVING “ GW ET “ ON UNIT 61 AT TIME STEP 15, STRESS PERIOD 3

 UBDSV3 SAVING “ UZF RECHARGE “ ON UNIT 61 AT TIME STEP 15, STRESS PERIOD 3

 UBDSV3 SAVING “ SURFACE LEAKAGE “ ON UNIT 61 AT TIME STEP 15, STRESS PERIOD 3

UNSATURATED ZONE PACKAGE VOLUMETRIC BUDGET FOR TIME STEP 15 STRESS PERIOD 3

 ------------------------------------------------------------------------------

 CUMULATIVE VOLUMES L**3 RATES FOR THIS TIME STEP L**3/T

 ------------------ ------------------------

 IN: IN:

 --- ---

 INFILTRATION = 423173664.0228 INFILTRATION = 84.7500

 OUT: OUT:

 ---- ----

 UZF ET = 640655272.4093 UZF ET = 108.5681

 UZF RECHARGE = 85558487.2589 UZF RECHARGE = 15.9795

 IN - OUT = -303040095.6455 IN - OUT = -39.7976

 STORAGE: STORAGE:

 -------- --------

 STORAGE CHANGE = -302943758.7172 STORAGE CHANGE = -39.7629

 PERCENT DISCREPANCY IS DIFFERENCE BETWEEN IN-OUT MINUS CHANGE IN STORAGE

 DIVIDED BY THE AVERAGE OF IN AND OUT TIMES 100

 PERCENT DISCREPANCY = 0.00 PERCENT DISCREPANCY = -0.03

Appendix 2 55
 VOLUMETRIC BUDGET FOR ENTIRE MODEL AT END OF TIME STEP 15 IN STRESS PERIOD 3

 ------------------------------------------------------------------------------

 CUMULATIVE VOLUMES L**3 RATES FOR THIS TIME STEP L**3/T

 ------------------ ------------------------

 IN: IN:

 --- ---

 STORAGE = 338021.0312 STORAGE = 0.2042

 CONSTANT HEAD = 0.0000 CONSTANT HEAD = 0.0000

 WELLS = 0.0000 WELLS = 0.0000

 HEAD DEP BOUNDS = 1648010.1250 HEAD DEP BOUNDS = 0.2083

 UZF RECHARGE = 85558496.0000 UZF RECHARGE = 15.9795

 GW ET = 0.0000 GW ET = 0.0000

 SURFACE LEAKAGE = 0.0000 SURFACE LEAKAGE = 0.0000

 STREAM LEAKAGE = 344701152.0000 STREAM LEAKAGE = 43.2896

 TOTAL IN = 432245664.0000 TOTAL IN = 59.6816

 OUT: OUT:

 ---- ----

 STORAGE = 70915512.0000 STORAGE = 21.0468

 CONSTANT HEAD = 0.0000 CONSTANT HEAD = 0.0000

 WELLS = 131400000.0000 WELLS = 10.0000

 HEAD DEP BOUNDS = 2065964.2500 HEAD DEP BOUNDS = 0.2630

 UZF RECHARGE = 0.0000 UZF RECHARGE = 0.0000

 GW ET = 126929936.0000 GW ET = 14.6589

 SURFACE LEAKAGE = 0.0000 SURFACE LEAKAGE = 0.0000

 STREAM LEAKAGE = 100745112.0000 STREAM LEAKAGE = 13.6927

 TOTAL OUT = 432056512.0000 TOTAL OUT = 59.6614

 IN - OUT = 189152.0000 IN - OUT = 2.0271E-02

 PERCENT DISCREPANCY = 0.04 PERCENT DISCREPANCY = 0.03

 TIME SUMMARY AT END OF TIME STEP 15 IN STRESS PERIOD 3

 SECONDS MINUTES HOURS DAYS YEARS

 -----------------------------------------------------------

 TIME STEP LENGTH 3.14103E+05 5235.0 87.251 3.6354 9.95332E-03

 STRESS PERIOD TIME 2.62800E+06 43800. 730.00 30.417 8.32763E-02

 TOTAL TIME 7.88400E+06 1.31400E+05 2190.0 ` 91.250 0.24983

...

 STRESS PERIOD NO. 12, LENGTH = 2628000. 

 -----------------------------------------------

 NUMBER OF TIME STEPS = 15

 MULTIPLIER FOR DELT = 1.100

 INITIAL TIME STEP SIZE = 82713.07 

...

56 Documentation of the Unsaturated-Zone Flow Package
 AREAL INFILTRATION RATE

 READING ON UNIT 19 WITH FORMAT: (10F4.0) 

 1 2 3 4 5 6 7 8 9 10

 ..........................................................................

 1 1E-09 2E-09 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 

 2 2E-09 1E-09 2E-09 0.000 0.000 8E-09 8E-09 5E-09 0.000 0.000 

 3 0.000 5E-09 1E-09 1E-09 8E-09 8E-09 5E-09 5E-09 0.000 0.000 

 4 8E-09 5E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 2E-09 0.000 

 5 8E-09 5E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09

 6 8E-09 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 1E-09 1E-09 5E-09

 7 8E-09 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 2E-09 2E-09 0.000 

 8 0.000 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 2E-09 5E-09 0.000 

 9 0.000 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 2E-09 5E-09 0.000 

 10 0.000 5E-09 5E-09 2E-09 2E-09 1E-09 1E-09 2E-09 5E-09 8E-09

 11 0.000 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 2E-09 5E-09 0.000 

 12 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 1E-09 2E-09 0.000 

 13 1E-09 1E-09 2E-09 2E-09 2E-09 2E-09 2E-09 1E-09 1E-09 0.000 

 14 0.000 2E-09 5E-09 5E-09 5E-09 5E-09 5E-09 2E-09 1E-09 0.000 

 15 0.000 0.000 8E-09 8E-09 8E-09 8E-09 0.000 0.000 0.000 0.000 

USING ET DEMAND RATE FROM PREVIOUS STRESS PERIOD. CURRENT PERIOD IS: 12 

USING ET EXTINCTION DEPTH FROM PREVIOUS STRESS PERIOD. CURRENT PERIOD IS: 12

USING EXTINCTION WATER CONTENT FROM PREVIOUS STRESS PERIOD. CURRENT PERIOD IS: 12

 SOLVING FOR HEAD

 7 ITERATIONS FOR TIME STEP 1 IN STRESS PERIOD 12

...

OUTPUT CONTROL FOR STRESS PERIOD 12 TIME STEP 15

 PRINT BUDGET

 PRINT HEAD FOR ALL LAYERS

 SAVE HEAD FOR ALL LAYERS

 SAVE BUDGET

Appendix 2 57
UNSATURATED ZONE PACKAGE VOLUMETRIC BUDGET FOR TIME STEP 15 STRESS PERIOD 12

 ------------------------------------------------------------------------------

 CUMULATIVE VOLUMES L**3 RATES FOR THIS TIME STEP L**3/T

 ------------------ ------------------------

 IN: IN:

 --- ---

 INFILTRATION = 1536788572.4216 INFILTRATION = 8.4750

 OUT: OUT:

 ---- ----

 UZF ET = 2130532183.6154 UZF ET = 29.5398

 UZF RECHARGE = 386682175.9374 UZF RECHARGE = 11.2488

 IN - OUT = -980425787.1312 IN - OUT = -32.3136

 STORAGE: STORAGE:

 -------- --------

 STORAGE CHANGE = -980298516.5008 STORAGE CHANGE = -32.4115

 PERCENT DISCREPANCY IS DIFFERENCE BETWEEN IN-OUT MINUS CHANGE IN STORAGE

 DIVIDED BY THE AVERAGE OF IN AND OUT TIMES 100

 PERCENT DISCREPANCY = 0.00 PERCENT DISCREPANCY = 0.24

...

 HEAD IN LAYER 1 AT END OF TIME STEP 15 IN STRESS PERIOD 12

 ---------------------------------------------------------------------------

 1 2 3 4 5 6 7 8 9 10

...................................................................................

 1 1092.37 1088.61 6999.00 6999.00 6999.00 6999.00 6999.00 6999.00 6999.00 6999.00

 2 1089.01 1085.68 1079.29 6999.00 6999.00 1073.62 1072.54 1070.97 6999.00 6999.00

 3 6999.00 1077.38 1074.52 1070.40 1070.89 1069.73 1068.19 1067.88 6999.00 6999.00

 4 1069.50 1068.75 1067.10 1063.67 1062.51 1061.22 1061.78 1064.46 1063.91 6999.00

 5 1065.39 1061.87 1059.39 1057.54 1057.25 1057.77 1057.60 1059.59 1063.01 1065.39

 6 1060.65 1056.02 1053.11 1051.72 1052.00 1053.31 1054.52 1056.79 1059.55 1063.65

 7 1056.94 1050.38 1047.40 1046.23 1046.83 1048.45 1050.78 1053.10 1054.36 6999.00

 8 6999.00 1042.92 1041.40 1040.73 1041.58 1043.41 1045.19 1047.53 1049.31 6999.00

 9 6999.00 1035.92 1035.34 1035.29 1036.46 1038.47 1039.35 1041.08 1043.88 6999.00

 10 6999.00 1028.78 1029.18 1030.18 1031.62 1033.35 1034.61 1036.48 1040.18 1045.29

 11 6999.00 1019.85 1021.63 1025.14 1027.16 1029.21 1030.68 1032.95 1036.39 6999.00

 12 1001.20 1008.88 1014.17 1019.07 1022.71 1025.93 1028.31 1031.01 1034.88 6999.00

 13 994.91 1000.19 1006.90 1013.08 1018.90 1024.47 1029.24 1034.27 1036.38 6999.00

 14 6999.00 1005.79 1011.36 1016.09 1021.14 1025.87 1030.11 1035.21 1039.48 6999.00

 15 6999.00 6999.00 1018.95 1022.40 1026.57 1029.77 6999.00 6999.00 6999.00 6999.00

 HEAD WILL BE SAVED ON UNIT 58 AT END OF TIME STEP 15, STRESS PERIOD 12

58 Documentation of the Unsaturated-Zone Flow Package
 VOLUMETRIC BUDGET FOR ENTIRE MODEL AT END OF TIME STEP 15 IN STRESS PERIOD 12

 ------------------------------------------------------------------------------

 CUMULATIVE VOLUMES L**3 RATES FOR THIS TIME STEP L**3/T

 ------------------ ------------------------

 IN: IN:

 --- ---

 STORAGE = 128218360.0000 STORAGE = 0.8641

 CONSTANT HEAD = 0.0000 CONSTANT HEAD = 0.0000

 WELLS = 0.0000 WELLS = 0.0000

 HEAD DEP BOUNDS = 6590192.5000 HEAD DEP BOUNDS = 0.2096

 UZF RECHARGE = 386682176.0000 UZF RECHARGE = 11.2488

 GW ET = 0.0000 GW ET = 0.0000

 SURFACE LEAKAGE = 0.0000 SURFACE LEAKAGE = 0.0000

 STREAM LEAKAGE = 1361713920.0000 STREAM LEAKAGE = 43.4029

 TOTAL IN = 1883204608.0000 TOTAL IN = 55.7253

 OUT: OUT:

 ---- ----

 STORAGE = 503634528.0000 STORAGE = 22.7426

 CONSTANT HEAD = 0.0000 CONSTANT HEAD = 0.0000

 WELLS = 394199968.0000 WELLS = 0.0000

 HEAD DEP BOUNDS = 8307182.0000 HEAD DEP BOUNDS = 0.2626

 UZF RECHARGE = 0.0000 UZF RECHARGE = 0.0000

 GW ET = 553775360.0000 GW ET = 19.6268

 SURFACE LEAKAGE = 48913.4141 SURFACE LEAKAGE = 0.0000

 STREAM LEAKAGE = 422406880.0000 STREAM LEAKAGE = 13.0608

 TOTAL OUT = 1882323968.0000 TOTAL OUT = 55.2894

 IN - OUT = 880640.0000 IN - OUT = 3.2516E-02

 PERCENT DISCREPANCY = 0.05 PERCENT DISCREPANCY = 0.06

 TIME SUMMARY AT END OF TIME STEP 15 IN STRESS PERIOD 12

 SECONDS MINUTES HOURS DAYS YEARS

 -----------------------------------------------------------

 TIME STEP LENGTH 3.14103E+05 5235.0 87.251 3.6354 9.95332E-03

 STRESS PERIOD TIME 2.62800E+06 43800. 730.00 30.417 8.32763E-02

 TOTAL TIME 3.15360E+07 5.25600E+05 8760.0 365.00 0.99932

Appendix 2 59
Following are the abridged contents of the UZFtest2.uzf1 output file for test simulation 2.

LOCATION OF SPECIFIED CELL FOR PRINTING VOLUMES IN UNSATURATED ZONE:

GAGE 1 ROW, COLUMN 3, 6 INITIAL LAYER ASSIGMENT 1

 GW UZ CUM. CUM. STORAGE SURF.

 LAYER TIME HEAD THICK. INFILT. RECHARGE TOTAL CHANGE LEAK.

 1 2.6280E+06 1.0692E+03 6.2800E+00 5.2560E+05 5.2560E+05 0.0000E+00 0.0000E+00 0.0

 1 2.7107E+06 1.0692E+03 6.2673E+00 6.5794E+05 5.5014E+05 4.0163E+06 6.4600E+04 0.0

 1 2.8017E+06 1.0692E+03 6.2543E+00 8.0352E+05 5.7573E+05 4.0889E+06 1.3717E+05 0.0

 1 2.9018E+06 1.0693E+03 6.2411E+00 9.6365E+05 6.0231E+05 4.1704E+06 2.1867E+05 0.0

 1 3.0119E+06 1.0693E+03 6.2280E+00 1.1398E+06 6.2958E+05 4.2621E+06 3.1040E+05 0.0

 1 3.1330E+06 1.0693E+03 6.2151E+00 1.3336E+06 6.5746E+05 4.3653E+06 4.1356E+05 0.0

 1 3.2662E+06 1.0693E+03 6.2025E+00 1.5467E+06 6.8581E+05 4.4812E+06 5.2950E+05 0.0

 1 3.4127E+06 1.0693E+03 6.1905E+00 1.7811E+06 7.1424E+05 4.6116E+06 6.5992E+05 0.0

 1 3.5739E+06 1.0693E+03 6.1790E+00 2.0390E+06 7.4267E+05 4.7581E+06 8.0640E+05 0.0

 1 3.7512E+06 1.0693E+03 6.1686E+00 2.3227E+06 7.7007E+05 4.9232E+06 9.7154E+05 0.0

 1 3.9462E+06 1.0693E+03 6.1594E+00 2.6348E+06 7.9663E+05 5.1086E+06 1.1569E+06 0.0

 1 4.1608E+06 1.0693E+03 6.1514E+00 2.9780E+06 8.2207E+05 5.3165E+06 1.3648E+06 0.0

 1 4.3968E+06 1.0694E+03 6.1446E+00 3.3556E+06 8.4609E+05 5.5492E+06 1.5975E+06 0.0

 1 4.6563E+06 1.0694E+03 6.1393E+00 3.7710E+06 8.6799E+05 5.8098E+06 1.8581E+06 0.0

 1 4.9419E+06 1.0694E+03 6.1353E+00 4.2278E+06 8.8781E+05 6.1009E+06 2.1492E+06 0.0

 1 5.2560E+06 1.0694E+03 6.1330E+00 4.7304E+06 9.0511E+05 6.4256E+06 2.4739E+06 0.0

 1 5.3387E+06 1.0694E+03 6.1328E+00 4.8958E+06 9.0892E+05 6.5450E+06 2.5933E+06 0.0

 1 5.4297E+06 1.0694E+03 6.1326E+00 5.0778E+06 9.1282E+05 6.6765E+06 2.7248E+06 0.0

 1 5.5298E+06 1.0694E+03 6.0730E+00 5.2780E+06 1.0730E+06 6.6659E+06 2.7142E+06 0.0

 1 5.6399E+06 1.0695E+03 6.0060E+00 5.4981E+06 1.2621E+06 6.6418E+06 2.6901E+06 0.0

 1 5.7610E+06 1.0696E+03 5.9388E+00 5.7403E+06 1.4641E+06 6.6221E+06 2.6704E+06 0.0

 1 5.8942E+06 1.0696E+03 5.8708E+00 6.0068E+06 1.6826E+06 6.6048E+06 2.6531E+06 0.0

...

 1 2.6363E+07 1.0699E+03 5.6163E+00 3.5232E+07 2.5589E+07 3.5513E+06 -4.0039E+05 0.0

 1 2.6454E+07 1.0699E+03 5.6203E+00 3.5250E+07 2.5598E+07 3.5173E+06 -4.3437E+05 0.0

 1 2.6554E+07 1.0699E+03 5.6241E+00 3.5270E+07 2.5610E+07 3.4793E+06 -4.7243E+05 0.0

 1 2.6664E+07 1.0699E+03 5.6283E+00 3.5292E+07 2.5622E+07 3.4374E+06 -5.1433E+05 0.0

 1 2.6785E+07 1.0699E+03 5.6326E+00 3.5316E+07 2.5635E+07 3.3911E+06 -5.6061E+05 0.0

 1 2.6918E+07 1.0699E+03 5.6370E+00 3.5343E+07 2.5651E+07 3.3399E+06 -6.1183E+05 0.0

 1 2.7065E+07 1.0699E+03 5.6415E+00 3.5372E+07 2.5668E+07 3.2834E+06 -6.6827E+05 0.0

 1 2.7226E+07 1.0699E+03 5.6459E+00 3.5404E+07 2.5686E+07 3.2210E+06 -7.3068E+05 0.0

 1 2.7403E+07 1.0698E+03 5.6509E+00 3.5440E+07 2.5706E+07 3.1534E+06 -7.9825E+05 0.0

 1 2.7598E+07 1.0698E+03 5.6571E+00 3.5479E+07 2.5724E+07 3.0820E+06 -8.6971E+05 0.0

 1 2.7813E+07 1.0698E+03 5.6646E+00 3.5522E+07 2.5741E+07 3.0069E+06 -9.4475E+05 0.0

 1 2.8049E+07 1.0698E+03 5.6733E+00 3.5569E+07 2.5756E+07 2.9278E+06 -1.0239E+06 0.0

 1 2.8308E+07 1.0698E+03 5.6832E+00 3.5621E+07 2.5768E+07 2.8444E+06 -1.1073E+06 0.0

 1 2.8594E+07 1.0698E+03 5.6939E+00 3.5678E+07 2.5779E+07 2.7557E+06 -1.1960E+06 0.0

 1 2.8908E+07 1.0698E+03 5.7052E+00 3.5741E+07 2.5787E+07 2.6609E+06 -1.2908E+06 0.0

 1 2.8991E+07 1.0698E+03 5.7083E+00 3.5757E+07 2.5789E+07 2.6365E+06 -1.3152E+06 0.0

 1 2.9082E+07 1.0698E+03 5.7116E+00 3.5776E+07 2.5790E+07 2.6097E+06 -1.3420E+06 0.0

 1 2.9182E+07 1.0698E+03 5.7151E+00 3.5796E+07 2.5792E+07 2.5804E+06 -1.3713E+06 0.0

 1 2.9292E+07 1.0698E+03 5.7189E+00 3.5818E+07 2.5793E+07 2.5484E+06 -1.4033E+06 0.0

 1 2.9413E+07 1.0698E+03 5.7230E+00 3.5842E+07 2.5795E+07 2.5134E+06 -1.4383E+06 0.0

 1 2.9546E+07 1.0698E+03 5.7272E+00 3.5868E+07 2.5796E+07 2.4750E+06 -1.4767E+06 0.0

 1 2.9693E+07 1.0698E+03 5.7318E+00 3.5898E+07 2.5798E+07 2.4330E+06 -1.5187E+06 0.0

 1 2.9854E+07 1.0698E+03 5.7365E+00 3.5930E+07 2.5799E+07 2.3870E+06 -1.5647E+06 0.0

 1 3.0031E+07 1.0698E+03 5.7417E+00 3.5965E+07 2.5800E+07 2.3368E+06 -1.6149E+06 0.0

 1 3.0226E+07 1.0698E+03 5.7470E+00 3.6004E+07 2.5800E+07 2.2818E+06 -1.6698E+06 0.0

 1 3.0441E+07 1.0697E+03 5.7525E+00 3.6047E+07 2.5801E+07 2.2215E+06 -1.7302E+06 0.0

 1 3.0677E+07 1.0697E+03 5.7579E+00 3.6095E+07 2.5801E+07 2.1552E+06 -1.7965E+06 0.0

 1 3.0936E+07 1.0697E+03 5.7635E+00 3.6146E+07 2.5801E+07 2.0825E+06 -1.8692E+06 0.0

 1 3.1222E+07 1.0697E+03 5.7691E+00 3.6204E+07 2.5801E+07 2.0024E+06 -1.9493E+06 0.0

 1 3.1536E+07 1.0697E+03 5.7746E+00 3.6266E+07 2.5801E+07 1.9142E+06 -2.0375E+06 0.0

60 Documentation of the Unsaturated-Zone Flow Package
Following are the abridged contents of the UZFtest2.uzf3 output file for test simulation 2.

LOCATION OF SPECIFIED CELL FOR PRINTING VOLUMES IN UNSATURATED ZONE:

GAGE 3 ROW, COLUMN 10, 5 INITIAL LAYER ASSIGMENT 1

 LAYER TIME GW HEAD UZ THICK. DEPTH WATER CONT.

...

 1 3.1536E+07 1031.6185 13.8814 0.34702721 0.13035826

 0.69405442 0.12027895

 1.04108167 0.12133145

 1.38810885 0.11728376

 1.73513603 0.11856612

 2.08216333 0.11264190

 2.42919040 0.11277605

 2.77621770 0.10000502

 3.12324476 0.10075631

 3.47027206 0.10229551

 3.81729937 0.10405234

 4.16432667 0.10543815

 4.51135349 0.10665397

 4.85838079 0.10000500

 5.20540810 0.10000503

 5.55243540 0.10000505

 5.89946270 0.10000508

 6.24648952 0.10000493

 6.59351683 0.10000496

 6.94054413 0.10000499

 7.28757143 0.10000502

 7.63459873 0.10000505

 7.98162556 0.10000508

 8.32865334 0.10000493

 8.67568016 0.10000496

 9.02270699 0.10000498

 9.36973476 0.10000501

 9.71676159 0.10000504

 10.06378937 0.10000506

 10.41081619 0.10000509

 10.75784302 0.10000512

 11.10487080 0.10000514

 11.45189762 0.10000517

 11.79892540 0.10000485

 12.14595222 0.10000488

 12.49297905 0.10000491

 12.84000683 0.10000493

 13.18703365 0.10000496

 13.53406143 0.10000499

 13.88143539 0.10000001

Appendix 2 61
Following are the abridged contents of the UZF2test2.uzf4 output file for test simulation 2.

UNSATURATED MASS BALANCE COMPONENTS FOR ENTIRE MODEL

 APPLIED ACTUAL SURFACE ET UZ STORAGE

 TIME INFILT. RUNOFF INFILT. LEAKAGE UZ GW CHANGE RECHARGE

2.62E+06 8.47E+00 0.0 8.47E+00 0.0 0.00E+00 2.04E+01 0.00E+00 8.47E+00

2.71E+06 6.78E+01 0.0 6.78E+01 0.0 1.24E+02 1.35E+01 -6.64E+01 9.70E+00

2.80E+06 6.78E+01 0.0 6.78E+01 0.0 1.24E+02 1.35E+01 -6.63E+01 9.61E+00

...

5.25E+06 6.78E+01 0.0 6.78E+01 0.0 1.24E+02 1.37E+01 -6.55E+01 9.04E+00

5.33E+06 8.47E+01 0.0 8.47E+01 0.0 1.24E+02 1.38E+01 -5.04E+01 1.09E+01

5.42E+06 8.47E+01 0.0 8.47E+01 0.0 1.24E+02 1.38E+01 -5.04E+01 1.09E+01

...

7.88E+06 8.47E+01 0.0 8.47E+01 0.0 1.08E+02 1.46E+01 -3.97E+01 1.59E+01

7.96E+06 1.69E+02 0.0 1.69E+02 0.0 1.04E+02 1.47E+01 4.93E+01 1.58E+01

8.05E+06 1.69E+02 0.0 1.69E+02 0.0 1.04E+02 1.47E+01 4.93E+01 1.58E+01

...

1.05E+07 1.69E+02 0.0 1.69E+02 0.0 9.15E+01 1.54E+01 5.22E+01 2.57E+01

1.05E+07 4.23E+01 0.0 4.23E+01 0.0 9.25E+01 1.55E+01 -7.44E+01 2.43E+01

1.06E+07 4.23E+01 0.0 4.23E+01 0.0 9.30E+01 1.55E+01 -7.47E+01 2.40E+01

...

1.31E+07 4.23E+01 0.0 4.23E+01 0.0 6.10E+01 1.84E+01 -2.66E+01 7.97E+00

1.32E+07 8.47E+01 0.0 8.47E+01 0.0 6.18E+01 1.83E+01 1.49E+01 8.04E+00

1.33E+07 8.47E+01 0.0 8.47E+01 0.0 6.21E+01 1.83E+01 1.44E+01 8.13E+00

...

1.57E+07 8.47E+01 0.0 8.47E+01 0.0 6.76E+01 1.81E+01 -1.65E+00 1.88E+01

1.58E+07 2.54E+01 0.0 2.54E+01 0.0 6.78E+01 1.81E+01 -6.00E+01 1.76E+01

1.59E+07 2.54E+01 0.0 2.54E+01 0.0 6.79E+01 1.81E+01 -6.01E+01 1.75E+01

...

1.83E+07 2.54E+01 0.0 2.54E+01 0.0 6.91E+01 1.80E+01 -5.28E+01 9.20E+00

1.84E+07 4.23E+01 0.0 4.23E+01 0.0 6.89E+01 1.80E+01 -3.53E+01 8.83E+00

1.85E+07 4.23E+01 0.0 4.23E+01 0.0 6.89E+01 1.80E+01 -3.57E+01 9.19E+00

...

2.10E+07 4.23E+01 0.0 4.23E+01 0.0 6.67E+01 1.80E+01 -3.53E+01 1.09E+01

2.11E+07 2.54E+01 0.0 2.54E+01 0.0 6.67E+01 1.80E+01 -5.07E+01 9.37E+00

2.11E+07 2.54E+01 0.0 2.54E+01 0.0 6.67E+01 1.80E+01 -5.06E+01 9.31E+00

...

2.36E+07 2.54E+01 0.0 2.54E+01 0.0 5.77E+01 1.82E+01 -4.34E+01 1.11E+01

2.37E+07 1.69E+01 0.0 1.69E+01 0.0 5.76E+01 1.82E+01 -5.36E+01 1.29E+01

2.38E+07 1.69E+01 0.0 1.69E+01 0.0 5.75E+01 1.82E+01 -5.34E+01 1.28E+01

...

2.62E+07 1.69E+01 0.0 1.69E+01 0.0 5.20E+01 1.86E+01 -4.71E+01 1.21E+01

2.63E+07 8.47E+00 0.0 8.47E+00 0.0 5.07E+01 1.86E+01 -5.21E+01 9.86E+00

2.64E+07 8.47E+00 0.0 8.47E+00 0.0 5.04E+01 1.87E+01 -5.17E+01 9.77E+00

...

3.09E+07 8.47E+00 0.0 8.47E+00 0.0 3.03E+01 1.95E+01 -3.09E+01 9.03E+00

3.12E+07 8.47E+00 0.0 8.47E+00 0.0 2.99E+01 1.96E+01 -3.15E+01 9.82E+00

3.15E+07 8.47E+00 0.0 8.47E+00 0.0 2.95E+01 1.96E+01 -3.24E+01 1.12E+01

62 Documentation of the Unsaturated-Zone Flow Package
Manuscript approved for publication, January 6, 2006

Prepared by the U.S. Geological Survey Publishing Staff, 

Publishing Service Center, Tacoma, WA 

Linda Rogers

Sharon Wahlstrom

Bobbie Jo Richey

Publishing Service Center, Salt Lake City, UT

Angelia M. Thacker

For more information concerning the research in this report, contact the 

Director, Nevada Water Science Center 

U.S. Geological Survey

2730 North Deer Run Road 

Carson City, NV 89701

http://nv.water.usgs.gov
TM 6-A19 Niswonger and others

Printed on recycled paper

Documentation of the Unsaturated-Zone (UZF1) Package for Modeling Unsaturated Flow 

Between the Land Surface and the Water Table with MODFLOW-2005