dissertation entitled

PREDICTING GROUNDWATER FLOW AND TRANSPORT USING
MICHIGAN'S STATEWIDE WELLOGIC DATABASE

presented by

Andreanne Simard

has been accepted towards fulfillment
of the requirements for the

Ph.D. degree in Civil and Environmental
Engineering

     

 

Major

firs]; Signature '
H a /

Date

MSU is an affinnative-acfion, equal-opportunity employer

 

 

| LJfTFvEllSle I

 

PLACE IN RETURN Box to remove this checkout from your record.
TO AVOID FINES return on or before date due.
MAY BE RECALLED with earlier due date if requested.

DATE DUE DATE DUE DATE DUE
200.8
D '9

 

PREDICTING GROUNDWATER FLOW AND TRANSPORT USING
MICHIGAN ’S STATEWIDE WELLOGIC DATABASE

By

Andreanne Simard

A DISSERTATION
Submitted to
Michigan State University
In partial fulfillment of the requirements
for the degree of
DOCTOR OF PHILOSOPHY
Department of Civil and Environmental Engineering

2007

 

ABSTRACT

PREDICTING GROUNDWATER FLOW AND TRANSPORT USING
MICHIGAN ’S STATEWIDE WELLOGIC DATABASE

By

Andreanne Simard

In this dissertation we present a new source of data and a statistically-based approach for
mapping long term mean and transient groundwater level and hydraulic conductivity
distributions. In particular, we propose integrating and making systematic use of a
massive amount of highly valuable water well data collected by drillers and accumulated
over the past several decades. Although these measurements can be relatively crude they
are available extensively and at a high resolution that increases with time. Many prior
investigations concluded that the data may be too noisy for most practical hydro geologic
applications, our systematic analysis shows that, if statistically processed, these data can
be extremely useful, enabling surprisingly accurate predictions of spatially-detailed and
temporally representative groundwater flow patterns and hydraulic conductivity
distributions. This finding may have major implications for groundwater research and

professional investigations in Michigan and beyond.

 

DEDICATION

TO MY FAMILY

iii

ACKNOWLEDGEMENTS

I would like to give my sincere thanks to Dr. Shu-Guang Li, my dissertation director and
my advisor, for all his advice and continuous support throughout my graduate studies at

Michigan State University.

I would also like to thank my Ph.D. committee members: Dr. David W. Hyndman, Dr.
Manta S. Phanikumar and Dr. Roger B. Wallace for their constructive advice and

suggestions throughout the course of my research.

This research was made possible through funding from the Michigan Department of
Environmental Quality (MDEQ). An advisory committee was formed within the MDEQ
to provided valuable feedback to this research. The Committee members included; Mr.
Elgar Brown, Dr. Bill Enslin, Mr. Brant Fisher, Mr. Steve Miller, and Mr. Scott Ross.
There continuous support for our research activities during the years is greatly

appreciated.

This research is a collaborative effort of different individuals. It involves developing,
implementing, testing, and applying an integrated statistical-geostatistical approach for
processing the Michigan statewide Wellogic database. The contributions of each are
described below.
0 Dr. Shu-Guang Li - Responsible for overall management and technical
supervision of the project, including conceptual design, algorithmic testing, and

code evaluation.

iv

 

0 Dr. Qun Liu - Responsible for the Interactive Groundwater (IGW) software
implementation. In particular, Dr. Liu has developed, implemented, tested,
debuged, and optimized the interface with the Groundwater Inventory and
Mapping Project (GWIM) database and their seamless integration with the IGW
environment.

0 Mr. Richard Mandle (MDEQ) - Responsible for testing the IGW interface with
respect to wellhead delineation. In addition, he has documented the IGW
interface, capabilities, functionalities and applications.

0 Ms. Andreanne Simard - Responsible for applying the integrated approach to
process the static water level and specific capacity data, predict groundwater flow
and transport in Michigan's aquifers, and delineate Michigan's Wellhead

protection areas.

 

PREFACE

One of the most significant challenges in groundwater investigations is data limitation.
The fact that the subsurface environment is inherently heterogeneous makes this data
limitation problem particularly acute. Among the various data often needed for
groundwater investigations the most critical components are water level and hydraulic
conductivity data. The collection of this data through traditional means (e.g., aquifer
tests, installing monitoring wells. . .etc) can be expensive. This is especially true when
the aquifer in question is highly heterogeneous and when data are needed over a large
area and/or the flow pattern is complex in response to natural heterogeneity, long term

pumping, surface water features, variable topography, recharge, and seasonal fluctuation.

In this research, we investigate a new source of data that has not been fully utilized in the
past to characterize groundwater systems. In particular, we propose integrating and
making systematic use of a/massive amount of highly valuable water well data collected
by drillers and accumulated over the past several decades. The collection of wells in this
database includes public water supply wells, irrigation wells, industrial wells, and

residential wells.
One of the key advantages of this database (e.g., the Wellogic database) is that the data it

contains is available at high-resolution, orders of magnitude higher than that of traditional

measured data. Most importantly, the data are cumulative and their quality improves

vi

 

with time. Of course, the disadvantage lies in that the quality of the data found in the

Wellogic database can be non-uniform, and sometimes highly uncertain.

The critical question facing us in this research is how to cope with the data uncertainties
to be able to extract the underlying “signals” from the highly noisy distributions. We
have be adopted a statistically-based approach for processing noisy water well
measurements to address this research issue. At the core of our statistical approach is:
(l) a sequence of filters that remove systematic data errors and statistical outliers, and (2)
an advanced geostatistical analysis including nonstationary vario gram modeling and
multiscale kriging that smoothes out small-scale random data noises, inconsistencies, and

subscale variabilities.

We have applied our statistically-based approach to map both long-term mean and
transient groundwater water level distributions, and aquifer hetero geneities in hydraulic
conductivity in a way that is spatially detailed across multiple scales (e. g., from county
scale to local site scale). We have combined both of these distributions to predict
groundwater systems. In order to illustrate the potential benefits of using this data source
in groundwater studies we have applied the Wellogic-based flow systems to (I) predict
contaminant transport migrations, (2) develop Wellhead Protection Areas (WHPAs), and

(3) predict flow systems in response to new pumping stresses.

We have evaluated these results using a systematic analysis with traditional measured

data. This analysis includes (1) comparing the water level distributions produced from

vii

using the Wellogic database to delineations produced from traditional data (existing
WHPAs), (2) comparing the water level distributions we obtain using the Wellogic
database with known contamination travel paths where both the source and the
contamination travel path have already been identified, and (3) evaluating the quality of
the temporal water level distributions that can be produced using this data by mapping
transient water level data at different study sites and comparing them to observed head

hydro graphs.

Our results indicate that the use of this data source can be potentially used to solve
complex groundwater problems. Since the Wellogic data are “free”, cumulative, and
available virtually “everywhere”, this finding has major implications for groundwater
research and professional investigations. In particular, the finding implies that if the water
level and hydraulic conductivity data can be systematically integrated we will be able to
accurately map the groundwater water level flow and velocity patterns in a way that is
spatially-detailed without collecting new measurements, drastically reducing the cost of
groundwater management investigations that require knowledge of the long term mean
groundwater flow system. Most importantly, our ability to characterize groundwater
systems can systematically increase with time as more water wells are drilled and the
Wellogic database is expanded. We recommend that the Michigan database be
continuously maintained and the systematic Michigan data integration effort be extended

nationwide and throughout the world.

viii

 

TABLE OF CONTENTS

LIST OF TABLES .............................................................................................................. xi
LIST OF FIGURES .......................................................................................................... xii
LIST OF ABBREVIATIONS ........................................................................................... xv
CHAPTER 1: INTRODUCTION ....................................................................................... 1
1. Motivation and Objectives ....................................................................................... l

2. The Michigan Wellogic Database ........................................................................... 4

3. Prior Investigations ................................................................................................ 10
CHAPTER 2: A STATISTICALLY BASED APPROACH ............................................ l3
1. Identifying Data Uncertainties ............................................................................... 13

2. Coping with Data Uncertainties ............................................................................ 15
2.1. Step 1: Addressing systematic error ............................................................. 16

2.2. Step 2: Removing statistical outliers ............................................................ 18

2. 2. 1 . Probabilistic analysis ............................................................................ 19

2.2.2. Residual analysis ................................................................................... 20

2.3. Step 3: Removing random errors .................................................................. 22
2.3.1. Variogram modeling .............................................................................. 22

2.3.2. Geostatistical interpolation ................................................................... 24

2.3.3. Computational issues ............................................................................. 26
CHAPTER 3: MAPPING WATER LEVEL DISTRIBUTIONS ..................................... 28
1. Mapping Long-Term Mean Water Level Distributions ........................................ 28
1.1. Comparisons with traditional measured data SWL maps .............................. 28

1.2. Comparisons with USGS monitoring wells ................................................... 35

1.3. Comparisons with delineations produced from traditional data .................... 38

1.4. Comparisons with known contamination travel paths ................................... 39

2. Mapping Large Scale Temporal Trends ................................................................ 43
2.1. Large-scale systematic trend ......................................................................... 43
2.1.1. Drawdown cone dynamics for the City of Lansing well field ................ 45

2.1.2. Groundwater drawdown dynamics for Monroe County, Michigan ...... 51
CHAPTER 4: MAPPING HYDRAULIC CONDUCTIVITY ......................................... 56
1. Specific Capacity Based Hydraulic Conductivity ................................................. 56
1.1. Analytical approach ....................................................................................... 56

1.2. Empirical approach ........................................................................................ 58

2. Lithology Based Hydraulic Conductivity .............................................................. 61

3. Evaluation .............................................................................................................. 63
CHAPTER 5: NUMERICAL MODELING OF NET DRAWDOWN ............................ 69
1. Mapping Groundwater Velocity Patterns .............................................................. 69

2. A Numerical Superposition Approach for Modeling of Net Well Drawdown ...... 7O

ix

 

CHAPTER 6: MAPPING WELLHEAD PROTECTION AREAS .................................. 76

1. Review of Traditional Delineation Methods ......................................................... 77
1.1. Fixed-radius method ...................................................................................... 77

1.2. Analytical method .......................................................................................... 78

1.3. Numerical method ......................................................................................... 79

2. Wellogic-Based WHPA Delineation ..................................................................... 81
2.1. Delineation based on mean advection ........................................................... 81

2.2. Delineation based on mean advection and dispersion ................................... 81
CHAPTER 7: SUMMARY .............................................................................................. 88
APPENDICES ................................................................................................................... 94
Appendix A ................................................................................................................... 95
Appendix B .................................................................................................................. 101
REFERENCES CITED ................................................................................................... 104

 

LIST OF TABLES

Table 1: Information collected by drillers during well construction. ................................. 5

xi

 

r

LIST OF FIGURES

Images in this dissertation are presented in color

Figure 1: A typical hydraulic conductivity variogram when the data is limited. ............... 3
Figure 2: Spatial distribution of water wells at state scale, county scale, and site scale. ....8

Figure 3: Spatial distribution and number of water wells (containing usable data) in
Michigan at different times .......................................................................................... 9

Figure 4: Interpolated static water level distributions by Hill-Rowley et al. (2003) in
Tyrone Township, Michigan. .................................................................................... 11

Figure 5: Interpolated static water level distribution by Fleis & Vanderbrink
Engineering, Inc. (2003) at the Augusta Creek site in Kalamazoo County, Michigan.

................................................................................................................................... 12
Figure 6: Spatial distribution of water wells in Michigan. Lefi: wells in glacial aquifer,
and Right: wells in bedrock aquifer ........................................................................... 17
Figure 7: Conceptual representation of a probability density distribution with an anomaly
indicating potential local outliers. ............................................................................. 19
Figure 8: Conceptual representation of the de-trending process. ..................................... 21
Figure 9: A typical variogram obtained from Wellogic data ............................................ 23
Figure 10: Kriging estimates (m) and kriging standard deviation (0) based on
observations at different locations with a nugget. ..................................................... 26
Figure 11: Illustration of the statistical outlier removal process ...................................... 30

Figure 12: Comparison of interpolated static water level distributions using measured
data (left) and Wellogic data (right) at two different study sites in Michigan .......... 32

Figure 13: Plot of predicted heads using Wellogic data versus traditional measured heads
at 2 selected study sites in Michigan. ........................................................................ 33

Figure 14: Long term mean comparison of the Wellogic data vs. the USGS monitoring
well data for 15 different sites in Michigan ............................................................... 35

Figure 15: Comparison of predicted and observed long-term mean water level
distributions as a function of time. ............................................................................ 37

xii

 

Figure 16: Predicted SWL distributions using statistically processed Wellogic data at 6-
wellhead protection study sites in Michigan. ............................................................ 40

Figure 17: Predicted static water level distribution using statistically processed Wellogic
data at the Roto-Finish Superfund Site located in Kalamazoo County, Michigan. ...41

Figure 18: Predicted static water level distribution using statistically processed Wellogic
data at the Village of Schoolcrafi in Kalamazoo, Michigan ...................................... 42

Figure 19: Illustration of a large-scale systematic decline over a 13-year time period. ...44

Figure 20: Groundwater withdrawals and population in Lansing, Michigan from 1910 to
1990 (Luukkonen, 1995). .......................................................................................... 46

Figure 2]: Predicted drawdown cone dynamics in response to the Lansing Wellfield
pumping ..................................................................................................................... 48

Figure 22: Comparison of predicted and observed water level distributions as a function
of time. ....................................................................................................................... 49

Figure 23: Groundwater drawdown dynamics for Monroe County, Michigan. The square
boxes in the upper left plot are locations of USGS monitoring wells ....................... 52

Figure 24: Comparison of predicted and observed water level distributions as a function
of time ........................................................................................................................ 53

Figure 25: Type I well locations in Michigan used to develop the Michigan specific
empirical equation to calculate specific capacity. ..................................................... 59

Figure 26: Correlation between aquifer test and Wellogic-based specific capacity at 300
aquifer test locations in Michigan. ............................................................................ 60

Figure 27: Specific capacity based and lithology based hydraulic conductivity (K) maps

at 4 sites in Michigan ................................................................................................. 65
Figure 28: Predicted velocity flow pattern for the City of Hudson, Michigan ................. 70
Figure 29: A simple example illustrating the numerical superposition approach. ............ 74

Figure 30: Predicting flow systems with response to pumping at 0, 200, and 2000 GPM.
................................................................................................................................... 75

Figure 31: An illustration of the conceptual difference in the WHPAs developed by
different methods ....................................................................................................... 80

xiii

 

Figure 32: WHPA comparison using (1) the Monte Carlo approach (top) and (2) with the
addition of dispersion (bottom) ................................................................................. 84

Figure 33: Comparison of Wellogic-based WHPA delineations (gray area) with existing

delineations based on measured data (black outlined polygon) at selected sites in
Michigan. ................................................................................................................... 85

xiv

 

LIST OF ABBREVIATIONS

GWIM — Groundwater Inventory and Mapping Project

K — Hydraulic conductivity

LTM — Long Term Mean

MDEQ — Michigan Department of Environmental Quality
SWL — Static Water Level

USEPA — United States Environmental Protection Agency
USGS — United States Geological Survey

WHPA - Wellhead Protection Area

XV

CHAPTER 1: INTRODUCTION

1. Motivation and Objectives

The most significant challenge in groundwater investigations is data limitation. The fact
that the subsurface environment is inherently heterogeneous makes this data limitation
problem particularly acute. Among the various data often needed for groundwater
investigations the most critical components are water level and hydraulic conductivity
data. The distribution of water level and hydraulic conductivity dictates the groundwater
speed and direction. Both of these data are needed for most hydrogeologic analyses such
as contaminant transport predictions, remediation design, wellhead delineation, and for
understanding the underlying hydrogeologic system and calibrating the conceptual

aquifer representation.

The collection of groundwater water level measurements through traditional means (e.g.,
installing and measuring monitoring wells) can be expensive, especially when data are
needed over a large area and/or the flow pattern is complex in response to natural
heterogeneity, long term pumping, surface water features, variable topography, recharge,
and seasonal fluctuations. The amount of water level data that is realistically collected is
typically limited which causes uncertainty in data interpretations. Even when there is a
large amount of data, the measurements may only represent a snapshot of the aquifer

condition at the time of data collection. This temporally unrepresentative data may not

be directly used in groundwater management investigations that require knowing long-
term system behavior. In areas where there has been a systematic change in water level
due to human activity, the lack of traditional measured water level data becomes much
more of an issue in groundwater investigations because a large amount of data, including

current and historical data, is needed to delineate spatial and temporal dynamics.

Collecting hydraulic conductivity measurements can also be very expensive because
these distributions are much more complex than water level distributions, and the amount
of data required to characterize their patterns is much larger. Obtaining such a large
amount of data is unfeasible due to measurement costs (e. g., aquifer tests). Although
aquifer tests provide effective parameters, they cannot be used to characterize spatial
variation. Recent research shows that hydraulic conductivity varies across multiple
scales, which is of great importance in transport applications. Typically, it is
prohibitively expensive to map all of the local details with traditional data and the amount
of hydraulic conductivity data typically used in groundwater investigations can only

reflect the local conditions.

The following figure presents a typical variogam that is often used to characterize the
statistical natural variability of hydraulic conductivity when the data is limited (Figure 1).
The variogram in Figure 1 is very noisy (very difficult to identify the spatial structure of
the data) and provides an uncertain structure, which makes the parameters needed to

characterize the variogram (nugget, variance, range) also very uncertain.

1.2

 

 

 

 

 

db
1.0..
d
0.0 i, 4“

o

I i
(is at "‘ II:
I: *
a 0.5.4: at an
> *i

,‘ i
M ——
-t * ‘ .3
l '3
I]. It ‘ *‘ 4* Resolved
.. VarIabIlIty .
. Sill
—————————————————————————————————————— *- —--——-—--—---—
‘fi * Nugget
4000 6000 0M0 10000

 

Separation distance, m

Figure l: A typical hydraulic conductivity variogram when the data is limited.

In this dissertation, we present a new source of data and a new approach for
characterizing groundwater flow systems. In particular, we integrate and make systematic
use of a massive amount of highly valuable data collected from Michigan's water wells
and accumulated over the past several decades. Although these measurements can be
relatively crude they are available extensively and at a resolution that increases with time.
Prior investigations concluded that the data may be too noisy for most practical
hydrogeologic applications, however, our systematic analysis shows that, if statistically
processed, the data can be extremely useful, enabling surprisingly accurate predictions of
spatially-detailed and temporally representative groundwater flow pattern and hydraulic
conductivity distributions. This finding may have major implications for groundwater

research and professional investigations in Michigan and beyond.

2. The Michigan Wellogic Database

Although we all recognize that the lack of sufficient data has often severely hampered our
ability to characterize Michigan’s groundwater, we are not making full use of available
existing information. Figure 2 shows approximately 500, 000 existing water wells in
Michigan — a massive source of potentially highly useful hydrogeologic information.
Although these wells are not monitoring wells, valuable data were collected, logged, and
filed with the Michigan Department of Environmental Quality (MDEQ) at the time when

the wells were drilled.

Table 1 shows key information contained in a typical well log. The “test rate” is the rate
that water is pumped from the well. The “pre-pumping depth to water” is the distance
between land surface and the level to which water has stabilized in the well with no
outside force acting upon it (the well is not pumping water). Subtracting the depth—to-
water from the land surface elevation results in the pre-pumping static water level (SWL).
Test drawdown refers to the difference between the pre-pumping SWL and the final
pumping level measured by the driller at the end of the well test. The ratio of the test rate
and test drawdown is referred to as “specific capacity” — a measure of well efficiency and
aquifer’s ability to transmit water. In general, wells with lesser drawdowns (higher
specific capacity) will be more efficient wells and the aquifer being pumped is more
permeable. Wells with higher or total drawdowns (lower specific capacity) are often very
inefficient, 10w yielding wells and the aquifer being pumped is less transmissive. The

"test time” is simply how long the driller pumped water during the test.

Table 1: Information collected by drillers during well construction.

 

 

 

 

Data Type Comments

Aquifer type Drift, Rock, or Unknown

Well Type Type of well: Heat pump, Household, Industrial, Irrigation,
Test well, Type 1, Type H, or Type 111 well

Well depth Distance between ground elevation and bottom of well (feet)

 

Screen location

Depth where the screen starts and ends (feet)

 

 

 

 

Well diameter The well casing diameter (inches)

Date the well was drilled The date the well was constructed (year, month, day)

Driller ID State of Michigan Water Well Dnllrng Contractor
Registration Number

Depth to water Distance between ground elevation and the pre-pumprng

static water level

 

Pump test type

Method used to develop the well: Air, Bailer, Plunger, Test
Pump, or Unknown

 

Test rate

Rate of water flow when the well was developed (GPM)

 

Test drawdown time

Duration of pumping when the well was developed (hours)

 

 

 

Eithologr cal layer Primary material encountered during well construction
escrrptrons

Them.“ Of each Measurement of the thickness of the stratum (feet)
lithologrcal layer

 

 

These data may be used to characterize Michigan’s groundwater systems. In particular,
the static water level (SWL) data can be used to map groundwater flow direction; the
specific capacity data can be used to estimate hydraulic conductivity; and the lithology

data can be used to map aquifer elevations and characterize aquifer heterogeneity.

Recognizing this potential, the MDEQ recently made a significant investment digitizing

the well logs throughout Michigan and created a statewide groundwater database. This

database is called Wellogic - an Intemet-based data entry program developed to provide
an easy method for water well drilling and pump installation contractors to submit water
well records (http://www.deq.state.mi.us/wellogic/main.html). Electronic well record
submittal satisfies state and county well record submittal requirements. Wellogic contains
approximately 425,000 water well records found within the State of Michigan since 1966,
and although it represents the best available data, it cannot be considered a complete
database of all the wells or well records in existence. Approximately 270,000 records are
actually usable after the inconsistent point locations and bad data values were removed. A
GIS data update utility was recently created to enable converting Wellogic data to a
general GIS format more accessible by other analysis and visualization programs

(GWIM, 2006).

The obvious advantage of these water well data lies in that they are now “free” and
available in virtually any geographic area in Michigan at a resolution orders of magnitude
higher than that of traditional measured data. Additionally, the number of water wells
increases steadily with time and the data accumulates. Over the past 40 years, an average
of 22,000 wells per year have been added to the data records. Figure 2 illustrates both the
extensive and detailed nature of the Wellogic data. Note even at a site-scale (e. g., an area
of a few square miles), there are hundreds to thousands of data points. Figure 3 shows
how the Wellogic data resolution has been rapidly increasing on a statewide scale over
the past four decades. This means that for sites that do not have enough data today for

certain applications, it may no longer be the case 5 or 10 years down the road.

Another significant feature of the water well data is that their quality improves with time.

More and more well drillers are beginning to take advantage of low cost global

positioning system (GPS) technology and reporting more accurate well locations. There

are 24 satellites positioned above the earth’s surface and at any given time; off-the shelf

GPS units are capable of linking to three or more of these and determining locations

within 100 feet or less. Such high precision locations greatly enhance our ability to use

the well report data to determine direction of ground water flow, head gradients, variation

of aquifer properties throughout Michigan.

Of course, the Wellogic data have their significant disadvantages. One potential Achille’s

heel lies in that their quality can be non-uniform, and sometimes highly uncertain.

Representing raw records directly entered by drillers and local health departments; these

“marginal data” were often collected under variable/ inconsistent conditions without

going through proper QA/QC procedures.

The questions that naturally arise are:

>

>

>

How can we use data that we suspect have a high degree of uncertainty?
How reliable are the Wellogic data?

Can we quantify data uncertainty?

Can the sheer quantity of this “marginal” data compensate for its lack of
quality?

Can lots of noisy data be even more useful than limited high quality data?

How can we maximize information extraction from noisy measurements?

 

 

8:5 N.~

  

 

     
    
     

 

m__o>> 003
5560 c9352; E 96 <
280 0:0

m__o>> mmw E
2:80 c2352..
280 3260

 

 

.258 26 new .3an >858

3 one o
...SHI

.x.
..

 

a

23m 88m 3 £33 .1855 no doesnwummu €0an “N 230E

w__m>> 000.000
2QO 99m

> .350 «gotmv 3 50332 5 A83 29%: 035.888 303 833 .«o hon—g: 05 nonznwammu 3.530 ”m 050E
..mo .oEE.

       

 

l1ii°

.82 8.9 .28. 82 82
.. ... \

 

_I 00000?

sue/VI lo JeqiunN

I 00000N

 

 

 

3. Prior Investigations

A number of researchers have attempted to utilize the Wellogic database for groundwater

research and professional investigations.

The recent Groundwater Inventory and Mapping Project (GWIM, 2006) created a
statewide “first-water” map using the SWL information in the Wellogic database. The
map enables predicting the approximate depth at which saturated water may be first
encountered in Michigan — a key parameter in characterizing groundwater vulnerability.
Preliminary feedback from the field shows that the predicted depth is reasonably
accurate. The GWIM report also estimates the hydraulic conductivity at Michigan’s water

wells using lithologies and land system information.

Reeves et al. (2003) recently utilized the Wellogic SWL data to evaluate the effects of
dewatering in Monroe County, Michigan. The results show that the predicted water level
is able to capture the significant drawdown effects caused by the pumping activities in the
region. Wayland et a1. (2002), Boutt et al.(2001), and Hoaglund et al. (2002) made use of
the lithologic information from Michigan’s water wells to aid the understanding and

development of conceptual representations for groundwater flow models.

Hill-Rowley et al. (2003) evaluated the Wellogic data more critically. In particular, they
collected a significant amount of traditional water levels from regular monitoring wells at
a site where there is also high-resolution Wellogic data. A systematic comparison of the

interpolated water level maps fiom the two different sources suggests that the Wellogic

10

SWLs may not be as effective for predicting the groundwater flow direction as the
general water table location (Figure 4). In fact, the comparative results show that the
Wellogic-based head gradient is significantly different from that based on the traditional
measured data and the difference actually increases with the number of Wellogic
measurements used. Hill-Rowley et al. (2003) state that the use of the Wellogic data
yields a “confused pattern of ground water flow” and the SWL map contains “numerous
closed contours, significant distortions, and obvious anomalous features” that are difficult
to reconcile with the site-specific hydrogeology. Hill-Rowley et a1. recommend randomly

sampling a subset of the data as one possible way to improve the quality of SWL

mapping.

 

 

 

 

 

 

 

o 1 Mile o—__i Mile
Figure 4: Interpolated static water level distributions by Hill-Rowley et al. (2003) in Tyrone Township,
Michigan Lefi is based on 857 Wellogic data points and right is based on 67 measured data points. Note
that the Wellogic-based head map is noisy and its corresponding hydraulic gradient is very different from
that based on measured head.

Fleis & Vanderbrink Engineering, Inc. (2003) recently attempted mapping Wellogic
SWL in the context of delineating a Wellhead Protection Area (WHPA) for the City of
Augusta, Michigan. The result again shows that the predicted local head gradient or flow
direction are problematic, even though the average water table location based on the
SWL may be reasonable (Figure 5). If we perform reverse particle tracking on this SWL

surface, the resulting WHPA will be qualitatively different from that obtained from

process-based groundwater modeling using traditional measured data.

 

 

 

 

 

 

 

 

 

 

 

 

Figure 5: Interpolated static water level distribution by F leis & Vanderbrink Engineering, Inc. (2003) at the
Augusta Creek site in Kalamazoo County, Michigan. Left: using Wellogic data. Right: using measured
data points. Note that the Wellogic map again contains numerous closed contours and the predicted
hydraulic gradient is problematic.

CHAPTER 2: A STATISTICALLY BASED APPROACH

In this research we revisit the Wellogic database and address the issues raised in Chapter
1 in a systematic fashion. In particular, we examine the viability and effectiveness of the
Wellogic database for mapping groundwater flow patterns and hydraulic conductivity
distributions, and present a statistically-based approach for processing noisy water well

measurements .

At the core of our statistical approach is: (1) a sequence of filters that remove systematic
data errors and statistical outliers, and (2) an advanced geostatistical analysis including
nonstationary variogram modeling and multiscale kriging that smoothes out small-scale

random data noises, inconsistencies, and subscale variabilities.

The following subsections discuss the methodology used in our statistically-based
procedure. A companion report documents in detail through a case study the step-by-step

procedure to extract and process Wellogic data using IGW (Mandle, 2006).

1. Identifying Data Uncertainties

To properly process the Wellogic data, we must first understand their special
characteristics and identify the sources of uncertainty, which have the greatest influence
on groundwater flow and contaminant transport. At most sites these include uncertainties

due to temporal variation, vertical variation, sub-well spacing variability, well—location

l3

inaccuracy, driller variability, measurement noise and recording errors. These different

sources of uncertainty are briefly discussed below.

Temporal variability. Wells in the Wellogic database were drilled at different times and
the data collected (if time-dependent) only reflect the condition at the time when the
wells were drilled. In this research, we focus on long term mean or large-scale temporal
trends since most groundwater management problems are associated with long time
scales (e.g., years, decades). We treat temporal variability around the large-scale trends as

randomly distributed from well to well.

Vertical variability. Wells in the Wellogic database were screened at different elevations
and possibly in different aquifers (e. g., glacial drift, bedrock) and the data collected (e.g.,
water level or specific capacity) may only reflect the condition at or near the screen
location. In this research we focus on vertically averaged conditions within a selected

aquifer. In this context, vertical variation around the mean is treated as random errors.

Small-scale horizontal variability. Although Wellogic data may potentially allow
mapping more spatial details than traditional measured data, variability at scales much
smaller than the well spacing cannot be resolved. This subscale variability is treated as

random errors.

Well location inaccuracy. The spatial locations of the wells in the Wellogic database are

approximate, especially for those constructed in early times and estimated based on

14

 

address matching. These “horizontal errors” may also translate into “vertical errors” in
the ground elevations, which may in turn propagate into the water level. We treat these

errors, due to inaccurate well location, as randomly-distributed from well to well.

Measurement errors. Lacking proper QA/QC procedures, the Wellogic data are often of
lower quality as compared to the traditional measured data. There can be systematic
inconsistencies, missing and bad values, recording errors and typographic errors, in

addition to the usual random measurement noises.

Driller variability. Wells in the Wellogic database were constructed by different drillers
and the data quality associated with each driller can be different. In this research, we treat

this “driller variability” as random errors.

It should be re-iterated that the term ‘random error’, mentioned above, not only includes
errors such as those caused by driller variability and measurement errors, but also

includes small scale fluctuations around the mean.

2. Coping with Data Uncertainties

To extract “signals” from noisy Wellogic measurements, we adopt a statistically-based

approach that consists of the following three steps:

15

2.1. Step 1: Addressing systematic error

We begin with addressing systematic errors, ensuring that data be extracted at the desired
spatial and temporal scale. We achieve this through the development of a number of

filters that are briefly explained below.

Aquifer code filter. This filter allows extracting data based on the “aquifer code”
information (e.g., drift, rock, or unknown) contained within the Wellogic database. This
filter is needed since the glacial drift and bedrock dynamics can be systematically
different and are best mapped separately. Figure 6 shows, respectively, the water well

distribution in the glacial drift and bedrock aquifers.

Temporal variability filter. This filter allows extracting data based on the “well
construction time” contained in the Wellogic database. This filter is useful for sites where
the flow condition exhibits a systematic temporal trend. This filter can also useful if one
is interested in mapping long term average seasonal dynamics. For example, one can map
the long term average dry (or wet) condition in an area by extracting and analyzing only

the dry (or wet) season data over the entire available data records.

16

0v one o
3.5“”

 

gobs—om x8600 5 £03

”Ema was .0033. Eon—w 5 2.03 ”we; dmeoaz 5 £03 333 00 .8020me Guam ”0 250E

 

l7

Well screen filter. This filter allows extracting data based on the “well screen location”
contained in the Wellogic database. This filter is useful when one is interested in
mapping the flow condition in a sub-aquifer layer or between a specified depth interval of

interest within a thick aquifer.

Driller ID filter. This filter allows extracting data based on the “state of Michigan driller
identification number” contained in the database. This filter is useful when the data
collected by a particular driller are deemed unreliable and must be excluded from the

analysis.

Bad data filter. This filter allows removing special coded numeric values, e.g., -999, 000,

888, 999, etc, representing missing or bad data in the Wellogic database.

2.2. Step 2: Removing statistical outliers

Our next step addresses more subtle errors such as local outliers that do not have an

explicit identifier in the Wellogic database.

Local outliers are observations that do not follow the statistical distribution of the bulk of
the data, and consequently may lead to erroneous results with respect to statistical
analysis. Local outliers can have several detrimental effects on the prediction including
the influence of neighboring values and effects on variogram modeling to be discussed in

the next subsection.

18

 

 

2.2.1. Probabilistic analysis

To detect these statistical outliers we first perform an exploratory probabilistic analysis
of the scatter data values. Potential anomalies or discontinuities in the probability
distribution may indicate the existence of data errors and local outliers. The data values at
which the probability distribution is discontinuous are candidates that may need to be
removed, although the actual removal must be implemented with professional judgment
and after consulting the actual well log to differentiate bad data from potential singular,

extreme heterogeneity/physical discontinuity. This is illustrated in Figure 7.

Probability Density Function

Probability

0.12
1i
0.11 i
0.10“
if
0.09't
4r
0.03:-
0.07"
0%“
0.05"
0.04"
0.03"
002‘:
0.01
0.00-

 

Outliers

"III... 7\:I 9

Variable Value

 

 

 

 

 

 

 

 

 

 

 

 

4 5

Figure 7: Conceptual representation of a probability density distribution with an anomaly indicating
potential local outliers.

l9

2. 2. 2. Residual analysis

We may also detect potential outliers by computing the global standard deviation
characterizing the data variability in the selected area of interest. If the residual for a
particular data value is found incompatible with or much larger than (e.g., more than 3
times) the overall standard deviation, this data value may potentially be an outlier.

In the case that the data being analyzed are nonstationary in space (e. g., water level and
hydraulic conductivity data), they must be detrended before a statistical residual analysis
can be performed (Figure 8). Since the overall trend in the data is not known a priori, and
can be very complex, we perform the outlier analysis based on a local moving window as

described below:

For a selected measurement,

1. Define a local “window” (e.g., in terms of nearest number of data points) around
the selected data point. This window should be large enough so that a meaningful
statistical analysis can be performed and small enough so that a simple function
(e.g., linear) can be used to adequately describe the local trend.

2. Perform simple linear regression to determine the local trend within the selected
window.

3. Detrend the data values within the local window and compute their standard
deviation.

4. Compare the residual of the selected measurement with the standard deviation

computed in step 3. If the residual is substantially larger than (c. g. by more than 3

20

times) the overall standard deviation, the selected measurement value may
potentially be an outlier.

5. Repeat steps 1 through 4 for next measurement.

“I Local linear trend within a
moving window around
point of interest

 

° I 0
0°. :0. ° °

..we. 30%?”
2: o .0 -- :9 o 0 0.:
2.5. ,5. " ‘:..
'0 0. ° °
8 r z
0 ‘ 8
2.3 o ’O O ; . : ‘uv: : .
<3 "3;." 1:" ' '5 If}.
D g .9 . 3
>~. c
I \

'7. Outlier \

I Outliers -—> '
A ' l l l l l l l
300 400 500 600 700 800

Distance, m

Figure 8: Conceptual representation of the de-trending process.

21

2.3. Step 3: Removing random errors

Finally, we address random data errors, inconsistencies, and subscale variabilities

through a systematic geostatistical analysis including variogram modeling and kriging.

2.3.1. Variogram modeling

A vario gram describes the relationship between measurements taken at some distance
apart. It provides a compact statistical characterization of “data errors” and the
underlying spatial correlation structure. A variogram is developed by graphing the
squared difference of pairs of samples separated by a given distance. Typically, the
variogram is curved and starts near zero and increases with the separation distance, but
for the Wellogic data, the variogram value at zero separation can be significantly higher
than zero. This is so because even for measurements taken from very closely “clustered”
wells, the observed values can be quite different, being collected at different times, at
different elevations, and by different drillers, and being influenced by local heterogeneity.
This nonzero variogram value at zero separation distance quantifies the overall data
uncertainty and is referred to as a nugget in geostatistics. This nugget effect must be

properly represented in a Wellogic-based analysis.

Figure 9 (right) illustrates key features of a typical variogram obtained for log
conductivity from the Wellogic data. The variogram increases to a plateau at larger
separation distances. The plateau value is referred to as the sill and reflects the overall

data variability including both data noises and natural variability. The separation distance

22

where the variogram reaches the sill is referred to as the range. If there is spatial
correlation in the data, the range indicates the distance that correlation is expected.
Figure 9 (left) also shows a typical Wellogic water level variogram. In this case, the
distribution is statistically nonstationary and the variogram increases indefinitely with

separation distance.

 

 

 

 

 

 

 

 

 

 

 

 

Static Water Level Hydraulic Conductivity
Variance Variance
100 1.4

l" 0
90: ‘1’" MW
”0 4i-
"), 1-0: v Resolved Variability

.I ------- l r ----------------------
60: 0.8T
50" 0 s” N

" - " u et .
40o 4- _ 99 em
30 _ o_4..

1. Range
20:“ 02..
10» I'

0 I Nugget 0.0 . v . : . at

0 . = = ‘r r = i = 0 5000 10000 15000
0 1000 2000 3000 4000
Separation distance (0). meters Separation distance (h), meters

Figure 9: A typical variogram obtained from Wellogic data. Left: SWL variogram. Right: hydraulic
conductivity variogram.

To actually use the experimental variogram for mapping purposes, we have to fit the
experimental vario gram with a theoretical model. Commonly used theoretical models
include an exponential model, a gaussian model, a power model, and a spherical model.
For large sites, data may potentially exhibit multiple scales of variations and there may be

a need to fit to the experimental variogram with a multiscale theoretical model or a linear

23

combination of the basic variograms with distinct correlation scales and variances. The

equation below provides an example of a multiscaled theoretical variogram:

uh) = is? (1 —e"’”"‘) (1)

Where n is the number of distinct scales, h is separation distance, M, 0,2 are respectively

the correlation scale and variance at the ith scale.

2.3.2. Geostatistical interpolation

Once the variogram is modeled, the data is interpolated using kriging or multiscale
kriging. Kriging is a geostatistical method used to interpolate spatial data to produce a
prediction surface. It estimates values at those locations, which have not been sampled.

The basic equation used in ordinary kriging is as follows:

FURY) = iW-f.‘ (2)

where m is the number of scatter points in the set, 1? is the value of the ith measurement

point, and w,- are weights assigned to f}.

Although Equation 2 is the same as that used in many traditional spatial interpolation
schemes, the weighting factors are chosen differently. Specifically, the weighting factors
in kriging are chosen based on the underlying statistical correlation structure as
represented by a site-specific variogram. Since there are new four equations and three
unknowns, a slack variable, 0, is added to the equation set. These weighting factors are

determined from the following equations:

24

ZWr70+fl= 710,- i=1,2,...n (3a)

i wi = 1 (3b)

where yg, is the variogram value representing spatial correlation between measurements 1'
and j, yio is the variogram value between measurement i and the location at which

interpolation is performed.

Figure 10 illustrates the role of the nugget effect and how it enables filtering out the
small-scale noises, subscale variations, and inconsistencies. Kriging is useful in capturing
the overall trend of the variable. When the nugget effect is modeled, kriging is no longer
an exact interpolator. As illustrated in Figure 10, the influence of the nugget parameter
will tend to produce smoother less variable kriging estimates that show an apparent

discontinuity at the data points.

25

Outliers

SWL

 

-3
Outliers / m o

 

Distance

Figure 10: Kriging estimates (m) and kriging standard deviation (0) based on observations at different
locations with a nugget. In this case, the local outliers are removed in our statistically-based filtering
procedure.

2.3.3. Computational issues

There are many important implementation issues that need to be considered when
applying kriging. One of the most important is computational cost. Solving Equation 3a
requires inverting a full m+1 by m+1 matrix for each grid point where the interpolation is
to be performed. This computation can be expensive when m is large as is almost always
the case when using the Wellogic data. In order to limit the size of the matrices involved,
the kriging calculation is often done with a moving neighborhood encompassing a
selected number of the closest points to the location where the unknown is being
calculated. The neighborhood is chosen to maximize the number of nearest points while

still being computationally efficient. Typically on the order of 100 nearest points is used.

26

It must be recognized that this approximation can lead to physically unrealistic results

when the number of nearest points used is too small.

27

CHAPTER 3: MAPPING WATER LEVEL DISTRIBUTIONS

1. Mapping Long-Term Mean Water Level Distributions

The new data source contains approximately 40 years of water level data, which makes it
ideal in providing long-term mean water level distributions. Evaluating the quality of the
SWL maps that can be produced using the Wellogic database was approached in four
phases: (1) comparisons with traditional measured data SWL maps, (2) comparisons
with United States Geological Survey (USGS) monitoring well data, (3) comparisons
with delineations produced from traditional data, and (4) comparisons with known

contamination travel paths.

1.1. Comparisons with traditional measured data SWL maps

In order to ensure that our evaluation is meaningful and effective, we must identify a
reference site(s) that can be used to establish the “ground truth”. After an extensive
literature review, we found two such sites in Michigan: (1) the Tyrone Township site (36
square miles) in southeastern Michigan (Hill—Rowley, et al., 2003) and (2) the Kalamazoo
site (116 square miles) in southwestern Michigan (Peerless Midwest, 2002). Both sites
were thoroughly characterized based on a significant amount of traditional measured
data. The geologic setting of both sites is typical of Michigan. Unconsolidated glacial

deposits overlying consolidated bedrock formations (Curran et al., 1981).

28

A total of 67 usable traditional SWL field measurements are available in the Tyrone
township site and 42 SWL field measurements are available at the Kalamazoo site.
Surface elevation points from major, hydraulically connected lakes and streams can also
be used as ground water elevations. Figure 12 (left) presents the measured static water
levels manually contoured by previous investigators (Hill-Rowley, et al., 2003; Peerless

Midwest, 2002).

The Wellogic lists 1937 wells completed in the glacial drift for the Tyrone Township site
and 2584 for the Kalamazoo site. Lithologic information taken from the well records for
both sites indicate that the glacial drift contains saturated sand and gravel layers,
separated by “randomly” distributed clay zones. Our evaluation process involves: (1)
processing and mapping wellogic-based SWL using our statistical approach; (2)
comparing the existing field derived SWL map with the Wellogic-based SWL map, and

(3) analyzing qualitatively and quantitatively the errors in the SWL predictions.

In order to illustrate the statistical outlier removal process, we have included a real case,
located approximately in the top ‘A of Tyrone Township in Livingston County, Michigan,
where the removal of statistical outliers was implemented. The initial probability
distribution of the Wellogic data set indicates the presence of a local outlier that does not
have an explicit identifier within the Wellogic database (Figure 11). Once the potential
presence of the local outlier was identified, the data went through a residual analysis
following the steps listed in Section 2.2.2. With the outlier removed, the remaining data

was processed through variogram analysis followed by a geostatistical interpolation to

29

obtain the SWL maps shown below. The location of the outlier removed was identified
and its actual well log was consulted. In this case, the outliers found (e.g., location 1 on

the top SWL map of Figure l 1), was in fact due to driller recording error within its

respective well log.

Probability Density Function
- Before outlier removal -

 

Probability
0.10
‘i

0.090
0.08"

v

0.07"

0

0'05" Potential

0

0.050 Outlier
0.04::
0.03;
0.021
(1012:
am?

 

 

 

 

 

 

 

 

 

 

2332402502802702032903w310320
Variable Value

 

 

 

 

 

Probability Density Function
Probability - After outlier removal -

0TB

 

0.07j»

MB .
0.050
0.ij
0:03
0.02::

Jr
0.014L J
0CD- .

260 270 280 23] 300 310

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Variable Value

Figure 11: Illustration of the statistical outlier removal process. The top row represents (1) the probability
distribution of the Wellogic data set and (2) the final SWL map betore the removal of statistical outliers.
The bottom row represents (1) the probability distribution of the Wellogic data set and (2) the final SWL

map a er the removal of statistical outliers.

30

The results presented in Figure 12 clearly show that when statistically processed, the
Wellogic SWL maps reasonably reflect the field conditions, capturing the most
important flow features at the test sites. Both measured and Wellogic SWL maps

show:

> The groundwater flow direction is generally from southeast to northwest at the
Tyrone township site. Water table elevations ranged from a high of about 1,010
feet in the southeast comer of the township to slightly less than 860 feet in the
northwest comer. Groundwater flow toward North Ore Creek on the west side of
the township was evident.

> For the Kalamazoo site, water table elevations ranged from a high of about 920
feet in the in the northern and southern end of the study area to slightly less than
800 feet in the middle area where the Kalamazoo river is located. Groundwater
flow toward the Kalamazoo River was evident.

> At both sites, the water table exhibits classic features such as curvature around

surface water bodies, and discharge to lakes and streams.

To make a more quantitative comparison between the map produced from field
measurements and the map produced from the Wellogic data, we also directly compared
the water levels at the field measurement locations. In particular, we plotted them against

each other and quantified the residuals statistically. The results are shown in Figure 13.

31

Tyrone Township. Michigan; 67 Measured Points

    
   
   

 

   

 

 

 

 

 

Kalamazoo. Michigan; 42 Measured Points
YJLEQ

  
 

 

 

 

 

   
 
  

 

 

  
 

 

 

5:"; I h! 1;,
.’ few & HIJ '5‘ , _ ,
9 ‘ 3' “~ r. i .. ,I p, - -_ - ./
' ' 5 ~. :: ‘ I. -> M ‘f
:7 . ‘ r . II:_\ t '1 l T
1'. . . , 20 . 9
. fir-i 1m 1 r l
. ‘ _‘ .“ .- P8001. /1I ‘
' -r ; .- “11.“..2-‘ . 5 grfigrn‘p 4
.. on. . « fié‘v-m g v v". (
soo ._¢_,- 5 , . ‘
.. __ H 4 d- ,, ’_ 920/
: ‘ 9.9 -. -. I! 5' ' I 15
"':‘“""’l_' i 7': , -9»99-~u5, 5-
.,-:.¢.....,.l 414 ...i a l _.
' ‘I‘"' I ' ‘ i I .
2.4.x"; "- 1 ~ .I' .
.f, -=- ' I ‘ .
.-n‘,..n— J
‘,. o I...
. 4 ' .
I "“1" i

 

 

 

-=
0 2 Mile

Figure 12: Comparison of interpolated static water level distributions using measured data (left) and
Wellogic data (right) at two different study sites in Michigan. The well locations are denoted as black dots
and the contours represent the head distribution in the area.

32

00030.2 E 8% Exam 0080.8 N E 0000: 0050008 .0002va 0:000.» 800 0.00:0? wEms 0000; 080.0000 00 SE ”2 03w.»

0 .mE0E050005. .0205. 0503 09.050 0.00.. .000;
000? 000 000 00s 00—? 000? 000
. _ . _ _ _ 02 _ _ _ _ .

mad 0 2.3.2.030 5:20:00 mad 0 20.0.0000 5:20:00

II 000

II 000

u ‘aseqerecj oifiouaM eur Suisn peugerqo S|9A91 JereM

 

II 000.
0:00.000800 .00.. .00; 0:00.0an0 .00.. .0003

3:000 00N0E0_mv. 3500 560:3... c. 0:0 QEmEsoh 0:05...

a .mE0E05000s. .0305. 050: 005050 000.. 0.03

000

000

II 000

fil 000w

 

II 00w?

1; ‘eseqereg ogBoueM our Bugsn peugerqo S|8A3'| 1818M

33

Clearly, the Wellogic and measured data are well correlated, with the relationship
following almost the 45-degree line, and the correlation coefficient being 0.98. The mean
error and the mean absolute error in the predicted water levels for the Tyrone Township
site are —7.0 fl and 9.0 ft and for the Kalamazoo site are —6.0 fl and 6.5 ft respectively. A
mean absolute error of 6.5 it indicates that, on average, the estimated SWL is within 6.5
ft of the measured value. The slight negative mean error for higher SWL values indicates

a minor bias toward underpredicting the SWL at higher elevations.

34

1.2. Comparisons with USGS monitoring wells

In order to make a more quantitative comparison between the long term mean (LTM)
water level data from field measurements and Wellogic SWL data, we have also directly
compared the long-term water levels measured at specific USGS monitoring well
locations to the Wellogic based LTM results at those specific USGS well locations. We
performed this type of analysis at 15 different USGS monitoring well locations across the
state where (1) there was enough Wellogic data to get a reasonable LTM SWL map and
(2) the USGS long-term hydrographs did not show signs of a systematic trend. The

results are plotted against each other and are shown in Figure 14.

1200 —1

1000 ——

800 —~

USGS Observation Well LTM Water Level, ft
1
b

 

 

600 I I I I I I
600 800 1000 1200
Wellogic LTM Water Level, ft

Figure 14: Long term mean comparison of the Wellogic data vs. the USGS monitoring well data for 15
different sites in Michigan.

35

Clearly, the Wellogic and USGS measured data are well correlated, with the relationship
following almost the 45-degree line. The reason for the over or underprediction in
Figure 14 can be attributed to the fact that the Wellogic wells in the immediate vicinity of
the USGS well location are screened at different elevations from that of the USGS
monitoring wells. This suggests that it is plausible that the groundwater flow direction is
systematically downward; however, additional study is needed to confirm this
conclusion. This explanation is consistent with our results which include; (1) when the
USGS long-term water levels are higher than the Wellogic based LTM results, the
Wellogic wells are deeper wells within the same aquifer (see Figure 15: Oakland County
Well 1 for an example hydrograph), (2) in areas where the USGS long-term water levels
are lower than the Wellogic based LTM results suggests that the USGS monitoring wells
are deeper wells within the same aquifer (see Figure 15: Oakland County Well 4 for an
example hydrograph) and (3) when both well types have similar well depths, the results

are very similar (see Figure 15: Oakland County Well 5 for an example hydrograph).

36

0.00

New...“ 200m 80% 250
_ _ . _ . _ . _

m
y 3., ? a

.08: mo 003003 0 00 maousflhmav _0>0. .5003 0008 800.000. 09/030 000 0806000 .0 08.000800 ”m. 05me

:a%

mmo

.I 08

1| mam

 

m .03 0.58 90:00
25 28:02, .9 900 .90.. .202, 000:

[Yam

moat.

r

_

0.00

moan: 5:. ..B wthw

_ _ . _ .

 

II ‘IGAG'I 1919M

 

 

0 :03 3:80 90:00

 

 

 

 

0.00
3.05 Sh 833 5:0 20$ $.03
98 . _ _ a . F . 89
T 03. M l 29
m
3
J
.l 1 .l
a
A
B.
l 98 nu « l R?
+
l 08 r. on?

25 00202, .0, 000 .23 502, 000:

F =25 5:80 20:00

2.... 0.00:0; .0> 900 _0>0... L0~0>> wOWD

1; ‘IeAa'I JelBM

37

1.3. Comparisons with delineations produced from traditional data

Although the detailed comparisons presented in the previous sections are truly
encouraging, we are not yet ready to generalize the findings to across the state. To ensure
that the Wellogic data can indeed be effectively used for mapping head patterns in
different hydrogeological settings in Michigan’s glacial aquifers, we performed
additional extensive testing in different parts throughout the state — focusing on the flow

direction.

In particular, we applied our statistical approach to a large numbers of sites that have an
existing WHPA delineated based on measured data and at which a clear flow direction
can be identified based on the orientation of the WHPA. These include sites of different
sizes, hydrogeologic characteristics, and wellogic data densities and distributions. We
overlaid existing WHPAs at the selected sites onto the corresponding Wellogic-based
SWL maps obtained from our statistical approach. The results are presented in Figure16.
Virtually in all cases, the wellogic-based flow direction matches that implied by the
WHPA orientation based on the traditional data. In some cases, the Wellogic-based flow
direction may deviate from that implied by the existing WHPA away fi'om the wells
being delineated, but this may not be the “fault” of the Wellogic data, rather, it is because
the direction derived from limited measured data around the wells may no longer be
accurate away from the wells when the regional flow is nonuniform. For example, at the
August Creek site, the flow direction based on the Wellogic data bends toward the well
connected lake to the north as expected while the flow based on limited measured data

around the wells is unidirectional and ignores the presence of the lake (see Figure 16).

38

1.4. Comparisons with known contamination travel paths

We also performed SWL mapping using the Wellogic data at two contamination sites
where the flow direction is known based on the plume orientation or the known
contaminant travel path. We overlaid the known contaminant migration pathways on the
Wellogic-based SWL maps obtained from our statistical approach. The results are

presented in Figure 17 and 18.

Again, the Wellogic-based flow directions are consistent with those implied by the

contaminant migration pathways.

39

 

  
 

I —. Ex: 2 ll \ l .L‘ J I
P x V, .._. ‘... -- .
II ‘ __ , .....

~ .

     
    

 

 

 

 

. Kalamazoo Conny: Au-usta‘WHPA

       
  
    

 

 

 
 
  

    
 

 

  

. s j
| \ '1-
» .. ..- - . 0 *“fi
'4: x \X‘. ”:13 “1 ”17—5.--
" ‘ "$060k; - 1 : M-
~ \‘. t A . , .. U.-.
‘ '. at»... '_ ‘ '.. ... -
1i. ——. f AVA \‘ _ . ‘épé‘l— : “a

 

 

Figure 16: Predicted SWL distributions using statistically processed Wellogic data at 6-wellhead
protection study sites in Michigan. The areas in gray are Wellhead Protection Areas (WHPAs) delineated
using traditional measured data. In all cases, existing WHPA orientations are consistent with Wellogic-
based flow directions.

40

$03800 .9300 uouufififluoo 8505— 05 0. 00: >00» «.030 BE. 0095.2 S2500
000050.02 E 08002 3% 003.0095 ammfiméaox 05 00 S00 0.3.33 00000005 3.0030080 wflms 0003506 .05. 5.03 0080 08200.5 u: 05m...—

Bsouu 802.5080 +

Sham 0030:3550 ¥

 

 

 

 

     

 

 

 

 

 

 

 

 

 

 

L
\a\ \\
.. . .
\ _._
.._
..
. .
... 22...?
:3 «E
x . mnmfll .
. .
x. . II. .. i.
. ll’ll Ono Ilia ..z , . . L.....ml. ,14//_.,1/.//.//

 

 

41

$03509 _0>0b 00000305000 0305— 05 mm 0:: >00»
x020 2F d0w30.2 .00N080.0M E cflsoonom .«0 0MOE> 05 00 800 030:0? 00000005 .3003350 MES: 00.03506 35. 0803 0580 0000.000...— “M: 0.50...—

o_.== n; o
.“I

 

0c._._0E00 0E3... £000.0050m v6 0m0_._> 'II

42

2. Mapping Large Scale Temporal Trends

Transient water level distributions can be obtained by using the water level data in this
database. The use of this database for transient analysis is most useful in predicting

large-scale systematic trends, which will be discussed in detail below.

2.1. Large-scale systematic trend

Transient evaluations become very pertinent in cases where there is a decline in the water
level due to human activities. This type of analysis is best applied to cases when there is
a systematic decline and where the time segments are large enough to provide a
reasonable transient data evaluation. Data segmentation filters are of great importance
when developing temporal groundwater flow distributions using the Wellogic database.
When evaluating a large-scale systematic decline, the Wellogic data should be carefully
filtered temporally with respect to the well construction date (located on the well record).
The results of this segmentation analysis greatly depend on the amount of available and
usable data in each segment. An increase in the number of time segments allows for

more temporal resolution.

Figure 19 illustrates a systematic decline in water level over a 13 year time period. In
this case it was determined that using a 2-year time segment was appropriate to obtain
enough data to be able to reproduce this systematic water level decline. The data density
will eventually become less of a problem as the database becomes more and more

populated with time.

43

 

600

 

'4» .
4:

I

.-

L: -

I .v
u-————-
—-—-—--

Mean water level
for selected time

I
I
I
I
I
I I I
I I I
O) I I I I
I; I I I
C; I I segment I
0 58° :General _.-: . . I
m 1 trend ',“ I I I
g) I ‘ I I I I
O I 2:. I I I
a I I—IO—I I I I
w . I I ‘30; I I I
«o—a | .~ . ‘ I I
(D 550 ‘ I I I ‘II ”I... I I I
3 . I I '0 ~ 5;. I I I
I: . I I I I I I
o—- 1‘ . I ' I. . . '
— 1 | I I -‘7 I I I
O.) .-' I I I "IT I I I
> I I "T ' I I I
3 : I I . I I
,_ 54° ”‘ I I I I :It I I
a) I, I I I I 3+- I I
w I I . I I I . I I
0 .: I firm I I F——'—|V a I
I I I I I
g , : segmem I I I I. I I
I I I I I ' I I -
I .2043 I I I ‘I‘I Il .;“
I I oi.—
520 I I I I I Htr.‘ - ‘I I
I I I I I .. 14.1.31:
. I I I I I I
I I I I I I
_' I I I I I

 

Year

Figure 19: Illustration of a large-scale systematic decline over a l3-year time period.

We applied our statistical approach to the Wellogic data at two selected test sites - the
Lansing wellfield site and the Monroe county dewatering site, where the water level is
known to exhibit systematic temporal trends in response to large-scale pumping
activities. In the following subsections we briefly discuss the temporal system dynamics
at these sites before presenting our Wellogic-based predictions and a comparison with the
USGS monitoring results. The results show that the four decades (1966 to 2006) of
extensive SWL data in Michigan enable accurate predictions of the long term water level

distributions at most geographic locations in the state (see Figures 21 through 24).

44

2.1.1. Drawdown cone dynamics for the City of Lansing well field

Groundwater is the primary source of water for domestic, commercial, and industrial use
within the Lansing area and most of the water pumped is from the deep productive
bedrock aquifer or the Saginaw formation. The overlying glaciofluvial deposits, formed
by continental glaciers, is withdrawn primarily by domestic wells. Historically, trends in
groundwater withdrawal rates have followed changes in population. In the City of

Lansing, the peak groundwater withdrawal occurred in 1980 (Luukkonen, 1995).

This significant and variable pumping stress creates a substantial drawdown cone of
depression that is continuously evolving over time underneath the Lansing area. Based on
USGS monitoring, the maximum drawdown in the area increased from approximately 50
feet in 1935 to 170 feet in 1967. The water level in the area began to recover the late
19805 in response to population decline in the region and as the Lansing Board of Water

and Light began to dig new wells away from the cone of depression.

In this research, we were interested in if the wellogic data can be used to characterize this
complex space-time drawdown dynamics. The results can be useful to assess potential
adverse impact by the Lansing wellfield operation, including induced downward
contaminant migration, water quality degradation, land subsidence, interference with the
Grand River and its associated lakes and wetlands, and degradation of aquatic

ecosystems.

45

 

Groundwater withdrawal
- I Population

Groundwater
Withdrawals (MGD)

  
 

14L;

_r

HUJJJD

120,000

100,000
80pm

60,000

Population

     

40pm

     

 

 

is; 20,000
1910 1920 1930

I940 I950 I960 I970 I980 1990

Year

Figure 20: Groundwater withdrawals and population in Lansing, Michigan from 1910 to 1990
(Luukkonen, 1995). The groundwater withdrawal rates in the City of Lansing have increased steadily
between 1910 and 1960, however, withdrawals since 1970 have been nearly constant and slowly
decreasing.

Our evaluation was conducted in several steps: (1) segmenting the water data over time in
decades, (2) preparing a static water level elevation map using the Wellogic data
available for each time segment, and (3) comparing the predicted and observed water
level as a function of time at selected USGS monitoring wells. The results are

summarized in Figures 21 and 22.

Figure 21 presents the mean spatial distribution of the head distribution in the bedrock
aquifer in five different time interval: 1966-1969, 1970-1979, 1980-1989, 1990-1999,
2000-2006. As expected, the Lansing wellfield operation has been exerting a significant
stress on the aquifer system, the water level dynamics are significantly altered, and a
substantial drawdown cone of depression is created underneath the Lansing area. The

deepest drawdown occurs in the 1970’s and the drawdown becomes more dispersed after

46

1990 where the water level begins to slowly recover at the pumping center, consistent

qualitatively with known observations and water use pattern.

Figure 22 presents a comparison of the water level as a function of time at a few selected
long-term USGS monitoring wells. The results show that the Wellogic-based prediction
is able to capture the general temporal trends in most cases and matches the water level

magnitude over time reasonably well.

There are some obvious differences between the Wellogic-based and measured
hydrographs. These differences can be attributed to the fact that the water wells are
screened at different elevations from that of the USGS monitoring wells and/or the fact
that the number of wells available for interpolation in different decades is different. The
errors caused by the inconsistency in well density is expected to decrease as more and

more data accumulate in the Wellogic database.

47

 

 

 
 

.73131tp910fi

I960 - I969
670 data points
‘ SWL range

 

 

. V

o 5' . < .
«1,... 5 .
. .pl -
‘ 1 .

I“. . rm
1663 data points
. SWL range
‘. 737 a to 909 ft
33E"

- 2001- 2006 4
I363 data points
SWL range

 

 

« . u e a, c
- ’ g a 766fito.9q4ft
o 3 no. ‘00. goo . /.
o . o g I“ 3 6
o . , o
‘ 0‘ .a ' o $0.6."
.00 1. N
o
.
‘ 820 ’ ,
‘0 . .0
' a ‘ ‘ .
o o o 9 I
. 1 o ‘ ..¢
0 ‘ N . i
f'
‘ 0
i . .
>00 ‘ o
- . ,.. ..‘ . . ‘ g .
I . 860 9" o ‘ , o
0 0 '
o. «a ‘ o
0 “ h " .. ob . go“
‘ . . ‘ . . 3. u
. ‘ ‘ u
c
..z.. .o . . ‘
o ‘ ‘ ‘
g . . . .. , 3.. . *
0
o o t g x a I g 34 a. o.
6 t 3 ,. ‘ .0... ‘
u o o 4
o t . . ‘a . “ 2900" - . ‘

 

48

1970 - 1979

  
  

, 956 data points
' SWL range

 

 

 

Figure 21: Predicted drawdown cone dynamics in
response to the Lansing Wellfreld pumping. The
square boxes in the upper left plot are locations of
USGS monitoring wells. The black dots represent
the water well distributions and the contours
represent the head distribution in the area.

Depth to water below land surface. ft

Quarantine Farm Well (Location 1)

 

20"
3:
3 40
a ’ z, A
W a
W
E “A
E 60‘ 't : ,"
2 _ .7
3 ‘ 3" A:
‘_ 80~ . "-14
.3 5 A .f ;
E l \ 4"?”
o w’
; RM
g loo
8
‘ Measured
120? "T”"r "7’ ' ‘V—"fi"‘fi'* ‘4' . "Y" r ".
”40 I960 1980 2000
Date
45 . . . .
t 1 Capital City Airport Well (Location 3)
0' J
0
g l
3 55_l -
s l
_ k A
E» l V
0
-° 2
E l -
E 657
9 l
5 l
a. .l
a l
Measured
75’L‘fi 7. m—r‘a *—*1 " ' ’ ' ‘i
1950 1970 1990
Date
40 1 Logan Well (Location 5)
l
l f“ A’ \A
80 ‘ g k
A

j A Wellogic
l ‘ Measured
1980
Date

2000

Depth to water below land surface. ft

0 l Maybel Street Well (Location 2)
'= l
3' l
‘2 20 .
3 fig? s.
g ‘ WV. 4'
2 A ” . of»
3 l "- f 4‘ I
2 4o. .3 V: ‘ A
3 l NA} .r'A
5 l
S l
0
°‘ 60
5 i
O.
8 l
‘ Measured
80 L ?‘#Y_TT—_7_“ 7’fi 7; ‘1— ‘ 'T’fi' ”—7—?
1940 1960 1980 2000
Date
0 '1 .
a: r Seymor Well (Location 4)
6 1 . 7
o A
g 4° l w ‘
3 l g ,
1: -
C
'9 30’“ i A
3 _
g A 05,!
n t
5 120. i i
g
.9.
5 160"
g AWcllogic
l ‘ Measured
2001 a a. .7. _ _ . 0%
1960 1980 2000
Date
20fi

 

 

Figure 22: Comparison of predicted and observed water level distributions as a function of time. (See
locations in Figure 21).

49

Depth to water below land surface. ft Depth to water below land surface. ft

Depth to water below land surface, ft

0 .1 Scott Park & Townsend Well (Location 7)

-l
l

20 ‘ t
. r ’1
f 7 A
40
'i i 5-
60 13
“AVA 31i-
W'il’ A
601
l AWeIlogic
i ‘ Measured
100’#'—‘ j#'mr*"*V' —"l 4 r* 1’ ‘—r
1960 1980 2000
Date

Fenner Arboretum Well (Location 9)

 

A Wellogic
' Marsured

 
 

 

7'—_“Y—‘T"W H "7"“ " "“7 7A A]
1960 I980 2000
Date
0 7 Meridian Township Well (Location 11)
A ‘3 xrr 23$.”
20 ,1
4o 7
60
80
A Wellogic
J ' Measured
100 t- 7 t—T_ a" .1 t -
1960 1980 2000
Date

 

a: 0 a P-5 Well (Location 8)
a 20 ’ ; TA
‘5 ‘ I?! A
E t.
E: ~. 3
8 A .
52 l M .‘
g 604. \ ‘2 7 , 1’ A
.9 f; r . ‘1;
g. 1 j! -;
3 60—, .

3 A Wellogic

‘ Measured
lmfl'y—I ’7" .' "'I'f‘W "'1 7T IVTT’W‘I‘Hjin' #V "7 ‘I
1930 1950 1970 1990 2010

Date

Spartan Village Well (Location 10)

Depth to water below land surface. ft

 

 

 

50

2.1.2. Groundwater drawdown dynamics for Monroe County, Michigan

During the past decade, groundwater levels have declined in bedrocks aquifers
throughout most of Monroe County. The declines are greatest in the northwestern
townships. From 1991 to 2001, ground-water levels declined 10 feet or more in 17 USGS
observation wells. The largest use is quarry-dewatering operations. The amount of
groundwater withdrawn by quarries has doubled during the past decade and is about 75

percent of the total ground-water use in the County (Nicholas et al., 2001 ).

Figure 23 presents the mean spatial distribution of the head distribution in the bedrock
aquifer in five different time intervals: 1991-1993, 1994-1996, 1997-2000, and 2001-
2003. As expected, there exists a number significant drawdown cones of depression in
Monroe county. Each one of them corresponds to a major quarry dewatering operation.
The deepest drawdown occurs in 2001-2003, consistent qualitatively with known

observations and water use patterns.

Figure 24 presents a comparison of the predicted water level using Wellogic as a function
of time at a few selected long-term USGS monitoring wells. The results again show that
the Wellogic-based predictions are able to capture the general temporal trends in most
cases and match the water level magnitude over time reasonably well. The difference
observed in some of the monitoring wells may be again attributed to the elevation
difference between the Wellogic and monitoring wells and the significant vertical head

variation that may be expected at this site during major pumping events.

51

 

 

 

 

 

 

 

 

«-~ . ',' 1994to1996
1991to1993 ' 213datapoints
342 data points . : €570-

   

 

 

 

 

 

SWI.renge541fltoSStm

SIM. range 535 ft to 653 ft . . ..

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1997to2000 . . 2001to2003
236 data points ‘ - ‘ 467 data points
SMrengeSflfttoSfl ’ “ : Mrsnge483ftto$83fl
l:—
0 6 miles

Figure 23: Groundwater drawdown dynamics for Monroe County, Michigan. The square boxes in the
upper left plot are locations of USGS monitoring wells. The black dots represent the water well
distributions and the contours represent the head distribution in the area.

52

Well G-3 (Location 1)

 

 

 

 

 

 

 

 

 

 

620 Well 64 (Location 2)

 

 

 

 

 

 

 

6007
T; 590- AWellogic g AWellogic
a w 3
a) '... -4": .. > A:
> 580 ‘: 4' O l " l we .r.' .
.3. ~ , . a "A"- * A
'33, A A a 580. F «- -
_ Q)
5 560 2
‘6 A ‘° 560'
E .. S
550, _ __ _ A -.f- m. _r A -_.I "Mr v 1——~r—— —r— ‘1 ' 'T ‘: " ”r* 'i“ 7* 1’ r— '* ’6’ fl *‘ ‘
- l I ‘ l 1991 1993 I995 I997 1999 2001 2003
1991 1993 1995 1997 1999 2001 2003
Date Date
640” . 5 -
Well G-5 (Location 3) _ 640j Well 6'6 (Location 4)
0)) . 7‘: “1 AWelIogic
it; 6201 AWellogrc g 6204 ‘ ‘ Measured
a) Measured g: r
m -_ o I
Q) .0 ‘- .
E A ‘ 3 l A «A
as 600-. «A § 600“; L-e— --
§ "" --. .E ' “T“ .-.,-.
E J A, '9 E A A
_ Q)
a) = — 580 a
E 5801 A 2 f
s. (U T'
2 E .-
E 560;
560...“ ..-A ..A AA, A _ AMA, ---.__.! 4‘ ‘rr r*— I an * *"*"' "Ma '1
‘ l l l 1991 I993 1995 1997 1999 2001 2003
I 991 1993 1995 1997 1999 2001 2003
Date Date
a; 640 A Well (3-7 (Location 5) 640 _: “(6119-8 (Location 6)
.92 i
g l J
g l AWeIlogic I g 620 ' AWeIIogic
g 6201 . Measured 6 Measured
% ' ' 0-... 0 A
-o-l J (D I . .- ..
8 i A 9 '
~0- l . j
s 6007 § 600‘; r.
3 f A .. -..--
E ..II “8 l A A7... ..
E2, 530} A ..’*-';»n .5 580; _ ..
(U I .6 .l
3 1 A E + A
j .2 560i
560; E:
'-_---_A-._]._-.-_AA-- -.-...---_-, . A-“ A“! .._-_-- ~- .,--
l f 1 , .
1991 I993 I995 I997 1.999 2001 2003 I991 1993 1995 1997 1999 2001 2003
Date Date

Figure 24: Comparison of predicted and observed water level distributions as a function of time. (See
locations in Figure 23).

53

Water level in feet above sea level

Water level in feet above sea level

Water level in feet above sea level

Well G-9 (Location 7)

 

 

 

 

 

 

 

 

 

 

 

640-..?
-l AWellogic
. 1 Measured
620 I A ——
l
l ——_
600—: A
1' A
1 A
580‘.1
4
5601
l
1”" 3" 1‘ 1” r "1' *‘1 A
1991 1993 1995 1997 1999 2001 2003
Date
650 .
Well G-12 (Location 9)
1; ' AWellogic
640 . .
- ~ Measured
630 A
A
6207
A
A
610 .....__,.._...22. AT A. , ALTA
I991 1993 1995 I997 I999 200] 2003
Date
650
1 Well G-20 (Location 11)
640
630- _
A A A
A
620
AWellogic
6101 A 2. .__. 1 M32579“?
I991 I993 I995 I997 1999 2001 2003
Date

 

Water level in feet above sea level

Water level in feet above sea level

Water level in feet above sea level

640“ Well 5.1] (Location 8)

1 A Wellogic
630
A " Measured

 

1l

5901-
2001 2003

1991 1993 1995' 199? 1999
Date

610] well 0-18 (Location 10)

 

 

 

 

 

600‘ A
l .
-.I _ A , " A -.1
590.1 A '
1
580 1
.-| AWellogic
1 . ’ Measured
570 3' _f’ T ""‘Y T T '_'_T __T_‘ I-'T——T fl
1991 I993 I995 I997 I999 200] 2003
Date
650'
Well G-2l (Location 12)
MOI A Wellogic
l . ’ Measured
630 A’
. “A
620
A
- A
610
600- I 1 1 1 I l :
199I I993 I995 I997 1999 2001 2003

Date

Water level in feet above sea level

Water level in feet above sea level

Water level in feet above sea level

 

 

66% Well G-22 (Location 13)
‘ AWellogic
' easur
6407 ’ :1; V
A
A .. a? ,
620. A A
600 ”"11""? ’ W ”" l"""7’ *
199! I993 I995 1997 1999 2001 2003
Date
1 Well G-24 (Location 15)
6501
<1 . AWellogic
‘ . ' " Measured
640 ‘ ‘
630al A . ‘
‘ ‘. -
’l
l A A
620"1‘
T' l‘”#?*r"re‘7' fi'l"j—‘ *fi
199! 1993 1995 1997 1999 2001 2003
Date
“’0 Well G-26 (Location 17)
650 A Wellogic
Measured
640 V ; ‘
‘22: ..
: . r ..,
630 i ‘ f;
"3" ‘ f“
_ A A
620’ h ‘ A
6’0"“ "1' l ' """WTV‘?
I991 1993 1995 I997 1999 2001 2003
Date

 

660* Well G-23 (Location 14)
g AWellogic
E ‘ Measured
8 640. _;
g3 -
o A -
a ~ e- A
a A , _
.2
E
a 620. A A
>
2
E
*5 _
3
600 .--L.‘ VA A A--. A r ..A‘
1991 I993 I995 1997 1999 2001 2003
Date
Well G-25 (Location 16)
_ 650-
a, l
5
E AWellogic
8 "Measured
"g 640< _
0
g ‘ ..
'5 7 A .
e . —.. A
5 630. A -
T, l
>
33
'3 l A A
E 620*
”’"r—T 171’ f T'fi‘ f 7 fii‘rr \
1991 1993 1995 1997 I999 200] 2003
Date
(:0 .
“ 1 Well G-27 (Location18)
a l
5 l
§ 640“ ‘ ‘
2 l A
O l A
tu l
16 6303 - . A
£1 l
.s ‘ A
E l
>
%’ will
0.)
(3‘3 l
‘
610 :l .1 Measured
I991 I993 I995 I997 1999 2001 2003

55

Date

CHAPTER 4: MAPPING HYDRAULIC CONDUCTIVITY

The Wellogic database can potentially provide us with 2 types of hydraulic conductivity
estimates, which are; (1) specific capacity based hydraulic conductivity and (2) lithology
based hydraulic conductivity. Both of these types of hydraulic conductivity estimates are

discussed below:

1. Specific Capacity Based Hydraulic Conductivity

There are analytical and empirical methods that can be used to estimate hydraulic
conductivity or transmissivity from specific-capacity data (for example, Thomasson and
others, 1960; Theis, 1963; Brown, 1963; Razack and Huntley, 1991; Huntley and others,
1992; El—Naqa, 1994; Mace, 1997). These techniques have been successfully used in the
Cretaceous sandstone aquifers of North Central Texas (Mace and others, 1994), the
Edwards aquifer (Hovorka and others, 1995, 1998; Mace, 1995), the Ogallala aquifer
(Myers, 1969; Mullican and others, 1997), and the Hill Country Trinity aquifer (Mace, in
prep). Prudic (1991) used specific-capacity data in his regional study of the Gulf Coast

regional aquifer systems.

1.1. Analytical approach

The Wellogic database includes specific capacity data for each well (a specific capacity

test is performed at time of installation of each well) and can be used to calculate

56

hydraulic conductivity. The analytical approach relates aquifer transmissivity to specific

capacity through the following logarithmic approximation of the classical Theis equation:

 

Tzzélog 2225“ _Q_ (4)
47: rwS s

where t is the pump test duration (time when the drawdown was measured), rw is the well
radius in the screened interval, and S is the storativity of the aquifer. This equation
assumes (1) a fully-penetrating well; (2) a homogeneous, isotropic porous media; (3)
negligible well loss; (4) and an effective radius equal to the radius of the production well

(Walton, 1970).

Since Wellogic contains data for test rate, drawdown, time, and well radius for most
water wells, we can estimate transmissivity at the well location if one knows the storage
coefficient. Typically S is on the order of 104 for a confined aquifer and 10" for an
unconfined aquifer, but there is no site-specific information on its values in the database.
Fortunately, the estimated T depends on the storage coefficient logarithmically and is not
sensitive to S when the argument in the logarithmic function is large (when r approaches
zero or is evaluated at well radius). Therefore we adopted these typical values for storage

coefficient in our estimation of T.

To estimate the hydraulic conductivity (K), the following simple relation can be used:

K = (5)

Dol'fl

57

where B is the aquifer thickness or effective aquifer thickness. Although most prior
research assume that the effective aquifer thickness B is to equal to the length of screen
interval for small capacity wells and the actual aquifer thickness for high capacity wells,
we let B be a linear function of the test pumping rate, varying linearly between the screen

interval and the entire aquifer thickness, as described by the following equation:

rbl when Q < Qmin

B:<bl+_g:gmi_n__(b2_bl) (6)
Qmax - Qmin
lb, when Q > me

 

where b; is the screen interval and b2 is the aquifer thickness. Note that Qmm is assigned a

value of 10 GPM and Qmax is assigned a value of 500 GPM.

Equations 4 through 6 enables computing hydraulic conductivity throughout Michigan at
water wells that contain information on test rate, drawdown, and time. Although this
equation is not explicit in K, it can be easily solved iteratively and the solution typically

converges in a few iterations.

1.2. Empirical approach

In this study, we expect a reasonable transmissivity spatial pattern due to the high data
resolution however; there is still a need to calibrate transmissivity in order to improve on
its magnitude. We have developed a Michigan specific empirical relationship by linearly
relating aquifer test-based log-transformed hydraulic conductivity to log-transformed
specific capacity calculated for the same well. To define an empirical relationship we

made use of the approximately 300 glacial drift water wells across the state that have

58

both aquifer test hydraulic conductivity and specific-capacity data, log-transformed
values of each parameter, plotted them against each other, and fit a line through the data
using least squares regression (Figure 26). The Type 1 Wells used in the development of

this equation are displayed in Figure 25.

 

_:Niles
0 1020 40

Figure 25: Type I well locations in Michigan used to develop the Michigan specific empirical equation to
calculate specific capacity.

59

 

”
'f
0

Log 10 Transmissivity

 

78% Confidence Interval

98.4% Confidence Interval

 

 

—

 

00 1 2 3 4 5

Log 10 Specific Capacity

Figure 26: Correlation between aquifer test and Wellogic-based specific capacity at 300 aquifer test
locations in Michigan.

The best-fit line through the data is:
0.58 (7)
T = 1936(2)

s
The relationship has an 80 percent prediction interval that spans approximately a half an
order of magnitude. The prediction interval means that we are about 80 percent confidant
that an estimate of conductivity for any given value of specific capacity is within a factor

of 5 of the estimate.

60

2. Lithology Based Hydraulic Conductivity

Estimation of the horizontal hydraulic conductivity using lithology provides an
evaluation of the entire usable Wellogic database lithologic data. The method used in
the State of Michigan Public Act 148 Groundwater Inventory and Map (GWIM) Project
(March 2006) involves assigning estimated hydraulic conductivity values to each
lithology reported on the well record (see Table B1 in Appendix B). The equivalent
hydraulic conductivities were estimated for the Wellogic records using an assigned value
for hydraulic conductivity for each lithology. This assigned value was determined using a
textbook range of values that had been adopted by the Source Water Assessment Program
(MDEQ, 2004). The mean of the log hydraulic conductivity values of the range was
assigned to each lithology. The ranges and assigned values are given in Appendix B

(Table B1).

In general, the equivalent horizontal hydraulic conductivity is controlled by the layer with
the highest hydraulic conductivity. For a system with horizontal layers, the horizontal

hydraulic conductivity is defined as (Freeze and Cherry, 1979):

_ 219.3,. (8)

where B is the total thickness of the unit, B, is the thickness of layer i, K, is the hydraulic
conductivity of layer i, n is the total number of layers in the unit and Kh is the equivalent

horizontal hydraulic conductivity (flow parallel to the layers).

61

In this analysis, the estimates of equivalent hydraulic conductivity were made using the
lithologies reported from the bottom of the well screen, or the top of rock for wells
completed in bedrock underlying the glacial deposits, up to the static water level reported
for the well. This choice limits the equivalent hydraulic conductivity estimate to the
material within the saturated thickness of the aquifer. Glacial deposits exist below the
level of well screens in most places in Michigan. Obviously, these lithologies are not
incorporated into the analysis because they are not present on the water-well record.
Only the used portion of the glacial deposits is evaluated in this analysis. This limitation
may be important in areas of the state where the glacial deposits are thick and most wells
only penetrate the top 100 to 200 ft. A benefit to this data is that it provides a great

vertical resolution (3-D variability) (GWIM Report, 2006).

62

3. Evaluation

Specific capacity based and lithology based hydraulic conductivity (K) results are
summarized in Figure 27 for 4 sites in Michigan. This figure compares both the lithology
and specific capacity based hydraulic conductivity probability distributions, variograms,

and maps for the City of Niles, Livingston County (2 sites), and Kent County.

Since the geostatistical part of our analysis implicitly assumes that the fluctuation around
the mean trend being processed is approximately normally distributed, we computed the
probability distribution of the lithology and specific capacity based hydraulic
conductivity values as part of our exploratory data analysis for the various test sites.
Figure 27 (row 1) presents the probability distribution function of log hydraulic
conductivity for different sites. These results show that log hydraulic conductivity for all

tests appears to be approximately normally distributed.

As an essential component of our hydraulic conductivity evaluation, we performed
variogram modeling at the various sites. Figure 27 (row 2) presents the experimental and
the fitted theoretical variograms for each site. The variograms show the experimental
values as dots and theoretical model as a line. The results reveal well-defined spatial
correlation structures for most sites. This in itself shows that the hydraulic conductivity
data are not dominated by spatially uncorrelated random noise. In particular, the
variograms show a decrease in variance for smaller separation distances indicating spatial

continuity. The variograms also have relatively large nuggets as expected, suggesting a

63

large amount of randomness due to local-scale heterogeneity and measurement errors,
inconsistencies and driller variability. The range, or the distance within which a
parameter is spatially correlated, is about 2,300 to 15,000 ft for hydraulic conductivity at
the sites tested. This means that hydraulic conductivity values are similar to other values
within about 2,300 to 15,000 ft. Theoretical semivariograms, exponential or gaussian
variograms with a nugget effect were visually fit to the experimental data. The fitted
theoretical variograms play an essential role in mapping the Spatial conductivity

distribution.

Finally, Figure 27 (row 3) shows the spatial distribution of lithology and specific capacity
based hydraulic conductivity distributions at selected sites. It should be noted that the
hydraulic conductivity estimated using the Wellogic specific capacity data are
significantly higher than those estimated by GWIM (2006) based on Wellogic lithologic
descriptions and slightly smaller than those estimated from standard aquifer tests. The

difference is probably due largely to the type and purpose of the well tested.

Aquifer tests are generally performed in the higher yielding municipal wells, while
Wellogic contains significant number of private wells and they do not require large yields
to supply a household and are usually completed when the desired yield is reached during
drilling. Consequently, private wells are usually screened in shallower water-bearing
zones and rarely penetrate the entire aquifer unit. The lithology-based K (GWIM, 2006)
was derived in the context of providing a conservative estimate of the aquifer yield in

Michigan.

64

Lithology K, City of Niles Histogram Specific Capacity K, City of Niles Histogram

 

 

 

 

 

 

 

 

 

 

    

Frequency Frequency
200
mu I I I
2 . 3 4 3 4 5 6
Var'ab‘e Va'ue Variable Value
Range = 7196 ft, Nugget = 0.34. Variance = 0.18 Range = 1.394 ft. Nugget = 0.29. Variance = 0.12
V . . an ...n a. ' I.
‘ I
0.5 0.4 ‘ . ‘ . I u 4- a
..I
0.;
0.3
0.3
0.2
0.2
0.1 0.1

 

 

 

 

 

 

U 2000 4003 5000 3000 10000 12000 14000 15000 0 11:00 21:00 3000 “In 5WD 6030 7003

Lag Distance, meters Lag Distance, meters

City Of Niles, Lithology Based K City of Niles, Specific Capacity Based K

169
n/day

370
ftlday

| 07
R1 day

   

20
fl/day
:— :—

6 miles 6 miles

Figure 27: Specific capacity based and lithology based hydraulic conductivity (K) maps at 4 sites in
Michigan. This figure Lithology Based K maps and Specific Capacity based K maps. Both types of maps
are plotted on a log scale. The black dots represent water well locations used to develop these maps.

65

Lithology K. Livingston County 1 Histogram Specific capacity K, Livingston County 1

 

 

 

 

 

 

 

Frequency Frequency Histogram
1m 14g
160
120
140
120 100
1m 00
90 so
80
40
w I I I I I I II
20
m I I II
c.ull - _ Illl ,. III ,.. . . -
1 2 3 4 2 3 4 5 8
Variable Value Variable Value

  
 
   

Range 6924 ft, Nugget = 0.45. Variance = 0.27 n n Range = 6686 ft. Nugget = 0. 45, Variance = 0.31

 

..ae-a.

  
 
   

 

 

 

 

 

‘1'...

0.4
0 3
0.2
0.1 0 2

0.1
wJ0 1000 2000 3&10 4000 5000 5000 3 30 1000 2C00 3000 4000 5030 60m

Lag Distance. meters Lag Distance, meters

Livingston County 1, Lithology Based K Livingston County 1, Specific Capacity Based K

96 548
ft/da ft/day

22 40
ft/day ftlday
:—
2 miles I:—
3 miles

66

 

Lithology K, Livingston County 2 Specific Capacity K, Livingston County 2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Frequency Hlstogram quency Hlstogram
1200 600
1000 500
900 400
600 300
400 200
u- c- 3 I .I I I . I II _
1 2 3 t ' 3 4 5 6
Variable Value Variable Value
Range = 1.194 ft, Nugget = 0.45. Variance = 0.2 Range = 2287 ft, Nugget = 0. 34. Variance = 0.15
05 II
0.4 I
I
0.3
0 3
0.2 0 2
0.1 01
“"0 2000 41300 5000 8000 10000 12000 14000 “Aug 51]] 10m 1500 20m 2500
Lag Distance, meters Lag Distance, meters
Livingston County 2. Lithology Based K Livingston County 2, Specific Capacity Based K

61 7
ft/day

 

13 56
ft/day ftlday
:— l:—
4.8 miles 4.8 miles

Lithology K, Kent County Histogram Specific Capacity K, Kent County Histogram

 

 

 

 

 

 

 

 

Frequency Frequency
5m
:10
‘m .
am 200
200
1m
.00 I
0_ a —--I 7., Will-
1 2 3 ‘ 2 3 4 5 5
Variable Value Variable Value
Range = 1.5e4 ft. Nugget = 0.53. Variance = 0.37 Range = 6792 ft, Nugget = 0. 36, Variance = 0.38
1.0
0 9 - n..-“ a...
0.8 ...Ifl ' -
0 7 ‘
0,5
0 5
0 4
0 3
0.2
0,2
01 01

 

 

 

 

 

0 211]] 4&10 5000 8000 10000 12000 14000 18000 0 2WD 4011) 6WD 8000 10000 121130 14000 15000

Lag Distance, meters Lag Distance, meters

Kent County, Lithology Based K Kent County Specific Capacity Based K

782
ft/day

 

1 .5 55
ft/day ft/day
:—
5.2 miles 5.2 miles

68

CHAPTER 5: NUMERICAL MODELING OF NET
DRAWDOWN

Recognizing that the SWL does not reflect the effect caused by pumping wells, we have
developed an approximate approach to model the net drawdown around the wells and

superimpose it onto the water level distributions.

1. Mapping Groundwater Velocity Patterns

The Wellogic database provides us with water level data and hydraulic conductivity data
(which includes historical pumping effects). Once the SWL and hydraulic conductivity

are mapped, we can compute the groundwater seepage velocity based on the Darcy’s law:

17(X,y,t)= K(x9y)j(x:yat) (9)

e

 

where K(x, y) is hydraulic Conductivity, J(x, y, t) is the hydraulic gradient, and ne is the
effective porosity. Since there are little site-specific data on effective porosity, we use

typical values ranging from 0.15 to 0.3.

Figure 28 presents results for the predicted velocity flow pattern for a selected site in

Michigan.

69

 

—=
0 1 mile

Figure 28: Predicted velocity flow pattern for the City of Hudson, Michigan.

2. A Numerical Superposition Approach for Modeling of Net Well
Drawdown

Since the SWL does not represent the effect of pumping, we present in this subsection a
numerical superposition approach to model the net drawdown around wells. The
fundamental basis of a superposition approach is that, for a linear or approximately linear

system, a complex problem can be decomposed into more simple sub-problems. The sum

70

of the solutions of the sub-problems will be the same as the solution to the whole, more

complex problem.

Since the net drawdown is defined by

s(x,y,t) = h0 - h(x, y, t) (10)
where ho is the static water level before including the effects of pumping well(s) being
evaluated, h is the post-pumping head after the wells are included. We can derive an

equation describing directly the net drawdown by subtracting the equation describing the

post-pumping condition from that describing the pre-pumping condition.

The governing equation describing the 2D depth averaged flow under pre-pumping

condition is:

 

S%=2[T%]+—6_ Taho +8 +Csw(hsw—h0)+Cdrain(hdrn-h0) (11)
5t 8x 6x 6y 6y

where T is transmissivity, S is elastic storage coefficient for a confined aquifer or specific
yield for an unconfined aquifer, a is head independent flux (e.g., recharge, background

pumping), and the last two terms represent head dependent surface water and drain

fluxes, C C

5w ’ drain

being surface water and drain leakance characterizing surface water-

groundwater connection (they are zero outside the surface water body or drain areas).

Equation 11 is subject to the following boundary and initial conditions:

ho = gl (x,y,t), (x, y) e an, (123)

ah 12b
—T-a"r;0'=g2(xayat)9(xay)ea§22 ( )

71

ho =f(x,y),t=0 (120)
where g1 is the head on boundary 6!), , g2 is the flux on boundary 6!), , and f (x, y) is

the initial head.

When the effect by the wells being delineated is taken into account, the governing

equation becomes:

6h 6 6h 6 ah
S—=-—T—+—T—++Ch—h+C.h—h+ .5_,;__
at ax( ax) @[ 6y] 8 sw( SW ) dram( drn ) ZQ, (x xw: y ym)
(13)

where Q,- is the pumping rate from the iLh well. The initial and boundary conditions
remain the same as those for pre-pumping condition, assuming that model boundaries are
selected such that the impact from the new pumping stresses do not propagate beyond the

boundaries.

Subtracting Equation 13 from Equation 11, we obtain the following equation describing

the net drawdown dynamics:

as_6

as 6 as
E—B;(T ]+—[TE]+CSW(O—S )+Cdrain(0_s )+;Qi6(x-xwl;y—ym)

E 6y
(14)

with the following “homogenized” initial and boundary conditions:

50 = 0, (x, y) e 60, (15a)

72

 

—T%=0,(x,y)e§flz (15b)
n

so =0,t=0 (15C)

Note Equation 14 is of the same form as the original flow equation (Equation 13) under
post-pumping condition, with zero initial and boundary conditions and surface water head
and all recharge and discharge terms removed, except for the pumping stresses being
evaluated. This creates an initial condition where there is no flux between the aquifer and
surface water features. Simulation of a new aquifer stress will induce flux from
represented surface water features in an amount that is equal to the depletion of rivers and
springs for the same stress in the full model. The results from this simulation represent
the impacts from the particular aquifer stress being evaluated in isolation of all other

recharge and discharge.

A simple example of an evaluation using numerical superposition would be an evaluation
of the impacts to river reaches due to pumping at a single well. Pumping at the well does
not affect any of the sources of recharge or discharge, which are not hydraulically
connected. The cone of depression from the pumping well will propagate radially from
the well until the resultant drawdown affects water levels near a river reach. At that time,
the pumping will result in a reduction of the river gain or increase in river loss. By
analyzing this stress using the numerical superposition model, all exchanges between the
river and aquifer will be due to the groundwater pumping being evaluated. Evaluation of

the same pumping well using the full model would be significantly more cumbersome.

73

It should be pointed out that the numerical superposition approach is, strictly speaking,
only valid for a linear, or confined aquifer model; unconfined aquifer model
representations are non-linear due to the fact that aquifer transmissivity changes as
aquifer water levels change. In most cases, the changes in water levels are very small

relative to the total saturated thickness, so these nonlinearities are considered negligible.

Figure 29 illustrates the basic concept of the numerical superposition approach.

 

 

SWL surface w/o pumping Calculated change in head due to SWL surface w_it[1 effects of
(Data) pumping (Model) pumping

Figure 29: A simple example illustrating the numerical superposition approach.

Figure 30 also illustrates the effect of the numerical superposition approach. In this
example, we have refined the water level surface to account for the effects of pumping in
which we have superimposed the net drawdown of a pumping well to the regional
groundwater flow map. This illustration represents the effects of pumping at 0, 200 and
2000 GPM. Note that there is very little effect shown between no pumping and pumping
at 200 GPM. The effects of pumping can be noticed by looking at the water level surface
in the immediate vicinity of the well or by looking at the distribution of the particles that
were released from the pumping well. Also note that pumping does not have a large

effect in areas of strong velocity flow fields, which are dominated by the regional system.

74

 

.No pumping\

\\

  

“

 

 

 

_, X“ ‘
\

      

 

 

 

 

 

 

A
Wells pumping at 200 GPM \

 

 

 

 

 

 

 

r .

  

 

 

 

 

 

 

' Wells pumping at 2000 GPM
" 110W \\\\\\‘

 

 

 

 

 

 

 

 

 

 

 

   

 

. x I
\ . Y-
- . '1 ~ -.
.\ ‘ - :37. V
, I :: . .
I ( Ix‘
' «"4“:
\a/
. h."
\‘l
\
pr. 11‘.
‘ I

 

0 125 mile

 

Figure 30: Predicting flow systems with response to pumping at 0, 200, and 2000 GPM. The dark gray
elongated polygons represent the WHPA that is developed using backwards particle tracking from the

pumping well.

75

CHAPTER 6: MAPPING WELLHEAD PROTECTION
AREAS

The long-term viability of Michigan’s groundwater resources can be threatened by
contaminated recharge water. This can give rise to increased treatment costs, and may
lead to the eventual abandonment of the resource altogether. Once an aquifer has become
contaminated, it is often extremely difficult, if not impossible, to clean up, owing to the
complex and tortuous flow pathways through the porous medium. It has been
increasingly recognized, therefore, that prevention is better than cure. This has resulted in
the development of various groundwater protection measures, as a means of regulating
activities potentially damaging to sub-surface water resources. One such approach is
well capture zones for groundwater source protection. This involves delineating areas
around a supply well, which can either define the region from which the well collects its

recharge water or the travel times for groundwater to reach the well.

The process for delineating well capture zones, however, can be expensive, especially
when the aquifer in question is strongly heterogeneous. At most sites in Michigan,
aquifer conditions and properties vary significantly, sometimes over surprisingly small
distance. For the most part, these conditions can only be observed at a few scattered wells
or boreholes. Even for the limited amount of data available they are often collected at a
particular time and do not necessarily reflect the long term mean condition that dictates

contaminant transport and pollution control measures.

76

All of these problems are familiar to the planners, managers and consultants responsible
for site assessments and resource management. Uncertainties about the way water and
contaminants move in the subsurface environment translates into problems for both risk
assessment and pollution control. If we do not properly characterize the groundwater
flow systems, especially the flow direction, we may greatly underestimate the risk to
human health from potential hazards. If we are unable to accurately map well capture

zones, it is less likely that wellhead protection will be successful.

1. Review of Traditional Delineation Methods

In this research, we make use of the Wellogic-based flow system to map Michigan’s
Wellhead Protection Areas (WHPAs). The United States Environmental Protection
Agency (USEPA, 1987) defines a WHPA as “the area surrounding a pumping well that
encompasses all areas or features that supply ground—water recharge to the well.” A
number of methods can be used to delineate WHPAs “footprints” for a variety of
hydrogeologic environments. These methods vary considerably in input data
requirements, difficulty of application and cost. In order of increasing technical

complexity, the general classes are as follows and are illustrated in Figure 31.

1.1. F ixed-radius method

The fixed-radius method simply consists of identifying the WHPA radius and drawing a
circle with that radius around the wellhead. The area inside this circle is the wellhead

protection area. Figure 31 includes an example of a circular wellhead protection area. The

77

method can be applied with no local knowledge of conditions that control flow and

transport in the vicinity of the wellhead.

Although this method is very easy to implement, it ignores, or vastly simplifies, the
hydraulic behavior of wells in their specific hydrogeologic setting. WHPAs identified by
this method may be either too large or too small, resulting in wellhead over-protection or
under-protection. If the radii are too large, more land use activities may require regulation

than necessary to effectively protect the selected wells.

Approximately 15% of the existing Michigan delineations utilized this method.

1.2. Analytical method

The analytical method involves a site-specific application of analytical equations

describing flow to a well.

The advantage of this analytical method is that it better incorporates local hydrogeologic
data and hydraulic behavior than the fixed radii method (see Figure 31). One
disadvantage of this method is that constant hydraulic parameters (such as hydraulic
conductivity and gradient) must be applied uniformly to the vicinity of the well when
these parameters may, in fact, vary over the region. Another disadvantage is the inability
to incorporate effects of nearby pumping wells. Compared to the fixed-radius method,

application of the analytical method is technically more involved.

78

Approximately 80% of the existing Michigan delineations utilized this method.

1.3. Numerical method

The numerical method is the most sophisticated among the different methods. This
method can incorporate the effects of complex flow fields, complex hydraulic property
zonation, boundary conditions and well interference through the use of numerical
computer codes which model groundwater flow. The benefit of this method is its ability
to better incorporate spatial variability in hydrogeologic parameters, multiple well
locations and irregular boundary conditions than analytical methods. Because of the
complexity of codes used to simulate flow and transport, WHPA delineation with this
method requires a high degree of technical expertise. This approach will tend to be the
most expensive but, if properly implemented, will yield the most realistic and precise

representation of a WHPA, especially in the more complex hydrogeologic settings.

Approximately 5 % of the existing Michigan delineations utilized this method.

79

 

 

      
  
      

 

l 9 ~‘. i

 

 

 

. Aer—f ~~~~~ i WHPA by n merical method:
i i l l Allows accounting for aquifer variability,
i E changing flow direction, but requires
ID ‘i ! significant site-specific data, very expensive.
'\ I
_. ‘ - a l l V l
WHPA by analytical method _ ' l l . 3 .. l-
\ Assumes uniform regional flow, ~ ‘l‘ a \9 . A-fi—VL l ‘ .3. i
“ rcqurres less data than numerical ‘2‘ I“ ‘i Q i l' 1
methods. less accurate away from the ° l l ' i {3" E ' fl” _ _
. . l . . e i" Viv ~1-
*— well. still expensrve. A. l l i l i ‘6‘ /
. y I”
l i x 1v
..‘4 ——-"—;‘/‘-

 

 
 
    

 

‘1‘
‘M

“ .

_ i

 

 

 

 

 

     
  

 

. L‘”--~T‘T"_—‘"’M’r—"~
i i i s i
l
l—___H : __ WHPA by fixed radius method: -'
i —__T\,A A lgnores site-specific regional flow
. ' and geology, very easy to ~—l
fl 1 implement.

1 Q1 I
l, x, l
i l
i‘ l
l

\\

 

 

 

Figure 31 : An illustration of the conceptual difference in the WHPAs developed by different methods.
The big arrows represent regional flow direction.

80

 

2. Wellogic-Based WHPA Delineation

2.1. Delineation based on mean advection

In this research the water level and hydraulic conductivity data obtained from the
Wellogic database were used to develop wellhead delineations. We advocate using
Wellogic data for WHPA delineation. The detailed water well data enables modeling the
effects of spatial variability, complex sources and sinks, and geological boundaries
without collecting new data. The approach consists of the following steps:

1. Processing and mapping static water level

2. Computing, processing, and mapping hydraulic conductivity using specific

capacity data
3. Modeling the effect of net drawdown using the numerical superposition approach
4. Computing and mapping groundwater seepage velocity using the Darcy’s Law

5. Performing reverse particle tracking

Coupled with our statistical approach, the Wellogic-based delineation can be readily

implemented.

2.2. Delineation based on mean advection and dispersion

Delineations based on mean advection tend to significantly underestimate the size of the
WHPAs. In this research, we have compared two methods to account for the effect of
unresolved small-scale variability to provide ensemble mean WHPAs which include; (1)
Monte Carlo simulations, and (2) the addition of dispersion. Monte Carlo simulations are

used to estimate the probability of an occurrence given the statistical parameters of an

81

aquifer. All of the possible realizations are simulated one after the other and the
probability of occurrence of any particular event is calculated sequentially for each
realization. This process requires a large amount of computational time, which makes it

expensive.

We also propose modeling the aggregated effect of dispersion caused by small scale
heterogeneity that cannot be resolved with the Wellogic data. This can be accomplished
by including a component of Brownian motion in the final particle-tracking step. The
resulting procedure will produce a WHPA that is wider due to transverse dispersion and
longer due to longitudinal dispersion. The longitudinal dispersivity can be estimated from

the following equation (Gelhar et al., 1992):

2

AL :me ’1 (16)

where xi is the unresolved scale of variability and can be assumed to be one tenth of the
correlation scale obtained from the conductivity variogram, 0,1,), K represents the
variance of subscale or local heterogeneity in log conductivity and can be conservatively
assumed to be equal to the nugget in the log conductivity variogram. Because Gelhar
original equation for transverse dispersivity is inaccurate, we have chosen to be more
conservative with respect to WHPA delineations. In this case, the transverse dispersivity
is assumed to be one third of the longitudinal dispersivity (EPA 1986). This method
provides equivalent WHPA envelops to account for small-scale variability (see Figure

32) and is much more cost effective because it can be performed within one realization.

82

Figure 33 illustrates the final results using our approach at 8 selected sites and compares
them with the corresponding existing delineations based on measured data. The
delineations are performed for lO-yr WHPAs both with and without including the effect
of dispersion. The areas in gray are WHPAs delineated using the well database and the

black outlined polygon is the WHPA delineated using traditional measured data.

The results show that our statistically-processed WHPA based on Wellogic data are
consistent with existing delineations based on measured data. When the effect of
dispersion is taken into account, our Wellogic based delineations are often fatter, longer,

or more conservative than existing delineations.

83

 

 

 

 

 

 

  

Monte Carlo:
300 Realizations
100 days

2 _
0an - 1
A = 20 m
At = 0.25 days

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

100 days
At = 0.25 days

Figure 32: WHPA comparison using (1) the Monte Carlo approach (top) and (2) with the addition of
dispersion (bottom). The colored zones in the WHPA developed with the Monte Carlo approach indicates
the probability of capturing contamination located anywhere between the outer bound of a zone and the
capture well within 10- years, and the grey area represents the lO-year WHPA developed with the effect of
dispersion. The heterogeneity (in hydraulic conductivity) in the Monte Carlo example is replaced by the
addition of dispersion in the bottom WHPA. The random field statistics of the Monte Carlo example (02 =
l, and k = 20) are replaced by the addition of dispersion (AL) computed using Gelhar’s equation (see
equation 12). Note: model dimensions are X-Length is 1000m and Y-Length is 750 m.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

84

Augusta Creek WHPA comparison Augusta Creek WHPA comparison
(no dispersion) (dispersion, AL=118m, AT=40m, ne=0.25)
' U a .
t ‘.

 

   

n“ V“ v m ,.22.2;
.; a7“ ‘

- ..

r: = "ma/j,”
/ j
4 .

7L...

\

t.

 

 

 

 

 

 

Durand WHPA comparison Durand WHPA comparison
(no dispersion) (dispersion, AL=70m, AT=23m, ne=0.25)

 

 

-= _=
0 2 Mile 0 2 Mile

Figure 33: Comparison of Wellogic-based WHPA delineations (gray area) with existing delineations based
on measured data (black outlined polygon) at selected sites in Michigan. The WHPA delineations on the
left represent those with no dispersion and those on the right represent WHPA delineations with the
addition of dispersion for each respective location

85

Lenawee County WHPA comparison
(no dispersion)

Lenawee County WHPA comparison
(dispersion AL=245m, AT=82m, ne=0. 25)

 

 
   

Avg-4M ~~v¢ /+i',""f f J“
V l V r 720
/

l
i 7’ A .

'/ //i." I 1740 xvf)‘ [if—i. 1.2.7 l

E/ ' .

      

 

, \t i l )
,---i:l:l;.£;”

 

 

 

 

 

 

Hudson WHPA comparison

1.5 Mile

Hudson WHPA comparison

 

no dis ersionL

 

 

 

(dispersion, AL=73m, AT=24m, ne=0.25)

 

 

 

0

86

 

n

M

2
law
m2
3%

am 0
m .
mm.
6
ma
HA
Wm
.mm
W:
cm.
am
mm
9%
lb
m

.k

M

2

0

 

 

 

Isabella County WHPA comparison
(no dispersion)

 

 

 

 

87

CHAPTER 7: SUMMARY

Most groundwater investigations are hampered by a lack of data with which to
characterize spatial variations in hydraulic conductivity and water level. However, data
from public water well records represent a vast yet under-utilized source of information.
One reason these data have been neglected is the general perception that they are of poor
quality. Indeed, short-duration pump testing, inaccurate recording of results in the field
and transcription errors in the public water well records may introduce noise in SWL and
specific capacity data sets. Nonetheless, we have demonstrated, using an integrated
statistical and geostatistical approach, that water well data are not overwhelmed by
random errors and can be used to produce meaningful maps of groundwater flow systems
with an accuracy sufficient for most groundwater investigations. In fact, statistical
Wellogic-based water level and hydraulic conductivity distributions have a number of

important advantages over traditional distributions that rely on limited measured data.

These advantages are:

0 The spatially-detailed Wellogic data, when statistically processed, allow
representing spatial variability in aquifer properties.

0 The spatially-detailed Wellogic data, when statistically processed, allows
modeling nonuniform regional flow reflecting the influence of nearby streams,
lakes, wetlands, pumping stresses, variable recharge, and geological variability
and boundaries. (e. g., in most existing delineations in Michigan, due to data

limitation, the regional mean flow is represented as being uniform and

88

unidirectional and the resulting WHPAs are usually less accurate away from the
well).

0 The 40-years of Wellogic data enable modeling long term mean regional flow
conditions and this long-term condition is what dictates most groundwater

protection measures.

Since the Wellogic data are “free”, cumulative, and available virtually “everywhere”
from historical to new water well logs, this finding has major implications for
groundwater research and professional investigations. In particular, the finding implies
that if the water level and hydraulic conductivity data can be systematically integrated we
will be able to accurately map the groundwater water level flow and velocity patterns in a
way that is spatially-detailed without collecting new measurements, drastically reducing
the cost of groundwater management investigations that require knowledge of the long
term mean groundwater flow system. Most importantly, our ability to characterize
groundwater systems can “steadily” improve with time as the Wellogic database expands

and the data quality improves.

Our systematic evaluations also reveal a number of limitations in using the Wellogic
database. These limitations are:

0 There are still sites where the Wellogic data may be too scattered or clustered to
provide a sufficiently accurate flow description. Although we may improve the
quality of Wellogic mapping by expanding the analysis area, “borrowing”
regional information to help local interpolation, this may not always be effective,

especially when local dynamics deviate significantly from the regional trend.

89

 

 

Errors caused by insufficient data should decrease with time as the Wellogic
database continues to expand and data quality improves.

Wellogic-based head distributions may become less accurate in areas
characterized by a very flat hydraulic gradient. Accuracy of Wellogic-based
predictions in a flat-gradient area should increase with time as the Wellogic
database expands and, especially, as the data quality improves.

Although our specific capacity-based hydraulic conductivity seems to be
generally consistent with the available lithologic descriptions, it is difficult to
verify how well the detailed spatial patterns reflect the true spatial distribution.
There are simply no measured data detailed enough to validate the hydraulic
conductivity distributions.

Since water wells are mostly screened at elevations where the aquifer material is
more permeable, the specific capacity-based hydraulic conductivity may be
systematically larger than the actual average conductivity for the entire drift
package. For the purpose of wellhead delineation, however, this is conservative.
The sufficiency of this estimate for other applications needs to be reassessed.
Although Wellogic data can be used to map transient flow, the accuracy of the
predicted transient pattern decreases as the scale of temporal variability decreases.
With the current data density, Wellogic-data can only be used to model large-
scale temporal trends (decadal trends, multi-year trend trends) in areas where the
Wellogic data is very detailed (e. g., in Livingston and Kalamazoo county). The
level of temporal details that can be modeled should increase with time as the

Wellogic database expands and the data quality improves.

90

 

0 Although our geostatistical approach enables smoothing out data errors and
noises, we felt sometimes it is too “aggressive”, killing some desired spatial
details. Future research should focus on how we can maximize the extraction of
information, including small-scale dynamics using, for example, a multiscale data
analysis technique. Our ability to map more detailed local systems should also

improve with time as the Wellogic expands and the data quality improves.

The Wellogic data and our statistical modeling approach can be potentially used in many
applications requiring knowledge of groundwater flow system. The following are a few
examples:

Contaminant fate and transport analysis. The Wellogic modeling system can be
potentially used to predict contaminant transport in Michigan’s aquifer systems. In
particular, the modeling system can be used to simulate the migration, diffusion,
dispersion, and attenuation of contaminants leaked from industrial sites, underground

storage tanks, landfills, agricultural fields, and other pollution sources.

Designing remediation capture systems. The Wellogic modeling system can be
potentially used to design contamination capture systems. In particular, the modeling
system can be used to evaluate cleanup options and optimize extraction rates, patterns,

and scheduling.

Evaluating the impact of aquifer use and stream depletion. The Wellogic modeling

system, coupled with the statewide GWIM database, can be potentially used to predict

91

 

how changes in land use and pumping stresses alter the base flow to Michigan’s streams,

lakes, and wetlands.

Evaluating well permitting and water conflicts. The Wellogic modeling system, coupled
with the statewide GWIM database, can be potentially used to evaluate well interference
(e. g., potential conflict between a high capacity well and nearby domestic wells) and

address water use conflicts.

Characterizing aquifer yield. The Wellogic modeling system, coupled with the statewide
GWIM database, can be potentially used to evaluate Michigan’s aquifer yield, accounting
for the influence of site-specific spatial variability, geological boundaries, and the

presence of streams, lakes, and other sources and sinks.

Mapping aquifer vulnerability/sensitivity. The Wellogic modeling system can be
potentially used to map aquifer “sensitivity”. Sensitivity is a term used to describe how
well an aquifer is protected from infiltrating contamination. A highly sensitive aquifer
would have little or no defense, whereas an aquifer with low sensitivity would be very
well protected. Wellogic can be used to quantify aquifer sensitivity as a function depth to
water, vertical permeability SWL, and recharge. This sensitivity map, coupled with a
WHPA map and a map of known and potential hazards, will provide valuable

information for wellhead protection.

92

 

Mapping three dimensional aquifer heterogeneity. Our analysis has so far focused on
two-dimensional depth averaged descriptions. In most cases, a contaminant plume is
localized and does not sample the entire aquifer thickness and a detailed 3-D
geostatistical description of aquifer heterogeneity is needed. Although such a 3-D
description is traditionally deemed impractical, the Wellogic database containing detailed

3-D lithologic information makes this possible.

The importance of water well reports in hydrogeological investigations can simply not be
overstated. Filling in those blanks doesn’t just satisfy some agency’s requirements; it also
provides hydrogeologists from the public and private sector valuable information
regarding Michigan’s ground water systems. Well drillers can provide important
contributions to the groundwater science by making careful observations and
measurements when recording that data on the well report. We strongly recommend that
the Michigan Wellogic database and related analysis tools be continuously maintained

and refined and the systematic Michigan data integration effort be extended nationwide.

93

APPENDICES

94

 

 

Appendix A
Specific capacity-based Hydraulic Conductivity (K) published Equations

95

 

Modified T heis Equation: (Driscoll)

The Modified Theis method is a commonly-used simplified variation of the Theis
method. Water well records also may report the results of test pumping the well to
demonstrate it has sufficient yield for the purpose drilled, and the information can be
used to estimate the specific capacity of the well. The use of specific capacity to estimate
aquifer yield was investigated by relating specific capacity to transmissivity. The
specific capacity of a well is the yield of the well in volume per time per unit of
drawdown (Driscoll, 198 6). In addition to characterizing the performance of a well,
specific capacity values reported on water well records may be used to estimate the
transmissivity of the aquifer near the well. The simplest relation between specific
capacity and transmissivity is obtained by considering the Theis solution describing

drawdown in an infinite confined aquifer (Freeze and Cherry, 1979).

S(r) = 2%— W(u) (1A)

where s(r) is the drawdown at a radial distance r from the pumping well, Q is the
pumping rate, T is aquifer transmissivity, W(u) is the well function in which 11 is a

dimensionless term that depends on aquifer characteristics, the radial distance, and time.

Rearranging equation 1A to arrive at an expression for transmissivity and to isolate
specific capacity (Q/sw), explicitly writing the equation for drawdown at the pumping

well, and including the definition of u yields,

96

 

 

 

 

“Hi-wt“)
sw 47: 4T t

where sW is the drawdown at the pumping well, rw is the radius of the pumping well, t is

time since pumping started, and S is the aquifer storativity.

Numerical testing has revealed that if the specific capacity is provided, then a simple
iterative approach can be used to determine transmissivity. To use equation 2A to
calculate aquifer transmissivity, aquifer storativity, radius of the pumping well,
drawdown at the pumping well at a given time, and discharge from the well must be
provided. Aquifer storativity is generally not known. Fortunately, storativity does not
vary over a large range once the aquifer is determined to be either confined or
unconfined, and, because this term appears in the well function, transmissivity is not a
strong function of storativity. If the well is known to be confined or unconfined, a
reasonable guess for storativity can be made to yield a reasonable approximation of

transmissivity.

Water-well records typically report the results of test pumping the well, and the
information requested on the Wellogic form includes test duration, test pumping rate,
drawdown at the well, well diameter, and test type. The test types include bailer, plunger,

airlifi or test pump tests.

A modified equation has been developed from the above equations to estimate the
potential transmissivity of a well (Driscoll, 1986). By re-arranging the terms of equation

2A, the specific capacity becomes the following:

97

 

 

 

Q T (3A)

 

 

If typical values are assumed for the variables in the log function of the equation such as t
=lday, r = 0.5ft, T = 30,000gpd/ft, and S = 1x10'3 for a confined aquifer and S = 7.5x10‘2
for an unconfined aquifer, the specific capacity of the aquifer is given by (Driscoll,

1986):

Q T (4A)

By assuming default values for all of the parameters in the Theis equation, it can be
shown that the value for the specific capacity multiplier for a confined aquifer is

approximately equal to 2000 and for an unconfined aquifer is equal to 1500.

Bradbury and Rothschild (1985)

Bradbury and Rothschild (1985) developed a computer program to estimate
transmissivity from specific capacity data. This method involves utilizing a method of
correcting the observed specific capacity for well losses and introduces a partial
penetration correction factor (Fetter, 1993). The following equation yields an estimate of
transmissivity, which is corrected for well loss and partial penetration and also

incorporates specific capacity data.

98

 

 

 

T = Q [rn(2'25T’]+ 2s ] (5A)
47t(s - SW) r: S p

where T is transmissivity (L2/'T), Q is the pumping rate (L3/T), sw is well loss (L), s is the
drawdown in the well (L), t is pumping time (T), S is storage coefficient (dimensionless),

rW is the well radius (L), and sp is the partial penetration factor.

Razack and Huntley (I991)

Razack and Huntley (1991) studied the relationship between specific capacity and
transmissivity in an alluvial groundwater basin in Morocco. They analyzed 215 data
pairs where transmissivity and specific capacity were known for a well. They were able

to find the following relationship:

 

T=15.3[ Q J (6A)

where T is transmissivity (mz/day), Q is the pumping rate (m3/day), and ho-h is the

drawdown (m).

99

Mace (I99 7)

Mace (1997) studied the relationship between specific capacity and transmissivity data
from 71 wells in the Karstic Edwards Aquifer of Texas (Fetter, 1993). He found the

following relationship:

 

1.03 7A
T=O.76[ Q J ( )

where T is transmissivity (mZ/day), Q is the pumping rate (m3/day), and ho-h is the

drawdown (m).

Well E fficien cy

This well efficiency factor only applies to the specific capacity data and only for the
methods based on the Theis nonequilibrium model such as the Theis, Bradbury and
Rothschild, and Modified Theis methods. The reported drawdown is multiplied by the
well efficiency factor. A Well efficiency factor of 1 represents a 100 percent efficient
well. Likewise a value of 0.6 represents a well that has 60 percent efficiency.
Multiplying the drawdown by the well efficiency factor results in a proportionally higher

specific capacity and estimated transmissivity.

100

Appendix B

101

«23:22

SE. mcuctml wtonoml £a_o£§u\ma.cE. 25m. 5.3:. gumbo—2:. a: Somme—:5—
Boo. 3:... fine. amass? at $6525.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

$03.23*.
oz oz 02 mo> oz mo> mm; .. . .5 v5. : mmwcmx
oz mo> mo> mo> mo> mo> wo> =55. xoucfiwmgflmol =2»:ng u:=2mto§<..c2>om.c.3§$as: manic.
mm> a; 2; m; 8> 8> «.82 E... . . .. . .. i: 205...
oz wo> mo> mo> oz mo> mo> EE 4 . . : .r: onmu.
<2 <2 <2 <2 <2 <2 02 E; . 1 : . ,. . .. i :95:

 

oz mo> mo> mm; 02 mo> mo> c_mo._moansmoncoao<Q>omngo~m>>dofiommEEd: mc0Nt<

         

E_o3\m=.,xm.o§wcco.23min: mxwm_<

mmocuom oo>>

4m: 2: 5 3% :03 museum“ Hm: 03$

102

. oEa w...“ .m ..ooccoo 50: 2300
SE .m__.o3<m32.2§>\9o\m3.. .3. 0.6.3.. E13333?

.m... m... :25... > ...m 3.95.... .33.:
own. 523on xouéw: 5. 23M. .23.. ESE. a:

\m3...o.m.m.m. 23.333395:

«59.2.0

 

xc .257 EC 2555.th

L

 

O.

 

 

 

 

 

 

 

 

LC Ct“? ..cc

 

xamm.>.o3OI=o>>ao.Mo:O\m3.cc.o.m5m.mou.32<s\\dz:

2.2m. EmL

 

 

 

 

 

 

 

oz mo> mo> mo> wo> mo> mo> 302
oz wo> mm> mo> mm> wo> mo> Eta . .o . ... . ... . ...E : t: mvm>mZ
mo> mm> wo> wo> wo> mo> wo> con. . : .. . a: «EVE... 5. 20223.95; mxwmbwz
98 93 .wam c3. 0.32.
.950“. 90 329.5 3.088 .96.. . $0.. :03 o... 0. $925 83 35m
050on 56>? Echo“. .Qm>> o_nm__m><

 

 

 

owmomfio

 

 

 

103

10.

11.

12.

13.

14.

REFERENCES CITED

. Aichele, SS. 2005. Water Resources in a Rapidly Growing Region-Oakland

County, Michigan. US. Geological Survey Open File Report 2005-1269.

Anderson, MP. and W.W. Woessner. 1992. Applied Groundwater Modeling:
Simulation of Flow and Advective Transport. Academic Press.

. Bradbury, K. R. and E. R. Rothschild. 1985. A computerized technique for

estimating the hydraulic conductivity of aquifers from specific capacity data.
Ground Water, vol. 23, no. 2, p. 240-246.

Batu, V. 2005. Applied Flow and Solute Transport Modeling in Aquifers. CRC.
Boutt, D. F., D.W. Hyndman, B.C. Pijanowski, and D.T. Long, 2001, Modeling
Impacts of Land Use on Groundwater and Surface Water Quality. Ground Water,
39 (1), 24-34.

Brown, R. H., 1963. Estimating the transmissivity of an artesian aquifer from the
specific capacity of a well. US. Geological Survey Water Supply Paper 1536-1, p.
336-338.

Charbeneau, R. J. 2000. Ground Water Hydraulics and Pollutant Transport.
Prentice-Hall, Old Tappan, N. J.

Chilés, J .P. and P. Delfiner. 1999. Geostatistics: Modeling Spatial Uncertainty. New
York: Wiley-Interscience, John Wiley & Sons, Inc., 695 p.

Curran, D., L.D. Leske, and L. Miller. 1981. Hydrogeologic Atlas of Michigan.
Western Michigan University Department of Geology, Kalamazoo, Michigan.

Davis, SN. and J. M. Roger De Wiest. 1966. Hydrogeology. Wiley: New York.
Davis, J .C. 1938. Statistics and Data Analysis in Geology. Wiley: New York.
Driscoll, F. 1986. Ground Water and Wells. Johnson Division: St. Paul, MN.
El-Naqa, A., 1994. Estimation of transmissivity from specific capacity data in
fractured carbonate rock aquifer, central Jordan. Environmental Geology, v. 23,

no. 1, p. 73-80.

Fetter, CW. 1993. Contaminant Hydrogeology. Macmillan, New York.

104

15.

16.

17.

18.

19.

20.

21.

22.

23.

24.

25.

26.

Firouzian, A. 1963. Hydrological studies of the Saginaw Formation in the Lansing,
Michigan area-1962. M.S. thesis: East Lansing, Michigan State University, 45 p.

F1eis& Vanderbrink Engineering, Inc. 2001. Village of Augusta Wellhead
Protection Plan.

Freeze, R. A. and J. A. Cherry. 1979. Groundwater. Englewood Cliffs, New Jersey,
Prentice-Hall, Inc, 604 p.

Gelhar, L.W. 1992. A critical Review of Data on Field-Scale Dispersion in
Aquifers. Water Resources Research, vol. 28, no. 7, p. 1955-1974.

Hill-Rowley, R., T. McClain, and M. Malone. 2003. Expanding a static water level
mapping methodology from a township to the watershed scale in Michigan. The
Great Lakes Geographer, Vol. 10, No.1.

Hill-Rowley, R., T. McClain, and M. Malone. 2003. Static Water Level Mapping in
East Central Michigan. Journal of the American Water Resources Association
(JAWRA) 39(1): 99-111.

Hoaglund 111, J .R., Huffman, G.C. Grannemann, N.G., 2002, Michigan basin
regional ground water flow discharge to three Great Lakes, Ground Water, 40(4), p.
390-406.

Holtschlag, D.J, C.L. Luukkonen, and J .R. Nicholas. 1996. Simulation of Ground-
Water Flow in the Saginaw Aquifer within Clinton, Eaton, and Ingharn Counties.
Michigan: US. Geological water-Supply Paper 2480, 49 p.

Hovorka, S. D., R. E. Mace, and E. W. Collins. 1995. Regional distribution of
permeability in the Edwards Aquifer. University of Texas at Austin, Bureau of
Economic Geology, Final contract report to the Edwards Underground Water
District under contract no. 93-17-FO, 126 p

Hovorka, S. D., R. E. Mace, and E. W. Collins. 1998. Permeability structure of the
Edwards aquifer, south Texas-Implications for aquifer management. University of

Texas at Austin, Bureau of Economic Geology, Report of Investigations No. 250,
55 p.

Huntley, D., R. Nommensen, and D. Steffey. 1992. The use of specific capacity to
assess transmissivity in fractured-rock aquifers. Ground Water, v. 30, no. 3, p. 396-
402.

Hyndman, D.W., M. J. Dybas, L. Fomey, R. Heine, T. Mayotte, M.S. Phanikumar,
G. Tatara, J. Tiedje, T. Voice, R. Wallace, D. Wiggert, X. Zhao and C. S. Criddle,

2000. Hydraulic Characterization and Design of a Full-Scale Biocurtain. Ground
Water, Vol. 38, No. 3, p. 462-474.

105

27.
28.

29.

30.

31.

32.

33.

34.

35.

36.

37.

38.

Isaaks, E.H. and R. M., Srivastava. 1990. An Introduction to Applied
Geostatistics. Oxford University Press. .

Kitanidis, P. K. 1997. Introduction to Geostatistics: Applications in Hydrogeology.
Cambridge University Press.

Kruseman, G. P., and N. A.,de Ridder. 1994. Analysis and evaluation of pumping
test data. International Institute for Land Reclamation and Improvement, The
Netherlands, 377 p.

Li, S.G. and Q. Liu. 2006. A Real-time, Computational Steering Environment for
Integrated Groundwater Modeling. Journal of Groundwater. 44(5), p. 75 8-63.

Li, S.G., Q. Liu, and S. Afshari. 2006. A Hierarchical Patch Dynamics Paradigm
(HPDP) for Modeling Complex Groundwater Systems Across Multiple Scales.
Environmental Modeling and Software. Vol. 21, Issue 5.

Li, S.G. and Q. Liu. 2006. Interactive Ground Water (IGW), Enviromnental
Modeling and Software. Vol. 21 , No. 3.

Li, S.G., and Q. Liu. 2003. Interactive Ground Water (IGW): An Innovative Digital
Laboratory For Groundwater Education and Research, Computer Applications in
Engineering Education 10(6), 31 p.

Liu, H, S. Shah, and W. J iang. 2004. On-line outlier detection and data cleaning.
Computers and Chemical Engineering, 28, p. 1635—1647.

Luukkonen, C.L., N.G. Grannemann, and DJ. Holtschlag. 1997. Ground-Water
Flow in the Saginaw Aquifer in the Vicinity of the North Lansing Well Field,

Lansing, Michigan—Part1, Simulations with a Regional Model. US. Geological
Survey Open-File Report 97-569, 13 p.

Luukkonen, C.L., N.G. Grannemann, and DJ. Holtschlag. 1997. Ground-Water
F low in the Saginaw Aquifer in the Vicinity of the North Lansing Well Field,
Lansing, Michigan--Part2, Simulations with a Regional Model Using a Reduced
Cell Size. US. Geological Survey Open-File Report 97-570, 25 p.

Luukkonen, CL. 1995. Ground-Water Withdrawals in Clinton, Eaton, and Ingham
Counties, Michigan. US. Geological Survey Fact Sheet 226-95, 2 p.

Mace, R. E., Nance, H. S., and Dutton, A. R., 1994. Geologic and hydrogeologic
framework of regional aquifers in the Twin Mountains, Paluxy, and Woodbine

Formations near the SSC site, North-Central Texas. Draft topical report submitted
to the TNRLC under contract no. IAC(92-93)-0301 and no. IAC 94-0108, 48 p.

106

39.

40.

41.

42.

43.

45.

46.

47.

48.

49.

50.

Mace, R. E., 1995. Geostatistical description of hydraulic properties in karst
aquifers, a case study in the Edwards Aquifer. In Charbeneau, R. J ., ed.,
Groundwater Management, Proceedings of the International Symposium: Water
Resources Engineering Division, American Society of Civil Engineers, American
Society of Civil Engineers, p. 193-198.

Mace, Robert E. 1997. Determination of Transmissivity from Specific Capacity
Tests in a Karst Aquifer. Groundwater: Vol. 35, No. 5, p. 738-742.

Mace, R. E. In review. Using specific-capacity data in hydrogeologic investigations.
UT.

Mandle, R. 2006. IGW WHPA Delineation Program, Users Manual v1.0. MDEQ.

MDEQ (Michigan Department of Environmental Quality), 2006. Wellogic System.
Available at http://www.deq.state.mi.us/wellogic/main.html. Accessed in April
2006.

Mullican, W. F., 111, N. D. Johns, and A. E. Fryar. 1997. Playas and recharge of the
Ogallala aquifer on the southern High Plains of Texas- an examination using
numerical techniques. University of Texas at Austin, Bureau of Economic Geology
Report of Investigations No. 242, 72 p.

Myers, B. N., 1969. Compilation of results of aquifer tests in Texas. Texas Water
Development Board Report 98, 532 p.

Nicholas, J. R., G. L Rowe, and J. R.Brannen. 1996. Hydrology, Water Quality,
and Effects of Drought in Monroe County, Michigan. Water-Resources
Investigations Report 94-4161. US Geological Survey.

Nicholas, J. R., S. P.Blurner, and R. M. McGowan. 2001. Comparison of
Hydrologic Data from Monroe County, Michigan, 1991-2001. US Geological
Survey Open File Report 01-498, Lansing, Michigan.

Peerless-Midwest, Incorporated. 2002. Kalamazoo Michigan: Water Pumping
Stations 25 and 39 Groundwater Flow Model and Capture Zone Delineations.

Merle GP and DC Wiggert. 2002Mechanics of Fluids, 3nd Ed. Prentice-Hall.

Prudic, D. E., 1991, Estimates of hydraulic conductivity from aquifer-test analyses
and specific-capacity data, Gulf Coast regional aquifer systems, south-central
United States. U. S. Geological Survey, Water-Resources Investigations, WRI 90-
4121, 38 p.

107

51.

52.

53.

54.

55.

56.

57.

58.

59.

60.

61.

62.

63.

Razack, M. and D. Huntley, 1991. Assessing transmissivity from specific capacity
in a large and heterogeneous alluvial aquifer. Ground Water, vol. 29, no. 6, p. 856-
861.

Reeves, H. W., K. V. Wright, and J. R. Nicholas. 2003. Hydrogeology and
Simulation of Regional Ground-Water-Level Declines in Monroe County,
Michigan. Water-Resources Investigations Report 03-4312. US Geological Survey.

State of Michigan Public Act 148 Groundwater Inventory and Map (GWIM)
Project. March 2006. Technical Report.

Theis, C.V., R.H. Brown, and RR. Meyer. 1963. Estimating the transmissibility of
aquifers from the specific capacity of wells, in Bentall, Ray, Methods of

determining permeability, transmissibility, and drawdown. US. Geological Survey
Water-Supply Paper 1536-1, p. 331—341.

Thomasson, H. J ., F.H. Olmstead, and ER. LeRoux. 1960. Geology, water
resources, and usable ground water storage capacity of part of Solano County, CA.
US. Geological Survey Water Supply Paper 1464, 693 p.

Using ArcGIS Geostatistical Analyst. 2003.

Vanlier, K.E., and ML. Wheeler. 1968. Analog simulation of ground-water
development of the Saginaw Formation, Lansing Metropolitan area, Michigan. Tri-
County Planning Commission, Lansing Ground-Water Report, 40 p.

Vanlier, K.E., and ML. Wheeler. 1968. Ground-water potential of the Saginaw
Formation in the Lansing metropolitan area, Michigan. US. Geological Survey
Open-File Report, unnumbered, 40 p.

Vanlier, K.E., W.W. Wood, and J .O. Brunett. 1973. Water-supply deve10pment and
management alternatives for Clinton, Eaton, and Ingham County, Michigan. US.
Geological Survey Water-Supply Paper 1969, 111 p.

Wang, HF. and MP. Anderson. 1982. Introduction to Groundwater Modeling:
Finite Difference and Finite Element Methods. Academic Press, San Diego.

Wayland, K.G, D.W. Hyndman, D.F. Boutt, B.C. Pijanowski, D.T. Long, 2002,
Modeling The Impact Of Historical Land Uses On Surface Water Quality Using

Ground Water Flow And Solute Transport Models, Lakes and Reservoirs, (7), p.
189-199.

Webster, R. Geostatistics for environmental scientists. 2001. John Wiley & Sons:
New York.

Williams & Works. 2000. Wellhead Protection Area Delineation: City of Clare.

108

64.

65.

Zheng, Chunmiao, and GD. Bennett. 2002. Applied Contaminant Transport
Modeling. Wiley Interscience.

Zhang, Dongxiao. 2001. Stochastic Methods for Flow in Porous Media: Coping
with Uncertainties. Academic Press.

109

   

Ill‘llllllllllllll