Special Publication SJ96SP5
WATER RESOURCE ALLOCATION AND QUALITY OPTIMIZATION MODELING
FINAL REPORT
July 5, 1995
Prepared for
Saint Johns River Water Management District
P.O. Box 1429
Palatka, Florida 321781429
Patrick C. Burger Graduate Research Assistant
Mark E. Holland Graduate Research Assistant
Kirk Hatfield, Ph.D. Principal Investigator
Donald W. Hearn, Ph.D. CoPrincipal Investigator
Wendy Graham, Ph.D. CoPrincipal Investigator
Department of Civil Engineering
University of Florida
Gainesville, Florida 32611
EXECUTIVE SUMMARY
PROJECT SCOPE
Due to the increasing demands placed on Florida's water resources, the state of
Florida adopted legislation in 1990 to improve water resource management and to direct
future growth through planning programs. This legislation requires each Water
Management District to completely evaluate their water needs and sources through the
year 2010 and delineate critical areas identified as water resource problems. Once
completed, Districts were expected to develop possible alternative water supply scenarios
which avoid adverse or otherwise unacceptable changes in the environment or the
availability of water.
This report presents systematic modeling methods for determining optimum water
supply strategies that satisfy various environmental and hydrological requirements. Five
site specific water resource allocation optimization models were developed for Volusia
County, Florida and were executed to investigate a variety of management objectives.
These optimization models incorporate both water quantity and quality aspects of water
resource management to determine optimum ground water allocation strategies which
satisfy future water service demands and minimize adverse environmental impacts at
specified areas. These areas include sensitive wetlands where projected water table
declines are predicted to induce a high level of vegetative harm and well fields where
excessive withdrawal is predicted to cause a degradation in water quality due to salt
water intrusion or upcoming.
The five optimization models were designed to elucidate water resource allocation
strategies for water service areas that would:
1) Satisfy the water demands of both municipal and agricultural water demands.
2) Explore development of both existing and proposed ground water supply
areas.
3) Select wastewater effluent as a feasible supply to supplement agricultural
demands.
Water quality aspects were incorporated in two models by constraining chloride
concentrations changes at wells while simultaneously minimizing the maximum
drawdown at sensitive wetland areas in one model and by minimizing relative chloride
concentration increases at wells in a second model.
The optimization models were a product of a research project that sought to
accomplish two goals. The first was to obtain a broad understanding of numerical
optimization modeling and its recent applications to ground water resource management.
The second goal was to demonstrate optimum resource allocation modeling at a selected
site that would satisfy specified environmental and hydrological requirements. In order
to achieve the established two goals and fulfill the scope of the project, the following
three objectives were accomplished:
1) Review of current literature of numerical optimization modeling with respect
to ground water resource management.
2) Construct site specific optimization models capable of generating water
resource allocation scenarios which satisfy projected demands and environmental
constraints for the year 2010.
3) Generate and summarize various optimum water resource allocation scenarios
to meet specific objectives.
LITERATURE REVIEW
To meet the first objective a literature review was conducted to investigate
applications of ground water resource management where tools of optimization were
applied to demonstrate optimum scenarios of resource allocation. Chapter 2 presents a
review of various water management models; a number of existing ground water
optimization models were summarized in tabular form. Based on knowledge acquired
through the literature review and the size and complexity of the problem, it was
determined that the most efficient combined simulation/optimization model would be one
developed using of the unit "response matrix" technique. This method consists
incorporating a matrix of coefficients which represent the response of ground water to
a specified change in withdrawal rate.
PROJECT SIMULATION MODELS
Chapter 3 presents details on the threedimensional ground water flow and solute
transport models used to generate the steadystate ground water system responses (i.e.,
predicted head values, solute concentrations, etc.) to various stresses (i.e., specified
pumping and recharge rates) in the study area; these system responses are generated as
needed to create response matrices that are later used to construct the larger optimization
models. The study area includes most of Volusia County where in order to alleviate the
problem of saltwater intrusion and satisfy water service demands, the general trend since
the 1950's has been to locate additional wells to the west toward central Volusia County.
The flow simulation model as originally acquired from the SJRWMD simulates the
aquifer systems of Volusia County using MODFLOW (McDonald and Harbaugh 1988).
This regional simulation model was revised to simulate flow in the study area. The
transport model for the same study area was also developed by the District using the
program DSTRAM (Huyakorn and Panday 1991).
WATER SUPPLY NEEDS AND SOURCES
Chapter 4 presents details on water supply needs and sources as currently viewed
by SJRWMD under the District Water Management Plan. Projections of future water
use for the year 2010 are presented that were based on historical trends, local
government comprehensive plans, and direct communication with both public and private
public supply utilities (SJRWMD, 1992). Major components of the water uses categories
include municipal and agricultural service area demands. Some of the larger projected
municipal demand increases include an increase in the Port Orange service area in excess
of 70 percent of the current use, an increase in Daytona Beach service area of 62
percent, and an increase of 400 percent over current uses in the Smyrna and Samsula
areas. The overall demand in the Volusia County area is projected to increased by
approximately 75 percent from 1988 to 2010. Also included in Chapter 4 are maps and
tables that depict wastewater treatment plant information used to incorporate water reuse
in the development of the ground water management models.
GROUND WATER OPTIMIZATION MODELS
Chapter 5 presents the five aquifer optimization models developed to investigate
optimal allocation of ground water to meet year 2010 demand in Volusia County project
area. These models were formulated to investigate future water allocation strategies
assuming feasible withdrawal scenarios must meet or exceed projected water service
areas demands and not exceed available water supplies. It was assumed with these
models that adverse environmental effects could be minimized at specific locations by
constraining pressure head changes (i.e., drawdown) to meet specified environmental
goals or standards. These models are essentially numerical optimization models
comprised of an objective function, decision variables, and constraints. Of the five
optimization models, three models dealt strictly with hydraulic constraints and two
models incorporated both hydraulic and water quality constraints. Objective functions
included the following: 1) minimizing the maximum drawdown value, 2) maximizing
the minimum pressure head, 3) minimizing average drawdown over all sensitive wetland
areas, 4) minimizing the maximum drawdown while limiting chloride concentrations, and
5) minimizing relative chloride concentration increases.
Data supplied by SJRWMD and gathered during the project were used to
formulate the allocation models (i.e., in the development of an objective function and
management constraints). This included well data specifying location, capacity,
withdrawal rates, corresponding service areas, and additional data identifying potential
well locations. Other data included information concerning effluent rates, capacities, and
locations of wastewater treatment plants that could supplement the demand of agricultural
areas. Water supply demand data consisted of municipal and agricultural water volumes
used in year 1988 and projected water quantity requirements for the year 2010.
The District also supplied environmental impact data specifying locations where
drawdown was expected to induce adverse effects on wetland vegetation. These target
locations were used to constrain the optimization model to identify resource allocation
scenarios that would achieve minimum environmental impacts.
Aquifer system response constraints were developed from steadystate unit
response matrices containing influence coefficients generated from ground water flow and
water quality simulation models described in Chapter 3. These steadystate unit response
matrices represent steadystate changes in pressure head and water chloride
concentrations with respect to pumpage. A ground water flow simulation model was
executed once for each existing and potential well location to determine the head
response at specific locations. A similar approach was taken with the solute transport
simulation model to determine chloride changes at wells due to pumpage. Water level
and water quality unit response matrices were developed for both municipal and
agricultural wells. These response matrices summarize the influence of each well upon
itself, upon all other wells, and also upon specified control points throughout Volusia
County.
The GAMS program (Brooke, Kendrick, and Meeraus 1989) was used to
formulate each allocation model into a linear program optimization model. A model
specific objective function and various constraints (i.e., water supply, demand, and
environmental constraints) constitute the basic framework of each resource allocation
model. All the models contain specific constraints that define allowable connections
between water supplies and demand areas. Three of the models incorporate data from
a unit response matrix created for various hydraulic management constraints, while two
models incorporates two response matrices needed to construct both hydraulic and water
quality management constraints. Once formulated the optimization models were executed
using the linear program solution algorithm available in the GAMS software.
Year 2010 water resource allocation strategies were examined with the
optimization models. Results of each model revealed where future well fields should be
developed and to what extent future and existing wells should be utilized. Aquifer
responses induced by the optimized allocation strategies were compared to those induced
by the projected year 2010 allocation strategy. Predicted pumping strategies identified
from each optimization model were then reviewed and tabulated. To verify optimal
resource allocation strategies identified by the five models the ground water flow and
transport simulation models were used to simulate the optimal pumping scenarios
generated. Predicted pressure heads, drawdowns, and chloride concentrations from the
simulation models were then compared to optimization model results.
SUMMARY AND CONCLUSIONS
Chapter 6 summarizes the overall research effort and presents several salient
conclusions. This work endeavor demonstrated the use of optimization modeling as a
valuable tool for the management of water resources. Optimization modeling elucidated
the potential for an abatement of adverse environmental impacts associated with declining
water table elevations and this was verified when the optimized allocation strategies were
simulated with ground water flow and transport models. It was also shown that the
minimum discharge constraint is an important determinant of the identified optimum
allocations. This constraint requires a minimum flow of at least 50 percent of the 1988
ground water pumping rate from all active wells. Several of the optimal allocation
scenarios incorporate new wells in areas of Port Orange West well field, Daytona Beach
West South Daytona water treatment plant well field, and Ormond Beach State Road
40 / Hudson well fields. In general, each model identified a scenario that decreased
ground water pumpage at the Daytona Beach West water service area and expanded
pumpage in areas east and west where the native vegetation is less sensitive to ground
water withdrawals.
Sensitivity analyses were performed with the optimization models to determine
where management strategies could be changed to induce improvements in ground water
levels. These analyses identified areas were the balance between supplies and demands
could be changed to decrease the environmental impacts of ground water withdrawals.
For example, it was shown that a ten percent decrease in the demand at the DBS
(Daytona Beach West South TP) water service area would induce a water level
improvement of 17 percent. These types of analyses give the water resource manager
valuable information on where conservation efforts and new well developments should
be pursued and let the manager quickly determine how the resource allocation scenario
given by each optimization model responds to changes in projected demands and
withdrawal limits.
TABLE OF CONTENTS
EXECUTIVE SUMMARY .................................. i
LIST OF TABLES .................... ................. xi
LIST OF FIGURES .................... ............... xiii
LIST OF VARIABLES ................... .............. xvii
1.0 INTRODUCTION ................... ................ 1
1.1 STUDY AREA BACKGROUND ...................... 4
1.2 TECHNICAL APPROACH .......................... 8
2.0 REVIEW OF CURRENT LITERATURE ................... 11
2.1 COMBINED SIMULATION/OPTIMIZATION MODELS ..... 11
2.1.1 Embedding Method ........................ 13
2.1.2 Unit Response Matrix Method ................. 19
2.2 SALTWATER INTRUSION MODELS ................. 32
2.3 INTERFACE PROGRAMS ........................ 36
2.4 SUMMARY ................... ............. 38
3.0 GROUND WATER FLOW SIMULATION MODEL ............. 44
3.1 REGIONAL FLOW SIMULATION MODEL ............. 47
3.2 FLOW SIMULATION MODEL FOR PROJECT AREA ...... 48
3.3 SOLUTE TRANSPORTATION MODEL ............... 54
4.0 WATER SUPPLY NEEDS AND SOURCES .................. 55
5.0 GROUND WATER OPTIMIZATION MODELS ............... 69
5.1 OPTIMIZATION MODEL DECISION VARIABLES ........ 70
5.2 OPTIMIZATION MODEL OBJECTIVE FUNCTIONS ....... 71
5.3 OPTIMIZATION MODEL CONSTRAINTS .............. 73
5.3.1 Aquifer Response Constraints .................. 73
5.3.2 Management Constraints ................. ... 77
5.3.3 Chloride Concentration Constraints ............. 80
5.3.4 Nonnegativity Constraints ................. ... 82
5.4 OPTIMIZATION MODEL FORMULATION ............. 83
5.4.1 General Optimization Model Formulation ......... 83
5.4.2 Specific Optimization Model Formulations ........ 86
5.5 PROCEDURE FOR IDENTIFYING STRATEGIES ......... 86
5.6 RESULTS AND DISCUSSION .... ................. 90
5.6.1 M odel 1 ... ....... ...... ....... ...... 90
5.6.2 Model 2 ....... ...................... 104
5.6.3 M odel3 .................. ...... ....... 118
5.6.4 Model 4 ....................... ....... 124
5.6.5 M odel 5 ............................... 138
5.6.6 General Observations ....................... 150
5.6.7 Sensitiviy Analysis ........................ 167
6.0 SUMMARY AND CONCLUSIONS ......................... 177
7.0 APPENDICES ..................... ....... .......... 181
7.1 PROGRAM CODES FOR DETERMINING RESPONSE
MATRICES ........... ....................... 182
7.2 OPTIMIZATION MODEL FORMULATION USING GAMS ... 190
7.2.1 Sets ............. .................... 190
7.2.2 Tables ..................... .......... 191
7.2.2.1 Source to Demand Link Tables ........ 191
7.2.2.2 Influence Coefficient Tables ............. 192
7.2.3 Parameters ........ ............ ......... 195
7.2.4 Variables ............ .... ............. 195
7.2.5 Equations ..... .......... ............... 196
7.2.6 GAMS Optimization Modeling ................ 196
7.3 GAMS OPTIMIZATION MODEL INPUT FILE ........... 197
8.0 REFERENCES ................... ................. 208
LIST OF TABLES
Tables page
2.1 Ground water optimization models .................... 41
4.1 Volusia County municipal well data for project area ............ 56
4.2 Volusia County agricultural well data for project area ............ 61
4.3 Year 1988 and projected year 2010 demand rates for municipal
water service areas ................... ............. 65
4.4 Year 1988 and projected Year 2010 demand rates for agricultural
water service areas ................... .............. 65
4.5 Volusia County wastewater treatment plant data for project area . 67
5.1 Specific optimization model formulations . . . ..... 86
5.2 Municipal water supply grid cell locations and withdrawal rates . 91
5.3 Agricultural water supply grid cell locations and withdrawal rates .. 95
5.4 Calculated pressure head and drawdown values while
minimizing maximum drawdown at sensitive wetland control
points (Results from initial formulation) . . . ...... ...... 97
5.5 Calculated pressure head and drawdown values while
minimizing maximum drawdown at sensitive wetland control
points (Results from revised model 1) ................ .... 101
5.6 Calculated pressure head and drawdown values while
maximizing minimum head at sensitive wetland control
points (Results from revised model 2) ..................... 119
5.7 Calculated pressure head and drawdown values while
minimizing average head at sensitive wetland control
points (Results from revised model 3) ................ .... 125
Tables page
5.8 Calculated pressure head and drawdown values while
maximizing minimum head at sensitive wetland control
points (Results from revised model 4) ..................... 130
5.9 Calculated chloride concentrations at municipal well grid cells
while minimizing maximum drawdown at wetland control points . 134
5.10 Calculated pressure head and drawdown values while
maximizing minimum head at sensitive wetland control
points (Results from revised model 4) ..................... 142
5.11 Calculated chloride concentrations at municipal well grid cells
while minimizing maximum drawdown at wetland control points .. 145
5.12 Pressure head, drawdown, and difference between projected and
optimized values at sensitive wetland control points ............. 154
5.13 Shadow prices for municipal well grid cells .................. 168
5.14 Shadow prices for agricultural well grid cells ................. 174
5.15 Shadow prices for wastewater treatment plants ................ 174
5.16 Shadow prices for municipal water service areas ............... 174
5.17 Shadow prices for agricultural water service areas .............. 175
LIST OF FIGURES
Figure page
1.1 Location of Volusia County, Florida within the St. Johns River Water
Management District .................................... 5
1.2 Hydrologic cross sections of Volusia County, Florida ............. 7
3.1 Map of Volusia County, Florida ......................... 49
3.2 Volusia County regional and subregional numerical model grid
discretizations ..................................... 50
3.3 Project Flow Model Cell Location Map ................... 52
4.1 Location of municipal water service areas within Project area ...... 60
4.2 Location of agricultural water service areas and wastewater treatment
plants with respect to numerical model ................. .... 63
4.3 Areas of potential impact to plant communities resulting from projected
changes in ground water withdrawals between 1988 and 2010 ....... 68
5.1 Example of nonlinear aquifer response with increased discharge rate at a
surficial aquifer control point ........................... 88
5.2 Correlation between simulation model results and optimization model 1
initial formulation results a) Pressure heads for year 2010;
b) Drawdowns from year 1988 to 2010 . . . ...... ...... 99
5.3 Correlation between simulation model results and revised optimization
model 1 results a) Pressure heads for year 2010;
b) Drawdowns from year 1988 to 2010 . . . ...... ....... 103
5.4 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for Port Orange East water service area . 105
5.5 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for Port Orange West water service area . 106
Figure
5.6 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for Daytona Beach water service area . ... 107
5.7 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for Spruce Creek water service area . ... 108
5.8 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for Holly Hill water service area . . .. 109
5.9 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for New Smyrna Beach water service area . 110
5.10 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for Ormond Beach water service area ........ 111
5.11 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for Tymber Creek & The Trails water service areas 112
5.12 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for agricultural water service areas 1, 3, and 9 .113
5.13 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for agricultural water service areas 2, 4, and 5 .114
5.14 Year 1988, projected year 2010, and model 1 optimized year 2010
withdrawal strategies for agricultural water service areas 6, 7, and 8 .115
5.15 Flow simulation model pressure head results in surficial aquifer
due to optimization model 1 withdrawal strategy for year 2010 . 116
5.16 Flow simulation model drawdown results in surficial aquifer
due to optimization model 1 withdrawal strategy for year 2010 . 117
5.17 Correlation between simulation model results and revised optimization
model 2 results a) Pressure heads for year 2010;
b) Drawdowns from year 1988 to 2010 .................... 121
5.18 Flow simulation model pressure head results in surficial aquifer
due to optimization model 2 withdrawal strategy for year 2010 . 122
page
5.19 Flow simulation model drawdown results in surficial aquifer
due to optimization model 2 withdrawal strategy for year 2010 ...... 123
5.20 Correlation between simulation model results and revised optimization
model 3 results a) Pressure heads for year 2010;
b) Drawdowns from year 1988 to 2010 .................... 127
5.21 Flow simulation model pressure head results in surficial aquifer
due to optimization model 3 withdrawal strategy for year 2010 . 128
5.22 Flow simulation model drawdown results in surficial aquifer
due to optimization model 3 withdrawal strategy for year 2010 . 129
5.23 Correlation between simulation model results and revised optimization
model 4 results a) Pressure heads for year 2010;
b) Drawdowns from year 1988 to 2010 .................... 133
5.24 Chloride concentration correlation between simulation model and
optimization model 4 ................... ............. 138
5.25 Flow simulation model pressure head results in surficial aquifer
due to optimization model 4 withdrawal strategy for year 2010 . 139
5.26 Flow simulation model drawdown results in surficial aquifer
due to optimization model 4 withdrawal strategy for year 2010 . 140
5.27 Chloride concentrations in upper Floridan aquifer due to optimization
model 4 water allocation strategy for year 2010 . . ..... 141
5.28 Correlation between simulation model results and revised optimization
model 5 results a) Pressure heads for year 2010;
b) Drawdowns from year 1988 to 2010. . . ........ 149
5.29 Chloride concentration correlation between simulation model and
optimization model 5 ................................ 150
5.30 Flow simulation model pressure head results in surficial aquifer
due to optimization model 5 withdrawal strategy for year 2010 . 151
5.31 Flow simulation model drawdown results in surficial aquifer
due to optimization model 5 withdrawal strategy for year 2010 . 152
page
Figure
5.32 Chloride concentrations in upper Floridan aquifer due to optimization
model 5 water allocation strategy for year 2010 . . ..... 153
5.33 Flow simulation model pressure head results in surficial aquifer
due to projected year 2010 withdrawal strategy . . ..... 158
5.34 Flow simulation model drawdown results in surficial aquifer
due to projected year 2010 withdrawal strategy . . ..... 159
5.35 Chloride concentrations in upper Floridan aquifer due to the
projected year 2010 water allocation strategy . . . ..... 160
5.36 Difference in induced pressure head between the projected year 2010
strategy and the strategy identified by model 1 (projectedoptimized) .162
5.37 Difference in induced pressure head between the projected year 2010
strategy and the strategy identified by model 2 (projectedoptimized) .163
5.38 Difference in induced pressure head between the projected year 2010
strategy and the strategy identified by model 3 (projectedoptimized) .164
5.39 Difference in induced pressure head between the projected year 2010
strategy and the strategy identified by model 4 (projectedoptimized) .165
5.40 Difference in induced pressure head between the projected year 2010
strategy and the strategy identified by model 5 (projectedoptimized) .166
LIST OF VARIABLES
ca, is the aquifer influence coefficient defining pressure head change at
each sensitive wetland area control point due to a change in discharge
rate at each municipal well grid cell i.
3,j is the aquifer influence coefficient defining pressure head change at
each sensitive wetland area control point due to a change in discharge
rate at each agricultural well n.
BOTELEVh is the bottom elevation of the surficial aquifer in feet from mean
sea level at each well grid cell control point h.
CA, is the capacity rate in cfd of each agricultural well grid cell i.
CCh is the chloride concentration in tenthousandths of a milligrams per liter
(mg/1) at each well grid cell control point h.
CCOh is the initial chloride concentration in mg/1 at well grid cell control
points h.
CIh is the increase in chloride concentration in tenthousandths of a
milligrams per liter (mg/1) at each well grid cell control point h.
CLh is the maximum chloride concentration limit in mg/1 of each well grid
cell control point h.
CM, is the capacity rate in cfd of each municipal well grid cell i.
CW. is the capacity rate in cfd of each wastewater treatment plant i.
DA, is the demand rate in cfd of each agricultural water service area o.
DDj is the drawdown at sensitive wetland control point.
xvii
DDPRIVj is the drawdown in feet at each sensitive wetland area control point
due to private wells not incorporated in the optimization process.
DDW, is the drawdown in millionths of a foot at each well grid cell h.
DMk is the demand rate in cfd of each municipal water service area k.
Yi,h is the aquifer influence coefficient defining pressure head change at well
grid cell control point h due to a change in discharge rate at each
municipal well grid cell i.
h defines all the well grid points incorporated in the optimization model.
HDj is the pressure head in millionths of a foot at each sensitive wetland
area control point.
HDWh is the pressure head in millionths of a foot at each well grid cell h.
HOj is the initial pressure head in feet at each sensitive wetland area control
point.
HWOh is the initial pressure head in feet at each well grid cell h.
i defines the 120 municipal well grid cells in the optimization model.
j defines the 100 grid cell points in the optimization model where there is
a high potential of harm to vegetation in sensitive wetland areas.
k defines the 14 Municipal water service areas in the optimization model.
m defines the five wastewater treatment plants in optimization model.
n defines the 28 agricultural well grid cells in the optimization model.
xviii
o defines the 9 different agricultural water service areas in the
optimization model.
40n,h is the aquifer influence coefficient defining chloride concentration
change at each well grid cell control point h due to a change in
discharge rate at each agricultural well n.
QAn,o is the discharge rate in cfd of each agricultural well grid cell n which
supplies each agricultural water service area o.
QAO, is the initial discharge rate in cfd of each agricultural well grid cell n.
QAT, is the total discharge rate in cfd of each agricultural well grid cell n.
QMik is the discharge rate in cubic feet per day (cfd) of each municipal well
grid cell i which supplies each municipal water service area k.
QMO, is the initial discharge rate in cfd of each municipal well grid cell i.
QMT, is the total discharge rate in cfd at each municipal well grid cell i.
QW,o is the effluent reuse rate in cfd of each wastewater treatment plant m
which supplements each agricultural water service area o.
QWTm is the total effluent reuse rate in cfd of each wastewater treatment
plant m.
RCIh is the relative chloride concentration increase in tenthousandths of a
miligram per liter (mg/1) at each well grid cell control point h.
SERVEi,k designates which municipal well grid cells i supply which municipal
xix
water service areas k.
SERVE2,,o designates which agricultural well grid cells n supply which
agricultural water service areas o.
SERVE3,,o designates which wastewater treatment plants m can supplement the
demand of which agricultural water service areas o.
TDDj is the drawdown in feet at each sensitive wetland area control point.
TDDWh is the drawdown in feet at each well grid cell h.
THDj is the total pressure head in feet at each sensitive wetland area control
point j.
THDWh is the total pressure head in feet at each well grid cell h.
On,h is the aquifer influence coefficient defining pressure head change at
well grid cell control point h due to a change in discharge rate at each
agricultural well grid cell n.
TCIh is the total increased chloride concentration, respectively, in mg/l at
each well grid cell control point h.
TCCh is the total chloride concentration in mg/1 at each well grid cell control
point h.
4i,h is the aquifer influence coefficient defining chloride concentration
change at each well grid cell control point h due to a change in
discharge rate at each municipal well grid cell i.
XX
CHAPTER 1
1.0 INTRODUCTION
In recent years, population growth, expanding agricultural and industrial activities,
and rapid urban development have significantly increased the demand for clean, fresh
water in the state of Florida. Within the boundaries of the St. Johns River Water
Management District (SJRWMD), demand for public water has increased 66% from 1975
to 1990 and is predicted to increase another 120% by the year 2010 (SJRWMD 1992).
This increased demand presents an ominous threat to ground water sources in terms of
both water quantity and quality. Along coastal areas of the District, increased
development has generated water quality problems related to high chloride concentrations
due to saltwater intrusion and localized upcoming. Sensitive wetland areas have also
been adversely affected by recent declines in the water table.
Due to the increasing demands placed on Florida's water resources, the state of
Florida adopted legislation in 1990 to improve water resource management and to direct
future growth through planning programs. This legislation requires each Water
Management District to completely evaluate their water needs and sources through the
year 2010 and delineate critical areas identified as water resource problems. Districts
are then expected to develop possible alternative water supply scenarios which avoid
adverse or otherwise unacceptable changes in the environment or the availability of
water.
2
A study of regional water supply needs and sources, could identify many possible
resource allocation plans that meet the needs of the District. The distribution of available
supplies to service areas is an allocation problem and can have several feasible solutions.
It is also likely that the number of plans will be too extensive to permit detailed
examination of each and every scenario. In order to identify a concise subset of water
allocation plans which best meet environmental and developmental goals of the District,
optimization modeling can be incorporated into the decision or plan elucidation process.
This project involves the development of a systematic method of determining
optimum water supply strategies that satisfy various environmental and hydrological
requirements. The purpose of this type of water supply strategy is to optimize the
pattern of water supply development and usage to meet projected needs. This resource
management problem requires the use of optimization modeling to identify desirable
scenarios of resource allocation; otherwise, resources may not be used in the most
effective and efficient manner. When environmental impacts are also incorporated, the
allocation problem expands to include identifying feasible scenarios that must also satisfy
environmental constraints (i.e., ground water quality standards, minimum water levels).
To balance projected needs against available sources, it is possible that the management
problem may become one of balancing projected development against adverse
environmental impacts.
The objective of this project is to accomplish two goals. The first is to obtain a
broad understanding of numerical optimization modeling and its recent applications to
ground water resource management. The second goal is to demonstrate optimum
resource allocation modeling at a selected site that would satisfy specified environmental
3
and hydrological requirements. The demonstration will give the SJRWMD essential
knowledge and experience to incorporate optimization technology into the District's
decisionmaking framework. In order to achieve the established two goals and fulfill the
scope of the project, the following three objectives were accomplished:
1) Review current literature of numerical optimization modeling with respect
to ground water resource management.
2) Construct site specific optimization models capable of generating water
resource allocation scenarios which satisfy projected demands and environmental
constraints for the year 2010.
3) Generate and summarize various optimum water resource allocation
scenarios which satisfy specific objectives.
To meet the first objective the literature review was limited to applications of
ground water resource management and only where tools of optimization were applied
to demonstrate optimum scenarios of resource allocation. Formulations of various water
management models were reviewed and a number of existing ground water optimization
models were summarized in tabular form. Based on knowledge acquired through the
literature review and the size and complexity of the problem, it was determined that the
most efficient simulation/optimization model would be one developed using the unit
"response matrix" technique. This method consists of incorporating a matrix of
coefficients which represent the response of ground water to a specified change in
withdrawal rate. To achieve the second and third objectives it was necessary to select
a study site within the district boundaries.
1.1 STUDY AREA BACKGROUND
Based on the extensive ground water flow and contaminant transport modeling
effort completed by the SJRWMD, four distinct regions within the district were identified
as candidates for project study. These areas included the Volusia County region, the
Wekiva River Basin region, the EastCentral Florida region, and the Geneva
Area/Seminole County region. The regional models and areas were evaluated with
respect to model accuracy and complexity, computer requirements and efficiency, data
availability, and applicability. Based on this review, A study area within Volusia County
Florida was selected as the site to apply the numerical optimization modeling techniques
of this project, shown in Figure 1.1.
The study area includes most of Volusia County, thus a general geographical and
hydrological characterization of this county is in order. Volusia County alone covers an
area of approximately 1,200 square miles in the eastcentral region of Florida. As shown
in Figure 1.1, the county is bounded by the Atlantic Ocean to the east, the St. Johns
River to the West, Flagler and Putnam Counties to the north, and Brevard County to the
south. Although the county consists of approximately 120 lakes larger than 5 acres in
size, ground water supplies are currently the sole means of meeting public water demands
(Knochenmus and Beard 1971, Kimrey 1990).
Most of the population of Volusia County occupies the region in and near the
coastal cities of Daytona Beach, New Smyrna Beach, and Ormond Beach. Extensive
development of these coastal areas has resulted in increased demands for ground water
and an encroachment of the freshwater/saltwater interface. In order to alleviate the
problem of saltwater intrusion and satisfy water service demands, the general trend since
PUTNAM
FLAER\ Project Area
MARION
VOLUSIA
LAKE
SEMINOLE
BREVARD
10 miles
0 5 kilometers
Figure 1.1
Location of Project in Volusia County, Florida within the
St. Johns River Water Management District
I
6
the 1950's has been to locate additional wells to the west towards central Volusia County.
The hydrogeologic system of Volusia County is made up of two aquifer systems,
the surficial and the Floridan aquifers, separated by an intermediate confining layer (See
Figure 1.2). The surficial aquifer is the uppermost formation, consists of silts, clays,
cemented shell, and quartz sands, and is considered to be unconfined. The water table
is usually at or near the surface in lowland and flatland regions and is generally a
suppressed image of the ground level in highland regions. Precipitation, lakes, wetlands,
and irrigation are the main sources of recharge into the surficial aquifer, which ranges
in thickness between 50 and 100 feet. Via the intermediate confining layer, leakage can
occur in and out of the aquifer depending on the difference in potentiometric head
between the surficial and Floridan aquifer systems (Tibbals 1990).
The Floridan aquifer system contains the Upper Floridan and Lower Floridan
aquifers and is divided by a middle confining unit. This confining layer of low
permeability is located completely within the Avon Park limestone geologic formation.
The lower portion of the Avon Park Formation along with the Oldsmar Formation
comprise the Lower Floridan aquifer, which consists mainly of saline water. The Upper
Floridan aquifer consists of the Ocala Limestone Formation and the upper portion of the
Avon Park Formation (Miller 1986, Tibbals 1990). Although the Upper Floridan aquifer
contains brackish water in the St. Johns River Valley, near the Atlantic coast, and north
of Volusia County line, this ground water supply is the main means of meeting public
water demands (Kimrey 1990). The entire thickness of the Floridan aquifer system
ranges from approximately 1,800 to 2,300 feet in Volusia County.
Separating the surficial aquifer from the Floridan aquifer throughout most of the
Section AA
100
50.
MSL
50
100
150
200
Feet
Feet
KEY
LSURFICIAL AQUIFER
16ONFINING UNIT
OFLORIDAN AQUIFER
Figure 1.2 Hydrologic cross sections of Volusia County, Florida.
(source: McKenzieArenberg, 1989)
8
county is the intermediate or upper confining layer. This layer consists of clay or silty
sand of the Miocene to Pleistocene age and has a thickness range of 0 to 60 feet. The
confining unit in the western portion of the county is thinner than the eastern portion and
is sometimes totally absent. Although the low permeability layer is leaky, it is capable
of confining the pressurized water in the Floridan aquifer system (Phelps 1990).
1.2 TECHNICAL APPROACH
Once the specific study site was selected, it was possible to initiate development
of the water resource allocation models. These models were essentially numerical
optimization models comprised of an objective function, decision variables, and
constraints. A total of five optimization models were formulated to examine water
resource allocations for year 2010, these included three models with hydraulic constraints
and two models that incorporate both hydraulic and water quality constraints. The five
model specific objective functions included: 1) minimizing the maximum drawdown
value, 2) maximizing the minimum pressure head, 3) minimizing average drawdown over
all sensitive wetland areas, 4) minimizing the maximum drawdown while limiting
chloride concentrations, 5) minimizing the maximum chloride concentration, and 6)
minimizing the maximum relative chloride concentration increase. From the literature
review, it was determined that the formulation of these water resource allocation models
would be predicated on the unit response matrix approach as described in the Literature
Review.
Data supplied by SJRWMD and gathered during the project were used to develop
the allocation models consisting of an objective function and management constraints.
9
The data supplied included well data specifying location, capacity, withdrawal rates,
corresponding service areas, and additional data identifying potential well locations.
Other data included information concerning effluent rates, capacities, and locations of
wastewater treatment plants that could supplement the demand of agricultural areas.
Water supply demand data consisted of municipal and agricultural water volumes used
in year 1988 and projected water quantity requirements for the year 2010.
The District also supplied environmental impact data specifying locations where
drawdown was expected to induce adverse effects on wetland vegetation. These target
locations were used to constrain the optimization models. Thus, various ground water
allocation scenarios could be examined that involve minimizing environmental impacts.
Aquifer system response constraints were incorporated in each optimization
model. These constraints were developed from steadystate unit response matrices
containing influence coefficients generated from ground water flow and solute transport
simulation models. These steadystate unit response matrices represent steadystate
changes in pressure head and water chloride concentrations with respect to pumpage.
A ground water flow simulation model was executed once for each existing and potential
well location to determine the head response at specific locations. A similar approach
was taken with the solute transport simulation model to determine chloride changes at
wells due to pumpage. Water level and water quality unit response matrices were
developed for both municipal and agricultural wells. These response matrices summarize
the influence of each well upon itself, upon all other wells, and also upon specified
control points throughout Volusia County.
Pumping strategies identified from each optimization model were reviewed and
10
tabulated. To verify these optimum water allocation strategies, the ground water flow
and the solute transport simulation models were used to simulate these pumping
strategies. Predicted pressure head, drawdown, and water quality changes generated from
the simulation models were then compared to predictions given by the five optimization
models to determine how well the optimization models emulated the simulation models.
When simulated and optimized model values for drawdown, pressure head, and chloride
concentration differed severely, response matrices were recreated using initial conditions
that better approximate the optimum pumping scenarios identified by each of the original
optimization models.
CHAPTER 2
2.0 REVIEW OF CURRENT LITERATURE
Over the last three decades, ground water management models have been
formulated using a variety of computational techniques to simulate and optimize the
management of aquifer systems. These models have been produced by applying
economic theory, heuristic or intuitive procedures, optimization algorithms, hydraulic
flow equations, and complex combined simulation/optimization algorithms. Although
each of these methods has its own advantages and disadvantages, the combined
simulation/optimization approach is found to be the most powerful since it is capable of
incorporating economic, physical, and policy considerations within a single model
environment. The literature review that follows is limited to only the combined
simulation/optimization approach and discusses the various methods and applications of
this modeling technique. Also reviewed are applications where optimization modeling
has been combined with saltwater intrusion models and computer programs which
facilitate model formulation and act as an interface between the simulation and
optimization models.
2.1 COMBINED SIMULATION/OPTIMIZATION MODELS
Models developed using the simulation/optimization approach contain both
simulation equations and operations research style optimization algorithms. The
simulation equations assure that the management model correctly emulates the aquifer
12
responses to external and internal fluxes. The optimization algorithms allow the water
management objective and restrictions to be specified as algebraic equations. The
combined models are then capable of identifying ground water withdrawal scenarios
which optimize given objective functions (Peralta and Willardson, 1992).
Ground water flow simulation models alone are only able to compute pressure
heads and fluxes resulting from specified boundary conditions and pumping rates. Using
the simulation model alone to determine an optimal pumping scenario can be tedious and
error prone process; the simulation model must be executed many times with different
pumping strategies in order to determine the best scenario. Unless all possible
combinations of pumping strategies are evaluated, the optimum pumping strategy is not
likely to be identified. The simulation/optimization models, however, directly compute
the optimum strategies under identified management objectives subject to specified system
constraints (i.e. assuring the heads and fluxes are constrained within desired limits).
Optimization methods which have been applied in the field of ground water
management include linear programming, nonlinear programming, mixedinteger linear
programming, quadratic programming, dynamic programming, differential dynamic
programming, and goal programming. While these applications differ greatly in their
mathematical basis, they may be grouped into general categories according to the method
by which ground water flow equations are incorporated into the problem formulation.
Ground water flow equations have been combined with mathematical optimization
algorithms through the use of various types of approaches, methods, or techniques.
Historically, most simulation/optimization models represent the relationship between
aquifer system response and stimuli by incorporating the "embedding" method or the
13
"unit response matrix" approach. For ground water management models, pressure heads
or fluid fluxes represent the system response and withdrawal rates represent the stimuli.
The following is a review of current literature on the embedding and unit response matrix
approaches.
2.1.1 Embedding Method
Models incorporating the "embedding" technique use numerical methods such as
finitedifferences or finiteelements to transform the partial differential ground water flow
equations into a set of linear algebraic equations. These numerical methods both involve
the use of a mesh or grid to discretize the aquifer system, with each grid cell specified
as having homogenous properties. Pressure heads and discharge (and/or recharge) rates
are then depicted by flow equations formulated at each cell node in terms of control or
decision variables. These equations are then incorporated into the optimization model
as constraints.
The embedding technique with respect to ground water management was first
proposed by Aguado and Remson (1974). Their initial work demonstrated that linear
programming could be used in conjunction with finitedifference approximations to study
the physical response of ground water systems under conditions which define optimal
resource management. The optimization models were developed to determine optimum
well production scenarios which produce a maximized sum total of hydraulic heads over
the area of interest. The method was demonstrated for both confined and unconfined
aquifer systems under onedimensional, twodimensional, steadystate, and transient flow
conditions. For the transient flow case, equations were developed for each time step
and then solved simultaneously. For the unconfined aquifer system, the nonlinear
14
equations were approximated by a succession of systems of linear difference equations.
Although the authors conclude that nonlinear flow equations and irregular
boundary or initial conditions can be modeled adequately, they concede that the resulting
linear programs are typically very large.
A sitespecific application of the embedding technique was presented by Aguado
et al. (1974) to determine an optimum aquifer dewatering strategy for a proposed dry
dock. A linear programming model incorporating steadystate, finitedifference
approximations of a twodimensional unconfined aquifer was developed in order to
predict the optimal number of wells, well locations, and corresponding withdrawal rates
required to produce and maintain desired ground water levels. The optimization model
formulation permits the identification of scenarios that minimize total pumping without
restricting the number or location of utilized wells. Model execution led to a solution
in which a maximum number of wells were positioned as close to the dewatering region
as possible. To verify the accuracy of the optimization model, aquifer response results
were compared to those predicted by unconstrained numerical models and also electrical
analog ground water models (Remson, Aguado, and Remson, 1974). The authors
concluded that a finer grid spacing in the management model would produce improved
accuracy with respect to aquifer response, thus an improved optimum pumping strategy.
Alley et al. (1976) demonstrated the embedding technique by using a finite
difference approximation of twodimensional transient flow. The optimization model
determined optimum withdrawal strategies depicting well location and pumping rates for
a 20 day time period. Since an existing withdrawal strategy was specified as an initial
condition, the optimization model predicted the rate of increased withdrawal (as opposed
15
to the absolute withdrawal rate). To optimize the transient condition, submodels were
created which divided the time period into four equal time steps. These models were
then executed sequentially to produce initial conditions for the subsequent time interval.
This stepwise approach was found to reduce size dimensions of the model, and thus
computation time, when compared to the lumped transient approach, which solves all
time steps simultaneously. However, this method has been criticized due to conflict
between longrange management goals and shortterm, sequential optimization (Gorelick
1983). The availability of future water resources is affected by decisions made in earlier
time periods.
The Galerkinfinite element method was incorporated by Willis and Newman
(1977) to simulate a hypothetical aquifer consisting of heterogeneous anisotropic porous
media. The twodimensional flow equations were embedded into an optimal control
ground water problem to determine well locations and pumping rates which would satisfy
water demands at minimum cost over a number of planning periods. Since the objective
function is nonlinear, the algorithm was developed to solve a sequence of linear
programming submodels until no improvement in the objective function could be
obtained. Additional constraints and "penalty" costs were added to the management
model to minimize deviations from the desired final system state. This modification
enabled the model to find a more equally distributed pumping pattern, thus minimizing
environmental impacts due to drawdown.
Yazicigil and Rasheeduddin (1987) extended the embedding approach to include
multiaquifer systems. Both a transient, singleobjective approach and a steadystate,
multipleobjective approach were presented. For both cases, a hypothetical, leaky
16
confined twolayer aquifer system was discretized to represent the physical system.
Quasithreedimensional finite difference approximations were formulated for each active
grid cell and embedded as constraints in the linear programming model.
For the transient case, the water management objective was specified to determine
optimal well locations and corresponding withdrawal rates which would minimize
drawdowns in the two aquifer system over a four season time period. Equations for all
four time steps were effectively solved simultaneously, but an attempt to sequentially
solve the four season problem generated an infeasible solution for the fourth time step.
These results supported earlier discoveries by Alley et al. (1976) which suggested that
the use of the embedding technique with sequential linear programming may result in
unacceptable long range management decisions.
To further demonstrate the management model, the authors applied the embedding
technique to steadystate, multiobjective analysis by replacing equations representing
transient responses and adding a second objective function to the original problem.
Constraint and weighing methods were then used to solve the resultant linear
programming problem and to develop graphical tradeoff curves between total withdrawal
rate and total hydraulic head.
Jones et al. (1987) applied a modified embedding approach by incorporating a
differential dynamic programming (DDP) algorithm to obtain optimal control of an
unsteady, unconfined aquifer. DDP is a successive approximation technique for solving
optimal control problems. Equations were incorporated which 1) represented the
simulation model of the aquifer system, 2) constrained the hydraulic head and pumping
parameters such as well capacity, drawdown, and supply demands, and 3) optimized
17
objective functions such as maximize the sum of heads or minimize yield. The method
was found to reduce the dimensionality problems associated with embedding hydraulic
equations, linearize the exponential growth in computer time with respect to stages and
time periods, and cause to solution to converge quadratically. The authors concluded
DDP to be a powerful tool for management of transient ground water hydraulics.
The embedding technique has not only been applied to ground water quantity, but
also to ground water quality. Willis (1979) demonstrated the use of ground water quality
constraints in linear programming by developing an embedded algorithm for contaminant
transport equations. The technique was applied to a ground water system utilized as both
a water supply resource and storage reservoir for waste water residuals. The Galerkin
method was used to transform the flow and transport differential equations, and the
resulting equations were embedded as constraints into the management model to simulate
the aquifer system. The model was developed to determine optimum pumping and waste
injection rates which satisfy water service demands while maintaining ground water
quality as mandated by current ground water standards. The analytical solutions to the
transport equations are nonlinear with respect to the decision variables (pumping and
injection rates). Therefore, the equations could not be embedded directly into a linear
programming model. In order to approximate this nonlinear response and to determine
optimal strategy, the management model was redefined to consist of two interdependent
linear submodels. Optimal pumping and injection strategies were first determined by a
flow submodel, then these results were input into a ground water quality submodel to
predict the maximum waste injection concentrations for each operational cycle.
The embedding technique was again applied to ground water quality management
18
by Shamir et al. (1984) for a coastal aquifer in Isreal. Aquifer response was simulated
by two sets of finitedifference equations which represent the flow of ground water and
the movement of the freshwater/saltwater interface. These equations were included as
constraints in the optimization model. However, in order to represent and constrain the
location of the freshwater/saltwater interface, a linearized model of saltwater intrusion
was incorporated to approximate the nonlinear response. Multiobjective linear
programming was incorporated to determine optimal annual operation of the coastal
aquifer.
Multiobjective optimization involved performing single objective linear
programming optimization and developing tradeoff functions for conflicting objectives
which are then utilized to obtain a final compromise solution. The management model
contained four objective functions which minimized 1) changes in water levels due to
pumping, 2) intrusion of the freshwater/saltwater interface, 3) chloride concentrations,
and 4) energy costs. Chloride mass balance equations and maximum upper concentration
limits were also specified for each aquifer cell location.
Gorelick et al. (1979) applied a numerical finitedifference approximation method to
solute transport convectivedispersive equations to determine chloride concentrations of
a transient pollutant source. Numerical approximation was utilized to generate a system
of linear equations which were then incorporated in a water quality optimization model
as part of the formulated constraint set. Water quality standards were specified in the
management model, and the model was executed to determine the maximum permissible
single source pollutant concentration. The method demonstrated the feasibility of
balancing water supply and waste disposal needs while satisfying water quality
requirements over long time frames.
As discussed, the embedding technique has been successfully applied to a wide
range of ground water management problems. The approach is capable of handling both
linear and nonlinear ground water hydraulics; steadystate and transient problems with
single or multipurpose objective functions. The primary advantage of the method is the
ability to emulate the physical processes of an aquifer system at a detailed level by
including ground water flow and/or solute transport equations for each grid cell directly
as constraints in the optimization model.
A disadvantage of the method is the dimensionality associated with large scale
management problems. The degree of dimensionality is a function of the number of
decision variables, aquifer grid cells, and time steps. Therefore, spatially large or multi
time step management models could potentially have hundreds or even thousands of
constraint equations. Attempts by Alley et al. (1976) and Yazicigil and Rasheedudd
(1987) to redefine multitime step embedded problems into a sequential set of smaller
problems have resulted in management decisions biased towards the shortterm. To
avoid these difficulties, the embedded technique should be limited to local aquifer
systems and steadystate management problems.
2.1.2 Unit Response Matrix Method
The "unit response matrix" method is utilized by simulation/optimization models
by incorporating the use of influence coefficients, superposition, and linear systems
theory. The simulation model is used initially to compute the system response to a unit
stimuli. For ground water management models, the response matrix approach is based
on the concept that the influence of discharging or recharging a single well on aquifer
20
drawdown at selected locations may be expressed as simple algebraic functions. The
individual influence functions are then combined using the principle of superposition to
obtain the aquifer response due to multiple wells. Therefore, a "response matrix" can
be developed which may be expressed algebraically as
m
dj = E [ aiJqi ] = [ aiJq + a(i1)Jq(i+) ".. amjqm
i=1
where
dj = drawdown at control point j
qi = pumpage at well i
aj = unit response function of well i on point j
i = 1,2,..., m
j = 1,2,..., n
m = number of wells
n = number of control points
The response matrix of influence coefficients, aj's, are generally determined from
analytical or numerical ground water simulation models. Since the response equations
need only be developed for selected points of interest, it is not necessary for equations
to be developed for each grid cell within the aquifer system. This characteristic
significantly reduces the dimensionality of the management problem when compared to
other simulation/optimization techniques.
One of the earliest pioneers to combine the response matrix approach with ground
water management, Deininger (1970) developed a management model which maximized
well field production by optimally selecting well locations and withdrawal rates. A
response matrix representing aquifer drawdown was developed from an analytical
solution to the Theis nonequilibrium equation and incorporated in the linear programming
model. Drawdown, pumping, and well characteristics were all specified as constraints
21
in the optimization model. Maximum drawdown limits were specified at well field
boundaries, and head versus discharge curves were approximated by piecewise
linearization to incorporate drawdown constraints at wells. Although head losses in well
casings were also linearized using a chordal approximation of the Manning formula, well
screen head losses were assumed to be constant. Deininger demonstrated the formulated
problem was easily solved by linear programming techniques.
Rosenwald and Green (1974) used the response matrix method in combination
with branch and bound mixedinteger linear programming to identify optimal well sites
in an "underground reservoir". Pressure coefficients, which represent potentiometric
drawdowns at wells due to changes in discharge rate, were determined by executing a
numerical simulation model. A set of constraint equations were then developed which
combined the pressure response matrix with maximum allowable drawdown limits. A
mixedinteger programming model was employed which assumed constant well flow rates
and involved a binary switching variable allowing each potential well location to be either
active or inactive.
One of the two cases used to demonstrate the technique involved the selection of
well locations in a hypothetical ground water aquifer. Enumeration of the alternatives
and comparison with simulation results verified that the response matrix approach based
upon linear superposition was applicable to ground water hydraulics. In the second case,
the authors applied the response matrix method to select well sites in a gas storage
reservoir. The response matrix initially developed erroneously estimated pressure heads
due to the nonlinear behavior of gas flow. Using various corrective techniques, attempts
were made to improve the pressure response equations. Although the results of the
22
optimization model were only slightly improved, the authors concluded that corrective
procedures offer some potential for improving the optimization of nonlinear systems with
response matrices.
Maddock (1972) discussed the use of transient algebraic technological functions
(response functions) for aquifers whose transmissivities vary with drawdown. These
algebraic functions were developed to relate seasonal pumping rates to drawdown levels
at specified wells, and the relationship was derived from an analytical solution of the
twodimensional transient ground water flow system. A hypothetical example is used to
demonstrate the method of combining response functions and quadratic programming to
determine optimal semiannual well withdrawal strategies.
Algebraic response functions were again developed by Maddock (1974) to
incorporate the optimization technique in nonlinear aquifer systems. Boussinesq's
equation for unsteady flow due to pumping in an unconfined aquifer was approximated
by an infinite power series. The total aquifer drawdown response function was then
expressed as the finite sum of an infinite power series in pumping values. The
approximation assumed fully penetrating wells, no vertical flow, and constant pumpage
over a single time horizon. The number of terms required for an accurate water
elevation estimate is dependent upon the ratio of drawdown to saturated thickness. To
demonstrate the methodology, the author formulated and executed the nonlinear
programming problem to determine least cost pumping distributions.
One of the first site specific applications of the response matrix technique was
presented by Heidari (1982) and involved determining optimum ground water allocation
strategies in Kansas. A twodimensional ground water simulation model developed for
23
the region was combined with the algebraic influence function proposed by Maddock
(1972,1974) to develop the response matrix. For simplicity, sixtyone hypothetical well
fields were created to represent total withdrawal from actual wells in the study area.
Once the response matrix was incorporated into the management model, constraint
equations were specified which limited drawdowns at each well to a fraction of the
saturated aquifer thickness.
The optimization model was formulated to maximize total pumpage under two
policy scenarios. With respect to the net appropriation (difference between appropriation
and recharge), maximum limits were specified for the first scenario but no limits were
specified for the second. The models were executed for five and ten year time periods
and each contained five time steps. Results revealed that removing the net appropriation
constraint created barely a noticeable increase in pumpage when the drawdown fraction
was allowed to be equal to or greater than 20%. The model also indicated that under
optimal conditions, only about 50% of the net appropriation could be satisfied over the
ten year time period.
Willis and Liu (1984) incorporated response functions to develop a biobjective
optimization model to allocate ground water to competing irrigation demands of the Yun
Lin basin in southwestern Taiwan. Objectives of the model were to maximize the sum
of hydraulic head and to minimize the total water deficit in the basin. The simulation
model, which was used to develop the response matrices, was developed using the
Galerkin finiteelement approximation method for ground water flow and was validated
using field data from over 350 monitoring wells. The response matrices were
incorporated as constraints and the model was initially executed as a steadystate linear
24
programming problem. Comparison of results from this initial steadystate optimization
to existing allocation strategies revealed that the total water deficit could be decreased
substantially without decreasing the sum total of hydraulic heads. By assigning weights
to the objective functions, solutions to the optimization problem were depicted in the
form of tradeoff curves to express the relationship between hydraulic head and total
water deficit.
The authors reformulated the optimization model to perform a transient analysis.
Piecewise linearization was used to incorporate time dependent boundary conditions for
the development of the response functions. Using the same objectives specified for the
steadystate condition, the model was executed and tradeoff curves were again developed
to depict the relationship between hydraulic head and water deficit. Results of the
transient analysis indicated that the steadystate formulation overestimated the reduction
in water deficit. However, the authors concluded that a significant water savings could
still be achieved with the use of the optimization model.
Danskin and Gorelick (1985) were some of the first to implement response
functions in the management of multilayer aquifer systems. The study area consisted
of a multilayer aquifer system connected to a surface water system, and the model was
formulated to evaluate the efficiency of flow between the two systems. Critical factors
controlling basin management decisions were also identified through execution of the
mixedinteger linear programming model. The response matrix was developed using a
quasithreedimensional finitedifference simulation model, and the influence coefficients
for the upper unconfined aquifer was linearized using an iterative approach. The
response functions were used to develop water elevation constraints and additional
25
limitations were specified with respect to surface water flow, vertical leakance, and water
service demands. Mass balance constraints were also incorporated to regulate the
concentration of dissolved solids.
Evaluation of historical basin management practices revealed that "the cost of
operation during the study period was twice that of optimal basin management."
However, the largest inefficiencies were localized and most activities were within 20%
of optimal. Sensitivity analyses indicated that the controlling factor in basin operations
involved surface water and ground water relationships.
Willis and Finney (1985) developed response equations using finitedifference
methods, quasilinearization, and matrix calculations. Nonlinear optimization for an
unconfined aquifer system was performed by incorporating a quasilinearization
optimization algorithm and projected Langrangean methods. The model was structured
as a discrete optimal control problem and determined the optimal pumping pattern while
satisfying water demands. Quasilinearization and MINOS (Modular InCore Nonlinear
Optimization System) algorithms were both found to be efficient for solving moderately
sized nonlinear ground water management models. However, Willis and Finney
concluded that large management problems could be solved more efficiently by applying
quasilinearization optimization.
Herrling and Heckele (1986) used both the embedding and response matrix
techniques to couple a finiteelement simulation model with an linear simplex
optimization model. Management of the nonlinear ground water system was performed
by optimizing well locations along with pumping and infiltration rates while satisfying
ecological constraints. These constraints consisted of sustaining ground water levels and
26
meeting contamination standards. The advantages of both coupling techniques were
incorporated. The embedding method involved the implicit consideration of the flow
model, and the influence function method reduced the amount of computer storage
required.
The response matrix technique was combined with a stochastic approach by Tung
(1986) for the management of ground water resources. The "chance constrained"
management model was formulated using linear systems theory to determine optimal
withdrawal rates in a well field subject to specified reliability constraints. Transient,
nonleaky drawdown response functions were developed from an analytical solution to the
CooperJacob equation while transmissivity and storativity values were treated as
independent random variables. Firstorder analysis was employed to estimate the
statistical characteristics of drawdown at each control point. Stochastic drawdown
response functions implicitly incorporating parameter uncertainty were then developed
as linear constraints. Drawdown limitations were stated at each control point which
required drawdowns to be less than a specified value times a reliability factor.
The management model was applied to a hypothetical confined aquifer where
optimal withdrawal rates were determined for three potential wells over three time
periods. While maximizing total well production, a sensitivity analysis was performed
on the model to investigate the effect of parameter uncertainty and of varying reliability
factors. Tung discovered that model results were insensitive to changes in storativity but
quite sensitive to changes in transmissivity and the specified reliability factor. A
postoptimal analysis indicated that the linear programming model produced acceptable
results in terms of complying with reliability requirements, but only when the
transmissivity uncertainty was small.
Lindner et al. (1988) applied the response matrix technique to a twoaquifer
system to determine the optimum ground water withdrawal while meeting established
environmental criteria. The twoaquifer system consisted of an upper unconfined and a
lower confined aquifer with a leaky intermediate layer. The environmental criteria
included meeting specified ecological conditions and ground water levels in the upper
aquifer. Though ground water was withdrawn from the confined aquifer only, the
unconfined aquifer was also affected due to leakage from the separating layer.
Galeati and Gambolati (1988) employed the response matrix technique to solve
a threedimensional aquifer dewatering problem. The water management model was
formulated to identify the optimal spatial distribution and corresponding well rates
required to maintain desired water levels during two planning periods. A three
dimensional finiteelement ground water flow model was developed and repeatedly
executed to develop a steadystate response for each abstraction well. The individual
influence coefficients were then combined using the principle of superposition, and
hydraulic head limits were specified to create a set of water elevation constraints. The
model was formulated to minimize total withdrawal while still dewatering the study area.
The authors assumed linear response since the dewatering area responds as a
hydrodynamically closed system having essentially no influence on the surrounding
unconfined aquifer.
Once the management model was executed and an optimal dewatering strategy
was prescribed, it was discovered that some wells which were deemed active during the
first planning period were deemed inactive during the second period, and vice versa.
28
Therefore, in order to avoid large installation costs and withdrawal rate nonuniformity
associated with the solution, the model was reformulated as a mixedinteger linear
programming model. Mixedinteger linear programming involved designating wells with
a withdrawal rate of either zero or some constant value. Additional constraints were
added which required wells be utilized in both time periods if they are utilized at all.
Although the solution to this reformulated problem reduced the number of wells by 36%,
the total pumping rate increased slightly. After performing yet another reformulation and
execution with respect to minimizing costs, the authors concluded that the intermediate
mixedinteger programming solution represented the best compromise between
withdrawal and installation costs.
Similar research was performed by Lall and Santini (1989) to extend the response
matrix approach to a nonlinear, multilayered, unconfined aquifer system. The Grinski
potential concept was utilized to develop a linear approximation of the nonlinear aquifer
system and was described as being analogous to the velocity potential for a confined
aquifer system. For this case, hydraulic head is directly related to the vertical depth and
hydraulic conductivity of each layer. The authors demonstrated that the steadystate
continuity equation is a linear function of Grinski potential when hydraulic heads were
converted to Grinski potentials. By incorporating the method of finite elements, the
aquifer responses were then modeled as linear functions of the Grinski potential. The
superposition principle was employed to sum the individual responses and obtain a
response matrix in terms of the Grinski potential. The use of the Grinski potential for
linear approximation was also shown to be applicable to transient ground water
hydraulics under certain conditions. The method was also demonstrated using three
different variations of the original dewatering problem.
Response functions and mixed integer programming were applied by Chau (1989) to
analyze pressure relief systems. Optimal well sites and discharge schedules were
determined while discharge was minimized and hydraulic heads were maintained with
respect to soil stability and ground water flow. Response functions were determined
from a simulation model and represented the drawdown as induced by discharge from
another well and time period. The effects of varying performance parameters which
represent aquifer characteristics, hydraulic head limits, well capacities, and discharge
elevation at well locations were evaluated. Two hypothetical examples were analyzed
to demonstrate the tradeoffs between system parameters and system performance.
Yazicigil (1990) incorporated the response function approach to determine optimal
planning and operating policies of a multiaquifer system in Eastern Saudi Arabia.
MODFLOW, a threedimensional finitedifference ground water flow model developed
by McDonald and Harbaugh (1988), was executed under transient conditions and used
to generate response coefficients for each well field. Since the flow model was
formulated to simulate an eight year planning period with monthly time steps, over
20,000 individual influence coefficients were developed. These response functions were
then combined with both linear and quadratic programming to determine optimal water
management strategies for the 52 well fields in the basin. Three different objective
functions were formulated which maximized withdrawals, minimized drawdowns, and
minimized pumping costs. Tradeoff costs which related withdrawal to drawdown,
aquifer dewatering, and pumping costs were also determined from results of the
optimization model.
The response matrix technique has not only been applied to management models
30
with respect to ground water quantity, but also ground water quality. Colarullo et al.
(1984) demonstrated this approach using a hypothetical unconfined aquifer which had
historically been used for surface waste disposal, but was soon to be developed as a
freshwater supply. A twodimensional ground water flow model was incorporated to
develop linear response equations which represented the aquifer response due to a change
in pumpage rate. Using this same method, response equations were also defined for
pumpage induced velocities and incorporated in the model as constraints to limit
contaminant flow. The optimization model was formulated to determine what quantity
of water could be removed to supplement water service demands and how interception
wells should be operated to avoid contamination of freshwater supplies. A nonlinear
optimization algorithm was utilized to identify optimal well discharges for supply and
interception wells. The authors verified previous assumptions of linearity for the
response matrices by executing the simulation model with the prescribed strategy.
Influence coefficients derived from an approximation of a solute transport model
were applied by Datta and Peralta (1985) to revise a quantity optimization model to
include quality. The authors applied the approach of Peralta and Killian (1985) to
optimize the potentiometric surface and identify the water use strategy required to
maintain the surface. Steadystate hydraulic stresses were determined from the
simulation model, and steadystate ground water concentrations were determined from
the solute transport model. The influence coefficients were determined based on the
hydraulic head levels required to meet quality limits, and then used in establishing new
hydraulic head constraints. The modified optimization model was created through the
use of constrained derivatives of the quadratic optimization model. The methodology
31
was found useful in determining concentrations in a subsystem of a larger model and/or
the influence of water quantity changes on water quality.
Similar to the effort presented by Tung (1986), Wagner and Gorelick (1987)
combined the response matrix approach with a stochastic optimization model for ground
water quality management. A twodimensional finite element flow and solute transport
model was utilized to develop the head and concentration response, respectively, as a
function of pumping and aquifer parameters. Simulation and multiple regression were
used to develop parameter estimates and minimize the differences between simulated and
observed values. In order to incorporate parameter uncertainty, the response equations
were transformed from deterministic constraints to probabilistic constraints through first
order, first and secondmoment analysis. The chance constrained nonlinear optimization
model was formulated to identify well locations and withdrawal rates for aquifer
remediation under specified parameter uncertainty.
The authors applied the management model to both steadystate and transient
conditions to determine which withdrawal strategy would satisfy water quality standards.
Results revealed that the prescribed location, number, and pumping rate of wells were
very sensitive to reliability level. Monte Carlo simulations verified that the true mean
concentration and the concentration predicted by firstorder analyses were nearly identical
and normally distributed.
Finney et al. (1992) incorporated response equations into a management model
for the control of saltwater intrusion in a multiaquifer system in Jakarta, Indonesia.
Hydraulic response equations were also used by Willis and Finney (1988) in the
development of a simulation/optimization model for the management of seawater
32
intrusion in YunLin ground water basin of southwestern Taiwan. These efforts along
with others are discussed in the following literature review of saltwater intrusion models.
2.2 SALTWATER INTRUSION MODELS
In many coastal areas, excessive abstraction of ground water has resulted in a
decline of freshwater potentiometric heads and therefore an intrusion of salt water.
When this encroachment may have adverse effects on water supply quality,
simulation/optimization models are applied to these coastal areas in order to minimize
the effect of saltwater intrusion on freshwater supplies. Simulation models which
predict the location of the freshwater/saltwater interface are combined with optimization
models to minimize encroachment of the interface while meeting pumping demands.
Some of the initial research involving the optimum exploitations of ground water
reserves located in coastal aquifers was performed by Cummings (1971) and Youngs
(1971). Cummings created a ground water management optimization model which dealt
with the economic impacts of aquifer exploitation and saltwater intrusion with respect
to time. The model focused on the interrelated problems between annual ground water
pumping rate and the annual cost of pumping. Although Cummings mentions how
hydrological simulation models and optimization models complement each other, most
of his research involved the economic side of ground water mining.
Youngs' (1971) research involved determining the optimum well conditions in
order to maximize freshwater output from coastal aquifers. By performing an analysis
of the horizontal seepage, Youngs determined that the maximum pumping rate of fresh
water could be calculated precisely and that the optimum conditions for well installation
33
can be found which will produce the maximum continuous pumping rate. Knowing that
pore water pressure is reduced near an operating well which causes an upcoming of saline
water, the analysis identifies the pumping rate required to raise the interface to the
bottom of the well. Youngs then determined that maximum withdrawal is obtained when
wells are positioned to a elevation equal to sea level and are pumped to maximum
capacity. The analysis is presented using two hypothetical examples.
The dynamic programming approach was applied by Nutbrown et al. (1975) in
a digital simulation/optimization model in order to effectively and efficiently manage a
coastal aquifer near Brighton, England. The goal of their research was to predict the
maximum yield of the aquifer while considering the limiting factor of saltwater
intrusion. The simulation model was created to describe the transient effects of natural
infiltration and abstraction of both fresh and saline ground water. The dynamic
programming method was combined with the simulation model to generate abstraction
regimes which would maximize the final storage of fresh water in the aquifer.
The approach was based on the methodology that at regular intervals of time, the
spatial abstraction distribution for succeeding intervals could be calculated on the basis
of existing water levels. The concept that pumping rate was proportional to the
magnitude of local ground water flow was assumed. The model was then executed under
average infiltration conditions until cyclic equilibrium was attained. Once validity had
been determined, the model was applied to various drought and recharge conditions to
determine the various optimum scenarios.
As mentioned earlier in the Embedding Method section of this literature review,
Shamir et al. (1984) developed an annual operating plan for a coastal aquifer in Israel.
34
A linearized model of saltwater intrusion was combined with a flow simulation model
to represent the movement of the saltwater interface. This simulation model was
incorporated into a multiobjective linear programming model to produce the management
model.
A planning model was developed by Willis and Finney (1988) to control salt
water intrusion and declining ground water levels in the Yun Lin ground water basin of
southwestern Taiwan. The simulation/optimization model applied hydraulic response
equations to relate the location and magnitude of ground water pumping and recharge to
the movement of the saltwater interface. Finitedifference methods were used in the
simulation model to approximate the aquifer's response to various management strategies.
The model was based on the following hydraulic assumptions: 1) Hydrodynamic
dispersion is negligible so the sharp interface theory is valid; 2) The Dupuit
approximation is valid throughout the aquifer system; 3) The aquifer base is
impermeable; 4) The hydraulic conductivity and storage coefficients are invariant with
depth; and 5) Leakage to or from the aquifer is at steadystate.
The solution algorithms were based on the influencecoefficient method combined
with quadratic programming and also the reducedgradient method combined with a
quasiNewton algorithm. The influencecoefficient algorithm is based on hydraulic
response equations of the ground water system, and the coefficients are applied to the
quadratic programming optimization problem of the management model. The quadratic
program determines the optimal direction vector which is then used to revise the current
solution using a gradientbased algorithm. The reducedgradient/quasiNewton algorithm
can also be used to determine optimal solutions of the planning model. This algorithm
35
is implemented using MINOS and is able to determine a feasible decent direction vector
at any iteration. The simulation/optimization model was created using historical ground
water data for the basin in order to develop optimal pumping and recharge schedules.
The authors made the following conclusions based on the results of the
optimization analysis: 1) Both algorithms produced stable and reliable solutions to the
saltwater management problem; 2) Greatly differing pumping and recharge schedules
produced essentially the same objective function values for the applied basin; 3) The
management problem is characterized by local optimality problems; and 4) Different
starting solutions for the algorithms produce different optimal schedules.
Finney et al. (1992) incorporated response equations into a management model
for the control of saltwater intrusion in a multiaquifer system in Jakarta, Indonesia.
Response functions relating pumpage to the location and movement of the saltwater
interface were developed from a finitedifference simulation model. Within the model,
multiple aquifer systems were linked through their recharge terms. The response
equations were then included as part of a nonlinear optimization model which sought to
minimize the total squared volume of salt water in each aquifer. Initial attempts to reach
a solution with the MINOS programming package resulted in solutions that were not truly
local optimum values. To correct the problem, Box's sequential search algorithm was
used. This resulted in a 20% improvement over the MINOS generated values. Finney
et al. concluded that the combined use of simulation and optimization was able to reduce
the magnitude of saltwater intrusion in the Jakarta basin by 6%.
2.3 INTERFACE PROGRAMS
Due to the increasing popularity of the U.S. Geological Survey Modular Three
Dimensional FiniteDifference Ground water Flow Model (MODFLOW) and of ground
water optimization models, recent efforts have been made by various groups to combine
these models via interface computer programs. AQMAN3D and MODMAN are two
such programs which simplify the interaction and iteration process between MODFLOW
and particular standard optimization programs (GAMS, MINOS, etc.).
Puig et al. (1992) developed AQMAN3D by modifying the twodimensional code
AQMAN originally developed by Lefkoff and Gorelick (1987). The revised version is
a mathematical programming system data set generator for aquifer management using
MODFLOW as its ground water flow simulation subroutine. The FORTRAN77
computer code can be used to formulate a variety of aquifer management models;
depending on the chosen objective function to be optimized and the constraints imposed
on hydraulic conditions. The program creates input files to be used by any standard
optimization program using Mathematical Programming System (MPS) input format.
When AQMAN3D is used with an optimization program, the optimum pumping
and/or recharge strategy can be determined while ground water hydraulic conditions are
maintained within specified limits. The applied management function may be linear or
nonlinear and restrictions can be applied to ground water heads, gradients, and/or
velocities. The program is limited to confined or quasiconfined aquifer systems and to
nonlinear sources/sinks in the ground water flow model.
The AQMAN3D aquifer management modeling process is initiated by
qualitatively and quantitatively conceptualizing the aquifer system by using MODFLOW.
37
The simulation flow model is then calibrated to steadystate or transient conditions and
then executed to simulate particular unmanaged aquifer conditions. The optimization
portion of the model is then initiated by developing a objective function and a system of
constraints which reflect the aquifer management condition. The AQMAN3D program
is now executed to interact with MODFLOW to produce the MPS input file containing
the objective function, constraints, and response coefficient matrix. The MPS file is
applied to the standard optimization model to determine the optimal scenario for the
specified management conditions. This optimal scenario should then be applied back into
the MODFLOW model to observe and verify the response of the aquifer system.
Another model which has been developed to perform the interface operations
between MODFLOW and a standard optimization program is recognized as MODMAN
(MODFLOW Management An Optimization Module for MODFLOW). GeoTrans, Inc.
(1990) developed MODMAN for use by the Southwest Florida Water Management
District to assist in determining optimum pumping scenarios in their region. MODMAN
is analogous to AQMAN3D in that it is a modified extension of AQMAN (Lefkoff and
Gorelick 1987) and accommodates threedimensional problems. However, unlike
AQMAN3D, MODMAN is actually a linked module or subroutine of a revised version
of MODFLOW and can be executed in two modes. The MODFLOW code was modified
by the creators of MODMAN in order to facilitate the process of determining the
response coefficient matrix.
When used together with standard optimization software, MODMAN determines
the optimum strategy depicting the location and rate of extraction and/or injection wells
while satisfying userspecified constraints. One of many objective functions may be
38
maximized or minimized depending on the goals of the ground water manager. Once the
optimal pumping/injection strategy has been determined, MODMAN is incorporated
again by acting as a postprocessor of the optimization output data. MODMAN
automatically inserts the optimal well rates into the MODFLOW input file, executes a
simulation based on the optimal well strategy, indicates which constraints are binding in
the optimization model, and warns the user if nonlinearities have significantly affected
the optimization process.
2.4 SUMMARY
The approach of combining simulation and optimization models was first used in
the development of management models in the late 1960's and early 1970's. This type
of an approach offers the water resource manager a means of considering physical,
policy, and environmental factors simultaneously. Generally, the combined approaches
can be divided into two separate categories according to the method in which ground
water flow equations are incorporated into mathematical optimization models.
In the "embedding" technique, numerical or analytical solutions to ground water
flow and/or transport equations are written for each control node and are subsequently
included as constraints in the optimization model. This method has been successfully
applied to a large number of ground water management problems and has been
demonstrated to handle both steadystate and transient conditions. This approach has also
been used with problems involving linear and nonlinear ground water hydraulics and
single and multipurpose objective functions. Although the ability to accurately simulate
physical processes of an aquifer system is a major advantage, there is a dimensionality
39
problem associated with the embedding approach when applied to large scale management
problems. Therefore, this technique is generally only applicable to small aquifer systems
and steadystate management problems.
The "response matrix" approach has also been demonstrated as a technique to
combine simulation and optimization models for ground water management. As opposed
to the embedding technique however, this method is capable of handling large scale
management problems. The response matrix is composed of coefficients which represent
the influence of specified wells on the aquifer system. The aquifer response in terms of
change in hydraulic head, or drawdown, is usually derived from a simulation model.
The response is then expressed as simple algebraic functions in the mathematical model
using the principle of superposition. Constraints are incorporated into the algebraic
functions to limit drawdown at specified control points. Because equations are developed
only for the selected control points, use of the response matrix method generally results
in a reduced problem.
Additional ground water optimization methods including differential dynamic
programming, multiobjective programming, and quadratic programming were reviewed.
Techniques such as combining embedding and response matrix approaches for a
simulation/optimization model and revising an existing quantity optimization model to
incorporate quality were also discussed.
An assortment of optimization models were reviewed which incorporated varying
methods of managing coastal aquifers while minimizing the adverse effects of saltwater
intrusion. Many of the models were combined with simulation models which applied
finitedifference methods to determine the aquifer's response to varying pumping
40
strategies. All of the models assumed that the effect of hydraulic dispersion between
fresh and salt water was negligible so that the "sharp interface" concept was acceptable.
In order to simplify the interaction process between the ground water simulation
model and the numerical programming software, interface programs have been developed
which link the ground water flow model MODFLOW with a standard optimization
program. Two such models which were reviewed are AQMAN3D and MODMAN.
Both models are extended modifications of AQMAN (a previously developed two
dimensional interface model), incorporate the response matrix approach, and
accommodate threedimensional problems. MODMAN differs from AQMAN3D in that
it is actually a module of MODFLOW, it can be applied in a second mode as a
postprocessor, and it is more "userfriendly".
Existing simulation/optimization models pertaining to ground water resource
allocation and management are displayed in tabular form in Table 2.1. The model name,
authorss, type, aquifer condition, optimization technique, objective function, and
constraints are shown for easy comparison.
Table 2.1 Ground water optimization models.
Author Type Equation / Solution Objective Function Constraints Optimization Model
[Year] [Dimension] Aquifer Condition Technique [Liner / Nonlinear] [Linear / Nonlinear] Technique Name
Ahlfeld Quantity / Saturated Confined Linear / Finite Minimize Total Pumping Rates, Linear VCON
[1988] Quality Steady / Unsteady Element Pumping Hydraulic Head, Programming
[2D] [Linear] Magnitude & Direction
of Groundwater
Velocity [Linear]
Ahlfeld & Quantity / Saturated Confined Linear / Finite Minimize Total Minimum Nonlinear GW2SEN
Pinder Quality Steady / Unsteady Element Volume of Pumping Concentrations at All Programming
[1988] [1 & 2D] [Linear] Nodes at Future Times
[Nonlinear]
Alley Quantity Unsteady Confined Linear / Finite Maximize Well Flow Rates, Linear AQMG
et al. [2D] Saturated Difference Hydraulic Heads Hydraulic Heads Programming
[1976] [Linear] [Linear]
Aly & Quantity Steady / Unsteady Theis Analytical Minimize Extraction, Hydraulic Gradients, Linear, Quadratic, or US/WELL
Peralta [2D] Confined / [Well Function] Injection Extraction or Injection Nonlinear
Unconfined Deterministic / [Linear / Nonlinear] Rates, Hydraulic Head Programming
Stochastic [Linear / Nonlinear]
Aquado & Quantity Steady Unconfined Linearized / Finite Minimize Total Hydraulic Heads, Linear OPAQD
Remson [2D] Saturated Difference Pumping Pumping Rates Programming
[1974] [Linear] [Linear]
Deninger Quantity Unsteady,Unconfined Analytical Maximize Well Drawdown, Linear SAWSS
[1970] Radial Saturated [Well Formula] Production Well Facility Programming
[Linear] [Linear]
Elango & Quantity Steady Confined Linear / Finite Maximize Total Flow Equilibrium, Linear
Rouve [2D] Saturated Element Pumpages Piezometric Heads, Programming
[1980] [Linear] Pumping Capacities
[Linear]
Table 2.1 continued
Author Type Equation / Solution Objective Function Constraints Optimization Model
[Year] [Dimension] Aquifer Condition Technique [Liner / Nonlinear] [Linear / Nonlinear] Technique Name
Gorelick Quality Steady State Linear / Finite Maximize Waste Water Quality Linear OLMWDF
& Remson [2D] Flow and Transport Difference Disposal [Linear] Programming
[1982] Confined Saturated [Linear]
Haimes & Quantity Aquifer / Stream Linearized / Maximize User's Water Requirements, Quadratic MGWSW
Dretzen [2D] System, Unsteady Cell Model Net Benefit Lift & Pumping limits, Programming
[1977] Unconfined [Quadratic] Capacity of Recharge
Saturated Facility [Linear]
Heidari Quantity Unsteady,Unconfined Linearized / Finite Maximize Pumping Pumping Demands, Linear LSTLP
[1982] [2D] Saturated Difference Rates over Time Drawdown Programming
[Linear] [Linear]
Koltermann Quantity Unsteady Confined Linear / Finite Maximize Hydraulic Groundwater Flow, Linear
[1983] [2D] Saturated Difference Head and Minimize Piezometric Heads, Programming
Water Transfer and Pumping Capacities
Recharge [Linear] [Linear]
Larson Quantity Steady Unconfined Linearized / Finite Maximize Steady Pumping Rates, Mixed Integer OTGWD
et al. [2D] Difference State Pumping Number of Wells, Linear
[1977] [Linear] Drawdown [Linear] Programming
Molz & Quantity Steady Confined Linear / Finite Minimize Total Hydraulic Heads, Linear HGCAQ
Bell [2D] Saturated Difference Pumping Head Gradients Programming
[1977] [Linear] [Linear]
Peralta Quantity / Confined Saturated Linearized Boussinesq, Minimize Changes in Steady/Unsteady Flow Linear Goal MODCON
et al. Quality Steady / Unsteady Linearized & Nonlinear Piezometric Head and Transport Water Programming
[1987] [2D] Transport/ Method for [Linear] Quality [Linear]
Characteristics
Peralta Quantity Steady / Unsteady Linear / Finite Minimize Withdrawal Head, Head Gradient, Linear US/REMAX
et al. [3D] Confined / Unconfined Difference and Recharge Rates Flow Velocities, Programming
[1992] Saturated/Unsaturated [Linear] Demand, Capacity &
Pumping [Linear]
Table 2.1 continued
Author Type Equation / Solution Objective Function Constraints Optimization Model
[Year] [Dimension] Aquifer Condition Technique [Liner / Nonlinear] [Linear / Nonlinear] Technique Name
Rosenwald Quantity UnsteadyUnconfined Linearized / Finite Minimize Total Production Demand, Mixed Integer OLWR
& Green [2D] Saturated Difference Pumping Number of Wells Linear
[1974] [Linear] [Linear] Programming
Shamir Quantity UnsteadyUnconfined Linearized / Finite Minimize Energy Pumping Demands, Linear OPOCAQ
et al. [2D] Saturated Difference Demands for Import/Export Fluxes, Programming
[1983] Pumping / Recharge Drawdown Position,
[Linear] Water Quality
Limits [Linear]
Takahashi Quantity Confined/Unconfined Linear / Finite Maximize Extraction Potentiomentric Head, Linear USUGWM
& Peralta [Quasi 3D] Steady Difference [Linear] Pumping Rate Programming
[1991] [Linear]
Willis Quantity / Confined Saturated Linear / Finite Maximize Lowest of Water Target, Linear PMMGWQ
[1979] Quality Unsteady Transport Element Waste Concentration Pumping & Injection, Programming
[2D] QuasiSteady Flow [Linear] Water Quality [Linear]
Willis Quantity Unsteady Confined/ Boussinesq Equation/ Maximize Sum of Agricultural Demand, Linear UARGWM
[1983] [2D] Unconfined Finite Element Heads, Minimze deficit Heads, Well Capacity Programming
Analytical [Linear] [Linear]
Yazicigil Quantity Unsteady Confined Linear / Finite Minimize Total Water Demands, Linear OPMRA
et al. [2D] Saturated Difference Pumping Maximum & Minimum Programming
[1987] [Linear] Pumping Rates,
Drawdowns [Linear]
Yazicigill & Quantity Saturated Confined Linear / Finite Maximize Total Water Demands, Linear OPMGW
Rashee [3D] Steady / Unsteady Difference Hydraulic Head Hydraulic head bounds, Programming
duddin [Linear] Maximize Pump Rates
[1987] [Linear]
Reference: ElKadi et al. (1991)
CHAPTER 3
3.0 GROUND WATER SIMULATION MODELS
From the literature review, it was determined that the most efficient general
formulation of the allocation models would be one containing aquifer response constraints
developed from steadystate unit response matrices created with ground water flow and
transport simulation models. Two and threedimensional finitedifference ground water
simulation models have been employed extensively for the last twenty years by the
United States Geological Survey and other ground water professionals. These numerical
simulation models apply particular boundary conditions to a spatially discretized aquifer
system in order to predict potentiometric heads, fluxes, and solute concentrations
throughout a specified region. These boundary conditions include the size, shape,
hydrogeologic framework, hydraulic characteristics, and particular known fluxes affecting
the aquifer system. Although there are various ground water flow simulation programs
available, MODFLOW was used in this project, because a simulation model had already
been developed for the project area. Similarly the solute transport model program
DSTRAM was incorporated into the project because it had been successfully implemented
in the project area.
MODFLOW is a modular threedimensional finitedifference ground water flow
program created by the United States Geological Survey (McDonald and Harbaugh
1988). The model is based on the following well known governing equations describing
44
the movement of an incompressible fluid through porous material:
a a h a h a ah ah
(K ) + (K ) + (K ) W= S (3.1)
x ax aY ay az zzaz sat
where K,, Ky, and K,, are values of hydraulic conductivity along the x, y, and z
coordinate axes, which are assumed to be parallel to the primary axes of hydraulic
conductivity (Lt'). H is the potentiometric head (L), W is a volumetric flux per unit
volume and represents sources and/or sinks of water (t'), S, is the specific storage of the
porous material (L'), and t is time (t).
Equation 3.1 approximates flow under nonequilibrium conditions in a
heterogeneous and anisotropic medium, assuming the principal axes of hydraulic
conductivity are aligned with the coordinate directions. Since analytical solutions to this
equation are rarely possible, the finitedifference numerical method is implemented to
procure approximate solutions. In this method, the continuous system described by the
partialdifferential equation is replaced by a finite set of discretized equations in time and
space. Continuous partial derivatives are replaced by discrete algebraic functions
describing the change in pressure head at the distinct points. A system of linear
algebraic difference equations results from this methodology, and its solution produces
values of pressure head at specific points and time.
MODFLOW is able to simulate aquifer systems with layers that are confined,
unconfined, or a combination of both. Each cell of every layer is specified as being
inactive (no flow), active (variable head), or having a constant head. Boundary
conditions such as specific flux, specific head, or a headdependant flux can be applied
46
in the model. The model is capable of simulating external flows such as discharge and
injection wells, rivers and streams, evapotranspiration, precipitation, and agricultural
drains. Because of its diversity and effectiveness, MODFLOW has become one of the
most popular ground water flow simulation models amongst hydrogeologic professionals.
DSTRAM developed by HydroGeologic Inc. is a threedimensional numerical
finite element program that simulates fluid flow and solute transport saturated porous
media. The code is capable of performing several types of analysis. These include
ground water flow analyses, trace concentration solute transport analyses, and density
dependent coupled flow and transport analyses. The code is based on the following
governing equation for threedimensional densitydependent flow and transport in an
aquifer system:
a I p gej) p i,j = 1,2,3 (3.2)
axi aPxi at
where p is fluid pressure, klj is the intrinsic permeability tensor, p is the fluid density,
1A is the dynamic viscosity, g is the gravitational acceleration, ej is the unit vector in the
upward vertical direction, and 4 is the porosity of the porous medium.
DSTRAM analysis can be performed in an areal plane, a vertical crosssection,
an axisymmetric configuration, or a fully threedimensional mode. Because of its special
design features, DSTRAM is capable of handling a wide range of complex three
dimensional, steadystate or transient, field problems and producing values of solute
concentration at specific points and times (Huyakorn and Panday, 1991).
3.1 REGIONAL FLOW SIMULATION MODEL
A regional model covering the Volusia County study area was acquired from
SJRWMD. The model was originally developed by Geraghty & Miller, Inc. (1991), but
several modifications were incorporated by SJRWMD to increase simulation accuracy
(Williams, 1993). The model as received is capable of simulating three different steady
state stress conditions with respect to time, redevelopment, year 1988, and year 2010.
The finitedifference grid of the regional model area consists of 86 columns, 91
rows, and 5 layers for a total of 39,130 cells. The model covers an area of
approximately 1,850 square miles and consists of planar cell spacings varying from 0.25
to 2.0 miles. The model is bounded by the St. John's river on the west and extends
approximately seven miles off the Atlantic Coast on the east.
The 5 layers of the model aquifer system consist of the surficial aquifer in layer
1, the Upper Floridan aquifer in layer 2, the middle semiconfining unit of the Floridan
aquifer in layer 3, and the Lower Floridan aquifer in layers 4 and 5. The upper
confining unit which separates the surficial and Upper Floridan aquifers are
interconnected through the use of leakance coefficients which represent the hydraulic
connection between the two layers. The surficial aquifer is modeled as unconfined and
the Floridan aquifer system is modeled as confined.
The boundary conditions used in the development of the regional model included
constant flux, constant head, and headdependant flux boundaries. Constant head
boundary conditions were applied to regions where impacts of pumping were assumed
to be negligible. These areas include the Atlantic Ocean, Halifax River, St. Johns River,
48
Blue Springs vicinity, and miscellaneous lakes within the county (See Figure 3.1).
Constant heads were also applied along the edge of the model near Blue Springs in the
surficial aquifer to allow flow from Lake County. Constant flux boundary conditions
were used in the model to represent the various discharge and recharge areas. The
discharge flow rates for each cell were calculated by summing the pumping rates from
all wells located within that cell for a specific layer.
Headdependant flux boundaries are often known as mixedtype boundary
conditions since they are essentially a combination of specified head and constant flux
conditions. They are applied to represent an unknown flux which is dependant on a
specified, hydraulically connected, and externallylocated head. Blue Spring, Ponce De
Leon Springs, Gemini Springs, and a variety of small creeks were modeled using head
dependant flux boundaries. The western and northern boundaries of the Floridan Aquifer
system were defined with a type of headdependant boundary called general head
boundaries in MODFLOW. These conditions provide adequate inflow to the model
while minimizing boundary effects due to varying pumping rates.
3.2 FLOW SIMULATION MODEL FOR PROJECT AREA
The five optimization models produced through the efforts of this project were
applied to the region defined by SJRWMD (See Figure 3.2). This project area covers
a smaller subarea of the above described regional flow simulation model. To facilitate
the development of this model, the Volusia County regional flow simulation model was
modified to create a MODFLOW subregional flow simulation model for the project area.
I RvA.rCa
EXPLANATION
= SURFACE WATER
***..** MAJOR DRAINAGE
  MIN DRAINAGE
BODIES
DMDE
nn fn
0
\
Jk
(from Rutledge, 1985o. fig 22)
09 SPRING
0 5 10 KILOMETERS
Figure 3.1
Map of Volusia County, Florida.
Fagler
Putnam
Lake
REGIONAL MODEL
FHI
PROJECT MODEL
volusia
0 5 miles
o0 5 kilomes I
Figure 3.2
Volusia County regional and project
numerical model grid discretizations.
51
The process of developing the project flow simulation model was initiated by
essentially cropping the northern, western, and southern portions of the regional model.
(The eastern border coincides with the larger regional model). This process involved
modifying all the regional model input files so they would correspond with the project
geometry. Computer programs written in FORTRAN code were developed to facilitate
this process. Removed from the regional model were 19 columns at the western border
and 20 and 11 rows from the northern and southern borders respectively. The location
of the project model with respect to the regional model is depicted in Figure 3.2. Figure
3.3 provides a reference for flow model cell locations.
The next step in the modification process involved applying the appropriate
boundary conditions to the project model in order to simulate the regional model as
closely as possible. Boundary conditions were applied using the general head boundary
option on the three boundaries of the model which do not coincide with regional
boundaries. This was performed by applying a constant external head value and a value
which represents the conductance between the external location and the boundary.
Effective conductance values between the regional and project borders were
calculated for each layer. Initially, the conductance of each cell located outside the
model was calculated in the direction perpendicular to the border. The effective
conductances were then calculated between the regional and project borders using the
individual cell conductances, or
SAmOH
Figure 3.3
Project Flow Model
Cell Location Reference Map
1 L[TITITl
0
..... :.. .
^ ^ ^ ^ ^ ^ ^ll IN 1111 1 11 1
1 1 1 1. (3.2)
C C1 C2 C.
where C is conductance, K is hydraulic conductivity, A is cross sectional area
perpendicular to flow, L is length of flow path, and n is the number of cells.
FORTRAN programs were developed to facilitate the process of gathering and
manipulating the geometric and hydraulic data for these calculations. Results from the
execution of the regional model were also used to determine the conductance and external
head values for the general head boundary conditions. These boundary conditions varied
depending on the simulated condition of the regional model redevelopmentn, year 1988,
or year 2010). For layers 2 through 5 of the regional model, general head boundary
conditions were applied to the northern and western borders. Therefore, external head
values used for the regional model were applied to the project model, and conductance
values applied to the project model were determined by combining the conductance
values used for the regional model with the conductance between the borders of the two
models. To determine the external head values to be applied at the southern border of
the project model, the regional model was executed for specific conditions to determine
head values along its southern border. This process assures that the influences of the
regional models on the smaller project model is incorporated in the simulation.
For the surficial aquifer (layer 1), there are often many bodies of water between
the regional and project borders which are modeled as constant heads (especially on the
western border). In these cases, the constant head values which were closest and
external to project border were used as the external head values for the project model.
54
The effective conductances were calculated in these cases between the project border and
the location of the external head. In cases where constant heads did not exist between
the regional and project borders, effective conductances were calculated between the two
borders as described previously. External head values for the project model were
determined by using the corresponding head values along the regional border.
3.3 SOLUTE TRANSPORT SIMULATION MODEL
In order to incorporate water quality constraints into several of the optimization
models a solute transport simulation model was needed. The DSTRAM program was
used for the project area was chosen for this application because SJRWMD had already
developed a model for another project which covers area common to the subregion
modeled with MODFLOW. The DSTRAM model incorporates an uniform grid of 30
columns by 51 row over the project area. This model uses 15 layers to simulate the
Upper Floridan aquifer, Lower Floridan aquifer, and the confining units in the system.
The majority of the information needed to create the DSTRAM project model was taken
from the regional MODFLOW model (Williams, 1994). This information includes
boundary conditions, the modeled region size and shape, hydrogeologic framework,
hydraulic characteristics, and known fluxes affecting the aquifer system. DSTRAM
simulation runs for the project were made by SJRWMD as needed during the
optimization model development.
CHAPTER 4
4.0 WATER SUPPLY NEEDS AND SOURCES
Water supply needs and sources are currently being evaluated by SJRWMD under
the District Water Management Plan. The needs and sources of Volusia County have
been recorded and projected by SJRWMD to the year 2010 for several use categories.
Projections of future water use are based on historical trends, local government
comprehensive plans, and direct communication with both public and private public
supply utilities (SJRWMD, 1992). Computer data files depicting municipal, agricultural,
and miscellaneous well information and water use needs have been obtained from the
District. Maps and tables have also been obtained which depict wastewater treatment
plant information. This data was used to incorporate water reuse in the development of
the ground water management model.
Municipal well data include the municipality and water treatment plant name, the
water treatment plant permitted capacity, existing and proposed wells which service these
areas, location of these wells in stateplanar and numericalgrid coordinates, pumping
rates for years 1988 and 1990, and the predicted pumping rate for year 2010. See Table
4.1. Within the project area, there are 97 existing and 74 proposed wells in 17 separate
well zones which supply water to 8 different municipalities or water service areas. Nine
different water treatment plants supply these service areas. These municipalities or water
service areas include Holly Hill, Port Orange, Spruce Creek, Daytona Beach, Ormond
Beach, Tymber Creek Utilities, The Trails, Inc., and a section of New Smyrna Beach,
55
Table 4.1 Volusia County municipal well data for Project area.
Project 2010
State Planar Flow Model Well 1988 1990 Projected Water Permitted
Municipal Coord. Location Cell Location Grid Cell Pumpage Pumpage Pumpage Treatment Capacity
Water Service Area X Y row col. Name (cu ft/mo) (cu ftmo) (cu ft/mo) Plant Name (mgd)
*PORT ORANGE *** Port Orange 6.21
WESTERN WELLFIELD
POW Proposed
POW Proposed
POW Proposed
457956
457862
457594
457505
457592
457946
456001
455911
455999
456175
453961
454491
454046
457948
454222
453874
453874
453874
456084
456084
456084
422871
420211
425531
422877
420218
425537
422865
420205
425525
457946
457946
458355
458355
458355
458672
458672
458672
1736744
1732805
1730684
1730381
1729068
1728664
1734120
1733009
1731999
1730888
1734022
1732708
1731597
1729674
1730486
1735638
1735638
1735638
1729070
1729070
1729070
1756899
1756905
1756894
1759929
1759935
1759924
1753869
1753874
1753864
1729674
1732805
1729068
1729674
1732805
1729068
1729674
1732805
MWELLI
MWELL9
MWELL5
MWELL5
MWELL4
MWELLIO
MWELL3
MWELL6
MWELL7
MWELLII
MWELL12
MWELL8
MWELL8
MWELL13
MWELL2
MWELL14
MWELL14
MWELL14
MWELL15
MWELL15
MWELLI5
MWELL16
MWELL17
MWELL18
MWELL19
MWELL20
MWELL21
MWELL22
MWELL23
MWELL24
MWELL25
MWELL26
MWELL27
MWELL28
MWELL29
MWELL30
MWELL31
MWELL32
1079040
1079040
1079040
1079040
1079040
1079040
1079040
1079040
1079040
1079040
1079040
1079040
1079040
1079040
1079040
1184841
1184841
1184841
1184841
1184841
1184841
1184841
1184841
1184841
1184841
1184841
1184841
1184841
1184841
1184841
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
961230
0
0
0
0
0
0
0
0
0
961230
961230
961230
961230
961230
961230
961230
961230
EASTERN WELLFIELD 491132 1745711 49 46 MWELL34 124503 136712 221822
491132 1746216 49 46 MWELL34 124503 136712 221822
491132 1746721 49 46 MWELL34 124503 136712 221822
491132 1746014 49 46 MWELL34 124503 136712 221822
490600 1745610 49 45 MWELL35 124503 136712 221822
491133 1748236 47 46 MWELL38 124503 136712 221822
491133 1748741 47 47 MWELL36 124503 136712 221822
491133 1749246 47 47 MWELL36 124503 136712 221822
490690 1749448 47 47 MWELL36 124503 136712 221822
490335 1749448 46 46 MWELL33 124503 136712 221822
489891 1749448 46 46 MWELL33 124503 136712 221822
489093 1749448 46 45 MWELL37 124503 136712 221822
488916 1749448 46 45 MWELL37 124503 136712 221822
Table 4.1 continued
Project 2010
State Planar Flow Model Well 1988 1990 Projected Water Permitted
Municipal Coord. Location Cell Location Grid Cell Pumpage Pumpage Pumpage Treatment Capacity
Water Service Area X Y row col. Name (cu fl/mo) (cu ft/mo) (cu ft/mo) Plant Name (mgd)
** DAYTONA BCH *** Daytona Beach 35.17
EASTERN WELLFIELD 478108 1765916 32 43 none 0 0 Marion Street
(7 wells inactive) 477221 1765513 32 42 none 0 0
476512 1765109 32 42 none 0 0
475182 1765110 32 41 none 0 0
474473 1764606 32 40 none 0 0
473852 1764303 32 40 none 0 0
472522 1763799 32 38 MWELL39 0 0
471990 1763496 32 38 MWELL39 2388368 2508913 1246680
471192 1763093 32 37 MWELL42 2388368 2508913 1246680
469862 1762589 32 36 MWELIA0 2388368 2508913 1246680
464897 1761482 31 32 MWEL41A 2388368 2508913 1246680
464188 1760978 31 32 MWELIAL 2388368 2508913 1246680
WESTERN WELLFIELD 462769 1760272 31 31 MWELL46 3796082 2092352 3110003 Ralph Brennan
461794 1761081 30 30 MWELLA7 3796082 2092352 3110003
458589 1747752 39 26 MWELL44 3796082 2092352 3110003
458145 1747551 39 25 MWELL52 3796082 2092352 3110003
457169 1747148 39 25 MWELL52 3796082 2092352 3110003
456371 1746846 39 24 MWELL43 3796082 2092352 3110003
455306 1746342 39 23 MWELIA5 3796082 2092352 3110003
454507 1746040 39 22 MWELLA8 3796082 2092352 3110003
454242 1746848 38 23 MWELLA9 3796082 2092352 3110003
453800 1747960 37 23 MWELL50 3796082 2092352 3110003
453269 1748769 36 22 MWELL51 3796082 2092352 3110003
WESTERN WELLFIELD 451854 1751800 34 22 MWELL53 0 2092352 3110003 South Daytona
(1988 construction) 451146 1752407 33 22 MWELL54 0 2092352 3110003
450969 1753115 33 22 MWELL54 0 2092352 3110003
451148 1753922 32 23 MWELL55 0 2092352 3110003
450972 1754933 31 23 MWELLS6 0 2092352 3110003
451062 1756044 31 23 MWELL56 0 2092352 3110003
449910 1756449 30 23 MWELL57 0 2092352 3110003
453265 1745638 38 21 MWELL58 0 2092352 3110003
452917 1750587 35 23 MWELL59 0 2092352 3110003
452297 1751093 34 23 MWELL60 0 2092352 3110003
DBW Proposed 451242 1757861 29 24 MWELL61 0 0 3110003
450977 1758872 29 24 MWELL61 0 0 3110003
450712 1759882 28 24 MWELL62 0 0 3110003
450448 1760893 27 24 MWELL63 0 0 3110003
450183 1761903 26 24 MWELL64 0 0 3110003
*SPRUCE CREEK *** 483674 1724907 56 34 MWELL66 658423 493538 326649 Spruce Creek 7.73
484295 1725513 56 35 MWELL65 658423 493538 326649
SC Proposed 484118 1725311 56 34 MWELL66 0 0 326649
483496 1724705 56 34 MWELL66 0 0 326649
484739 1725815 56 35 MWELL65 0 0 326649
0** HOLLY HILL *** I I I Holly Hill 1.43
EASTERN WELLFIELD 484940 1784396 22 53 MWELL68 69192 66195 79694
485117 1784093 22 54 MWELL67 69192 66195 79694
485294 1783689 22 53 MWELL68 69192 66195 79694
485472 1784800 21 54 MWELL69 69192 66195 79694
484409 1784700 21 53 MWELL70 69192 66195 79694
482991 1784094 21 52 MWELL71 69192 66195 79694
WESTERN WELLFIELD 461629 1773101 22 34 MWELL72 593074 567314 597705
461896 1773606 22 34 MWELL72 593074 567314 597705
462251 1774009 22 35 MWELL73 593074 567314 597705
462694 1774211 22 35 MWELL73 593074 567314 597705
461632 1775424 21 35 MWELL74 593074 567314 597705
462163 1774918 21 35 MWELL74 593074 567314 597705
463402 1773705 22 35 MWELL73 593074 567314 597705
HHW Proposed 462783 1775019 21 35 MWELL74 0 0 597705
Table 4.1 continued
Project 2010
State Planar Flow Model Well 1988 1990 Projected Water Permitted
Municipal Coord. Location Cell Location Grid Cell Pumpage Pumpage Pumpage Treatment Capacity
Water Service Area X Y row col. Name (cu ft/mo) (cu fl/mo) (cu fl/mo) Plant Name (mgd)
*NEW SMYRNA BCH New Smyrna 3.60
SAMSULA WELLFIELD 471322 1701078 60 20 MWELL75 1207921 1290193 1521308 Beach
470701 1701483 60 20 MWELL75 1207921 1290193 1521308
SR44 WF Proposed 461913 1702702 59 16 MWELL76 0 0 1521308
461380 1702804 59 15 MWELL77 0 0 1521308
462358 1704318 58 16 MWELL78 0 0 1521308
461470 1704117 58 16 MWELL78 0 0 1521308
461915 1704722 58 16 MWELL78 0 0 1521308
461560 1705127 58 16 1MWELL78 0 0 1521308
ORMOND BCH *** Ormond Beach 8.00
DIVISION AVE 479276 1795409 12 53 MWELL80 917447 645276 659368
WELLFIELD 479365 1795308 12 53 MWELL8O 917447 645276 659368
478833 1793692 13 52 MWELL81 917447 645276 659368
478479 1794197 13 52 MWELL81 917447 645276 659368
478302 1794399 13 52 MWELL81 917447 645276 659368
478124 1794096 13 52 MWELL81 917447 645276 659368
476176 1793794 12 50 MWELL79 917447 645276 659368
475644 1793290 13 50 MWELL83 917447 645276 659368
475024 1793088 13 49 MWELL84 917447 645276 659368
478480 1796116 12 53 MWELL80 917447 645276 659368
478125 1795914 12 52 MWELL82 917447 645276 659368
475822 1794704 12 50 MWELL79 917447 645276 659368
SR 40 WELLFIELD 467319 1793397 10 44 MWELL86 2582442 645276 659368
464038 1789966 11 41 MWELL87 2582442 645276 659368
467940 1794912 9 45 MWELL88 0 645276 659368
468294 1794305 10 45 MWELL89 0 645276 659368
462089 1789564 11 39 MWELL85 2582442 645276 659368
HUDSON WELLFIELD 450488 1791092 6 32 MWELL90 0 645276 659368
450489 1792002 6 32 MWELL90 0 645276 659368
450490 1792911 5 32 MWELL91 0 645276 659368
450491 1793719 4 33 MWELL92 0 645276 659368
450492 1794628 4 33 MWELL92 0 645276 659368
450494 1795537 3 33 MWELL93 0 645276 659368
450495 1796446 3 33 MWELL93 0 645276 659368
451553 1792909 5 33 MWELL94 0 645276 659368
451997 1793818 5 34 MWELL95 0 645276 659368
451998 1794727 4 34 MWELL96 0 645276 659368
452001 1797151 3 35 MWELL97 0 645276 659368
449516 1792912 5 32 MWELL91 0 645276 659368
448719 1792913 4 31 MWELL98 0 645276 659368
OB Proposed 419697 1764986 15 6 MWELL99 0 0 659368
CENTRAL RECHARGE 420758 1764075 16 7 MWELL100 0 0 659368
WELLFIELD 421643 1763265 16 7 MWELL100 0 0 659368
421995 1761951 17 7 MWELL101 0 0 659368
422080 1760436 18 7 MWELL102 0 0 659368
422166 1758820 20 6 MWELL103 0 0 659368
422163 1757708 20 6 MWELL103 0 0 659368
422160 1756193 21 5 MWELL104 0 0 659368
422159 1755385 22 5 MWELL105 0 0 659368
423769 1762452 18 8 MWELL106 0 0 659368
424833 1762652 18 9 MWELL107 0 0 659368
425719 1762549 18 10 MWELL108 0 0 659368
426872 1762749 18 11 MWELL109 0 0 659368
425544 1763762 17 10 MWELL110 0 0 659368
425280 1764671 17 10 MWELLI10 0 0 659368
420039 1759430 19 5 MWELLI1 0 0 659368
419152 1759129 18 4 MWELL112 0 0 659368
418000 1759333 18 3 MWELL113 0 0 659368
417820 1757920 19 3 MWELL114 0 0 659368
417639 1756405 20 2 MWELL115 0 0 659368
417548 1755395 21 2 MWELL116 1 0 0 659368
OBProposed 441343 1777571 13 23 MWELL117 0 0 659368
RIMARIDGE 441168 1778784 12 24 MWELL118 0 0 659368
**TYMBER CREEK UTIL *** 459436 1792900 8 39 MWELL119 372103 409982 728944 Tymber Creek 160000
THE TRAILS INC. '** 453145 1790887 7 34 MWELLI20 1299484 1299484 2355224 ???? ????
Total 100230505 105037423 174947036
(Figure 4.1).
Agricultural water service data include both agricultural and golf course well data.
This data depicts well locations in stateplanar and numericalgrid coordinates, water
application types, and pumping rates for year 1990. See Table 4.2. These rates have
been predicted to remain approximately constant through the year 2010 (Geraghty &
Miller 1992). According to data obtained from SJRWMD, there are 80 agricultural and
golf course wells throughout the Volusia County subregion. The ground water
management model was formulated to supply a portion of the ground water needs with
effluent from wastewater treatment plants. The agricultural areas and wastewater
treatment plants are shown in Figure 4.2.
The miscellaneous well data includes the well description, the well location in
stateplanar and numericalgrid coordinates, and the 1988 and 1990 pumping rates. The
Volusia County subregion consists of 9 miscellaneous wells which supply the Florida
Mining and Materials Department, the Florida Department of Education, and the Tomoka
Correctional Facility. Additional data depicting location and discharge rates of private
wells within the region for years 1988 and 2010 are also included. Data were supplied
by SJRWMD in a form to facilitate execution of the simulation model and subsequent
incorporation into the optimization model. However, these nonmunicipal wells were not
optimized during this study.
In order for the simulation/optimization model to predict a feasible water
management strategy, water service demands must be incorporated into the model.
Municipal and agricultural area water supply demands for the year 2010 were calculated
by taking the sum of all projected well discharge rates which supplied a specific water
HOLLY HILLS
TYMBER CREEK
& THE TRAILS
DAYTONA BEACH
NEW SMYRNA BEACH
o 5mile
0 5kitmeras
SPrivate areas not managed
Figure 4.1
Location of Municipal Water Service Areas
within Project Area.
F'I
Table 4.2 Volusia County agricultural well data project area.
Project 1990 2010
State Planar Flow Model Estimated Projected
Agricultural Application Coord. Location Cell Location Well Grid Pumpage Pumpage
Service Area Number X Y row col Cell Name (cu ft/mo) (cu ft/mo)
AGAREA1 21270396AN 495388 1736822 54 46 AGWELL2 55355 55355
21270396AN 496541 1736519 54 46 AGWELL2 55355 55355
21270396AN 495388 1737024 54 46 AGWELL2 55355 55355
21270396AN 496452 1737327 54 47 AGWELL1 55355 55355
21270396AN 496452 1736923 54 47 AGWELL1 55355 55355
AGAREA2 21270269AN 479539 1789550 16 51 AGWELL4 191139 191139
21270269AN 481045 1789953 17 53 AGWELL5 191139 191139
21270269AN 482639 1790256 17 54 AGWELL3 191139 191139
AGAREA3 21270279AU 505496 1758941 45 58 AGWELL6 14255 14255
21270279AU 505674 1757931 45 58 AGWELL6 14255 14255
21270279AU 506117 1758032 45 58 AGWELL6 14255 14255
21270279AU 505674 1758638 45 58 AGWELL6 14255 14255
21270279AU 505496 1758941 45 58 AGWELL6 14255 14255
21270279AU 505851 1758436 45 58 AGWELL6 14255 14255
21270279AU 505762 1758739 45 58 AGWELL6 14255 14255
21270279AU 505940 1758234 45 58 AGWELL6 14255 14255
21270279AU 505851 1758436 45 58 AGWELL6 14255 14255
21270279AU 505585 1758840 45 58 AGWELL6 14255 14255
21270279AU 505674 1758739 45 58 AGWELL6 14255 14255
21270279AU 505585 1757931 45 58 AGWELL6 14255 14255
21270279AU 506738 1756921 46 58 AGWELL7 14255 14255
21270279AU 506206 1757830 46 58 AGWELL7 14255 14255
21270279AU 506294 1757729 46 58 AGWELL7 14255 14255
21270279AU 506383 1757527 46 58 AGWELL7 14255 14255
21270279AU 506206 1756719 46 58 AGWELL7 14255 14255
21270279AU 506472 1757325 46 58 AGWELL7 14255 14255
21270279AU 506826 1757022 46 58 AGWELL7 14255 14255
21270279AU 505940 1756820 46 58 AGWELL7 14255 14255
21270279AU 506560 1757123 46 58 AGWELL7 14255 14255
21270279AU 506383 1757628 46 58 AGWELL7 14255 14255
21270279AU 506206 1757224 46 58 AGWELL7 14255 14255
21270279AU 506383 1757325 46 58 AGWELL7 14255 14255
21270279AU 506117 1756618 46 58 AGWELL7 14255 14255
AGAREA4 21270565ANV 484418 1806921 6 58 AGWELL12 92010 92010
21270565ANV 484506 1807022 6 59 AGWELL11 92010 92010
21270565ANV 485480 1805708 7 59 AGWELL10 92010 92010
21270565ANV 485745 1804799 8 59 AGWELL8 92010 92010
21270565ANV 485833 1804395 8 59 AGWELL8 92010 92010
21270565ANV 485834 1804698 8 59 AGWELL8 92010 92010
21270565ANV 485833 1803183 9 58 AGWELL9 92010 92010
AGAREA5 21270647AUS 469188 1804304 3 49 AGWELL16 193220 193220
21270647AUS 468566 1802688 4 48 AGWELL14 193220 193220
21270647AUS 470692 1803596 4 50 AGWELL13 193220 193220
21270647AUS 470071 1802182 5 49 AGWELL15 193220 193220
Table 4.2 continued
Project 1990 2010
State Planar Flow Model Estimated Projected
Agricultural Application Coord. Location Cell Location Well Grid Pumpage Pumpage
Service Area Number X Y row col Cell Name (cu ft/mo) (cu ft/mo)
AGAREA6 21270147AU 483584 1722786 56 33 AGWELL17 302712 302712
21270147AU 485536 1723189 56 35 AGWELL18 302712 302712
21270147AU 485536 1723189 56 35 AGWELL18 302712 302712
21270147AU 484118 1726321 56 35 AGWELL18 302712 302712
AGAREA7 21270236AN 465995 1701284 59 17 AGWELL19 857472 857472
21270236AN 468304 1701787 59 18 AGWELL20 857472 857472
21270236AN 466614 1698456 60 17 AGWELL21 857472 857472
AGAREA8 21270237AN 487399 1719754 57 35 AGWELL22 141258 141258
21270237AN 488907 1717128 58 35 AGWELL24 141258 141258
21270237AN 488907 1717128 58 35 AGWELL24 141258 141258
21270237AN 488907 1717128 58 35 AGWELL24 141258 141258
21270237AN 488907 1717128 58 35 AGWELL24 141258 141258
21270237AN 488907 1717128 58 35 AGWELL24 141258 141258
21270237AN 488907 1717128 58 35 AGWELL24 141258 141258
21270237AN 490416 1719652 58 37 AGWELL23 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
21270237AN 490681 1715511 59 36 AGWELL25 141258 141258
AGAREA9 21270085AN 483427 1768641 32 47 AGWELL26 24536 24536
21270085AN 483782 1768843 32 48 AGWELL28 24536 24536
21270085AN 483516 1768742 32 48 AGWELL28 24536 24536
21270085AN 484047 1768943 32 48 AGWELL28 24536 24536
21270085AN 483427 1768944 32 48 AGWELL28 24536 24536
21270085AN 483604 1769247 32 48 AGWELL28 24536 24536
21270085AN 484668 1770155 31 49 AGWELL27 24536 24536
Total 9686206 9686206
Area 4
ea \
Area 5 1'
I
, \
Daytona Beach Bethune Pt. WIP
Area 9
Il \Area 3
F1 /
/ Port Orange WrP
// Areal
/ olusia Co. Spruce Cr. WIP
SArea 6 (r 7
L j
Area 8
Area 7 r
0 ,5 miles
0 5 kilanem rs
Figure 4.2
Agricultural Areas and Wastewater Treatment Plants
with respect to Numerical Model Grid Location
Hd
\
X
64
service area. Municipal and agricultural service area demands are summarized in Tables
4.3 and 4.4, respectively, for years 1988 and 2010.
Municipal water service areas were generally defined by individual well fields
and/or water treatment plants which supply ground water to the municipalities. Proposed
well fields which were not in service before 1991 were not considered as individual water
service areas. They were incorporated in the optimization model as possible well sites
and their year 2010 withdrawal rates contributed to projected demand. However, these
locations were not specified as demand areas. The Central Recharge and Rima Ridge
well fields are two such proposed areas which will supply the Ormond Beach
municipality by the year 2010. Therefore, the year 2010 projected demand of these areas
were equally divided between the State Road 40 and Hudson well fields. These well
fields were selected because of the proximity to the proposed Central and Rima well
fields.
Defining water service areas as individual well fields is only one method of
designating demand areas. The definition of a water service area within the optimization
model can easily be revised depending on the objectives of the water resource manager.
This includes defining each municipality, or the entire Volusia County region for that
matter, as a single water service area.
Effluent from wastewater treatment plants can often be used to replace ground
water currently being used to irrigate agricultural areas and golf courses. The
wastewater treatment plant data includes tables depicting facility names and locations,
current permitted capacities and mean flow rates, and predicted permitted capacities and
mean flow rates for the year 2010. There are 9 wastewater treatment facilities in the
Table 4.3 Year 1988 and projected year 2010 demand rates for
municipal water service areas.
Projected
Year 1988 Year 2010
Well Demand Demand
ID Water Service Area Name Rate (cfd) Rate (cfd)
POW Port Orange West Wellfield 532127 916458
POE Port Orange East Wellfield 53211 94807
DBM Daytona Beach East Marion TP 392609 204934
DBB Daytona Bch West Brennan TP 1372832 1124716
DBS Daytona Beach West South TP 0 1533702
SCK Spruce Creek 43294 53695
HHE Holly Hill East Wellfield 13650 15720
HHW Holly Hill West Wellfield 136489 157205
NSB Smyrna Beach / Samsula 79425 400124
OBD Ormond Beach Division Ave WF 361952 260135
OB4 Ormond Beach State Rd. 40 WF 254706 346848
OBH Ormond Beach Hudson Wellfield 0 541950
TCU Tymber Creek Utilities 12234 24298
TTI The Trails, Inc. 42723 78507
Total 3295252 5753099
Table 4.4 Year 1988 and projected year 2010 demand rates for
agricultural water service areas.
Projected
Year 1988 Year 2010
Well Demand Demand
ID Agricultural Application Number Rate (cfd) Rate (cfd)
AGAREA1 AG0396AN 9100 9100
AGAREA2 AG0269AN 18852 18852
AGAREA3 AG0279AU 11717 11717
AGAREA4 AG0565ANV 21175 21175
AGAREAS AG0647AUS 25408 25408
AGAREA6 AG0147AU 39808 39808
AGAREA7 AG0236AN 84573 84573
AGAREA8 AG0237AN 102170 102170
AGAREA9 AG0085AN 5645 5645
Total 318448 318448
66
Volusia County subregion with permitted capacities of 0.1 mgd or greater. The
wastewater and agricultural data was reviewed and compared to determine which
agricultural areas could supplement their ground water demand with wastewater plant
effluent. From this review and based on proximity to agricultural areas, it was
determined that 5 of the 9 treatment plants could feasibly supplement the agricultural
demands. Constraints for the ground water optimization model were then formulated
based on the results of the review and comparison. See Table 4.5.
The needs and sources evaluation of Volusia County by SJRWMD includes a
water resources impact assessment. This assessment involves identifying potential
problem areas when the year 2010 projected allocation strategy is implemented. One
such problem area involves the possibility of hydrophytes being harmed as a result of
declines in water table. Based on the evaluation, several areas of Volusia County have
a high potential for vegetative harm if proposed strategies were incorporated.
In order to reduce the possibility of harm to wetland vegetation, areas depicted
as having high potential for harm were incorporated into the optimization model as
control points. These control points are locations where pressure head and/or drawdown
values were constrained or optimized. To reduce the number of control points and
therefore the size of the optimization model, 100 out of approximately 500 points were
selected throughout the high potential harm areas. Figure 4.3 displays the control point
locations along with areas of high, medium, and low potential harm.
Table 4.5 Volusia County subregional model wastewater treatment plant data.
Coordinate Numerical Grid
Location Cell Location Permit Capacity
Supply ID Plant Name latitude longitude Row Col. (mgd) (cfd) Wastewater Reuse Service Area
WASTE1 Port Orange WWTP 290812 805949 51 53 12.00 1604400 AGAREAl and/or AGAREA3
WASTE2 Holly Hill WWTP 291426 810240 21 54 2.40 320880 AGAREA2
WASTE3 Daytona BeachBethune Pt. WWTP 291205 810031 31 55 12.00 1604400 AGAREA3 and/or AGAREA9
WASTE4 Ormond Beach WWTP 291720 810426 6 52 6.00 802200 AGAREA4 and/or AGAREA5
WASTES Volusia Co. Spruce Creek WWTP 290443 810318 55 36 0.35 46795 AGAREA6 and/or AGAREA8
 I
I 
I ,
Daytona Beach WSA
I " A I
I
Port Orange WSA
 1%
I ,l .
 Low Harm
SModerate Harm
SHigh Harm
SControl Point
5 mits
I I l I
I I I I I I
0 5 kemdeas
Figure 4.3
Areas of potential impacts to plant communities resulting from
projected changes in ground water withdrawals between 1988 and 2010.
4
CHAPTER 5
5.0 GROUND WATER OPTIMIZATION MODELS
Five aquifer optimization models were developed to investigate optimal allocation
of ground water to meet year 2010 demand in the project area. These models were
developed to investigate future water allocation strategies assuming that feasible
withdrawal scenarios meet or exceed projected water service area demands and do not
exceed available water resource supplies. It was assumed with these models that adverse
environmental effects could be minimized at specific locations by constraining pressure
head changes (i.e. drawdown) to meet specified environmental goals or standards. For
example, one hundred control point locations were chosen at which ground water levels
changes were constrained. These points were in areas where native vegetation could be
harmed by declines in the surficial aquifer due to pumping. In addition to the ground
water level constraints two optimization models developed that incorporate a set of
constraints to allocate water in a manner that preserves ground water quality.
It is important to note the optimization models were developed using data
generated from both numerical flow and transport simulation models (e.g., information
describing aquifer responses to changing stresses such as pumping). The modeling grid
from the flow simulation model was incorporated in the formulations of all the
optimization models (see Chapter 3). Information on needs and sources was also
70
included in the formulation of each the optimization model (see Chapter 4). For
example, elemental discharge rates and pressure heads given by each optimization model
correspond to elemental cumulative discharges (from wells located in a grid cell) and
elemental average pressure heads in associated cells defined in the flow simulation
model. In addition, the term "well grid cells" refers to numerical model cells where one
or more wells are located, and is used herein to reflect the fact that the optimization
model identifies the cumulative well flows in each grid cell (in contrast to individual well
flows).
5.1 OPTIMIZATION MODEL DECISION VARIABLES
Each model is comprised of combinations of several groups of decision variables.
One group defines steadystate drawdowns and pressure heads at specified control points.
These are DDj, the drawdown at sensitive wetland control point; DDWh, the drawdown
at each well grid cell h; HDj, the pressure head at sensitive wetland control point; and
HDWh, the pressure head at each well grid cell h.
Another group defines cumulative and elemental flows from well grid cells used
to meet service area demands. These are QM,k, the discharge rate of each municipal
well grid cell i which supplies each municipal water service area k; QA,o, the discharge
rate of each agricultural well grid cell n which supplies each agricultural water service
area o; QW,,,, the effluent reuse rate each wastewater treatment plant m which
supplements each agricultural water service area o; QMT1, the total discharge rate at each
municipal well grid cell i; QAT,, the total discharge rate of each agricultural well grid
cell n; and QWT,, the total effluent reuse rate of each wastewater treatment plant m.
Additional decision variables used in one optimization model that incorporate
water quality constraints define steadystate chloride concentrations and changes in
concentrations at well points. These are CIh, the increase in chloride concentration at
each well grid cell control point h and CCh, the chloride concentration at each well grid
cell control point h.
5.2 OPTIMIZATION MODEL OBJECTIVE FUNCTIONS
Five optimization model formulations were identified under the assumption that
it was desirable to determine feasible ground water allocation in which minimum aquifer
system responses are maximized or maximum aquifer system responses are minimized.
For example, when the objective is to minimize the maximum drawdown then the value
of the objective function must be less than or equal to all drawdowns in the management
area. However, if the objective is to maximize a minimum pressure head, the objective
function must have a value greater than or equal to all the pressure heads in the
management area. In the optimization models presented, objective functions will appear
as statements specifying that the value of the objective function S is to be maximize or
minimize.
The objective function of the first model is to minimize the maximum drawdown
at all sensitive wetland control points (i.e. minimize S). For this model, the following
constraint was used to define S:
S DD, for all j (5.1)
where DDi is the drawdown at sensitive wetland control point.
The objective of the second model is to maximize the minimum pressure head at
all sensitive control points (i.e. maximize S). The appertain constraint associated with
this objective is:
S HDj for all j (5.2)
where HDj is the pressure head at sensitive wetland control point.
The objective of the third model is to minimize the average drawdown of all
sensitive wetland control points (i.e. minimize S). This model objective, which is equal
to maximizing average pressure heads, requires the following constraint to define S:
P
SDD (5.3)
S t J for all j
P
where p is the number of sensitive wetland control points.
The objective of the fourth model is to minimize the maximum drawdown while
constraining concentration levels. The objective function of this model is formulated the
same as the first model, equation 5.1. The difference between model one and model four
is the addition of water quality constraints to model four. These additional constraints
function to elucidate pumping strategies which satisfying specified water quality
standards.
Finally, the objective of the fifth model is to minimize the maximum relative
chloride concentration increase. This objective function is written as follows:
S RCLh (5.4)
where RCLh is the relative chloride concentration increase at all well grid cells h.
5.3 OPTIMIZATION MODEL CONSTRAINTS
Other than model specific objective functions and their appurtenant constraint
equations (i.e., Equations 5.1, 5.2, 5.3, and 5.4) the six optimization model share a large
common block of constraint equations. This block of equations includes aquifer response
constraints, management constraints, and nonnegativity constraints. The common block
and water quality constraints were used in different combinations to formulate the five
different optimization models with unique goals to identify optimal allocation strategies.
The following sections discuss these constraint in greater detail.
5.3.1 Aquifer Response Constraints
The optimization models were developed to allow pressure head and/or drawdown
to be constrained or optimized. Drawdown constraints at the specified control points
were developed using influence coefficients that describe pressure head changes at each
control point created by ground water pumpage at each well grid cell. The following
general drawdown constraint for a control point includes a linear combination of aquifer
responses to the municipal, agricultural, and private wells.
DD, = aj (QMT QMOi) + E P, (QAT QA0)
n (5.5)
+ DDPRIVW (10'6) for all i, j, and n
where DDj is the drawdown at sensitive wetland control point j, aj is the aquifer
influence coefficient defining pressure head change at each sensitive wetland area control
point due to a change in discharge rate at each municipal well grid cell i, QMTi is the
total discharge rate in cfd at each municipal well grid cell i, QMO, is the initial discharge
rate in cfd of each municipal well grid cell i, f,n is the aquifer influence coefficient
defining pressure head change at each sensitive wetland area control point j due to a
change in discharge rate at each agricultural well n, QATn is the total discharge rate in
cfd of each agricultural well grid cell n, QAO, is the initial discharge rate in cfd of each
agricultural well grid cell n, and DDPRIVj is the drawdown in feet at each sensitive
wetland area control point j due to private wells not incorporated in the optimization
process.
QMT, values were calculated by summing the discharge rates at a specific well
grid cell over all municipal service areas for which it supplies. QATn values were
calculated by summing the discharge rates at a specific well grid cell over all agricultural
service areas for which it supplies. QMO, and QAO, values were determined using the
year 1988 water allocation strategy. QMO, and QAO, were used also in calculating the
drawdown from years 1988 to 2010 at the control points. DDPRIVj is the drawdown
from year 1988 to 2010 under steadystate conditions. These values were calculated by
increasing the discharge of private wells alone in the simulation model and were used in
calculating the total drawdown at the control points. To facilitate the linear programming
75
solver, the original influence coefficients were multiplied by 10+6 to give resulting
coefficients values on the order of one. As a result, calculated DDj drawdown values
were divided by 10+6 to produce the decision variable TDDj, which has units of feet:
DD.
TDD, for all j (5.6)
10+6
where TDDj is the total drawdown at each sensitive wetland area control point.
Pressure head at the sensitive wetland control points, HD,, are determined by a
constraint equation subtracting the drawdown at control point j from the initial pressure
head:
HDj = HOj (10+6) DDi for all j (5.7)
where HDj is the pressure head in millionths of a foot at each sensitive wetland area
control point and HOj is the initial pressure head in feet at each sensitive wetland area
control point.
Initial pressure heads at sensitive wetland area control points, HOj, were
determined by executing the simulation model at year 1988 conditions. Again, a second
decision variable was created to describe pressure head in feet. The total pressure head,
THDj, values were calculated by dividing HDj by 10+6:
HD.
THD. J (5.8)
D 10+6
where THDj is the total pressure head in feet at each sensitive wetland area control point
j.
Similar to the aquifer influence coefficient matrices use to create constraint
76
equations 5.6 and 5.8 for sensitive wetland control points, influence coefficient matrices
were also developed for constraints expressing the aquifer response at each well grid
locations due to pumpage within a well grid cell and at every other well grid cell. The
following constraint equations define the pressure head at well grid cell h:
DDWh = Y, (QMTi QMO,) + E O, (QAT, QAOn) (5.9)
i n
HDWh = HWOh (106) DDWh for all h (5.10)
where DDWh is the drawdown in millionths of a foot at each well grid cell h, 7,h* is the
aquifer influence coefficient defining pressure head change at well grid cell control point
h due to a change in discharge rate at each municipal well grid cell i, Oh is the aquifer
influence coefficient defining pressure head change at well grid cell control point h due
to a change in discharge rate at each agricultural well grid cell n, HDWh is the pressure
head in millionths of a foot at each well grid cell h, and HWOh is the initial pressure
head in feet at each well grid cell h.
Initial pressure heads at well grid cells, HWOh, were determined by executing the
simulation model at year 1988 conditions and were used in calculating the year 2010
pressure head at the well grid cells. Again, i,h and On,h were multiplied by one million
to obtain values on the order of one. THDWh and TDDWh values were calculated by
dividing HDWh and DDWh by 10+6:
TDDWh = DDWh (5.11)
10+6
HDW
THDWh = h for all h (5.12)
10+6
where TDDWh is the drawdown in feet at each well grid cell h and THDWh is the total
pressure head in feet at each well grid cell h.
An additional constraint was incorporated into the optimization model to preclude
pumpage that would dewater the surficial aquifer. The following constraint was
implemented to ensure water table elevations do not decrease below a level of one foot
above the bottom of the surficial aquifer. This constraint enables the pressure heads to
be constrained in order to avoid drying out of well grid cells:
HDWh > BOTELEVh + 1.0 for all h (5.13)
where BOTELEVh is the bottom elevation of the surficial aquifer in feet from mean sea
level at each well grid cell control point h. BOTELEVh values were obtained from the
input file of the simulation model.
5.3.2 Management Constraints
Management constraints used in these optimization models define the capacity
of available resources, the demand for available resources, and the source to demand
links. The first set of constraints specify limits on capacities associated with the
production of water from aquifer systems. These following constraints were incorporated
in the model to limit the maximum withdrawal rates at the well grid cells:
QMT, = QMi CM, (5.14)
k
QAT, = QAn,o CAn (5.15)
0
QWT = E QWr,,o CWm for all i, k, m, n, and o (5.16)
where QM,k is the discharge rate in cubic feet per day (cfd) of each municipal well grid
cell i which supplies each municipal water service area k, CM, is the total capacity rate
in cfd of each municipal well grid cell i, QAn,o is the discharge rate in cfd of each
agricultural well grid cell n which supplies each agricultural water service area o, CA,
is the total capacity rate in cfd of each agricultural well grid cell n, QWT, is the total
effluent reuse rate in cfd of each wastewater treatment plant m, QW,,o is the effluent
reuse rate in cfd of each wastewater treatment plant m which supplements each
agricultural water service area o, and CWm is the capacity rate in cfd of each wastewater
treatment plant i.
Municipal well grid cell capacities, CM, (for all i), were set at 600,000 cfd. CA,
values for each agricultural well grid cell were set at the service area demand for which
the cell supplied. CWm limits reflect the available wastewater effluent which could be
used to supplement the agricultural demand.
Minimum withdrawal rates on municipal well grid cell were also incorporated into
the optimization model to prevent the shutting off of existing wells. This is achieved by
placing lower discharge limits on the well grid cells. These constraints were formulated
to require minimum flows that equal a percentage of the 1988 withdrawal rates. From
a review of the projected year 2010 water allocation strategy, discharges at existing well
79
grid cells were allowed to decrease to approximately 50 percent of the 1988 distribution
discharge rate estimates. Thus, the following constraints incorporated in the optimization
models, which requires year 2010 discharges at municipal well grid cells to be greater
than or equal to half of the year 1988 rates:
QMT, = QMi 0.50 (QMO) (5.17)
k
In order to be a feasible water allocation strategy, the strategy must meet or exceed the
demands of the service areas. The following demand constraints ensure that water needs
of municipal and agricultural service areas are satisfied:
SQMj ; DMk (5.18)
SQAn,o + QWmo DAo for all i, k, m, n, and o (5.19)
n m
where DMk is the demand rate in cfd of each municipal water service area k and DA, is
the demand rate in cfd of each agricultural water service area o.
Municipal demand, DMk, and agricultural demand, DAo, were calculated using
projected year 2010 discharge rates (See Chapter 4). These demands ensure that the
model identifies discharge from well grid cells that meet future demands of the water
service areas. It may be seen from constraint equations 5.19 that agricultural well grid
cells are not constrained to lower limits as long as demand can be satisfied with
wastewater effluent.
Also, since every water supply area does not (and can not) supply every water
service area, the optimization model was constructed to link specific water sources with
specific demand areas. The following constraints were formulated to specify which well
