WATER RESOURCES
research center
Publication No. 27
A Physical Model for Prediction and Control of Saltwater
Intrusion in the Floridan Aquifer
By
Bent A. Christensen
and
Andrew J. Evans, Jr.
Department of Civil & Coastal Engineering
University of Florida
Gainesville
'I.
96L
ui :' '. . .
w"
II a~3'~
L,
r
5L

.t~ _.I.,
,
..i . .~
;,a

UNIVERSITY OF FLORIDA
. 1
A PHYSICAL MODEL FOR PREDICTION AND CONTROL OF SALTWATER
INTRUSION IN THE FLORIDAN AQUIFER
By
Bent A. Christensen
and
Andrew J. Evans, Jr.
PUBLICATION NO. 27
FLORIDA WATER RESOURCES RESEARCH CENTER
RESEARCH PROJECT TECHNICAL COMPLETION REPORT
OWRR Project Number A022FLA
Annual Allotment Agreement Number
143100013809
Report Submitted: January 29, 1974
The work upon which this report is based was supported in part
by funds provided by the United States Department of the
Interior, Office of Water Resources Research as
Authorized under the Water Resources
Research Act of 1964
ACKNOWLEDGEMENTS
Appreciation is extended to the Office of Water
Resources Research, United States Department of the
Interior, for financial support of this project, and to
Dr. W. H. Morgan for his administrative assistance, and
Ms. Mary Robinson for accounting assistance. The par
ticipation of Mr. Floyd Combs in all aspects of the
project was especially valuable.
The efforts of Mr. C. L. White, Mr. William
Whitehead and the staff of the Mechanical Engineering
machine shop in the construction of the model are grate
fully acknowledged. The assistance of Mr. Richard Sweet
in the operation of the model and the patient attention
of Ms. Cheryl Combs in typing the report deserve special
thanks.
TABLE OF CONTENTS
Page
ACKNOWLEDGEMENTS.................................. ii
LIST OF TABLES .................................. .. v
LIST OF FIGURES..................................... vi
KEY TO SYMBOLS OR ABBREVIATIONS..S.................. viii
ABSTRACT........................................... X
CHAPTER
I. INTRODUCTION AND STATEMENT OF PROBLEM ........ 1
Topography...... ....................... 2
Western Highlands..................... 2
Marianna Lowlands..................... 2
Tallahassee Hills..................... 4
Central Highlands .................... 4
Coastal Lowlands...................... 4
Climate................................ 5
Geology.............................. .. 5
II. MODELS, NUMERICAL AND PHYSICAL................ 8
A: NUMERICAL METHODS.......................... 8
Method of Finite Differences............ 8
Method of Finite Elements................ 9
Relaxation Methods...................... 9
B: PHYSICAL MODELS................... ....... 10
The Sandbox Model........................ 11
The HeleShaw Analog.................... 12
Electric Analog......................... 12
The Continuous Electric Analog........ 12
The Discrete Electric Analog.......... 12
The Ion Motion Analog................. 13
The Membrane Analog..................... 13
Summary............. ........ ......... 13
III. THE HELESHAW MODEL........................... 16
Viscous Flow Analog..................... 16
Scaling ................................. 21
Time................................ 23
Anisotropy ........ .... ............ .. 24
Leakage........................ ... .. 29
Storativity........................... 31
Discharge............................. 32
Accretion.................... ...... 33
Volume................................. 33
iii
TABLE OF CONTENTS (Continued)
Page
IV. SITE SELECTION AND PROTOTYPE GEOLOGY
AND HYDROLOGY................................. 35
Site Selection........................... 35
Prototype Geology ........................ 37
Prototype Hydrology..................... 39
V. DESIGN, CONSTRUCTION AND OPERATION OF MODEL... 40
Design ................................... 40
Prototype............................... 40
Model................ ................... 40
Construction............................ 50
Frame.................................. 50
Plexiglas Plates and Manifolds.......... 52
SaltWater System ..................... 58
FreshWater System, General............. 64
FreshWater System, Accretion........... 64
FreshWater System, Wells............... 64
FreshWater System, Flow Meters......... 64
Operation ................................. 71
VI. RESULTS, CONCLUSIONS AND RECOMMENDATIONS...... 72
Results...................... ............ 72
Conclusions and Recommendations........... 86
LIST OF TABLES
TABLE Page
2.1 APPLICABILITY OF MODELS AND ANALOGS........... 14
5.1 PROTOTYPE PARAMETERS.......................... 41
5.2 MODEL PARAMETERS .............................. 46
5.3 SIMILARITY RATIOS.............................. 47
LIST OF FIGURES
FIGURE Page
1.1 TOPOGRAPHIC DIVISIONS OF FLORIDA.............. 3
1.2 MEAN ANNUAL PRECIPITATION..................... 6
3.1 FREE BODY FLOW DIAGRAM HELESHAW MODEL........ 18
3.2 PLAN SECTION OF ANISOTROPIC GROOVED ZONE
IN HELESHAW MODEL .......................... 26
3.3 HEADLOSS IN GROOVED ANISOTROPIC ZONE.......... 26
3.4 PLAN SECTION OF LEAKY ZONE IN HELESHAW MODEL. 30
3.5 STORATIVITY .... ................... ........ .. 30
4.1 REGIONAL AREA OF PROTOTYPE.................... 36
5.1 GHYBENHERZBERG INTERFACE MODEL............... 42
5.2 CRADLE AND CRADLE DOLLY....................... 51
5.3 CRADLE ROTATION................................. 51
5.4 FRAME AND MODEL SETUP ........................ 53
5.5 AIR HOSE AND VALVE ARRANGEMENT................ 53
5.6 STUB SHAFT AND PILLOW BLOCK ARRANGEMENT....... 54
5.7 BACKUP AIR SUPPLY............................ 54
5.8 INTERNAL SUPPORT AND SEALING SYSTEM........... 55
5.9 MODEL MOUNTING SYSTEM......................... 56
5.10 MODEL BACK AND FRONT PLATES................... 56
5.11 MODEL BACK PLATE, DETAIL...................... 57
5.12 MODEL BACK PLATE, DETAIL...................... 57
5.13 FRONT PLATE WITH ACCRETION MANIFOLDS.......... 59
5.14 BACK AND FRONT PLATE CLAMPUP................. 59
5.15 ACCRETION MANIFOLDS........................... 60
5.16 ACCRETION MANIFOLDS, DETAIL................... 60
5.17 FLUID SUPPLY NETWORK SCHEMATIC................ 61
5.18 CONNECTIONS BETWEEN MODEL AND FLUID
SUPPLY SYSTEM... ................ .... ....... 62
5.19 SALTWATER RESERVOIR AND PUMP................. 62
LIST OF FIGURES (Continued)
5.20 SALTWATER CONSTANT HEAD TANK................
5.21 FRESHWATER SUPPLY SYSTEM...................
5.22 FRESHWATER RESERVOIR AND ACCRETION PUMP.....
5.23 WELL SUPPLY MANIFOLD AND PUMP................
5.24 WELL SUPPLY MANIFOLD AND ACCRETION
SUPPLY MANIFOLD............................
5.25 OPPOSITE VIEW OF FIGURE 5.24.................
5.26 FLOW METER BANK ..............................
5.27 FLOW METER DETAIL...........................
5.28 FLOW METER TO MODEL CONNECTIONS..............
5.29 FLOW METER PRESSURE SENSING LINES............
5.30 FLOW METER SWITCHING DEVICE AND
PRESSURE TRANSDUCERS.......................
5.31 CARRIER DEMODULATOR AND STRIP CHART RECORDER.
INTERFACE
INTERFACE
INTERFACE
INTERFACE
INTERFACE
INTERFACE
INTERFACE
INTERFACE
INTERFACE
INTERFACE
INTERFACE
INTERFACE
INTERFACE
LOCATION,
LOCATION,
LOCATION,
LOCATION,
LOCATION,
LOCATION,
LOCATION,
LOCATION,
LOCATION,
LOCATION,
LOCATION,
LOCATION,
LOCATION,
tm = 0 Min...............
tm = 32 Min..............
t = 48 Min..............
t = 62 Min..............
tm = 79 Min..............
t = 92 Min.............
m = 102 M....
t = 107 Min.............
t = 117 Min.............
m
t = 132 Min.............
t = 132 Min.............
m
tm = 147 Min.............
tm = 162 Min.............
tm = 172 Min.............
vii
6.1
6.2
6.3
6.4
6.5
6.6
6.7
6.8
6.9
6.10
6.11
6.12
6.13
Page
KEY TO SYMBOLS OR ABBREVIATIONS
SYMBOLS
A = area, Fourier coefficient
b = width of model and prototype
B = body force, Fourier coefficient
f = subscript denoting freshwater, Floridan aquifer
g = acceleration due to gravity
h =height of water
j = summation limit
k = intrinsic permeability
K = hydraulic conductivity
1 = distance between storativity tubes, subscript
denoting leaky layer
L = distance center to center of groove in
anisotropic zone, flow meter tube length
m = subscript denoting model
n = effective porosity
e
p = pressure, subscript denoting prototype
q = specific discharge
Q = total flow
r = subscript denoting ratio
R = accretion
R = Reynolds' number
S = specific storage
o
t = time
T = transmissivity
U = volume
V = velocity
x = horizontal direction parallel to test section
y = horizontal direction perpendicular to test
section
z = vertical direction
viii
KEY TO SYMBOLS OR ABBREVIATIONS (Continued)
1,2,3 = subscripts denoting zone 1, zone 2 and zone 3
a = width adjustment factor
y = unit weight
X = length adjustment factor
= absolute viscosity
5 = geometric parameter in anisotropic zone
p = mass density
u = kinematic viscosity
= potentiometric head
= velocity potential, potential
ABBREVIATIONS
D = substantial derivative
Dt
V2 = Laplace operator
A = difference operator
OC = Centigrade
cfs = cubic feet per second
oF = Fahrenheit
fps = feet per second
g/cm3 = grams per cubic centimeter
gpd = gallons per day
Hg = Mercury
I.D. = inner diameter
msl = mean sea level
O.D. = outer diameter
pcf = pounds per cubic feet
psid = pounds per square inch differential
psig = pounds per square inch gage
RPM = revolutions per minute
ABSTRACT
Continuing development of the coastline zone in
the middle Gulf area of Florida is increasing the demand
for ground water supplies, and in turn increasing the
probability of saltwater intrusion. Methods must be
developed to make longrange predictions on the effects of
increased demands on the Floridan aquifer.
A HeleShaw model is a physical model which fits
the requirements for longrange planning. It is well
suited to handling anisotropic aquifers, difficult boundary
conditions and can simulate years of field conditions in
minutes of model time.
The section selected for study lies in a line from
the Gulf coast nar Tarnnn Anrinag to a Doint near Dade
City and passes through the EldridgeWilde well field. The
EldridgeWilde well field is the major water producer for
pinelas County. This region has experienced several years
of dry weather, and pumping has lowered the water levels in
the aquifer by a significant amount. This loss of fresh
water head is certain to induce saltwater intrusion.
A HeleShaw model has been built for this section, and
all pertinent geological and hydrological features of the
area are included. Steadystate characteristics of the
aquifer system have been considered. In particular, the
longterm effects due to pumping and artificial recharge
were examined.
CHAPTER I
INTRODUCTION AND STATEMENT OF PROBLEM
Florida, with the possible exception of California,
is the fastest growing state of the United States. The
rapid influx of people since World War II has greatly in
creased the demands for land and water. In the past, there
has been an almost total lack of widerange planning for
the uses of these resources. Even fewer investigations
have been made into the consequences of their rapid and un
ordered development. Recently, water supplies have had to
be rationed in South Florida. Overall, land and wetlands
required for fish and wildlife have so diminished that, in
some instances, there has been a marked decrease in their
numbers. It seems reasonable to conclude that in some areas
of the state, land and water resources cannot support much
larger populations with current locally available supplies
without almost irrecoverable damage to the ground water
system in tne rorm or saltwater encroachment.
With the growing affluence of the American people,
and the availability of economically priced airconditioning
units, it can be expected that even more people will leave
the colder northern climates for the southern and western
states. Florida can expect to receive more than its share
of the migration. Frequently now, environmental protection
groups are making forecasts of impending doom. At worst,
their predictions may come true and people are beginning to
look at all growth with a jaundiced eye.
It is doubtful, however, that growth and development
can be stopped. The history of man indicates a continual
effort to better his life style, his private environment.
There is little doubt that this has sometimes caused a deg
radation of other portions of his world. Unless the cessa
tion of all growth and development is acceptable, new ways
must be found of forecasting, or predicting, the results of
all growth so as to combat possible undesirable results.
Consequences of all growth must be known, even of those re
sulting when pragmatic shortterm solutions are used.
Hopefully, the remainder of this report will present
a modeling method which will be useful in forecasting the
results of pumpage and use of ground water in our coastal
zones so that we may better plan their usage. But first,
a little background on Florida.
i
F
Topography
Florida lies between latitudes 24o40' and 3100'
North, and longitudes 80000' and 87040' West, and is the
most southerly unit of the continental United States. In
its southernmost extension it is less than 10 of latitude
north of the Tropic of Cancer.
Florida is bounded on the east by the Atlantic
Ocean; on the south by the Straits of Florida and the Gulf
of Mexico; on the west by the Gulf of Mexico and the state
of Alabama; and on the north by Alabama and Georgia. The
shape of the state in relation to the remainder of the
United States suggests two distinctive parts: the Floridan
panhandle and the peninsula of Florida. The panhandle is
a strip roughly 225 miles long that stretches in an east
west direction. The peninsula is a southsoutheast exten
sion at approximate bearing S 170 E. From the northern
boundary of the state to the tip, not including the chain
of keys, the peninsula is approximately 415 miles long and
includes 2/3's of the land mass of the entire state. Its
coastline, some 1350+ miles long, is the longest with the
exception of Alaska. No place in the interior of Florida
is more than 60 miles from either the Gulf or the Atlantic
coast.
Cooke (1945) divided the terrain of Florida into
five sections, Figure 1; the Western Highlands, the Marianna
Lowlands, the Tallahassee Hills, and a narrow band of
Coastal Lowlands, which comprise the panhandle, and the
Central Highlands and Coastal Lowlands, which comprise the
peninsula. The topography of each is described briefly.
Western Highlands
Extending eastward from the Perdido River (the west
ern boundary of Florida) to the Apalachicola River, the
northern part of this region near the Alabama stateline is
not much higher than 300 feet. It is considered to be
hilly when compared to the broad gently rolling southern
parts of this region which drop to 100 feet elevation as one
approaches the coastal lowlands. The highest elevation in
the state, 345 feet, is found in this region in the north
west corner of Walton County. The western highlands are
underlain with sand of the Pliocene Citronelle Formation.
The steepness of the bankslopes at the headwaters of the
many streams is the most unique physiographic characteristic
of this section.
Marianna Lowlands
This roughly triangularshaped region of Holmes,
Jackson and Washington Counties, with somewhat smaller con
tributions from Bay and Calhoun Counties lies between the
Tallahassee hills and the western highlands. It is difficult
FIGURE 1.1
TOPOGRAPHIC DIVISIONS OF FLORIDA
Western Highlands
Marianna Lowlands
Tallahassee Hills
Central Highlands
Coastal Lowlands
I
86
0 50 100 150
Scale of Miles
50 Fathoms
After:
Cooke (1939)
to recognize this area of gently rolling hills as lowlands.
Cooke (1945) attributes the lower elevations to the solu
bility and consequent degradation of the underlying lime
stone. This area is one of the two in the state where the
Ocala Formation is exposed to the surface and the only area
of the state where the Marianna limestone, the soft white
limestone of the Oligocene Age, is found exposed. The
region is dotted with sinks, sinkhole lakes and springs.
Tallahassee Hills
From the Apalachicola River east to the Withlacoochee
River, the Tallahassee hills extends along the Georgia
Florida border and is only 25 miles in width. The western
section is a nearly level plateau some 300 feet above mean
sea level. The remainder consists of rolling hills carved
out of the Citronelle Formation. In addition to this, a
red clayey sand and Fuller's earth of the Hawthorn Formation
are found in this area. This is a fertile farming region.
Central Highlands
The central highlands forms the backbone of the
Floridan peninsula and extends from the Georgia line between
the Withlacoochee and St. Mary's Rivers southsoutheastward
some 250 miles into Glades County west of Lake Okeechobee.
This region is highly diversified. It includes high swampy
plains, hills and innumerous lakes. Soils are sandy. Much
of them were derived from Pleistocene (Ice Age) marine ter
races. However, a distinguishable amount comes from the
Miocene Hawthorn and Pliocene Citronelle Formations. The
lakes and sinks which dot the entire area indicate the
presence of limestone below the surface. Elevations of this
region average just slightly more than 150 feet; however,
they vary from less than 100 feet to approximately 300 feet.
Coastal Lowlands
The coastal lowlands, or coastal plains as it is
sometimes called, borders the entire 1350 mile Florida
coastline. Flanking on both sides of the central highlands,
the coastal lowlands is widest just south of Lake Placid
and narrowest between the western border and the Chocta
whatchee Bay just south of the western highlands. The ele
vations everywhere within this region are less than 100 feet.
The soil for the most part is sandy except in the Everglades
and Big Cypress Swamp locales where Pliocene limestone, muck
and peat prevail near the surface. The keys, which extend
some 100 miles into the Straits of Florida, are mostly sandy
oolitic limestone like that of the mainland; however, some
limestone with coral heads is found. The islands seldom
reach 15 feet elevation. The entire region is generally
flat, typical of recently deposited material with little or
no erosion.
Climate
The sea surface temperatures east and west of Florida
average, respectively, 780 and 770 Fahrenheit. Water tem
peratures range from 740 to 830 Fahrenheit in the east and
70 to 840 Fahrenheit in the west. The coldest month in
both cases is February; the warmest month is likewise August.
The relative homogeneous distribution of sea temperature,
the lack of high relief and the peninsula shape of Florida
contribute greatly to its climate. The temperature is every
where subtropical. Mean annual average temperature in the
north is 68* Fahrenheit, and in the southern tip 750 Fahren
heit.
The tradewinds, which shift from northern Florida
to southern Florida and back semiannually, bring a mildly
monsoon effect to Florida. In November, the tradewinds are
at their southernmost extension and Florida's climate is
controled by frontal, or cyclonic, activity moving in from
the continental United States. Rainfall during this period
is of low intensity and longer duration.
Beginning in early May, the tradewinds move north,
again bringing with them the moist warm air of the Atlantic.
The cyclonic activity is greatly reduced over the state and
convectional instability begins to become established. June
through September is known as the rainy season in Florida.
The thunderstorms of this period are intense and very spa
cially varied. They usually occur during the hottest part
of the day and only on rare occasions last longer than two
hours. About 60 percent of the total average annual rain
fall occurs during this period, Figure 2. The mean average
rainfall of Florida is in the neighborhood of 53 inches.
It varies from 38 to 40 inches in the lower keys to over
65 inches in the southeast corner of the peninsula and the
western portion of the panhandle. Most of the interior,
that is the central highlands, receive approximately the
mean annual average.
Geology
The Floridan peninsula and the offshore submerged
lands above 50 fathoms, which Vaughn (1910) called the
Floridian Plateau, have existed for several million years.
The region has not been subject to violent earth movement
and, consequently, there has been a gentle doing resulting
in the formation of an oval arch above the basement rock.
The rock of the core underlying the plateau is hypothesized
to be preCambrian; however, no drill has penetrated the
core. The oldest rocks penetrated, to date, are a quartz
ite found at about 4500 feet below the surface in Marion
County. The borehole encountered another 1680 feet of
FIGURE 1.2
MEAN ANNUAL PRECIPITATION
(In Inches)
SOver 66
S62 66
58 62
L:I 54 58
1 50 54
S46 50
SLess than 46
46\
0 50 100 150
Scale of Miles
Note: Precipitation normals
compiled from records pub
lished in Climatological
Data: Florida Section, Dec,
1971.
__
quartzite before drilling was suspended. This metamorphized
rock, believed to be a continuation of the Piedmont region
of Georgia, was assigned by Cooke (1945) to the Pennsylvanian
period. The arch above the metamorphized basement, composed
of almost pure porous limestone, is known as the Lake City,
Avon Park and Ocala Formations. Dated in the Eocene period,
the Ocala Formation has an estimated maximum thickness of
360 feet. It is found at, or above, mean sea level through
out northeast and north central Florida and is this section's
principal aquifer. In southern Florida, in the vicinity of
the Everglades, the Ocala is found at depths approaching
1200 feet. The Lake City and Avon Park limestones found be
low the Ocala are the principal aquifer used by agricultural
interests in central and south central Florida and are known
locally as the Floridan aquifer.
Above the Eocene series are the formations of the
Oligocene epoch. These are represented by the Marianna
limestone and the Byran limestones found and mined in the
Marianna lowlands of the northern part of the state, and the
Suwanne limestone found over the Ocala Formation as far
south as Hillsborough County.
The next higher formations are those of the Miocene
epoch. These are well represented by the Tampa limestone
of the early Miocene which are found above the Suwanne and
Ocala limestone in south Florida, the Chipola and Shoal
River Formations of the Alum Bluff group found in northwest
and north central Florida, the Hawthorn Formation and the
Duplin marls. The latter three formations, Alum Bluff,
Hawthorn and Duplin, chiefly are sands, clays and marls
that form a confining layer over the Eocene and Oligocene
limestones.
The Hawthorn, with the possible exception of the
Ocala, is the most extensive formation within the state.
It occurs at, or near, the surface in most of north Florida.
It overlies the Tampa limestone formation in Hillsborough
County, and is, itself, overlain by the Duplin marls and
younger deposits in the south central and southern parts of
the state.
The surface material of most of the coastal lowlands
are of the Pliocene, Pleistocene and Recent periods. The
most widely distributed are the sands formed along the old
shorelines of previous ocean levels. Cooke (1945) defines
seven of these marine terraces. Some small deposits of
coquina, oolite, coral reef limestone and freshwater marls
are found among these deposits.
CHAPTER II
MODELS, NUMERICAL AND PHYSICAL
The purpose of this chapter is to enumerate some
of the more widely used modeling techniques in ground
water flow, along with a brief description of each. The
reader is referred to Bear (1972) for additional informa
tion and references.
A: NUMERICAL METHODS
Numerical methods are used in many cases where the
partial differential equations governing flow through
porous media cannot be solved exactly. Various techniques
have been developed for obtaining numerical solutions.
Method of Finite Differences
The method of finite differences is one such tech
nique. The first step is to replace the differential
equations by algebraic finite difference equations. These
difference equations are relationships among values of the
dependent variable at neighboring points of the applicable
coordinate space.
The resulting series of simultaneous equations is
solved numerically and gives values of the dependent vari
ables at a predetermined number of discrete or "grid"
points throughout the region of investigation.
If the exact solution of the difference equations is
called D, the exact solution of the differential equation
is called S, and the numerical solution of the difference
equation is called N, two quantities of interest may be
defined. They are the truncation error, IS DI, and the
roundoff error, ID NI. In order for the solution to
converge, it is necessary that IS DI > 0 everywhere in
the solution domain. The stability requirement is such
that ID NI => 0 everywhere in the solution domain. The
general problem is to find N so that IS NJ is smaller
than some predetermined error. Noting that (S N) =
(S D) + (D N), it is seen that the total error is com
posed of the truncation error and the roundoff error.
The arbitrary form selected for the finite difference equa
tion leads to the truncation error. This error is fre
quently the major part of the total error.
The actual computation proceeds by one of two schemes.
They are the explicit, or forwardintime, scheme and the
implicit, or backintime, scheme. The explicit scheme is
simpler but more time consuming than the implicit scheme due
to the stability constraint. The implicit scheme is more
efficient but requires a more complicated program as compared
to the explicit scheme.
Method of Finite Elements
The finite element technique employs a funtional asso
ciated with the partial differential equation, as opposed to
the finite difference method which is based on a finite dif
ference analog of the partial differential equation. A cor
respondence which assigns a real number to each function or
curve belonging to some class is termed a functional.
The calculus of variations is employed to minimize
the partial differential equation under consideration. This
is done by satisfying a set of associated equations called
the Euler equations. Thus, one seeks the functional for
which the governing equations are the Euler equations and
proceeds to solve the minimization problem directly, rather
than solving the differential equation.
The procedure is continued by partitioning the flow
field into elements, formulating the variational functional
within each element and taking derivatives with respect to
the dependent variables at all nodes of the elements. The
equations of all the elements are then collected. The bound
ary condition is expressed in terms of nodal values and in
corporated into the equations. The equations are then solved.
Relaxation Methods
This method may be applied to steadystate problems
which are adequately described by the Laplace or Poisson
equations. The process involves obtaining steadily improved
approximations of the solution of simultaneous algebraic
difference equations.
The first step of the procedure is to replace the
continuous flow domain under investigation by a square or
rectangular grid system. The governing differential equation
is also replaced by corresponding difference equations.
Next, a residual, say R is defined corresponding to the
point o on the grid. Ro represents the amount by which the
equation is in error at that point. If all values of the
equation are correct, R will be zero everywhere. In the
initial computational step, values are assigned at all grid
points and in general, the initial residuals will not be
zero everywhere. The process now consists in adjusting
values at each point so that eventually all residuals ap
proach zero, to at least some required accuracy.
The reduction of residuals is achieved by a "relaxa
tion pattern" which is repeated at different grid points so
as to gradually spread the residuals and reduce their value.
B: PHYSICAL MODELS
As implied in section A, direct analytical solutions
are frequently inadequate or impractical for engineering
application. In many cases, the analytical solutions which
are found are difficult to interpret in a physical context.
In an attempt to circumvent some of the shortcomings of a
purely mathematical approach, model and analog methods are
frequently employed. The analog may be considered as a
single purpose computer which has been designed and built
for a given problem.
Modeling, then, is the technique of reproducing the
behavior of a phenomenon on a different and more convenient
scale. In modeling, two systems are considered: the proto
type, or system under investigation, and the analog system.
These systems are analogous if the characteristic equations
describing their dynamic and kinematic behavior are similar
in form. This occurs only if there is a onetoone corre
spondence between elements of the two systems. A direct
analogy is a relationship between two systems in which cor
responding elements are related to each other in a similar
manner.
A model is an analog which has the same dimensions
as the prototype, and in which every prototype element is
reproduced, differing only in size. An analog is based on
the analogy between systems belonging to entirely different
physical categories. Similarity is recognized in an analog
by two characteristics: (1) for each dependent variable
and its derivatives in the equations describing one system,
there corresponds a variable with corresponding derivatives
in the second system's equations, and (2) independent vari
ables and associated derivatives are related to each other
in the same manner in the two sets of equations. The ana
logy stems from the fact that the characteristic equations
in both systems represent the same principles of conserva
tion and transport that govern physical phenomena. It is
possible to develop analogs without referring to the mathe
matical formulation; an approach which is particularly ad
vantageous when the mathematical expressions are excessively
complicated or are unknown.
Analogs may be classed as either discrete or continu
ous with respect to space variables. In both cases, time
remains a continuous independent variable.
The need for complete information concerning the flow
field of a prototype system is obvious, and no method of
solution can bypass this requirement. However, in many
practical cases involving complicated geology and boundary
conditions, it is usually sufficient to base the initial
construction of the analog on available data and on rough
estimates of missing data. The analog is then calibrated
by reproducing in it the known past history of the proto
type. This is done by adjusting various analog components
until a satisfactory fit is obtained between the analog's
response and the response actually observed in the proto
type. Once the analog reproduces past history reliably,
and within a required range of accuracy, it may be used to
predict the prototype's response to planned future operations.
The Sandbox Model
A reduced scale representation of a natural porous
medium domain is known as a sandbox model, or a seepage
tank model. Inasmuch as both prototype and model involve
flow through porous media, it is a true model.
A sandbox model is composed of a rigid, watertight
container, a porous matrix filler (sand, glass beads, or
crushed glass), one or several fluids, a fluid supply system
and measuring devices. The box geometry corresponds to that
of the investigated flow domain, the most common shapes
being rectangular, radial and columnar. For onedimensional
flow problems, the sand column is the most common experimen
tal tool. Transparent material is preferred for the box
construction, especially when more than one liquid may be
present and a dye tracer is to be used. Porosity and permea
bility variations in the prototype may be simulated by vary
ing the corresponding properties of the material used as a
porous matrix in the model according to the appropriate scal
ing rules. The porous matrix may be anisotropic. In order
to measure piezometric heads and underpressures, piezometers
and tensiometers may be inserted into the flow domain of
the model.
Wall effects are often eliminated by gluing sand
grains to the walls of the box. This effect can also be re
duced by making the porous matrix sufficiently large in the
direction normal to the wall. Inlets and outlets in the
walls connected to fixedlevel reservoirs or to pumps are
used to simulate the proper boundary and initial conditions
of the prototype.
Water is usually used in models which simulate ground
water aquifers, although liquids of a higher viscosity may
be used to achieve a more suitable time scale.
The sandbox model is used extensively because of its
special features which permit studies of phenomena related
to the microscopic structure of the medium such as: hydro
dynamic dispersion, unsaturated flow, miscible and immiscible
displacement, simultaneous flow of two or more liquids at
different relative saturations, fingering, wettability and
capillary pressure. The capillary fringe in a sandbox model
is disproportionately larger than the corresponding capillary
rise in the prototype, and for this reason the sandbox model
is usually used to simulate flow under confined rather than
phreatic conditions.
The HeleShaw Analog
The HeleShaw or viscous flow analog is based on the
similarity between the differential equations governing two
dimensional, saturated flow in a porous medium and those
describing the flow of a viscous liquid in a narrow space
between two parallel planes. In practice, the planes are
transparent plates, and the plates are usually mounted in a
vertical or horizontal orientation.
The vertical HeleShaw analog was selected for this
study because it is more appropriate for the prototype system
under investigation. Also, it is not possible to model a
free ground water table or percolation in a horizontal model.
A detailed description of this analog is presented
in Chapter III of this report.
Electric Analog
Three types of electric analogs are powerful tools
in the study of flow through porous media. They are the
continuous electric analog, the discrete electric analog and
the ion motion analog.
The Continuous Electric Analog
This analog takes two forms: the electrolytic tank
and the conducting paper analogs. The analogy rests on the
similarity between the differential equations that govern
the flow of a homogeneous fluid through a porous medium, and
those governing the flow of electricity through conducting
materials.
In particular, Darcy's law for flow in a porous
medium and Ohm's law for the flow of an electric current in
a conductor may be compared. Also, the continuity equation
for an incompressible fluid flowing through a rigid porous
medium may be compared with the equation for the steady flow
of electricity in a conductor. One concludes from this com
parison that any problem of steady flow of an incompressible
fluid having a potential may be simulated by the flow of
electric current in an analog.
The Discrete Electric Analog
This analog also takes two forms: the resistance
network analog for steady flow, and the resistancecapacitance
network for unsteady flow.
In this analog, electric circuit elements are concen
trated in the network's node points to simulate the properties
of portions of the continuous prototype field around them.
The unknown potentials are the solution of the problem, and
they can only be obtained for those points which correspond
to the nodes of the analog network. The discrete electric
analog is based on the finitedifference approximation of
the equations to be solved; therefore, the errors involved
in the discrete representation are the same as those occur
ing in this approximation.
The electric resistor corresponds to the resistance
of soil to flow through it, and capacitors are used at the
nodes to simulate storage capacity of the prototype.
The Ion Motion Analog
This analogy uses the fact that the velocity of ions
in an electrolytic solution under the action of a DC voltage
gradient is analogous to the average velocity of fluid par
ticles under imposed potential gradients in a porous medium.
In this case, both electric and elastic storativities are
neglected. The primary advantage of the ion motion analogy
is that, in addition to the usual potential distribution,
it permits a direct visual observation of the movement of an
interface separating two immiscible fluids. In ground water
interface problems where gravity is involved, this analog
cannot be used. Scaling for the analog is based on the sim
ilarity between Darcy's law and Ohm's law governing the ion
motion in an electrolytic solution.
Physically the analog consists of an electrolytic
tank having the same geometry as the investigated flow domain.
Inflow and outflow boundaries are simulated by positive and
negative electrodes, and two and threedimensional flow
domains may be investigated.
The Membrane Analog
The membrane analog consists of a thin rubber sheet,
stretched uniformly in all directions and clamped to a flat
plane frame. The achievement of equilibrium of various
forces and stresses in the membrane (caused by distorting
the frame or transversal loads) leads to the Laplace equation
and the Poisson equation. The analogy is based on the simi
larity between these two equations and the corresponding
equations that describe the flow in the prototype.
This method is applicable mainly to cases of steady
twodimensional flow involving complicated boundary geometry
and point sources and sinks within the flow field.
Summary
Following Bear (1972), Table 2.1 is presented as a
summary of the models and analogs discussed in section B of
this chapter. In section A, the numerical methods discussed
are most likely to be carried out on a digital computer. It
is important for the investigator to examine both the cost
and the applicability of these various numerical and physical
11111
TABLE 2.1
APPLICABILITY OF MODELS AND ANALOGS
Characteristic
Sandbox Model
HeleShaw Analog
Vertical Horizontal
Electric Analogs
Electrolytic RC Network
Dimensions of field
Steady or unsteady flow
Simulation of phreatic surface
Simulation of capillary
fringe and capillary pressure
Simulation of elastic
storage
Simulation of anisotropic
media'
Simulation of medium
inhomogeneity
Simulation of leaky formations
Simulation of accretion
Flow of two liquids with
an abrupt interface
Simultaneous flow of two
immiscible fluids
Hydrodynamic dispersion
Observation of streanlines
two or three two
both both
yes1 yes'
yes yes
yes, for two yes
dimensions
yes yes
kxyes yes
yes yes
yes yes
yes yes
approximately yes
yes no
yes no
yes, for two yes
dimensions, near
transparent walls
for three dimen
sions
two
both
no
no
yes
yes
k xk
x y
yes5
yes
yes
yes
r.o
no
yes
two or three
steady
yes2
no
yes, for two
dimensions
yes
yes
yess
yes, for two
dimensions
no6
no
no
no
two or three
both
no'
no
two(horizontal)
steady
no
no
yes
kx ky
yes
two(horizontal)
steady
no
no
no
yes
kxfky
yes
no
yes
no
Subject to restrictions because of the presence of a capillary fringe.
By trial and error for steady flow.
By trial and error for steady flow, or, as an approximation, for relatively small phreatic surface fluctuations.
By scale distortion in all cases, except for the RC network and sometimes the HeleShaw analog where the hydraulic
conductivity of the analog can be made anisotropic.
With certain constraints.
For a stationary interface by trial and error.
Ion Motion
Membrane
Analog
~al~m~~slr*RnTi~anrr~li~ I
15
methods to his particular case. An analog is usually pre
ferred to a digital solution when the accuracy and/or
amount of field data is small. In many simple cases, the
analog is likely to be less expensive than a digital com
puter; whereas, for large regions or unsteady three
dimensional problems, the computer may be less expensive.
The HeleShaw model also has definite advantages
when demonstration of the saltwater intrusion phenomenon
to a public body, or other laymen involved in political
decisionmaking, is considered. This type of model allows
for direct observation of the phenomenon without the
numerical interpretations used in the computer models.
CHAPTER III
THE HELESHAW MODEL
The viscous flow analog, more commonly referred to
as a parallelplate or HeleShaw model, was first used by
H. S. HeleShaw (1897, 1898, 1899) to demonstrate two
dimensional potential flow of fluid around a ship's hull
and other variously shaped objects. The analog is based
on the similarity of the differential equations which
describe twodimensional laminar flow, or potential flow
for that matter, of a viscous fluid between two closely
spaced parallel plates; and those equations which describe
the field of flow below the phreatic surface of ground
water, namely Darcy's law:
q= 2 ; qz K (3.1la;b)
qx x ax z z 3z
where:
qx, q = Darcy velocity of specific discharge in the
xdirection and zdirection, respectively.
K K = hydraulic conductivity in the xdirection
S and zdirection, respectively.
x = horizontal direction (major flow direction).
z = vertical direction.
= potentiometric head
and, by use of the conservation of mass principle, the
Laplace equation:
2 + a = 0 (3.2)
3x2 z2
Viscous Flow Analog
To demonstrate the analogy of model and prototype,
the equations of motion and continuity for laminar flow of
a viscous fluid between two closely spaced parallel plates
will be developed and then compared to equations 3.1a, 3.1b
and 3.2.
CHAPTER III
THE HELESHAW MODEL
The viscous flow analog, more commonly referred to
as a parallelplate or HeleShaw model, was first used by
H. S. HeleShaw (1897, 1898, 1899) to demonstrate two
dimensional potential flow of fluid around a ship's hull
and other variously shaped objects. The analog is based
on the similarity of the differential equations which
describe twodimensional laminar flow, or potential flow
for that matter, of a viscous fluid between two closely
spaced parallel plates; and those equations which describe
the field of flow below the phreatic surface of ground
water, namely Darcy's law:
q K ; q =K (3.1a;b)
where:
qx, qz = Darcy velocity of specific discharge in the
xdirection and zdirection, respectively.
Kx, Kz = hydraulic conductivity in the xdirection
and zdirection, respectively.
x = horizontal direction (major flow direction).
z = vertical direction.
= potentiometric head
and, by use of the conservation of mass principle, the
Laplace equation:
1 + l = 0 (3.2)
[ x2 z2
Viscous Flow Analog
To demonstrate the analogy of model and prototype,
the equations of motion and continuity for laminar flow of
a viscous fluid between two closely spaced parallel plates
will be developed and then compared to equations 3.1a, 3.1b
and 3.2.
Consider a viscous incompressible fluid flowing
ever so slowly between two parallel plates which are spaced
such that the Reynolds' number, R based on the interspace
width is less than 500 (Aravin and Numerov, 1965). In
Cartesian coordinates, the general NavierStokes equations,
i.e., the equations of motion, are:
DV
S= B + 1 { + 2V (3.3a)
Dt x p 5x x
DV
= B + 1 {+ V2V (3.3b)
Dt y p ay y
DV B + 1 k + pV2V (3.3c)
Dt z p z
where:
substantial derivative = + V + V
Dt at x 5 y 
+ V
z 9z
V2 = Laplace operator = + +
3x2 ay2 az2
V V V = velocity in the x, y and zdirections,
respectively.
B x B B = body forces in the x, y and
S zdirections respectively.
p = pressure.
p = density of the fluid.
p = absolute viscosity of the fluid.
y = horizontal direction (minor flow
direction).
t = time.
Referring to the freebody diagram of the idealized
flow regime shown in Figure 3.1, if no slip conditions
(adherence to the walls) of the fluid particles are assumed
in the molecules closest to the walls of the parallel plates,
it is easily seen that the velocity gradient in the y
direction is much larger than the velocity gradient in either
the x or zdirections. Thus, the first and second order
partial derivatives taken with respect to both x and z may
be neglected when compared to those taken in the ydirection.
Secondly, because of the very low velocities ("creeping"
motion) the inertia terms, that is the terms on the left side
of equations 3.3, are very small when compared to the viscous
_a__
1"~n """""""~""~ ~ ~ L~""
x I xm
pVz+ (pV
pV x, 
pV
___ z
FIGURE 3.1
FREE BODY FLOW DIAGRAM
HELESHAW MODEL
19
terms, those on the right side of equations 3.3, and may
be neglected. Thirdly, because of the restriction to two
dimensions, the velocity in the ydirection is taken to
be zero; consequently, all rates of change of velocity in
the ydirection must be zero. Finally, the only non
cancelable body force acting on the fluid is gravity which
acts only in the vertical. Mathematically, Bx = (gz)
= 0; = (gz) = 0; and Bz (gz) g = 32.17
ft./sec.2. Incorporating all of the above arguments and
values into equations 3.3, the equation of motion becomes:
32V
+ V x = 0 (3.4a)
2x ~y2
a = 0 (3.4b)
a2V
pg + = 0 (3.4c)
2 y2
Defining the potentiometric head, or potential, ( = z + p/y,
where y equals the unit weight of the fluid, and taking the
partial derivative with respect to x, y and z, the following
results are obtained after multiplying through by the unit
weight of the fluid:
Y = a = y; 1} = (3.5a;b;c)
ax ax ay ay az 3Z
Introducing these relationships into equations 3.4 and
dividing through by the unit weight yields:
S= P 2V (3.6a)
ax 7 ay2
= 0 (3.6b)
a2V
3b= z (3.6c)
r Y ay2
It is evident from equation 3.6b that the potentiometric
head is constant in the ydirection. It is possible
then to integrate the first and third equations of
equations 3.6 with respect to y. After separating
variables and integrating once, using the boundary condition
3V aV
y= 0, ax = 0 and z = 0, the following equations are
ay dy
obtained:
aV
y 1 x (3.7a)
Sx Y ay
av
y z (3.7b)
3z y ay
Integrating, once again, using the second boundary condition
y = b/2, Vx = 0 and Vz = 0 (no slip) the above becomes,
after solving for the respective velocities:
= {y2 (3.8a)
x 21 T x
V = {y2 b (3.8b)
z 211 4y
Note, where b is the spacing between the plates, that if
a potential = {y } is defined, equations 3.8
can be written:
V = V (3.9a;b)
x ax z z
4 is what Shames (1962) calls a velocity potential. It is
dependent only on y. Integrating the velocity profiles
established by equations 3.8 between the limits of + b/2,
and dividing by b, the directional specific discharges are
obtained:
+b/2 +b/2
Sf dy 1 y 9 y b2y b2 1
x b x b V2p I x 3 4 a12 ? x
b/2 b/2
(3.10a)
+b/2 +b/2 2
q f Vzdy = Y3 b2 Yi
b b 2 z 4 z
b/2 b/2
(3.10b)
r~
21
It is obvious that for a model of constant spacing b, the
b2
quantity 2 does not vary in either the x or zdirection.
Defining the model hydraulic conductivity as K = K =
2 Y xm zm
2 equations 3.10a and 3.10b become:
q = K ; q = K (3.11a;b)
xm Dx z zm 3z
which, of course, is analogous to equations 3.1.
Consider, now, the twodimensional continuity equa
tion for flow between parallel plates:
x + 0 (3.12)
The specific discharge, or Darcy velocity, is re
lated to the velocity by the vector equation n q = V,
where ne is the effective porosity of the flow media. In
the model, ne equals 1. From the analogy Vx = qx; Vz = qz
and substituting the relationship obtained from equation
3.11 into equation 3.12:
mK x} K = 0 (3.13)
or dividing by Km and recalling that for a model Km
K =K :
xm zm
i2_ + a2l = 0 (3.14)
Dx2 Dz2
which is clearly analogous to equation 3.2.
The similarity of equations 3.1 and 3.11 and equa
tions 3.2 and 3.14 establish the analogy.
Scaling
The twodimensional equation along the free surface,
or water table, of an anisotropic porous media given by
Bear (1972) is:
K 2 } + K ({pp = nep (3.15)
xp x zp z z epat
P P P P
where the subscript p denotes the prototype. For a Hele
Shaw model, using the subscript m, the equation can be
written as:
K m! 2
Kxm {x2
m
In
zm m 2
m
m nm tm
zM em t
m m
Introducing the similitude ratios, denoted by the
subscript r, of the corresponding parameter of model and
prototype:
K
S xm
K 
xr K
xp
x
m
x = m ;
r x
p
P
p
K
K =zm
K z
zr K
zp
z
z m
r z
P
n
n em
er n
ep
(3.17a;b)
(3.17c;d)
(3.17e;f)
(3.17g)
t
t = m
r t
P
and substituting these relationships into equation 3.15,
the following is obtained:
RK (x_ /x) K (7z T/r) (7 ) n
Kxr m r zr /r m r em
xr m r zr m r m r) mer
(tm/tr)
(tm/tr)
(3.18)
The ratios of model to prototype quantities are constant
and can be removed from behind the differential; therefore,
equation 3.18 can be rearranged in the following way:
r m 2K
Kxr2 K fax m
Kxr r m
Zr 2 r 2m 2
zm Izr r 2m
Zr (mI tr n m
lzr r m n nem tm
(3.16)
(3.19)
P
23
CQmparing the equations 3.16 and 3.19, it is evident that, if
the equations are identical, the following must be true:
x 2 z2 z t
1 r r r r (3.20)
K K ~ K n (3.20)
xr r zr r zr r er r
Solving the third equality for zr, the following important
relationship is found:
Zr = r (3.21)
The second equality, after crossmultiplying, yields:
x 2 K
z K
{} : =xr (3.22)
r zr
Recalling the definitions of K and K (equations 3.17a;b)
xr zr
and remembering that K = K in an isotropic model, the
xm zm
above equation can be rewritten:
x 2 Kxm/K K
{r} r K/xp = z (3.23)
K /K K
r zm zp xp
The ratio of Zp is called the ratio or degree of anisotropy
K
xp
of the prototype.
Time
Using the fourth equality of equation 3.20, the time
ratio of the model and prototype is established:
n z n x
S er r er r (3.24a;b)
r KR tr K
zr xr r
Substituting the vertical ratio, zr, for the potentiometric
head established by equation 3.20 and the similitude ratios
of time, hydraulic conductivity and porosity into equation
3.24b results in:
t n K x
m em E xp r (3.25)
t n K z
p ep xm r
The effective porosity of an isotropic model, nem, is unity.
The hydraulic conductivity of the model was defined pre
viously as 2, thus the time scale for the model is finally
written as:
K x2
t 12 u xp _r t (3.26)
m g nep b2 zr P
g ep b2 rr
Agisotropy
The HeleShaw model is normally isotropic. This is
because of the nonvariance of the spacing of the parallel
plates. There are, however, two methods for simulating
anisotropy in a model. Equation 3.23 gives a clue as to
the first possibility of simulating anisotropy:
x 2 K /K K
{,} xm xp z= (3.23)
2 K /K K
r zm zp xp
Since Kxm = K the x or z ratio can be adjusted so that
ie model's hydraulic conductivities are kept equal. This
is usually done by choosing a suitable horizontal ratio.
Knowing the prototype parameters, a vertical scale for the
model is computed so that the aforementioned conductivities
Share kept equal, demonstrating:
z K 1/2 x
z _X (3.27)
r z xK 
p zp p
Solving for z :
K 1/2 X
z = { 12 xm z (3.28)
m K x p
zp p
Unfortunately, the geometric distortion method is adequate
for modeling only one ratio of anisotropy. If there is a
second aquifer, within the prototype which has a different
vertical or horizontal hydraulic conductivity, the second
aquifer cannot be correctly simulated; unless, of course,
the second aquifer's ratio of anisotropy is the same as
the ratio of the first. This restriction would severely
limit the use of the HeleShaw analog in modeling of re
gional ground water problems unless another method were
available to correct the ratio of anisotropy.
PolubarinovaKochina (1962) suggests using a grooved
plate within the model to correct the ratio of anisotropy
of the second flow zone. The plate may be grooved in any
PI__
one of several methods. It matters little whether a grooved
plate is sandwiched between the parallel plates, or if rec
tangular bars are attached to the front or back plate. The
degree of anisotropy of the second aquifer and the amount
of geometric distortion used to model the first flow zone
determines the directions the grooves, or bars, are placed;
however, the grooves are normally placed horizontally or
vertically. Collins and Gelhar (1970) have developed the
conductivity equations for the flow zone in which Polubarinova
Kochina's grooved plate is used. The analysis assumes one
dimensional flow and can be used equally well with either
vertical or horizontal orientation of the grooves.
Following Collins and Gelhar (1970), consider flow in
a grooved portion of a model. For simplicity, assume Figure
3.2 is a plan section of the grooved zone. Assuming such,
the horizontal direction then corresponds to the xdirection
and the grooves, which are vertical, lie in the zdirection.
Area 1 is associated with the wider spacing of length ab.
Area 2 is associated with the narrower spacing of length b.
Since flow area 1 is the much larger of the two areas, most
of the frictional head loss occurring through the total
length L is developed in flow area 2 which has length (1X)L.
Lambda, X, is a length correction factor. Referring to
Figure 3.3, the potentiometric gradient across area 2 is:
2 2
2 (3.29)
3x (lX)L
For high values of a:
3 = "2 (3.30)
dx L L
Ac2 3 22
but, from equation 2.29, = (1X) so that:
Lx
S(l) (3.31)
7x ( X)
Applying Darcy's law to area 2:
a2
q K 2 (3.32)
x x, Tx
and, substituting the previous expression for :
K
q x2 3_ (3.33)
x T(1A 3x
Rctancguiar
SBar  7 Back Plate
_L___IL 1 XI)L Front Plate
FIGURE 3.2
PLAN SECTION OF ANISOTROPIC GROOVED ZONE
IN HELESHAW MODEL
x or z
FIGURE 3.3
HEADLOSS IN GROOVED ANISOTROPIC ZONE
r
27
The effective hydraulic conductivity in the xdirection
then is:
K x 2 b2 (3.34)
xm (1X) 12u(1X)
Consider vertical flow through the grooved zone
illustrated in Figure 3.2. In particular, consider flow
downward through areas 1 and 2. The total discharge
through these areas can be written as the sum of the dis
charge through each, that is, Qzm = Q + Q2. Applying
Darcy's law for the total discharge, Q:
Q = K (atbL) (3.35)
Q1 = Kzr az
Q = K (1X)Lb (3.36)
2 Z2 3z
Adding Q1 and Q2:
Q = {K (aX) + K (lX)} bL (3.37)
zm zi Z2 az
For flow area 2 it is not unreasonable to assume that the
frictional forces in the fluid boundary to either side of
area 2 are negligible. Therefore, the horizontal conduc
tivity in this area is the same as defined by the earlier
analysis, that is:
K 2 1 = b2 (3.38)
z2 12 p 12 U
where u is the kinematic viscosity. Furthermore, if
b/ab << 1, the flow in area 1 can be assumed roughly equi
valent to flow through a rectangular hole. According to
Rouse (1959), the equation of motion through a rectangular
crosssection of length XL and width ab is given by:
Vz= x(x XL) + sin x
z 2p dz j=l
{A cosh i + Bj sinh '} (3.39)
where:
A _2y(AL)2 ( j (3.40)
S2y(X) (cos jrl) (3.40)
A3 j3ff3
28
and,
cosh jb 1
B.=A. j (3.41)
R3 3 Isinh XLb
XL
By integrating V over the area and dividing by the total
z
area XLab, the mean velocity is given by:
= (ab)2 (3.42)
z p 12 az
where:
1216 fab2 { {cos(jrI) 1}2tanh
=1Tj 4j
(3.43)
Since the terms of the infinite series decrease as j5,
only the first term of the series need be considered and
retained, so that:
192 ab iXL
s =1 {2} tanh 2b (3.44)
From equation 3.42, the equivalent hydraulic conductivity
in area 1 is given by:
K = j (Cb)2 45)
zi ( 12 3.45)
Finally, introducing the values found for K and K
into equation 3.37: Z Z2
Q  (b (aX) + b (1X)} bL (3.46)
zm 1 1212
or,
QO  b ('XE + 1 X)} bL (3.47)
zm p 12 az
from which it's seen that the effective vertical hydraulic
conductivity is given by:
S+ 1 ) (3b248)
zm = Ty (alx + 1 X) (3.48)
29
after defining Qzm/bL = Vz, where V = q is the effective
vertical specific discharge. Equations 3.34 and 3.48 give
the second method available to correct the hydraulic con
ductivity of a model so that it can simulate the true
ratios of anisotropy found in the prototype.
Leakage
An aquiclude can be defined as a soil stratifica
tion in which the hydraulic conductivities are zero. In
certain geohydrologic problems, it is convenient to assume
such conditions. However, in reality few soil masses are
truly impervious. The degree of perviousness in a stratum
is referred to as leakance and it is generally assumed that
the direction of flow is only vertical. There is no hori
zontal flow, that is,
K = 0 (3.49)
xp
Bear, et al., (1968) suggests the use of vertical
slots to model such a semipervious layer. To accomplish
this, the spacing between the parallel plates of the Hele
Shaw analog is filled with a slotted middle plate. See
Figure 3.4.
The analysis to determine the effective vertical
hydraulic conductivity of a model's leaky layer closely
parallels that for an anisotropic grooved zone. Again,
following Collins (1970), Darcy's law for flow through a
vertical slot is:
Q = K XLab (3.50)
z z K z
The effective specific discharge through the slot found by
integrating the Rouse equation (equation 3.39) is the same
as equation 3.42 from which is found the hydraulic conduc
tivity:
K = (ab)2 (3.51)
z p 12
and, introducing the above into equation 3.50:
Q = (b)2 SLab (3.52)
z i 12 dz
Again, the effective specific discharge, or mean velocity,
is equal to:
Sz 5 b (3.53)
z bL 12
Back Plate
Flowspace
Solid ab
I I
I I
Front Plate
FIGURE 3.4
PLAN SECTION OF LEAKY ZONE
IN HELESHAW MODEL
bm
Eli
z^
_zz
Area, Am
 lm m
END SECTION
ELEVATION
FIGURE 3.5
STORATIVITY
r
Front
Plate
Back
Plate

i
1
so that the effective hydraulic conductivity of a leaky
layer in the model is:
K a3~X b (3.54)
zm P 12
Storativity
While the problem of storage has not been completely
solved, it has, in general, been neglected by most researchers.
Bear (1960) suggests that discrete tubes attached to either
the front or back plate and connected to the aquifer be used
to model the specific storage of a confined aquifer. For a
nonisotropic aquifer, the righthand side of equation 3.13
is not zero, but, in fact, equals the specific storage, So,
times the rate of change of the potentiometric head. Rewrit
ing the twodimensional equation 3.13 for both model and
prototype to include the above is as follows:
K  + K  = S 3 (3.55)
xp ax2 zp 9z2 op at
P P
a2,m 2+m a m
K 2 + K = S o (3.56)
xm ax2 zm 8z2 m
m m
Defining a ratio of storativity:
S
S = om (3.57)
or S
op
it follows from inspection that,
z K z
K r_ = zr = S r (3.58)
xr 2 z or t
x r r
or that,
K t S
S zr r om (3.59)
or 2 S
z op
r
Referring to Figure 3.5, the storage represented by
the model in the discrete length 1m is equal to:
A
m
Som = b 1 (3.60)
om b m
mmm
where A is the crosssectional area of the storativity
m
tube. Introducing the above into equation 3.59 and solving
for A :
m
t
A =blz S K r
m mmm op zr 2 (3.61)
Discharge
The discharge scales are obtained from Darcy's law.
Written for both prototype and model with the usual sub
scripting, these are in the xdirection:
Q K b z (3.62)
xp = axp x pp
P
and,
am
Qxm= K bm z (3.63)
xm xm ax mm
m
Dividing equation 3.63 by equation 3.62 and recalling the
definitions for the various parameters'ratios, it follows
that:
z z
Q = Kxr r br = K b r (3.64)
r r
Similarly, in the zdirection,
x
Q r K b K b x (3.65)
zr zr r r z zr r r
r
Solving equation 3.22 for the hydraulic conductivity in the
xdirection:
x 2
Kx =K zr } (3.66)
xr zr zr
and, substituting this result into equation 3.64, it follows
that:
X 2
Q =K {} b = K b x (3.67)
xr zr z r x zr r r
r r
or,
Qx r Q = Q (3.68)
Accretion
Accretion, R, is the rate at which a net quantity
(percipitation and surface inflow minus evapotranspiration,
runoff, etc.) of liquid is taken into the flow system at
the phreatic surface. It is measured as a volume per unit
horizontal area per unit time, that is:
Qr
Rr = r (3.69)
r r
From equations 3.64 or 3.65, it follows that:
K z 2
R = xr r =K (3.70)
r 2 zr
x
r
Volume
On occasion, volume, U, is of some importance. The
volume scale follows directly from continuity, that is:
Ur = Qtr (3.71)
Substituting the values found from equations 3.65 and 3.24a
for Q and tr respectively, the above equation becomes:
z
r b n x z (3.72)
zr
Ur = b x n er K = n r r z(3.72)
As inferred earlier by this section's opening sentence, the
volume scale is usually neglected; however, in the case of
free surface water bodies, lakes, rivers, etc., if the vol
ume exchange of liquid is of interest and has to be modeled,
the volume scale requires an additional restriction. In
the following analysis, the bar above the width dimension
indicates the free water surface of a river, lake, ocean,
or such.
In the portion of the model simulating the body of
water, the spacing of the model is increased to maintain
hydrostatic pressure distributions within the model. The
narrower spacing of the model is, of course, a measure of
the hydraulic conductivity of the aquifer. In the proto
type, however, the width of the open water and the aquifer
are equal and this leads to the following (Bear, 1960) for
the model and prototype, respectively:
U = n b x z (3.73a)
r er r r r
U = n b x z (3.73b)
r er r r r
The same volume ratio must be applicableto both the narrow
and the enlarged interspace; therefore, U = U It fol
lows that:
n = n b (3.74)
er r er r
but,
n
n em = 1 (3.75)
er
ep
so,
S ner b (3.76)
Note that for an anisotropic media, nem does not necessarily
equal one.
CHAPTER IV
.SITE SELECTION AND PROTOTYPE GEOLOGY AND HYDROLOGY
Site Selection
The site selected for this study is the middle Gulf
area of Florida. This region has a rapidly expanding popu
lation with a corresponding growth in water demand. The in
creased pumping to satisfy this demand also increases the
likelihood of saltwater intrusion, and, in fact, a number
of municipal supply wells in the coastal zone have been shut
down in recent years due to chloride contamination.
One of the major water supply systems in this region
is the Pinellas County Water System, and this study is cen
tered around the EldridgeWilde well field of this system.
The location of EldridgeWilde in relation to several of the
population centers of this region is shown in Figure 4.1.
It is about 8 miles east of the Gulf of Mexico and encom
passes an area in the northeast corner of Pinellas County,
at the intersection of the boundaries of Pinellas, Hills
borough and Pasco Counties.
In 1970 (Black, Crow and Eidsness, Inc., 1970), the
waterworks facilities at EldridgeWilde included: sixtyone
water wells, over 11 miles of rawwater collection piping,
water treatment facilities consisting of aeration and chemi
cal treatment, including chlorination and fluoridation, and
high service pumping units.
All wells are open hole and penetrate the Floridan
aquifer at depths from 140 to 809 feet below ground surface,
averaging 354 feet. The design capacity of the field at the
present i c Q mi n gallons per day, although the
maximum allowable oumacre has been set by the Southwest
Florida War Manaement District at 28 million gallons per
day on the average with a maximum day of 44 million gallons
per day.
In selecting the prototype location within the site
area, two characteristics of the vertical HeleShaw analog
must be considered. The first characteristic is that there
can be no general flow normal to the parallel walls of the
model. This means that the flow from one end of the model
to the other is streamline flow. The second characteristic
is that the ends of the model are finite.. Therefore, the
prototype must be along a streamline in the flow domain and
have boundary conditions which are "infinite" reservoirs or
water divides.
35
i ri
The prototype selected meets the above requirements
and includes the point of interest, i.e., EldridgeWilde
well field. The center line of the prototype is shown in
Figure 4.1 as the unbroken line passing through Eldridge
Wilde in a southwest to northeast direction. The dotted
contours in the figure define the potentiometric surface of
the Floridan aquifer in feet above mean sea level as of
May, 1971. They were obtained from a map publication en
titled "Potentiometric Surface of Floridan Aquifer South
west Florida Water Management District, May, 1971" prepared
by the U. S. Geological Survey in cooperation with the
Southwest Florida Water Management District and the Bureau
of Geology, Florida Department of Natural Resources. Now,
in a flow field, streamlines are perpendicular to potentio
metric lines. As can be seen from the figure, the proto
type orientation reasonably satisfies the streamline re
quirement. The prototype is terminated on the southwestern
end at the 15 feet depth contour in the Gulf of Mexico, and
it is assumed that this satisfies the infinite reservoir
boundary condition. The northeastern terminus is located
in the center of the 80 feet contour, southwest of Dade City.
This location satisfies the water divide boundary condition.
The area in the vicinity of the 80 feet contour is known as
the Pasco High. The overall length of the prototype is 36
statute miles. The width of the prototype is taken to be
3.5 statute miles. This dimension is sufficient to include
the cone of depression caused by pumping in EldridgeWilde
well field, and is based on the results of a study by Mr.
Evans (employing a numerical model) for Black, Crow and
Eidsness, Inc. The land surface contours of the prototype
are obtained from U. S. Geological Survey topographic maps.
The bottom boundary of the prototype is taken to be the
base of the Lake City Formation, with depths being deter
mined from available well logs of wells in the prototype
vicinity. The maximum depth from highest land surface to
deepest point is 1340 feet.
Prototype Geology
Stewart (1968) identifies eight formations as being
of interest in terms of water production in the prototype
area. They are in descending order, the Undifferentiated
Deposits, Tampa Limestone, Suwannee Limestone, Crystal River
Formation, Williston Formation, Inglis Formation, Avon Park
Limestone and Lake City Limestone. Underlying the Lake City
Limestone is the Oldsmar Limestone which is not used as a
source of water at present.
The Undifferentiated Deposits are interbedded sand,
silt and clay of PostMiocene age and range in thickness
from zero near the Pasco High to 60 feet in the Eldridge
Wilde well field. The thickest deposits are in northeast
Pinellas County around the north end of Lake Tarpon where
sand dunes, as much as 40 feet high, overlie alternating
layers of clay, thin limestone beds and sand greater than
70 feet thick.
The Tampa Limestone is a hard, dense, sandy, white
to lighttan or yellowishtan fossiliferous limestone of
Miocene age. This limestone is near the surface in the
area of the Pasco High and about 80 feet below land surface
at the EldridgeWilde well field. At EldridgeWilde, the
thickness varies erractically from about 20 to 240 feet.
The Tampa Limestone is a poor to fair producer of water.
The Suwannee Limestone is a soft to hard, nodular
or granular, fossiliferous white to tan limestone of Oligo
cene age and is about 200 feet thick. The Suwannee and
Tampa Limestones are the major water producers for wells
in the area.
The Crystal River, Williston and Inglis Formations
comprise the Ocala Group of late Eocene age. The Crystal
River and Williston Formations are lithologically similar
units of white to cream, porous, soft, coquinoid limestone
and are generally poor producers of water. The Inglis
Formation is a hard, cream to brown to gray fossiliferous
limestone and is generally a good producer of water.
The Avon Park and Lake City Limestones are litho
logically similar units of soft to hard, cream to brown,
fossiliferous limestone with beds of dolomitic limestone
and some gypsum. Both formations are good producers of
poor quality water.
The Oldsmar Limestone is a fragmental dolomitic
limestone with lenses of chert, thin shale beds and some
gypsum.
In this report, two formations are considered, the
Undifferentiated Deposits and the Floridan aquifer. The
Floridan aquifer is considered to contain all formations
from the Tampa to, and including, the Lake City Limestone.
The transmissivitv of the Floridan aauifer ranges
from about 165.Onn .n 55n.000 gallons per day per foot,
and the coefficient of storage ranges from about 0.0005 to
0.0015. The coefficient of leakage is approximately 0.0015
gallons per day per cubic foot.
Based on ground water discharge and water levels,
the estimated recharge (leakance) to the Floridan aquifer
was computed s.9tfwart 19 9 o 1 F bTvhnll* 1 m 11nn gaI lnns
per day. Based on aquifer test data, the estimated recharge
for a 250 square mile area was 90 million gallons per da.
The Undifferentiated Deposits act as a confining
layer, and the Floridan aquifer is thus under artesian
conditions.
Prototype Hydrology
The surface waters of the area consist of many
lakes and few streams. Because of the flat topography,
little water runs off into streams, and swampy wetlands
are numerous. Most rainfall evaporates or is transpired
by plants.
The Floridan aquifer is recharged through the
Undifferentiated Deposits by surface and ground water
derived from local rainfall. Many millions of gallons of
water are also admitted to the aquifer by numerous sink
holes in the region. Water levels in the Floridan aquifer
respond to rainfall since this is the recharge source.
This response is not immediate, but usually fluctuates with
the wet and dry seasons. Water levels in wells which are
not directly affected by local pumping show yearly lows in
the dry season, April and May, and yearly highs during the
wet season, late summer or early fall, (Black, Crow and
Eidsness, Inc., 1970).
The aquifer recharge has been estimated (Black,
Crow and Eidsness, Inc., 1970) from available data and the
use of the following formula:
Aquifer Recharge = P + SWI + GWI ET R GWO
The basin area is 575 square miles, changes in
storage are assumed zero, and evapotranspiration is assumed
to be 75% of the precipitation. The applicable values are
listed below in million gallons per day:
P = Precipitation = + 1492
SWI = Surface Water Inflow = + 0
GWI = Ground Water Inflow = + 0
ET = Evapotranspiration = 1119
R = Runoff = 218
GWO = Ground Water Outflow = 37
Aquifer Recharge = 118
This value is in reasonable agreement with previous
reported values (Stewart, 1968).
CHAPTER V
DESIGN, CONSTRUCTION AND OPERATION OF MODEL
Design
Prototype
The selection of the prototype area was discussed
in Chapter IV. Table 5.1 is a summary of the prototype
characteristics. The leaky layer is synonymous with the
undifferentiated deposits. The top and bottom of the
Floridan aquifer were determined by straightline extra
polation from available well logs.
The hydraulic data are within the reported range
of values and are the result of a trial and error process
to stay within the range and still produce a reasonable
model.
Model
The purpose of the HeleShaw analog in this study
is to model saltwater intrusion. Before discussing the
model design, it seems appropriate at this point to pro
vide some background about saltwater intrusion. Water,
in general, whether it be surface water or ground water,
is continually migrating towards the sea, where an equi
librium, or moving freshwater/saltwater interface, is
established. The two fluids are messible, but because of
the difference in densities and the very low velocities,
the interface is formed. Across the interface, the salin
ity varies from that of the fresh ground water to that of
the ocean. The transition zone, as it is called, is due
to hydrodynamic dispersion and, although it is anything
but abrupt, it is usually assumed to be. The interface
then is generally selected to occur at some measured
electric conductivity or salt (chloride) concentration.
The earliest investigations of saltwater encroach
ment were made by BadonGhyben (1888) in Holland and
Herzberg (1901) in Germany. Working independently, both
investigated the equilibrium relationships between the
shape and position of the freshwater/saltwater interface.
Figure 5.1 shows a coastal phreatic aquifer and the Ghyben
Herzberg interface model. BadonGhyben and Herzberg
assumed static equilibrium and a hydrostatic pressure dis
tribution in the fresh ground water and stationary saline
ground water near the interface.
TABLE 5.1
PROTOTYPE PARAMETERS
Parameters
a. Geometric
x (ft.)
Zp (ft.)
b (ft.)
P
Floridan Aquifer
190,080
1,340
18,480
Leaky Layer
190,080
55
18,480
b. Hydraulic
Txp (gpd/ft.)
Tzp (gpd/ft.)
K
xp
K
zp
Leakance (gpd/ft.3)
Kzp (gpd/ft.2)
Kxp (gpd/ft.2)
S
vp @ 770 F (ft.2/sec.)
g = 32.2 ft./sec.2
161.4
196.7
0.00158
0.965 (105
0.965 (10
0.0015
0.09
0
0.965 (105
0.965 (10
225,000
184,426
1.218
FIGURE 5.1
GHYBENHERZBERG INTERFACE MODEL
Considering a point P on the interface, and
choosing mean sea level as the datum, the pressure at
point P is:
P = h y (5.1)
p ss
where:
h = vertical distance from mean sea level to
s point P.
s = unit weight of sea water.
This pressure may, also, be expressed by:
Pp = (hf + hs) Yf (5.2)
in which hf equals the vertical distance from mean sea
level to the phreatic line at the location of P and Yf
equals the unit weight of freshwater. Equating the
preceding two equations:
hsys = hfYf + hsYf (5.3)
and, rearranging,
h Ys hsYf = hfYf (5.4)
Solving for h :
h = h (5.5)
s Y f f
Introducing y = pg, where p is the density, factoring and
canceling out the gravity term, the GhybenHerzberg rela
tion is found:
Pf
h p= hf (5.6)
s p P f
s f
for a saltwater density of 1.025 lb sec.2/ft." and a
freshwater density of 1.000 lb sec.2/ft.4, the quantity
Pf/(Ps Pf) = 40. The implications of equation 5.6 are
rather dramatic. For instance, for every foot of fresh
water above the tum th fet of freshwater
below the datum. More importantly however, consider the
effects of lowering the phreatic surface. For every one
foot drop of the water table, the interface raises 40
teeth. it must be remembered that the above analysis
assumes static conditions. This, in fact, is not always
the case. The position of the interface is a function of
dynamic conditions rather than static. Even so, in cases
where flow is quasihorizontal, i.e., the equipotential
lines are nearly vertical, equation 5.6 is valid.
Many investigators have incorporated dynamic forces
into the analysis of the stationary interface. Hubbert
(1940) was able to ascertain a more accurate determination
of the shape of the interface near the coast line. He
assumed that at the interface the tangential velocity was
zero in the saltwater, but increases with horizontal dis
tance in the freshwater as the coast line is approached.
This then is the cause for the interface to tilt upwards
as the sea is approached and the greater depths found than
those estimated by the GhybenHerzberg relationship.
Hubbert showed that the GhybenHerzberg equation holds be
tween points on the water table and the interface along an
equipotential line, rather than along a vertical plane.
R. E. Glover (1959) modeled an infinitely deep
coastal aquifer by assuming no flow in the saltwater
region, a horizontal water table and a horizontal seepage
face located seaward of the coast line. He found an exact
solution for the shape of the wedge, giving the following
relationship:
z2= 2qx + q2 (5.7)
S I K 2 P SPf 2
K{I s f K2 Sf
Pf Pf
where x and z are the horizontal and vertical directions,
respectively, q is the seepage rate per unit width and K
is the hydraulic conductivity. De Wiest (1962) using
complex variables and a velocity potential of D = K x
(p pf)/Pf derived the same equation.
Bear and Dagan (1964), using the Dupuit assumptions
and the GhybenHerzberg equation, developed the approximate
shape of the interface for a shallow aquifer of constant
depth.
All of the above investigated the equilibrium posi
tion of the saltwater/freshwater interface.
If there is a change in the freshwater flow
regime, a transition period is caused during which the
interface moves to a new point of equilibrium. The non
linear boundary conditions along the interface make the
solution for the shape and position of the transient inter
face all but impossible except for the simpliest geometries.
Bear and Dagan (1964b), as well as other investigators,
have used the Dupuit assumptions to approximate the rate
of movement of an interface in a confined aquifer. Fol
lowing PolubarinovaKochina's (1962) suggestion, they
assumed quasisteady flow and were able to approximate
the interface shape and position for both a receding
motion and landward motion of the interface.
Characteristically, the solutions obtained by
investigators to date have all had simple geometries and
involved simplifying assumptions, some of which have had
little resemblance to actual conditions. Therein lies
the advantages of a HeleShaw analog, complex geometries
and boundary conditions can be modeled with relative ease.
In order to satisfy additional similitude require
ments for the flow of two liquids with an abrupt interface,
as in this study, and to provide a suitable time ratio,
two liquid silicone fluids were chosen to be used in the
model. Dow Corning Corporation Series 200 silicone fluid
was used to model freshwater. Series 510 silicone fluid
from the same company was used to simulate saltwater.
The 200 Series fluid and the 510 Series fluid have densities
of 0.977 gm/cm3 and 1.00 gm/cm3, respectively, at 250 C.
Both fluids have absolute viscosities of 500 centistokes at
250 C. Dow Corning 200 fluid is a clear dimethyl siloxane
which is characterized by oxidation resistance, a relatively
flat viscositytemperature slope and low vapor pressure.
Dow Corning 510 fluid is a clear phenylmethyl polysiloxane
which also has a relatively flat viscositytemperature slope.
In order to locate and follow the interface movement, the
denser fluid was dyed blue and the lighter fluid dyed orange.
The dimensions of the model were selected so that a
unit of reasonable size would be produced. These dimensions
are 11.75 feet long, 0.5 inch wide inside and 2 feet deep.
The parameters x r z and b are therefore set, and the
r r r
result is a distorted model. This distortion requires the
use of a slatted inner zone as discussed in Chapter III.
Tables 5.2 and 5.3 list the applicable parameters.
The following analysis is shown for the design of the
model Floridan aquifer, leaky layer and storage coefficient:
a) Slatted Anisotropic Zone for the Floridan Aquifer.
1) Compute kxm/kzm from Eq. 3.22
x 2 K
z K
{r = xr (3.22)
r zr
Noting that k = K (5.8)
g
TABLE 5.2
MODEL PARAMETERS
Parameters Floridan Aquifer
a. Geometric
xm (ft.) 11.75
z (ft.) 2.0
b (ft.) 0.0417
a 16.50
x 0.8097
b (in.) 0.030
L (in.) 1.235
E 0.6905
ab (in.) 0.495
XL (in.) 1.00
(1 X)L (in.) 0.235
b. Hydraulic
k (in.2) 3.94 x 10
xm
k (in.2) .1884
zm
k m/k .00209
1
Sm (ft.) 0.0369
3
v @ 770 F (ft.2/sec.) 5.382 x 103
m 32.2 ft./sec
gm = 32.2 ft./sec.2
Leaky Layer
11.75
0.082
0.0417
.25
0.3846
0.500
1.235
0.8349
0.125
0.475
0.760
0
1.05 x 104
0
5.382 x 10
5.382 x 10
TABLE 5.3
SIMILARITY RATIOS
Parameters
a. Geometric
x
r
z
r
b
r
b. Hydraulic
k
xr
k
zr
kzr/kxr
S
or
Vr @ 770 F
R
r
Q
t
r
Floridan Aquifer
s
6.182 x 10
3
1.493 x 10
6
2.255 x 106
3 x 104
1.748 x 107
582.8
2.677 x 104
557.7
3.136 x 104
6
4.372 x 106
6
1.903 x 106
1
Leaky Layer
5
6.182 x 10
3
1.493 x 103
6
2.255 x 106

1.748 x 107
557.7
3.136 x 104
6
4.372 x 10
6
1.903 x 106
1
48
Then (3.22) becomes,
xr 2 kx k
f^.kzp (5.9)
r xp zm
and finally,
kxm k x 2 r 2
xm x r2 6,182 x 10
t = (1.22){682 x l0 = 0.00209
k k z =s
zm zp r 1.493 x 10
This is the value which the slatted zone must produce.
2) Select dimensions for the slatted zone as shown
in Figure 3.2, as follows:
i) hold ab = 0.5 inch
ii) select b
iii) find a
iv) select (1 A)L
v) select XL
vi) find L
vii) find X
b2 1
3) Compute k = 1 ' (3.34)
xm 12 (1 A)
b2
4) Compute kzm = 12 (aS + 1 X) (3.48)
5) Compute kxm/kzm and compare to the results of
step (1).
6) Repeat the process until the result of step (5)
equals the result of step (1).
b) Slotted Leaky Layer. Refer to Figure 3.4 for this
section.
kxlm = 0 (5.10)
b2
kzlm = b2 a3EX (3.54)
Also, continuity of flow between the leaky layer and
the Floridan aquifer requires that:
Qlr = Qfr (5.11)
Also,
r = Qzr = Qxr = Kzrbr (3.67) & (3.68)
Therefore,
K br = K zfrbrx
zlr r r zfr r r
(5.12)
so,
Klr = Kfr (5.13)
or,
klr = kzfr (5.14)
Therefore,
kzlm =kzlp kfr (5.15)
Now rewrite Eq. 3.54 as,
12
zim b2
1) set b = 0.5 inch
2) a = constant, since b is set and kzlm can be
computed from previous information.
3) select ab and find a
4) select AL & L and find X
5) compute
6) compute a 3X
7) compare the result of step (6) to the result of
step (2).
8) repeat the process until the result of step (6)
equals the result of step (2).
9) make sure the physical size of this layer is
compatible to the slatted layer, especially in terms of
slot spacing, i.e., blockage.
c) Storage Coefficient Manifolds. Refer to Figure 3.5
for this section.
1) take average storage coefficient and average depth
in prototype to compute S
op
2) select convenient time ratios, in this case
1 minute = 1 year, and compute Sor from Eq. 3.59.
3) compute S
om
4) Eq. 3.60 is now used to compute Am for various m's
with zm averaged over m.
In this case the model was apportioned into five
zones with one manifold per zone.
Construction
Because the HeleShaw analog is capable of modeling
complex geometries and boundary conditions, it is desirable
that it be as adaptable to as many different prototype
geometries and hydraulic parameters as possible. This
would facilitate model construction and investigations of
many different areas in the state of Florida where salt
water intrusion is, or in the future might be, a problem.
A reduction in cost of investigation would also be achieved
if many of the parts were reusable.
A list of general specifications would then be as
follows:
1) The HeleShaw model should be housed in a frame
in which it can be easily installed and removed.
2) The front and back plates with interior model
parts should not be permanently sealed together.
3) The front and back plates should be as adapt
able as possible to different situations.
4) The model should have as few opaque parts as
possible.
5) The model should be mobile.
Frame
As shown in Figure 5.2, the frame is composed of
two assemblies; the cradle and the cradle dolly. The
function of the cradle is to support and orient the Plexi
glas plates. It also contains the inflatable neoprene
hose which seals the plates. It is fabricated of 21/2"
x 2" x 3/8" steel angles which are welded into a channel
21/2" high by 4" wide. The channel has the shape of an
elongated rectangular "U" which is 321/2" high and 147"
long.
The function of the cradle dolly is to support the
the cradle and provide mobility. The cradle dolly is
constructed of 31/2" x 31/2" x 1/4" steel tubing. The
length of the center tube is 14 feet, the short cross
pieces at the ends are 2 feet long. The cradle load is
transmitted to the cradle dolly through two stub shafts
and pillow blocks which are mounted on pedestals at each
end of the center tube. The center line of the stub
shafts is aligned on the cradle so that the cradle may
be easily rotated as shown in Figure 5.3. This rotational
ability allows easy insertion or removal of the model into
or out of the cradle. It also allows a convenient orienta
tion to be selected for working on the model.
The height of the stub shaft center line from the
floor is 37". When the cradle is rotated into a horizontal
position, its internal supports correspond to the height of
the table upon which the model is built and assembled.
This allows the model to be slid from the table and into
4sbg~~ hV::B
mil
'iii vi
FIGURE 5.2
CRADLE AND CRADLE DOLLY
.
I?;** *, *
the cradle. The cradle is then rotated into the vertical
position and the whole assembly rolled into the laboratory
for testing. Figure 5.4 shows the frame and model in a
completed setup.
Figure 5.5 shows the air hose and valve arrangement
for pressurizing the hose seal around the model.
Figure 5.6 shows the stub shaft and pillow block
mounting arrangement. Each pillow block has a jam screw
which is tightened against the stub shaft to hold the
designed orientation of the cradle. Spreader arms are
added to the cross member at the bottom of the cradle
dolly to provide stability, anchoring and leveling capa
bilities. The spreader arms increase the cross member
length to 48".
Figure 5.7 shows the backup pressurizing system
consisting of SCUBA bottle, regulator and air hose.
Figure 5.8 shows a typical crosssection of the
internal support and sealing system. The Plexiglas plates
are supported on pairs of angle brackets and fastened into
the cradle with cap screws and bolts. There are nine
pairs of angles in the bottom of the cradle and two pairs
in each end. The sides of the cradle are drilled slightly
oversize to permit adjustment in the angle brackets and
Plexiglas plates. The neoprene sealing hose is 1" O.D.
with 1/8" walls. It is contained in an aluminum channel
which provides lateral support and proper alignment with
the Plexiglas plates. A pressure of 25 psig was found
sufficient to swell the hose against the inner edges of
the plates and provide positive sealing. Figure 5.9 shows
a portion of the front plate of the model with support
angles and bolts located in the cradle channel.
Plexiglas Plates and Manifolds
The front and back plates of the model were fabri
cated from 5' x 6' x 1/2" sheets of G grade Plexiglas.
The sheets were sawn in half to obtain two pieces 21/2'
x 6' x 1/2". These pieces were then glued together to
form a single plate 21/2' x 12' x 1/2". The bottom
corners of each such plate were cut to a 4" radius to
accommodate the sealing hose. Figure 5.10 shows, from
left to right, the back plate with internal strips and
the front plate. They are on the construction table men
tioned previously. The plates are reinforced around
their edges with a Plexiglas strip and a plate over the
center butt joint.
Figure 5.11 and Figure 5.12 show the detail of the
model leaky layer and Floridan aquifer. The 0.030" gap
in the Floridan aquifer is maintained by seventynine
stainless steel tabs which are glued to the strips at
regular intervals. In Figure 5.11, (turn the page side
ways), the base of the aquifer is shown as the tapered
FIGURE 5.4
FRAME AND MODEL SETUP
z
O
LAn
0
H
FIGURE 5.6
STUB SHAFT AND PILLOW BLOCK
ARRANGEMENT
j i i!'
FIGURE 5.7
BACKUP AIR SUPPLY
FIGURE 5.8
INTERNAL SUPPORT AND SEALING SYSTEM
Back Plate
3/8 x 24
Socket Head
Cap Screw
I\I
3/8 x 24
Threaded
Rod
6 x 32
Machine Screw
1/4 x 28
Socket Head
Cap Screw
Neoprene Hose
H
mow"
Vre
^&
*
FIGURE 5.10
MODEL BACK AND FRONT PLATES
~b~m,~*u~~~;~~';
' "
:iu'
i:raP. h'' P
~a..
*2'
AMSWAL
INFO"p
$1
9, >
FIGURE 5.12
MODEL BACK PLATE, DETAIL
cl 'e ~ ~~C'~
~ r
i C.~) ~ r
M ~e
FRONT PLATE
FIGURE 5.13
WITH ACCRETION MANIFOLDS
lAt
Ael
r
:1
FIGURE 5.14
BACK AND FRONT PLATE CLAMP UP
o I e
i^^rttlHea
'
"
~,
..
.. "
iPp~n /
FIGURE 5.15
ACCRETION MANIFOLDS
Cr
# f
41
~, i, ia
a a 1
FIGURE 5.16
ACCRETION MANIFOLDS,
DETAIL
Constant Head
Tank
G u Gage
"A" % Accretion
P Adjustable Pressure
ByPass Valve
SShutOff Valve.
Well Pump Accretion Pump
FIGURE 5.17
FLUID SUPPLY NETWORK SCHEMATIC
O
t3
H
H
O0
toc
H
ri
H
FIGURE 5.19
SALTWATER RESERVOIR AND PUMP
FIGURE 5.21
FRESHWATER SUPPLY SYSTEM
FIGURE 5.20
SALTWATER CONSTANT HEAD TANK
4dd=
this tank. Fluid overflowing the center tube of the tank
returns to the reservoir. The tank is connected to the
Gulf of Mexico end of the model by three tubes and valves
as shown in Figure 5.20. The tank's elevation may be
raised or lowered in order to set the proper sea level in
the model.
FreshWater System, General
Figure 5.21 shows the three subsystems of the fresh
water system. At the right top of the figure is the fresh
water reservoir. The accretion pump is at the lower right
of the figure. The well pump is at the lower left, and the
flow meters are mounted on the peg board at the upper left
of the figure.
FreshWater System, Accretion
Figure 5.22 shows the freshwater reservoir and the
accretion pump unit. This unit consists of a gear pump,
a variable speed, reversible motor and an adjustable
pressure bypass valve. The unit supplies fluid at a pre
selected pressure to a manifold mounted on the flow meter
board. This manifold can be seen in Figure 5.23 in the
top center of the figure.
FreshWater System, Wells
Figure 5.23 shows the well pump unit and its supply
manifold. The pump unit consists of a gear pump, a vari
able speed, reversible motor and a sensitive compound gage.
This unit pulls a slight vacuum on the manifold in order
to remove fluid from the model at well locations. Fluid
removed from the model is returned to the freshwater
reservoir. In addition, a line from the reservoir dis
charge is kept open to the pump to prevent operating in a
"starved" condition.
FreshWater System, Flow Meters
Figure 5.24 shows, from left to right, the well
supply manifold and the accretion supply manifold with
their connections to the bottoms of the flow meters. There
are twelve flow meters for each manifold. Figure 5.25 is
the opposite view. The two manifolds are connected by a
tube which is only used during startup or shutdown oper
ations; the tube is most clearly seen in Figure 5.24.
Figure 5.26 is an overall view of the twentyfour
flow meter bank. Meters one through twelve supply accre
tion manifolds and recharge wells; meters thirteen through
twentyfour handle well pumping.
Figure 5.27 is a detail view of several flow meters.
The operation of the meters is based on Poiseuille's rela
tionship for laminar flow in a tube:
Q Ah (5.16)
where d equals tube diameter and Ah equals head loss.
65
FIGURE 5.22
FRESHWATER RESERVOIR AND ACCRETION PUMP
FIGURE 5.23
WELL SUPPLY MANIFOLD AND PUMP
FIGURE 5.24
WELL SUPPLY MANIFOLD AND ACCRETION SUPPLY MANIFOLD
r. s
Ih l'd H ^n?1[f
'1,
gtiiipm 4( tiiiti
Btt'T^T ~sI
o*Il4 : T T
FIGURE 5.25
OPPOSITE VIEW OF FIGURE 5.24
r~,: "s~
fj t
lilt' i
67
FIGURE 5.26
FLOW METER BANK
4 ,
,: ::;, ... "
,. 'i : ,, ,, ,
.v cy:, ':'
.: '.. ." 41 : .. I,'&."
~a~p9
For a given tube and fluid, the flow rate and pressure
drop in the tube are proportional. Therefore, if the
pressure drop in a tube is measured, the flow rate may be
determined. The flow meter consists of a 125/8" length
of 1/8" brass pipe with a nominal I.D. of 0.278". The
pipe is threaded to accept a hose barb at the bottom end
and a regulating needle valve at the top end. Pressure
taps are installed in the pipe and are 7.67" apart. This
spacing is consistent with the maximum flow rate required
of the meter and the pressure sensing capability of the
transducers. The pipe is mounted to the peg board with
two Plexiglas blocks. A brass tee from the needle valve
carries an air bleed at the top and a hose barb for hook
up to the model. Figure 5.28 shows some of the connec
tions to the model. The cluster of tubes in the center
of the model is connected to ports simulating Eldridge
Wilde well field. The longer tubes going over the top of
the model are connected to accretion manifolds. The two
tubes extending out of the figure to the left are for
recharge wells near Lake Tarpon. There is an additional
well field near the Pasco High which is not shown and
which was not hooked up.
The direction of flow for accretion flow meters
is from bottom to top. The directions for well pumping
flow meters is from top to bottom.
Figure 5.29 shows the back of the flow meter
mounting board and the pressure sensing lines from each
meter. These lines originate at the pressure taps on the
meters as shown in Figure 5.27 and terminate at a rotary
valve switching device, seen in the center of Figure 5.29.
A more detailed look at this device, and the two pressure
transducers that are connected to it,is shown in Figure
5.30. The switching device consists of four fluid switch
wafers coaxially mounted in a control unit. Each wafer
(Scanivalve, Inc., Model 60 Fluid Switch Wafer, Scan Co.
#W1260/1P12T) has twelve input ports and one output port.
Therefore, two wafers can switch twelve pairs of input
lines from twelve flow meters into one differential pres
sure transducer. Four wafers and two transducers are re
quired for twentyfour flow meters. The control unit has
twelve clickstop positions, and each position connects
two flow meters and two transducers. The pressure trans
ducers used are Celesco #P7D Differential Pressure Trans
ducer, 1/2 PSID.
The output of the pressure transducers is sensed
and conditioned by a carrier demodulator unit, Celesco
Model CD12, which sends an output signal to a Hewlett
Packard dual channel strip chart recorder. Figure 5.31
shows the carrier demodulator positioned above the
recorder.
69
.
N
~ "N
~;6 I N
. .~h1~
FIGURE 5.30
FLOW METER SWITCHING DEVICE AND PRESSURE TRANSDUCERS
p*
I
FIGURE 5.31
CARRIER DEMODULATOR AND STRIP CHART RECORDER
~':(iY~ar
II~
'p"*~8~
.s~S~s8~BSP~
Rather than measure pressure drop across the flow
meters and then computing the flow rate, the meters were
directly calibrated for deflection in inches on the strip
chart recorder as a function of measured flow rate. The
carrier demodulator and recorder were set to provide zero
reading at zero flow rate and near maximum reading at max
imum flow rate.
Each flow meter was operated at various flow rates
as determined by measurement with a stopwatch and gradu
ated cylinder. The recorder reading and fluid temperature
were noted at each point. All data were corrected to
250 C and fitted to a straightline by the method of least
squares. The correlation coefficient for twentyfour
meters ranged from 0.997 to 1.000.
Operation
The accretion pump is used to fill the accretion
supply manifold, well supply manifold (through crossover
tube), flow meters and pressure sensing lines. Air bleeds
are located throughout the system to remove trapped air.
To fill the model, the saltwater pump is started
and the constant head tank filled and adjusted. The valves
on the accretion flow meters are opened and the valves for
the saltwater opened (slightly). The inflow rates are
adjusted to bring the salt and freshwater levels up to
gether. Once the model is near being full and the inter
face formed, the pressure bypass valve on the accretion
pump, pump RPM and the accretion flow meters are set to
provide the flow rates required for the test under con
sideration. The saltwater valves are fully opened. The
well pump is adjusted by seeing that the well supply mani
fold is closed to the pump and that the pump RPM is suffi
cient to maintain a vacuum of = 4" Hg. when recirculating
fluid from the freshwater reservoir. The valve to the
well supply manifold is then opened and the well flow
meters adjusted. The fluid temperature is monitored
throughout the test so that flow rates may be corrected.
Excess accretion is removed from the model through various
drain ports.
When the model is shut down between tests, it may
be drained by gravity or reversing the pumps, or the two
fluids may be left in the model where they become horizon
tally stratified. Starting up again after either proce
dure is a tedious process and the selection of shut down
is left to the subjective evaluation of the operator.
CHAPTER VI
RESULTS, CONCLUSIONS AND RECOMMENDATIONS
Results
The objectives of this project were to (1) develop,
(2) construct, and (3) calibrate a HeleShaw model for a
section of Pinellas County, Florida. Objectives (1) and
(2) have been met. Objective (3) has been met from a qual
itative point of view, but operational features in the
model, in the area of the Pasco High, have presented a prob
lem in obtaining quantitative calibrations. These features
are presently being modified to overcome the difficulties;
this is discussed later, and the Hydraulic Laboratory will
then proceed to a quantitative examination of saltwater
intrusion in this area.
Figures 6.1 through 6.13 show the freshwater/salt
water interface as it proceeds from an equilibrium position
through successive stages of penetration and upcoming due
to pumpkin n EldridqeWilde well field, and then through
further stages as freshwater recharging in the coastal
zone is carried ou.
Figure 6.1 shows the position of the wedge at a steady
state condition. The wedge is correctly located at the coast
and in the approximate present location. EldridgeWilde is
the cluster of dots in the center of the picture, it is
directly above the landward toe of the wedge. The wedge and
the Gulf of Mexico are the dark region at the right end of
the model. The freshwater is the light region.
The time ratio is approximately 1 minute of model time
= 1 year of prototype time. The sequence beginning with
Figure 6.1 and proceeding through Figure 6.1 represents about
l00 years of pumping at low rates with low accretion rates
set to maintain an equivalent pumping water level of 40 feet.
These figures show the progressive penetration of the wedge
and upconing under the well field
Figure 6.7 is the beginning of freshwater recharge
in the coastal zone The initial location of the upper end
of the wedge in this zone is accentuated by the dark line on
the interface as shown in Figure 6.8. The recharge points
are the two dots located above this line. A recharge period
of about 70 years is shown in Figures 6.7 through 6.13.
a,
0
ai
. I^
* f I
FIGURE 6.1
LOCATION, tm
B
f "~4
.'
6* 9 4
9 2.
yWp p
= 0 Min.
INTERFACE
9
9
I, l'
INTERFACE
FIGURE 6.2
LOCATION, tm = 32 Min.
*; I
0
__
f
4k*jhx>aAA
INTERFACE
FIGURE 6.3
LOCATION, tm
= 48 Min.
FIGURE 6.4
INTERFACE LOCATION, tm
= 62 Min.
FIGURE 6.5
INTERFACE LOCATION, tm
= 79 Min.
;^ n*a r r~rrl~rr~~r^_L__xr__lr ; ^ r_.
A
FIGURE 6.6
INTERFACE LOCATION, tm
= 92 Min.
FIGURE 6.7
INTERFACE LOCATION, tm
= 102 Min.
FIGURE 6.8
INTERFACE LOCATION, tm
= 107 Min.
FIGURE 6.9
INTERFACE LOCATION, tm
= 117 Min.
FIGURE 6.10
INTERFACE LOCATION, tm
mn
= 132 Min.
,I b
S~
1 *L*
a I ar
jsS
N930Am
SrAJ7rf
FIGURE 6.11
INTERFACE LOCATION, tm
= 147 Min.
. I
a
IN
*'4
H
z
13
C]
tx21
OH
C:C
0
rt )
II
H
0r
H.
W
The ultimate effect of the recharging is shown in
Figure 6.13 in which the wedge has been depressed and moved
seaward in the area under the recharge wells; and, in the
area under Lake Tarpon, the wedge has been forced upward.
The dark lines on the model were added to bring out the wedge
location at the beginning and end of the recharge period.
Conclusions and Recommendations
As can be seen from the preceding figures, the model
operation is satisfactory from a qualitative point of view.
The problem in obtaining quantitative results at this time
involves maintaining the correct water level in the aquifer
in the area of the Pasco High. During the run from which
the pictures were taken, the water level in this area was
much too low and the figures would therefore indicate more
penetration and upcoming than would be likely.
It is felt that the change in wedge shape due to
recharging is a likely result, and that recharge wells would
be better located at points farther inland. They might be
located in Lake Tarpon or just to the east of the lake. If
the recharge water is to be treated wastewater, the proxim
ity of the recharge wells to the EldridgeWilde well field
would be of great importance.
In order to maintain the correct water level in the
aquifer, an overflow manifold system will be installed on
the back plate of the model. This system should allow the
continuation of the quantitative study in regard to long
term pumping in Eldridge well field, the most advantageous
location of recharge wells and the effect of a new well field
in the Pasco High area. The use of a tracer dye in the re
charge wells would facilitate the evaluation of their effect
on the wedge and on EldridgeWilde well field.
Objectives (1) and (2) are met in that the model is
designed for ease of operation, quick change adjustments,
ability to handle silicone oil and other such features as
discussed in chapter V.
BIBLIOGRAPHY
Aravin, V. I., and Numerov, S. N. Theory of Motion of
Liquids and Gases in Undeformable Porous Media.
Trans. A. Moscona. Jerusalem: Israel Program
for Scientific Translation, 1965.
BadonGhyben, W. "Nota in Verban met de Voorgenomen
Putboring nabij Amersterdam", Tijdschrift van
net Koninklijk Instituyt van Ingenieurs. The
Hague: 1888.
Bear, Jacob. "Scales of Viscous Analogy Models for
Ground Water Studies", Journal of the Hydraulics
Division, ASCE. Vol. 86, No. HY2, (February,
160), 1123.
Dynamics of Fluids in Porous Media. New York:
American Elsevier Publishing Company, Inc., 1972.
Bear, Jacob, and Dagan G. "Some Exact Solutions of
Interface Problems by Means of the Hodograph
Method", Journal of Geophysical Research. Vol. 69,
No. 8, (1964b), 15631572.
"The Unsteady Interface Below a Coastal Collector",
The Transition Zone Between Fresh and Salt Waters
in a Coastal Aquifer, Progress Report #3.
Technion., Haifa: 1964.
Bear, Jacob, Zaslavsky, D., and Irmay, S. Physical
Principles of Water Percolation and Seepage.
Paris: UNESCO, 1968.
Black, Crow and Eidsness, Inc. "Water Resources
Investigation for the Pinellas County Water System,
Pinellas County, Florida", Gainesville, Florida, 1970.
Collins, M. A., and Gelhar, L. W. Ground Water Hydrology
of the Long Island Aquifer System, Hydrodynamics
Laboratory, Report No. 122. M. I. T., Cambridge
Massachusetts: 1970.
Cooke, C. W. "Scenery of Florida Interpreted by a
Geologist", State Geological Survey, Bulletin
No. 17. Tallahassee: 1939.
"Geology of Florida", Florida Geological Survey,
Bulletin No. 29. Tallahassee: 1945.
BIBLIOGRAPHY (Continued)
De Wiest, R. J. M. "Free Surface Flow in Homogeneous
Porous Medium", ASCE Transactions. Vol. 127,
(1962), 10451089.
Glover, R. E. "The Pattern of FreshWater Flow in a
Coastal Aquifer", Journal of Geophysical Research.
Vol. 64, No. 4, (1959), 457459.
HeleShaw, H. S. "Experiments on the Nature of Surface
Resistance in Pipes and on Ships", Transactions
Institute Naval Architects. Vol. 39, (1897),
145156.
"Experiments on the Nature of Surface Resistance
of Water and Streamline Motion under Certain
Experimental Conditions", Transactions Institute
Naval Architects. Vol. 40, (1898), 2146.
"Streamline Motion of a Viscous Film", Rep. of
the British Association for the Advancement of
Science, 68th Meeting. Vol. 136, (1899), 136142.
Herzberg, A. "Die Wasserversorgung einiger Nordsee bader",
Journal Gasbeleuchtung und Wasserversorgung.
Vol. 44, Munich: (1901) 815819, 824844.
Hubbert, K. M. "The Theory of Ground Water Motion",
Journal of Geology. Vol. 48, No. 8, (1940),
785944.
PolubarinovaKochina. Theory of Ground Water Movement.
Trans. J. M. R. De Wiest. Princeton: Princeton
University Press, 1962.
Rouse, H. Advanced Mechanics of Fluids. New York:
John Wiley & Sons, Inc., 1959.
Shames, I. H. Mechanics of Fluids. New York: McGraw
Hill Book Company, Inc., 1962.
Stewart, J. W. "Hydrologic Effects of Pumping From the
Floridan Aquifer in Northwest Hillsborough,
Northeast Pinellas, and Southwest Pasco Counties,
Florida", U. S. Geological Survey, Open File
Report, Tallahassee, Florida: 1968.
U. S. Weather Bureau. Climatological Data for the United
States by Section, Florida Section. Washington,
D. C: U. S. Government Printing Office, 1971.
Vaughan, T. W. "A Contribution to the Geologic History
of the Floridan Plateau", Papers from the Trotugas
Laboratory, No. 4. Carnegie Institute, Publication
133, (1910) 99185.
`f 
1 1
