r

NEVADA

\

CALIFORNIA \

l

NEW MEXICO \

U. S. GEOLOGICAL SURVEY PROFESSIONAL PAPER 1406-D
AVAILABILITY OF BOOKS AND MAPS OF THE U.S. GEOLOGICAL SURVEY

Instructions on ordering publications of the U.S. Geological Survey, along with prices of the last offerings, are given in the currentyear issues of the monthly catalog "New Publications of the U.S. Geological Survey." Prices of available U.S. Geological Survey publications released prior to the current year are listed in the most recent annual "Price and Availability List." Publications that may be listed in 

various U.S. Geological Survey catalogs (see back inside cover) but not listed in the most recent annual "Price and Availability List" may 

be no longer available.

Order U.S. Geological Survey publications by mail or over the counter from the offices given below.

BY MAIL OVER THE COUNTER

Books

Professional Papers, Bulletins. Water-Supply Papers, Techniques of Water-Resources Investigations, Circulars, publications 

of general interest (such as leaflets, pamphlets, booklets), single 

copies of Earthquakes & Volcanoes, Preliminary Determination of 

Epicenters, and some miscellaneous reports, including some of the 

foregoing scries that have gone out of print at the Superintendent 

of Documents, are obtainable by mail from

U.S. Geological Survey, Information Services 

Box 25286, Federal Center, Denver, CO 80225

Subscriptions to periodicals (Earthquakes & Volcanoes and 

Preliminary Determination of Epicenters) can be obtained ONLY 

from the

Superintendent of Documents

Government Printing Office

Washington, DC 20402

(Check or money order must be payable to Superintendent of 

Documents.;

Maps

For maps, address mail orders to

U.S. Geological Survey, Information Services 

Box 25286, Federal Center, Denver, CO 80225

Books and Maps

Books and maps of the U.S. Geological Survey are available 

over the counter at the following U.S. Geological Survey Earth 

Science Information Centers (ESIC), all of which are authorized 

agents of the Superintendent of Documents:

 ANCHORAGE, Alaska Rm. 101, 4230 University Dr.

 LAKEWOOD, Colorado- -Federal Center, Bldg. 810

 MENLO PARK, California Bldg. 3. Rm. 3128. 345 

Middlcfield Rd.

 RESTON, Virginia USGS National Center, Rm. 1C402, 

12201 Sunrise Valley Dr.

 SALT LAKE CITY, Utah Federal Bldg., Rm. 8105, 125 

South State St.

 SPOKANE, Washington U.S. Post Office Bldg.. Rm. 135. 

West 904 Riverside Ave.

 WASHINGTON, D.C. Main Interior Bldg., Rm. 2650, 18th 

and C Sts., NW.

Maps Only

Maps may be purchased over the counter at the following 

U.S. Geological Survey offices:

 ROLLA, Missouri 1400 Independence Rd

 STENNIS SPACE CENTER, Mississippi Bldg. 3101
Simulation of Ground-Water Flow

in Alluvial Basins in South-Central Arizona

and Parts of Adjacent States

flyT.W. ANDERSON and GEOFFREY W. FREETHEY

REGIONAL AQUIFER-SYSTEM ANALYSISSOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

U.S. GEOLOGICAL SURVEY PROFESSIONAL PAPER 1406-D

UNITED STATES GOVERNMENT PRINTING OFFICE, WASHINGTON
U.S. DEPARTMENT OF THE INTERIOR 

BRUCE BABBITT, Secretary

U.S. GEOLOGICAL SURVEY 

GORDON P. EATON, Director

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.

Library of Congress Cataloging in Publication Data

Anderson, T.W. (Thomas Warren), 1940-

Simulation of ground-water flow in alluvial basins in south-central Arizona and parts of adjacent States / T.W. Anderson and

Geoffrey W. Freethey. 

p. cm. (Regional aquifer-system analysis southwest Alluvial basins, Arizona and adjacent States) (U.S. Geological Survey

professional paper; 1406-D) 

Includes bibliographical references. 

Supt. of Docs, no.: 119.16:1406-D 

1. Groundwater flow Arizona Computer simulation. I. Freethey, Geoffrey W II. Title. III. Series. IV. Series: U.S. Geological

Survey professional paper; 1406-D. 

GB1197.7.A54 1995

551.49'09791 dc20 94-31125

CIP

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

Box 25286, Federal Center, Denver, CO 80225-0425
FOREWORD

THE REGIONAL AQUIFER-SYSTEM ANALYSIS PROGRAM

The Regional Aquifer-System Analysis (RASA) program was started in 1978 

after a congressional mandate to develop quantitative appraisals of the major 

ground-water systems of the United States. The RASA program represents a systematic effort to study a number of the Nation's most important aquifer systems 

which, in aggregate, underlie much of the country and which represent important 

components of the Nation's total water supply. In general, the boundaries of these 

studies are identified by the hydrologic extent of each system, and accordingly, 

transcend the political subdivisions to which investigations have often arbitrarily 

been limited in the past. The broad objective for each study is to assemble geologic, hydrologic, and geochemical information, to analyze and develop an understanding of the system, and to develop predictive capabilities that will contribute 

to the effective management of the system. The use of computer simulation is an 

important element of the RASA studies, both to develop an understanding of the 

natural, undisturbed hydrologic system, and of any changes brought about by 

human activities, as well as to provide a means of predicting the regional effects of 

future pumping or other stresses.

The final interpretive results of the RASA program are presented in a series of 

U.S. Geological Survey Professional Papers that describe the geology, hydrology, 

and geochemistry of each regional aquifer system. Each study within the RASA 

program is assigned a single Professional Paper number and where the volume of 

interpretive material warrants, separate topical chapters that consider the principal elements of the investigation may be published. The series of RASA interpretive reports begins with Professional Paper 1400 and thereafter will continue in 

numerical sequence as the interpretive products of subsequent studies become 

available.

Gordon P. Eaton 

Director

III

CONTENTS

Abstract...................................................................................................... Dl

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

Purpose and scope......................................................................... 2

Approach......................................................................................... 2

Description of the study area........................................................ 3

Geohydrologic setting.............................................................................. 3

Description of sedimentary units................................................. 5

Predevelopment hydrologic conditions....................................... 6

Development of water resources ................................................. 7

Basin categories........................................................................................ 7

Ground-water modeling........................................................................... 9

Modeling approach......................................................................... 11

Selected basin models.................................................................... 12

Upper San Pedro basin........................................................ 12

Benson basin......................................................................... 16

Lower San Pedro basin........................................................ 20

Willcox basin......................................................................... 24

Douglas basin........................................................................ 25

Tucson basin.......................................................................... 29

Avra Valley............................................................................. 35

Waterman Wash basin.......................................................... 38

Ground-water modeling Continued

Selected basin models Continued

Harquahala Plain................................................................... D39

McMullen Valley.................................................................... 44

Parker Valley......................................................................... 45

Mohave basin........................................................................ 48

Basin-group models........................................................................ 49

Southeast basins................................................................... 49

Central and west basins....................................................... 54

Results from specific-basin and basin-group sensitivity analysis...... 57

Relation of basin category and model input............................... 57

Information adequacy and transferability................................... 64

Hydrologic boundaries........................................................ 65

Flow components................................................................. 66

Hydrologic factors....... ......................................................... 67

Hydrologic insights gained from modeling southwest alluvial

basins.......................................................................................... 70

General modeling considerations........................................................... 73

Summary.................................................................................................... 74

Selected references.................................................................................. 76

ILLUSTRATIONS

Page

FIGURE 1. Map showing Southwest Alluvial Basins study area........................................................................................................................................... D4

2. Sketch showing generalized vertical sequence of sediments that fill most basins......................................................................................... 5

3. Sketch showing typical predevelopment inflow and outflow mechanisms in basins of the study area...................................................... 6

4. Graph showing annual ground-water pumpage, 1915-80................................................................................................................................... 7

5. Map showing areas of water-level decline from the beginning of development through about 1980 .......................................................... 8

6. Graphs showing water-level declines in developed basins ............................................................................................................................... 9

7. Map showing categories of basins and basins for which ground-water models were developed................................................................ 10

8, 9. Maps and graphs showing location of basin, estimated and simulated water-budget components and pumpage, ground-water 

levels, model grid and boundary, and observed and simulated water levels for 

8. Upper San Pedro basrn.................................................,.................................^ 14

9. Benson basin.................................................................................................................................................................................................... 18

10. Maps and graphs showing location of basin, estimated and simulated water-budget components and pumpage, ground-water

levels, model grid and boundary, simulated changes in water-budget components, and observed and simulated water levels, 

lower San Pedro basin ................................................................................ 22

11,12. Maps and graphs showing location of basin, estimated and simulated water-budget components and pumpage, ground-water 

levels, model grid and boundary, and observed and simulated water levels for 

11. Wfflcox basin.......................................................................... 26

12. Douglas basin........................................................ 30

13. Maps and graphs showing location of basin, estimated and simulated pumpage, ground-water levels, model grid and boundary,

and observed and simulated water-level declines, Tucson basin.............................................................................................................. 32

14. Maps and graphs showing location of basin, estimated and simulated pumpage, measured and simulated ground-water levels,

model grid and boundary, and observed and simulated water-level declines, Avra Valley.................................................................... 36

15-17. Maps and graphs showing location of basin, estimated and simulated pumpage, ground-water levels, model grid and boundary, 

and observed and simulated water levels for 

15. Waterman Wash basin..................................................................................................................................................................................... 40

16. Harquahala Plain......................................................................................................................................................... 42

17. McMullenVaUey.................................................... 46
VI CONTENTS

Page

FIGURES 18,19. Maps and graphs showing location of basin, estimated and simulated water-budget components, model grid and boundary, 

and observed and simulated water levels for 

18. Parker Valley.............................................................................................................................................. D50

19. Mohave basin............................................................................................................................................................................................... 52

20-25. Sketches showing20. Typical geographic and geohydrologic features in the southeast basins............................................................................................. 55

21. Model grid, boundary conditions, and data input for steady-state basin-group model of the southeast basins............................ 56

22. Effects of modifications of mountain-front recharge, interlayer-leakance values, lower-layer transmissivity, and basin

width on water-level profiles of the steady-state basin-group model of the southeast basins................................................. 58

23. Typical geographic and geohydrologic features in the central and west basins................................................................................. 60

24. Model grid, boundary conditions, and data input for steady-state basin-group model of the central and west basins................ 61

25. Variations in data used and results of analysis of basin shape and recharge location and distribution in basin-group model

of the central and west basins.......................................................................................................................................................... 62

26-32. Graphs showing26. Effect of lower-layer boundary location on simulated drawdown in McMullen Valley..................................................................... 66

27. Variation in components of loss from river flow, evapotranspiration, and drain flow resulting from changes in the leakance

factor for drains and rivers in the Parker Valley model................................................................................................................. 67

28. Variation in net loss from river flow in Mohave basin and Parker Valley related to change in average evapotranspiration

rate ....................................................................................................................................................................^ 67

29. Sensitivity of model-calculated heads in the Benson and Douglas basin models to variations in interlayer-leakance value...... 69

30. Maximum simulated drawdown using a variable area! distribution of specific yield and a uniform specific yield

representing the average value for the models of upper San Pedro basin and Douglas basin................................................. 69

31. Response of the basin-group model of the central and west basins to variations in areal distributions and values of

specific yield ............................................................................................................................................................................^ 71

32. Water-level declines resulting from variable areal distributions of pumping from the lower layer and presence or absence

of a fine-grained facies in the basin-group model of the central and west basins..................................................................... 72

TABLES

Page

TABLE 1. Range in values of geohydrologic properties in southeast basins and central and west basins and values used in basin-group

models............................................................................................................................^ D54

2. Relative importance of hydrologic boundaries, flow components, and geohydrologic factors by basin group for modeling steadystate and transient conditions.....................................................................................................................................^ 64

3. Geohydrologic factors evaluated as part of the modeling process ..................................................................................................................... 65

4. Sensitivity of model-calculated heads and flow quantities to changes in values of flow components and hydrologic factors in

southeast basin model......................................................................................................................................................................................... 68
CONTENTS VII

CONVERSION FACTORS AND VERTICAL DATUM

Multiply inch-pound unit

inch (in.)

foot (ft)

mile (mi)

acre

square mile (mi2)

gallon per minute (gal/min)

foot per second (ft/s)

foot squared per day (ft2/s)

foot per mile (ft/mi)

cubic foot per second (ft^/s)

cubic foot per second per square mile

[(ft3/s)/mi2]

acre-foot (acre-ft)

By

25.4

0.3048

1.609

0.004047

2.590

0.06309

0.3048

0.0929

0.1894

0.02832

0.01093

0.001233

To obtain metric unit

millimeter

meter

kilometer

square kilometer

square kilometer

liter per second

meter per second

meter squared per day

meter per kilometer

cubic meter per second

cubic meter per second

per square kilometer

cubic hectometer

Sea Level: In this report "sea level" refers to the National Geodetic Vertical Datum of 1929 

a geodetic datum derived from a general adjustment of the first-order level nets of both the 

United States and Canada, formerly called Sea Level Datum of 1929.

REGIONAL AQUIFER-SYSTEM ANALYSIS-SOUTHWEST ALLUVIAL BASINS, 

ARIZONA AND ADJACENT STATES

SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL 

BASINS IN SOUTH-CENTRAL ARIZONA AND PARTS OF

ADJACENT STATES

By T.W. ANDERSON and GEOFFREY W. FREETHEY

ABSTRACT

Numerical modeling was used to examine the characteristics of 

ground-water flow systems in the alluvial basins of south-central Arizona and parts of adjacent States and to evaluate similarities and differences among basins. The 72 alluvial basins of the study area were 

categorized into five groups that have similar geohydrologic characteristics and similar requirements for model development. Basin categorization is based on the relative importance of various components of the 

geologic and hydrologic systems. The principal geologic differences 

among groups are the presence and extent of fine-grained sediments in 

the basin center and the presence of saturated stream alluvium. These 

units represent the least and most productive parts of the aquifer, 

respectively. The principal hydrologic differences among groups are 

magnitude of the component of mountain-front recharge and presence 

of a perennial stream in the basin. These groups were termed southeast, 

central, west, Colorado River, and highland.

Ground-water flow models of 12 specific basins were developed to 

evaluate the utility of transferring hydrologic concepts and property values within and among groups of basins. In addition, two models representing generalized basin systems were developed that incorporated 

average physical and hydrologic characteristics of the southeast basin 

group in one model and of the central and west basin groups in the 

other. These basin-group models were used, along with the 

specific-basin models, to examine the degree of model sensitivity to 

changes in basin shape and the values and areal distributions of hydrologic properties.

Each model of a specific basin was based on available field data and 

previous conceptual or numerical models. Quantitative determinations 

of model-data requirements were used where field data were available. 

Transferability of information among basins of a particular group was 

tested and generally was successful. The properties that were most 

important to the successful development of basin models were summarized in groups including hydrologic boundaries, flow components, and 

hydrologic properties.

Systematic model-sensitivity analyses were done on all models to 

qualitatively evaluate how uncertainties in hydrologic property values 

and distributions might affect results of the simulations. Sensitivity analyses indicated that simulations of basins in the same group are affected 

most by the same changes in model input. The degree of sensitivity of a 

model to a specific property depends on the geologic setting of a basin,

the magnitudes of predevelopment downvalley flow, and the degree of 

development that has taken place.

Changes in hydraulic conductivity, natural recharge quantities, and 

extent and integrity of confining layers affected computed water levels 

mainly in models of southeast basins. Changes in specific yield and location and quantity of pumping affected computed water levels the most in 

models of central and west basins in which storage depletion is common. Changes in quantity of water moving between streams and the 

aquifer, evapotranspiration, and the return quantity of excess applied 

irrigation water affected computed water levels the most in models of 

the Colorado River and probably of the highland basins in which the 

hydrologic system is dominated by a perennial stream.

INTRODUCTION

The Southwest Alluvial Basins, Regional AquiferSystems Analysis (Swab/RASA), is one of a series of 

projects designed to systematically evaluate the Nation's 

major ground-water systems (Bennett, 1979). The Southwest Alluvial Basins study began in October 1978 and was 

one of 29 systems included in the RASA Program (see 

Foreword).

The Swab/RASA study area (fig. 1) includes about 

82,000 mi2 and consists of 72 basins, most of which are 

virtually independent hydrologic systems. The physical 

character of the basins is diverse, although similarities 

exist in geohydrologic character. Tectonic activity, which 

formed the Basin and Range physiographic province 

(Fenneman, 1931), was areally variable. Rocks that bound 

the individual basins differ areally and are the source of 

detrital material that subsequently filled the basins. 

Deposition of the alluvial material that constitutes the 

aquifer(s) probably was unique in each basin.

The basins of the study area are considered a regional 

system for several reasons:

Dl
D2 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

 Several basins are linked by natural hydraulic connection in which outflow from one basin represents inflow 

to a downgradient basin.

 The hydrologic response to development of the water 

resources and water-use characteristics generally have 

been similar.

 Ground-water resources of the basins represent a common element in the economy of the entire region. 

This study integrated and explored some of the concepts 

used to evaluate individual ground-water systems and to 

compare the individual basin characteristics and hydrologic properties.

The overall purpose of the Swab/RASA study was to 

describe and define the hydrologic systems in the alluvial 

basins. The study included quantitative estimates and 

qualitative descriptions of the flow systems before and 

after development. Numerical models were developed to 

simulate ground-water flow in specific individual basins 

and in basin systems that represent individual groups of 

basins. The study emphasized ground-water flow simulation and was designed to 

1. Simulate ground-water flow in a few representative 

basins of the 72 basins in the study area.

2. Analyze findings related to geologic and hydrologic 

patterns that exist throughout the study area.

3. Document similarities and differences in basin and 

regional geohydrology to support the use of information transfer from basins that were modeled to those 

that were not modeled.

PURPOSE AND SCOPE

The findings of the Swab/RASA study are described in 

U.S. Geological Survey Professional Papers 1406-A 

through 1406-D, including this report, Professional Paper 

1406-D. This report is closely related to Professional 

Paper 1406-B (Anderson and others, 1992), which 

describes the geohydrology and water resources of the 

Swab/RASA study area.

This report describes the modeling approach used in 

the study, the design and calibration of the ground-water 

models, and the results from specific-basin and 

basin-group models. The models were developedby using 

geologic and hydrologic knowledge of the systems as presented by Anderson and others (1992) and in other previous studies. The emphasis of the report is on 

documentation of the knowledge gained from the study of 

geologic and hydrologic similarities and differences 

among the basins of the study area. The models provide a 

means of testing the general understanding of the 

hydrologic-system function of the basins. The report also 

includes a general discussion of the modeling techniques 

used to simulate the response of a system or part of a system to a specific stress.

The transferability of information and the limitations 

of interpretations are discussed, and all facts and assumptions on which the interpretations depend are presented. 

Possible use of the models for predictive and management tools is also discussed, although the models developed as part of this study generally were not used to 

predict possible future conditions.

APPROACH

The general approach to this study was to analyze 

available data and develop numerical models of selected 

basins for which data were adequate for history matching. 

The study rationale assumes that some system functions 

are common and that selected information is transferable 

among basins or groups of basins. Information transfer 

was used as a means of quantitatively defining properties 

of the hydrologic systems and was evaluated through the 

use of numerical models. The importance of various 

hydrologic properties and flow components that define 

the ground-water systems also was evaluated for predevelopment and development conditions. The informationtransfer approach was used to reduce the need for extensive data collection and to reduce disparities in hydrologic-system definition that result from different political 

jurisdictions and economies.

For studying individual basins as hydrologic units, the 

study area was divided into 72 basins. The boundaries 

were selected to represent surface-water divides, zones of 

minimal intefbasin connection, or in some instances, arbitrary boundaries through areas of minimum ground-water 

development. As defined, the basins serve as the smallest 

practical unit for analysis and numerical modeling. Each 

basin can be modeled individually with minimal boundary 

effects because contiguous areas were not included; 

basins can be readily combined to model or analyze larger 

areas.

The 72 basins were categorized on the basis of physiographic, geologic, and hydrologic properties and water 

chemistry to develop and evaluate information transfer 

that could be used to define and describe both the individual basin systems and the regional system. The effects of 

poorly known system properties were explored by using a 

range of values in the model and evaluating the effect of 

each range of values on model performance. The results 

indicated basin conditions for which detailed and accurate values of hydrologic properties are needed and, conversely, conditions for which general approximations of 

these properties are acceptable.

The basins of the study area were grouped initially into 

categories on the basis of discernible surface features 

climate, runoff, geology, altitude, shape, and general 

physiography that were thought to have an influence on 

the ground-water hydrology of the basin (Anderson and
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D3

others, 1992). A typical basin that had sufficient data 

available for model construction was then selected from 

each group, and a numerical model of the ground-water 

flow system was developed. Subsequently, the information gained from this model was transferred to one pr several other basins for which fewer data were available. 

Models were then developed for those basins. If the information remained acceptable for use in the model when 

transferred to*the basins with fewer data, then greater 

confidence could be placed in its applicability to basins 

that are less developed or undeveloped and for which virtually no data are available. Tests of the adequacy of data 

to define hydrologic-system properties and sensitivity of 

model response to variations in property values were 

made on all models to relate those results to the types and 

magnitude of data-collection programs. Two basin-group 

models were developed by using average values of physical and hydrologic properties of the southeast basins and 

of the central and west basins. These basin-group models 

were used to test sensitivity of hydrologic response to 

variations and generalizations of properties used in the 

models.

Data used in the numerical models consisted of information that is available for most developed basins in the 

Swab/RASA study area. Recharge rates were estimated 

from average annual precipitation. Values and areal distributions of transmissivity, hydraulic conductivity, and specific yield were estimated from aquifer-test results, 

specific-capacity values, or examination of descriptions 

of samples in drillers' logs. Flow-net analyses were used 

to estimate distributions of hydraulic conductivity, transmissivity, and mountain-front recharge for initial use in 

the models. Alternative values and areal distributions of 

the hydrologic properties were tested in the ground-water 

flow models through sensitivity analyses.

DESCRIPTION OF THE STUDY AREA

The study area encompasses about 82,000 mi2 in southern and central Arizona, western New Mexico, and small 

parts of southeastern California and Nevada (fig. 1). The 

study area includes the surface-water drainage of the Colorado River between Hoover Dam and the international 

boundary near Yuma. The San Bernardino and Douglas 

basins in southeastern Arizona and parts of several basins 

in southwestern Arizona that are cut by the international 

boundary also are included. These basins drain into Mexico but are included in the study because of similarities in 

physical, climatic, and land-use characteristics to those of 

basins within the Colorado River drainage.

The area includes approximately equal parts of highly 

productive alluvial aquifers in the basins and nearly 

impermeable consolidated rocks in the surrounding 

mountains that function as boundaries to the aquifer

systems. The aquifers are composed mostly of alluvial 

deposits and occur in the basins formed by down-faulted 

blocks between mountain masses, which represent the 

up-faulted blocks. An estimated 900 million acre-ft of 

recoverable ground water is stored in the upper 1,200 ft of 

alluvial deposits that fill the basins (Freethey and Anderson, 1986). The volume of ground water in storage, occurrence and rate of movement of ground water, and 

productivity of wells have a wide range from basin to 

basin.

Physiography, climate, and geohydrology vary widely 

throughout the study area. Land-surface altitude ranges 

from about 150 ft above sea level near Yuma, Arizona, to 

more than 10,000 ft in the mountainous regions of eastern 

Arizona and western New Mexico. Average annual precipitation ranges from less than 3 in. to more than 30 in. and 

generally decreases from the higher eastern basins to the 

lower western basins. Many streams in the northeastern 

part of the study area are perennial and typically recharge 

the alluvial aquifers over which they flow. Only the Colorado River presently contains sufficient flow to reach the 

dry southwesternmost part of the area. The Gila River 

was perennial through its entire reach in the early 1900's 

but has become ephemeral in places because of upstream 

storage and diversion of surface water and increased 

losses of streamflow resulting from ground-water 

development.

GEOHYDROLOGIC SETTING

The basins of the study area were formed as a result of 

the Basin and Range structural disturbance, which 

occurred 15 to 10 million years ago (Scarborough and 

Peirce, 1978). Movement along high-angle normal faults 

near the present basin edges resulted in the formation of a 

series of generally north- to northwest-trending basins 

and mountain ranges. Basin subsidence occurred at varying rates throughout the study area and was accompanied 

by deposition of locally derived sediment within the internally drained basins. Deposition continued after the 

establishment of external drainage in most basins. The 

mountains consist of igneous, metamorphic, and sedimentary rocks. The bedrock of the mountains is virtually 

impermeable and forms the side and bottom boundaries 

of the alluvial aquifers.

The basins are filled with alluvial deposits that range 

from less than 2,000 ft to more than 10,000 ft in thickness 

(Oppenheimer and Sumner, 1981). These deposits form 

the major aquifers of the area and store large quantities of 

ground water. Basin subsidence and sediment deposition 

occurred at different rates throughout the area, and as a 

result, the thickness, areal extent, and grain size of the 

alluvial deposits are highly variable. Sediments range in
D4 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

115° 114° 113 110° 109° 108°

50 100 MILES

UNITED STATES

BASIN AND RANGE

PHYSIOGRAPHIC

PROVINCE

50 100 KILOMETERS

STUDY AREA

INDEX MAP 

FIGURE 1. Southwest Alluvial Basins study area

EXPLANATION

Alluvial Deposits Locally may 

include evaporite deposits and 

volcanic rocks

Bedrock of the mountains 

Basin boundary 

Boundary of study area
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D5

POORLY TO WELL-DEFINED TRANSITION 

FROM LOWER TO UPPER BASIN FILL

Greater than 60 percent fines 

" n " percent fines

Greater than 

n percent fines

Not to scale

EXPLANATION

Stream alluvium

Basin fill Shows a facies change 

from coarse grained to fine grained

Pre-basin and range sedimentary and 

volcanic rocks

Bedrock of the mountains

llf* Normal fault Arrows show 

I relative vertical movement

FIGURE 2. Generalized vertical sequence of sediments that fill most basins.

grain size from clay to gravel and can be unconsolidated 

to highly consolidated in a single basin.

Basins are interconnected in a dendritic pattern, much 

like the surface drainage, and form an integrated regionalflow system. Individual basin-aquifer systems serve 

mainly as ground-water reservoirs with only small areas 

of interconnection and a small quantity of flow from 

basins of higher altitude to those of lower altitude.

DESCRIPTION OF SEDIMENTARY UNITS

A similar vertical sequence of sedimentary units is 

present in most basins and consists of pre-Basin and 

Range rocks, lower and upper basin fill, and stream alluvium (fig. 2). Each unit represents a different depositional 

environment, which is reflected in the lithology. The units 

are distinguishable on the basis of structural relation, 

degree of consolidation and deformation, source area of 

clasts, geologic age, and water-bearing characteristics. 

Although each basin is unique, general spatial patterns of

textural composition, thickness, and extent of the units 

have been documented among groups of basins (Anderson, 1986; Freethey and others, 1986).

Pre-Basin and Range rocks consist of moderately to 

highly consolidated continental sedimentary deposits that 

range in grain size from silt to gravel. The rocks were 

deposited before the Basin and Range structural disturbance and were faulted and tilted at the same time as the 

underlying bedrock. Basin fill was deposited in the basins 

formed by the Basin and Range disturbance. Lower basin 

fill is moderately to highly consolidated, includes extensive fine-grained material near the center of basins, and 

contains evaporites that range from disseminated gypsum 

to massive halite deposits. The unit locally contains interbedded volcanic material consisting of basalt flows and 

tuff. Lower basin fill was deposited in topographically 

closed basins, and upper basin fill was deposited during a 

transition period from a closed to an integrated drainage. 

Upper basin fill is less consolidated, generally is thinner,
D6 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

UNDERFLOW IN

INFLOW OUTFLOW

^ '^

STREAM FLOW

2

UNDERFLOW OUT

FIGURE 3. Typical predevelopment inflow and outflow mechanisms in basins of the study area.

and contains fewer fine-grained sediments than the lower 

basin fill. Stream alluvium was deposited contemporaneously with the establishment of the present-day 

surface-drainage system. Stream alluvium ranges in texture from clay to boulders and is unconsolidated. Sediments that fill the basins generally are interconnected 

hydraulically and form a single aquifer system within each 

basin.

Upper and lower basin fill exhibit gradational patterns 

in grain-size distribution. The deposits become coarser 

grained in the direction of the source areas, which are 

principally the mountains that bound the basins. The 

change in percentage of fine-grained material from the 

source area to the basin center can be highly variable 

from basin to basin within equivalent units (Anderson and 

others, 1992). The sediment is continental in origin except 

the marine Bouse Formation in the Colorado River region 

(Metzger, 1968).

PREDEVELOPMENT HYDROLOGIC CONDITIONS

Ground water generally occurs under unconfined conditions in the aquifers of the study area. In several basins 

where extensive fine-grained deposits overlie the principal water-bearing unit, ground water occurs under confined conditions. Hydraulic heads in the basin fill range 

from at or slightly above land surface near perennial 

streams to more than 600 ft below land surface near the 

mountain fronts in several basins. Within the study area, 

an estimated 900 million acre-ft of ground water was 

stored in the upper 1,200 ft of alluvial deposits before 

development (Freethey and Anderson, 1986). The quantity of water in storage is hundreds to thousands of times 

greater than the quantity of water that moves into and out 

of the aquifers annually.

Ground-water inflow to the aquifers occurs mainly 

through three mechanisms: (1) infiltration of streamflow 

along the major streams, (2) infiltration along the mountain fronts, and (3) underflow from adjacent basins 

(fig. 3). Direct infiltration of precipitation that falls on the 

basin floor is negligible. The total average annual 

recharge to all the basins in the study area is estimated to 

be 2.5 million acre-ft for predevelopment conditions 

(Anderson and others, 1992). Recharge varies from year 

to year and from place to place within a basin. Annual 

recharge probably ranges from near zero in dry years to 

several times the long-term annual average in abnormally 

wet years. In general, water enters the aquifer at the 

upstream end and along the mountain fronts and flows 

toward the basin center and the outlet at the downstream 

end. Specific flow characteristics in a basin or for a group 

of basins depend on the hydrologic properties of the aquifer and on the total annual downvalley flow through the 

aquifer.

Ground-water outflow from the aquifers occurs 

through three basic mechanisms: (1) discharge to springs 

and streams, (2) evaporation and transpiration from the 

water-table zone, and (3) underflow to adjacent basins 

(fig. 3). The quantity of natural discharge from each basin 

probably was nearly constant from year to year before 

development. The influence of wet and dry years was 

dampened by the large quantity of ground water stored in 

the basin. Average inflow to and outflow from the aquifers 

before development were assumed to be equal at about 

2.5 million acre-ft/yr, and no long-term change in storage 

occurred.

Predevelopment water-level data indicate the general 

direction of ground-water flow through the basins 

(Freethey and Anderson, 1986). The data indicate that the 

basins are linked to various degrees and form an integrated regional-flow system. The shapes of the water-level
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D7

FIGURE 4. Annual ground-water pumpage, 1915-80.

contours generally indicate areas of inflow and outflow 

(Anderson and others, 1992). Basins that receive a large 

quantity of mountain-front recharge have contours that 

are U-shaped and nearly parallel to the mountain fronts 

(Anderson and others, 1992, fig. 15A). Basins that receive 

little or no mountain-front recharge have contours that 

generally are straight and nearly perpendicular to the 

mountain fronts (Anderson and others, 1992, fig. 15E). 

Both extremes of contour shape, as well as intermediate 

shapes, exist within the study area and indicate the 

extreme areal variability in recharge magnitude. Significant variations in water-level gradients also occur in the 

basins and generally indicate variable physical and hydrologic properties of the aquifer material.

DEVELOPMENT OF WATER RESOURCES

Development of water resources in the study area has 

been principally for agriculture. Use of the water 

resources has resulted in depletion of base flow of formerly perennial streams in several basins and overdraft of 

ground water in most developed basins. Ground-water 

pumping began in the late 1800's to augment the highly 

variable surface-water supplies. Development began at 

different times in various basins. As early as 1920, development was mostly limited to east and west Salt River 

Valley, Eloy basin, and Stanfield basin (fig. 1). In 1940, 

nearly 1 million acre-ft of ground water was withdrawn in 

east and west Salt River Valley; another 0.5 million acre-ft 

was withdrawn in Eloy, Stanfield, Tucson, and Yuma 

basins. The other basins were either undeveloped or 

slightly developed until after 1940. Ground-water pumpage was 1.7 million acre-ft in 1942 and increased to 3.8 

million acre-ft in 1952. During 1950-80, an average of 4.8 

million acre-ft/yr was pumped (fig. 4), most of which was 

pumped from 16 basins. From the beginning of development through 1980, about 184 million acre-ft of ground 

water was pumped from the basins (U.S. Geological Survey, 1982). This pumpage does not represent a total loss 

from storage because some water recharges the systems 

through leakage from irrigated land and from canals. 

Added recharge was induced by increased streamflow 

infiltration, and some water that previously discharged 

through evapotranspiration by natural vegetation has 

been captured. On the basis of the volume of sediments 

dewatered, about 50 percent of the total pumpage is estimated to be depletion of ground water in storage.

From predevelopment time through 1980, groundwater levels declined from less than 50 ft in areas of 

minor development to more than 200 ft in the major developed areas (fig. 5). A maximum decline of more than 

450 ft occurred in the Stanfield basin between 1940 and 

1980. Average annual rates of decline have ranged from a 

few tenths of a foot per year to more than 10 ft/yr (fig. 6). 

The magnitude of ground-water depletion varies from 

basin to basin and reflects the influence of the geohydrologic environment, as well as the magnitude and duration 

of ground-water withdrawals and stream diversions. In 

the basins along the Colorado River in which irrigated 

agriculture is extensive and surface water is a major 

source of supply, a few feet of water-level rise has 

resulted from the infiltration of excess applied irrigation 

water.

BASIN CATEGORIES

The alluvial basins of Arizona can be-grouped into five 

categories on the basis of hydrologic and geologic similarities (Anderson, 1984; Anderson and others, 1992). The 

most important hydrologic factor used in categorizing the 

basins is the total downvalley flow, which represents the
D8 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

11 O"7o 5° 114° 113° 112° UTAH 111°

j' "T 1 ARIZONA '

*y 

* ! *f

110° 109° 10J """ I

36

35

34

- i - > i _____ 

(

C

_______ L_ 

)

I 

) 50

50

I

100 MILES

I

I 

100 KILOMETERS

EXPLANATION

Water-level decline, in feet Alluvial Deposits Locally may include 

evaporite deposits and volcanic rocks

Bedrock of the mountains 

Basin boundary 

Boundary of study area

FIGURE 5. Areas of water-level decline from the beginning of development through about 1980.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D9

LU 

O

tr

Z) 

CO

tr

LU

o

100

200

-i 300

O

_j

LU 

CD

LU 400

LU

500

600

700

STANFIELD BASIN 

T. 6 S., R. 2 E., Sec. 1

o 10 01

o 

oo

O)

LU

O 100

<

LL

tr

Z) 

CO

o

_l 

LU 

CD

ILU 

LU 

LL

O

200

300

SAN SIMON BASIN 

T. 12S., R. 28 E., Sec. 17

100

CL

LU 200 

Q

300

DOUGLAS BASIN 

T. 21 S., R. 26 E., Sec. 24

o

<001

o r~- 

01

o 

oo 01

FIGURE 6. Water-level declines in developed basins.

total renewable water resource available on an annual 

basis before development. The principal geologic property used in categorizing the basins is the character of the 

sediments that make up the aquifer system (Pool, 1984). 

Sediment character controls the movement of ground 

water. These two factors downvalley flow and basinsediment characteristic control the manner in which the 

aquifer system will respond to a given stress. On the basis 

of these two factors, the basins of the study area were 

grouped into five categories southeast, central, west, 

Colorado River, and highland (fig. 7). Generally, basins of 

the west category have the smallest downvalley flow, and 

basins of the Colorado River category have the largest 

downvalley flow. The principal aquifer system in each category generally consists of the following: (1) stream alluvium along the basin axis and adjacent areas in the 

southeast basins, (2) upper basin fill and stream alluvium 

in central basins, (3) lower basin fill near the basin center 

in west basins, (4) stream alluvium in Colorado River 

basins, and (5) stream alluvium that overlies preCenozoic sedimentary rocks in highland basins. The general geographic grouping of basins of a category reflects

the interrelation of hydrology, lithology, and physiographic factors. The response of a ground-water system 

to the stresses of development can be related to the basin 

category; therefore, the response of a ground-water 

system to a proposed development is qualitatively 

predictable.

GROUND-WATER MODELING

Two types of numerical ground-water flow models 

were used to analyze the hydrologic systems within the 

basins. One type included models of selected basins 

(fig. 7). The goal in developing this type was to simulate 

the hydrologic system of specific basins and, in doing so, 

establish reasonable ranges for hydrologic properties and 

water-budget components in each modeled basin. The 

other type comprised two basin-group models that simulated basin hydrology for three different geohydrologic 

environments of the study area. The goal in developing 

these basin-group models was to analyze how flow models for these three environments were affected by 

changes in values of individual properties, providing the
DIG REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

115° 114 110° 109° 108°

PART OF 

SOUTHEAST

50 100 MILES

_l

1

2

3

4

Mohave Basin 

Parker Valley 

McMullen Valley 

Harquahala Plain 

Boundary of study area

5

6

7 

Q

 I | 

50 100 KILOMETERS

EXPLANATION

Modeled basin

Waterman Wash Basin 

Avra Valley 

Tucson Basin 

Lower San Pedro Basin 

Basin boundary

y Benson Basin

t U Upper San Pedro Basin

| "I Willcox Basin

1 ^ Douglas Basin

- - -" Boundary of general 

categories of basins

FIGURE 7. Categories of basins and basins for which ground-water models were developed.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS Dll

changes were within previously established or reasonable 

ranges. The basin-group models eliminated site-specific 

problems commonly encountered in the specific-basin 

models, and an appraisal of the relative importance of 

various properties in simulating basin hydrology of different categories could be obtained.

Twelve specific-basin models were developed, and history matched to a predevelopment or another appropriate 

steady-state condition. Ten of the basin models were evaluated further by using available transient or development 

conditions. Transient conditions were not simulated in 

the Parker Valley and Mohave basin models because of 

the rapid change to a new steady-state condition in 

response to development in the Colorado River basins. 

For these two basins, two or more steady-state conditions 

were simulated in which the magnitude of flow components and head distributions differed from one steadystate condition to another.

Basin-group models used average values of physical 

and hydrologic properties that were based on actual values in the basins of the particular group being represented. These models were used to evaluate the 

significance of particular physical and hydrologic properties and flow components on modeling results and to systematically explore model sensitivity to variations in data 

values used in the models.

MODELING APPROACH

Basins selected for simulation encompassed the ranges 

in conditions and hydrologic-property values that exist in 

the study area. All modeling was accomplished by using 

the three-dimensional finite-difference model described 

by Trescott (1975) and Trescott and Larson (1976). A 

quasi-three-dimensional approach was used in which the 

aquifer system is simulated as two layers. The two layers 

represent upper and lower segments of the aquifer that 

are connected hydraulically. To allow vertical flow 

between simulated layers, this vertical hydraulic connection is represented by a leakance coefficient, which was 

either calculated by the model or assigned to each cell. 

Data entry to the models included the geographic extent 

and thickness of the two layers of the aquifer; waterbearing properties of the layers, such as horizontal and 

vertical hydraulic conductivity; specific yield and storage 

coefficient; ground-water recharge and discharge; and 

pumpage. The upper layer of each modeled basin was 

simulated by using hydraulic-conductivity values to allow 

transmissivity to change as water levels changed because 

of stress imposed on the system. The lower layer was simulated by using transmissivity values, which did not 

change during transient simulations. The decrease in 

transmissivity caused by lowering water levels below the 

top of the lower layer was not adjusted because the

magnitude of change was not considered significant, 

although this was not tested as part of a model analysis. 

In a few simulations, ground-water withdrawals resulted 

in sufficient water-level declines that a secondary 

specific-yield value had to be used for the simulated lower 

layer rather than the initial storage-coefficient value. 

Results of the simulations consist of hydraulic head at 

specified locations, depths, and times; a mass balance 

that represents the summation of flow into and out of the 

simulated system; and simulated change in storage.

Each basin model was developed by using available 

data and known geologic and hydrologic concepts developed by previous investigators. Information transfer also 

was used in basins where the geohydrologic system was 

poorly defined. Measured and reported water levels in 

wells were used to develop potentiometric-surface configurations. Aquifer-test results yielded values for transmissivity and hydraulic conductivity, and analyses of various 

types of well logs were used to define aquifer lithology, 

saturated thickness, and aquifer extent in some areas. 

Flow-net analyses were used to estimate the quantity and 

distribution of recharge along aquifer boundaries where 

reliable water-level configurations were available.

For development of models for those basins in which 

data were sparse, values of hydrologic properties and 

mountain-front recharge were obtained by transferring 

values from geohydrologically similar basins and testing 

these values during model calibration. Empirical methods 

of estimating transmissivity from specific capacity (Theis 

and others, 1963) and of estimating transmissivity and 

specific yield from drillers' logs (Davis and others, 1959; 

Kisser and Haimson, 1981; Evans and Haimson, 1982) 

were used to supplement quantitative data where possible. For basins in which water-level data near mountain 

fronts are unavailable, mountain-front recharge was estimated by using the equation: log Qrech = -1.40 + 0.98 logP, 

where P is volume of annual precipitation in excess of 8 

in./yr. Use and limitation of the equation are explained in 

Anderson and others (1992).

Final development of steady-state models was accomplished through a trial-and-error process of adjusting initial values used in the models within an acceptable range. 

The model was assumed to be an adequate, although not 

unique, representation of the actual system when the 

standard error between field-measured and modelcalculated water levels was decreased to an acceptably 

small value and the model-calculated ground-water budget components agreed within an acceptable range with 

estimates from previous studies. The acceptability of the 

steady-state model depended on the adequacy of field 

estimates and area! distribution of water-level data.

Transient models were developed by using the steadystate model and were used to simulate water-level and 

flow-component changes for various simulated pumping
D12 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

periods. During the history-matching process, values for 

poorly known properties, such as specific yield, were varied within previously established limits of acceptability; 

values for the better known properties, such as pumpage, 

generally were not varied. This process provided a means 

of obtaining representative transient models and also 

demonstrated the effect of possible errors in estimating a 

poorly known property on simulation results.

Examination of model sensitivity to variations in property values used as input data was carried out in specificbasin models and basin-group models for steady-state and 

transient conditions. Properties were tested for sensitivity 

by using the simulations developed for basins in the category in which that property was thought to be most significant. A few properties, including transmissivity, specific 

yield, and recharge, were tested by using simulations 

developed for all categories of basins. Basin-group models were used mainly to investigate the effect of property 

generalization on model results.

SELECTED BASIN MODELS

UPPER SAN PEDRO BASIN 

GEOHYDROLOGIC CONDITIONS

Upper San Pedro basin includes 1,650 mi2 , of which 

about 950 mi2 is in the United States. The basin is in a 

north-trending structural trough in southeast Arizona 

(fig. 7). The sediments that compose the aquifer system 

include pre-Basin and Range rocks, lower and upper 

basin fill, and stream alluvium. The pre-Basin and Range 

rocks occur at depth in most of the basin and are exposed 

only near the mountain fronts; they are not considered an 

important part of the aquifer (Freethey, 1982, p. 7). The 

lower basin fill is the principal water-bearing unit and is 

hydraulically connected to the upper basin fill, which is 

saturated only in places. Stream alluvium is limited in 

areal extent and occurs along the flood plain of the San 

Pedro River and major tributaries.

Average annual precipitation ranges from 12 in. in the 

valley to more than 25 in. in the mountains that surround 

the basin (University of Arizona, 1965; Sellers and Hill, 

1974). Recharge to the aquifer occurs as infiltration along 

the mountain fronts and along the major stream channels. 

Ground water also enters the area from the south as 

underflow from Mexico. Ground water generally moves 

from the basin margin toward the central axis of the basin 

where most ground water is discharged along gaming 

reaches of the stream or lost to evapotranspiration by 

riparian vegetation. A small quantity of water leaves the 

basin as underflow. Estimates for inflow and outflow are 

given in figure 8.

The hydrologic system of the upper San Pedro basin is 

typical of basins of the southeast group. The upper San

Pedro basin was one of the initial basins used for groundwater simulation. The model was used to evaluate sensitivity to generalizations of hydrologic properties. Knowledge about model sensitivity gained in this effort was 

transferred to models of other basins in the southeast 

group. Inflow to the ground-water system is principally 

through mountain-front recharge; however, a significant 

quantity of water is recharged through infiltration of surface flow. Under predevelopment steady-state conditions, 

outflow occurred through evapotranspiration and discharge to surface flow of the San Pedro River. Basin fill 

and stream alluvium are stratigraphically complex. Water 

levels in wells drilled in the basin-fill sediments indicate 

that ground water in the aquifer is confined in places; 

however, ground water generally occurs under unconfined conditions. Only the upper more permeable part of 

the aquifer has been drilled and is used for ground-water 

withdrawal (Freethey, 1982).

Ground-water resources were developed initially for 

agricultural purposes. The first substantial pumping 

occurred in the 1940's, and in 1977 about 13,200 acre-ft of 

water was withdrawn (fig. 8). About 70 percent of the 

pumpage was for agricultural use; 14 percent was for 

industrial use; and 16 percent was for public-supply, 

domestic, and livestock uses (Konieczki, 1980). As a 

result of increased ground-water withdrawal, surface 

flow has decreased, underflow entering the area at the 

international boundary has increased, and discharge 

through evapotranspiration has decreased (Freethey, 

1982).

MODEL CHARACTERISTICS

The conceptual model of the hydrologic system of the 

upper San Pedro basin was based on available data and 

interpretations from previous studies (Bryan and others, 

1934; Heindl, 1952b; Brown and others, 1966; Roeske and 

Werrell, 1973). The numerical model was developed to 

simulate basinwide ground-water flow and to test the conceptual model. The numerical model developed and documented by Freethey (1982) used two layers to simulate a 

single aquifer system. The upper layer represents the 

upper 1,000 ft of the aquifer for which field data are 

available; the lower layer represents the aquifer deeper 

than 1,000 ft for which few data are available. Hydrologic 

properties representing the upper and lower model layers 

are given in figure 8. Quantities of recharge and values of 

hydrologic properties were adjusted within acceptable 

limits during the modeling process to improve the comparisons between field data and model results (Freethey, 

1982, p. 49).
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D13

The following model characteristics were adopted to 

simulate the hydrology of the upper San Pedro basin 

(Freethey, 1982, p. 15-16).

 A finite-difference grid with variable cell size was used 

to provide high resolution in areas where data density 

was high or where large variations in aquifer properties 

or stresses occurred.

 Two layers represented (1) the basin fill less than 1,000 

ft below land surface for which data were available and 

(2) the basin fill deeper than 1,000 ft below land surface for which few data were available. This arbitrary 

separation is the approximate maximum depth to 

which wells were drilled. Although the upper 1,000 ft of 

basin fill is of primary interest, the possibility of 

ground-water movement to or from the lower layer 

cannot be ignored. Hydrologic properties for this lower 

layer are unknown; however, values were selected in 

an assumed reasonable range to reflect what is known 

about the geology of the lower layer and to simulate 

the structural shape of the basin.

 The upper layer was simulated as, an unconfined 

aquifer; the lower layer was simulated as a confined 

aquifer.

 Vertical connection between layers was calculated by 

the model from the assigned hydrologic properties of 

each layer. The connection is represented quantitatively by the vertical leakance, which is the vertical 

hydraulic conductivity (K) of the fine-grained sediment divided by its thickness (b).

 Most recharge simulated by the model occurred in the 

upper layer. Specified heads in the lower layer allowed 

a small quantity of recharge directly into the lower 

layer.

 Interaction takes place between perennial streams and 

the upper layer. Stream leakance, which represents the 

hydraulic property that allows or restricts vertical 

movement of water between the stream and the aquifer, is constant.

 Evapotranspiration discharges water from the upper 

layer and was simulated by a linear relation between a 

maximum evapotranspiration rate and a depth to water 

of 10ft where evapotranspiration was assumed to 

cease.

The upper San Pedro basin model was used to evaluate 

the effect of using uniform representative values for 

hydrologic properties on model results, hereinafter 

referred to as generalizing. Freethey (1982, p. 16) lists the 

alternatives that were explored to evaluate model sensitivity to values and areal distributions of boundary 

recharge, aquifer hydraulic conductivity, and aquifer 

storage.

 Boundary recharge was evaluated by using two 

alternative modes: (1) distributed uniformly along 

mountain fronts to represent a situation of minimum

data availability and (2) distributed on the basis of sitespecific data and a flow-net analysis. 

Aquifer hydraulic conductivity was analyzed by comparing model results using (1) uniform values in three 

geohydrologically similar subareas of the basin and (2) 

an areal distribution on the basis of meager aquifer-test 

and specific-capacity data and a flow-net analysis. 

Aquifer storage was examined by comparing three 

model simulations using (1) a variable distribution of 

specific yield estimated from drillers' logs, (2) a uniform value for specific yield based on few data in the 

basin, and (3) the uniform specific-yield value used in 

the previous simulation and uniform values for boundary recharge and aquifer hydraulic conductivity. 

Generalized boundary recharge, hydraulic conductivity, and specific yield were used together to represent a 

crude approximation of the hydrologic system that 

might be developed where few data are available. The 

results were compared with the final calibrated model 

that used all the available information.

RESULTS OF SIMULATIONS

After steady-state model history matching, several 

model runs were made by using an area-weighted arithmetic average value to replace a single data array used in 

the model. Values for mountain-front recharge, hydraulic 

conductivity, and specific yield were independently varied 

in this manner to determine the effects on model results. 

These results were then compared with results of the 

steady-state simulation that most accurately matched 

measured water levels. Evapotranspiration rate and 

extinction depth, stream leakance, and ratios for vertical 

anisotropy of hydraulic conductivity remained unchanged 

during this process.

Initially, mountain-front recharge was distributed 

around the model perimeter on the basis of a flow-net 

analysis, and recharge varied from cell to cell (nongeneralized). Subsequently, a uniform average value was specified at each recharge cell (generalized). Generalizing 

mountain-front recharge resulted in simulated water levels that were more than 25 ft lower than water levels in 

the history-matched steady-state model over about 17 percent of the modeled area. The areas of greatest difference 

coincided with those areas of steep hydraulic gradient 

near the mountain fronts. The greatest water-level differences correspond with the greatest differences between 

the quantities of nongeneralized and generalized 

recharge. Water levels that are higher than those of the 

history-matched steady-state simulation are caused by the 

generalized recharge quantity that is greater than the nongeneralized recharge quantity. Water-level increases typically were less than 20 ft; however, the area influenced is
D14 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

ESTIMATED AND SIMULATED WATER-BUDGET COMPONENTS 

[Values are in acre-feet per year; dashes indicate no data. Estimated values modified from Freethey (1982, table 2)]

Transient state

Steady state

Components

Inflow:

Storage ...................

Total ...................

Outflow:

Underflow to Benson basin . . .

Pumpage .................

Total ..................

Estimated

10,000-11,000

800- 4,300

700- 3,500

0

1,900-14,300

0 

3,600-12,300

0

5,500 - 26,600

Simulated

8,300 

800 

7,400 

0

16,500

8,300 

400 

7,800 

0

16,500

1968

Estimated

10,000- 11,000

10,000-11,000 

13,500

13.500

Simulated

8,300 

1,200 

7,500 

6,000

23,000

5,400 

300 

6,500 

10,800

23,000

1977

Estimated

10,000- 11,000

10,000- 11,000 

13,200

13,200

Simulated

8,300 

1,400 

7,700 

5,600

23,000

5,900 

300 

6,200 

10,600

23,000

20 1

SIMULATED PUMPING PERIODS

_^^_

3

| Estimated annual pumpage 

D Simulated annual

pumpage (10

pumping periods)

-\25

50 MILES

0 50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

UPPER SAN PEDRO BASIN (SHADED)

115

125

135

80

90

WELL A

WELLB

4,245

100

1940

4,235

UJ

>

UJ

_l 

<

UJ 

CO

LJLJ

>

O 

4,225 5

4,180 m

4,170

1950 1960 1970

4,160

1980

FIGURE 8. Location of basin, estimated and simulated water-budget components and pumpage, ground-water levels, model gnd and boundary, and 

observed and simulated water levels, upper San Pedro basin.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D15

110°30'

0°00'

0 | | | 5 MILES 

0 5 KILOMETERS

Base from U.S. Geological Survey 

State base map, 1:500,000, 1974

5 KILOMETERS

EXPLANATION

Alluvial deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating groundwater flow

Cell representing

Constant head 

Evapotranspiration 

Mountain-front recharge 

River-aquifer connection

--4200-

 4250-

 B

Geologic contact

Boundary of ground-water area

Boundary of modeled area 

Lower and upper layers

Water-level contour, 1977-78 

Shows altitude of water level. 

Contour interval 100 feet. Datum is 

sea level

Observed (Konieczki, 1980) 

Simulated (Freethey, 1982) 

Well for which hydrograph is shown

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Range 

Average

Specific yield of upper layer

Storage coefficient of lower layer

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

0.1 - 18.3 

4.1

9 - 2,307 

684

0.08

10-5

5x1 Om7- 2x1 0'9

FIGURE 8. Continued.
D16 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

limited to a narrow zone along the basin edge where 

recharge occurs, gradients are large, and close matching 

of field- and model-head values are difficult to attain. 

These results indicate that generalizing recharge along a 

basin boundary does not greatly affect overall modeling 

results but can cause nonrepresentative water levels to be 

simulated in site-specific areas near the mountain fronts.

Specifying the same value of hydraulic conductivity 

throughout the basin resulted in water levels that were 

more than 25ft lower than simulated water levels in 

about 53 percent of the modeled area. The model has a 

greater response to the generalization of hydraulic conductivity than it does to generalization of recharge quantity because the generalized hydraulic-conductivity value 

represents the arithmetic average of a data set that has a 

large range in values. The average value in the model can 

differ markedly from the field value, and large head differences will result, particularly in the basin center where 

field values of hydraulic conductivity are low. A conclusion of this analysis is that use of generalized hydraulicconductivity values in a basin model is transferable but is 

not considered a desirable approach. The effectiveness in 

simulating an actual hydrologic system by using a generalized hydraulic-conductivity value depends on the variability of the hydraulic conductivity and the accuracy of the 

field information on which the generalized value is based.

A transient analysis of the model was used to evaluate 

model sensitivity to a uniform specific-yield value. A uniform specific-yield value (0.08) resulted in from 4 ft less 

to 1 ft more drawdown in the final pumping period of the 

simulation and in virtually the same volume of storage 

depletion in the aquifer. Generalizing recharge, hydraulic 

conductivity, and specific yield resulted in drawdown in 

the final pumping period that ranged from 8 ft less to 7 ft 

more than that in the simulation using nongeneralized 

values. The pumpage that represented a removal of water 

from storage increased from 53 to 58 percent; other 

flow components remained within previously estimated 

ranges.

The numerical model developed for the upper San 

Pedro basin simulated the hydrologic system to a degree 

of accuracy acceptable for most purposes related to 

water-resources appraisals. The accuracy, or how closely 

the simulation represents the actual hydrologic system, is 

a function of the dynamics and variability of local hydrologic conditions and the magnitude of pumpage relative to 

other water-budget components. By use of current knowledge and definition of the system, the numerical model 

can simulate the key hydrologic processes presumed to 

be taking place. The model can be used to simulate areal 

water levels that approximate the field data and to simulate ground-water-budget values that approximate values 

estimated in previous investigations (fig. 8). These results 

demonstrate that the model is a reasonable, though not

unique, simulation of the real system. Predictive capabilities are limited to a general assessment of changes in 

inflow, outflow, and storage values and changes in 

regional-flow directions as a result of pumping. The purpose of the model was not to analyze site-specific 

ground-water conditions or to predict water-level changes 

in individual wells. Reliability of the predictive capability 

of the model also depends on future changes in the hydrologic system and how well the changes can be incorporated into numerical representations.

BENSON BASIN 

GEOHYDROLOGIC CONDITIONS

Benson basin includes about 900 mi2 in southeastern 

Arizona and is in the same north-trending structural 

trough as the upper San Pedro basin (fig. 9). Benson basin 

is downstream from the upper San Pedro basin; outflow 

from the upper San Pedro basin is inflow to Benson basin. 

About half the surface area is bedrock of the mountains 

that forms the east and west sides of the basin. The basin 

is about 32 mi long, and sedimentary deposits that make 

up the valley floor average 15 mi wide. The aquifer system 

is composed of the same units as in the upper San Pedro 

basin, and basin fill and stream alluvium are stratigraphically complex. In the Benson basin, however, the lower 

basin fill contains a gypsiferous, fine-grained facies as 

much as 800 ft thick near the basin center. This facies 

forms the confining unit for ground water in underlying 

lenses and layers of coarse-grained material. Water from 

wells that penetrate the coarse-grained material may flow 

at the land surface. On the perimeter of the fine-grained 

facies and in the overlying sediments, ground water 

occurs under unconfined conditions. Vertical leakage 

upward probably occurs through the confining unit and is 

a source of inflow to the upper, unconfined part of the 

aquifer.

Average annual precipitation ranges from less than 12 

in. in the valley to more than 25 in. in the mountains that 

bound the basin (University of Arizona, 1965). Recharge 

to the ground-water system takes place at the mountain 

fronts and along the San Pedro River; ground-water movement is from the mountain fronts toward the central axis 

of the basin. Most discharge occurs as evapotranspiration 

in the zone of riparian vegetation along the San Pedro 

River. A small quantity of discharge occurs as underflow 

at the north end of the basin. Estimated and simulated 

inflow and outflow are given in figure 9.

Ground-water development began in the early 1900's, 

and pumpage is estimated to have been less than 3,000 

acre-ft in 1940. In 1967, ground-water withdrawal was 

12,400 acre-ft; withdrawal increased to nearly 40,000 

acre-ft in 1976 and was 35,000 acre-ft in 1977 (fig. 9).
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D17

Pumping of ground water in Benson basin resulted in a 

decrease in the number of flowing artesian wells and a 

decrease in the discharge from wells that continue to 

flow. Reaches of the San Pedro River that were perennial 

now flow only intermittently.

MODEL CHARACTERISTICS

A conceptual model of the ground-water system of 

Benson basin was formulated on the basis of results of 

previous studies (Bryan and others, 1934; Roeske and 

Werrell, 1973; Konieczki, 1980) and through analysis of 

available geologic and hydrologic data Pertinent data for 

Benson basin are sparse, and the approach and many of 

the data used in the model were transferred from the 

results and modeling approach used in the upper San 

Pedro basin (Freethey, 1982).

The numerical model of the Benson basin had the following characteristics:

 A variable grid size for higher resolution along the central axis of the basin (fig. 9).

 Two layers representing (1) the stream alluvium that 

extends along the axis of the basin but does not extend 

to the mountain boundaries and in which ground water 

is unconfined and (2) the lower layer that extends 

throughout the basin and in which ground water is confined in an area that underlies the basin-center 

fine-grained facies and unconfined outside this area 

extending to the mountain fronts.

 A simulated "quasi" third layer connects the two aquifers and simulates vertical flow through the 

fine-grained unit. This layer is represented quantitatively by a vertical-leakance value, K/b, which is the 

vertical hydraulic conductivity divided by the thickness of the fine-grained unit.

 Recharge along the mountain fronts is to the lower 

layer, which is the only active layer adjacent to the 

mountains.

 Interaction takes place between the San Pedro River 

and the stream alluvium. Stream leakance, which represents the hydraulic property that allows or restricts 

vertical movement of water between the stream and 

the aquifer, is constant.

 Ground water discharges by evapotranspiration from 

the upper layer.

 Pumping was simulated from both layers. Irrigation 

pumpage from the upper layer was reduced by 30 percent to account for return of excess applied irrigation. 

Thirty percent of the irrigation pumpage from the 

lower layer was recharged to the upper layer. 

The finite-difference grid for the Benson basin model 

consisted of 19 columns and 30 rows making 570 rectangular cells in each layer that were from 0.36 to 1.0 mi2 in

area (fig. 9). The upper layer was represented by 91 active 

cells, and the lower layer had 367 active cells. Mountainfront recharge was applied to the lower layer at a uniform 

rate per mile of mountain front. This rate differed for each 

mountain range on the basis of altitude of the mountains, 

which is related to the volume of annual precipitation in 

the mountains. Underflow into the basin at Fairbank and 

out of the basin at the Narrows was simulated by 

constant-head cells in both model layers. Upward leakage 

from lower to upper layers was controlled by the verticalleakance value.

RESULTS OF SIMULATIONS

To evaluate the numerical-model simulation, steadystate and transient model runs were made. The steadystate model was used to determine if the water-budget 

components and hydrologic properties of the aquifer that 

were estimated on the basis of existing information functioned together as a reasonable representation of the 

hydrologic system. The transient model was used to evaluate the reliability and uniqueness of the steady-state 

model. The steady-state model for Benson basin verified 

that all the proposed hydrologic concepts were compatible in the simulation and that values of simulated flow 

components approximated the values of flow components of the conceptual model (fig. 9). Outflow from 

the upper San Pedro basin was used as inflow to the 

Benson basin. Simulated water levels for the predevelopment period were the same general configuration as measured water levels. The results of simulated transient 

conditions indicated that evapotranspiration, streamaquifer exchange, and upward leakage from the lower 

layer decreased as a result of pumping and that as much 

as 10,000 acre-ft of water was removed from storage in 

the upper aquifer through 1977. On the basis of the 

approximate areal extent, saturated thickness, and specific yield of this unconfined aquifer, 10,000 acre-ft represents only about 1 to 2 percent of the available ground 

water in storage. Water-level declines are greater than 20 

ft only locally near pumping wells.

The transient model simulated historical pumping conditions for three time periods from predevelopment 

(arbitrarily selected as 1915) through 1945, 1946 through 

1967, and 1968 through 1977 (fig. 9). Water-level declines 

have occurred primarily in the area underlying the flood 

plain of the San Pedro River where most pumping has 

occurred. Simulated head declines at the end of 1945 typically were less than 5 ft in the unconfined aquifer and less 

than 10 ft in the confined aquifer. Field evidence to corroborate these head declines consists only of early observations about the decreasing discharge of flowing 

artesian wells (Roeske and Werrell, 1973). Simulated heaci
D18 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

ESTIMATED AND SIMULATED 

WATER-BUDGET COMPONENTS

[Values are in acre-feet per year. Estimated values modified from 

Heindl (1952b) and unpublished records of the U.S. Geological Survey]

Steady stale

Components Estimated Simulated

Inflow:

Mountain-front recharge ......... 8,100 - 10,500 9.400

River losses .................. 200 - 9,700 1,400

Underflow from upper basin ...... 400_________600

Total ...................... 8,700 - 20,600 11,400

Outflow:

River gains ................ .. 900 - 5,900 4,000

Underflow to lower basin ........ 100- 200 300

Evapotranspiration ............. 5,100 - 16,700____7,100

Total ...................... 6,100-22,800 11,400

SIMULATED PUMPING PERIODS

COH 

33 UJ

50

40

30

10

| Estimated annual pumpage 

D Simulated annual pumpage 

(three pumping periods)

NOTE: The first simulated 

pumping period began in 

1915

70

60 ?i

s§ 50 gm

2"

40 . m

0.0-

30

115 109°

20

10 55

0

108°

<u50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

BENSON BASIN (SHADED)

10

20

Well A

Unconfined aquifer

3,648

3,638

3,628

3,618

3,695

- 3,685

- 3,675

- 3,665

1940 1950 1960 1970 1980

3,655

FIGURE 9. Location of basin, estimated and simulated water-budget components and pumpage, ground-water levels, model grid and boundary, and 

observed and simulated water levels, Benson basin.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS

now

D19

Base from U.S. Geological Survey 

State base map, 1:500,000, 1974

5 KILOMETERS

EXPLANATION

Alluvial deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating groundwater flow

Cell representing

I Constant head 

Ejl Evapotranspiration 

@ Mountain-front recharge 

0 River-aquifer connection 

FIGURE 9. Continued.

-- 3500

 3500- 

 B

Geologic contact

Boundary of ground-water area

Boundary of modeled area 

Lower layer 

Upper layer

Water-level contour, 1977-78 

Shows altitude of water level. 

Contour interval 100 feet. Datum is 

sea level

Observed (Konieczki, 1980) 

Simulated

Well for which hydrograph is shown

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Range 

Average

Specific yield 

Lower layer 

Upper layer

Storage coefficient of lower layer

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

2-45 

17.1

11 -4,445 

832

0.114 

0.125

1C'4

10-7-10'11
D20 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

decline at the end of 1967 typically ranged from 5 to 15 ft 

in the unconfined aquifer (hydrograph for well A, fig. 9) 

and from 10 to 25 ft in the confined aquifer (hydrograph 

for well B, fig. 9). At the end of 1977, head declines in 

most of the area ranged from 10 to 40 ft in the unconfined 

aquifer and 15 to 50 ft in the confined aquifer. A local area 

near the center of the modeled area had simulated head 

declines of 50 to 100 ft for both aquifers. This head 

decline was not documented by field data. The discrepancy could be caused entirely or in part by the method 

used to distribute simulated pumping, by no-flow boundaries in the upper layer of the model, by incorrect hydrologic properties assigned to the area, or by a combination 

of these factors.

LOWER SAN PEDRO BASIN

GEOHYDROLOGIC CONDITIONS

Lower San Pedro basin is the northernmost extension 

of the narrow north- to northwest-trending structural 

trough that forms the upper San Pedro and Benson 

basins. The lower San Pedro basin includes about 1,700 

mi2 in southeastern Arizona (fig. 10). About 40 percent of 

the surface area is basin fill, and 60 percent is bedrock of 

the mountains. The alluvial basin is about 75 mi long and 

averages less than 10 mi wide. Aquifer lithology in the 

lower San Pedro basin is similar to that of the upper San 

Pedro and Benson basins. The ground-water system consists of an upper unconfined aquifer and a lower confined 

aquifer separated by an intermediate fine-grained unit. 

The upper aquifer consists mainly of stream alluvium 

underlying the flood plain of the San Pedro River and is of 

limited areal extent. The lower aquifer consists of thin 

sand layers underlying a fine-grained unit.

Average annual precipitation ranges from less than 12 

in. on the valley floor to more than 30 in. in the mountains 

that bound the basin (University of Arizona, 1965). 

Ground water moves from areas of recharge along the 

mountain fronts and along stream channels to areas of 

discharge along the central axis of the basin. Discharge as 

underflow from the upstream Benson basin enters the 

lower San Pedro basin at the south end. Some ground 

water moves vertically from the deep confined aquifer to 

the upper unconfined aquifer as leakage through the finegrained unit. Most discharge occurs as evapotranspiration 

in the zone of riparian vegetation in the flood plain of the 

San Pedro River or as base flow to the river. Heindl 

(1952a) estimated that 35,000 acre-ft of water discharged 

from the basin through evapotranspiration annually.

Ground water is the primary source of water used in 

the lower San Pedro basin. Many ditches diverted streamflow to irrigate small acreages prior to 1940 and probably 

still do when streamflow is available. After World War II, 

however, ground-water withdrawal increased markedly to 

about 20,000 acre-ft in 1951 (Heindl, 1952a) and increased 

gradually to about 23,000 acre-ft in 1960 (fig. 10). From 

1965 through 1976, pumpage increased to more than 

50,000 acre-ft/yr (U.S. Geological Survey, 1985) and has 

remained near 50,000 acre-ft/yr since 1976. Most pumpage 

for irrigation is from the stream deposits along the San 

Pedro River in which ground water is unconfined. Most 

pumpage for industrial and public-supply uses is from the 

deep aquifer system in which ground water is confined. 

Before 1953, all ground-water withdrawal was for irrigation; however, by 1976, 53 percent was for uses other than 

irrigation. A summary of the ground-water pumpage history in the lower San Pedro basin is included in figure 10.

MODEL CHARACTERISTICS

The conceptual model of the hydrologic system of the 

lower San Pedro basin is similar to those of the upper San 

Pedro and Benson basins. The results of previous studies 

provided the basis for developing the initial model; these 

studies include Bryan and others (1934), Heindl (1952a), 

Page (1963), Roeske and Werrell (1973), and Jones (1980). 

The numerical model was used principally to test the 

transferability of information on mountain-front recharge, 

river-aquifer interaction, evapotranspiration, transmissivity, hydraulic conductivity, and vertical interlayer leakance. The characteristics of the developed models of the 

upper San Pedro and Benson basins were used as initial 

property values for the model of the lower San Pedro 

basin.

The numerical model developed for the lower San 

Pedro basin had the following characteristics and 

assumptions:

 A variable column width is specified to enhance resolution along the axis of the basin where discharge to the 

stream and evapotranspiration occurs. Row width is a 

uniform 1 mi (see fig. 10).

 Two layers representing (1) the stream alluvium along 

the flood plain of the San Pedro River in which ground 

water is unconfined and (2) the basin-fill sediment 

throughout the basin extending to the mountain fronts. 

Ground water is confined in the basin-fill layer underlying the basin-center fine-grained facies and unconfined 

outside this area extending to the mountain fronts.

 The fine-grained facies separating the stream alluvium 

from the underlying confined aquifer is represented 

quantitatively by the vertical-leakance value, which is
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D21

the vertical hydraulic conductivity of the fine-grained 

facies divided by its thickness.

 Mountain-front recharge was applied to the lower layer 

only.

 Stream-aquifer interchange takes place and is controlled by the stream elevation, stream-bottom elevation, stream-leakance value, and head in the upper 

layer.

 Evapotranspiration was simulated as a discharge of 

ground water from the upper aquifer and controlled by 

a linear relation between the maximum evapotranspiration rate and a depth to water of 10 ft, where evapotranspiration was assumed to cease.

 Transient simulation used two pumping periods in 

which withdrawals from both layers were represented. 

Irrigation pumpage, which was simulated as withdrawal from the upper layer, was reduced by 30 percent to account for return of excess applied irrigation 

water.

The finite-difference grid for the lower San Pedro basin 

model (fig. 10) consisted of 18 columns and 69 rows making 1,242 rectangular cells in each layer from 0.4 to 1 mi2 

in area. The upper layer included 248 active cells, and the 

lower layer included 865 active cells. Mountain-front 

recharge was introduced into the lower layer as a constant flow at a rate that varied depending on the altitude 

of the adjacent mountain range and average annual precipitation in the mountains. Underflow into the basin at 

the Narrows (south end) and out of the basin at Winkelman (north end) was simulated by constant-head cells; 

underflow into the basin at the mouth of Aravaipa Canyon 

was simulated as a constant flow. Vertical leakage 

through the fine-grained confining unit was controlled by 

an array representing leakance values (K'/b~) and ranged 

from 10~9 to 10~6 (ft/s)/ft, depending on the thickness and 

character of the fine-grained facies.

RESULTS OF SIMULATIONS

The steady-state model was adjusted through a trialand-error process to simulate the general shape of the 

estimated water-level configuration and to obtain flowcomponent values that were reasonable approximations 

of estimated water-budget values. Because correlative 

field data were sparse, however, a detailed definition of 

the hydrologic system could not be established. Matching 

measured water levels to within appropriate accuracy criteria was not done because actual predevelopment potentiometric levels for the confined aquifer were not 

available. After an acceptable steady-state simulation was 

obtained, the model was stressed by simulated pumping. 

The response of the transient model was compared with 

available field evidence to determine qualitatively the

representation of the actual ground-water system by the 

model.

Steady-state simulation produced a hydraulic-head 

configuration similar to the estimated predevelopment 

potentiometric surface, although differences of as much 

as +50 ft were present in places near the mountain fronts. 

Water-budget outflow components generally were smaller 

than estimated outflow components in previous investigations; however, inflow components generally matched 

estimated inflow components (fig. 10). The difference is 

in the balance of inflow with outflow. Previous estimated 

inflow-outflow values were not balanced, whereas the 

model inflow-outflow values were balanced. The model 

can be adjusted to simulate larger or smaller inflow and 

outflow. For example, inflow and outflow components 

can be increased by a factor of two by increasing 

mountain-front recharge, transmissivity, vertical leakance, and hydraulic conductivity, with little or no resulting change in water-level configuration. The model inflow, 

however, would then be twice the inflow estimated on the 

basis of field data. The model is not sufficiently sensitive 

to variations in water-budget components to define 

flow components or aquifer characteristics with more 

precision.

The transient simulation represented two pumping 

periods. The termination of these periods represents the 

ends of 1967 and of 1978, which are times for which 

water-level maps are available (Roeske and Werrell, 1973; 

Jones, 1980). Differences in the water-level configuration 

shown on the two maps are small; however, water-level 

measurements in a deep well near Mammoth indicated a 

water-level decline of about 20 ft between 1967 and 1978. 

Model results indicate an average water-level decline of 

about 4 ft in the lower layer for the same period and localized declines of as much as 150 ft. The reason for this 

large difference between field and simulated conditions is 

unknown but could be because data are insufficient to 

define local geohydrologic conditions and pumping distributions. Between 1967 and 1978, net pumpage from the 

aquifer approximately doubled. Model results indicated 

that storage depletion increased from about 1,700 acre-ft 

in 1967 to 8,250 acre-ft in 1978 (fig. 10). The model results 

also indicate that natural discharge by evapotranspiration 

decreased from nearly 30,000 acre-ft/yr in the early 1940's 

to less than 23,000 acre-ft in 1967 and to about 16,000 

acre-ft in 1978 as a result of water-level declines. Model 

results also indicate that base flow supplied by aquifer 

discharge decreased from more than 22,000 acre-ft/yr in 

the early 1940's to about 10,000 acre-ft in 1967 and to less 

than 6,000 acre-ft in 1978. The head decreases in both 

parts of the aquifer system were accompanied by a 

decrease in the net upward flow from the lower aquifer to 

the upper aquifer.
D22 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

ESTIMATED AND SIMULATED 

WATER-BUDGET COMPONENTS

[Values are in acre-feet per year: dashes indicate no data. Estimated 

values modified from Heindl (1952a) and unpublished records of 

the U.S. Geological Survey]

Steady state

Components Estimated Simulated

Inflow:

Mountain-front recharge ......... 25,000 - 27,500 26,400

River losses .................. 2,200 900

Underflow from Benson basin..... 1,200_______1,700

Total ..................... 28,400 - 30,900

Outflow:

River gains ................... 21,700

Underflow to Gila River ......... 2,700

Evapotranspiration ............. 12,700-38,000

Springs ...................... ___3,200_______ 

Total ...................... 40,300-65,600 29,100

29,000

11,400

2,100

15,600

60

SIMULATED PUMPING PERIODS

1

| Estimated annual pumpage

D Simulated annual pumpage 

(two pumping periods)

0 50 MILES

0 50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

LOWER SAN PEDRO BASIN (SHADED)

DU

H CC 

LU LU

111 LL

LL ^ 40

LU O

cc. **:

o o

< iLL ^

O O 20

tO LL

Q 2 z - 

< 

to

g o

I

z £

 LL

> 3

H 0 20

z 5

< O

=3 CC 

0 ^

5 0 40

O !j -i LL

LL KO

cr\

-

Z

o

" H

0

Z

- o 

o

LU ""<""

 1-

to

- 1 >-

Q

- rjj '"H"

to

-

Mountain-front recharge

-

-

Storage

0. ,. depletion Streamflow \T

infiltration x V^ 

\ Underflow in _ ^^- "

^s^- -. -a- -

'--.._ Underflow out *. ""

Discharge "^»^^ . ""

to river .I'^-C "*

(base flow) -^ ^ . " *^ % -.«*

 . ^ * "*" x » « ""

p "^ « *^ *

.---- """""" % \

Evapotranspiration \ 

Net \

pumpage\

\

i i , i . i ,

1940 1950 1960 1970 1980

SIMULATED CHANGES IN WATER-BUDGET COMPONENTS

AS A RESULT OF INCREASED GROUND-WATER 

WITHDRAWAL IN LOWER SAN PEDRO BASIN, 1940-78

2,552

1950 1960 1970

1,902

1980

FIGURE 10. Location of basin, estimated and simulated water-budget components and pumpage, ground-water levels, model gnd and boundary, 

simulated changes in water-budget components, and observed and simulated water levels, lower San Pedro basin.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D23

EXPLANATION

Base from U.S. Geological Survey 

State base map, 1:500,000, 1974

5 MILES j

5 KILOMETERS

Alluvial deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating groundwater flow

Cell representing

Constant head 

Evapotranspiration 

Mountain-front recharge 

River-aquifer connection 

Vertical leakance

2000---

-2000 

 B

Geologic contact

Boundary of ground-water area

Boundary of modeled area 

Lower layer 

Upper layer

Water-level contour, 1978-79 

Shows altitude of water level. 

Contour interval 250 feet. Datum 

is sea level

Observed (Jones, 1980) 

Simulated

Well for which hydrograph is shown

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Range 

Average

Specific yield 

Lower layer 

Upper layer

Storage coefficient of lower layer

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

16-32 

31.6

67 - 5,346 

947

0.086 

0.146

1C'4

10-6. irj-9

FIGURE 10. Continued.
D24 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

WILLCOX BASIN 

GEOHYDROLOGIC CONDITIONS

Willcox basin is a topographically closed basin in the 

southeastern part of the study area (fig. 7). Although 

small quantities of ground water may have discharged to 

adjacent basins before development, the basin is internally drained, and the basin center is a large, barren playa. 

The bedrock of the mountains crops out in about 

one-third of the 1,660-mi2 basin, bounds the basin, and 

forms its physical and hydrologic boundaries. The 

remaining two-thirds of the basin area is underlain by 

alluvial sediments that form the aquifer system.

The sediments that fill Willcox basin consist of preBasin and Range rocks, lower and upper basin fill, and 

stream alluvium. Pre-Basin and Range rocks occur at 

depth in the basin and are not penetrated by wells. Lower 

and upper basin fill form the aquifer system. Coarsegrained facies of these units occur near the mountain 

fronts, and very fine grained facies, composed primarily 

of silt and clay, are present at the basin center underlying 

and forming the playa. Near the playa, the sediments are 

heterogeneous and consist of interfingering fine and 

coarse facies. Stream alluvium consists primarily of sands 

and gravels that underlie the many washes that enter the 

basin from the mountains.

Average precipitation ranges from less than 12 in./yr on 

part of the valley floor to more than 35 in./yr in the surrounding mountains (University of Arizona, 1965). 

Ground water originates as precipitation falling in the 

mountains. Recharge is a result of infiltration of surface 

runoff into the coarse sediments along the mountain 

fronts. Before development, ground water moved from 

the mountain fronts toward the playa where discharge 

occurred through spring flow, evaporation, and transpiration (fig. 11).

Development of ground-water resources in Willcox 

basin paralleled that of other basins that have agriculturally oriented uses. Ground-water withdrawal generally 

was less than 5,000 acre-ft/yr until the end of World War II 

(U.S. Geological Survey, 1982). Withdrawal increased 

gradually to about 40,000 acre-ft in 1952; however, in 1953, 

ground-water withdrawal more than doubled to about 

94,000 acre-ft. Pumping increased steadily after 1953 to 

about 300,000 acre-ft in 1967 (fig. 11) (U.S. Geological Survey, 1982). Since 1967, rates of ground-water withdrawal 

generally have remained stable. Although exact quantities 

are not well documented, the quantity of ground water 

used for domestic, livestock, public supply, and industrial 

purposes is small compared to that used for irrigation.

MODEL CHARACTERISTICS

A numerical model of Willcox basin was developed to 

examine the hydrologic similarities of this basin and others of the southeast group. Only a small number of data 

were available to provide the initial values of hydrologic 

properties for use in the model. Previous studies by Meinzer and Kelton (1913), Coates (1952), Brown and others 

(1963), Brown and others (1966), and Brown and Schumann (1969), provided the basis for model development. 

One objective of developing this model was to determine 

if the relation between precipitation and mountain-front 

recharge that was used in the model for the upper San 

Pedro basin remained valid when used in the model for 

Willcox basin.

The numerical model was developed to simulate 

steady-state and transient conditions and had the following characteristics:

 A uniform row and column width of 1 mi.

 Two layers representing (1) the aquifer system in the 

upper 1,000 ft of basin fill from which most withdrawal 

occurs and (2) the aquifer system below 1,000 ft for 

which few field data were available (fig. 11).

 The two aquifer layers were connected hydraulically 

by a leakance value to allow vertical flow. Confined 

ground-water conditions, prevalent at the edge of the 

playa on the east side, could not be adequately represented because of the model scale. Representation of 

these areas as unconfined probably caused actual 

drawdown to exceed simulated drawdown in those 

areas.

 Perimeter recharge applied to the upper model layer 

generally was uniform within the influence of each 

individual range of mountains. Simulated recharge 

from the Chiricahua Mountains on the southwestern 

perimeter of the basin was distributed over a large area 

of the model on the basis of chemical characteristics of 

the ground water (EN. Robertson, hydrologist, U.S. 

Geological Survey, oral commun., 1981).

 Drain cells were used to represent spring discharge at 

four points at the margin of the playa. Simulated spring 

discharge was dependent on the difference between 

head in the aquifer and altitude of land surface.

 In areas where depth to water was less than 25 ft 

(Meinzer and Kelton, 1913), evaporation from the playa 

and transpiration by mesquite and grass immediately 

adjacent to the playa were simulated with a single 

maximum evapotranspiration rate and extinction 

depth of 25 ft.

 Four periods of pumping were simulated. Irrigation 

return flow was not incorporated, although it probably 

does occur and could have been simulated to improve 

the match of measured and simulated water levels. The 

simulated pumping periods are a compromise between
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D25

periods of generally uniform pumping and periods for 

which pumpage and water-level records are available. 

The finite-difference grid used for simulating the 

hydrologic system of Willcox basin consisted of 55 rows 

and 24 columns (fig. 11). A total of 711 cells were active in 

each of the two model layers. Pumpage data for the transient runs were distributed on the basis of well location 

and date of drilling because pumpage data for individual 

wells were not available.

RESULTS OF SIMULATIONS

The steady-state (predevelopment) model was 

matched to the water-level configuration of 1910 (Meinzer 

and Kelton, 1913). The transmissivity distribution 

obtained for the upper layer by using transmissivity estimated by Brown and Schumann (1969) worked well in the 

model and, for the most part, was not altered. A recharge 

rate of 47,000 acre-ft/yr along the mountain fronts provided the most reasonable match of hydraulic properties 

and water-level configurations. Estimates of mountainfront recharge from earlier investigations ranged from 

20,000 to 75,000 acre-ft/yr.

The transient model was calibrated for 1910-74 with 

intermediate pumping periods ending in 1946, 1952, and 

1962. Simulated water levels for the first 50 years of the 

transient period were higher than measured water levels, 

but in the last 10 to 15 years of the transient period, simulated water levels were lower than measured water levels. 

This difference may indicate that a return of excess 

applied irrigation water to the aquifer may have occurred 

after a significant time lag or that transmissivity and specific yield used in the simulation change significantly with 

depth. Measured water levels in the Kansas Settlement 

area (fig. 11) were consistently lower than simulated 

water levels. The difference may be caused by the extensive interfingering of coarse material extending from the 

basin perimeter toward the playa and fine-grained lacustrine deposits of the playa. The interfingering, in turn, 

causes water-level response to pumping to be more characteristic of a confined aquifer system rather than an 

unconfined aquifer system as was simulated. The deep, 

steep-sided cone of depression also suggests a confinedaquifer system response. In addition, much of the pumping in the Kansas Settlement area occurs close to the 

playa. The very fine grained lacustrine material that 

underlies the playa could effectively create a boundary 

condition that results in increased drawdown as the cone 

of depression intersects this poorly permeable unit.

The regional geohydrologic character of the Willcox 

basin is similar to that of other basins in southeastern Arizona. Recharge and discharge mechanisms are similar 

except that no underflow occurs into or out of adjacent

basins. The distribution of pumping in the basin, however, 

is not similar because of a clay deposit at the basin center. 

The clay deposit extends from land surface to below the 

present depths of pumping wells. Boundary conditions 

created by this clay deposit may increase the rate of storage depletion, and delayed drainage may be a long-term 

consideration.

The quantity of mountain-front recharge into the basin 

was evaluated by using the steady-state model. The initial 

value simulated, about 40,000 acre-ft/yr, was based on the 

same relation between total rainfall greater than 8 in. and 

recharge used for the upper San Pedro basin and other 

basins in the study area (Anderson and others, 1992). The 

transmissivity of the upper model layer, which was simulated as the product of saturated thickness and hydraulic 

conductivity, was deemed to be one of the most reliable 

property values and was not altered during simulation. To 

accommodate different recharge to the system, transmissivity of the lower model layer was altered. Values of simulated mountain-front recharge tested in the steady-state 

model were about 24,000, 47,000, and 74,000 acre-ft/yr. 

After an initial model run, the largest value was not tested 

further because the transmissivity of the lower unit 

required for that annual rate of recharge was considered 

unreasonably large. The two smaller recharge values 

were tested further in the transient-model analysis. The 

most reasonable match between water-level-decline rates 

measured in observation wells and water-level-decline 

rates simulated was obtained by using the recharge rate 

of 47,000 acre-ft/yr.

DOUGLAS BASIN 

GEOHYDROLOGIC CONDITIONS

Douglas basin in southeastern Arizona includes about 

1,440 mi2 in the United States and Mexico (fig. 7). About 

840 mi2 is underlain by basin fill; the remaining area is 

bedrock of the mountains. Douglas basin is geohydrologically similar to the upper San Pedro basin. The principal 

aquifer is composed largely of fine-grained sediments in 

both the upper and lower basin fill, and lithologic differences between the units are difficult to distinguish. The 

sediments grade from coarse near the mountain fronts to 

as much as 80-percent silt- and clay-size particles near the 

basin center. Ground water generally is unconfined in the 

basin; however, ground water is confined in places, particularly at depth in the southern part of the basin. The 

model of Douglas basin served as a means of testing and 

verifying the applicability of information transfer and testing the empirical relation between precipitation and 

mountain-front recharge. Previous studies by Meinzer and 

Kelton (1913), Coates and Cushman (1955), White and
D26 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

ESTIMATED AND SIMULATED 

WATER-BUDGET COMPONENTS

[Values are in acre-feet per year; dashes indicate no data. Estimated 

values derived from information in Meinzer and Kelton (1913), Coates 

(1952), Brown and others (1963), and Mann and others (1978)]

Steady-state

Components Estimated Simulated

Inflow:

Mountain-front recharge. 20,800 - 75,000 47,000

Outflow: 

Evaporation 

Springs .. 

Total ..

75,000 46,300

700

75,000 47,000

SIMULATED PUMPING PERIODS _s^_ 

2

PUMPAGE, IN THOUSANDS 

OF ACRE-FEET

ro tco -» 

0 C

3 OC

| Estimated a 

D Simulated a 

,_ (four pumj

NOTE: The firs 

pumping peri 

1910

........ Ji 111 III

nnual pum 

nnual pum 

aing perioc

t simulatec 

od began i

I II:

page 

page 

s)

j __. .

n

-

 :

=> ro .(*CO 01-^ Oo

SIMULATED PUMPAGE, IN 

CUBIC SECOND FEET PER

115° 113° 109° 108°

50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

WILLCOX BASIN (SHADED)

200

250

300

350

50

I 10°

150

200

250

50

100

150

200

250

1940

WELL A

WELLB

WELLC

1950 1960 1970

4,270

4,220

4,170

4,120

4,200

LLJ

4,150 >

4,100

4,050

4,050

1980

4,000 -

LU 

Q

3,950 t 

4,250 <

4,200

4,150

4,100

FIGURE 11. Location of basin, estimated and simulated water-budget components and pumpage, ground-water levels, model grid and boundary, and 

observed and simulated water levels, Willcox basin.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D27

110°00'

109°30'

Base from U.S. Geological Survey 

State base map, 1:500,000, 1974

5 MILES 

 

5 KILOMETERS

EXPLANATION

n

Alluvia] deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating groundwater flow

Cell representing

Mountain-front recharge

Evapotranspiration

Spring

-- 4140 - 

 4140 

 B

Geologic contact

Boundary of ground-water area

Boundary of modeled area 

Upper and lower layers

Water-level contour, spring 1975 

Shows altitude of water level. Contour 

interval 20 and 40 feet. Datum is sea 

level

Observed (from Mann and 

others, 1978)

Simulated 

Well for which hydrograph is shown

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Range 

Average

Specific yield of upper layer

Storage coefficient of lower layer

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

0.4 - 47 

8.0

400 - 54,980 

9,109

0.098

1C'4

10'7- 10"9

FIGURE 11. Continued.
D28 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

Childers (1967), and Mann and English (1980) provided 

the basis for model development.

Precipitation on the valley floor is less than 12 in./yr; 

however, the surrounding mountains receive more than 

20 in./yr, and parts of the Chiricahua Mountains to the 

east receive as much as 35 in./yr (University of Arizona, 

1965). Ground water moves from the basin margins to the 

basin axis and then toward the south end of the basin. 

The origin of ground water at the basin margins is surface 

runoff from the mountains that infiltrates the coarse sediments of the mountain fronts. Ground water discharges 

from the basin as base flow in Whitewater Draw and as 

underflow into Mexico. Before development, ground 

water also was discharged by evapotranspiration.

Ground water is withdrawn for irrigation, industrial, 

and public-supply uses in the Douglas basin. Withdrawal 

began about 1910 but was less than 5,000 acre-ft/yr until 

the close of World War II. At that time, about half the 

pumpage was for other than agricultural purposes. Beginning about 1950, irrigation became the major use of 

ground water. In 1974, pumpage was 138,000 acre-ft 

(fig. 12), and about 96 percent was used for irrigation 

(U.S. Geological Survey, 1982). Pumpage subsequently 

decreased to 80,000 acre-ft in 1978 (U.S. Geological Survey, 1982); however, the system continues to be overdrafted that is, more water is withdrawn than recharged 

back into the system.

MODEL CHARACTERISTICS

The numerical model developed for the Douglas basin 

had the following characteristics:

 Uniform row and column widths of 1 mi.

 Two layers representing (1) the upper part of the aquifer that was penetrated by wells and (2) the lower part 

of the system as yet impenetrated by wells (fig. 12).

 A vertical leakance connecting the two layers represents the vertical anisotropy in hydraulic conductivity 

within the aquifer system, which results hi head differential with depth.

 Perimeter recharge was applied to the upper layer of 

the model to simulate infiltration of runoff from surrounding mountains and underflow in the northeastern 

part of the basin, north of the Swisshelm Mountains.

 River cells simulate the discharge of ground water 

from the aquifer along the lower reach of Whitewater 

Draw.

 Discharge of ground water by evapotranspiration from 

the upper part of the aquifer in the area occurs where 

depth to water was less than 15 ft in 1910 (Meinzer 

and Kelton, 1913). The model uses a linear relation

between a maximum evapotranspiration rate and the 

depth at which evapotranspiration ceases. 

 Four periods of simulated pumping represent a generalized pattern of the actual pumping history; simulated 

pumping used an average rate applied over several 

years. Simulated irrigation pumpage was reduced by 5 

to 15 percent to represent a small quantity of recharge 

from excess applied irrigation.

RESULTS OF SIMULATIONS

Development of the steady-state model involved 

matching the model-calculated water-level configuration 

for the upper model layer with the measured water-level 

configuration for 1910 (Meinzer and Kelton, 1913). Initial 

values of horizontal and vertical conductivity were 

altered until the average absolute difference between 

model-calculated and measured heads reached an optimum low value of about 3 ft. The largest differences were 

near the basin margin where only a few water-level measurements were available. Less than 2 percent of the 

active cells of the upper layer indicated a head difference 

of more than 10 ft.

The validity of the steady-state solution was tested further by transient simulation involving four pumping periods: (1) predevelopment (1910) through 1946, (2) 1947 

through 1951, (3) 1952 through 1965, and (4) 1966 through 

1978. Information on water-level configurations was available for the winter periods of 1951, 1965, and 1978, and 

some water-budget information was available for 1951 

and 1965 for use in the history-matching process (fig. 12). 

The simulated water level for the 1951 pumping period 

matched well with observed water levels for 1951. However, because water-level declines at the south end of the 

basin may be affected by the occurrence of confined conditions locally in the aquifer system, simulated water levels for pumping periods ending in 1965 and 1978 (fig. 12) 

show more decline than was recorded in the northern and 

northwestern parts of the basin. This difference indicates 

that either hydrologic properties such as transmissivity, 

specific yield, and recharge values were too small; that 

actual pumpage was less than the quantity simulated; or 

that a combination of some of these factors occurred. 

More hydrologic data are needed to confirm these 

possibilities.

Initial values used in the Douglas basin model for rates 

of evapotranspiration, river leakance, and interlayer leakance were transferred from the model of upper San 

Pedro basin. Subsequently, only minor changes were 

needed in these values to reflect different basin geometry. 

The general relation between precipitation greater than 8 

in. and mountain-front recharge developed as part of this 

regional study (Anderson and others, 1992) was used to
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D29

determine the initial values for mountain-front recharge. 

Use of this relation provided adequate values, although no 

refinement or further definition could be accomplished 

because of the lack of detailed data about hydrologic 

properties and the exact shape of the water-level contours in the recharge areas.

TUCSON BASIN 

GEOHYDROLOGIC CONDITIONS

Tucson basin includes about 2,870 mi2 in southeastern 

Arizona and northern Mexico. The basin consists of about 

equal proportions of alluvial material and bedrock of the 

mountains. The length of the basin from the international 

boundary on the south to the outflow, about 20 mi northeast of Tucson, is about 85 mi, but only the northernmost 

65 mi was considered in this study (fig. 7). The numerical 

model included an area about 65 mi long and from 5 to 20 

mi wide (fig. 13). Average precipitation is about 12 in./yr 

on the basin floor and may exceed 30 in./yr in the highest 

parts of the surrounding mountains.

The basin occupies a structural trough that trends 

north to northwest. The lithology of the sediments that fill 

the trough are well defined (Pashley, 1966; Davidson, 

1973; Anderson, 1987; Anderson, 1988). Three units have 

been defined Pantano Formation, Tinaja beds, and Fort 

Lowell Formation that are nearly coincident with the 

pre-Basin and Range, lower basin fill, and upper basin fill, 

respectively. The Pantano Formation (Finnell, 1970) is at 

depth in the basin except where it outcrops on the basin 

perimeter. More than 8,000 ft of basin-fill sediments, 

which include the Tinaja beds, Fort Lowell Formation, 

and stream alluvium, are known to overlie the pre-Basin 

and Range rocks along the central axis of the basin 

(Anderson, 1987). The sedimentary units are hydraulically 

connected and form a single aquifer system in the basin.

Ground water generally is under unconfined conditions although head variations with depth, which reflect 

local ground-water confinement, are associated with the 

presence of fine-grained facies. Recharge occurs by infiltration of streamflow along the main stream channels and 

along mountain fronts. Before development, ground 

water moved from the areas of recharge toward the center and the downstream end of the basin, where it discharged by evapotranspiration or as surface flow in the 

Santa Cruz River. Small quantities of water also recharged 

and discharged from the aquifer system as underflow. 

After development of the ground-water resources, capture of discharge occurred gradually until only minor 

quantities of natural discharge remained. The principal 

discharge mechanism that has continued to occur is 

underflow at the downgradient end of the basin. As development increased, more ground water was withdrawn

from storage within the aquifer. The regional groundwater flow system has been disrupted by pumping, and 

most of the ground-water movement in the basin is now 

toward pumping wells and their associated cones of 

depression.

On the basis of archeological evidence, use of the available surface-water resources in the Tucson basin for irrigation of cropland predates the arrival of the Spanish in 

the 1500's (Harshbarger and others, 1966). Ground-water 

withdrawal began prior to 1900, but the quantity of 

ground-water withdrawal did not become significant until 

after 1920. Pumpage in 1940 was 62,000 acre-ft (Anderson, 

1972) and averaged about 280,000 acre-ft/yr in the 

mid-1970's (U.S. Geological Survey, 1985).

Water levels declined as much as 150 ft between 1940 

and 1978 near the center of two cones of depression in the 

basin (fig. 13). Rates of water-level decline were as great 

as 7 ftVyr in the mid-1970's (Murphy and Hedley, 1984) in 

the area of greatest withdrawals. Because of large-scale 

withdrawals and associated water-level decline, reaches 

of the Santa Cruz River that formerly contained perennial 

flow now contain flow only in direct response to runoff 

events.

MODEL CHARACTERISTICS

The numerical model of the Tucson basin simulated 

steady-state (predevelopment) and transient (pumping) 

conditions. The quasi-three-dimensional model was based 

on an earlier two-dimensional electrical-analog model 

(Anderson, 1972). The steady-state model was used to 

refine values and areal distribution of hydrologic properties and to evaluate the rates of ground-water movement 

within and between the two simulated layers of the aquifer. The model also was used to explore alternative conceptual ground-water models. The transient model was 

used to determine if the hydrologic properties used in the 

steady-state model realistically represented and simulated system response to pumping stress.

The steady-state and transient models developed for 

the Tucson basin have the following characteristics:

 A finite-difference grid with uniform cell size of 1 mi2 

except at the south end where an expanding row width 

was used because of sparse data (fig. 13).

 Two layers representing (1) the Fort Lowell Formation, 

which is equivalent to an upper basin-fill unit, and 

stream alluvium and (2) the Tinaja beds, which are a 

lower basin-fill equivalent, and part of the Pantano Formation, which are pre-Basin and Range rocks (Davidson, 1973).

 A vertical-leakance value was used to allow vertical 

ground-water flow between layers. Vertical anisotropy 

of hydraulic conductivity in each principal layer is used
D30 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

ESTIMATED AND SIMULATED WATER-BUDGET COMPONENTS 

[Values are in acre-feet per year; dashes indicate no data. Estimated values modified from 

Coates and Cushman (1955), White and Childers (1967), and unpublished records of the 

U.S. Geological Survey]

1951

Components

Inflow: 

Mountain-front recharge . . .

Storage ................

Outflow:

Pumpage ..............

Total ................

Estimated

20,000 

pnnn

1 30,000

co nnn

1 ,400 

1 ,400 

8,000- 13,000 

200

41 ,000

52,000 - 57,000

Simulated

15,500 

21,600

37,100

0 

1,700 

7,100

28,500

37,300

1965

Estimated

20,000 

270,000

90,000

0 

0 

0

90,000

90,000

Simulated

15,500 

44,300

59,800

0 

1,100 

2,000 

200 

57,300

60,600

^Calculated by using specific yield of 8 percent. 

Calculated by using specific yield of 20 percent.

SIMULATED PUMPING PERIODS

150 

O

z

100

z or 

 o

au- 50

| Estimated annual pumpage

Q Simulated annual pumpage

(four pumping periods)

NOTE: The first simulated 

pumping period began in 

1910

0 MM,,.jfe> _,,..:,.

200

150 S§

sec

100 ^ °-°-

3t

PLU 

<: Li50 ^O

115° 113 109° 108°

50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

DOUGLAS BASIN (SHADED)

50

50

Q 75

Z

UJ 

CO

IL 100 UJ 

LL.

125

75

100

125

150

1940

WELL A

3,950

WELLB

3,925

4,105

4,080

UJ

4,055 ai

4,030

WELLC

HI

Q

D

4,005 t

4,207 5

4,182

4,157

4,132

1950 1960 1970 1980

FIGURE 12. Location of basin, estimated and simulated water-budget components and pumpage, ground-water levels, model grid and boundary, and 

observed and simulated water levels, Douglas basin.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D31

109°45'

0 5 MILES 

0 5 KILOMETERS

05 MILES 

0 5 KILOMETERS

Base from U.S. Geological Survey 

State base map, 1:500,000, 1974

EXPLANATION

D

Alluvial deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating groundwater flow

Cell representing

Constant head 

Evapotranspiration 

Mountain-front recharge 

River-aquifer connection

 4050-

 4050-

 B

Geologic contact

Boundary of ground-water area

Boundary of modeled area 

Lower layer 

Upper layer

Water-level contour, 1978 

Shows altitude of water level. 

Contour interval 50 feet. Datum 

is sea level

Observed (Mann and English, 1980) 

Simulated

Well for which hydrograph is shown

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Range 

Average

Specific yield of upper layer

Storage coefficient of lower layer

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

1 -78 

7.8

86-1,996 

1,008

0.14

10-5

10-7-1Q-10

FIGURE 12. Continued.
D32 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

2,467

300

SIMULATED PUMPING PERIODS 

12345

- 400

w 

Q

Z 

< Wlrjiij O m

^ i-ii 

ztr 

 o

LU"< 

O LL 2°

^

ID 

Q.

200 -

| Estimated annual pumpage

D Simulated annual pumpage 

(five pumping periods)

100 -

115° 113' 109° 108°

0 50 MILES

0 50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

TUCSON BASIN (SHADED)

250

50

100

150

WELLC

- 2,417

2,227

2,037

1,987

1920 1930 1940 1950 1960 1970 1980

1,937

FIGURE 13. Location of basin, estimated and simulated pumpage, ground-water levels, model grid and boundary, and observed and simulated waterlevel declines, Tucson basin.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D33

COLUMNS

10 12 14 16 18 20 22 24 26 28

.10°45'

5 MILES ' 

5 KILOMETERS

0 5 MILES 

0 5 KILOMETERS

Base from U.S. Geological Survey 

State base map, 1:500,000, 1974

EXPLANATION

D

Alluvial deposits

Bedrock of the 

mountains

Grid of model cells 

used in simulating 

ground-water flow

Cell representing 

steady-state 

boundary conditions

Constant head 

Ground-water pumpage 

Mountain-front recharge 

Streamflow recharge

 Geologic contact

 Boundary of ground-water area

Boundary of modeled area 

 Lower layer

 Upper layer

Line of equal water-level decline,

1940-78 Interval 25 and 50 feet

 - - 50 - - - Observed

 50 Simulated

 B Well for which hydrograph is shown

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Range 

Average

Specific yield of upper layer

Storage coefficient of lower layer

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

3.4-716 

86.6

18-40,300 

5,444

0.133

10'4

5x1 0'4 - 4x1 0'6

FIGUEE 13. Continued.
D34 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

to calculate the harmonic mean of the vertical leakance at the start of the simulation. The leakance value 

does not change during the simulation but converts to 

zero if the upper layer is dewatered. 

The steady-state model has the following additional 

characteristics:

 Recharge to the lower layer along the mountain fronts 

and along stream channels traversing sediments of the 

Tinaja beds or Pantano Formation, except hi the area 

adjacent to the Tucson Mountains where no recharge is 

simulated.

 Recharge to the upper layer along the stream alluvium 

of the main stream channels.

 Simulated discharge from the lower layer that represents the quantity and area! distribution of pumping hi 

1940.

 Underflow simulated by constant heads hi both layers 

as inflow at the south end and northeast end of the 

basin.

 Underflow simulated by constant heads hi both layers 

as outflow near Rillito at the northwest end of the 

basin.

The transient model had these additional characteristics:

 A superposition approach was used so that only 

changes hi the ground-water system from steady-state 

conditions were simulated.

 A simulated surface of uniform altitude representing 

the head hi each of the two layers at the beginning of 

the transient period.

 Discharge from simulated pumping wells representing 

only the increase in quantity of pumpage since 1940.

 Recharge hi areas where the quantity of pumpage had 

been reduced since 1940 and where natural discharge 

was captured after 1940.

 Pumping only from the lower layer of the model.

Conditions in 1940 were assumed to represent steadystate conditions of inflow equal to outflow with no longterm change in storage (Anderson, 1972). Although some 

development had occurred in the northwestern part of the 

basin before 1940, conditions were assumed to have 

reached a new equilibrium condition by this time. The 

transient phase was simulated by five pumping periods 

1940^7, 1948-55, 1956-65, 1966-70, and 1971-78. The 

total water-level decline during 1940-78 (fig. 13) was used 

for matching water-level declines simulated by the model. 

For the superposition approach used hi the transient analysis, flow quantities used in the model consisted of 

changes hi flow components that occurred after 1940, and 

model output consisted of changes in water levels since 

1940.

RESULTS OF SIMULATIONS

The electrical-analog model developed by Anderson 

(1972) was used hi development of the quasi-three-dimensional numerical model of the Tucson basin. Because of 

differences hi computational procedures, more complex 

configuration, and time constraints, the simulated 1940 

water levels did not match measured water levels as 

closely as water levels generated by the electrical-analog 

model. The greatest difference between the observed and 

simulated water levels occurred near the basin margins 

where few water-level data are available (fig. 13). 

Improved matches between measured and simulated 

head values could have been achieved by selectively 

increasing recharge or decreasing transmissivity, depending on the range of values considered reasonable for these 

properties. Total recharge used hi the numerical model 

was 5,000 acre-ft/yr less than that used .hi the electricalanalog model, and increasing this flow component 

appears to be a viable option hi improving the match of 

water levels. Values of hydrologic properties used in the 

numerical model were about the same as those used in 

the electrical-analog model.

Water-level-decline maps and hydrographs of specific 

wells were used in the history-matching process. The distribution and magnitude of simulated head decline were 

consistently similar to actual measured declines. Simulated declines in the later pumping periods were less than 

measured declines, which indicated that simulated pumpage was too small, specific yield was too large, additional 

discharge was not being accounted for hi the model, or a 

combination of these factors. Model-generated declines 

are much more uniform than those based on field data 

because of the averaging and simplification of a highly 

variable, complex hydrologic system that are necessary in 

the model. A major limitation of the superposition 

approach used hi this basin model was the nonlinearity of 

the saturated thickness-transmissivity relation. Although 

water-level declines represented only a few to no more 

than 10-percent change in saturated thickness, changes in 

transmissivity of about 45 percent occurred in a few areas 

between 1940 and 1978. This large change is of special 

concern hi areas where highly transmissive surficial 

deposits have been dewatered. The superposition 

approach is not recommended for future use because of 

this limitation hi accuracy of simulation. An alternative 

approach might be to include more model layers and for 

each layer to represent a reasonably homogeneous segment of the aquifer.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D35

AVRA VALLEY 

GEOHYDROLOGIC CONDITIONS

Avra Valley occupies a north-trending structural basin 

in south-central Arizona (fig. 7). The basin includes about 

520 mi2, of which about 75 percent is basin-fill sediments, 

and 25 percent is bedrock of the mountains. Avra Valley is 

a broad, gently sloping plain that drains to the north and 

northwest. The valley is bounded on the east and west by 

mountains composed of nearly impermeable consolidated 

rocks. Average precipitation is about 10 in./yr on the valley floor and about 12 in./yr in the mountains (University 

of Arizona, 1965).

The sediments that fill Avra Valley are hydraulically 

connected and form a single aquifer system. Upper basin 

fill is as much as 700 ft thick in the northern part of the 

valley. Only about 400 ft of these sediments were saturated before development. Upper basin fill thins to less 

than 200 ft thick in the southern part of the basin. This 

part of the unit is above the water table and therefore is 

not part of the aquifer system. The lower basin fill 

includes extensive fine-grained material in the northern 

part of the valley and coarse-grained material in the 

south.

Before development, the ground-water system of Avra 

Valley consisted of a large storage reservoir having small 

quantities of recharge and discharge. Recharge occurs as 

underflow from the south; a small quantity of mountainfront recharge may enter the system from the east and 

west. Additional underflow enters the basin from the 

adjacent Tucson basin to the east Before development, 

the movement of ground water generally was from south 

to north paralleling the direction of surface drainage and 

discharge of ground water occurred as underflow to the 

northwest. Development has resulted in local cones of 

depression that have altered the original ground-water 

flow pattern.

The hydrologic system of Avra Valley is typical of those 

of the central group of basins (fig. 7) but does not have 

perennial streams or a large recharge from streamflow 

infiltration. Significant ground-water development began 

in the early 1940's for agricultural use. Withdrawals were 

estimated to be 12,000 acre-ft in 1940 (fig. 14) and 

159,000 acre-ft in 1975 (U.S. Geological Survey, 1985). 

Annual withdrawals decreased as agricultural lands were 

purchased and left unplanted by the City of Tucson, 

which purchased the land to obtain the rights to pump 

ground water underlying the property. Because withdrawal exceeded recharge, water levels declined as much 

as 150ft between 1940 and 1985 (Cuff and Anderson, 

1987). Maximum rates of water-level decline were about 

10 ft/yr in the early 1960's (fig. 14).

MODEL CHARACTERISTICS

The numerical model for Avra Valley was developed for 

steady-state and transient analysis and was based on previous work by Moosburner (1972). The model had the following characteristics:

 A finite-difference grid with uniform cell size of 1 mi2. 

The model grid consisted of 21 columns and 44 rows 

(fig. 14). The upper layer had 420 active cells, and the 

lower layer had 368 active cells.

 Two layers connected by a vertical-leakance value 

were used to simulate the aquifer system. In the northern part of the valley, the upper layer represents the 

upper basin fill, and the lower layer represents the 

lower basin fill. In the south and central parts of the 

valley, the break between model layers is arbitrary, and 

both layers represent lower basin fill. The verticalleakance value was calculated by the model as a function of the ratio of vertical to horizontal hydraulic conductivity for each layer.

 Perimeter recharge was considered negligible in the 

steady-state and transient simulations, especially in 

contrast to the large quantity of pumpage and extensive depletion of water in storage during 1940-77.

 The steady-state model used constant-head boundaries 

for both layers in areas of inflow and outflow. The transient model used constant-flow boundaries at the same 

locations. The quantity of flow in the transient model 

was equal to the flow that was simulated in the steadystate (predevelopment) model.

 The pumping history was simulated in six pumping 

periods 1940-50, 1951-55,1956-60, 1961-64,1965-72, 

and 1973-77 (fig. 14).

RESULTS OF SIMULATIONS

The steady-state model represents a conversion of the 

two-dimensional electrical-analog model of Moosburner 

(1972) to a quasi-three-dimensional numerical model. 

Water-level altitudes from the steady-state simulation generally matched measured water-level altitudes within 5 ft 

(fig 14). Larger differences were noted in the southernmost part of the area where data are sparse.

Measured water levels for 1964 and 1977 were used for 

the history-matching process. Simulated water levels at 

the end of 1964 compared closely with water levels of 

Moosburner (1972) and measured water levels. Simulated 

water levels for the end of 1977 resulted in a poor match 

with measured water levels, and significantly greater 

magnitude and rate of drawdown were simulated by the 

model in the north-central part of the area. These results 

indicated that the model could be improved by an 

increase in recharge, an increase in specific yield, a 

decrease in pumping, or some combination of these
D36 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

150

200

SIMULATED PUMPING PERIODS

j^_ 

3

05 h- 

=>UJ

zcr 100 

 o 

ul< 

(5 L-

| Estimated annual pumpage

D Simulated annual pumpage 

(six pumping periods)

uJM

250

si

200 ulg

£05 

20C 150 =)UJ

Q. DD IujW 

trw 100 5"- 

^

*%

50 W 0

115 109° 108°

50 MILES

_l

50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

AVRA VALLEY (SHADED)

250

350

150

Simulated WELL A

Measured With 

recharge

No recharge**

1,755

1,655

250

350

250

WELLB

Measured

\

With 

V «, \, recharge \-*-^ -I

No 

recharge'

1,555

1,819

1,719

350 -

450

1,619

1,940

- 1,840

1940 1950 1960 1970 1980

1,740

FIGURE 14. Location of basin, estimated and simulated pumpage, measured and simulated ground-water levels, model grid and boundary, and 

observed and simulated water-level declines, Avra Valley.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D37

0 5 MILES

k . 0 5 KILOMETERS

Base from U.S. Geological Survey 

State base map, 1:500.000,1974

5 KILOMETERS

EXPLANATION

Alluvial deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating groundwater flow

Cell representing

Inflow

Outflow 

Geologic contact

 50---

-50- 

 B

Boundary of ground-water area

Boundary of modeled area 

Upper and lower layers

Line of equal water-level decline, 

1940-78 Interval 25 and 50 feet

Observed (based on unpublished 

data from City of Tucson, 

University of Arizona, and 

U.S. Geological Survey)

Simulated 

Well for which hydrograph is shown

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Range 

Average

Specific yield of upper layer

Storage coefficient of lower layer

Vertical leakance between layers, feet 

per second per foot 

Range

IN MODEL

VALUE

0.6-110 

26.6

66 - 4,360 

1,530

0.131

io-4

9x10'4 -4x10'6

FIGURE 14. Continued,
D38 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

factors. The extent of the area in which matches between 

simulated and measured water levels were poor was similar to the known extent of perched ground water. The correspondence in the two areas probably indicates that 

recharge is entering the regional ground-water system as 

a result of drainage from the perched ground water. The 

perched ground water probably originates, at least in part, 

from excess applied irrigation water.

The model was further used to develop quantitative 

information on the magnitude and rate of the possible 

recharge. In the area in which recharge was indicated, 

depths to water were 200 to 300 ft in the early to mid1960's, which is the period when simulated and measured 

water levels began to diverge (hydrographs for wells A 

and B, fig. 14). In the area to the south where depths to 

water were greater than 350 ft, the measured and simulated water levels did not diverge (hydrograph for well C, 

fig. 14). Where recharge from deep percolation of excess 

applied irrigation water apparently occurred, the volumetric difference caused by the divergence of measured 

and simulated water levels for 1940-77 was about 500,000 

acre-ft. A transit time of 15 to 20 years was indicated for 

the excess applied irrigation water to move through the 

250-ft-thick unsaturated zone to reach the saturated zone 

(Anderson, 1983, p. 38).

WATERMAN WASH BASIN 

GEOHYDROLOGIC CONDITIONS

Waterman Wash basin is in south-central Arizona about 

25 mi southwest of Phoenix (fig. 7) and is typical of the 

basins of the west group. The basin is drained by Waterman Wash, an ephemeral stream tributary to the Gila 

River. The basin includes about 470 mi2, of which about 75 

percent is alluvial sediments, and 25 percent is nearly 

impermeable bedrock of the mountains. The alluvial 

basin is about 25 mi long and 15 mi wide. Average precipitation ranges from less than 7 in./yr on the valley floor to 

more than 10 in./yr in the mountains (University of 

Arizona, 1965).

The sediments that fill Waterman Wash basin include 

both upper and lower basin-fill units. Pre-Basin and 

Range rocks may occur at depth in the basin; however, 

data are not available to verify this speculation. The lithologic units are interconnected hydraulically and form a 

single aquifer system. Ground water occurs under unconfined conditions in the basin.

Before development, ground-water movement in 

Waterman Wash basin was minimal because recharge to 

the ground-water system was small. Ground water moved 

from the southeast to the northwest, where it discharged 

as underflow and possibly by evapotranspiration.

Natural recharge to the aquifer of Waterman Wash 

basin has been estimated to be 1,500 to 2,500 acre-ft/yr 

(Halpenny and others, 1952, p. 155). Withdrawal of ground 

water in Waterman Wash basin began as early as 1940 but 

probably did not exceed natural recharge until about 

1950. During the early 1950's, the amount of cultivated 

acreage increased from virtually nothing to about 15,000 

acres, and ground-water withdrawal increased from less 

than 5,000 acre-ft in 1950 to 60,000 acre-ft in 1960 (fig. 15). 

Since 1960, the quantity of pumpage has fluctuated 

between 45,000 and 72,000 acre-ft/yr (U.S. Geological Survey, 1982). Water levels have declined more than 170 ft in 

the area of greatest withdrawal.

MODEL CHARACTERISTICS

The numerical model developed to simulate the hydrologic system of the Waterman Wash basin was based on a 

previous model by Matlock (1981). Modification of Matlock's single-layer model to a multilayered model enabled 

analysis of vertical components of ground-water flow. A 

steady-state simulation was developed to assess the 

previous concepts of the natural hydrologic conditions, 

the significance of the quantity and distribution of natural 

recharge, and the significance of interlayer vertical flow 

under steady-state conditions. A transient simulation 

using two periods of pumping 1951-60 and 1961-80 

was used to analyze the significance of hydrologic properties during the development period (fig. 15).

The characteristics of the model were as follows:

 A finite-difference grid with uniform cell spacing of 1 

mi in both horizontal directions. The model grid consisted of 14 columns and 27 rows, of which 209 cells 

were active.

 Two layers representing (1) the upper basin fill, which 

consists of poorly sorted gravel, sand, silt, and clay, 

and (2) the lower basin fill, which consists of moderately indurated sand and gravel and a small percentage 

of fine-grained material.

 A vertical-leakance value connected the two layers to 

represent the vertical anisotropy within the aquifer 

system.

 Discharge as underflow was simulated in the area 

where the channel of Waterman Wash exits the basin. 

Recharge as underflow was simulated at two cells at 

the south end of the model array.

 Simulated pumpage was from the upper layer only.

RESULTS OF SIMULATIONS

Although the quasi-three-dimensional steady-state 

model developed as part of this study was nearly identical 

with the two-dimensional model developed by Matlock
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D39

(1981), one significant difference was the simulation of 

natural recharge. In Matlock's model, natural recharge of 

2,000 acre-ft/yr was distributed around the model perimeter; in this model, 2,000 acre-ft/yr was simulated as entering through two cells at the south end of the area. The 

differences in water-level configuration are minor and 

indicate that the effect of natural recharge on the system 

is minor and that the spatial distribution of recharge is not 

an important model consideration. This situation probably reflects the combination of a small quantity of 

recharge entering part of the aquifer having a large 

hydraulic conductivity. Also, few wells are available in 

which to measure water levels in the mountain-front area. 

Thus, a small change in hydraulic gradient, if it occurs, 

cannot be well defined. These factors result in water-level 

contours being interpreted as generally straight lines 

across the basin and nearly perpendicular to the mountain fronts.

The transient model in this study used 1980 as the final 

stress period, whereas Matlock's model used 1975 for the 

final stress period. Transient results generally were satisfactory, although the results of the first pumping period 

do not simulate actual water-level conditions as well as 

those of the second period. The conclusion was that additional time increments of pumping should have been used 

to better simulate the temporal and spatial changes in 

pumping stress.

HARQUAHALA PLAIN 

GEOHYDROLOGIC CONDITIONS

Harquahala Plain is a northwest-trending structural 

basin in west-central Arizona about 70 mi west of Phoenix 

and is typical of the west group of basins (fig. 7). The 

basin includes about 760 mi2 , of which about 30 percent is 

consolidated rocks, and 70 percent is alluvial material. 

The area is arid, and average precipitation is about 7.5 in./ 

yr (University of Arizona, 1965).

Basin fill in Harquahala Plain consists of upper and 

lower units. Pre-Basin and Range rocks are presumed to 

occur at depth in the basin but have not been penetrated 

by wells. The upper basin-fill unit includes a heterogenous mixture of coarse- and fine-grained sediments in the 

southeastern part of the basin and an areally extensive 

fine-grained facies that increases in thickness from the 

southeast toward the northwest. On the basis of well logs, 

the fine-grained facies is more than 1,300 ft thick near the 

center of the basin. The lower basin-fill unit includes 

fine-grained material near the basin center and unconsolidated silt, sand, and gravel and conglomerate at the southeast end of the basin. This material, plus the coarser 

material of the upper basin-fill unit, composes the principal water-bearing unit in the basin. All the sediments are

interconnected hydraulically and form a single aquifer 

system.

Before development, ground-water movement in the 

aquifer was from the northwest to the southeast (fig. 16); 

the water-level gradient throughout much of the central 

and southeastern parts of the basin was about 2 ft/mi. 

Natural recharge to and discharge from the aquifer system before development were small. Metzger (1957, p. 32) 

postulated that 90 percent of natural recharge was from 

infiltrating streamflow and that the total was less than "a 

few thousand acre-feet annually." Some water enters at 

the northwest end of the basin as underflow from 

McMullen Valley. Before development, discharge 

occurred as either underflow or as a combination of 

underflow and evapotranspiration at the southeast end of 

the basin.

Irrigation wells were first drilled in Harquahala Plain in 

the late 1930's, but pumpage was small until the early 

1950's when the first large-capacity wells were developed. 

Water-level declines of 3 to 4 ft/yr occurred during the 

mid-1950's. Withdrawals increased markedly from 30,000 

acre-ft in 1955, to 125,000 acre-ft in 1960, and to 200,000 

acre-ft/yr during 1962-65 (fig. 16). Since 1965, pumpage 

has diminished steadily to about 100,000 acre-ft/yr during 

1978-80 (fig. 16). During 1957-63, water-level declines 

ranged from 60 to 160 ft in areas of ground-water withdrawals. By 1980, total water-level decline exceeded 

300 ft in much of the southeastern part of the basin (Graf, 

1980).

MODEL CHARACTERISTICS

The numerical model for Harquahala Plain included 

both steady-state and transient analyses and had the following characteristics:

 A variable-sized finite-difference grid was used. The 

model grid consists of 38 rows and 17 columns and 

includes 298 active cells in each of two layers.

 Two layers were used to represent (1) upper basin fill, 

which consists of a mixture of fine- and coarse-grained 

sediments at the southeast end and fine-grained sediments at the center of the basin, and (2) lower basin 

fill, which includes conglomerate and silt, sand, and 

gravel at the southeast end and fine-grained sediments 

in the basin center.

 Interconnection between the two layers was simulated 

by a model-calculated vertical-leakance value, which 

was calculated on the basis of the harmonic mean of 

the vertical hydraulic conductivities of the units 

divided by one-half the thickness of the units.

 No mountain-front recharge was simulated. For the 

steady-state simulation, only underflow into the basin 

from the northwest and underflow out of the basin at
D40 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

SIMULATED PUMPING PERIODS 250

75

CO

Q

z 

<

CO I- cn

3 LU 50

 ^LU 

Z CC ~O 

LU"<

O Li1

25

| Estimated annual pumpage

Q Simulated annual I

pumpage (two iJ 

pumping periods) III

100

115° 113 109° 108°

0 50 MILES 

0 50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

WATERMAN WASH BASIN (SHADED)

350

300

400

200

3°°

400

100

200

300

400

WELL A

932

WELLB

WELLC

WELLD

832

936

836

877

777

677

890

790

690

590

1940 1950 1960 1970 1980

FIGURE 15. Location of basin, estimated and simulated pumpage, ground-water levels, model grid and boundary, and observed and simulated water 

levels, Waterman Wash basin.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D41

0 | 5 MILES 

5 KILOMETERS

Base from U.S. Geological Survey 

State base map, 1:500,000,1974

5 MILES

5 KILOMETERS

EXPLANATION

E3

Alluvial deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating groundwater flow

Cell representing

Inflow

Outflow 

Geologic contact

---750--

 750 

 B

Boundary of ground-water area

Boundary of modeled area

Upper and lower layers

Water-level contour, 1980-81 

Shows altitude of water level. 

Contour interval 25 feet. Datum 

is sea level

Observed (R.S. Stulik, Arizona 

Department of Water Resources, 

written commun., 1981)

Simulated 

Well for which h'/drograph is shown

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Single uniform value

Specific yield of upper layer

Storage coefficient of lower layer

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

8-20 

14.6

7,000

0.117

1C'4

10'4 -4x10-5

FIGURE 15. Continued.
D42

SIMULATED PUMPING PERIODS

COh3UJ 

O LLi

ZOC o

O u- 

<O 

D.

£3U

200

150

100

50

n

| Estimated annual pumpage (based on 

power usage)

D Estimated annual pumpage (be 

^ irrigated acreage)

. D Simulated annual

pumpage (eight 

pumping periods)

"t

t J

ll

...... if 1

ased on

In - 1 In (HI ,

; --

1

...lilllL,

300

200

0. Q100 5^5 

§

CO O

115° 113 109° 108°

0 50 MILES

0 50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

HARQUAHALA PLAIN (SHADED)

200

300

400

< 500

600

200

300

400

500

600

WELL A

Measured

Simulated 

variable 

^/ specific

yield 

k \ 

t \ 

\ 

A

Simulated

small specific

yield

WELL B

Measured

Simulated 

variable 

» /'specific 

/ yield

Simulated

small specific

yield

1,000

900

800

700

600

1,006

906

806

706

606

1940 1950 1960 1970 1980

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Range 

Average

Specific yield of upper layer

Storage coefficient of lower layer

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

0.1 - 43 

2.2

1,000-41,000 

13,900

0.05-0.15

io-4

10-6-10-7

FIGURE 16. Location of basin, estimated and simulated pumpage, ground-water levels, model grid and boundary, and observed and simulated water 

levels, Harquahala Plain.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D43

Base from U.S. Geological Survey 

State base map, 1:500,000, 1974

EXPLANATION

Alluvial deposits

Bedrock of the 

mountains

2 4

Mill Grid showing model 

2 MM! cells used in 

simulating 

ground-water 

flow 

FIGURE 16. Continued.

Cell representing

Inflow 

Outflow 

Geologic contact

Boundary of groundwater area

Boundary of modeled 

area Upper and lower 

layers

 750 

 750 

Water-level contour, 1980 

Shows altitude of water level. 

Contour interval 20 and 50 feet. 

Datum is sea level

Observed (Graf, 1980) 

Simulated

' B Well for which hydrograph is shown
D44 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

the southeast end were included. These flows were 

simulated by use of constant-head cells. 

 All simulated pumping was from the lower layer. The 

pumping history was simulated in eight pumping 

periods 1953-56,1957-58,1959-60,1961-65,1966-69, 

1970-73,1974-76, and 1977-SO.

RESULTS OF SIMULATIONS

Few data are available on the magnitude and distribution of hydrologic properties of the sediments and predevelopment water-level conditions; thus, the developed 

steady-state model is of uncertain reliability. The steadystate simulation generally was insensitive to change in 

recharge location and quantity within a range thought to 

be practical for the basin. On the basis of the steady-state 

model results, natural recharge to and discharge from the 

basin were estimated to be about 1,200 acre-ft/yr.

Transient-model results indicate that the basin is a typical storage-depletion system in which pumpage is balanced almost entirely by removal of water from storage. 

The small rate of natural recharge had no discernible 

effect on the rate of water-level decline and the mining of 

ground water. Under these conditions, specific yield is a 

major controlling property in the simulation but is the 

poorest known. The effect of poorly known specific yield 

on model results was evaluated by making two transientmodel runs using different distributions and values of specific yield. The values were assigned on the basis of 

extent of fine-grained sediments and were based on values determined by indirect means for other basins in the 

study area. The first model run used relatively small 

specific-yield values of 0.05 and 0.10, which corresponded 

to areas of fine-grained facies and to areas of coarser, heterogeneous sediments, respectively. The model response 

at two locations is shown as "small specific yield" on 

hydrographs for wells A and B in figure 16. A second 

model run used specific-yield values of 0.05, 0.10, 0.13, 

and 0.15 to represent various areas of decreased content 

of fine-grained sediments. The model response for this 

scenario is shown as "variable specific yield" on hydrographs for wells A and B in figure 16. The model-generated response indicates that use of the variable specificyield values provided good comparisons of water-leveldecline rates during 1950-65. The divergence of field and 

model data starting about 1965 cannot be explained by 

altering the specific yield and is thought to indicate 

recharge from excess applied irrigation water, although 

other factors may also be involved.

MCMULLEN VALLEY 

GEOHYDROLOGIC CONDITIONS

McMullen Valley is a roughly rectangular, northeasttrending basin in west-central Arizona (fig. 17). The basin 

includes about 720 mi2, of which about 25 percent is 

nearly impermeable consolidated rocks of the mountains, 

and the remaining 75 percent is permeable alluvial material. The climate is arid, and average annual precipitation 

is less than 10 in. (Sellers and Hill, 1974).

Basin fill, which is the major aquifer in the valley, consists of upper and lower units. The upper basin fill is 

unsaturated in most of the basin. The small part of upper 

basin fill that was saturated before development subsequently has been dewatered. In the central part of the 

basin, lower basin fill includes coarse-grained sediments 

overlain by as much as 1,000 ft of fine-grained sediments. 

The lithology of the basin was documented by Pool (1987, 

p. 14-18).

Ground-water movement through the basin was from 

northeast to southwest before development. Natural 

recharge, although small, occurs in the northern part of 

the basin and possibly at the base of the peripheral mountains as mountain-front recharge. Discharge occurred as 

underflow into Harquahala Plain and, in part, possibly as 

evapotranspiration.

Ground-water development altered the natural 

ground-water flow system after 1952. During 1954-59, 50 

irrigation wells were completed in the valley (Kam, 1964). 

Irrigated area increased from 11,000 acres in 1958 to 

16,600 acres in 1965 (Briggs, 1969). By 1980, irrigated 

acreage had increased to 34,200 acres (Remick, 1981). As 

shown in figure 17, pumpage increased from 5,000 acre-ft 

in 1952 to 139,000 acre-ft in 1975 and has decreased 

slightly since 1975 (U.S. Geological Survey, 1985). Maximum water-level declines of more than 250 and 150 ft 

occurred in the southwest and northeast parts, respectively, of the basin. Rates of decline ranged from 5 to as 

much as 15 ft/yr during the early 1970's (hydrograph for 

well A, fig. 17) in the principal cones of depression.

MODEL CHARACTERISTICS

The numerical model for McMullen Valley included 

steady-state and transient analyses and were documented 

by Pool (1987). The steady-state analysis was used to 

evaluate areal distribution and magnitude of aquifer properties and predevelopment recharge and discharge. The 

transient simulation was used to analyze the relation of 

pumpage and water-level decline. The numerical model 

had the following characteristics: 

 A finite-difference grid of variable cell size was

designed to allow more detailed evaluation in areas of
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D45

steep water-level gradients and known geohydrologic 

boundaries. The model grid consisted of 37 rows and 

29 columns.

 Two layers were used to represent (1) upper and lower 

basin fill of which the lower basin fill is the principal 

aquifer in most of the basin and (2) pre-Basin and 

Range rocks. The layers simulate different area! 

extents because of the physical and hydrologic boundaries of the aquifer.

 Hydraulic connection between the two layers was simulated by a vertical-leakance value.

 Natural inflow and outflow were simulated as constant-head cells during steady-state calibration. The 

flow quantity from the constant-head cells simulated in 

the steady-state model was transferred to the transient 

simulation and held constant.

 Pumping was simulated in four periods 1951-57, 

1958-65,1966-72, and 1973-80.

RESULTS or SIMULATIONS

Water levels in the steady-state model were matched 

with the pre-1952 water-level configuration (Kam, 1964), 

although definition of the predevelopment system was 

poor because of the lack of data. The transient analysis 

included the period 1952-80. The history-matching procedure resulted in a good comparison of simulated and 

field-measured water-level contours, with little change in 

the general distribution of hydrologic properties. Transient simulation required adjustments to vertical leakance, specific yield, and the location of an impermeable 

boundary as documented by Pool (1987), all of which 

were undefined or poorly defined by available field data 

The northeastern part of the basin was simulated as a typical storage-depletion system. Stratigraphic complexities 

and the presence of buried consolidated rocks along the 

south boundary made simulation of the southwestern part 

of the system difficult. A fine-grained facies that inhibits 

vertical movement of ground water also affected the size 

and shape of the simulated cone of depression. Simulations indicate that development in McMullen Valley has 

affected water levels throughout the valley. In some areas, 

water-level declines have been accelerated because the 

cones of depression have intersected impermeable 

boundaiies. Model results also indicated that vertical leakance increases in the northwestern part of the valley; the 

fine-grained facies may be absent in this part.

PARKER VALLEY 

GEOHYDROLOGIC CONDITIONS

Parker Valley includes about 230 mi2 of the Colorado 

River flood plain in western Arizona and is typical of the 

basins in the Colorado River group (fig. 7). The climate is 

arid, and average annual precipitation on the flood plain 

is less than 4 in. (Hely and Peck, 1964, pi. 3).

In Parker Valley, the sediments include the Miocene(?) 

fanglomerate (Metzger and others, 1973); the overlying 

Bouse Formation, which is equivalent to upper basin fill; 

older alluvial deposits, which are almost completely 

eroded; and stream alluvium, which makes up the flood 

plain of the Colorado River and fills channels cut into the 

Bouse Formation and older alluvium (Metzger and others, 

1973; Tucci, 1982). The stream alluvium and at least the 

upper part of the Bouse Formation are connected hydraulically and constitute the principal aquifer in Parker Valley. Ground water occurs under unconfined conditions in 

the aquifer.

The hydrologic system of Parker Valley is dominated 

by the Colorado River, which controls the quantity of 

inflow to and outflow from the aquifer adjacent to the 

river. In response to extensive agricultural development, 

diversion and application of surface water to irrigated 

cropland have resulted in shallower depths to ground 

water because of recharge from excess applied irrigation 

water. Discharge from the ground-water system occurs as 

evapotranspiration from crops and riparian vegetation, as 

leakage to drains used to prevent waterlogging, and as 

discharge to the Colorado River.

MODEL CHARACTERISTICS

The numerical model for Parker Valley was used to 

evaluate present knowledge and concepts of the 

ground-water system and to simulate the area! system 

(Tucci, 1982). Model results were initially assumed to be 

transferable to other similar basins along the lower Colorado River; this assumption was proved to be only partly 

valid on the basis of a model of the Mohave basin, which 

is about 30 mi north of Parker Valley. The numerical 

model for Parker Valley had the following characteristics:

 Variable spacing was used in the finite-difference grid 

to allow increased resolution in simulation of the 

extent and alignment of the Colorado River.

 Two layers were used to represent (1) stream alluvium 

and (2) the upper zone of the Bouse Formation 

(Metzger, 1968).

 Flow between layers was controlled by a leakance 

value to provide a vertical connection and was calculated by the model on the basis of the vertical hydraulic 

conductivity and thickness of the model layers.
D46 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

SIMULATED PUMPING PERIODS

IOU

in

Q

z

ID LU 100

z tr 

 o

UJ< 

O LL

<0 50

D.

^

D.

n

| Estimated annual pumpage

Q Simulated annual pumpage 

(four pumping periods)

, ......Liillll

'- 200

ll

150 LUQ

^> C/5 

ID LU

__J Q

50 I §

J. 1 1 L U 1 1 ! L A

109° 108°

0 50 MILES

0 50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

McMULLEN VALLEY (SHADED)

50

150

250

350

150

lu 250

LLJ

350

450

350

450

550

1940

WELL A

WELLB

WELLC

1950 1960 1970

1,816

1,716

1,616

LJJ 

> 

LIJ 1,516 -1 

<

1,776 UJ w

LIJ

> 

O

CQ 

<

1,676 tii

LIJ

1,576

1,476

1,849

1,749

1980

1,649

FIGURE 17. Location of basin, estimated and simulated pumpage, ground-water levels, model grid and boundary, and observed and simulated water 

levels, McMullen Valley.
Alluvial deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating groundwater flow

Cell representing

Inflow 

Outflow 

Geologic contact

Boundary of groundwater area

SIMULATION OF GROUND- WATER FLOW IN ALLUVIAL BASINS

113°30'

D47

0 5 MILES 

0 5 KILOMETERS

EXPLANATION

Boundary of modeled area 

Lower layer 

Upper layer

Water-level contour, December 

1980 (Pool, 1987) Shows altitude 

of water level. Contour interval 

50 feet. Datum is sea level

 -1500 - - Observed

-7500 Simulated

 B Well for which hydrograph is shown

HYDRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer, in 

feet per day 

Range 

Average

Transmissivity of lower layer, in feet 

squared per day 

Range 

Average

Specific yield of upper layer

Storage coefficient of lower layer

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

1 -20 

3.5

100-30,000 

6,450

0.10

10-4

10'e -10-11

FIGURE 17. Continued.
D48 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

 Uniform hydraulic conductivity was used throughout 

the upper layer.

 Uniform transmissivity was used throughout the lower 

layer.

 River-aquifer connection was simulated in the upper 

layer as leakage through a confining layer. Leakage to 

drains was similarly simulated (Tucci, 1982).

 Evapotranspiration by native vegetation was simulated

in the upper layer, and the maximum effective depth of

evapotranspiration was simulated as 20 ft.

The finite-difference grid consisted of 30 rows and 52

columns (fig. 18). A row of constant-head cells was used

at the downgradient end of the model, and constant-flow

cells were used at areas of tributary underflow into the

basin. River and drain cells were treated similarly by

using an altitude for the water level in the river or drain

and a leakance value that controls the interconnection of

the surface-water and ground-water systems.

RESULTS OF SIMULATIONS

Steady-state model simulations were made for 

ground-water conditions in 1940^1, the mid-1960's, and 

1980 (Tucci, 1982, p. 38). In 1940^1, the Colorado River 

was a losing stream throughout Parker Valley and represented the major source of recharge to the aquifer; evapotranspiration by native vegetation on the flood plain was 

the dominant discharge mechanism. Model results for 

conditions in the mid-1960's indicate that the Colorado 

River was a gaining stream in the northern part of the valley as a result of ground-water-level rises caused by infiltration of applied irrigation water. Recharge to the aquifer 

from infiltration of surface water was decreased as a 

result of (1) shallower ground-water levels and (2) 

smaller gradients away from the river along reaches 

where it was still a losing stream. Model results for conditions in 1980 indicated a further decrease in surface-water 

infiltration associated with continued ground-water level 

rises.

Model results generally compared favorably with available field data. The model was most sensitive to changes 

in the average evapotranspiration rate. The model was 

less sensitive to changes in river- and drain-leakance values. The hydrologic system in Parker Valley can be adequately simulated by a simple two-dimensional approach 

because of the magnitude and dominance of the nearsurface hydrologic interactions.

MOHAVE BASIN 

GEOHYDROLOGIC CONDITIONS

Mohave basin is along the flood plain of the Colorado 

River in western Arizona (fig. 7). The modeled area 

includes 103 mi2 near the center of the basin (fig. 19). The 

climate is arid, and average annual precipitation is less 

than 6 in. (Hely and Peck, 1964, pi. 3).

The geohydrologic setting of Mohave basin is similar to 

that of Parker Valley and other basins transected by the 

lower Colorado River. The aquifer is composed of stream 

alluvium of the Colorado River and upper basin fill, which 

includes older alluvium and the underlying upper zone of 

the Bouse Formation (Metzger, 1968; Metzger and Loeltz, 

1973). As in Parker Valley, the stream alluvium occupies 

channels cut into the older alluvium and the Bouse Formation. The older alluvium occurs only laterally adjacent 

to the stream alluvium at the edge of the present flood 

plain of the Colorado River. The ground-water system of 

Mohave basin is dominated by infiltration of surface 

water from the Colorado River, which is the principal 

source of recharge. Evapotranspiration by crops and 

riparian vegetation is the principal discharge.

MODEL CHARACTERISTICS

The numerical model for Mohave basin was developed 

to test the transferability of modeling concepts and hydrologic properties from the model of Parker Valley. The 

model results were matched to steady-state conditions in 

the 1960's and had the following characteristics:

 Variable spacing was used in the finite-difference grid 

to provide additional resolution in simulation of the 

extent and alignment of the Colorado River.

 Two layers were used to represent (1) stream alluvium 

and the adjacent older alluvium and (2) the upper zone 

of the Bouse Formation (Metzger, 1968).

 Flow between layers was controlled by a leakance 

value to provide a vertical connection and was calculated by the model on the .basis of vertical hydraulic 

conductivity and thickness of the model layers.

 Uniform hydraulic conductivity was used throughout 

the upper layer.

 Uniform transmissivity was used throughout the lower 

layer.

 The river-aquifer connection was simulated in the 

upper layer as leakage through a confining layer as in 

the model for Parker Valley.

 Constant-head cells were used to simulate groundwater flow along the southeast edge of the model and 

to simulate the open-water surface of Topock Marsh 

adjacent to the Colorado River in the southern part of 

the area (fig. 19).
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D49

 Evapotranspiration by native vegetation was simulated 

in the upper layer. The maximum effective depth of 

evapotranspiration was simulated as 20 ft. 

The finite-difference grid consisted of 52 rows and 23 

columns (fig. 19). The outer boundary of the model consisted of no-flow cells except in the southeast corner of 

the model, where tributary ground-water inflow along 

Sacramento Wash was simulated by constant-head cells. 

As an initial approach to model development, hydraulic 

conductivity of the upper layer, transmissivity of the 

lower layer, river leakance, and the evapotranspiration 

and depth-to-water relation were transferred directly 

from the calibrated model for Parker Valley.

RESULTS OF SIMULATIONS

Transfer of hydrologic-property values from the model 

of Parker Valley to the model of Mohave basin was only 

partially successful. Although the geohydrologic settings 

are similar in both basins, sufficient differences exist that 

make direct transfer of some property values inappropriate. The most significant difference was the need to 

decrease hydraulic conductivity of the upper unit where 

the older alluvium occurred laterally adjacent to the 

stream alluvium. The older alluvium is considered to be 

much finer grained than the stream alluvium, on the basis 

of field observations and analysis of drillers' logs. The 

hydraulic conductivity was decreased from 311 to 6 ft/d in 

areas underlain by the older alluvium outside the flood 

plain of the Colorado River; this decrease resulted in 

much improved simulated water-level gradients. Hydraulic conductivity of the stream alluvium was maintained at 

311 ft/d as used in the Parker Valley model. Also, transmissivity of the lower model unit equivalent to the 

upper part of the Bouse Formation was decreased by 

about one-third to adjust for the unit being thinner and 

having a higher silt and clay content than in Parker Valley. 

The only other modification to model properties was in 

the relation between evapotranspiration and depth to 

water. On the basis of known differences in type and density of riparian vegetation, the maximum evapotranspiration rate was decreased about 30 percent.

A water budget and a map of head distribution were 

used as means of evaluating the degree to which the 

model simulated the actual hydrologic system. Individual 

water-budget components based on field estimates were 

compared with model-calculated components (fig. 19). 

The small imbalance in the budget is because of intrinsic 

inaccuracies in estimating some of the individual budget 

values. Within the flood plain, where water levels were 

best known, model-calculated water levels were within 10 

ft of measured levels. Information transfer from a model 

of a hydrologically similar basin provided initial estimates

of aquifer properties and thus decreased the time necessary to calibrate the numerical model. The entire process 

served to prove the adaptability of information transfer 

developing the conceptual model and the numerical 

model to basins of this group.

BASIN-GROUP MODELS

Two basin-group models one representing the southeast basins and the other representing the central and 

west basins were used to explore model sensitivity to 

use of generalized values of geohydrologic properties 

(table 1). Models that used average values were compared 

with models that used various alternative combinations of 

boundary, hydrologic, and flow properties.

Because of the ready transferability of property values 

from the Parker Valley model to the Mohave basin model, 

no further basin-group modeling of the Colorado River 

group was undertaken. In addition, owing to the paucity 

of data, no models were developed for the basins of the 

highland group. These basins are thought to have distinct 

similarities to the Colorado River group in which coarse, 

unconsolidated sediments make up the shallow alluvial 

aquifer, and the hydrology is dominated by the presence 

of perennial surface flow. Significant differences do exist, 

however, in size of the basins and climatic characteristics.

SOUTHEAST BASINS

GEOHYDROLOGIC CONDITIONS REPRESENTED 

AND MODEL CHARACTERISTICS

The basin-group model of the southeast basins consists 

of two layers. The upper layer represents an unconfined 

aquifer 5 mi wide and extending the length of the basin. 

The lower layer represents the underlying basin-fill aquifer and extends the entire length and width of the 15- by 

50-mi basin (fig. 20). Water in this layer is simulated as 

occurring under confined conditions in the part underlying the upper layer and under unconfined conditions in 

the rest of the aquifer. The model was developed by using 

a uniform finite-difference grid that has a cell size of 1 mi2, 

includes 17 rows and 52 columns, and represents the 

average physical size of the basins in the southeast group. 

The upper layer was. represented by 250 active cells, and 

the lower layer by 690 active cells. Transmissivity values 

for the lower layer ranged from 300 ft2/d at the outer edge 

of the simulated area to 7,000 ft2/d near the center. Transmissivity of the upper layer was initially 4,000 ft2/d. Flow 

components that are typical of the southeast basins were 

included; initial model data and water-level configuration 

are shown in figure 21.
D50 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

ESTIMATED AND SIMULATED WATER-BUDGET COMPONENTS 

[Values are in acre-feet per year. Modified from Tucci, 1982, table 1]

1940-41

Components

Inflow:

Outflow:

Outflow to Palo Verde

Vallow

Estimated

12,000 

213,000 

12,000

237,000

235,000 

2,000

237,000

Simulated

12,000 

253,000 

12,000

277,000

271 ,000 

1,000

272,000

Mid-1960'S

Estimated

225,000 

92,000 

12,000

329,000

163,000

3,000 

172,000

338,000

Simulated

223,000 

97,000 

12,000

332,000

161,000

3,000 

168,000

332,000

1980

Estimated

236,000 

50,000 

12,000

298,000

71,000

3,000 

224,000

298,000

Simulated

235,000 

82,000 

12,000

329,000

120,000

2,000 

206,000

328,000

HYPRAULIC CHARACTERISTICS USED

CHARACTERISTIC

Hydraulic conductivity of upper layer , in 

feet per day 

Single uniform value

Transmissivity of lower layer, in feet 

squared per day 

Single uniform value

Vertical leakance between layers, in feet 

per second per foot 

Range

IN MODEL

VALUE

311

8,980

5x1 0'6 - 7x1 0'6

109° 108°

50 MILES -H

50 KILOMETERS

1940-

INDEX MAP SHOWING LOCATION OF 

PARKER VALLEY (SHADED)

FIGURE 18. Location of basin, estimated and simulated water-budget components, model grid and boundary, and observed and simulated water 

levels, Parker Valley.

Mid-1960's
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D51

COLUMNS

10 12 14 16 18 20 22 24 26 28 30

5 MILES

5 KILOMETERS

0 5 MILES

I I'M'I I' I

D 5 KILOMETERS

Base from U.S. Geological Survey 

State base maps, 1:500,000, 

Arizona 1974 and California 1968

Alluvial deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating 

ground-water flow

H

EXPLANATION

Cell representing

Constant head inflow

Constant head outflow

Drain (mid-1960s and!980 simulation)

River-aquifer connection

Geologic contact 

Boundary of flood plain 

Boundary of ground-water area

-295-

295-

Boundary of modeled area 

Upper and lower layers

Water-level contour, 1940-41, 

mid-1960's, and 1980 (Tucci, 

1982) Shows altitude of water 

level. Contour interval 5 and 

10 feet. Datum is sea level

Observed 

Simulated

FIGURE 18. Continued.
D52 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

ESTIMATED AND SIMULATED WATER-BUDGET COMPONENTS

[Values are in acre-feet per year. Estimated values modified

from Metzger and Loeltz (1973)].

Components

1960's

Estimated Simulated

Inflow: 

River losses. .......

Tributary inflow and 

unmeasured runoff.

Total ............

Outflow:

Evapotranspiration 

Pumpage. ........

Total..........

150,000

5,000

155,000

147,000

11,000

158,000

155,000

8,000

163,000

150,000

11,000

161,000

115' 109° 108°

0 50 MILES

0 50 KILOMETERS

INDEX MAP SHOWING LOCATION OF 

MOHAVE BASIN (SHADED)

HYDRAULIC CHARACTERISTICS USED IN MODEL

CHARACTERISTIC VALUE

Hydraulic conductivity of upper layer, in

feet per day

Range 6 - 311 

Average 120

Transmissivity of lower layer, in feet 

squared per day

Single uniform value 2,590

Vertical leakance between layers, in feet 

per second per foot

Single uniform value 6.8x10"6

FIGURE 19. Location of basin, estimated and simulated water-budget components, model grid and boundary, and observed>and simulated water 

levels, Mohave basin.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D53

COLUMNS

2 4 6 8 10 12 14 16 18 20

5 MILES

5 KILOMETERS

0 i i i i 5 MILES 

0 5 KILOMETERS

Base from U.S. Geological Survey 

State base maps, 1:1500,000, 

Arizona 1974, California 1968, 

Nevada 1964

EXPLANATION

Alluvial deposits

Bedrock of the 

mountains

Grid showing model 

cells used in 

simulating 

ground-water flow

Cell representing

Constant head Ground-water 

inflow, upper and lower 

layers

Constant head Open-water 

surface, upper layer only

River-aquifer connection 

Geologic contact

-485--

-500 

Boundary of ground-water area

Boundary of modeled area 

Upper and lower layers

Water-level contour, 1960's 

Shows altitude of water level. Contour 

interval 5 and 25 feet. Datum 

is sea level

Observed (Metzger and Loeltz, 1973) 

Simulated

FIGURE 19. Continued.
D54 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

TABLE 1. Range in values of geohydrologic properties in southeast basins and central and west basins and values used in basin-group models

I indicates no data]

Southeast basins Central and west basins

Basin characteristics

Maximum and

minimum values for Basin-group 

hydrologic properly or model 

basin characteristics

Maximum and

minimum values for Basin-group 

hydrologic property or model 

basin characteristics

Basin-fill width, in miles............................................... 7-21 15 7-27 12

Basin-fill length, in miles.............................................. 30-80 50 17-68 30

Length-to-width ratio.................................................... 1.5-8.0 3.33 1.0-8.5 2.5

Axial water-level gradient, in feet per mile................ 3^5 15 1-18 3

Mountain-front recharge, in acre-feet per year.......... 5,000-50,000 22,600 50-5,000

Saturated thickness, in feet.......................................... 100-1,000 200-1,000 100-1,000

Hydraulic conductivity (K), in feet per day............... 1-50 20 3^0

Transmissivity (T), in feet squared per day................ 100-50,000 300-7,000 2,000-15,000

Specific yield(s)............................................................. .03-.25 .20 .05-.18

Vertical-leakance rate (TK), in feet per second 10~4-iO~n 10~5-iO~n 10~4-KTn

per foot.......................................................................

Riverbed-leakance rate (RC), in feet per second i

per foot....................................................................... TOO X ^ - 2 

Maximum evapotranspiration rate (QET), in feet

per year at 100-percent density................................ 3-10 6 

Effective depth at which evapotranspiration

ceases, in feet............................................................. 0-50 10 

Vegetation density, in percent...................................... 20-80 50 

Cell coverage by phreatophytes, in percent............... 10-90 50 

Average annual pumpage, 1960-80, in acre-feet

per year....................................................................... 4,000-250,000 50,000-260,000 4,000-138,000 63,700

1,000-2,000 

200-700

9.5

2,000-6,600 

.12

10-8-io-n

RESULTS OF SIMULATIONS

The basin-group model of the southeast basins was 

used to investigate the relation of water-level contour 

shape to the controlling properties. Water-level contours 

in the southeast basins typically are U-shaped, are concave opposite the direction of ground-water flow, and 

qualitatively indicate the quantity of recharge entering a 

basin at the mountain fronts. The concave shape ranges 

from a broad U-shape in Douglas and Willcox basins to a 

V-shape in the lower San Pedro basin. Initial property values used in the basin-group model of the southeast basins 

were modified systematically to determine which properties exerted the greatest influence on water-level configuration. Comparison of water-level profiles at a particular 

section across the model (fig. 22) illustrates the effect of 

changing basin shape, transmissivity distribution, vertical 

leakance between layers, or increased recharge.

Changing the transmissivity distribution had the greatest effect on the shape of the simulated water-level 

profile, and altering the vertical-leakance value had the 

greatest effect on the computed head value. The model 

results indicate that areas within a basin where waterlevel gradients are steep generally near mountain

fronts have either hydraulic conductivities or saturated 

thicknesses that are on the order of 10 times smaller than 

areas where water-level gradients are gentle. Shallow 

depth to a pedimentlike surface is the most likely cause 

for small saturated thickness; hydraulic conductivity may 

be small because of the presence of consolidated rocks 

indigenous to the nearby mountains or indurated 

pre-Basin and Range rocks.

CENTRAL AND WEST BASINS

GEOHYDROLOGIC CONDITIONS REPRESENTED 

AND MODEL CHARACTERISTICS

The central and west basins, which are hydrologically 

similar, were analyzed by use of a single basin-group 

model. The principal differences in basins of the two 

groups are size, depth, and vertical extent of the finegrained facies. Also, because of larger and more significant inflow and outflow components in the central basins 

than in the west basins, the central-basin system response 

to development is delayed slightly as a result of capture of 

natural discharge. Subsequent to total capture of this 

flow, the response of the systems generally is the same.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D55

Change in pattern implies gradation in A' 

unconsolidated material from coarse 

grained to fine grained (mudstone and 

£ evaporate deposits) toward the center,.

PRE-BASIN AND 

RANGE SEDIMENTS

UPPER BASIN FILL , , , ; 

Confining unit , _

EXPLANATION

Stream alluvium

Alluvial deposits

Bedrock of the 

mounatins

Direction of groundwater flow

Evapotranspiration

FIGURE 20. Typical geographic and geohydrologic features in the southeast basins.
D56 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

COLUMNS 

12345678 910111213141516171819202122232425262728293031323334353637383940414243444546474849505152

1

2

3

4

5

6

7 

W 8

> Q O 9

tr 10 

11

12

13

14

15

16

17

o 

o

0

o

o 

o

*~

o 

o 

o

CM

.

o 

o

0

m

o 

o

0

^

0

o

0

m

o 

o 

o

(O

p ; 300

800

3.000 '

5.000

7.000

ol I 

o 1 [ 5.000 |

 * ! , I J . I

7.000

5.000

3.000

800

300

0

o 

o

i^

o 

o

(O

0 

0

m

o 

o

^

0

o

CO

1

o

0

CM

o

o

0

r

O

o

-H

I

...

o

o

'

TRANSMISSIVITY DISTRIBUTION OF THE LOWER LAYER 

^~~^~ Boundary of active model area 

5,000 Transmissivity Value in feet squared per day

MODEL BOUNDARY CONDITIONS

H River-aquifer connection cells, upper layer 

Ki Evapotranspiration cells, upper layer

H Constant head cells, upper and lower layer

[ 1 Constant flow cells to represent recharge, 

lower layer

WATER-LEVEL CONFIGURATION AT END OF STEADY-STATE SIMULATION

 1,400- WATER-LEVEL CONTOUR Shows altitude of simulated water level. Heavier line represents 

upper layer, finer line represents lower layer. Contour interval 100 feet. Datum is sea level

MODEL-INPUT CHARACTERISTICS

Maximum evapotranspiration rate (QET)

Depth at which evapotranspiration effectively ceases (ETDIST)

Interlayer leakance (TK)

Riverbed leakance (PC)

Hydraulic conductivity of the upper layer (K)

Saturated thickness of the upper layer (STRT-BOTTOM)

Grid cell dimensions (AX and AY)

1.5x1 CT8 feet per second

10 feet

1CT7 feet per second per foot

10~8 feet per second per foot

20 feet per day

200 feet

5,280 feet x 5,280 feet

FIGURE 21. Model grid, boundary conditions, and data input for steady-state basin-group model of the southeast basins.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS P57

The basin-group model that represents the central and 

west basins includes one layer for those analyses that 

involve basin geometry and two layers for analyses that 

involve pumping stress-related reactions. The model, 

which was developed to represent a 12- by 30-mi area, had 

a uniform finite-difference grid having 360 active cells and 

included 14 rows and 32 columns (fig. 23). The singlelayer simulation was developed to analyze the character 

of saturated lower-basin fill under unconfined conditions. 

A two-layer simulation represented saturated upper and 

lower basin fill separated by fine-grained sediments near 

the center of the basin. Evapotranspiration and riveraquifer exchange were not simulated because these components are negligible in the central basins and are not 

typical of the west basins. Flow components for the central and west basins can include mountain-front recharge, 

underflow into and out of the basin, and streamflow infiltration along the major basin drainage. The basin-group 

model representing these two groups was developed by 

using only the flow components of mountain-front 

recharge and underflow out of the basin because that configuration is most typical (fig. 24).

RESULTS OF SIMULATIONS

The basin-group model that represents the central and 

west basins was used to analyze the significance of basin 

shape and recharge location on water-level gradients. 

Ranges and location of recharge are the principal hydrologic differences between the central and west basins. In 

the west basins, recharge is small on the order of a few 

hundred to about 5,000 acre-ft/yr. In the central basins, 

recharge is large generally ranging from about 5,000 to 

about 50,000 acre-ft/yr and occurs as a combination of 

underflow and mountain-front recharge.

Four basin shapes were considered (1) a straightsided channel with open ends, (2) a slope-sided trough 

with open ends, (3) a scoop with one end constricted and 

one end open, and (4) a long bowl with both ends 

constricted (fig. 25). Initial flow through each basin shape 

was established by constant-head cells at each end of the 

model grid and an initial water-level gradient of 3 ft/mi. 

Model-calculated flow through each basin was then 

applied as a uniformly distributed mountain-front 

recharge and as an intermittently spaced recharge around 

the basin perimeter. Results indicate that the water levels 

are practically insensitive to basin shape. Recharge location causes a significant change in water-level configuration. Where recharge is distributed around the margin of 

the basin rather than as underflow at the upgradient end, 

water-level gradients change from uniform to variable 

throughout the basin. Differences in water levels resulting 

either from uniformly distributed recharge or intermittently spaced recharge were negligible. Intermittently 

spaced recharge depicted in the model probably is more 

realistic and simulates recharge at intermittently spaced 

points where streams emerge from the mountains.

RESULTS FROM SPECIFIC-BASIN AND 

BASIN-GROUP SENSITIVITY ANALYSIS

Numerical models were developed for 12 of the 72 alluvial basins in south-central Arizona and for 2 idealized 

alluvial basins representing a typical southeast basin and 

a typical central-west basin. The models were analyzed to 

determine (1) the relation between the basin category and 

the importance of various model-input parameters and (2) 

the adequacy and transferability of information used to 

develop basin models.

RELATION OF BASIN CATEGORY AND MODEL INPUT

A controlled examination of model response to 

changes in property values can be used to show the significance of each property and the data that define that 

property. The significance of a particular property 

depends on the hydrologic and geologic characteristics of 

a basin. The basins were grouped on the basis of these 

characteristics.

If large changes in simulation results take place with 

only minor changes in a tested property, the model reliability is highly dependent on the correct value for that 

property. In contrast, if large changes in property values 

result in insignificant changes in simulation results, the 

model is insensitive to that property, and an estimate of 

the actual value will suffice for use in the model. Knowledge of the model sensitivity to changes in each property 

may allow streamlining and redirection of data-collection 

activities. Also, the validity of transferring hydrologic 

boundaries, flow components, and hydrologic factor values from one basin to another can be qualitatively judged 

where alternative model results can be compared.

The properties that are significant in development of 

models of the basins are divided into three groups 

hydrologic boundaries, flow components, and hydrologic 

factors. Components of each property are listed in table 2 

with a subjective evaluation of their importance in steadystate and transient model development. The hydrologicboundaries property represents the physical constraints 

on aquifer extent. The flow-components property 

includes recharge to and discharge from the aquifer. The 

hydrologic-factors property controls the movement and 

storage of water within the aquifer.

Selected properties are more important in steady-state 

simulations than in transient simulations. For steady-state 

conditions, water-level configuration and hydraulic 

conductivity are essential for model development. Flow
D58 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

1,900

1,800

1,700

LU 

CO 

LU

§ 1 -600

CO

1,500

1,400

1,300

1,200

87654321012345678

DISTANCE FROM BASIN AXIS, IN MILES

EXPLANATION

CASE 1

Input characteristics for steady-state 

model (used as standard for comparison)

Hydraulic conductivity (K=20) of upper layer, in 

feet per day

Transmissivity (5,000) of lower layer, in feet 

squared per day

Interlayer leakance (TK) =10~7 feet per second 

per foot

Mountain-front recharge (Qm)=26 cubic feet 

per second

 Simulated water-level profile

CASE 2

Increased recharge

Mountain-front recharge = 52 cubic feet per 

second (All other parameters same as in 

CASE 1)

Simulated water-level profile

FIGURE 22. Effects of modifications of mountain-front recharge, interlayer-leakance values, lower-layer transmissivity, and 

basin width on water-level profiles of the steady-state basin-group model of the southeast basins.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D59

CASES Modified interlayer leakance

A. Interlayer leakance=10~H feet per second per 

foot (All other parameters same as in CASE 1)

Simulated water-level profile

B. Interlayer leakance =10 feet per second per 

foot (All other parameters same as in CASE 1)

Simulated water-level profile

CASE 4

\

K=20 Modified lower-layer transmissivity

A. Transmissivity as shown at left. (All other 

parameters same as in CASE 1)

Simulated water-level profile

\

K=20

TK

o

o

in

B. Transmissivity as shown at left (All other 

parameters same as in CASE 1)

mm mm mm mm Simulated water-level profile

CASES

9 miles 

3 miles

K=2Q

TK

Decreased basin width

Basin width as shown at left (All other 

parameters same as in CASE 1)

Simulated water-level profile

FIGURE 22. Continued.
D60 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

Change in pattern implies gradation in 

unconsolidated material from coarse 

grained to fine grained (mudstone and 

evaporate deposits) toward the center

UPPER BASIN FILL 

LOWER BASIN FILL

PRE-BASIN AND 

RANGE SEDIMENTS

EXPLANATION 

Alluvial deposits

Bedrock of the 

mountains

Direction of groundwater flow

FIGURE 23. Typical geographic and geohydrologic features in the central and west basins.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D61

COLUMNS 

12345678 91011121314151617181920212223242526272829303132

500

1.000

1.500

2.000

| 2,500 |

2.000

1.500

1.000

500

1,000

SATURATED THICKNESS 

Boundary of active model area 

Saturated thickness of model area Value in feet

u

n

LJ LJ LJ u u

E

E

n n n n n

MODEL BOUNDARY CONDITIONS 

| Constant head cells 

HI Constant flow cells to represent recharge

o

CO

o

1,030

WATER-LEVEL CONFIGURATION AT 

END OF STEADY-STATE SIMULATION

Water-level contour Shows altitude of simulated 

water level. Contour interval 100 feet. Datum is sea level

MODEL-INPUT CHARACTERISTICS

Hydraulic conductivity of the upper layer (K) 

Grid cell dimensions (AX and&Y)

20 feet per day 

5,280 feet x 5,280 feet

FIGURE 24. Model grid, boundary conditions, and data input for steady-state basingroup model of the central and west basins.
D62 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

STRAIGHT-SIDED CHANNEL

Total flow=2,078 acre feet per year

TROUGH

Total flony=f,332 acre feet per year

SCOOP

Total flow=1,250 

acre feet per year

BOWL

Total How=1,190 

^cre feet per year

UNDERFLOW ONLY

COLUMNS 

2345678 91011121314151617181920212223242526272829303132

1

2

3

4

5

CO 6

5 7

0 8 01 q

y 

1011

12

13

14

c

C

c

3 C

3 T

3 C

3 C

e

3 C

3 C M C

3 C

3 C

0 *

3 «

3 C 

t I

3 C

3 C

n c

3 C

3 C D r

3 C

3 C

- 0

3 C

3 C

0 C

3 C

3 

D

3

COLUMNS

CASE 1

12345678 91011121314151617181920212223242526272829303132

1

2

3

4

5

CO 6 

§ 7

0 ft

cc 9 

10

11

12

13

14

c 

c 

e

9 C

3 i 

3 C

3 C

e

3 C "*

3 C M C 

3 C "*

3 C 

0 * 

3 C "*

3 C

t H 

3 C "*

3 C

n c

3 C "*

3 C D r

3 C "*

3 C

« e

3 C **

3 < 

§ t

**

3

n

3

CASE 4 

COLUMNS

2345678 91011121314151617181920212223242526272829303132

1

2

3

4

5

CO 6 

§7 

0 8

10

11

12

13

14

c 

c 

c

3 C

3 T

3 C

3 C

e

3 C

3 C M C 

3 C

3 «

1 «'3 (

t I

3 <

3 C

n t

3 <

3 <

o r

3 <

I

ll 3 O ic 

-- 00 

3 0 

I

/

I

COLUMNS

CASE 7

12345678 91011121314151617181920212223242526272829303132

1

2

3

4

5

CO 6

JM^ 9 

10

11

12

13

14

c 

c 

c

I I ° ? 30 C

I

3 C M C 

3 C

3 C 

O * 

3 C

3 C

t u

3 C

3 C

>

3 <

1 <r

\ 3 0 C

«. 00 C 

3 0 C

J

3

n

3

CASE 10

WATER-LEVEL CHANGE STATISTICS, IN

NOTE: Case 1

Case: 1 2

Maximum 

difference

Minimum 

difference

Average 

difference

Average absolute 

difference

- 27.9

- -13.2

3.3

9.2

3

28.6

-13.5

4.1

9.4

FEET

is the standard for comparison

4

0.9

-2.0

-0.6

0.9

5

31.3

-9.5

6.3

9.2

6

31.8

-9.9

7.1

9.7

7

5.3

0.07

1.7

1.7

8

33.2

-9.0

8.3

10.4

9

34.0

-9.3

9.3

11.0

10

4.6

-4.3

-0.4

2.4

11

26.3

-18.8

0.1

11.0

12

27.2

-19.1

1.1

11.0

FIGURE 25. Variations in data used and results of analysis of basin shape and recharge location and distribution in basin-group model of the central 

and west basins.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D63

MOUNTAIN-FRONT RECHARGE 

UNIFORMLY DISTRIBUTED

COLUMNS 

12345678 91011121314151617181920212223242526272829303132

1

2

3

4

5

co 6 

§ 7

1

2

3

4

5

co 6 

§ 7 

0 8 

^9

10

11

12

13

14

1

2

3

4

5

co 6 

§7

1

2

3

4

5

co 6 

§ 7

CASE 2 

COLUMNS

12345678 91011121314151617181920212223242526272829303132

_L_Z

CASE 5 

COLUMNS

12345678 91011121314151617181920212223242526272829303132

-r 

CASE 8 

COLUMNS

2345678 91011121314151617181920212223242526272829303132

MOUNTAIN-FRONT RECHARGE ONLY WHERE 

STREAM CHANNELS EXIT MOUNTAINS

COLUMNS 

2345678 91011121314151617181920212223242526272829303132

1

2

3

4

5

co 6 

§7

m

o o

if) U>

o o

_HL

CASES 

COLUMNS

2345678 91011121314151617181920212223242526272829303132

m/ m

CASE 6 

COLUMNS 

2345678 91011121314151617181920212223242526272829303132

1

2

3

4

5

co 6 

§7°l

10

11

12

13

14

1

2

3

4

5

co 6 

§ 7

8|

10

11

12

13

14

id \ E3

I

CASE 9 

COLUMNS

2345678 91011121314151617181920212223242526272829303132

E3

i- CM CO

|M| / fil

CASE 11 CASE 12

BASIN CHARACTERISTICS

Basin length 

Basin width

30 miles 

12 miles

Maximum saturated thickness 700 feet

Hydraulic conductivity 

Initial water-level gradient 

Minimum saturated thickness

9.5 feet per day 

3 feet per mile 

200 feet

EXPLANATION

1,030 Water-level contour Shows altitude of water level

calculated by the model. Contour interval 10 feet. Datum 

is sea level

I I Constant-head cell 

H Constant-flow cell

FIGURE 25. Continued.
D64 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

TABLE 2. Relative importance ofhydrologic boundaries, flow components, and geohydrologic factors by basin group for

modeling steady-state and transient conditions

[1, important; 2, moderately important; 3, not important; indicates property does not exist in this category]

Basin group

Steady-state conditions

Hydrologic 

boundaries

c

configuratio

2

1

bottom >r 

S

q*

)oundaries

3

1

5

Flow components

§0

13 6

B 3

H

CD

recharg ront 

Js

o

inrflow 

p

rflow out

 a

^

infiltration

I

8

£

discharge iv 

0

<2

i>

oa

anspiration

S

§*

H

Hydrologic factors

conductivity

u

si

X

anisotropy

H

fi

£

leakance ti 

|(-i

g 

«

Southeast

Central

West

Colorado River

1

1

1

2

1

1

2

2

3

3

3

3

2

2

3

1

2

3

3

2

3

3

3

2

3

3

3

2

2

2

1

2

2

2

1

1

1

2

1

1

1

2

2

1

1

2

3

2

2

2

1

Transient conditions

Basin group

Hydrologic 

boundaries

§

Water-level configurat bottom Aquifer boundaries Lateral 

Flow components

1recharge Areal 

00

Mountain-front rechar Underflow in Underflow out

2infiltration Stream 

00 

CO

0

CQ

Drain discharge Evapotranspiration

Pumping

Hydrologic factors

ft

Hydraulic conductivi Vertical anisotropy River/drain leakancc Storage coefficient

Southeast 

Central 

West 

Colorado River

133 2 233222212121

122 2 333233312231

122 3 33333 313231

233 2 333111123313

'Includes recharge from excess applied irrigation water in some basins during recent periods. 

2Includes recharge from canals and laterals in a few basins.

component values require varying degrees of reliability 

owing to the different significance of the components as 

related to the hydrologic categories. For example, 

although mountain-front recharge is important in the 

southeast basins because it represents the largest single 

predevelopment recharge source, in the west basins 

mountain-front recharge is practically nonexistent and 

therefore is not important.

For transient simulations, reliable water-level-change 

information, pumpage, and specific yield or storage coefficient are the most important data in the history-matching process, especially for those basins in which 

significant storage depletion has occurred. The relative 

significance of other boundary conditions and hydrologic 

factors is related to the individual setting of the geohydrologic system. Flow components also have differing importance depending on the individual basin setting; however, 

the magnitude and distribution of pumping commonly are 

critical items. This is particularly true in the central and 

west basins where large-scale pumping has resulted in a 

large-scale overdraft of ground-water resources.

INFORMATION ADEQUACY AND TRANSFERABJLTTY

Model-sensitivity analyses provided a qualitative 

assessment of the adequacy of values and areal distributions of properties used in the models. In general, the 

more sensitive a model is to a particular property, the 

greater the data requirements will be. Similarly, transferability of property values and distributions is related to
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D65

TABLE 3. Geohydrologic factors evaluated as part of the modeling process

[ indicates this factor was evaluated in the model]

Factors

Southeast

^c

c§ XI

2

1

CO

»^<

 -1

.5

X! 

1=O 

en

I

.5

c§

XI

2 o

Sanwer 

>_4 9

c

aX!

O

i

c

1XI

!

Central

.5

o

en

o

Sf 1

 i

West

.a

XI

1

aterman

jE

ri

B

1

SMcMull

Colorado River

c 

1

XI

I1

s

c

1XI

fc-.

Basin group

a .5

Southea

'i

X! 

 M

cn

1

ntral anc

cu

O

Hydrologic boundaries

Water-level configuration

Aquifer bottom

Lateral boundaries

Flow components

Area! recharge

Mountain-front recharge

Underflow in

Underflow out

Stream infiltration

Base-flow discharge

Drain discharge

Evapotranspiration

Pumping

Hydrologic factors

Hydraulic conductivity

Vertical anisotropy

River/drain leakance

Storage coefficient

model sensitivity. In the situation where a model is insensitive to a property, it is directly transferable to all the 

basins of a group. The models developed as part of this 

study were used to explore the effects of various geohydrologic properties on model results. The properties that 

were specifically addressed in each model are indicated 

in table 3.

HYDROLOGIC BOUNDARIES

The most significant hydrologic boundary was the 

water-level configuration, which represents the upper 

boundary of the ground-water flow system. Accurate 

area! and vertical head-distribution data are vital in understanding the flow system. In general, because of similarities in hydrologic properties of sediments that make up

the aquifers and the existence of regional similarities or 

trends in magnitude of some components of inflow to the 

basins, hydraulic gradients are similar in basins of a particular group. Moreover, the general shape of the 

water-level contours can be reliably predicted on the 

basis of basin categorization (see section "Results of 

Simulation" for the basin-group models).

Explicit definition of the physical limits of the aquifer 

in the typical geohydrologic setting of the study area, 

including the bottom and sides of the aquifer system, is 

least important for steady-state modeling. For transient 

modeling, lateral boundaries of the aquifer in close 

proximity within several miles of pumping centers 

need to be well defined because of boundary effects on 

the size and shape of a cone of depression. The model for 

McMullen Valley was used to evaluate the sensitivity of
D66 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

40

HI

o

CD

ZJ 

0_

LU 

CO 

< 

LU 

DC 

O

01 

O

DC

30

20

10

5,000 10,000 15,000

DISTANCE OF IMPERMEABLE BOUNDARY 

FROM PUMPING NODE, IN FEET

FIGURE 26. Effect of lower-layer boundary location on simulated 

drawdown in McMullen Valley (from Pool, 1987).

simulated drawdown near an impermeable bedrock 

boundary (Pool, 1987). The location of the bedrock 

boundary was varied from one-quarter mile to several 

miles from the pumping cell. In this storage-depletion system, part of the well field is within a few miles of the lateral boundary of the basin. The size and shape of the cone 

of depression are influenced greatly by boundary effects. 

Drawdown increases significantly if the pumping center is 

within 7,000 ft of the aquifer boundary (fig. 26). The distance is unique to the specific situation and will change if 

aquifer character, geometry, or magnitude and duration of 

pumping are changed.

FLOW COMPONENTS

Pumpage estimates are available in most areas or can 

be made on the basis of irrigated area. Although pumpage 

is the most well-documented flow component in most 

basins, it still is poorly known. Data documenting other 

flow quantities are, for practical purposes, nonexistent, 

and values must be estimated by indirect means. The 

other flow components, however, generally are much 

smaller than pumpage, so that errors in these other values

have less effect on results of simulations than do errors in 

pumpage estimates. Modeling is useful in determining the 

reasonableness of flow-quantity estimates. Large errors in 

estimated quantities are possible because the typical 

water-budget method of estimation allows errors to accumulate in the unknown quantities.

In the southeast basins, mountain-front recharge represents the principal inflow; however, interbasin underflow 

and stream-aquifer exchange can be of nearly equal 

importance. Evapotranspiration is the main discharge 

component, and model results are sensitive to the evapotranspiration rates; however, underflow and stream base 

flow are important in several basins. Varying the evapotranspiration rate within a reasonable range of values in 

the model of the upper San Pedro basin (Freethey, 1982) 

changed the net discharge by about 40 percent and 

resulted in an imbalance of the other flow components.

Results from the upper San Pedro basin model indicated that generalizing mountain-front recharge was feasible (Freethey, 1982) but not recommended if sitespecific water-level matches are a goal. In models for five 

other southeast and central basins, the use of generalized 

mountain-front-recharge values presented no undue difficulty in calibration.

Other significant flow components were investigated 

during model development for specific basins. Results of 

simulations of Avra Valley indicate that return of excess 

applied irrigation water may be about 30,000 acre-ft/yr 

since the early to mid-1960's (Anderson, 1983). Results 

from nine other models also indicate that return of excess 

irrigation water to the principal basin aquifer probably is 

a significant flow component. Few field data exist with 

which to verify irrigation return flow, however, and only 

indirect evidence is available with which to quantify the 

flow component.

The models of the Colorado River basins were sensitive to the rate of river leakage to the principal aquifer and 

the rate of consumptive use by riparian vegetation. The 

magnitude of these two components overwhelms other 

much smaller components. The value of riverbed leakance controls the quantity of infiltration from surface 

flow. Drains have been installed in a few basins, and the 

drain leakance is not as great a control on the water budget as is the riverbed leakance (fig. 27), probably because 

of the much smaller areal extent of drains. The groundwater systems of the Colorado River basins are limited in 

areal extent, and modeling on a regional scale as in this 

study is useful mainly in examining the magnitude of the 

water-budget components. Balancing the water budget of 

a basin so that all flow components are within a reasonable range of values was the best method of calibration 

rather than relying on head values and configurations.

Models of the Colorado River basins were most sensitive to changes in the evapotranspiration rate (Tucci,
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D67

tr 200

LLJ

tr

LLJ 

0.

150

LU 

LU 

U_

LLJ 

CC 

O 

<

LL 

O

CO 

Q

^ 100

O

X

LU

tr

Q

LU

I

1

CO

50

 Drain value, adjusted

 River value, adjusted

 River flow

A Evapotranspiration

 Drain flow

i i i i i ml i i i i i ml

0.01 0.1 1.0 10 100

LEAKANCE, IN FEET PER DAY PER FOOT

FIGURE 27. Variation in components of loss from river flow, 

evapotranspiration, and drain flow resulting from changes in the 

leakance factor for drains and rivers in the Parker Valley model (from 

Tucci, 1982).

1982). Increasing the average evapotranspiration rate by 

1 (acre-ft/acre)/yr caused about the same increase in river 

infiltration in Mohave basin as in Parker Valley about 

60,000 acre-fVyr (fig. 28). Because of model sensitivity to 

changes in evapotranspiration, changes in vegetation type 

and density are important considerations for similar 

basins in the Colorado River area Sensitivity to evapotranspiration also is expected for models of other alluvial 

aquifers traversed by a perennial stream, such as some 

basins of the highland and southeast groups.

HYDROLOGIC FACTORS

For a steady-state model, the most significant hydrologic factor of a basin is the hydraulic conductivity; for a 

transient model, specific yield or storage coefficient is the 

most significant property. Vertical anisotropy is important 

in basins containing an extensive fine-grained unit that 

functions as a leaky confining unit between two layers of 

the aquifer system.

Hydraulic conductivity and transmissivity arrays were 

developed for most basin models by using specificcapacity data and flow-net analyses. Hydrologic factors

tr

CC >

!£tr

^ LLJ tr Qik

200

100

31LLJ ^ 

QZ 0~

AVERAGE EVAPOTRANSPIRATION RATE, 

IN ACRE-FEET PER ACRE PER YEAR

FIGURE 28. Variation in net loss from river flow in Mohave basin and 

Parker Valley related to change in average evapotranspiration rate.

were generalized to analyze the adequacy of using a constant rather than a variable distribution of values for areas 

where no field data existed. Results from this analysis in 

the upper San Pedro basin indicated that the model was 

sensitive to generalized aquifer hydraulic conductivity 

(Freethey, 1982). Use of a single average value for hydraulic conductivity in the model of the Benson basin proved 

unsatisfactory because of the large variability in the 

potentiometric-surface gradient created by variations in 

the hydraulic conductivity. The greatest disparity 

occurred at the basin margins where water levels were 

more than 150 ft lower than those in the nongeneralized 

simulation; water levels in the flood-plain area were about 

30 ft higher than those in the nongeneralized simulation. 

In the same test, simulated water-budget values also fell 

outside the acceptable range established from field data 

and previous studies.

The basin-group model of the southeast basins was 

used in a sensitivity analysis of system response to variations in values of hydraulic conductivity of the upper 

layer, transmissivity of the lower layer, interlayer leakance, and riverbed leakance. The results of varying these 

properties in terms of average head change and percentage change in flow quantities allow a qualitative comparison of the general significance of each property (table 4).
TABLE 4. Sensitivity of model-calculated heads and flow quantities to changes in values of flow components and hydrologic factors in southe

[* indicates that reasonable range is exceeded; indicates that multiplication factor was not tested]

Multiplication 

factor

Flow components

Maximum 

evapotranspiration rate

Average absolute 

head change, 

in feet

Lower 

layer

Upper 

layer

Change 

in 

quantity 

of flow, 

in percent

Riverbed leakance

Average absolute 

head change, 

in feet

Lower 

layer

Upper 

layer

Change 

in 

quantity 

of flow, 

in percent

Mountain-front recharge

Average absolute 

head change, 

in feet

Lower 

layer

Upper 

layer

Change 

in 

quantity 

of flow, 

in percent

Hydrologic factor

Hydraulic conductivity, 

upper layer

Average absolute 

head change, 

in feet

Lower 

layer

Upper 

layer

Change 

in 

quantity 

of flow, 

in percent

Transmissivity, 

lower layer

Average absolute 

head change, 

in feet

Lower 

layer

Upper 

layer

Change 

in 

quantity 

of flow

in percen

;rease

w

o> 

Q

, 0.0001 

^ .001

.01

.1

.2

.25

.3

.5

1

2

3

4

4.5

1 5

1 1 10

£ T 100

* 

*

*

*

*

*

4.98

0

4.07

4.85

*

*

*

* 

*

*

*

*

*

4.72

0

4.02

4.81

*

*

*

* * 

* *

* *

* 2.61

* _

* _

-19.93 

 

0 0

 

5.76 

 

5.89 

* _

* 3.09

* *

* * * 

* * *

* * *

2.91 -22.25 *

_ _ *

_ _ *

_ _ *

 14.33

000

 26.60

_ _ *

_ _ *

_ _ *

_ _ *

3.40 41.88 *

* * *

* 

*

*

*

*

*

*

3.83

0

5.51

*

*

*

*

*

*

* 

*

*

*

*

*

*

-28.44

0

52.42

*

*

*

*

*

*

* 

*

*

*

3.13

1.39

0

1.35

2.57

*

*

* 

*

*

*

3.07

1.33

0

1.23

2.25

*

*

* 

*

*

*

-8.40

-5.11

0

9.69

37.25

*

*

* 

*

*

*

*

62.32

21.64

0

16.24

41.93

*

*

*

*

* 

*

*

*

*

3.07

2.17

0

7.80

30.20

*

*

*

*

* 

* 

* 

* 

* 

1.13 

 

-.35 

0

5.78 

 

19.05 

* 

* 

* 

* 
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D69

250

240

E co

LU <

at 

§1UJ CO 

I Z

_ LJJ ? co

UJ CC

oo

z u,

UJ z

cc ?

Ul =: u. cc

It CO

Sg eo

 " ~ ~ ~ Douglas basin model 

Benson basin model 

Layer 1 

Layer 2

DC 

LJJ

0.01 0.1 1.0 10 100 

RELATIVE VALUE OF LEAKANCE BETWEEN LAYERS

FIGURE 29. Sensitivity of model-calculated heads in the Benson and 

Douglas basin models to variations in interlayer-leakance value.

The value used for vertical interlayer leakance was found 

to have almost no effect on simulation results until the 

value was less than 10"8 (ft/s)/ft. Head differential 

between layers increases significantly when smaller 

values of vertical leakance are used. Interlayer-leakance 

values were varied through several orders of magnitude in 

the models for Benson and Douglas basins. On the basis 

of well-log data, the basin-center fine-grained facies is 

more extensive and thicker and serves more as a confining bed in Benson basin than in Douglas basin. Accordingly, model sensitivity to variations in interlayer leakance 

is greater in the model for Benson basin than in the model 

for Douglas basin (fig. 29). Both model responses indicate 

an increasing sensitivity as the value of interlayer leakance is decreased. The conclusion is that a model is 

most sensitive to interlayer leakance where the finegrained sediments form a massive confining bed and the 

vertical hydraulic conductivity is small.

In transient models of the upper San Pedro and Douglas basins, a uniform average specific yield was used in

20

40

UPPER SAN PEDRO 

BASIN MODEL

B. DOUGLAS

BASIN MODEL

200

10

YEARS OF PUMPING

FIGURE 30. Maximum simulated drawdown using a variable area! 

distribution of specific yield and a uniform specific yield 

representing the average value for the models of upper San Pedro 

basin and Douglas basin.

several sensitivity runs. Results indicate that using a uniform average value of specific yield does not change 

model results significantly compared with simulations 

using a moderately variable array of specific yield. As the 

degree of development increases, however, the difference 

resulting from the two arrays increases. If the deposits 

within the basin area affected by pumping have a large 

areal difference in grain size, an average specific-yield 

value can result in significant differences in simulated 

water levels.

Comparison of maximum drawdowns using average 

and variable specific yield in the upper San Pedro model 

(fig. 3(L4), in which all other values used in the model are 

fixed, indicates that using a uniform average value results
D70 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

in a slightly smaller maximum drawdown because the 

average value is slightly larger than the variable value in 

the cone of depression. In the model of the Douglas basin, 

a similar analysis using an average specific yield resulted 

in a larger maximum drawdown (fig. 3QB). In this area, 

the actual specific yield is smaller than the basinwide 

average. Use of a basinwide uniform value of specific 

yield provides satisfactory results for a general model calibration, but for site-specific analyses, especially where 

the aquifer material is heterogeneous, areally variable 

specific-yield values are needed.

The effect of varying hydrologic properties also was 

examined by imposing pumping stress on the basin-group 

model of the central and west basins in conjunction with 

changing specific yield and by altering the location of the 

area where vertical leakance is small because of the presence of a confining layer. First, the basin-group model 

was used to examine variations in specific yield and 

pumping locations and quantity. The second approach 

introduced a fine-grained layer in the basin-group model 

to examine the resulting effects on drawdown. Location 

and distribution of pumping were varied to approximate 

conditions that exist in three of the modeled basins of the 

west group McMullen Valley, Harquahala Plain, and 

Waterman Wash basin.

In the single-layer model, a simulated pumping rate of 

1 (ft3/s)/mi2 was used over an area of 88 mi2 in the central 

part of the modeled area for a 20-year duration. The 

resulting water-level declines using six different specificyield distributions indicate the model sensitivity to the 

absolute value and area! distribution of specific yield 

(fig. 31). Using a constant value of 0.12 for specific yield 

results in virtually the same water-level decline as using a 

random spatial distribution of specific-yield values of 

0.06, 0.12, and 0.18. Water-level decline is proportional to 

the average value of specific yield in these storagedepletion systems. Actual ground-water withdrawal in 

several west basins is concentrated near the outflow end 

of the basin or in more than one pumping center. The 

effect of concentrating the pumping is shown by waterlevel declines on an axial profile through the basin-group 

two-layer model (fig. 32). The increase in decline because 

of the presence of a fine-grained facies is also shown. The 

magnitude of the increase in water-level decline depends 

on the quantity of pumpage in areas where the finegrained unit is present. Redistribution of pumpage to 

within 2 mi of the basin boundary increased water-level 

declines less than 10 percent.

HYDROLOGIC INSIGHTS GAINED FROM 

MODELING SOUTHWEST ALLUVIAL BASINS

The degree of ground-water development represented 

in the 72 alluvial-basin aquifers in south-central Arizona

ranges from virtually pristine to being dewatered beyond 

the point that results in land subsidence. Development in 

many basins began prior to comprehensive data-collection efforts, thus precluding conceptualization of predevelopment hydrologic conditions on the basis of fact. 

Results of modeling 12 alluvial basins in various stages of 

development aided in more clearly understanding hydrologic processes occurring before pumping began. Modeling more than one basin in four of the groups has shown 

that basins within each of these four groups of basins 

have similar hydrologic systems and demonstrate a similar response to pumping stress. Although no basins in the 

highland group were modeled because of a lack of information, their geographic and hydrographic environments 

suggest that their hydrologic systems are similar to those 

of the Colorado River basins; thus, insights about Colorado River basins gained from modeling results are also 

considered to apply to the highland basins.

Results of simulations of hydrologic systems of 

the southeast group indicated that (1) mountain-front 

recharge probably occurs nonuniformly along the margins of the basins; (2) hydrologic factors such as transmissivity, hydraulic conductivity, specific yield, and confined 

storage coefficient are laterally variable; (3) the distribution of storage coefficient or specific yield is not as critical as the average value of these properties in determining 

the amount of water-level decline caused by pumping; 

(4) accurate areal distribution and values of transmissivity or hydraulic conductivity are critical in determining 

the configuration and gradient of the potentiometric surface; (5) variably confined aquifer conditions near basin 

margins probably are a result of interlayering of coarse 

sediments from the mountains and fine sediments along 

the axis of basins originating from lake deposits; and 

(6) pumping causes discharge by evapotranspiration and 

stream gains to decrease, even if pumping is from substantial depth in the basin fill.

The model of Willcox basin (closed drainage) indicated 

that extensive lacustrine deposits in the center of a basin 

can create a boundary condition that will accelerate 

water-level declines in nearby pumping wells. The model 

of the Tucson basin indicated that postdevelopment 

water-level declines can result in a decrease in transmissivity of the principal deep aquifer and that these 

decreases are not proportional to the decrease in saturated thickness because the upper basin fill is coarser 

grained and more permeable than the lower basin fill. As 

the upper basin fill is dewatered, that part of the aquifer 

having larger hydraulic conductivity is eliminated, leaving 

only the underlying finer-grained basin fill to transmit 

water to wells.

Models of the hydrologic systems of the central and 

southeast basin groups indicated that (1) excess applied 

irrigation water can eventually percolate down to
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D71

40

LJJ 

LJJ

uT 80

o

LJJ 

D

tr

LJJ

I

120

160

200

6 8 10 12 14 16 18 20 22 24 26 

DISTANCE FROM THE DOWNSTREAM END OF THE BASIN, IN MILES

28 30

Model characteristics

Transmissivity distribution as in figure 22, Case 1

General basin shape

EXPLANATION

CASE1

Model grid and boundary conditions

Specific yield variations

Uniform values 

throughout model area

Water-level decline,

Specific yield=0.06 

Specific yield=0.12 

Specific yield=0.18 

Variable values

I' 1*1 1*1 ! ] 1*1

HI HI HI HI

E

( 

Location 

of profile 

^ shown 

above

Specific yield randomly assigned 

to cells; values 0.06, 0.12, or 

0.18; average value 0.12

CASE 2

 Constant-head cell 

13 Constand-flow cell

D Pumping 1 cubic foot per

second per cell for 20 years

1 0.06

0.12

I

0.18

CASES

0.06

0 12

0.18

Variable specific yield

concentric areas with

low value at basin center

(as shown in sketch at left)

 Water-level 

decline

Variable specific yield

with low value at 

downstream end of basin 

(as shown in sketch at left)

mm mm. m^ Water-level 

decline

FIGURE 31. Response of the basin-group model of the central and west basins to variations in area! distributions and values of 

specific yield.
D72 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

UJ

50

100

UJ 150 Q

_i 

UJ

UJ_J

£ 200

I

250

300

25 50

PERCENT PUMPING WHERE 

FINE-GRAINED UNIT PRESENT

D 2 4 6 8 10 12 14 16 18 20 22 24 

DISTANCE FROM THE DOWNSTREAM END OF THE BASIN, IN MILES

EXPLANATION

Pumping distribution model scenarios 

Model characteristics

Transmissivity distribution as in figure 22, Case 1 

Specific yield of upper layer=0.12

Interlayer leakance (7"K)=10~' feet per second per foot

26 28

General 

basin 

shape

CASE 1 Pumping from lower layer uniformly distribution over 88 square miles 

at rate of 1 cubic foot per second per square mile for 20 years

(1A) 

Water-level decline along basin axis, no fine-grained

unit 

Water-level decline along basin axis with fine-grained

unit at basin center, interlayer leakance=10~11 feet

per second per foot

Pumping from lower layer uniformly distributed over 44 square miles at rate of 2 cubic 

feet per second per square mile near downstream end of basin for 20 years

(2) = Water-level decline along basin axis, no fine-grained

unit

Water-level decline along basin axis with fine-grained 

unit at basin center, interlayer leakance=10~H feet 

per second per foot

(2A) ==== =

Pumping from lower layer uniformly distributed over 44 square miles in two pumping 

centers at rate of 2 cubic feet per second per square mile for 20 years

Qj^ Water-level decline along basin axis, no fine-grained 

_ unit

Water-level decline along basin axis with fine-grained

unite at basin center, interlayer leakance=10-11 feet

per second per foot

(3A) -- 

30

FIGURE 32. Water-level declines resulting from variable areal distributions of pumping from the lower layer and presence or absence of 

a fine-grained fades in the basin-group model of the central and west basins.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D73

recharge the aquifer even in basins where the aquifer is 

several hundred feet below land surface and (2) the time 

required for this water to move through the unsaturated 

zone may be several decades.

Models of hydrologic systems in the central and west 

basin groups indicated that (1) the direction and rate of 

ground-water movement prior to development are 

defined primarily by the location of recharge sources, 

even though the magnitude of this recharge is small; 

(2) after development in the basin, the direction and rate 

of ground-water movement are governed primarily by the 

location of pumped wells and the amount of withdrawal 

occurring at each well or well field; and (3) because 

recharge is small and storage depletion is large, waterlevel declines near pumped wells can be accelerated if the 

cone of depression intersects buried bedrock along the 

margins of the basin.

Models of hydrologic systems in the Colorado River 

basin group indicated that near-surface hydrologic processes, such as (1) stream-aquifer water exchange, 

(2) evapotranspiration prior to development, (3) drainaquifer water exchange, and (4) recharge from irrigation 

after development, are the primary controlling components of ground-water movement in these basins.

GENERAL MODELING CONSIDERATIONS

The four principal items that are required in the development of a numerical model are as follows:

 Water-level maps that represent the predevelopment 

period and subsequent development periods, so that 

field-determined potential distribution can be compared to the simulation results;

 Definition of the aquifer boundaries, including information to define the size and shape of the basin;

 Area! distribution of aquifer properties, including 

hydraulic conductivity or transmissivity and storage 

characteristics; and

 Definition of quantities and areal distributions of the

various flow components of the water budget of a

basin.

This study developed information and understanding 

of these four essential items throughout the study area. 

The steady-state water-level map for a particular basin 

can be of fair to poor quality depending on the rapidity of 

the initial ground-water development. A representation of 

the steady-state flow system throughout the study area 

was developed as part of this study (Freethey and Anderson, 1986). The process of developing this map included 

the use of all available data and information transfer to 

extrapolate water-level shapes and gradients to areas of 

minimum data.

The steady-state water-level map was used in conjunction with knowledge of the physical boundaries of the

basins to develop an areawide synthesis of boundary-flow 

conditions. Flow-net analyses of individual basins were 

used to develop and refine estimates of location and quantities of boundary flows. Through a process of balancing 

and adjusting the flow components of individual basins 

and the entire study area, a general synthesis of inflow 

and outflow quantities was developed (Freethey and 

Anderson, 1986).

Aquifer properties are a function of the sediments that 

make up the aquifers, which are geologically complex in 

the basins of the study area In developing a model, the 

complexities of the aquifer are generalized through the 

use of average values of hydrologic properties to simulate 

a small block of the system. In this study, that concept has 

been extended to look at generalizations within groups of 

basins. The models represent another step in the development of a complete understanding of the geohydrologic 

systems of the basins. The systems were simulated in 

quasi-three-dimensional form. Previous models had simulated the basins as two-dimensional systems, although 

evidence documents the existence of head variation with 

depth. In most instances, however, the head variations are 

small and spatially variable. Head variations are thought 

to be local and associated with the presence of finegrained facies and faults (Davidson, 1973). In the southeast basins, however, areally extensive confined groundwater conditions have been documented. Major depositional or structural discontinuities in a basin must be considered in development of a three-dimensional model 

because of differences in hydrologic properties. Location 

and extent of the fine- and coarse-grained sediments that 

make up the aquifers have been described in a general 

manner by Freethey and others (1986). These data should 

allow the development of first-approximation threedimensional models of all basins in the study area.

The values of hydrologic factors for use in a model can 

be significantly refined if more field data are available. 

Two issues need to be considered, however (1) the 

accuracy of existing estimates of property values and (2) 

the sensitivity of the model to variations in those values. 

A particular property may not be of equal importance in 

all ground-water flow models. Ranking of the properties 

by importance is dependent on the geohydrologic environment that defines each aquifer system. Properties that 

exhibit the greatest influence on steady-state groundwater flow hydraulic conductivity, natural recharge 

rates, and extent and integrity of confining units are 

most significant in the southeast basins. Properties that 

exhibit the greatest influence on transient ground-water 

flow specific yield, pumping rates, and well location 

are most significant in the central and west basins. Properties that exhibit the greatest influence on shallow watertable conditions evaporation and transpiration, streamaquifer exchange, and irrigation-return flow are most
D74 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

important in those basins in which streamflow is 

perennial.

Information needed to adequately simulate a hydrologic system that has been stressed by pumping or landuse changes caused by development will differ between 

basin groups. When ground-water withdrawals in a basin 

become greater than the natural downvalley flow and 

withdrawal is sufficiently widespread, duplication of 

stress-strain relations, in terms of matching pumpage and 

water-level decline, requires the documentation of actual 

withdrawal rates and locations of pumped wells. In a 

ground-water overdraft environment such as exists in the 

more intensively developed basins, a model may have a 

generally coarse grid with the smallest cells being 1 mi2 

because the data to simulate the dominant stress, which 

is pumping, do not exist in greater detail. Grid sizes of less 

than 1-mi width were used in the models of Parker Valley, 

Mohave basin, upper and lower San Pedro basins, and 

Benson basin, in which the small spacing was required to 

simulate the river alignment. In the models for Harquahala Plain and McMullen Valley, smaller grid spacing was 

used where a pumping center was close to a poorly permeable boundary and image effects greatly influenced the 

size and shape of the cone of depression. In most basins, 

available data did not warrant the use of small grid 

spacing.

Models of the southeast basins need to be multilayered 

to allow for vertical flow through leaky confining layers. 

Significant head differences may exist where fine- and 

coarse-grained deposits are interlayered, and the value of 

interlayer leakance may be smaller than 10~10 (ft/s)/ft. History matching in steady-state and transient simulations is 

greatly affected by the steep hydraulic gradient at the 

basin perimeter; this effect may indicate a vertical component of flow that is not being simulated or the presence of 

relatively thin, steeply inclined, permeable deposits that 

transmit mountain-front recharge into the basin-fill aquifer. Because evapotranspiration is a principal discharge 

mechanism in the southeast basins, the mathematical 

relation that simulates the evapotranspiration function is 

important and will have a great effect on the historymatching process.

Models of the central and west basins need to be multilayered to allow adequate simulation of vertical-flow components. Although one-layer models have been developed 

and calibrated, their limitations must be recognized. History matching for transient simulations is most sensitive 

to pumping rates and specific yield. Variable recharge, 

delayed drainage, and depth-dependent variations in 

hydrologic properties are of major importance in transient simulations because of their variation with time.

Simple one-layer models of the Colorado River basins 

can provide satisfactory results for purposes such as 

those of this study. The process of matching simulated

and measured water levels is sensitive to the evapotranspiration rates and riverbed-leakance value. For these 

basins, water-budget estimates are of equal or greater 

value than a water-table map in providing the data necessary for final model development. The alluvial aquifers in 

the highland basins are considered to be scaled-down versions of the Colorado River group; therefore, properties 

critical to model development should be the same. Specific basin models of this group were not developed 

because of insufficient data, Evapotranspiration and 

riverbed leakance are presumed to be the principal factors to which the model history-matching process will be 

sensitive.

Simulation of the stream alluvium is important in some 

basins, especially in the southeast basins and the Colorado River basins. This upper layer needs to be simulated 

in models of some basins because it is integral to the 

hydrologic processes of evapotranspiration, stream-aquifer exchange, and vertical movement of ground water to 

or from the deeper aquifer being pumped. Transient history matching for a simulation where storage depletion is 

important, such as those of the west basin group, depends 

on accurate values for pumpage and specific yield. Without reliable values for pumpage and specific yield, the 

range of possible matches between simulated and measured water levels is large, and the resulting simulation 

may provide questionable results. Other factors that 

result in time-dependent changes in flow quantities will 

be important in the simulation; these factors include 

recharge from excess applied irrigation water, delayed 

drainage from fine-grained sediments, and the quantity of 

water derived from compaction of fine-grained material.

SUMMARY

Numerical modeling was used as a tool to examine and 

evaluate regional similarities and differences in groundwater flow systems within the alluvial basins that occupy 

82,000 mi2 of south-central Arizona and parts of adjacent 

States. The 72 alluvial basins of the study area represent 

virtually independent hydrologic units but can be categorized into five groups that have similar geohydrologic 

characteristics and similar requirements for model development. The groups southeast, central, west, Colorado 

River, and highland are composed mainly of adjacent 

basins that have similar climatic characteristics and altitude ranges. Models of 12 specific basins and 2 basin 

groups were developed to evaluate the current understanding of the basin flow systems and to explore the 

utility and applicability of general information-transfer 

techniques.

The basins are filled with variable thicknesses of alluvial deposits that constitute the aquifer system of the 

basin. A general sequence of sediments is present in the
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D75

basins and includes pre-Basin and Range rocks, lower and 

upper basin-fill units, and stream alluvium. The sediments 

that fill the basins generally are interconnected hydraulically. The part of the aquifer system that includes the 

principal water-bearing unit in a basin depends on the 

geohydrologic setting, which also influences basin 

categorization.

The basins are interconnected in a dendritic pattern 

much like the surface-drainage system and form an integrated regional-flow system. Individual basin-aquifer systems serve mainly as reservoirs for the storage of ground 

water with only small areas of interconnection and a 

small quantity of flow that occurs from basins of higher 

altitude to those of lower altitude.

Basin categorization is based on the relative importance of various components of the geologic and hydrologic systems. The principal geologic differences between 

groups are the presence and extent of fine-grained sediments in the basin center and the presence of saturated 

stream alluvium. These units represent the least and most 

productive parts of the aquifer, respectively. The principal 

hydrologic differences between groups are (1) magnitude 

of the component of mountain-front recharge and (2) 

presence of a perennial stream in the basin.

Ground-water modeling provides a means of improving 

the understanding of the geohydrology by allowing evaluation of alternative property values and distributions. 

Representative basins from each category were analyzed 

by using numerical models to test the concepts of categorization. Each model of a specific basin was developed on 

the basis of available data and previous conceptual or 

numerical models. Initial property values used in the 

models were based on available field information or, 

where no field data were available, were transferred from 

a basin or basins within the same geohydrologic category 

for which data were available.

Of the 12 basin models developed, all 12 were history 

matched to predevelopment steady-state conditions, and 

10 were history matched to transient development conditions. The history-matching procedure for steady-state 

models depended on the most reliable data available. 

Models of basins having plentiful predevelopment waterlevel information were history matched by reducing the 

standard error between observed and model-calculated 

water levels to an acceptably small value. Modelcalculated ground-water budget values were compared 

with corresponding budget values from previous field 

investigations where the measured values were available.

Hydrologic-property values derived from steady-state 

simulations were used subsequently in the transient models. Transient models were used to simulate water-level 

and flow-component changes for various periods of 

pumping stress; they also served as a means of evaluating 

the consistency of property values determined from

steady-state history matching. The least-known values, 

such as specific yield, were varied within acceptable limits; the best-known values, such as pumpage, were not 

varied.

Among the specific topics investigated by groundwater simulation were (1) magnitude and distribution of 

mountain-front recharge, (2) spatial distribution of 

hydraulic conductivity and specific yield of the aquifer 

material, (3) relation between pumping and capture of 

natural discharge, and (4) recharge from deep percolation 

of excess applied irrigation water. Both specific-basin and 

basin-group models were used to analyze the significance 

of these topics to the local hydrologic setting and their 

effect on system response to development.

Because mountain-front recharge is a major water-budget component in the southeast basins, these models 

were used to test the general relation between precipitation and recharge that resulted from the regional waterbudget analysis of all basins. The specific-basin models 

and basin-group model for the southeast basins were further used to analyze the effect of using a uniform spatial 

distribution of mountain-front recharge. The empirical 

relation between precipitation and recharge proved adequate for estimating the magnitude of mountain-front 

recharge, but because recharge does not occur uniformly 

along the basin margin, this approach to spatially distributing the recharge is not recommended. The quantity of 

mountain-front recharge in the west basins is small; therefore, models are insensitive to the exact magnitude and 

spatial distribution used in the simulation.

All basin models can be used to evaluate the effect of 

poorly defined spatial distributions of hydraulic conductivity and specific yield. Simulation results indicate that 

aquifer properties are highly variable laterally and that 

limited generalization is acceptable for modeling purposes. The relatively large magnitude of water-budget 

components in the southeast and Colorado River basins 

results in significant long-term capture of natural discharge by pumping. The components of evapotranspiration and discharge to springs and streams are most 

affected by capture.

The need to include recharge from deep percolation of 

excess applied irrigation water in simulations is indicated 

in 10 of the 12 specific-basin models. Such recharge 

becomes important after a lag time of a few years to several decades following initiation of irrigation. The overall 

result has been a change in the stress-strain relation over 

time; a decreased rate of decline has been documented in 

recent time in most basins, although the pumping rate has 

remained nearly uniform.

Basins of the southeast group for which ground-water 

flow models were developed are upper San Pedro, Benson, lower San Pedro, Willcox, and Douglas basins. These 

basin aquifers were simulated as two layers connected
D76 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

through a vertical-leakance value that represented the 

effects of anisotropy of the basin-center fine-grained 

facies. The most important properties used in the models 

of the southeast group include specific yield, interlayer 

leakance, the relation between evapotranspiration and 

depth to water, and pumpage quantity and location.

Basins of the central group for which ground-water 

flow models were developed are the Tucson basin and 

Avra Valley. Models of the west basins group were developed for Waterman Wash basin, Harquahala Plain, and 

McMullen Valley. These aquifer systems are most effectively simulated by using multiple layers. Recharge and 

discharge mechanisms are more important in the central 

group than in the west group. Storage depletion is common in both basin groups as a result of development, 

although some capture of discharge occurs in the central 

basins during the early phase of development.

Basins of the Colorado River group for which groundwater flow models were developed are the Parker Valley 

and Mohave basin. Because the hydrologic system of 

these basins is dominated by streamflow in the Colorado 

River, hydrologic properties related to shallow watertable conditions are most important. These include riverbed leakance, the relation between evapotranspiration 

and depth to the water table, and return of excess applied 

irrigation water. Adequate ground-water flow models of 

these basins can be developed by using a two-dimensional 

approach. Because of insufficient data, no specific-basin 

models were developed for the highland group, although 

these basins are thought to be scaled-down versions of 

the Colorado River basins.

Systematic model-sensitivity analyses were done on all 

models to qualitatively evaluate uncertainty in the definition of values of hydrologic properties and the effect on 

model results. A direct relation exists between basin category and sensitivity of the developed model to selected 

hydrologic properties. The specific properties depend on 

the degree of development, the geologic setting, and the 

magnitude of predevelopment downvalley flow. The 

transfer of information among basins of a particular 

group proved generally successful. The properties that 

are important to adequate modeling of the basins can be 

summarized into three groups hydrologic boundaries, 

flow components, and hydrologic factors. Properties in 

the hydrologic-boundaries group represent the physical 

limits of the aquifer system, those in the flow-components 

group include recharge to and discharge from the aquifer, 

and those in the hydrologic-factors group control the 

movement and storage of water in the aquifer.

SELECTED REFERENCES

Anderson, S.R., 1987, Cenozoic stratigraphy and geologic history of the 

Tucson basin, Pima County, Arizona: U.S. Geological Survey WaterResources Investigations Report 87-4190, 20 p.

 1988, Potential for aquifer compaction, land subsidence, and 

earth fissures in the Tucson basin, Pima County, Arizona: U.S. Geological Survey Hydrologic Investigations Atlas HA-713, 3 sheets.

Anderson, T.W, 1972, Electrical-analog analysis of the hydrologic system, Tucson basin, southeastern Arizona: U.S. Geological Survey 

Water-Supply Paper 1939-C, 34 p.

 1983, Implication of deep percolation to ground-water systems in 

south-central Arizona based on numerical-model studies, in Deep 

Percolation Symposium, Scottsdale, Arizona, October 26, 1982, 

Proceedings: Arizona Department of Water Resources Report 4, 

p. 30-40.

 1984, Southwest Alluvial Basins, RASA study An overview, in 

Replogle, J.A., and Renard, K.G., eds., Water today and tomorrow: 

Flagstaff, Arizona, American Society of Civil Engineers Specialty 

Conference Proceedings, p. 606-613.

 1986, Study in southern and central Arizona and parts of adjacent

States, in Sun, Ren Jen, ed., Regional Aquifer-System Analysis 

Program of the U.S. Geological Survey Summary of projects, 

1978-84: U.S. Geological Survey Circular 1002, p. 116-131.

Anderson, T.W, Freethey, G.W, and Tucci, Patrick, 1992, Geohydrology 

and water resources of alluvial basins in south-central Arizona and 

parts of adjacent States: U.S. Geological Survey Professional Paper 

1406-B, 67 p.

Andrews, D.A., 1937, Ground water in Avra-Altar Valley, Arizona: U.S. 

Geological Survey Water-Supply Paper 796-E, p. 163-180.

Bennett, G.D., 1979, Regional ground water systems analysis: Water 

Spectrum, Fall 1979, p. 36-42.

Briggs, PC., 1969, Ground-water conditions in McMullen Valley, Maricopa, Yuma, and Yavapai Counties, Arizona: Arizona State Land 

Department Water-Resources Report 40, 31 p.

Brown, S.G., 1976, Preliminary maps showing ground-water resources in 

the lower Colorado River region, Arizona, Nevada, New Mexico, 

and Utah: U.S. Geological Survey Hydrologic Investigations Atlas 

HA-542, 3 sheets.

Brown, S.G., Davidson, E.S., Kister, L.R., and Thomsen, B.W, 1966, Water 

resources of Fort Huachuca Military Reservation, southeastern Arizona: U.S. Geological Survey Water-Supply Paper 1819-D, 57 p.

Brown, S.G., and Schumann, H.H., 1969, Geohydrology and water utilization in the Willcox basin, Graham and Cochise Counties, Arizona: 

U.S. Geological Survey Water-Supply Paper 1859-F, 32 p.

Brown, S.G., Schumann, H.H., Kister, L.R., and Johnson, P W, 1963, Basic 

ground-water data of the Willcox basin, Graham and Cochise Counties, Arizona: Arizona State Land Department Water-Resources 

Report 14, 93 p.

Bryan, Kirk, Smith, G.E.P, and Waring, G.A., 1934, Ground-water supplies and irrigation in San Pedro Valley, Arizona: U.S. Geological 

Survey open-file report, 167 p.

Coates, D.R., 1952, Douglas basin, Cochise County, in Halpenny and others, eds., Ground water in the Gila River basin and adjacent areas, 

Arizona A summary: U.S. Geological Survey open-file report, 

p. 187-194.

Coates, D.R., and Cushman, R.L., 1955, Geology and ground-water 

resources of the Douglas basin, Arizona, with a section on Chemical quality of the ground water, by J.L. Hatchett: U.S. Geological 

Survey Water-Supply Paper 1354, 56 p.

Cuff, M.K., and Anderson, S.R., 1987, Ground-water conditions in Avra 

Valley, Pima and Pinal Counties, Arizona 1985: U.S. Geological 

Survey Water-Resources Investigations Report 87-4192, scale 

1:125,000, 3 sheets.
SIMULATION OF GROUND-WATER FLOW IN ALLUVIAL BASINS D77

Davidson, E.S., 1973, Geohydrology and water resources of the Tucson 

basin, Arizona: U.S. Geological Survey Water-Supply Paper 1939-E, 

81 p.

Davis, G.H., Green, J.H., Olmsted, F.H., and Brown, D.W, 1959, Groundwater conditions and storage capacity in the San Joaquin Valley, 

California: U.S. Geological Survey Water-Supply Paper 1469, 287 p.

Evans, L.G., and Haimson, J.S., 1982, Swab/RASA aquifer parameter 

study: Arizona Department of Water Resources, Final report for 

U.S. Geological Survey contract 14-08-0001-18268, 21 p.

Fenneman, N.M., 1931, Physiography of western United States: McGrawHill, 534 p.

Finnell, T.L, 1970, Pantano Formation, in U.S. Geological Survey, 

Changes in stratigraphic nomenclature, 1968: U.S. Geological Survey Bulletin 1294-A, p. A35-A36.

Freethey, G.W, 1982, Hydrologic analysis of the upper San Pedro basin 

from the Mexico-United States international boundary to Fairbank, 

Arizona: U.S. Geological Survey Open-File Report 82-752, 64 p.

Freethey, G.W, and Anderson, T.W., 1986, Predevelopment hydrologic 

conditions in the alluvial basins of Arizona and adjacent parts of 

California and New Mexico: U.S. Geological Survey Hydrologic 

Investigations Atlas HA-664, 3 sheets.

Freethey, G.W., Pool, D.R., Anderson, T.W, and Tucci, Patrick, 1986, 

Description and generalized distribution of aquifer materials in the 

alluvial basins of Arizona and adjacent parts of California and New 

Mexico: U.S. Geological Survey Hydrologic Investigations Atlas 

HA-663, 4 sheets.

Graf, C.G., 1980, Maps showing ground-water conditions in the Harquahala Plains area, Maricopa and Yuma Counties, Arizona 1980: Arizona Department of Water Resources Hydrologic Map Series 

Report 1, scale 1:125,000, 3 sheets.

Halpenny, L.C., and others, 1952, Ground water in the Gila River basin 

and adjacent areas, Arizona A summary: U.S. Geological Survey 

open-file report, 224 p.

Harshbarger, J.W., Lewis, D.D., Skibitzke, H.E., Heckler, W.L., and Kister, 

L.R., 1966, Arizona water U.S. Geological Survey Water-Supply 

Paper 1648, 85 p.

Heindl, L.A., 1952a, Lower San Pedro basin, in Halpenny and others, 

eds., Ground water in the Gila River basin and adjacent areas, 

Arizona A summary: U.S. Geological Survey open-file report, 

p. 87-99.

 1952b, Upper San Pedro basin, in Halpenny and others, eds., 

Ground water in the Gila River basin and adjacent areas, Arizona 

A summary: U.S. Geological Survey open-file report, p. 69-4J6.

Hely, A.G., and Peck, E.L., 1964, Precipitation, runoff and water loss in 

the lower Colorado River-Salton Sea area: U.S. Geological Survey 

Professional Paper 486-B, 16 p.

Jones, S.C., 1980, Maps showing ground-water conditions in the lower 

San Pedro basin area, Pinal, Cochise, Pima, and Graham Counties, 

Arizona 1979: U.S. Geological Survey Open-File Report 80-954, 

scale 1:125,000, 2 sheets.

Kam, William, 1964, Geology and ground-water resources of McMullen 

Valley, Maricopa, Yavapai, and Yuma Counties, Arizona: U.S. Geological Survey Water-Supply Paper 1665, 64 p.

Kisser, K.G., and Haimson, J.S., 1981, Estimates of aquifer characteristics using drillers' logs, in Proceedings of 1981 meeting of American Water Resources Association and Hydrology Section: ArizonaNevada Academy of Sciences, May 1-2, 1981, v. 11, p. 112-116.

Kister, L.R., Brown, S.G., Schumann, H.H., and Johnson, PW, 1966, Maps 

showing fluoride content and salinity of ground water in the 

Willcox basin, Graham and Cochise Counties, Arizona: U.S. Geological Survey Hydrologic Investigations Atlas HA-214, 6 p.

Konieczki, A.D., 1980, Maps showing ground-water conditions in the 

upper San Pedro basin area, Pima, Santa Cruz, and Cocliise Counties, Arizona 1978: U.S. Geological Survey Water-Resources Investigations Report 80-1192, scale 1:125,000, 2 sheets.

Mann, L.J., and English, C.S., 1980, Maps showing ground-water conditions in the Douglas basin area, Cochise County, Arizona 1978: 

U.S. Geological Survey Open-File Report 80-700, scale 1:125,000, 3 

sheets.

Mann, L.J., White, N.D., and Wilson, R.P., 1978, Maps showing groundwater conditions in the Willcox area, Cochise and Graham Counties, Arizona 1975: U.S. Geological Survey Water-Resources Investigations Report 78-60, scale 1:125,000, 3 sheets.

Matlock, D.T., 1981, Simulative models for the analysis of ground-water 

flow in Vekol Valley, the Waterman Wash area, and the Bosque area, 

Maricopa and Pinal Counties, Arizona: U.S. Geological Survey 

Open-File Report 82-77, 77 p.

Meinzer, O.E., and Kelton, F.C., 1913, Geology and water resources of 

Sulphur Spring Valley, Arizona, with a section on Agriculture, 

by R.H. Forbes: U.S. Geological Survey Water-Supply Paper 320, 

231 p.

Metzger, D.G., 1957, Geology and ground-water resources of the Harquahala Plain area, Maricopa and Yuma Counties, Arizona: Arizona 

State Land Department Water-Resources Report 3, 40 p.

 1968, The Bouse Formation (Pliocene) of the Parker-BlytheCibola area, Arizona and California, in Geological Survey Research 

1968: U.S. Geological Survey Professional Paper 600-D, p. 126-136.

Metzger, D.G., and Loeltz, O.J., 1973, Geohydrology of the Needles area, 

Arizona, California, and Nevada: U.S. Geological Survey Professional Paper 486-J, 54 p.

Metzger, D.G., Loeltz, O.J., and Irelan, Burdge, 1973, Geohydrology of the 

Parker-Blythe-Cibola area, Arizona and California: U.S. Geological 

Survey Professional Paper 486-G, 130 p.

Moosburner, Otto, 1972, Analysis of the ground-water system by 

electrical-analog model, Avra Valley, Pima and Pinal Counties, 

Arizona: U.S. Geological Survey Hydrologic Investigations Atlas 

HA-215, 2 sheets.

Murphy, B.A., and Hedley, J.D., 1984, Maps showing ground-water conditions in the upper Santa Cruz basin area, Pima, Santa Cruz, Pinal, 

and Cochise Counties, Arizona 1982: Arizona Department of 

Water Resources Hydrologic Map Series Report 11, scale 1:125,000, 

3 sheets.

Oppenheimer, J.M., and Sumner, J.S., 1981, Gravity modeling of the 

basins in the Basin and Range province, Arizona: Arizona Geological Society Digest, v. 13, p. 111-115.

Page, H.G., 1963, Water regimen of the inner valley of the San Pedro 

River near Mammoth, Arizona (a pilot study): U.S. Geological Survey Water-Supply Paper 1669-1,22 p.

Pashley, E.F, Jr., 1966, Structure and stratigraphy of the central, northern, and eastern parts of the Tucson basin, Arizona: U.S. Geological 

Survey open-file report, 273 p.

Pool, D.R., 1984, Aquifer geology of alluvial basins of Arizona, in Replogle, J.A., and Renard, K.G., eds., Water today and tomorrow: 

Flagstaff, Arizona, American Society of Civil Engineers Specialty 

Conference Proceedings, p. 683-690.

 1987, Geohydrology of McMullen Valley, west-central Arizona: 

U.S. Geological Survey Water-Resources Investigations Report 87- 

4140, 51 p.

Remick, W.H., 1981, Maps showing ground-water conditions in the 

McMullen Valley area, Maricopa, Yavapai, and Yuma Counties, 

Arizona 1981: Arizona Department of Water Resources Hydrologic Map Series Report 6, scale 1:125,000,3 sheets.

Roeske, R.H., and Werrell, W.L., 1973, Hydrologic conditions in the San 

Pedro River valley, Arizona: Arizona Water Commission Bulletin 4, 

76 p.

Scarborough, R.B., and Peirce, H.W, 1978, Late Cenozoic basins of Arizona, in Callender, J.F., Wilt, J.C., and Clemmons, R.E., eds, Land of 

Cochise: New Mexico Geological Society Guidebook, 29th Field 

Conference, p. 253-259.
D78 REGIONAL AQUIFER-SYSTEM ANALYSIS SOUTHWEST ALLUVIAL BASINS, ARIZONA AND ADJACENT STATES

Sellers, W.D., and Hill, R.H., eds., 1974, Arizona climate 1931-1972: Tucson, University of Arizona Press, 616 p,

Shafiqullah, M., Damon, RE., Lynch, D.J., Reynolds, S.J., Rehrig, W.A., 

and Raymond, R.H., 1980, K-Ar geochronology and geologic history 

of southwestern Arizona and adjacent areas: Arizona Geological 

Society Digest, v. 12, p. 201-260.

Smith, G.E.R, 1938, The physiography of Arizona valleys and the occurrence of groundwater: University of Arizona, Agricultural Experiment Station Technical Bulletin 77, 91 p.

Theis, C.V., Brown, R.H., and Meyer, R.R., 1963, Estimating the transmissibility of aquifers from the specific capacity of wells, in Bentall, 

Ray, comp., Methods of determining permeability, transmissibility, 

and drawdown: U.S. Geological Survey Water-Supply Paper 1536-1, 

p. 331-340.

Trescott, PC., 1975, Documentation of finite-difference model for simulation of three-dimensional ground-water flow: U.S. Geological Survey Open-File Report 75-438, 32 p.

Trescott, PC., and Larson, S.P, 1976, Documentation of finite-difference 

model for simulation of three-dimensional ground-water flow (supplement to Open-File Report 75-438): U.S. Geological Survey OpenFile Report 76-591, 21 p.

Tucci, Patrick, 1982, Use of three-dimensional model for the analysis of 

the ground-water flow system in Parker Valley, Arizona and California: U.S. Geological Survey Open-File Report 82-1006, 54 p.

University of Arizona, 1965, Normal annual precipitation Normal MaySeptember precipitation 1931-1960, State of Arizona- Tucson, 

University of Arizona map.

U.S. Geological Survey, 1982, Annual summary of ground-water conditions in Arizona, spring 1980 to spring 1981: U.S. Geological Survey 

Open-File Report 82-368, scale 1:1,000,000, 2 sheets.

 1985, Annual summary of ground-water conditions in Arizona, 

spring 1983 to spring 1984: U.S. Geological Survey Open-File 

Report 85-410, scale 1:1,000,000, 2 sheets.

White, N.D., and Childers, Dallas, 1967, Hydrologic conditions in the 

Douglas basin, Cochise County, Arizona: Arizona State Land 

Department Water-Resources Report 30, 26 p.

White, N.D., Stulik, R.S., Morse, E.K., and others, 1963, Ground-water 

conditions in the Rainbow Valley and Waterman Wash areas, Maricopa and Final Counties, Arizona: U.S. Geological Survey 

Water-Supply Paper 1669-F, 50 p.

White, N.D., Mattock, W.G., and Schwalen, H.C., 1966, An appraisal of the 

ground-water resources of Avra and Altar Valleys, Pima County, 

Arizona: Arizona State Land Department Water-Resources Report 

25, 66 p.

Wilkins, D.W., 1986, Geohydrology of the Southwest alluvial basins, 

Regional Aquifer-System Analysis, parts of Colorado, New Mexico, 

and Texas: U.S. Geological Survey Water-Resources Investigations 

Report 84-4224, 61 p.


SELECTED SERIES OF U.S. GEOLOGICAL SURVEY PUBLICATIONS

Periodicals

Larthqnakes & Volcanoes < : S M I d hiit'i 'iithh , 

Preliminary Determination of Lpkenurs IISSIK-"! iMdiii:'-

Technical Books; and Reports

Professional Papers .'IK- nunniv ..ompiehcn,! * e ^:,.n\-,\.

1C ports ot Wkk and LsUHg IliK N s{ .Hid L !> I 'i »l 1 atlCC tO pt O! C'->- . >V ,

scientists 'lEiil esiL'incers Iniiuku aie icpotl- on 'he UsiiK- ''

K'MUncC -HidleS illlii '<'. MpOLliaphn. ll\ dK .. ^U and gcoloy.,

iiivestisjaiMii- ! hev .sh') ii\ hide n'll. clio'is d rri iled ( ap. : 

adcli'.\ssing ddiereiti i'-pcctso! i MP.L:^ scieiintu topu

Bulletins corn un. M_,t it.iani data s;.-i intciprfUUon1 - the:! '. 

ol laMing scientific inter ; hu r ,IK- ivi.ualiv pxiic i-niitul n> v ;\ 

or gcViMaphK i.! -uaue Map Pt.H. s-,,011 ti I'ap- >s |'!H \ mchkK i 1 

IC'SlllN DI Jt-'SOCK C Sll.di 's uild o! VOiOgk and tOpi'L'iaptlk I'i'.v i

Lotions i*. well ,,-, volh-c'ion- o 1 ;,!UKI [;,:[icr- tetiii \l to .1 ^pt.' 1 n^

tOpK

Uater-Supph 1'apers a-c «.'iinfn'C!ivitM\i- iepoits tiur <n--^ f n

i Imestigations Maps are geologic maps on topogtaphic 

; : ;; ! " >. a! vatioiis scale 1- show ins.' be chock 01 surficial

-ii i;u,Taph\. and siiikHiral leialions in certain eoai- 

. i e , r

and («a^ hnestigations Charts show sttatigiaphie intor

 ,vT!.;i!i <M! and L\IS lields and otliei ateas having pciios Field Studies Maps are multicolor or blaekmans on topographic 01 pl.inimetnc bases for (juadran- 

' jiii,.r aieas at ' anoiis -,ealcs, Prc-1^71 maps show bed- 

! o_'.\ in relation ro specific mining or inineral deposit 

p«-q i l > ' j map-- aie piimaiiK black-and-white maps on 

i let! 1 - sut. ii as envuo.smenlal studies 01 wilderness ininH\(ir<»Iogic Imestigalions Atlases are multicolored ot 

r 'i-o v. hiti maps oil i.'pographn. or pianimetrie bases preVj -, uie lange ot geonxdrologic data of both regular and

.i.:; :tc is. pi incipal scale is 1 ;2 l-.OOO, and regional studies are 

>,i ("Hi ^.ilc- 01 smallei

\\atci

Circulars prcsctit jai'iinisiirsri-,.' udiiifii.unni ot m.pi.ri ; ;i, 

scienlslR- u.loniiatu)!! oi v\ >uj niipula; n:LTcs! in a forma! desk;*, . i 

tc<r di^iihuiion at no LOSI 10 the puh'io IniorinalMii is nsijull\ o; 

sluin tcmi interest

\\ater- Resources im estimations Reports au1 papcis ut ai; 

inlcipiclivc nature made axaiLinle k 1 the piih-lk nuisid.- (he ioini.il 

I S(iS pnhiienlioiii sencs. Copies aie lepiodtued on lecjuesi unj;ke 

toimal I'SCiS pubiieations, and the\ ait ;.iKo a\ailab!e tot pnhik 

inspection at depositones u'^liealed in I'SCjS ^ataiogs.

Open-File Reports include unpublished manuscript tvpoiiv 

map.s. and other mateiial that ate made available tor puhliL eon<«il 

tation at depositories 1 liey ate a iionpermaneM! iorm ot pubh. a 

tion that ma\ be uied in othet publications as Mimces of 

intoi (nation.

Maps

Geologic Quadrangle Maps ate multiculor geologic map^ 

on topoiiraphic bases in 7.5- 01 1 .^-minute c|uadrangkj formats 

(scales mainls 1:24.0(10 or l.d2,SOO) showing bedrock, suilieuil. 

or engineeiing geolog\. Maps gencia!l\ include brief texts: some 

maps include structure and columnar sections only

Geophysical Investigations Maps are on topographic or 

planimetric bases at various scales. thc\ shou results of sui\e\s 

using geophysical techniques, such as gtavity, magnetic, seismic. 

or radioactivity \\hich icflecl subsurface structures thai ate ot economic or geologic sigmlicance. Main maps include correlations 

with the geology

Miscellaneous Imestigations Series Maps are on planimeitic or topographic bases of tegular and irregular areas ai various 

scales; they piesent a wide variety ot toimat and subiect mailer. 

The series also includes 7 5-minute quadtangle photogeologic 

maps on planimetric bases that show geolog\ as interpreted from 

aenas photographs. Series also includes maps of Mats and tin 

Moon

Catalogs

 vnaancni catalogs, as well as some others, giving compre- 

!'c;'-v i.snnuN ot ! S Geological Siu\ey publications aie a\r ailat-lt ui.dei ;tie conditions indicated below I'rom the I'S. 

(^,:(.hiL',iLai Sii;ve\, Inlotmation Suvices. Bov 252KO. Federal 

( ; titci. I >en\er. C'O S0225 (See iatesi Price and Availabditv List.)

"Publications of the Geological Survey, 1879-1961" may 

i\ i UH iuised K mail and o\et the counter in papeiback book foim 

I'id ;u a >ct oi nuciotiehe.

'"Publications of the (Geological Survey, 1962-1970" may 

!v pi'tehased bv mail and ovei the counter in paperback book form 

 iiki n> a =.ei *>i microtiche.

'Publications of the L'.S. Geological Survey, 1971-1981" 

m..\ be ['iiichased b\ mail and over the counter in paperback book 

!onr ,iu.i volumes, publiealions listing and index) and as a set of 

ii;k loik'he

Supplements to! 19S2, 1983. 1984, 19H5. 198(v and for subsequent Nears since the last permanent catalog may be purchased 

b\ ma;! and ovei the counter in paperback book form.

Stale catalogs, "List of U.S. Geological Survey Geologic 

and \\ater-Supply Reports and Maps For (State)," may be purchased by mail and over the counter in paperback booklet form 

only

"Price and Availability List of LT .S. Geological Survey 

Publications," issued annually is available i'ree of charge in 

paperback booklet form only

Selected copies of a monthly catalog "New Publications of 

the V.S. Geological Sur\ey" are available free ol charge by mail 

01 mav he obtained ovei the counter in paperback booklet form 

oniv Those wishing a free subscription to the monthh catalog 

"New Publications of the L'.S. Geological Survey" should write to 

the I S. Geological Survey. 582 National Center, Re.ston. VA 

220^2.

Note- - Prices of Government publications listed in older catalogs, announcements, and publications mav be incorrect. Therelore, the prices charged may differ from the prices in catalogs, 

announcements, and publications.