Analysis of storm water seepage basins in peninsular Florida

MISSING IMAGE

Material Information

Title:
Analysis of storm water seepage basins in peninsular Florida
Series Title:
Florida Water Resources Research Center Publication Number 39
Physical Description:
Book
Creator:
Rubin, Hillel
Glass, John P.
Hunt, Anthony A.
Publisher:
University of Florida
Place of Publication:
Gainesville, Fla.
Publication Date:

Notes

Abstract:
Rapid urbanization, currently proceeding in Florida, has resulted in significant problems with regard to both flood control and pollution abatement. The objective of the reported study was to search for the design procedures that may improve efficiency, safety and adequacy of drainage systems in peninsular Florida. Special attention was given to possibilities of recharging groundwater aquifers with excess storm water. Through such a system, partial solution to problems of inadequate potable water supply in some areas can be achieved simply as a by-product of flood control systems. In Chapter 2 of this report a conceptual design framework is presented. In developing this framework, a variety of disciplines involved in the solution of urban drainage systems in peninsular Florida are considered. Comparatively high levels of groundwater in many locations require consideration of two factors relating to the ability of seepage ponds to divert surface water to the aquifer. The engineer should study the ability of the pond to seep water and he should analyze the response of the groundwater table to the seeping water. With respect to these two factors, analyses and solutions adapted for engineering application are presented in Chapter 3 of this report. The continuing reduction in the availability of potable water in coastal zones of Florida has prompted several communities to design recharging systems utilizing treated storm water as well as effluents. Chapters 4 to 6 of this report concern flow conditions in the Floridan Aquifer that should be considered in connection with this subject. Methods by which the particular phenomena associated with flow conditions in the aquifer can be evaluated are presented. These methods cover several spectral expansion approaches as well as a complete numerical approach.

Record Information

Source Institution:
University of Florida Institutional Repository
Holding Location:
University of Florida
Rights Management:
All rights reserved by the source institution and holding location.
System ID:
AA00001513:00001


This item is only available as the following downloads:


Full Text



2#~r
AL.


Publication No. 39


Analysis of Storm Water Seepage Basins
In Peninsular Florida


By


Hillel Rubin, John P. Glass and Anthony A. Hunt


Department of Civil Engineering
University of Florida
Gainesville


* ~


WIWI


9*1w EWA


77








TABLE OF CONTENTS


ACKNOWLEDGEMENTS . . . . . . iii

ABSTRACT . . . . . . iv

CHAPTER 1: :E!E: pL INTRODUCTION . . . . 1
S'- -.:- of the Study . . . . . 1
Objectives . . ... . . . 2
Methodology . . . . . . 2
Research Program . . . . . 3

CHAPTER 2: FACTORS AFFECTING THE DESIGN OF URBAN DRAINAGE
SYSTEMS IN PENINSULAR FLORIDA . . . . 5
Introduction . . 7. .. .. .. 5
Legal Principles . . . . . 7
Federal, State and Local Regulations . . . 9
Project Feasibility . . . . 14
Preliminary and Detailed Design . . . 19
Conclusions . . . . . . 24
Appendix Flow Diagram for Urban Drainage Design . .. 26
References . . . . . . 32

CHV'TER 3: ANALYSIS OFTRANSIENT GROUNDWATER FLOW FROM SEEPAGE
- i' . . . . . . . 33
Introduction . . . . . . 33
Calculation of Unsteady Seepage from a Pond . . 33
Calculation of the Response of the Groundwater Table . 39
Discussion and Conclusions . . . . 45
Notation . . . . . . 47
References . . . . . . 49

CHAPTER 4: THERMAL CONVECTION IN A CAVERNOUS AQUIFER . 50
Introduction . . . . . . 50
Basic Equations . . . . . 50
Linear Stability Analysis . . . . 55
Finite Amplitude Disturbances and Nonlinear Stability
Analysis . . . . . . 58
Results and Discussion . . . . 65
Conclusions . . . . .. 71
Appendix Expressions for Coefficients of Heat Advection,
Friction and Heat Dispersion Spectra . . . 72
Notation . . . . . . 73
References . . . . . . 77

CHAPTER 5: SEMINUMERICAL APPROACH FOR THE MATHEMATICAL MODELING
OF SINGLY DISPERSIVE CONVECTION IN GROUNDWATERS . . 78
Introduction . . . . . . 78
Basic Equations . . . . . 79
The Flow Field Stability . . . . 82
Analysis of the Steady State Convection . . 83
Numerical Calculations . . . . 86









Conclusions . . . . . 89
Appendix I Expressions for the Coefficients of the
Spectral Functions . . . . . 90
Notation . . . . . . 91
References . . . . . . 94

CHAPTER 6: NUMERICAL SIMULATION OF SINGLY DISPERSIVE CONVECTION
IN GROUNDWATERS . . . . . . 96
Introduction . . . . . 96
Basic Equations . . . . . 96
The Flow Field Stability . . . . 100
Numerical Calculation of the Steady State Convection .. 102
Numerical Results and Discussion . . . 108
Conclusions . . . . . . 114
Notation . . . . . . 115
References . . . . . . 118

CHAPTER 7: SUMMARY AND CONCLUSIONS . . . 120










ACKNOWLEDGEMENTS


This report summarizes the research results obtained in 1975 and 1976

in the Ai,,L .!lic Laboratory of the University of Florida. The research

project was :,r.-:v r",. by the Office of Water Research and Technology (OWRT)

through the Florida Water Resources Research Center.

The investigators are grateful to Dr. William H. Morgan, Director of

the Florida Water Resources Research Center for his help through all

phases of this investigation.

Dr. Sylvester Petryk initiated the study and served as its principal

investigator through October, 1975. His initiative and participation in

this research are greatly appreciated.

Chapters 2 and 3 of this report are mainly based on the Master of

Science theses of A. A. Hunt and J. P. Glass. The investigators are very

grateful to Dr. B. A. Christensen who served as the chairman of their

respective graduate study committees.

Dr. W. Huber served as co-investigator in this research and partici-

pated in the graduate study committees. His activity is very much

appreciated.










ABSTRACT


Rapid urbanization, currently proceeding in Florida, has resulted

in significant problems with regard to both flood control and pollution

abatement.

The objective of the reported study was to search for the design

procedures that may improve efficiency, safety and adequacy of drainage

systems in peninsular Florida. Special attention was given to possibil-

ities of recharging groundwater aquifers with excess storm water.

Through such a system, partial solution to problems of inadequate potable

water supply in some areas can be achieved simply as a by-product of

flood control systems.

In Chapter 2 of this report a conceptual design framework is presented.

In developing this framework, a variety of disciplines involved in the

solution of urban drainage systems in peninsular Florida are considered.

Comparatively high levels of groundwater in many locations require

consideration of two factors relating to the ability of seepage ponds

to divert surface water to the aquifer. The engineer should study the

ability of the pond to seep water and he should analyze the response of

the groundwater table to the seeping water. With respect to these two

factors, analyses and solutions adapted for engineering application are

presented in Chapter 3 of this report.

The continuing reduction in the availability of potable water in

coastal zones of Florida has prompted several communities to design

recharging systems utilizing treated storm water as well as effluents.

Chapters 4 to 6 of this report concern flow conditions in the Floridan










Aquifer that should be considered in connection with this subject. Methods

by which the particular -.T :. associated with flow conditions in the

aquifer can be evaluated are presented. These methods cover several

s.. tral ;...ion ,;- -.:aches as well as a complete numerical approach.










CHAPTER 1

GENERAL INTRODUCTION


SCOPE OF THE STUDY

Since ancient times, people throughout the world have had to cope

with periodic floods and inundations of lands and communities. Notably

pragmatic solutions to flood control were developed in various locations.

These solutions were dependent on human imagination and ingenuity applied

through available resources to a broad variety of situations and conditions.

Even for a relatively limited area such as the United States, we cannot

imagine a singular methodology applicable in every situation. Effective

flood control depends not only on resources and materials provided by

'Mother Nature' but, also on man's attitude toward the application of

these resources. The very definition of the problem depends on the

attitude of the people and their understanding.

Peninsular Florida is a unique system in many respects. In the

early part of this century, flood control projects were begun in earnest

over an area frequently exposed to devastating hurricanes and extensive

flooding. Since World War II, tremendous changes, inherent in the rapid

development of the state, brought about improved techniques and better

understanding of flood control. Presently, the state of Florida faces

simultaneous problems of excessive storm water and limited water supplies

particularly in coastal communities. Rapid growth and the influx of new

people also resulted in reduced water quality, thus adding another factor

to water management and flood control. This project has been initiated

with the idea of simultaneously improving flood control techniques and

enriching groundwater resources in peninsular Florida.










UbULLiiVtc3

The investigators found it necessary to direct their efforts toward

three principle objectives:

(a) General management and development of conceptual design

framework.

(b) Development of a simplified analysis for the evaluation of

seepage ponds in diverting stormwater to groundwater storage.

(c) Development of approaches in the analysis of migration of

contaminants in groundwater due to natural conditions as

well as situations induced by artificial seepage.


METHODOLOGY

With respect to the listed objectives, Chapter 2 of this report

concerns the development of the conceptual design framework as related

to urban drainage systems in central and southern Florida.

Chapter 3 of this report concerns management of water quantities.

This effort in the investigation involved determining the ability of

seepage ponds to divert collected stormwater to the groundwater aquifer.

Emphasis was given to development of simplified methods that can easily

be applied by drainage design personnel.

Chapters 4 through 6 concern mathematical methods that can be used

for the analysis of flow conditions in an aquifer similar to the Floridan

Aquifer. Basic models are suggested and a variety of mathematical

methods are checked from the point of view of applicability, efficiency

and accuracy. These chapters form an introduction to analyses of water

quality problems that have yet to be resolved in Florida.










. 9-r'--i PROGRAM

An initial research pr.:.'. was outlined by Dr. Sylvester Petryk in

1974 when he submitted the research proposal. The study, as envisioned

in that f .::;.1, would consist of the following five tasks

1) Review of existing literature

2) Modeling of precipitation, ,..i,.;,7, and storm water

seepage basins

3) Economic analysis and optimization of design

4) Design procedure

5) .>:. mental measurements.

The research conducted during the past two years has been concerned

with all of these areas. However, the program was modified with respect

to the emphasis and the degree of effort devoted to each task.

A review of the existing literature indicated that our efforts should

be directed more toward improvement of the complete design framework and

project .1.'.e':e.i, rather than involvement with particular design techniques

or procedures. Chapter 2 of this report covers this part of our activity,

which is partially related to each of tasks 3, 4, and 5 mentioned above.

The design of seepage ponds in Florida often presents special problems

because of high *.,..,, eter tables. These problems are the subject of

r--'".. 3 and are related to task 2 mentioned above.

Problems of water quality have become increasingly important in

Florida as well as in other parts of the country. One of our primary

concerns in this regard was the effect that injection of impure water

might have on the quality of Florida's groundwaters. Chapters 4 to 6 of

this report are concerned with this problem, which is related to task 2

mentioned above.










A series of field studies and measurements was conducted as suggested

in task 5. Good results were obtained but the job was more or less

routine in nature and we did not find it necessary to include them in this

report. Our main objective in the field study program was to acquaint

ourselves with local problems and to advise local people with respect to

these subjects.










'. i OTER 2

F.', AFFECT'- THE DESIGN OF .Li DRAINAGE
TEH. IN PENI[. AR FLORIDA




demand for housing in Florida continues to climb as Americans

seek a place in the sun. .: irement and tourism are rapidly displaci.

agriculture in developing and utilizing an attractive natural environment.

'-.th, since the end of World War II, has been phenomenal. Land devel.r-

ment and home construction have become an integral part of the economy.

These activities draw heavily on Florida's natural resources which are

not limitless and, in many instances, are key elements in a sensitive

environmental system. The demands of a progressive economy cannot be

inn:?-1, however, these demands can be adjusted to provide optimal use

of a limited system of supply. Concerted efforts are being made to

7...-ide both understanding of a complex and dynamic system, and a reason-

able balance of supply and demand.

Water, one of Florida's most abundant resources, is a critical factor

in this dynamic s. :t-r Early management concentrated on drainage and

fl;,-. control to such an extent that damage to the system resulted. The

past twenty years has seen a shift in management emphasis as water short-

.;-. and environmental damage became more apparent. Contemporary manage-

ment practices are ...,.jing rapid change. The idea of water as a

scarce commodity has ::. ,tEd basin-wide regulation of consumptive use

and natural flow. .-:-.'adation of water quality has led to more stringent

pollution control laws and the development of improved abatement techni,.: :-.o










All of these efforts have recognized and addressed urban development

as a leading cause of problems in maintaining water quality and managing

flow. Urban storm runoff has been found to be a source of water pollu-

tion that is of equal or greater magnitude when compared with any other

identifiable source. Home building and land developmen'L require drainage

systems that accelerate runoff, consequently, creating or aggravating

flood prone situations. Homes, streets, parking lots and commercial build-

ings retard or prevent infiltration of rainwater and necessary recharge of

groundwater aquifers.

In view of the current activity in water management and pollution

control, the drainage engineer must be aware of the many facets of contem-

porary design and he must apply sound and systematic methods in development

of his design. The natural environment of Florida presents a unique set

of circumstances. Laws, regulations, procedures and techniques, although

rooted in historic precedent and established standards, have all been

developed with some consideration for these circumstances. The drainage

engineer should have a thorough knowledge, not only of applied techniques,

but of social, environmental and political impacts of development. He

should understand recent changes in legal and governmental philosophies.

The engineer needs to know sources of information, plan requirements, and

the procedural aspects of government regulation and approval. Most

importantly, systematic method in design synthesis permits the engineer to

organize and evaluate his data, identify and augment weak points and

effectively manage the design process.

The objective of the present study is to review and suggest a

conceptual design framework as related to the following topics:










1. Basic legal principles associated with drai..ay;., groundwater,

land use and water courses.

2. Highlights and ..*, ..-:' of legislative and regul =i. ," .

on the federal, state and local level.

3. Elements of the feasibility study as the first *-i.:,s in the

design --. .

4. Considerations and techniques applied to preliminary and

detailed design of drainage systems.

5. Integration of the various design factors by means of a flow

diagram.

Finally, this discussion is oriented toward design factors particular to

- .insular Florida.


L5. PRINCIPLES

runoff is of legal concern for a variety of reasons. 'n,:..,-

done concentrated ".'rc, either as a floodwave passing downstream or

a backwater flooding of upstream lands, may result in injury to the

:'--,ted parties. Alteration of flow directions, capacities or other

drainage characteristics can and often does lead to environmental damage,

e.g.,lowered water tables, water pollution, damage to vegetation, erosion

and accretion. -,..*?ty damage from flooding of homes and businesses may

be the direct result of poorly managed runoff and improper land use.

Two basic principles of law concerning disposal of surface waters are

'the civil law rule' and 'the common enemy rule'. Under the civil law

rule, the upland or dominant owner has an easement on the downstream owner

for .;.:- of surface runoff in its natural manner. The common enemy

rule stipulates that the servient or lower owner may take measures neces-

7










sary to keep these waters off of his land. Generally, Florida courts have

followed the civil law rule modified by 'reasonable use'. For example,

the general rule regarding drainage into a natural watercourse states that

a riparian owner may discharge surface runoff without regard to either the

'civil law rule' or the 'common enemy rule'. However, this right is

subject to three limitations which have been imposed by the courts in

varying degrees. These limitations include:

1. Drainage must be reasonable.

2. Drainage must not come from outside the natural basin.

3. The natural capacity of the stream must not be exceeded.

The concept of reasonable use has been applied in cases involving both

land use and water rights. In several instances, the courts have ruled

that land use, adversely affecting surface or groundwater flows, incurs

no liability on the owner as long as his use is reasonable and legitimate

[2, 3]. Rulings have also been made in light of the riparian doctrine,

requiring reasonably unimpaired or undiminished flow.

In conclusion, the engineer must consider possible liabilities and

legal consequences in his design of drainage works. Stormwater systems

designed to parallel natural discharges from a site are desirable from

both an engineering and a legal viewpoint. Further, conflicts between

water rights and legitimate land use have not been resolved in an

entirely consistent manner. It is apparent that the body of law governing

these activities and rights is not well defined and does not truly

recognize the interrelation of land use and natural flow.











p'":!' STATE AND LOCAL FF ,' ..ATIONS

Florida's rapid :th has produced he.:,' and often conflicting demands

on existing water resources. The geology of Florida '... .ides extensive

groundwater sources that have been utilized in satisfying these demands.

.. -lesale acceleration and channelization of runoff and failure to provide

te recharge to ,-.s::rlying aquifers combined with concentrated with-

drawal from these aquifers has led, in several instances, to serious

problems that cannot be resolved on a single situational basis. Total basin

' ,-it has become the imperative solution. Comprehensive legislation

at all levels of government has been enacted to provide a framework for

.,- rm t. Agencies with permitting and other regulatory powers have

been created to affect solutions to apparent conflicts on a regional basis.

The Federal Water Pollution Control Act of 1948 began what has

become a massive effort to improve the quality of our national water

resources. .-e::;t'legislative efforts authorized federal assistance

in research and development of methods of controlling pollution. Most

recently, the F.W.P.C.A., Amendments of 1972 placed stronger fmirqpc,p is on

storm runoff as a source of pollution. Section 208 of this act is

directed toward area-wide planning and management of both 'point source'

and ': .r-,'.cint source' discharges. Stormwater runoff, a non-point source,

is to be evaluated on the basis of extensive data collection and monitor-

ing in the field. Methods for reducing pollutant loading are to be

developdF at the community level and, subsequently integrated into a

basin-wide management plan [1].

Current technology in this area has not established a total 'cause

and ( '-:ct' relationship between runoff and stream quality. Loading rates










are highly variable and depend on land use, storm duration and intensity,

and seasonal weather patterns. This relationship or a reasonable approxi-

mation is singular in its importance with regard to effective modeling and

comprehensive understanding of the problem.

A second program of peaticular interest to the drainage engineer

stems from the National Flood Insurance Act of 1968 and the Flood Disaster

Protection Act of 1973. Through this program, federally subsidized flood

insurance is available to homes and businesses located in flood-prone

communities. A prerequisite to qualification requires that designated

communities adopt restrictive land use ordinances for areas subject to

flooding on a frequency of once in one-hundred years (equivalent to a one-

percent chance in any one year) [4].

The drainage engineer should recognize three significant points in

these programs. First, the federal government has initiated comprehensive

efforts to maintain and improve the quality of natural waters through

strengthened programs at the state and local levels. Secondly, the

flood insurance program is a substantial legislative effort directed

toward better land use practices and flood control measures. Finally,

these efforts actively involve the community and its citizens by promot-

ing improved local ordinances and related decision making.

During 1975. Florida's environmental agencies were reorganized and

consolidated into two major agencies, the Department of Environmental

Regulation and the Department of Natural Resources. Additionally, the

state's water management districts were redefined and were invested with

broadened duties and powers. The Department of Environmental Regulation

is responsible for enforcing pollution control laws and maintaining or










-i.: oving water quality. Permitting of 'i:iii;> treatment plants is a duty

of this =-.:/. In granting permits for point sources such as ..--

treatment plants, the D.E.R. will review and actively consider stormwater

control, treatment and di,;::...,l as part of the permit request [7]. D~.',

highlights of r- 2 guidelines include [7]:

1. T.::..t of stormwater dis(cd'.:.e- on the receiving waters

will be viewed considering designated use of the waters,

practical considerations and cost-effectiveness.

2. Design of stormwater management systems should include

a variety of techniques in reducing impact. Sod filters,

vegetated buffers, sediment traps, primary considerations

in design.

3. =... ,ion basins should be considered as essential in

residential, business, industrial and highway devel..pa~erit..

Interim drainage systems, serving construction sites,

should be anticipated and included.

The Florida Water Resources Act of 1972 empowers the Department of

Natural Resources to accomplish the conservation, protection, management

and control of the waters of the state. In effect, this act makes all

waters in the state subject to regulation [6]. This act further

directed the environmental agencies to formulate a state-wide water use

plan. Elements of the plan are being prepared at the district level and

will be comprehensive in scope and objectives. Other activities of the

Department of ii. 1 e: uces through the water management districts

include permitting for .:c.,- npive use by all users, constructing and

maintaining lands and works incident to management activities, permitting

and re, nation of wells, and management and storage of surface waters.










It is readily apparent that the state is assuming the role of the

riparian owner in many respects. Maintaining undiminished water quality

and quantity is now a valid public concern as opposed to a strictly

riparian right. In addition, the proposed water use plan will include

recommendations as to flood plain zoning and protection from flood

hazards. Thus, the state is assuming an ever increasing responsibility

over land use management and water resources which, traditionally, have

been inherent in private land ownership.

Subdivision regulations and zoning ordinances form development

regulations and standards applicable at the community level. Further,

these regulations outline procedures in obtaining community approval for

a proposed development. They provide an orderly exchange or basis of

communication between the engineer and the community. Frequent discus-

sions during the approval process, enable the engineer and community

planners and officials to effectively demonstrate their needs and desires.

Individual citizens may express their opinions during public hearings

before the local governing body as a matter of state law.

Subdivision regulations vary considerably across the state. Commun-

ity resources and extent of development tend to impose practical limits

on regulatory programs and requirements. Obviously, the most sophisticated

regulations available are relatively ineffective without adequate staffing

and community support. Federal and state efforts are aimed at supplement-

ing community programs and preventing abusive and costly development such

as has occurred in the past.

Development standards tend to fall into two broad categories. The

first is concerned with conditions or programs unique to a community.










Such tors as flood control .- :..:Tfer recharge and salt water intrusion

me,: ire regulations or procedures written specifically -', these

conditions. For example, -, ,.: County has recognized the need for con-

servi natural recha,-: areas. Their -.ulations outline :.-..ial measures

maintaining natural high-infiltration in these areas. "Typical methods

include the use of filtered rer c':,.; wells, bottomless inlets, :. :,. ted

pipe, grading to retard .: artificial ::.:... basins, swales in

street ri.,;':.-of-way, and utilization of natural ,. ...,lation areas" [5].

The second ca'-- of regulations covers more detailed requirements.

Minimum pipe size and material specification, street and lot i.'.ing,

allowable overland flow distances, design methodology, allowable velocities

in sales and ditches, drai.!.-.:e right-. -way i -,-.irements, and storm water

details are t- ical items included in this category. Serviceability and

ease maintenance are important considerations in these -!..-:ifications

relatiJ :.- ly to the 'hardware' of drainage. Routine maintenance e.

dra-" ..- systems and streets is a -.-.jor cost item in municipal '...:- ts.

Flooded streets and overgrown ditches are common complaints .., tax-

,-.~' s and -:2. .. is made to minimize such situations in writi .-.;

these regulations.

In summary, legal principles ,... : .ing surface flow and .'....... ter

are derived from ri- ... associated with private ownership. Conflicts

arising from alteration of natural flow patterns or from consumptive uses

have been resolved considering reasonable use of land and water. It is
--:-.rent that these principles are not adequate to resolve *.:.o iicts

consider. cont.- :--y demands. As a result, regional or basin-wide

, :.:.:. of water resources is being affected t:-..'.- combined 1.-,:isla-










tive efforts on all levels of government. Preemption of private water

rights appears to be justified in that the burden of flood damages and

optimization of consumptive uses rests with the entire community. Finally,

lawmakers have broad powers in preserving and protecting natural resources

in the public interest. Thus, considerable legislation has been enacted

to improve land use, relieve flood hazards and reduce pollution of

natural waters.


PROJECT FEASIBILITY

Project initiation usually begins with a discussion including the

client and key staff ,I-.-. -ionals. Subsequent efforts depend to a great

extent on several factors. Among these are the intended scope of the

project, potential market, capital availability, anticipated scheduling of

planning, engineering, construction and sales, legal considerations and

immediate technical needs. Certainly, forty acres in an established

residential area will have significantly different needs than ten-thousand

acres of undeveloped and agricultural lands. The experience of the client

or developer will also influence project initiation. An experienced

developer often will have established project economics, scope and

approximate timing prior to discussion with design professionals. Project

feasibility may, however, include multi-disciplinary studies and analyses

of the previously mentioned factors.

One of the most fruitful approaches to evaluation and design is

organization and coordination of all related efforts. Project planning

draws each factor into the design process as it becomes pertinent. A

systematic approach to project management helps to ensure optimization of










technical i and ;:. :..tion and, that the quality of the : r.;c will

be ':-.. ....ily .. ...',ssional.

i. ::. design as an -i Yineering .: ., is well suited to an

or ized ... -- : ._.ign .. ,thesis and coordination of -:i .. in drain-

S-ineeri., is : .,!cally demonstrated in the flow diagram as shown

in the '. -,,.x .-J from Woodson [8]. Flow di c.:ams are an .- -.

and proven :.... -.t tool. Their inherent flexibility makes them

adaptable to a variety of situations and design conditions. i. included

flow di am,presented in the Appendix, is divided into three distinct

phases allegingg the balance of this discussion. I,.'- phases include

-"-.'ect i- ability, preliminary design and detailed design.

Intensive data collection is the first step following project initi-

ation. Gathering information on existing drainage patterns and both

current and future land'use follows several avenues. A thorough under-

standing of the existing drainage situation is essential not only in

determining required improvements but, also in demonstrating the effect-

iveness of ,*'r-,.~c~ improvements. Standard sources of information that

will readily provide data on drainage, land use and flood p.;,-stial should

first be invest ...-. .*".L:g these sources are:

1. U. S. Geological Survey for topographic ',-.-:, flood studies,

groundwater studies, geologic investigations, well :....,

stream and lake gage records and general hydrologic information.

2. ";,.'icipal or County Engineering and Planning for subdivision

regulations, zoning ordinances, aerial mapping, established

drai .r--. networks and other i '.-<.on directly related to the

site.










3. Soil Conservation Service for soil surveys, drainage

techniques proven in the area, recommended best use of

the site and extensive technical expertise in drainage

and hydrology.

1. State age,: f-k for re.,ia:y informationn and design

guidelines in stormwater management. Specific informa-

tion as to local environmental considerations and

current practices.

5. Federal agencies for flood insurance and flood plain

information, possible federal permits, design standards

for qualification under federal loan programs and many

other technical services.

The reconnaissance survey of the site of a proposed development is

particularly important in that information resulting from visual inspec-

tion is of primary importance and is generally available from no other

source. Initial data may be summarized and plotted or noted on a site

plan to facilitate observations and comparisons during the survey. Veri-

fication of land use, drainage patterns and improvements, both on the

site and in the surrounding catchment, is an important objective. Other

factors to be considered during reconn' are:

1. Attractive or beneficial features on thesite that will

contribute to aesthetics and environmental quality.

2. Hydraulic constrictions and control points, ponds, lakes,

streams, springs and natural depressions collecting run-

off.











3. High water marks on trees and drai:.,-.-, structures, water

-- elevations and groundwater or water table eleva-

tions.

4. Soil t .- and area distribution. Soil -.,. ..its such as

'c muck and that must be considered in design

and construction.



receiving waters ..--bI': le quality of current .".

and receive -waters.

6. Identification of chronic drainage ::,.l':,.ems and existing

i,-.,....:,'A,.: ,- both on the site and in the immediate catch-

ment.

'.- i..:- or following field investigation, the engineer should outline require-

ments ,- additional ;... ,aphic surveys, soils investigations, stream .

ing, water quality evaluations, environmental surveys and similar e-.'. .

and techni -:.: for field investigations are '.,-er....:. on site

conditions and the extent of proposed development. Detailed field notes

and -;.-rrnlementary "',,Co -.. I are common to all surveys and provide

documentation to :. :. subsequent j'.-m.!i:., public presentations and

technical discussions.

Evaluation and incorporation of collected information involves simu-

lation of existing drainage conditions, potential for improvement, delin-

eation of drai....?-- systems and land use planning. 'r..!!';' calculations

or :. ...- imations should be made in evaluating the hydraulic and hydrologic

characteristics of the existing system and the ..;.* improvements.

drai ....- engineer should work very closely with those involved in land









planning. Land use, of course, will have multiple effects on any proposed

drainage system. Compatibility of land use and soil characteristics is

of particular importance in Florida. Planners must be aware of probable

water table elevations and seasonal fluctuations. Much of Florida con-

sists of fine sands-with small percentages of silt in the surficial layers.

Often, a semi-permeable layer of consolidated fines and organic material

underlies surface sands. Infiltration of stormwater is severely retarded

and lateral or surficial flow may be very slow, particularly on the

marine terraces and in the 'piney flat woods'. The drainage engineer

should evaluate cost-effective measures in reducing high water tables

and flat gradients and, he should advise land planners of his findings.

Land planning should be coordinated to reflect drainage methods and

alternatives best suited for the site, required rights-of-way, pond loca-

tions, flood elevations, fill sources and quantities, recharge sites and

natural features to be preserved. Every effort should be made at this

time to promote an aesthetically pleasing drainage network. An open

ditch or floodway is patently obvious and prominent in appearance. The

additional costs associated with a landscaped and meandered channel are

more than justified considering heightened sales potential and consumer

ideals.

The preliminary engineering report or feasibility study reflects a

combined effort in engineering, planning and project economics. Form

and content for such a study are variable and, generally, depend on the

scope of the project and desires of the client. It cannot be over-

emphasized that such a report is preliminary in nature and that subse-

quent efforts may contradict stated conclusions, cost projections and










1 .,ity estimates. Misunderstandings are frequent in this c.: t and

the engineer should use care and discretion in making this point clear.

[.: lations, -,,:,-mitting and o.,.-..- ei t coordination are often disegrJ-..3

at this stage of the process. The feasibility study should include a

Y'.L* :ji! discussion of all aspects of the permitting and
J L_-1 ,o

Preliminary engineering design provides the client with valuable

i '.fen-tion -.. project scheduling and funding. The feasibility report

serves the developer in anticipating construction costs, cash flow,

marketing and sales. Finally, the engineer and client can selectively

evaluate alternatives and determine preliminary design objectives.


PREL -IV.i.Y AND.DETAILED DESIGN

Preliminary design should be approached considering overall site

i logy. Simulation techniques commonly used include the rational

method, unit "--ii ,;, synthetic hydrographs, regression equations and

more complex computer modeling techniques such as those developed by the

r*... of Engineers and the Soil Conservation Service. The object of all

of these methods is to systematically quantify runoff and to route these


T :;.ically, the outfall point is identified together with available

records or some estimation of appropriate water surface elevation.

The flat gradients, so common to much of Florida, necessitate a thorough

investigation and a reasonable approximation of downstream water surfaces

corresponding to the events being simulated. Collector systems are

usually a combination of open channels and storage basins that depend on










both physical ground slopes and hydraulic head in maintaining positive

flow. Central Florida's highlands may require grade stabilization

structures where slopes are relatively steep. Peak velocities in unlined

channels generally should not exceed two feet per second in Florida's

sandy soils. Less frequent storm events may be permitted to generate

somewhat higher velocities. Storage-elevation curves and stage-discharge

relationships may be approximated in routing surface runoff through the

system. Routing begins with generation of water quantities at the up-

stream end of the system and proceeds down the slope to the receiving

waters. Peak discharges and water surface elevations, derived from the

routing process, may then be used to calculate water surface profiles.

Most systems generate peak flows over a relatively short period of time

so that only peaks may be considered in calculating the profiles.

Subcritical flow profiles must begin at the receiving waters and, based

on the appropriate surface elevation, proceed upstream. Adjustments are

made in channel geometry, storage, control structures and overland flow

times until close agreement is achieved between the water surface profiles

and the elevations generated by routing. This procedure begins in the

preliminary design phase and is continually refined through the end of

the design sequence. Many factors that are elemental to simulation will

have considerable variation depending on the event being considered. Rain-

fall intensities, for example, will generate system responses markedly

different for various events. It is advisable to simulate storms of

several frequencies and durations even though regulations may specify

design based on one or two standard events. Practical considerations,

such as allowable design time, must be considered in managing this

procedure.










Internal drainage systems direct runoff to the collectors. The

rational method is almost universally used to design local drainage such

as swales, culverts, inlet size and location, street slopes and lot -i'd.

ing. Home building and related land development require simultaneous

evaluation of earthwork and drainage. Generally, surface slopes may be

adjusted to minimize drainage improvements. Cutting and filling required

to achieve ideal slopes may, however, result in excessive costs in earth

moving. Once again, several iterations of site engineering may be

required to achieve optimal conditions in detailed design.

Incorporation of detention-retention ponds in the drainage system is

desirable to attenuate increased discharges resulting from development.

Pondi .. is a recommended method in reducing pollutant loading. Common

practice is to provide sufficient storage within the system for the

volume of runoff resulting from a given storm event. These events or

design rainfalls generally require enough storage so that runoff from

high ,.*. .,ny events, such as afternoon thunderstorms, is totally

retained or detained for an extended period. Silts and other suspended

solids have a chance to settle, evaporation and infiltration effectively

reduce the volume of stored runoff, biological action and dilution

reduce the concentration of pollutants, filtering devices prevent dis-

charge of surface borne contaminants such as oil and grease, and trash

and debris are prevented from entering the receiving waters. Natural

storage areas, including cypress ponds and mangrove swamps, may be

utilized in developing detention. These areas have the advantage of acting

as filters for runoff, thus improving discharge quality. Care must be

used to maintain minimum water levels in these areas and to prevent










damaging pollutants from entering such a system that utilizes natural

vegetation. Other considerations necessary to retention-detention ponds

are:

1. Storage capacity recovery through evaporation, infiltration,

structural devices regulating outflow, modified underdrains

and pumping.

2. Siltation or clogging in the pond. Two-stage ponds should

be provided where sediment loads will be high. Periodic

cleanout of settling basins in the first stage will be

required. A subsurface berm or equivalent structure should

separate the settling basin from the second stage.

3. Maintenance, safety, vegetation, normal depth, fill needs,

aesthetics, emergency overflow, and drawdown capability are

among the additional considerations.

Systems of interconnected lakes and ponds are frequently employed in

Florida to provide fill material, aesthetic appeal and recreational

benefits. In addition, extensive open systems tend to lower surrounding

water tables in chronically wet areas.

All drainage systems should be viewed considering the effects of

the maximum probable storm. Surface drainage systems may be inundated

in some areas causing streets and low lying areas to act as open channels

or ponding areas. Utility services may be interrupted. High water in

the outfall system may result in flood conditions for an extended period.

In this regard, auxiliary floodways may be desirable, road elevations

should be sufficiently high to allow limited continuous access, floor

elevations near points of concentration should be higher than arbitrary










minimums and off-site structures and embankments should be evaluated as

flood hazards.

In ? :-.:t, preliminary design involves evaluation K. a selected

development plan. .- .ign efforts should reflect sufficient detail such

that .' limited alternatives remain together with final details neces-

sary to construction. The flow diagram,- .'.:' ted in the Appendix, outlines

elements of this phase through client and agency approval. Presentation

of the preliminary plan C.I...,,- with supporting data, becomes part of a

comprehensive review by the client and staff. These plans may then be

submitted to local agencies for preliminary approval and required zoning
*", .:,o Review and discussion with state and federal agencies should

be conducted on an informal basis. Results of these interviews and

evaluations should be summarized for incorporation into detailed design.

:-e third and final phase, detailed design, completes this sequence.

Construction plans, specifications and detailed cost estimates are primary

objectives. :.-'isions to the preliminary design should be solicited from

staff, client, agencies and potential construction contractors. The

design -r.:-... must provide a thorough check of every detail and should

coordinate with other members of the staff in locating and resolving

possible conflicts. Almost invariably, insufficient clearances between

sanitary sewerage, storm mains and water supply will become evident and

necessitate : *e!...*. Experienced contractors will often point out plan

discrepancies, ,*.,.ity errors and unfavorable site conditions during bid

0=.ip- ation or contract negotiation. During this phase, the engineer

should coordinate closely with permitting and approving agencies in

..-iding '- ting data for selected design features. Unnecessary and










costly changes, based on arbitrary interpretation or judgmental differ-

ences, may be avoided through foresight and preparation. All parties

should be informed of required revisions and all plans should be noted

accordingly.


CONCLUSIONS

Drainage design has become considerably more than estimating runoff

and selecting culverts. The design professional is obliged to assume a

broad viewpoint in developing his program and directing design efforts.

In summary of this discussion, the following conclusions are presented:

1. The engineer should be familiar with basic legal principles

related to drainage and how they have been applied in

resolving conflicts and water rights litigation.

2. The engineer must be aware of rapid changes taking place

with regard to water rights and resource management.

Federal and state legislation in this area has shifted

rights and responsibilities away from the landowner and

toward the public entity.

3. Local regulations and ordinances are being re-written or

supplemented in order to effect better land use practices

and flood protection measures.

4. Substantial efforts are being made to develop and implement

improved pollution abatement techniques. Urban runoff, as

a prime source, may be improved through a variety of

techniques presently available to the designer.

5. Legislative mandates have been very direct. The engineering

profession should take precautions in protecting their design










pr*.-.H*;vves and, at the same time, align their efforts

with these mandates.

6. The feasibility study is, possibly, the most if-.o t..Iit

phase in design sy.:': sis. This study develops a set of

conditions influencing the entire project.

7. Florida's environment and .di['I! ..' on well */, f,2er

water resources emphasizes the efforts of the drainage

engineer. He must actively participate in several areas

of project development.

8. "-.i.ition-retention systems have become a necessity in

designing improvements. Ponding or storage improves

water quality, reduces flood potential and provides

increased aquifer recharge.

9. Management of the design process and attainment of

technical efficiency may be aided through the use of a

flow diagram presented in the Appendix. Considerable

flexibility allows this technique to be adapted to

most projects.

Finally, the attitude of the engineer must be receptive to change and inno-

vation. Bending to accommodate is not the answer. Embracing and assertive

leadership is a trademark of the professional.











APPENDIX
FLOW DIAGRAM FOR URBAN DRAINAGE DESIGN

Phase I: Project Feasibility




Client Desires

Staff Discussion

Project Scope -Project Initiation -

Project Economics

Time Schedule
Outline First Efforts






Is there sufficient ? No
information to proceed?

Yes

Existing Land Use

Soil Survey

Topographic Surveys

Historical Records

Interviews Data Collection ---- ---

Future Land Use Studies

Laws and Regulations

Reconnaissance Survey --











C Establish Preliminary
DDesign Objectivesl o




Are pre1imif"i.t design -? No
objectives consistent and
reasonable?

Yes Phase I
---------------------------------C o m p l e t e---



Phase II: Preliminary Design



Development of Internal
-. inage System

Development of Collector
and Outfall System

Identification of
Critical Points

Development of Hydrologic
Model

Development of Hydraulic
Model Development of
Preliminary Design ---
System Evaluation
Through Models

Identification of
Limited Alternatives

Water Quality Evaluation
Prepare Preliminary
Design Plans












Have all sources -
been tried?





Delineation of Existing
Flow Patterns

Identification of
Constraints

Simulation of Existing
Drainage Conditions

Proposed Drainage
Schemes

Proposed Land Use







Are proposed land use
and drainage plans workable?





Cost Projections

Quantity Estimates

Anticipated Permits
and Approvals

Staff Review and
Coordination

Owner Review and
Approval


Presentation of
Feasibility Study











Have preliminary design
objectives been realized?





Design Features

Limited Alternatives


Special Problems and
Solutions

Cost Estimate

Additional Data Needs


Presentation of
Preliminary Design


Obtain Client and
Agency Approval


Are final design features
and objectives approved
by client and agencies?


Phase II
Complete


Phase III: Detailed Desiqn


Incorporation of Testing
Results and Survey Data

Safety Considerations

Cost Optimization

Coordination of Design
and Drafting











Design and Specification
of Special Features

Preparation of Detailed
Cost Estimate

Preparation of General
Specifications

Checks and Revisions of
Calculations and Drawings

Project Plans, Estimates
and Specifications Review


Development of
Detailed Design


Issue Plans Reflecting
Final Design
E#


Is the final design complete ?
and acceptable to the client?


Yes


Permit Applications
and Procedures

Local Governmental
Approval

Contractor or
Bidder Review

Final Staff and
Client Review


Revision of
Detailed Design


Obtain Final Approval
by All Parties


i









Are plans, specifications,
estimates, and permits
complete? Are Ep!.rovals
in receipt?


Phase III
Complete


Symbols Key


Flow Direction:



Specific Activity or Information:



General Activity:



Executive Action:



Decision:


.^^ --^ -BBgffi.^ --l l....DO. ,r.


I









REFERENCES


1. "Executive Summary of Section 208 Program for Designated Area -
Federal Water Pollution Control Act Amendments of 1972", United
States Environmental Protection Agency, Washington, D. C., 1974.

2. Linsley, R. and Franzini, J., Water Resources Engineering, Second
Edition, MIcGaw-Hill, New York, 1964.

3. Maloney, F., Plager, S., and Baldwin, F., Jr., Water Law and
Administration, The Florida Experience, University of Florida
Press, Gainesville, 1968.

4. "National Flood Insurance Program", United States Department of
Housing and Urban Development, Washington, D. C., 1974.

5. "Orange County Subdivision Regulations", Orange County, Orlando,
1974.

6. Tebeau, C., "South Florida Water Management" and "Environments of
South Florida: Past and Present", Miami Geological Society,
Miami, 1974.

7. Thabaraj, G., "Regulatory Aspects of Storm Runoff Control",
Proceedings of the Storm-Water Management Workshop, Florida Techno-
logical University, Orlando, February, 1975.

8. Woodson, T., Introduction to Engineering Design, McGraw-Hill,
New York, 1966.










CHAPTER 3

ANALYSIS OF TRANSIENT GROUNDWATER FLOW FROM SEEPAGE PONDS


INTRODUCTION

The ability of a seepage pond to divert storm runoff into the ground-

water aquifer depends on two processes: (a) the rate of seepage from the

pond, and (b) the reaction of the groundwater table.

---The-acstenr- of-a seepage-pond-in-dispos-ng-of-storm-watertis-a

transient phenomenon. The inflow to the pond, specified as a runoff hydro-

graph, causes variations in the depth of ponded water. Outflow is consid-

ered here to be entirely by infiltration, which varies not only because

of changes in pond depth but also due to variations in such soil proper-

ties as hydraulic conductivity and storativity. The growth and decay of

the groundwater mound in the underlying aquifer is also time dependent.

These facts are well known among scientists but are usually not taken

into consideration in the design of a seepage pond. The objectives of

this study are to develop methods by which the designer will be able to

estimate seepage pond effectiveness by considering both the rate of

seepage from the pond and the variation in the groundwater table while

taking into account the transient nature of the problem.


CALCULATION OF UNSTEADY SEEPAGE FROM A POND

The adequacy of a seepage pond is evaluated by a storage routing

procedure, which is basically an account of the inflow, outflow, and

change of stored volume over successive discrete time increments. The

inflow is represented by the runoff hydrograph for the design storm on

the area to be drained. For the purposes of the present exposition it










will be considered a given function of time. The outflow, on the other

hand, is dependent on time, depth of ponding, and the properties of the

soil.

Consider an initially dry seepage pond constructed in unsaturated

soil. As it fills, the water begins to escape from it by vertical unsat-

urated flow, or infiltration. The most convenient way to describe this

flow is by the formula of Green and Ampt [1911] as modified by Bouwer

[1969]. It has been shown that this approach, which was long thought to

be purely empirical, is soundly based on physical principles [Morel-Seytoux &

Khanji, 1974] and gives very good answers [Whistler and Bouwer, 1970].

Green and Ampt based their derivation on a simplified model of infil-

tration which treats the soil as a bundle of vertical capillary tubes.

The vertical hydraulic conductivity and moisture content of the unsatur-

ated flow are considered constant, as is the capillary suction potential

of the advancing wetting front. Applying Darcy's law to this idealized

flow, the infiltration rate is described as



H + L + 'n (1)
w = K (1)


where w = infiltration rate

Kt = hyaraulic conductivity of the transmission zone

H = depth of ponded water

ln = capillary suction potential

L = depth of penetration of the wetting front.










.i,- rate of advance of the wetting front is


dL w (2)
dt f


where t = time

f = the volumetric fraction of fillable pore space.

The equation of continuity applied to the pond yields

dv dH T f d (3)
dt A = I f(AL) (3)


where v = stored volume in the pond

I = inflow

A = area of pond surface

Af = area involved in infiltration.

Introducing equation 3 into 1 and 2 we obtain, after differentiation, the

following non-linear differential equation

S2 dL dAf
dt2(L2/2) + (fAf/A-l)-- + fL/A dt- I/A = 0 (4)
dt

which can be solved numerically.

However, in most practical cases the rate of change of H is about

an order of magnitude smaller than the rate of change of L. In such cases

an integration of equation 2 between ti and ti+1 provides a simple

discrete presentation of the problem, as follows

i+1 i+1l
f dt = w dL.
ti L f









The integration leads to


ti+ -t = f/Kt[Li+1-Li-(H +*n))n( H ++L++p )]. (5)


In dimensionless form this can be expressed as


AtK AL i' Li+i/Li + r. i/L
fLi L. Li 1 + r /L(6)


where rF = Hi + n

L. = value of L at time t.

L i+1 = value of L at time t.i+

At = ti+1 t

AL = Li L.
H = ()(H. + H i+), the mean value of H.

A chart that can be applied for engineering design purposes is shown

in Figure 1. The curves are based on the solution to equation 6.

The application of the discrete approach of equation 6 was

checked against experimental data reported by Weaver and Kuthy [1975].

In this experiment a seepage pond was constructed and filled with water

at a controlled rate. The variation of the pond volume and area vs. depth

is shown in Figure 2. The experimenters listed soil test data which led

to the following values of soil parameters.


Kt = 1.2 ft/hr f = 0.2 n = 0.5 ft.


According to the reported measurements, the rate of change of H was

5 to 10 times smaller than the rate of change of L. This would indicate












0.5 -- ---_
_= L____ __ _





r=o.5
S 0.05


0.0 =0.1

0.01 r = 0.05 1 =..

0.005
r=0.02 --- _

r=0.01 -;_:

0.001 r=o.005 : r=0.001 ___ __
0.001 0.005 0.001 0.005 0.01 0.05 0.1 0.5 1.0
(AtK i)/(fLi)

Figure 1. Solutions to Equation 6 fo: Li+ /L = 1.2
i+1







AVERA-.E A PE (thousand square feet)

3 4 5


5 10 15 20 25 30 35

POND VOLUME (thousand cubic feet)

Figure 2. Depth vs. Average Area and Volume for Test Pond










that the method utilizing -..i',,tion 6 should give :J.oo results. The

,.ii-." ison between the .-.., : '.l data and the theoretical ," fiction

is shown in Figure 3. The .,_ .n:en. is, indeed, quite good.


CA,.. ACTION OF THE ,'.: ..;: OF THE GROUNDWATER TABLE

The reaction of the water table aquifer to vertical infiltration

a ... g pond is characterized by the growth of a groundwater mound. For

analysis we consider a circular pond constructed in homogeneous and iso-

tropic soil. The aqui:..-- is underlaid by a horizontal impervious layer

and initially has a horizontal free surface. When the wetting front of

vertical infiltration reaches this free surface the mound begins to form.

From this time the infiltration is assumed to continue at a constant

rate.




































1 2 3 4 5 6 7 8 9 10 11


10

9

8

7

6 u

5 -
L-
LD
A im


12 13 14 15


TIME (hours)


Figure 3. Comparison of Experimental and Calculated Storaqe Routinq










In a i'ndri cal coordinate L .-:.centered on the ..e :.:- -,. as
.. in Fi ,.. 4 the axis :..'i, .c flow in the saturated mound is described



h3 K K s- h h. + + (7)
at n 3r nr ar n

where Ks = saturated J J.'.-.:.ic conductivity

n = .-,- ific yield

w = rate of vertical infiltration.

w and n are -... .tions of radial distance, r. The value of specific

yield is :t. -tly reduced in the zone of infiltration under the pond.
on the .- i.:-. mental evidence of Bodman and Coleman [1943] it seems

reasonable to assign a s.i:.'-fic yield to the zone of infiltration with a

value of of the specific yield elsewhere. The rate of infiltration,

w, is zero when r is greater than the pond radius, R .
ing that h = S + a, equation 7 can be rearranged:


3S Ks +S Ks [S)2 s 9S w
tS = S + a) + 2 + r S + a + (8)
at n n ar rn ar n

To make the .--..tion dimensionless, introduce the following dimension-

less variables:

t = tks a/R, r= r/R., S' = S/a

Substitution of these variables along with some manipulation and omitting

the .* 'mes we obtain:

.a -2(S2/2) + 1 + r 2) + +
4 (S2/2) + 3 rr (S'1) + + W/%) (9)







































Figure 4. Groundwater Mound










In this form, the non-linear terms are distinctly separated from the
linear terms. They are treated differently in the numerical analysis.


S = 3 0 at r = 0
(10)
and

S = 0 at r =

The initial condition is:

S = 0 at t = 0 (11)

To solve the problem by finite differences, we must set up an
appropriate y;'"e time grid. The time increments are designated by i's
and the space increments by j's. The time derivative of S is approximated
by a ,.r,,:rd difference.


a i Si, j/At (12)
A i+lj 1,j

The non-linear spatial derivatives are represented by central difference




r (S2/2 S iS2 /(4Ar) (13)


4 (S2/2) S2i,j+l 2S,j + S2,j+l/2Ar2 (14)

The linear space derivatives are the most influential terms in the
equation. If the solution were to become numerically unstable, it would









probably be because of these terms. Therefore, they are approximated by

the mean of the finite difference representations on the (i + 1)th and the

(i)th time rows. This is the implicit method of solution developed by

Crank and Nicolson [1947] which has good stability and convergence charac-

teristics:

aS 1 ( )/(2r)
r 2 I Si+l ,j+l i+l,j-_ )/(2Ar)
+ (S i,+- Si,jl)/(2Ar)] (15)



r-4 I l Si+l,j+l 2S + Si+lj-_l)/Ar2

+ (Si,+l 2Sij + Si,jl)/Ar2 (16)


When these finite difference approximations are substituted into

equation 9, a tridiagonal system of equations is generated which can be

solved by Gaussian elimination.

The computer program developed from these finite difference operators

was tested with several combinations of grid mesh ratio and grid size.

It was found to give stable results with a ratio of At/Ar2 < 0.11. The

mesh size chosen was Ar = 0.03. To simulate the boundary condition at

infinity, the calculations were carried out to r = 25 at which point S

was required to equal zero always.

For comparison with the linearized analytical solution of Hantush

[1967], the finite difference program was run with constant specific

yield, n. The results of this run are shown by the dashed lines in











Figure 5. They are practically identical to the results obtained by the

method of Hantush. ,c n was allowed to vary as a function of r, the

results were .-. much different, as shown by the solid lines in Fi..L 5.


D? .ION AND -I.,LI-:% ONS

analysis of the hydraulic .,* .ation of a ::. .....d. should be

conducted in two phases. The first phase is concerned with the rate at

which water seeps out of the -e.:Ci by vertical infiltration through its

bottom. The equation for vertical unsaturated flow devel..:o. by Green

andi Amt [1911] can be used to describe the outflow in a s:. -.. routing

-:,,;.s. The storage routing is simple enough to be done by hand when

the soil is assumed to be homogeneous. The Green and Ampt equation can

also be adapted for use in soil where the hydraulic conductivity varies

monotonically with depth [Bouwer, 1969] which is a very common

occurrence. This calculation is about as simple as the s,:.:.e routing

e...we but the combination of variable ponded dJei and variable

!--w!.aulic conductivity would be cumbersome enough to make the use of a

,-o,,.i-.r or a programmable calculator desirable.

The second phase of the analysis deals with the response of the water

table to recharge from the seepage pond. The determination of adequate

S--.:,_ity of the groundwater aquifer requires use of the design charts for

each specific case. The design chart itself is based on the numerical

solution of the non-linear differential equation.

The i ,~. ":'..:- of the finite difference solution of the i.ev.,.,ter

mound is that it makes possible the consideration of an axisymmetric

variation r- the specific yield. It would also be possible to include

the effects of axisymmetric variation in other soil properties.











1.14


.. ... .3 P 0 .3
'J ... 0 8 ... ......
1 12 -~~ ~.. .

W 2. P = 0.15....
.... .... ...- -.

1.10






... ...... ...

S .02 .05 .. .2 .. 1 2 5







c-a ... .. ...




1.02 __' ^' '.'" '_ _.*J ^ -'. -___________* -- -Variable n
.02 .05 .1 .2 .5 1 2 5 10

DIMENSIONLESS TIME, t'


Figure 5. Mound Height at r' = 0, Comparison of Constant and Variable n:










It should be noted that the saturated ,-.aulic conductivities

to in the two phases the anal ,is are, ...- ...:., not the same. .' V1-

tration is concerned with vertical flow while the ,l .,i..:-.ter mound is

mainly i --enced horizontal flow. ,: difference between the saturated

'. .... ic conductivities in the two directions usually stems -... hori-

zontal la -"i the soil.


NOTATION

a initial thickness -- -ffer EL]

A area FL 2

A f area of i IStration [L ]

f volumetric ', tion of fillable .;. ; ...

h thickness of aquifer EL]

H : ". of ,-. -.. -' water EL]

I i :ow to E&.,,. [L3]

Ks .- -aulic conductivity of saturated soil [L/T]

Kt '.K... aulic conductivity of the transmission zone in unsaturated
flow [L/T]

L .-:th of .--tration of the wetting -.,. [L]

n specific yield

P dimensionless recharge intensity, wR /K a2

r radial distance EL]

r -'.. .ionless radial distance, r/R0

R0 radius of :-.. pond [L]

S rise of groundwater mound above water table [L]

So dimensionless mound height, S/a

t time [T]










t dimensionless time, tK a/nR

v stored volume in pond EL3]

w infiltration rate [L/T]

r On + H [L]

In capillary suction potential at field capacity EL]













., G B and lean E A. ":",sture and .:,:.. Conditions
. ,:- .: .rd Entry Water into Soils", Soil Science Society of
:__. -.-_ Vol. 8, 3, pp. 116-122.


, nfiltration Water into '.. .,;form Soil" '. -I.=-. ..-
.. .. :.., ^ .. L._,2 Irrigation and : .'
1. '** pp '.: .


Di isr, H


'rank, J. and ,:colson, P., "A Practical Method for Numerical

Heat Conduction T..L-", .- inqs Cambridge Phil(.- ;i.cal Society,
Vol. : -.17, pp. 50-67.

Green, W. H. and R- G. A., "Studies on Soil --K.,ics. I. The
Flow of Air and Water Th'..- Soils", Journal Agricultural Science,
Vol. 4, I pp. 1-24.

Hantush, M. S., r'---'th and Decay of Groundwater-Mounds in Response
to Uniform Percolation", Water *..L'.rces Research, Vol. 3, No. 1, 1967
'-234.

I:1- *toux, W. J. and -.;r,;i, J., "Derivation of an Equation of
Infiltration, Water ..i.rces :1.. -..rch, Vol. 10, No. 4, '-'74,
)p -'. :

Weaver, R. J. and Kuthy, R. A., Field Evaluation of a a-e.-. -ge Basin,
.:. York .'te Department of Tr-,: .. ,ation, Engineering Research
and Devel .'" Bureau, :.:.,-.-ch '.'-:. 26 1975.

Whistler, F. D. and Bouwer, H., "Comparison of v;, : for Calculating
Vertical :i,.:.. and Infiltration for Soils", Journal of Hydrol.._,
Vol. 10, '!. 1, 1970, pp. 1-19.










CHAPTER 4

THERMAL CONVECTION IN A CAVERNOUS AQUIFER


INTRODUCTION

In a previous article [Rubin, 1976], referred as I in this paper,

the author analyzed instability criteria related to onset of thermohaline

convection in an aquifer whose properties are similar to those of the

Boulder Zone of the Floridan Aquifer. Such an aquifer is characterized

by extremely large pore size and transmissivity leading to very intensive

solute and heat dispersion as well as to invalidity of the laminar Darcy

law even when flow velocities are extremely small.

In I it was found that for moderate Reynolds numbers the doubly

diffusive convection can be approximated by the singly diffusive convec-

tion, for in such cases mechanical dispersion is larger than molecular

solute and heat diffusion.

The objective of this study is to analyze transport phenomena in the

cavernous aquifer subjected to singly diffusive convection.


BASIC EQUATIONS

The analysis is related to a flow field similar to the Boulder Zone

(the deep regions) of the Floridan Aquifer. We assume that density

gradients are induced by a single component only, referred to as temper-

ature. In the cavernous strata, turbulent effects, as well as mechanical

heat dispersion, are induced by even extremely slow fluid motions. In I

the basic equations related to the calculation and approximated by the

Boussinesq approach were presented as follows:











0
x,


+ pgn + (1 + b)ui = 0


I -_+ u T (E
S axi. ij axj


p = pO [1 a(T TO


u. = velocity vector

p = :" _.sure

p = fluid density

K = ,. i,,--.-,!,ility

= porosity

u = viscosity

T = i.-'.cient of thermal

a'. ..i:;icient of thermal


:h.n:. ion


p oj = densi, and t-mn-rature of reference

= a .-.Y:,cient J.k-ined by

pC + sCs(1 )
pC


:. ,,4 a brief literature survey *:-x.ined in I it was suggested

that the friction function, b, can be approximated by:


b = 0.014 _-


In an isotropic medium the di:-.-- -.ion tensor, Ei.. can be ,. -:...-

as a sum of the i ,I,-:ic molecular diffusivity and the second order


where


(4)









symmetric mechanical dispersion tensor, as in


E = K6ij + Et4 (7)


where E4. = E* 6.. + (E* E*) uiuj/U2 = (E K) 6..
13 t _J k. t it

+ (E Et) ui u /U2 (8)


Here subscripts t and k refer to transversal and longitudinal

components respectively.

An assumption of singly diffusive convection is justified for a

multi-component system for moderate Reynolds numbers. A model suitable

for the description of the mechanical dispersion tensor in such cases

was suggested by Saffman [1959]. According to this model


t s + 2 Re (9)
v 2(s + 1)(s + 3) -) (9)


E* 2
S_ (s + 1) Re (10)
v 2(1 s)(s + 2)(s + 3) ()

where s is a power coefficient describing the dependence between the

velocity and the pressure drop (s varies between unity, for laminar flow,

and half, for turbulent flow).

The flow field model considered in this study while unperturbed

conditions prevail, is described in Figure 1. This is a saturated porous

layer of infinite horizontal extent bounded by two impermeable planes

located at bottom and top of the aquifer. Temperatures on bottom and

top of the porous layer are To and T -AT, respectively. Through the




















o ~r -7/ / / 77
^^^^JZ.J^^^^^^^^^Z


TEMPERATUF2

*OF ILE E


IMPERMEABLE

BO'JDARI


Fi,.,: 1. Schematical description of unperturbed conditions


'


Uo


. y









porous layer the fluid flows uniformly in the longitudinal, x, direction.

y and z are transversal and vertical coordinates respectively.

The flow field variables can be nondimensionalized as follows:


x = (x ut / )/d


u1 = uid/E
1 i t

t = tEt/d2


(p pogz)K xu0
vp E (1 + b ) Et


T' = (T T )/ T


U' = Ud/Et


E = E /Et
Ii ii t


(11)


where .i is a unit vector in the x direction. Substituting the dimension-

less variables of (11) in (1)-(3) and omitting the primes we obtain


-P-- RTn. + (1 + B)u = 0


+t i ax. ax. ij ax.
1 1 J


(12)


(13)


Here R is the Rayleigh number defined by


R = ATKd
v (Et(I + b)


(14)


The power s in (9) and (10) according to I can be calculated

through


s n(Re)
Zn(1 + b) + kn(Re)


(15)


However, (15) can be utilized only when Re 2.77 if b = 0.014 Re.










For ...- turbed conditions (12) and (13) yield


u = 0 (16)

= 0 (17)

T = -z (18)


p = P Rz2/2 (19)

E = E t E E = 1 Eij(i j) = 0 (20)


LT!C:- STABILITY ANALYSIS

Stability criteria of the flow field are determined through the

linear stability analysis. The flow field is subjected to small distur-

bances in the velocity (u, v, w), temperature (6), dispersion tensor

(e ij), friction function (W), friction power coefficient (c) and

pressure. Disturbances are very small and through the linear stability

analysis second order terms are negligible.

For the linear stability analysis the module of the velocity vector

in the perturbed flow field is approximated by


U = [(u d/Et + u)2 + 2 + w2]0.5 u 0d/Et + u (21)

Substituting (21) in (15), (9) and (10) we obtain expressions for the

principal components of the dispersion tensor. By applying (8), all

:,.nr..-nits of the dispersion tensor in the perturbed flow field are

obtained.










It is convenient to express the velocity components by utilizing a

scalar function Q as follows



2 22
u a v 2 w = -2 0 (22)


where v= 2 + ~ (23)
ax ay


Introducing flow field perturbations in (12) and (13), neglecting

second order terms and eliminating the pressure perturbation we obtain


v (Re + v2 ) = 0 (24)


y + V + 2 (25)
t ax2 ay az 1 2 axaz2

2 2 2 2
tt 3x- y z


where v + 2 + -
Dx2 @Y2 Z2
ax ay az

E E E*
S____ t
c u d T = u d (26)
0 0

Assuming vanishing values of 0 and Q on top and bottom of the

aquifer, these disturbances can be expanded by the following normal modes

(e,a) = (oe,3,) sin (Trz)exp [i(a x + ayy) + (or + io)t] (27)


where 0e and 2V, are constants, ax and ay are the wave number components.

For point of stability (or = 0) substitution of (27) in (24) and

(25) yields the following secular equation











a + 7T

a2 ia(cla2 Cr-2)

(E/Et )(a x/a V)2 +
X = 2
(ax/a ) + 1


R
2 2
Xa + Tr + iyW= 0


The minimal value of R satisfying (28) is the critical Rayleigh

r. -yields the following criteria of point of inst


ability


X = 1

a


ax 0

o -


ay = a


W = 0


(30)


convection cells are two dimensional rolls whose

parallel to the unperturbed velocity vector. Superposition

,-.,,-turbed velocity and the convection velocity leads to a

field.


axes are

of the

helical flow


Convection currents conducted in two dimensional rolls were also

obtained for the ordinary B6nard problem [Malkus and Veronis, 1958;
:i.luter et al., 1965] as well as for free convection in a porous layer

while mechanical dispersion effects are negligible [Straus, 1974].

However, in those cases only nonlinear stability analysis associated with

stability analysis of the steady convection motion yield such a result,

whereas, in our study the linear stability analysis indicated that

phenomenon. In our case, calculations concerning the anisotropy of

the di ':;;.'; on led to the conclusion that convection cells should be

two dimensional rolls.


(28)


(29)


numbe










Convection sets out in planes where the effective Rayleigh number

attains maximal values, namely, where the coefficients of hydrodynamic

dispersion attain minimal values.

Inertial effects associated with the invalidity of the laminar Darcy

law introduce the friction function, b, in the expression f'r R but do

not affect the linear stability analysis and predictions.


FINITE AMPLITUDE DISTURBANCES AND NONLINEAR STABILITY ANALYSIS

The effect of the convection motion on transport processes through

the aquifer can be predicted through the solution of the nonlinear equations

of motion and heat transport related to supercritical conditions.

The convection motion is two dimensional, therefore, the velocity

components can be expressed by the stream function as follows


n w =-- (31)
3- y az y


Substituting the finite amplitude disturbances in (12), (13) and elimin-

ating the pressure perturbation we obtain


R -+ V2+ = B(i,8) (32)


-Y +e 2 a- = H(p,O) + D(Eiz) + F(cij,e) (33)


where B and H are the friction and the heat advection spectra, D and F

are two parts of the heat dispersion spectrum, all of which are of the

form


B(ip,) ( i -) (34)
1 1x
M 1 x









H(Q,O) = ; ( 7 (35)


)- (e ) (!)
(iz x iz

F j -i x (37)
i i jx

Here g and cij are finite amplitude disturbances in the friction
funtio ant thpe disptson tensor respectively.
As 1.,.. as the convection velocity is smaller than the unperturbed
velocity the absolute value of the velocity in the flow field can be
.apco..imated by

U [(u d/Et)2 + V2]0.5 = (od/Et)(1 + x) (38)


where V2 v2 + w2 (39)

E t 2 V t 2 V t 4 V2
x = ( ) [1 ( ) 4- + (_) 8-
0 0 0
Et 6V6
-5 )6 V + + ] (40)


,-.L(:'lying (40), (6) and (15) we obtain

= [b/(1 + b)]Ax (41)

(1 + b)[zn(1 + b)] b[zn(Re)] (42)
(1 + b)[En(Re)]rn[(1 + b)Re]

Introducing (42) in (9), (10) and the dimensionless form of the
dispersion tensor we obtain after minor approximations









t 1
S. [A + s ( +
13 E s+2


d2

(Uod)2


E*
+ E -(3 -


1 1 3
si + s + 3 ij


1)(1 + ) +E t
t


tE
r


1
{- s
1 s


s)] 2s + 5 E
s] (s + 2)(s + 3) (-E*-
t,


- 1)} (1 -


As long as the Rayleigh number is not very high, x is very small and can

be approximated by the first term of (40), leading to the following

approximations


B = ~1V2

= (1V2


1.. = & V2
13 1


(44)

(45)

(46)


ii + 2 u u.
l.) 2 1 j


where


b
1 2(1 + b)


2
E ts (1 +
1 2(u d)2 (1 +

E*E
& + 1
2(u 0d)2 1


E- t
2 (u d)


b)[Ekn(1 + b)]- b[n Re)]
b)[Tn(Re)],n[(1 + b)Re]


1 1 1 E
s + 2 s + 1 s + 3 E


E 9
t 1)
t


x)}u u


(43)


Et 2
(u
0


(47)










If ( ") and ( .) are

conditions


,e = 0 at z = 0,1


aect to the '"..lowi-, ;.-'.. : -.,::.. boundary



( ,,


then the system (-) and ( .) can be solved by means of a set of trun-

cated ei tions '-:- -.si the finite ".-*litude disturbances.



as follows


= 1 q si 'pay) sin(qrz)
p 5 p 3q


p=0 p,q
q=1


cos(.-..y) sin(qmTz)


11 4 4)


(~7)


The calculation can be s .!-lified by using the complex variable

presentationn ..- sin(:,.y) and cos(pay) leading to


p


(5))




( ..)


T ei 'Y sin(-:.,,z)




e sin(qrz)






1q Pq


1'-
,,q 2 p1q


(53)


Z


00
a
I


;ded that


= 'V


p,q


7,q









The solution of the system (32) and (33) through the Fourier series
expansion means the determination of the Fourier series coefficients which

can be done by power series expansion. Such a method was first used by

Kuo [1961] for the analysis of the ordinary Benard convection. Palm et al.

[1972] and btobi, [1975] through different modifications applied this

approach for the analysis of free convection with no dispersion and

turbulence in porous media. Rubin and Christensen [1975] suggested some

guidelines for the utilization of such an approach when analyzing insta-

bilities induced by salinity gradients in a saturated porous layer. The

advantage of this method lies in its simplicity.

According to (30) the inception of convection motion is of the

marginal instability of exchange, namely, steady convection follows point

of instability. Therefore, the first term in (33) vanishes when steady

convection is attained. In such conditions the coefficients pq and

e of the Fourier series expansion can be expressed through a power

series expansion as follows:

N (n) (n) n(55)
(p,q p,q n=1 pq (55)

where n is a small parameter defined by

n = [(R-R )/ ]0'5 (56)


The Rayleigh number is also expanded by a finite power series as

follows:

R = R + Rs n2j (57)









where = RS/(1 n S = N/2 (58)

flow r-i d perturbations like v, w, g, and .ij as well as
..tion spectra H, B, D, and F can be also -. .. i in double .-...'ier

and series F ions for the coefficients H(n) Bn) D(n) and
,C o ,p,q' pq D p,q
F(n) are .-.ted in the ....dix.
pq
If the anal ..s is L:o-i..;:., for N = 1, it reduces to the linear
.I.> *.i l "' i c .. 1 ,.

number as given in (30).

There is a .,.: 'tive relationship between increases in '.l,,eigh
numbers and wave numbers. H;v.vr, taking the assumption that under
itical conditions the wave number remains constant does not signifi-

cantly a .*.-.. -Jictions of transport phenomena for quite a wide range
of -".,,leigh numbers [Straus, 1974]. Such an assumption is not required
f.. the method used here but considerably simplifies the analysis.

Substitut' .,:i the series p.,,.ions in (32) and (33) we obtain

nO (n-2i) 2 2 2 (n) +(n)
R pir e + Res pi- e(n 2i)' + i (p +q )pn + B = 0 (59)
pq os i=l pq ,q pq

2 2 2)(n) H(n) + D(n) + F(n) = 0 (60)
pq p,q p,q p,q p,q


The '*..-. ,ions q) and e are generated by superpositions of trigonometric
4.0-...ons, For n = 1 the only nonvanishing coefficients are'T1j. Therefore,

only -.:. Tcients Y (n) and e(n) with even values of 1pj+q have nonvanishing
p,q pq
values. .-.. the values of the subscripts Ijp and q are smaller or
equal to n. Coefficients with other subscripts vanish. According to (60).








(n) (H(n) + D(n)+ ,F(n) /Ti22) (61)
o,q o,q o,q o,q
For p00, (59) and (60) yield

(n) [(P2 2) -4] + Rs (n-2i) + B()/(npR )
p,q 2 2 R i=1 p,q p,q o
4pr(p ( q) 0

(H(n) + D(n) + F,)/[pq2)] = 0 (62)
p,q p,q p,q
For p=q=l
e(n) + E 0(n-2i) + B(n+2)'/R 2(H(n+2) + D(n+2) + F(n+2))/R (63)
1,I i 1 1,1 l /(os)-2 1 1 1,1 s (63)


Through (61), (62) and (63) all the coefficients ,(n) and O(n) are
p,q pq
determined. According to (63) the coefficients T1 n) (n+l) and e0n+2)
must be determined simultaneously. However, by simple arrangements e(n)
can be expressed explicitly. Such a procedure avoids any trial process.
Calculations of the mean horizontal temperature and the Nusselt
number followed the determination of the series coefficients. These para-
meters are given by


T = -z +Z (n) [sin(qfTz)]nn (64)
n=l ql1 o,q
N n-1
Nu N1 n r (n) + (il ) sq T(n-i)
Nu 1 no,q 1 q1= pq=- q,s=1 pq -ps
n=l q=l1

=3 n-1 i-1 = qs (j) 'ij) (n-i)
*F .E Z- E qs T T ea (65)
1 i= 2 j=i p,k=- q,k,s=l p,q k, IP+kI,s


Calculation of Nu determines the convergence of the method and the
termination of the series expansion. Through the calculation, presented in
the next section, we took N = 6, 8, 10, 12, 14 and 16 according to the
variation in the Nusselt number. If Nu varied by less than 2% as N was
increased from N to N + 2 then the expansion was terminated.










F- .! I 'J .ION

.. rding to (44)-(47) convection effects are not determined only 0L.

the .leigh number but also ". the ', I:1ds (Re) and Prandtl (Pr) numbers,

as well as :. the ratio between the *.'. .v: layer thickness and the charac-

teristi: ; K size (, .'.3.) whole .: 1 .:is .:,. a.- is cable

only when mechanical dirvs.sion is at least *..'..rable with the molecular

dii ion, s relationship is determined by the magnitude of Re and

It seems that if Pr ? 1, which is reasonable for practical pur..,..:

of hydrol ., [Somerton, '"-A], then dispersion effects are significant

even f.a.. Re Y 3. However, if the Prandtl number is smaller than unity

(if the solid fraction is a -.' conductor) then dispersion effects

become si '..ficant at hi.h-.' Reynolds numbers.

The : t of Re, Pr and d/d on the convection :w.n 1.j-.:.. is intro-

duced in the analysis through S1, &1 and &2 (and &3 = a. + 2) which

determine 2,---.ts of turbulence and dispersion induced by the convection

motion. Fi .,.- 2 demonstrates :-;,,.c. in these ,:c',: clients due to

variations in *, Pr and d/d when the porosity is 0.4. The *:- t; of Pr

vanishes ".. high Oldss numbers as in such cases mechanical di *..-..ion
_- .ts tr'- -, esses more than the molecular diffusion. .i,.- the

':- .ds number increases, according to model (9) and (10), the mechanical

di- ,-..ion tensor becomes more and more isotropic. This phenomenon leads

to a reduction in the value of Ia. Turbulent Jf5 :':. become more and

more si -ificant when the =-- .lds number increases, leading to a positive

relationship between increases in the Reynolds number and the coefficient

51'










According to (63) turbulent effects introduced through the friction

spectrum B amplify dispersion effects. However as 2 decreases when Re

increases, it was .*. i through the calculation that the net effect of

an increase in Re led to a minor reduction in the influence of the

spectra B, D and F.

According to Figure 2 the main parameter that determines the effect

of the spectra B, D and F on the convection phenomenon is d/d The

coefficients &a' &2 and I are aii,,:.. inversely proportional to the

square of this parameter. It should be mentioned that the absolute

value of the coefficient r, is usually much smaller than al. Hence

changes in s due to the convection motion are almost negligible. We may

usually neglect terms d:e!- ,ing on this coefficient in expressions (47).

As an example we present in Figures 3 and 4 descriptions of mean hori-

zontal temperature and Nusselt numbers for an aquifer when p = 0.4,

Re = 3 and Pr = 1.

Figure 3 demonstrates profiles of mean horizontal temperature for

various Rayleigh numbers and values of d/d As explained before, a

decrease in d/d increases the effects of dispersion and turbulence

induced by the convection motion. The convection motion increases the

value of the effective dispersion coefficient. However, this effect leads

to a reduction in values of temperature gradients at top and bottom of the

porous layer as demonstrated in Figure 3.






01




-



SFigure 2. -.cription of
-,^L /0.5 S B, .5,- and a2 vs. Re

5 various values of
-T.- and d/d (1=0.4).











0.5






-... ---- -. -- -- / _- -A r 4


3 .


I




























-1.0 -0.8 -0.6


1.0


0.8


0.6


0.4


0.2


-0.4 -0.2 0.0


T


Figure 3. Mean horizontal temperature profiles for
various values of R/R and d/d (=O.4, Pr=l,
Re=3).











At high Rayleigh numbers a thick region in the center of the r....

layer should achieve a nearly isothermal state in the mean. However, as

.. ted in Fi,,-. 3 such conditions lead to a positive temperature

gradient or reversal of temperatures. Such a phenomenon was also ident-

ified in the ordinary Benard problem [Kuo, 1961; Veronis, 1966]. Veronis

[1966] tried to explain the origin of such a strange phenomenon. However,

better choice of wave numbers diminishes this effect. Figure 3 also

demonstrates the creation of boundary layers on top and bottom of the

aquifer leading to invalidity of the continuum approach even for moderate

7.'leigh numbers if d/d is not very large.

Figure 4 presents variations of Nusselt number with Rayleigh number

for various values of d/d According to this presentation, the net

effect of the mechanical dispersion induced by the convection motion

leads to a reduction in the heat transport through the aquifer. Disper-

sion and turbulence induced by the convection motion act as a stabilizing

mechanism in the flow field (ironical interpretation).

However, this mechanism is more complicated than just a stabilization.

In the calculation it was found that the phenomenon was associated with

increased values of the higher modes of the series expansion, leading to

a reduction in the convergence of the series expansions.

Through the calculation, the maximal values of X and V were continu-

ously checked in order to examine the validity of approximations (44)-

(47) and to follow .!.. -,c in the Reynolds and the Peclet numbers. It








Nu


6-


5


Description of Nu vs.
R/R for various values
of d/d (4=0.4, Pr=1,
Re=3).


3


2-

I-


-pO

- ----p /P=1


) 2 4 6 8
) 2 4 6 8


Figure 4.


10
R/Ro










was .. :,. that i .' the -., R/R!0 10, T_, was much smaller T,:,; unity

which u;.. ifies utilization of (44)-(47) and yields only minor

in the "_-. Ids number.


1 IONS

The anisotropic character of the dispersion tensor leads to convec-

tion motions conducted in two dimensional rolls whose axes are parallel

to the ..,.bed velocity vector. By choosing a new definition :,

the Ray1l-..'( number, turbulence and dispersion effects can be introduced

in the linear stability analysis with no substantial complications in

the calculations and results.

Finite amplitude analysis for homogeneous boundary conditions can

be conducted by Fourier series and power series expansions.

The significance of mechanical dispersion and turbulence induced by

the convection motion .:ps ic, on Pr, Re and mainly on d/d P An increase

in Re diminishes the effect of Pr and reduces effects of dispersion and

turbulence associated with the convection motion. An increase in -'1d

reduces significantly e-~.': s of dispersion and turbulence due to the

convection motion. For d/d 102 these effects practically vanish.

--'. r, for smaller values of d/d these effects lead to a reduction in


The si. ly diffusive analysis presented in this article can be

.'.,lied when density gradients are induced by a single component or

when ..,!--. conditions and effective dispersion coefficients are

identical f. all components in a multi,,ri.one,! t system. The latter

conditions *....: ... .'! -ly satisfied when mechanical dispersion e:: s

are 1,.. :- than : molecular diffusion in the aquifer.








Expressions for coefficients of heat
advection, friction and heat dispersion spectra


T (i) p (n-i) (v=1,2;3)] + k[s (n-i) (v=2;1,3)]}
k, |jp-kks v V p-k, ,
(66)
Z(i) 2[ 2 (n i)
k, l +3)k(p-k)&2(p-k)2+2&i 2 (n-i) (v=1,2;3)]}
S(67)
(67)


s3= -q


(68)

(69)


Olp-k ,sv (v=2;1,3)= O p-kj,s2-e p-ki,sl -p-k1,s3


" (j) (i-j)
k,m=-m k,k m,h
t,h=1


{zh[4k(p-k-m)+(p-k-m)2 +3 s 2](n-i) (v=l

+[3k2hSv-Y2m(p=-k-m)][(ni) (v=3,4,5;1


(n) .4 n-l i-1 c,
pi,q- 16 i=2 j=1 k,m=-cc
9,h=1


(j)
k,9


,2,3,4;5,6,7)]

,2,6,7)]}


(i-j)
m,h


{th[(1 + &3)k(p-k-m)+&3(p-k-m) 2+&s 2][en-i) ,s(v=5,6,7;1,2,3,4)]

+m[2&1k2-2 2)(p-k-m)+&1k(p-k-m)2+&3ks2][O(n-i) (v=3,4,5;1,2,6,7)]


+h[&2k2+2 2k(p-k-m)-212-(& +& 3)mk][sE e(n-i) (v=1,3,5,7;2,4,6)]}

In (70) and (71) the subscript sv may obtain seven different values as

follows:


s3=q-z+h


s4=q+z-h


s2=q+(+h

s7=A-h-q


APPENDIX


2 n-1
2 1=1

3 n-1
-7 i=1


H(n)
Ip ,q


where



where


CO
k -
9,=1

k=-co


s =q-P,


B(n)
p,q


4 n-i i-I
1 4 i=2 j=1


(70)


(71)


s5= +h-q


s2=q+k


s, =q-y-h

s6=h-k-q


(72)











NOTATION

a

ax a y
a 0

b

B
B (n)
B B (
psq pq

ClC2
C

Cs
d

d
p
D

D D(n)
pq p,q
Et,E.
I3 i.j

E*, E

E*,E
EtEt

F

F F (n)
p,q p,q
g
H

H ,H(n)
p,q p,q
K



ni


wave number

(o-' ./i,!ts of wave number

critical wave number

friction function

friction spectrum

c.L:-f:.i-ents in -series-expansions for-B

coefficients defined in (26)

specific heat of fluid

specific heat of solid

porous layer thickness

characteristic pore size

part of heat dispersion spectrum

coefficients in series expansions for D

heat dispersion tensors (mechanical and hydrodynamical respec-
tively)

longitudinal heat dispersion coefficients

transversal heat dispersion coefficients

part of heat dispersion spectrum

coefficients in series expansions for F

gravity acceleration

heat advection spectrum

coefficients in series expansions for H

permeability

unit vector in the longitudinal direction

unit vector in the vertical direction









N total number of terms in the series expansion

Nu Nusselt number

p pressure

pO pressure at the coordinates origin

Pr Prandtl number (v/K)

R Rayleigh number defined in (14)

R critical Rayleigh number

Ros parameter defined in (58)

Re Reynolds number (pUd /v)

s power coefficient

S = N/2

t time

T temperature

T temperature at z = 0

u longitudinal velocity perturbation

ui velocity vector

u unperturbed velocity

U module of velocity vector

v lateral velocity perturbation

V module of velocity perturbation

w vertical velocity perturbation

xi coordinate

x,y,z coordinate system

coefficient of volumetric thermal expansion

1' 2' 3 coefficients defined in (47), (c3 = a1 + a2).
6 perturbation in the friction function











,1 coefficient defined in (47)

y parameter defined in (5)
6.. Kronecker's delta
13
AT difference in temperature between bottom and top of the porous


Eij dispersion tensor perturbation

n small parameter defined in (56)

6 ..:mre. ture perturbation

06 constant A. .ined in (27)
(n)
0e e e( ; i>*:,'cients in series expansions for e

K thermal diffusivity of saturated porous medium

x parameter defined in (40)

viscosity

v kinematic viscosity

C perturbation in s

C1i ,-.:.-ficient defined in (47)

p fluid density

p density at z = 0

ps solid density

or parameter expressing amplification of small disturbances

Sporosity

x parameter defined in (29)

9 stream :iuc.tion

T p ?T(n) coefficients in series expansions for i
pq, p,q' p,q
om parameter expressing oscillations










scalar function defined in (22)

1 constant defined in (27)










REFERENCES


Kuo, H. L., "Solution of the Non-linear Equations of Cellular
Convection and Heat Transport", J. Fluid Mech., 10, 1961, pp. 611-
634.

' Ikus, W. V. R., and Veronis, G., "Finite Amplitude Cellular
Convection", J. Fluid Mech., 4, 1958, pp. 225-260.

Palm E., Weber, J. E., and Kevernold, 0., "On Steady Convection
in a Porous Medium",J. Fluid Mech., 54(1), 1972, pp. 153-161.

.'r n,, H., ',;. the Analysi- o*f Ce' ,ular Convection in Porous Media",
Int. J. Heat & :t.,Trans., 18, 1975, pp. 1483-1486.

Rubin, H., "Onset of Thermohaline Convection in a Cavernous Aquifer",
to be published in Water Resourc. Res., 1976.

Rubin, H, and Christensen, B. A., "Convection Currents Associated
With Hydrodynamic Dispersion in a Porous Medium", 16th Congress of
-.-, Sao Paul, Brazil, 1975.

.*f'',- P. G., "A Theory of Dispersion in a Porous Medium", J.
Fluid *!_., 6(3), 1959, pp. 321-349.

SchlUter, A., Lortz, D., and Busse, F., "On the Stability of St-.dv
Finite Amplitude Convection", J. Fluid Mech., 23(1), 1965, pp. 129-
144.

Somerton, H. W., "Some Thermal Characteristics of Porous Rocks",
J. Petrol] Tech., 10, Note 2008, 1958, pp. 61-65.

Straus, J. M., "Large Amplitude Convection in Porous Media", J.
Fluid Mech., 64(1), 1974, pp. 51-63.

Veronis, G., "Large Amplitude Benard Convection", J. Fluid Mech.,
26(1), 1966 pp. 49-68.








CHAPTER 5

SEMINUMERICAL APPROACH FOR THE MATHEMATICAL
MODELING OF SINGLY DISPERSIVE CONVECTION
IN GROUNDWATERS


INTRODUCTION

Considerable interest is now focused in the State of Florida, as well as

in other locations in the United States, for the possible utilization of deep

saline aquifers for waste disposal. Vernon (16) delineates the properties of

the deep zones of the Floridan aquifer which makes this stratum available

for such application. Henry and Kohout (2) mention the fact that in a thick

system, like the Floridan aquifer, the effect of geothermal activity should

be considered, too. One of these authors has postulated in previous articles

(3,4) that geothermal activity induces groundwater circulation in the Floridan

aquifer. Singly diffusive convection is the convection motion induced by a

single dissolved component (e.g. temperature or salinity) in a fluid layer.

The hydrodynamics of this phenomenon in a saturated porous media has been

studied while, in most instances, assuming that the fluid is initially at rest

(7,8). However, groundwaters are generally subject to hydraulic gradients

leading to slow, effectively horizontal flow of the subsurface water. This

movement through a formation, similar to the Boulder Zone in the deep saline

region of the Floridan aquifer (1,4), leads to intensive mechanical dispersion

of heat and soluted materials in the aquifer, as well as inertial effects,

as demonstrated by invalidity of the laminar Darcy law. In such a system

molecular diffusion effects are usually less significant than mechanical

dispersion. Convection phenomenon under such conditions may therefore, be

called dispersive convection. A system is subject to singly dispersive convec-

tion if one of the following criteria is satisfied: a) density gradients are








introduced by a si: le component; or b) all molecular diffusivities of the

dissolved components are much smaller than the mechanical dispersion, and all

of them have the same boundary conditions. The objective of this article is

to present a rather simple method by which flow conditions in such an aquifer

can be simulated.



The basic ,:-....:,tions applied for the analysis are the equations of

Sitnuit-. motion, and ,-'.persion, subject to the Boussinesq a.:-:.ima.i,.:.

you = 0 . . . . . . . (1)

vp + pg + (1+b)u = 0 . . . . . (2)

y +- + O vT = v7 (E VT) . . . . . (3)

The .-.,-rfcient y appearing in Equation 3 is defined by:

Y = [pCwV + PsC )] . . . . . (4)

In a brief literature survey, presented in a previous study(10), it was

s5'..,= ted that the friction function, b, appearing in Equation 2 can be

approximated by:

b = 0.014 Re . . . . . . (5)

It is assumed that the fluid density depends linearly on the dissolved

....- t, which is the temperature, as follows:

p = po[I- (T-T )] . . . . . . (6)














zf
LmA /T To -AT///////////
ZZ//


T///////T


TEMPERATURE
PROFILE \


IMPERMEABLE
BOUNDARIES
I


7777
x


Schematical description of the aquifer with no
convection motion.


0


Figure 1.


9 9 0


y








Sow .eld model considered prior to the inception of the convection

motion, as presented in Figure 1, is a cavernous aquifer consisting of a porous

layer of infinite horizontal extent. It is bounded by two impermeable planes

on which the -,, ..ture is constant. Through the porous layer the fluid flows

uni.. -ly in the 1 ... tudinal x direction. The transversal and vertical

coordinates are y, and z, respectively.

A movi-.; coordinate system is applied with the velocity u /y in the

x direction, and .' ider

d, AT, e et/d, d2/et, et(1+b)/K . . . . (7)

as characteristic length, temperature, dispersion coefficient, velocity, time,

and :-, respectively.

In such a manner the following dimensionless basic equations are obtained

(in .. :. :,ions all variables are dimensionless):

v u = 0 . . . . . . . . (8)

vp RTn + (1+3)u = 0 . . . . . . (9)

@T 4
y + u VT = V7 (!ovT) . . . . . .(10)

the Rayleigh number, R, and the variable a, are defined as follows:
v t (ATKd


In an isotropic medium the dimensionless dispersion tensor can be



E = E t + (E Et)u u/U . . . . . (12)

As long as there is no convection motion Equations 8-10 yield:
u = 0 = 0 T = -z p = P Rz /2

S =E =E =1 E =x E = Eyz =0 . . .(13)
xx E yy zz xy xz yz .(.. )








THE FLOW FIELD STABILITY

The stability of flow conditions presented by Equations 13 can be

determined by analyzing the growth of small disturbances in the aquifer.

These are disturbances in the velocity, u, temperature, o, dispersion tensor, E

and pressure, P.

Second order terms depending on these disturbances are negligible.

It is convenient to express the velocity vector through the scalar vari-

able o as follows:

U = where = nx vx . . . . (14)

By substituting the small disturbances in the basic equations, and applying

the boundary conditions

o'v20 = 0 at z = 0, 1 . . . . . .. .(15)

(where v2 =

we obtain:


ylV2 RV2 2 V2d R[(E -1)/u ][V(V )]*k

+ R(aEt/au )[V(Q.-1)] . . . (16)

where
2 2 2
7v = v-V 7d = (E -1)P: VV + V . . . .. (17)

By expand-i in the foll-winq normal mode,

[sin(7rz)] exp[i(axx ayy) + (ol+ia2)t]. . . . (18)

we obtain for the point of instability (a1 = 0, and minimum value of R)

ax = 0 ay = a G2 = 0 ao = R = 4T2 (19)

Hence the instability is of the marginal of exchange type. Convection

cells are two dimensional rolls whose axes are parallel to the unperturbed

velocity vector. Superposition of the unperturbed motion and the convection

motion yields a helical flow field.








_ .IS OF THE STEADY C-,i ..TCTION

In the .- ous section it was -, .--;" that convection cells are two

dimensional rolls. It is convenient to :.'; .: the velocity by the stream

.,,tion
-> 4-
u = x 7 = -[ n, vs] . . . (20)

where the square brackets s ,-, lize the box product.

Substituting Equation 20 in the basic qur.*tions, and eliminating the pressure

perturbation, we obtain:

R[~ ?, VT] + V2 = B( . . . . . (21)

DO 2
-Y-1t + ve [n, v ] = H( 9,e) + D(c*) + F(, ). . . (22)

where T is the finite amplitude disturbance in the dispersion tensor.
.iables B, H, D and F are nonlinear terms defined as follows:

!, 8) = -v.(ev ) . . . . . . (23)

H( 0, e) = -[ v v ] . . . . . (24)

D(to ) = v-(Tn) . . . . . . (25)

F(e, e) = -v ( voe) . . . . . . (26)

As long as convection velocity is smaller than the unperturbed velocity

the absolute value of the flow field velocity can be approximated by:

U = (1 + A 1/2X2 + 1/2x3 . . . . (27)

where

S = V2/(2u ) V2 = u . . . . (28)

If x << 1 terms depending on high orders of A can be neglected. In

such conditions, the expressions for g and the dispersion tensor can be

. :. .-. i mated by:

S= 2 . . . . . . (29)
e = alV I + e2uu + a3(ut + zu) + a4V 2t. . . . ( )







where

1 = b/[2(l+b)u ] . . . . . . (31)

aI = (1/2Uq)(aE t/au) a2 = (E,-l)/u2

a3 = (E-)/o 4 = (-E)/u + (1/2u )[(E -Et )/u ] . (32)

The system of the differential equations 21 and 22 is subject to the

following boundary conditions:

4, e = 0 at z = 0,1 . . . . . . (33)

Such boundary conditions are simple according to the definitions presented

by Orszag (9). In such cases, accurate simulation of incompressible flows

can be obtained by spectral methods.

Assuming that + and e are periodic in the horizontal direction, these

variables can be presented by sets of truncated eigenfunctions as suggested

by Veronis in similar studies (17,18).

0 = ? ^ sin(pay)sin(qrz) . . . . (34)
p,q=l p,q

e = o e cos(pay)sin(qrz). . . . . .. (35)
p=0 pq
q=1

The calculation can be simplified by using the complex variable presentation

of sin(pay) and cos(pay) leading to

= -i i p [sin(qfrz)] exp(ipay). . . . (36)
0=-W p,q


o = 2 e [sin(qrz)] exp(ipay) . . . . (37)
p=-o
q=1

provided that






Pages
Missing
or
Unavailable








The convergence of the expression for the Nusselt number also determines

the truncation of the Fourier series expansion (the value of N).

NUMERICAL CPACLJ!ATIQNS

Experiments concerning heat transfer characteristics of porous rocks

(e.g. 5, 6) showed that the phenomenon of heat dispersion due to the fluid

movement in the stratum is very similar in nature to the characteristics of

solute dispersion in porous media. We may adopt, then, models of mechanical

dispersion, which are available in the literature, for the quantitative

evaluation of the effect of the convection phenomenon on the intensity of

transport processes in the aquifer.

Saffman (11) suggested the following expressions for the coefficients

of dispersion.

et 1 s + 2 Re
e 2Is + 2 Re) . . . . (42)



ej + (s + I)2 Re
v Pr 2(l-s)(s+2)(s+3) ( 43

where s is a power coefficient describing the dependence between the flow

velocity and the pressure gradient. In a previous article (10) it was suggested

that this coefficient can be calculated through the following expression

s = in(Re)/{tn[(1+b)Re]} . . . . . (44)

Equation 4+ :an R= armlied 'r T > 2.77. Through the calculation, it was

assume that Pr = 1, which seems to be a reasonable value for limestone and

dolomite aquifers (13).

The convergence of the numerical integration was very moderate. Several

approaches were applied to speed convergence of the calculations. These were:

(a) Each calculation was conducted in two steps, in the first step we took








d (which leads to :;- '=0), and, ..:' obtaining results :'- such

conditions the actual value d/d was introduced into the ~n ram; (b) .- -Its

obtained with lower values N 1 ., !/or R were used as initial -7.-.;ities

higher values -: these parameters. As a criterion for *.' state, Straus

(14l5) ted. in such calculations, to take

{ (' /dt)/e }max < 10-4 .()

in our calculation, it was "::-.j that such a criterion is ... conservative.

11y when ,1 a criterion of one hundred times less conservative,

vari. ,ons less than 1% were obtained in the Nusselt Iw .-, T. :.;, the
., ..- process, deca :ng oscillations were detected in the value .. the

7. t number.

Fi ,. 2 presents variations of the J;:;-.lt number with -:.leigh number

various values :, d/d P These are the maximal values of Nu obtained when

S i the value of the wave number (12). The wave number increased with

increase -: values of '/R o
nationss in the wave number induced minor -,h ... in Nu ,:.'- a ..'- ,.

ylei number

S... to Fi .. *. 2, mechanical di -..ion, due to the convection motion,

leads to a reduction in the intensity of transport processes through the aqui.'r.

calculation indicated that in small values of d/d the magnitude :" the

di .ion cents increases, when convection occurs, leading to a reduction

in the :- .itude of the ....-, :, gradients on top and bottom -'! the .q .-., .

This reduction in the i**..'... ture gradients affects transport ,.' :: more

than the increase in the value of the dispersion .-,. icients. It was also

found that lower values of d/d give rise to the higher modes of the Fourier

series expansions and reduce the *:: .,-. :-, of the calculation. For very







6-


4 --





d/d- //
220
/ ---d/dp= 20
---- d/dp= 10



0I ---- I I
0 2 4 6 8 10
R/Ro

Figure 2. Description of Nusselt number vs. R/R0 for various
values of d/dp (=0.4; Pr=1; Re=3).









small values of d/d the analysis fails to follow the physical -...r -.,.

even at moderate Rayleigh numbers, as the thickness of the boundary layers

developed on top and bottom of the aquifer is of the same order of magnitude

as the characteristic pore size. In such cases the continuum ;n-,,.:;.

applied th ..' the analysis is not satisfied. Through the calculation, the

maximal value of A was continuously checked in order to examine the validity

of the ,, -....,tions presented in Equations 27-32. It was found that for

,I. <10,- the maximal val.,: XH x was much smaller than unity.


CONCLUSI'.

Singly dispersive convection in aquifers can be analyzed by expanding the

flow field disturbances in eigenfunctions.

The ani:-.;, ..,..ic character of the mechanical dispersion determines the

plane in which convection motions are conducted. According to the analysis,

convection cells are two dimensional rolls whose axes are parallel to the

.' i:.- :ed velocity vector.

Analysis of steady state convection can be conducted by transforming the

equations of motion and continuity to a set of general first order dif,.:realal

equations that can be integrated through available subroutines.

The effect of mechanical dispersion and turbulence, induced by the convection

motion, depends mainly on the ratio between the porous layer thickness and

the characteristic pore size; Prandtl and Reynolds numbers have less significant

effects on the physical phenomenon.

The analysis of the singly dispersive convection presented in this study

can be :-il-ied when density gradients are induced by a single component, or when

boundary conditionL and effective dispersion coefficients are identical for

all components in a mul'icomponent system.







APPENDIX I EXPRESSIONS FOR THE COEFFICIENTS OF THE SPECTRAL FUNCTIONS


N-1 N-Ik|
H = H_ pq = (a/2) H z
Pq -Pq k=1-N s=1


+ k[w v p-kl, w (v=2;1,3)]} . . .. . ....... (46)


D q = D ,q
psq -pAq


2 N-1 N-IkI
= (Tra /2) Z E
k=1-N s=1


sT ks{[ (ctl +cx3)k(p-k)


+ 2ua (72/a2)w2 u2(p-k)2] p-k,w (v=1,2;3)]} . . ... (47)

where


wj = q-s


w2 = q+s


w3 = s-q


eIp-klWv (v=2;1,3) = e p-ki' p-klw2-e p kl 3 *.. ...... (48)


Bpq = -B
Pq -PAq


2 2 N-1
= (2a2/4) N
k,m=1 -N


N-jIk N-Im|
E h=
s=1 h=1


{sh[4k(p-k-m)+(p-k-m)2 + 3(2 /a2)w2 ][ kmw (v=1,2,3,4;5,6,7)]
v p-k-m,w ,
V


p,q -p,q


2 2 N-1 N-Ikl N-Iml
= (ra /4) z zE
k,m=1-N s=l h=l


{sh[(al+a3)k(p-k-m)+a3(p-k-m)2+al(T 2/a2)W2][eipkm l,wv(v=5,6,7;1,2,3,4)]

+ m[(p-k-m)(2 k2a22- 2s2) + 2k(a2 s2)(p-k- )2 + 23kw 2][ p -k-m2,Wy

(v=3,4,5;1,2,6,7)] + h[a2k2 +2a2k(p-k-m)-2as2 (72/a2) (al+a3)mk]


V='


wI = q-s-h

w5 = s+h-q


.,6)]} . . . . (50)


w2 = q+s+h

W6 = h-s-q


w3 = q-s+h


w4 = q+s-h


W7 = s-h-q . . . (51)


k,s m,h


Tk,s Tm,h


where


T ks fs(p-k)[E) lp-ki, w (v=1,2;3)]
v













a = wave number;

ao = critical wave number;

ax, a = wave number components;
b = friction function defined in Equation 5;

B = friction spectrum defined in Equation 23;

B =-coefficients -in the -Fourier series expanded for B-
q,q
Cs, C = specific heats of soil and water, respectively;

d = porous layer thickness;

d = characteristic pore size;

D = part of heat dispersion spectrum defined in Equation 25;

D = coefficients in the Fourier series expanded for D;
p,q
e et = longitudinal and lateral dispersion coefficients, respectively;

E Et = dimensionless dispersion coefficients (dispersion coefficients
divided by et existing prior to convection conditions);

E = dispersion tensor;

F = part of heat dispersion spectrum defined in Equation 26;

F = coefficients in the Fourier series expanded for F;
p,q
g = gravity acceleration;

H = heat advection spectrum defined in Equation 24;

H = coefficients in the Fourier series expanded for H;
P,q
I = unit matrix;

K = permeability;

9 = unit vector in the longitudinal direction;

n = unit vector in the vertical direction;

N = truncation .*ii:r;

Nu = Nusselt number;

p = pressure;


'^ ",! t h..I ,








Po = pressure at z = 0;
Pr = Prandtl number (=v/K);

R = Rayleigh number defined in Equation 11;

R = critical Rayleigh number;

Re = Reynolds number (=pUd /v);

s = power coefficient;

t = time;

T = temperature

u = velocity vector;

u = longitudinal velocity existing prior to convection conditions;

U = module of velocity vector;

V = module of convection velocity;

x, y, z = coordinates;

a = coefficient of volumetric thermal expansion;

i(i=1, ..,4) = coefficients defined in Equation 32;

8 = perturbation in friction function;

a1 = coefficient defined in Equation 32;

y = parameter defined in Equation 4;

= operator defined in Equation 14;

AT = temperature difference between bottom and top of the aquifer;

e = dispersion tensor perturbation;

e = 'mpe. !Ae p, LcuroL en;

^q,, pq = coefficients in the Fourier series expanded for o;

K = thermal diffusivity of saturated porous media;

x = parameter defined in Equation 28;

p = viscosity;










= kinematic viscosity;

= fluid density;

= fluid density at z = 0;

= solid density

= parameters expressing growth and oscillation of disturbances;


-p,q, Tp,q


porosity;

stream function;

coefficients in the Fourier series expanded for {;

scalar variable defining velocity perturbations;

constant defined in Equation 18;


P

p

ps

01, C02








REFERENCES


1. Burke, R. G., "Sun Whips Florida's Boulder Zone", Oil and Gas Journal,
Vol. 65, No. 7, 1967, pp. 126-127.

2. Henry, H. R., and Kohout, F. A., "Circulation Pattern of Saline Ground-
water Affected by Geothermal Heating as Related to Waste Disposal",
Proceedings of the 1st Symposium on Underground Waste Management and
Environmental Implications, T. D. Cook (ed.), American Association
of Petroleum Geologists, Tulsa, Oklahoma, 1972, pp. 202-221.

3. Kohout, F. A., "A Hypothesis Concerning Cyclic Flow of Salt Water Related
to Geothermal Heating in the Floridan Aquifer, New York Academy of Science
Transactions, Vol. 28, No. 2, 1965, pp. 249-271.

4. Kohout, F. A., "Groundwater Flow and the Geothermal Regime of the Floridan
Plateau", Symposium Geological History of the Gulf of Mexico Caribbean
Antillen Basin, Gulf Coast Association of Geological Societies Transactions,
Vol. 17, 1967, pp. 339-354.

5. Kunii, D., and Smith, J. M., "Heat Transfer Characteristics of Porous
Rock", American Institute of Chemical Engineering Journal, Vol. 6, No. 1,
1960, pp. 71-78.

6. Kunii, D., and Smith, J. M., "Heat Transfer Characteristics of Porous Rocks:
II. Thermal Conductivities of Unconsolidated Particles with Flowing Fluids",
American Institute of Chemical Engineering Journal, Vol. 7, No. 1, 1961,
pp. 29-34.

7. Lapwood, E. R., "Convection of a Fluid in Porous Media", Proceedings of
Cambridge Philosophical Society, Vol. 44, 1948, pp. 508-521T

8. Nield, D. A., "Onset of Thermohaline Convection in a Porous Medium", Water
Resources Research, Vol. 4, No. 3, 1968, pp. 553-560.

9. Orszag, S. A., "Numerical Simulation of Incompressible Flows Within Simple
Boundaries; Accuracy", Journal of Fluid Mechanics, Vol. 49, Pt. 1, 1971,
pp. 75-112.

10. Rubin, H., "Onset of Thermohaline Convection in a Cavernous Aquifer", Water
Resources :search, Vol. 1?, No. 2, 1976, pp. 141-147.

11. >-i. ,*, P. G., "A Theory of ;ispersion in a Porous Medium", Journal of
Fluid Mechanics, Vol. 6, Pt. 3, 1959, pp. 321-349.

12. Schluter, A., Lortz, D., and Busse, F., "On the Stability of Steady Finite
Amplitude Convection", Journal of Fluid Mechanics, Vol. 23, Pt. 1, 1965,
pp. 129-144.

13. Somerton, H. W., "Some Thermal Characteristics of Porous Rocks", Journal
of Petroleum Technology, Vol. 10, Note 2008, 1958, pp. 61-65.




Full Text

PAGE 1

WATER IiRESOURCES researc center Publication No. 39 Analysis of Storm Water Seepage Basins In Peninsular Florida By Hillel Rubin, John P. Glass and Anthony A. Hunt Department of Civil Engineering University of Florida Gainesville UNIVERSITY OF FLORIDA

PAGE 2

ACKNOWLEDGEMENTS. ABSTRACT .... TABLE OF CONTENTS CHAPTER 1: GENERAL INTRODUCTION Scope of the Study Objectives Methodology ........ Research Program . CHAPTER 2: FACTORS AFFECTING THE DESIGN OF URBAN DRAINAGE SYSTEMS IN PENINSULAR FLORIDA . Introduction ........... Legal Principles ......... Federal, State and Local Regulations Project Feasibility ...... Preliminary and Detailed Design. Conclusions. . . .. .... Appendix -Flow Diagram for Urban Drainage Design References .... . CHAPTER 3: ANALYSIS OF TRANSIENT GROUNDWATER FLOW FROM SEEPAGE PONDS ......................... Introduction ................... Calculation of Unsteady Seepage from a Pond .... Calculation of the Response of the Groundwater Table Discussion and Conclusions Notation ................... References . . . . iii iv 1 1 2 2 3 5 5 7 9 14 19 24 26 32 33 33 33 39 45 47 49 CHAPTER 4: THERMAL CONVECTION IN A CAVERNOUS AQUIFER 50 Introduction . . . . 50 Basic Equations . . . . 50 Linear Stability Analysis. . . . 55 Finite AmplitudeDisturbances and Nonlinear Stability Analysis . . . . 58 Results and Discussion . . . 65 Conclusions. . . . . . .. 71 Appendix Expressions for Coefficients of Heat Advection, Friction and Heat Dispersion Spectra ............ 72 Notation . . . . .. . .. 73 References .. . . . . . .. 77 CHAPTER 5: SEMINUMERICAL APPROACH FOR THE MATHEMATICAL MODELING OF SINGLY DISPERSIVE CONVECTION IN GROUNDWATERS 78 Introduction . . . 78 Basic Equations . . . 79 The Flow Field Stability . 82 Analysis of the Steady State Convection 83 Numerical Calculations .. . 86 i

PAGE 3

Conclusions ................ Appendix I Expressions for the Coefficients Spectral Functions Notation ................. References . . . '.' of the CHAPTER 6: NUMERICAL SIMULATION OF SINGLY DISPERSIVE CONVECTION 89 90 91 94 IN GROUNDWATERS . . . . .. 96 _'-_, _. '._' ,_--" __ '-! _____ -,-_._-,--_ 96 Basic Equattons_ .. _. __ ... ... _.. .. __ The Flow Field Stability. . . . 100 Numerical Calculation of the Steady State Convection 102 Numerical Results and Discussion 108 Conclusions .. 114 Notation . . . 115 References . 118 CHAPTER 7: SUMMARY AND CONCLUSIONS 120 i i

PAGE 4

ACKNOWLEDGEMENTS This report summarizes the research results obtained in 1975 and 1976 in the Hydraulic Laboratory of the University of Florida. The research project was sponsored by the Office of Water Research and Technology (OWRT) through the Florida Water Resources Research Center. The investigators are grateful to Dr. William H. Morgan, Director of the Florida Water Resources Research Center for his help through all phases of this investigation. Dr. Sylvester Petryk initiated the study and served as its principal investigator through October, 1975. His initiative and participation in this research are greatly appreciated. Chapters 2 and 3 of this report are mainly based on the Master of Science theses of A. A. Hunt and J. P. Glass. The investigators are very grateful to Dr. B. A. Christensen who served as the chairman of their respective graduate study committees. Dr. W. Huber served as co-investigator in this research and participated in the graduate study committees. His activity is very much appreciated. ii i

PAGE 5

ABSTRACT Rapid urbanization, currently proceeding in Florida, has resulted in significant problems with regard to both flood control and pollution abatement. The objective of the reported study was to search for the design procedures that may improve efficiency, safety and adequacy of drainage systems in peninsular Florida. Special attention was given to possibilities of recharging groundwater aquifers with excess storm water. Through such a system, partial solution to problems of inadequate potable water supply in some areas can be achieved simply as a by-product of flood control systems. In Chapter 2 of this report a conceptual design framework is presented. In developing this framework, a variety of disciplines involved in the solution of urban drainage systems in peninsular Florida are considered. Comparatively high levels of groundwater in many locations require consideration of two factors relating to the ability of seepage ponds to divert surface water to the aquifer. The engineer should study the ability of,the pond to seep water and he, should analyze the response of the groundwater table to the seepi,ng water. With respect to these two factors, analyses and solutions adapted for engineering application are presented in 3 of this report. The continuing reduction in the availability of potable water in coastal zones of Florida has prompted several communities to design recharging systems utilizing treated storm water as well as effluents. Chapters 4 to 6 of this report concern flow conditions in the Floridan iv

PAGE 6

Aquifer that should be considered in connection with this subject. Methods by which the particular phenomena associated with flow conditions in the aquifer can be evaluated are presented. These methods cover several spectral expansion approaches as well as a complete numerical approach. v

PAGE 7

SCOPE OF THE STUDY CHAPTER 1 GENERAL INTRODUCTION Since ancient times, people throughout the world have had to cope with periodic floods and inundations of lands and communities. Notably pragmatic solutions to flood control were developed in various locations. so utions were dependent on human imagination and ingenuity applied through available resources to a broad variety of situations and conditions. Even for a relatively limited area such as the United States, we cannot imagine a singular methodology applicable in every situation. Effective flood control depends not only on resources and materials provided by 'Mother Nature but, also on man's attitude toward the application of these resources. The very definition of the problem depends on the attitude of the people and their understanding. Peninsular Florida is a unique system in many respects. In the early part of this century, flood control projects were begun in earnest over an area frequently exposed to devastating hurricanes and extensive flooding. Since World War II, tremendous changes, inherent in the rapid development of the state, brought about improved techniques and better understanding of flood control. Presently, the state of Florida faces simultaneous problems of excessive storm water and limited water supplies particularly in coastal communities. Rapid growth and the influx of new people also resulted in reduced water quality, thus adding another factor to water management and flood control. This project has been initiated with the idea of simultaneously improving flood control techniques and enriching groundwater resources in peninsular Florida. 1

PAGE 8

iJb0t.L. 11 V t.:) The investigators found it necessary to direct their efforts toward three principle objectives: (a) General management and development of conceptual design framework. (b) Development of a simplified analysis for the evaluation of seepage ponds in diverting stormwater to groundwater storage. (c) Development of approaches in the analysis of migration of contaminants in groundwater due to natural conditions as well as situations induced by artificial seepage. METHODOLOGY With respect to the listed bbjectives, Chapter 2 of this report concerns the development of the conceptual design framework as related to urban drainage systems in central and southern Florida. Chapter 3 of this report concerns management of water quantities. This effort in the investigation involved determining the ability of seepage ponds to divert collected stormwater to the groundwater aquifer. Emphasis was given to development of simplified methods that can easily be applied by drainage design personnel. Chapters 4 through 6 concern mathematical methods that can be used for the analysis flow conditions in an aquifer similar to the Floridan Aquifer. Basic models are suggested and a variety of mathematical methods are checked from the point of view of applicability, efficiency and accuracy. These chapters form an introduction to analyses of water quality problems that have yet to be resolved in Florida.

PAGE 9

RESEARCH PROGRAM An initial research program was outlined by Dr. Sylvester Petryk in 1974 when he submitted the research proposal. The study, as envisioned in that proposal, would consist of the following five tasks 1) Review of existing literature 2) Modeling of precipitation, runoff, and storm water 3) Economic analysis and optimization of design 4) Design procedure 5) Experimental measurements. The research conducted during the past two years has been concerned with all of these areas. However, the program was modified with respect to the emphasis and the degree of effort devoted to each task. A review of the existing literature indicated that our efforts should be directed more toward improvement of the complete design framework and project management rather than involvement with particular design techniques or procedures. Chapter 2 of this report covers this part of our activity, which is partially related to each of tasks 3, 4, and 5 mentioned above. The design of seepage ponds in Florida often presents special problems because of high groundwater tables. These problems are the subject of Chapter 3 and are related to task 2 mentioned above. Problems of water quality have become increasingly important in Florida as well as in other parts of the country. One of our primary concerns in this regard was the effect that injection of impure water might have on the quality of Florida's groundwaters. Chapters 4 to 6 of this report are concerned with this problem. which is related to task 2 mentioned above. 3

PAGE 10

A series of field studies and measurements was conducted as suggested in task 5. Good results were obtained but the job was more or less routine in nature and we did not find it necessary to include them in this report. Our main objective in the field study program was to acquaint ourse-l ves with local protrlems and to-advise local p-eoplewith respect to these subjects. 4

PAGE 11

CHAPTER 2 FACTORS AFFECTING THE DESIGN OF URBAN DRAINAGE SYSTEMS IN PENINSULAR FLORIDA INTRODUCTION The demand for housing in Florida continues to climb as Americans seek a place in the sun. Retirement and tourism are rapidly displacing agriculture in developing and utilizing an attractive natural environment. Growth, since the end of World War II, has been phenomenal. Land development and home construction have become an integral part of the economy. These activities draw heavily on Florida's natural resources which are not limitless and, in many instances, are key elements in a sensitive environmental system. The demands of a progressive economy cannot be ignored, however, these demands can be adjusted to provide optimal use of a limited of supply. Concerted efforts are being made to provide both of a complex and dynamic system, and a reasonable balance of supply and demand. Water, one of Florida's most abundant resources, is a critical factor in this dynamic system. Early management concentrated on drainage and flood control to such an extent that damage to the system resulted. The past twenty years has seen a shift in management emphasis as water short-ages and environmental damage became more apparent. Contemporary manage-ment practices are undergoing rapid change. The idea of water as a scarce commodity has prompted basin-wide regulation of consumptive use and natural flow. Degradation of water quality has led to more stringent pollution control laws and the development of improved abatement techniques. 5

PAGE 12

All of these efforts have recognized and addressed urban development as a leading cause of problems in maintaining water quality and managing flow. Urban storm runoff has been found to be a source of water pollu tion that is of equal or greater magnitude when compared with any other --------------------------------identifiable source. Home building and_tanddevelopment require dr-ainage systems that accelerate runoff, consequently, creating or aggravating flood prone situations. Homes, streets, parking lots and commercial build ings retard or prevent infiltration of rainwater and necessary recharge of groundwater aquifers. In view of the current activity in water management and pollution control, the drainage engineer must be aware of the many facets of con tem-porary design and he must apply sound and systematic methods in development of his design. The natural environment of Florida presents a unique set of circumstances. Laws, regulations, procedures and techniques, although rooted in historic precedent and established standards, have all been developed with some consideration for these circumstances. The drainage engineer should have a thorough knowledge, not only of applied techniques, but of social, environmental and political impacts of development. He should understand recent changes in legal and governmental philosophies. The engineer needs to know sources of information, plan requirements, and the procedural aspects of government regulation and approval. Most importantly, systematic method in design synthesis permits the engineer to organize and evaluate his data, identify and augment weak points and effectively manage the design process. The objective of the present study is to review and suggest a conceptual design framework as related to the following topics: 6

PAGE 13

1. Basic legal principles associated with drainage, groundwater, land use and water courses. 2. Highlights and purposes of legislative and regulatory efforts on the federal, state and local level. 3. Elements of the feasibility study as the first phase in the design process. 4. Considerations and techniques applied to preliminary and detailed design of drainage systems. 5. Integration of the various design factors by means of a flow diagram. Finally, this discussion is oriented toward design factors particular to peninsular Florida. LEGAL PRINCIPLES Surface runoff is of legal concern for a variety of reasons. Damage done by concentrated runoff, either as a floodwave passing downstream or a backwater flooding of upstream lands, may result in injury to the affected parties. Alteration of flow directions, capacities or other drainage characteristics can and often does lead to environmental damage, e.g., lowered water tables, water pollution, damage to vegetation, erosion and accretion. Property damage from flooding of homes and businesses may be the direct result of poorly managed runoff and improper land use. Two basic principles of law concerning disposal of surface waters are 'the civil law rule and 'the common enemy rule' Under the civil law rule, the upland or dominant owner has an easement on the downstream owner for passage of surface runoff in its natural manner. The common enemy rule stipulates that the servient or lower owner may take measures neces-7

PAGE 14

sary to keep these waters off of hi s 1 and. Generally, Flori da courts have followed the civil law rule modified by 'reasonable use'. For example, the general rule regarding drainage into a natural watercourse states that a riparian owner may discharge surface runoff without regard to either the the_ rule'. right is subject to three limitations which have been imposed by the courts in varying degrees. These limitations include: 1. Drainage must be reasonable. 2. Drainage must not come from outside the natural basin. 3. The natural capacity of the stream must not be exceeded. The concept of reasonable use has been applied in cases involving both land use and water rights. In several instances) the courts have ruled that land use, adversely affecting surface or groundwater flows, incurs no liability on the owner as long as his use is reasonable and legitimate [2, 3J. Rulings have also been made in light of the riparian doctrine, requiring reasonably unimpaired or undiminished flow. In conclusion, the engineer must consider possible liabilities and legal consequences in his design of drainage works. Stormwatersystems designed to parallel natural discharges from a site are desirable from both an engineering and a legal viewpoint. Further. conflicts between water rights and legitimate land use have not been resolved in an entirely consistent manner. It is apparent that the body of law governing these activities and rights is not well defined and does not truly recognize the interrelation of land use and natural flow. 8

PAGE 15

FEDERAL. STATE AND LOCAL REGULATIONS Florida's rapid growth has produced heavy and often conflicting demands on existing water resources. The geology of Florida provides extensive groundwater sources that have been utilized in satisfying these demands. Wholesale acceleration and channelization of runoff and failure to provide adequate recharge to supplying aquifers combined with concentrated withdrawal from these aquifers has led, in several instances, to serious problems that cannot be resolved on a single situational basis. Total basin management has become the imperative solution. Comprehensive legislation at all levels of government has been enacted to provide a framework for management. Agencies with permitting and other regulatory powers have been created to affect solutions to apparent conflicts on a regional basis. The Federal Water Pollution Control Act of 1948 began what has become a massive effort to improve the quality of our national water resources. Subsequent'legislative efforts authorized federal assistance in research and development of methods of controlling pollution. Most recently, the F.W.P.C.A., Amendments of 1972 placed stronger emphasis on storm runoff as a source of pollution. Section 208 of this act is directed toward area-wide planning and management of both 'point source' and 'non-point source' discharges. Stormwater runoff, a non-point source, is to be evaluated on the basis of extensive data collection and monitoring in the field. Methods for reducing pollutant loading are to be developed at the community level and, subsequently integrated into a basin-wide management plan [lJ. Current technology in this area has not established a total 'cause and effect' relationship between runoff and stream quality. Loading rates 9

PAGE 16

are highly variable and depend on land use, storm duration and intensity, and seasonal weather patterns. This relationship or a reasonable approxi mation is singular in its importance with regard to effective modeling and comprehensive understanding of the problem. A second ar intel"e-st tothe drai-nage engineer stems from the National Flood Insurance Act of 1968 and the Flood Disaster Protection Act of 1973. Through this program, federally subsidized flood insurance is available to homes and businesses located in flood-prone communities. A prerequisite to qualification requires that designated communities adopt restrictive land use ordinances for areas subject to flooding on a frequency of once in one-hundred years (equivalent to a one percent chance in anyone year) [4J. The drainage engineer should recognize three significant points in these programs. First, the federal government has initiated comprehensive efforts to maintain and improve the quality of natural waters through strengthened programs at the state and local levels. Secondly, the flood insurance program is a substantial legislative effort directed toward better land use practices and flood control measures. Finally, these efforts actively involve the community and its citizens by promot ing improved local ordinances and related decision During 1975: Florida1s environmental agencies were reorganized and consolidated into two major agencies, the Department of Environmental Regulation and the Department of Natural Resources. Additionally, the state1s water management districts were redefined and were invested with broadened duties and powers. The Department of Environmental Regulation is responsible for enforcing pollution control laws and maintaining or 10

PAGE 17

improving water quality. Permitting of sewage treatment plants is a duty of this agency. In granting permits for point sources such as sewage treatment plants, the D.E.R. will review and actively consider stormwater control, treatment and disposal as part of the permit request [7J. Other highlights of agency guidelines include [7J: 1. Impact of stormwater discharges on the receiving waters will be viewed considering designated use of the waters, practical considerations and cost-effectiveness. 2. Design of stormwater management systems should include a variety of techniques in reducing impact. Sod filters, vegetated buffers, sediment traps, primary considerations in design. 3. Detention basins should be considered as essential in residential, business, industrial and highway development. Interim drainage systems, serving construction sites, should be anticipated and included. The Florida Water Resources Act of 1972 empowers the Department of Natural Resources to accomplish the conservation, protection, management and control of the waters of the state. In effect, this act makes all waters in the state subject to regulation [6J. This act further directed the environmental agencies to formulate a state-wide water use plan. Elements of the plan are being prepared at the district level and will be comprehensive in scope and objectives. Other activities of the Department of Natural Resources through the water management districts include permitting for consumptive use by all users, constructing and maintaining lands and works incident to management activities, permitting and regulation of wells. and management and storage of surface waters. 11

PAGE 18

It is readily apparent that the state is assuming the role of the riparian owner in many respects. Maintaining undiminished water quality and quantity is now a valid public concern as opposed to a strictly riparian right. In addition, the proposed water use plan will include -recomrnendat fonsast6 flood pTa-in z6nfng and prole-eli on lromfl 60a--hazards. Thus, the state is assuming an ever increasing responsibility over land use management and water resources which, traditionally, have been inherent in private land ownership. Subdivision regulations and zoning ordinances form development regulations and standards applicable at the community level. Further, these regulations outline procedures in obtaining community approval for a proposed development. They provide an orderly exchange or basis of communication between the engineer and the community. Frequent discus sions during the approval process, enable the engineer and community planners and officials to effectively demonstrate their needs and desires. Individual citizens may express their opinions during public hearings before the local governing body as a matter of state law. Subdivision regulations vary considerably across the state. Community resources and extent of development tend to impose practical limits on regulatory programs and requirements. Obviously, the most sophisticated regulations available relatively ineffective without adequate staffing and community support. Federal and state efforts are aimed at supplement ing community programs and preventing abusive and costly development such as has occurred in the past. Development standards tend to fall into two broad categories. The first is concerned with conditions or programs unique to a community. 12

PAGE 19

Such factors as flood control, aquifer recharge and salt water intrusion may require regulations or procedures written specifically for these conditions. For example, Orange County has recognized the need for con serving natural recharge areas. Their regulations outline special measures for maintaining natural high-infiltration in these areas. "Typical methods include the use of recharge wells, bottomless inlets, perforated pipe, grading to retard runoff, artificial seepage basins, swales in street rights-of-way, and utilization of natural percolation areas" [5J. The second category of regulations covers more detailed requirements. Minimum pipe size and material specification, street and lot grading, allowable overland flow distances, design methodology, allowable velocities in swales and ditches, drainage right-of-way requirements, and storm water details are typical items included in this category. Serviceability and ease of maintenance are important considerations in these specifications relating primarily to the 'hardware' of drainage. Routine maintenance of drainage systems and streets is a major cost item in municipal budgets. Flooded streets and overgrown ditches are common complaints from tax payers and every effort is made to minimize such situations in writing these regulations. In summary, legal principles governing surface flow and groundwater are derived from rights associated with private ownership. Conflicts arising from alteration of natural flow patterns or from consumptive uses have been resolved reasonable use of land and water. It is apparent that these principles are not adequate to resolve conflicts considering contemporary demands. As a result, regional or basin-wide management of water resources is being affected through combined legisla-13

PAGE 20

tive efforts on all levels of government. Preemption of private water rights appears to be justified in that the burden of flood damages and optimization of consumptive uses rests with the entire community. Finally, lawmakers have broad powers in preserving and protecting natural resources -i n-the-llub-lt-c s icrt-ibn l1a s been enacte a to improve land use, relieve flood hazards and reduce pollution of natural waters. PROJECT FEASIBILITY Project initiation usually begins with a discussion including the client and key staff professionals. Subsequent efforts depend to a great extent on several factors. Among these are the intended scope of the project, potential market, capital availability, anticipated scheduling of planning, engineering, construction and sales, legal considerations and immediate technical needs. Certainly, forty acres in an established residential area will have significantly different needs than ten-thousand acres of undeveloped and agricultural lands. The experience of the client or developer will also influence project initiation. An experienced developer often will have established project economics, scope and approximate timing prior to discussion with design professionals. Project feasibility may, however, include multi-disciplinary studies and analyses of the previ mentioned factors. One of the most fruitful approaches to evaluation and design is organization and coordination of all related efforts. Project planning draws each factor into the design process as it becomes pertinent. A systematic approach to project management helps to ensure optimization of 14

PAGE 21

technical effort and production and, that the quality of the effort will be thoroughly professional. Drainage design as an engineering process, is well suited to an organized approach. Design synthesis and coordination of effort in drainage engineering is graphically demonstrated in the flow diagram as shown in the Appendix adapted from Woodson [8J. Flow diagrams are an accepted and proven management tool. Their inherent flexibility makes them adaptable to a variety of situations and design conditions. The included flow diagram, presented in the Appendix, is divided into three distinct phases paralleling the balance of this discussion. These phases include project feasibility, preliminary des"ign and detailed design. Intensive data collection is the first step following project initiation. Gathering information on existing drainage patterns and both current and future land' use follows several avenues. A thorough under standing of the existing drainage situation is essential not only in determining required improvements but, also in demonstrating the effectiveness of proposed improvements. Standard sources of information that will readily provide data on drainage, land use and flood potential should first be investigated. Among these sources are: 1. U. S. Geological Survey for topographic maps, flood studies, groundwater studies, geologic investigations, well surveys, stream and lake gage records and general hydrologic information. 2. Municipal or County Engineering and Planning for subdivision regulations, zoning ordinances, aerial mapping, established "drainage networks and other information directly related to the site. 15

PAGE 22

3. Soil Service for soil surveys, drainage techniques proven in the area, recommended best use of the site and extensive technical expertise in drainage and hydrology. -guidelines in stormwater management. Specific informa-tion as to local environmental considerations and current practices. 5. Federal agencies for flood insurance and floodplain information, possible federal permits, design standards for qualification under federal loan programs and many other technical services. The reconnaissance survey of the site of a proposed development is particularly important in that information resulting from visual inspec tion is of primary importance and is generally available from no other source. Initial data may be summarized and plotted or noted on a site pl an to facil itate observati ons and compari sons duri ng the survey. Verification of land use, drainage patterns and improvements, both on the site and in the surrounding catchment, is an important objective. Other factors to be considered during I recon I are: 1. Attrar.tive or beneficial features on the*site that will contribute to aesthetics and environmental quality. 2. Hydraulic constrictions and control points, ponds, lakes, streams, springs and natural depressions collecting run-off. 16

PAGE 23

3. High water marks on trees and drainage structures, water surface elevations a groundwater or water table elevations. 4. Soil and areal distribution. Soil deposits such as ic muck and peat that must be considered in design and construction. ----receiving waters. Probable quality of current runoff and recei waters. 6. Identi cation of chronic drainage problems and existing improvements both on the site and in the immediate catch ment. During or following field investigation. the engineer should outline requirements for additional topographic surveys, soils investigations. stream gag ing, water quality environmental surveys and similar efforts. Methods and techniques for field investigations are dependent on site condi ons and the extent of proposed development. Detailed field notes and supplementary photographs are common to all surveys and provide documentation to support subsequent judgments, public presentations and technical discussions. Evaluation and incorporation of collected information involves simulation of existing drainage conditions, potential for delineation of drainage systems and land use planning. 'Rough' calculations or approximations should be made in evaluating the hydraulic and hydrologic characteristics of the sting system and the proposed improvements. The drai engineer should work very closely with those involved in land 17

PAGE 24

planning. land use, of course, will have multiple effects on any proposed dra i nage system. Compati bil ity of 1 and use and soil character; sti cs is of particular importance in Florida. Planners must be aware of probable water table elevations and seasonal fluctuations. Much of Florida cons;-sts of fine sands with-small -percentages-of silt in the surf-icial layers. Often, a semi-permeable layer of consolidated fines and organic material underlies surface sands. Infiltration of stormwater is severly retarded and 1 atera 1 or surfi ci a 1 flow may be very slow, parti cul arly on the marine terraces and in the 'piney flat woods' The drainage engineer should evaluate cost-effective measures in reducing high water tables and flat gradients and, he should advise land planners of his findings. Land planning should be coordinated to reflect drainage methods and alternatives best suited for the site, required rights-of-way, pond locations, flood elevations, fill sources and quantities, recharge sites and natural features to be preserved. Every effort should be made at this time to promote an aesthetically pleasing drainage network. An open ditch or floodway is patently obvious and prominent in appearance. The additional costs associated a landscaped and meandered channel are more than justified considering heightened sales potential and consumer idea 1 s. The preliminary engineering report or feasibility study reflects a combined effort in engineering, planning and project economics. Form and content for such a study are variable and. generally, depend on the scope of the project and desires of the client. It cannot be over emphasized that such a report is preliminary in nature and that subse quent efforts may contradict stated conclusions, cost projections and 18

PAGE 25

quantity estimates. Misunderstandings are frequent in this respect and the engineer should use care and discretion in making this point clear. Regulations, permitting and government coordination are often disregarded at this stage of the process. The feasibility study should include a thorough discussion of all aspects of the permitting and approval process. Preliminary engineering design provides the client with valuable information for project scheduling and funding. The feasibility report serves the developer in anticipating construction costs, cash flow, marketing and sales. F1nally, the engineer and client can selectively evaluate alternatives and determine preliminary design objectives. PRELIMINARY AND. DETAILED DESIGN Preliminary design should be approached considering overall site hydrology. Simulation techniques commonly used include the rational method, unit hydrographs, synthetic hydrographs, regression equations and more complex computer modeling techniques such as those developed by the Corps of Engineers and the Soil Conservation Service. The object of all of these methods is to systematically quantify runoff and to route these quantities through the system. Typically, the outfall point is identified together with available gage records or some estimation of appropriate water surface elevation. The flat gradients, so common to much of Florida, necessitate a thorough investigation and a reasonable approximation of downstream water surfaces corresponding to the events being simulated. Collector systems are usually a combination of open channels and storage basins that depend on 19

PAGE 26

both physical ground slopes and hydraulic head in maintaining positive flow. Central Floridais highlands may require grade stabilization structures where slopes are relatively steep. Peak velocities in unlined channels generally should not exceed two feet per second in Florida's -sandy soils. Less frequent storm events maybe perm; tted to generate somewhat higher velocities. Storage-elevation curves and stage-discharge relationships may be approximated in routing surface runoff through the system. Routing begins with generation of water quantities at the upstream end of the system and proceeds down the slope to the receiving waters. Peak discharges and water surface elevations, derived from the routing process, may then be used to calculate water surface profiles. Most systems generate peak flows over a relatively short period of time so that only peaks may be considered in calculating the profiles. Subcritical flow profiles must begin at the receiving waters and, based on the appropriate surface elevation, proceed Adjustments are made in channel geometry, storage, control structures and overland flow times until close agreement is achieved between the water surface profiles and the elevations generated by routing. This procedure begins in the preliminary design phase and is continually refined through the end of the design sequence. Many factors that are elemental to simulation will have considerable variation depending on the event being considered. Rain-fall intensities, for example, will generate system responses markedly different for various events. It is advisable to simulate storms of several frequencies and durations even though regulations may specify design based on one or two standard events. Practical considerations, such as allowable design time. must be considered in managing this procedure. 20

PAGE 27

Internal drainage systems direct runoff to the collectors. The rational method is almost universally used to design local drainage such as swales, culverts, inlet size and location, street slopes and lot grad ing. Home building and related land development require simultaneous evaluation of earthwork and drainage. Generally, surface slopes may be adjusted to minimize drainage improvements. Cutting and filling required to achieve ideal slopes may, however, result in excessive costs in earth moving. Once again, several iterations of site engineering may be required to achieve optimal conditions in detailed design. Incorporation of detention-retention ponds in the drainage system is desirable to attenuate increased discharges resulting from development. Ponding is a recommended method in reducing pollutant loading. Common practice is to provide sufficient storage within the system for the volume of runoff resulting from a given storm event. These events or design rainfalls generally require enough storage so that runoff from high frequency events, such as afternoon thunderstorms, is totally retained or detained for an extended period. Silts and other suspended solids have a chance to settle, evaporation and infiltration effectively reduce the volume of stored runoff, biological action and dilution reduce the concentration of pollutants, filtering devices prevent discharge of surface borne contaminants such as oil and grease, and trash and debris are prevented from entering the receiving waters. Natural storage areas, including cypress ponds and mangrove swamps, may be utilized in developing detention. These areas have the advantage of acting as filters for runoff, thus improving discharge quality. Care must be used to maintain minimum water levels in these areas and to prevent 21

PAGE 28

damaging pollutants from entering such a system that utilizes natural vegetation. Other considerations necessary to retention-detention ponds are: 1. Storage capacity recovery through evaporation; infiltration, outfJow,-modifi-edunde-Y'drainsand pumping. 2. Siltation or clogging in the pond. Two-stage ponds should be provided where sediment loads will be high. Periodic cleanout of settling basins in the first stage will be required. A subsurface berm or equivalent structure should separate the settling basin from the second stage. 3. Maintenance, safety, vegetation, normal depth, fill needs, aesthetics, emergency overflow, and drawdown capability are among the additional considerations. Systems of interconnected lakes and ponds are frequently employed in Florida to provide fill material. aesthetic appeal and recreational benefits. In addition, extensive open systems tend to lower surrounding water tables in chronically wet areas. All drainage systems should be viewed considering the effects of the maximum probable storm. Surface drainage systems may be inundated in some areas cat1sing streets and low lying areas to act as open channels or ponding areas. Utility services may be interrupted. High water in the outfall system may result in flood conditions for an extended period. In this regard, auxiliary floodways may be desirable. road elevations should be sufficiently high to allow limited continuous access, floor elevations near points of concentration should be higher than arbitrary 22

PAGE 29

minimums and off-site structures and embankments should be evaluated as flood hazards. In effect, preliminary design involves evaluation of a selected development plan. Design efforts should reflect sufficient detail such that very limited alternatives remain together with final details neces sary to construction. The flow diagram,presented in the Appendix,outlines elements of this phase through client and agency approval. Presentation of the preliminary plan together with supporting data, becomes part of a comprehensive review by the client and staff. These plans may then be submitted to local agencies for preliminary approval and required zoning changes. Review and discussion with state and federal agencies should be conducted on an informal basis. Results of these interviews and evaluations should be summarized for incorporation into detailed design. The third and final phase, detailed design, completes this sequence. Construction plans, specifications and detailed cost estimates are primary objectives. Revisions to-the preliminary design should be solicited from staff, client, agencies and potential construction contractors. The design manager must provide a thorough check of every detail and should coordinate with other members of the staff in locating and resolving possible conflicts. Almost invariably, insufficient clearances between sanitary sewerage, storm mains and water supply will become evident and necessitate change. Experienced contractors will often point out plan discrepancies, quantity errors and unfavorable site conditions during bid preparation or contract negotiation. During this phase, the engineer should coordinate closely with permitting and approving agencies in providing supporting data for selected design features. Unnecessary and 23

PAGE 30

costly changes, based on arbitrary interpretation or judgmental differences, may be avoided through foresight and preparation. All parties should be informed of required revisions and all plans should be noted accordingly. CONCLUSIONS Drainage design has become considerably more than estimating runoff and selecting culverts. The design professional is obliged to assume a broad viewpoint in developing his program and directing design efforts. In summary of this discussion, the following conclusions are presented: 1. The engineer should be familiar with basic legal principles related to drainage and how they have been applied in resolving conflicts and water rights litigation. 2. The engineer must be. aware of rapid changes taking place with regard to water rights and resQurcemanagement. Federal and state legislation in this area has shifted rights and responsibilities away from the landowner and toward the public entity. 3. Local regulations and ordinances are being re-writtenor supplemented in order to effect better land use. practices and flood protection measures. 4. efforts are being made to develop and implement improved pollution abatement techniques. Urban runoff, as a prime source, may be improved through a variety of techniques presently available to the designer. 5. Legislative mandates have been very direct. The engineering profession should take precautions in protecting their design 24

PAGE 31

prerogatives and, at the same time, align their efforts with these mandates. 6. The feasibility study is, possibly, the most important phase in design synthesis. This study develops a set of conditions influencing the entire project. 7. Florida's environment and dependence on well managed water resources emphasizes the efforts of the drainage engineer. He must actively participate in several areas of project development. 8. Detention-retention systems have become a necessity in designing improvements. Ponding or storage improves water quality, reduces flood potential and provides increased aquifer recharge. 9. Management of the design process and attainment of technical efficiency may be aided through the use of a flow diagram presented in the Appendix. Considerable flexibility allows this technique to be adapted to most projects. Finally, the attitude of the engineer must be receptive to change and inno vation. Bending to accomodate is not the answer. Embracing and assertive leadership is a trademark of the professional. 25

PAGE 32

APPENDIX FLOW DIAGRAM FOR URBAN DRAINAGE DESIGN Phase I: Project Feasibility Client Desires Staff Discussion Project Scope Project Economics Time Schedule Is there sufficient information to proceed? Existing Land Use Soil Survey Topographic Surveys Historical Records Interviews Future Land Use Studies Laws and Regulations Reconnaissance Survey Project Initiation Outline First Efforts Yes Data Collection Summarize and Catalog Data 26 No

PAGE 33

Are preliminary design objectives consistent and reasonable? Development of Internal Drainage System Development of Collector and Outfall System Identification of Critical Points Development of Hydrologic Model Development of Hydraulic Model System Evaluation Through Models Identification of Limited Alternatives Water Quality Evaluation 27 Establish Preliminary Design Objectives __ No Yes Phase I Phase II: Preliminary Design Development of Preliminary Design Prepare Preliminary Design Plans

PAGE 34

Have all sources been tried? Delineation of Existing Flow Patterns Identification of Constraints Simulation of Existing Drainage Conditions Proposed Drainage Schemes Proposed Land Use Are proposed land use and drainage plans workable? Cost Projections Quantity Estimates Anticipated Permits and Approvals Staff Review and Coordination Owner Review and Approval t >-..... ---No Yes Evaluation and Incorporation Develop Alternative Site Plans 28 >-...... ---No Yes Presentation of Feasibility Study

PAGE 35

Have preliminary design objectives been realized? Design Features Limited Alternatives Special Problems and Solutions Cost Estimate Ad9itional Data Needs Are final design features and objectives approved by client and agencies? Incorporation of Testing Results and Survey Data Safety Considerations Cost Optimization Coordination of Design and Drafting ">---... --No Yes Presentation of Preliminary Design Obtain Client and Agency Approval >-----1 ..... --No Yes t Phase II Complete Phase III: Detailed Design 29

PAGE 36

Design and Specification of Special Features Preparation of Detailed Cost Estimate Preparation of General Spec i-fi cations -. Checks and Revisions of Calculations and Drawings Project Plans, Estimates and Specifications Review Development of Detailed Design Issue Plans Reflecting Final Design I s the fi na 1 des i gn complete --..... and acceptable to the client? Permit Applications and Procedures Local Governmental Approval Contractor or Bidder Review Final Staff and Client Review 30 Yes Revision of Detailed Design Obtain Final Approval by All Parties i I I

PAGE 37

Are plans, specifications, estimates, and permits complete? Are approvals in receipt? >-------No Yes End Design Sequence Phase III Complete SymbOl s Key Flow Direction: Specific Activity or Information: General Activity: Executive Action: Decision: 31

PAGE 38

REFERENCES 1. "Executive Summary of Section 208 Program for Designated Area Federal Water Pollution Control Act Amendments of 197211, United States Environmental Protection Agency, Washington, D. C., 1974. 2. Linsley. R. and Franzin;, J., Water Resources Engineering, Second Edition-,nMcGraw-HilL, New York, -1964. _. _._. 3. Maloney, F., Plager, S., and Baldwin, F., Jr., Water Law and Administration, The Florida Experience, University of Florida Press, Gainesville, 1968, 4. IINational Flood Insurance Program", United States Department of Housing and Urban Development, Washington, D. C., 1974. 5. "Orange County Subdivision Regulations", Orange County, Orlando, 1974. 6. Tebeau, C., IISouth Florida Water Management" and IIEnvironments of South Florida: Past and Presentll, Miami Geological Society, Miami, 1974. 7. Thabaraj, G., "Regulatory Aspects of Storm Runoff Control", Proceedings of the Storm-Water Management Workshop, Florida Technological University, Orlando, February, 1975. 8. Woodson, T., Introduction to Engineering Design, McGraw-Hill, New York, 1966. 32

PAGE 39

CHAPTER 3 ANALYSIS OF TRANSIENT GROUNDWATER FLOW FROM SEEPAGE PONDS INTRODUCTION The ability of a seepage pond to divert storm runoff into the ground water aquifer depends on two processes: (a) the rate of seepage from the pond, and (b) the reaction of the groundwater table. transient phenomenon. The inflow to the pond. specified as a runoff hydro graph, causes variations in the depth of ponded water. Outflow is consid ered here to be entirely by infiltration, which varies not only because of changes in pond depth but also due to variations in such soil properties as hydraulic conductivity and storativity. The growth and decay of the groundwater mound in the underlying aquifer is also time dependent. These facts are well known among scientists but are usually not taken into consideration in the design of a seepage pond. The objectives of this study are to develop methods by which the designer will be able to estimate seepage pond effectiveness by considering both the rate of seepage from the pond and the variation in the groundwater table while taking into account the transient nature of the problem. CALCULATION OF UNSTEADY SEEPAGE FROM A POND The adequacy of a seepage pond is evaluated by a storage routing procedure, which is basically an account of the inflow. outflow, and change of stored volume over successive discrete time increments. The inflow is represented by the runoff hydrograph for the design storm on the area to be drained. For the purposes of the present exposition it 33

PAGE 40

will be considered a given function of time. The outflow, on the other hand, is dependent on time, depth of ponding, and the properties of the soi 1. Consider an initially dry seepage pond constructed in unsaturated -soi1. As it fills, the water begins to escape from it by vertical unsaturated flow, or infiltration. The most convenient way to describe this flow is by the formula of Green and Ampt [1911J as modified by Bouwer [1969J. It has been shown that this which was long thought to be purely empirical, is soundly based on physical principles [Morel-Seytoux & Khanji, 1974J and gives very good answers [Whistler and Bouwer, 1970J. Green and Ampt based their derivation on a simplified model of infiltration which treats the soil as a bundle of vertical capillary tubes. The vertical hydraulic conductivity and moisture content of the unsatur ated flow are considered constant, as is the capillary suction potential of the advanc i ng wett i ng front. App lyi ng Darcy I slaw to th is idea 1 i zed I flow, the infiltration rate is described as where W ::: K H + L + 1/In t L w = infiltration rate K t = hyarau1ic conductivity of the transmission zone H = depth of ponded water 1/In = capillary suction potential l = depth of penetration of the wetting front. 34 ( 1)

PAGE 41

The rate of advance of the wetting front is (2) where t = time f = the volumetric fraction of fillable pore space. The equation of continuity applied to the pond yields dv = AdH = I -fddt{AfL) dt dt (3) where v = stored volume in the pond I = inflow A = area of pond surface Af = area involved in infiltration. Introducing equation 3 into 1 and 2 we obtain, after differentiation, the following non-linear differential equation d 2 2 dL dAf dt2(L /2) + (fAf/A-l)dt + fL/A dt I/A = 0 (4) which can be solved numerically. However, in most practical cases the rate of change of H is about an order of magnitude smaller than the rate of change of L. In such cases an integration of equation 2 between ti and ti+l provides a simple discrete presentation of the problem, as follows t'+l J 1 t, 1 35

PAGE 42

The integration leads to (5) In dimensionless form this can be expressed as f1tKt AL r'+l L'+l/L + r+1/L. __ = _0. w( 1 1 1 "2 1). fL L l 1 + r. +1 / L 1 1 1 1 "2 1 (6) where r '+1 = 1 "2 H "+1 + 1 "2 l/In L. = value of l at time t. 1 1 li+l = value of L at time ti+l f1t = ti+l -t. 1 f1l = li+l -l. 1 H '+1 = 1 "2 (H. 1 + Hi+l)' the mean value of H. A chart that can be applied for engineering design purposes is shown in Figure 1. The curves are based on the solution to equation 6. The application of the approach of equation 6 was checked against experimental data reported by Weaver and Kuthy [1975J. In this experiment a seepage pond was constructed and filled with water at a controlled rate. The variation of the pond volume and area vs. depth is shown in Figure 2. The experimenters listed soil test data which led to the follow;n9 values of soil parameters. K t = 1. 2 ft/hr f = 0.2 l/In = 0.5 ft. According to the reported measurements, the rate of change of H was 5 to 10 times smaller than the rate of change of l. This would indicate 36

PAGE 43

L'lL c:1 w -....J 0.5 0.1 0.05 T 1 I 1 11111 I I I T F III I I I J I [Il ++ I l r 11l I11I1 I I11I1 I I III TIV J I II I i-III I 1 I I 111111 j= I 1111111 1 1 J#f118 I I I II -E r=0.5 _L, I r=0.2 LJ r=O.l Mill-i-' :1 /-e--Wf1JW' I l..k--:=::::' l,_ A/17" -+--HI I I /1 d' I r-tD III' l I R=R -iii I I fffii I III111 -l-d 0.01 Ir=o.o) A' IA.a2' .----1 1L1N +-r=O.02 I' f 0.005 I I I I iHF::fzbi711 + ill I I III H 0.001 I r =0. 0 1 I r=0.005 0.001 0.005 0.001 0.005 0.01 (L'ltK; )/(fL;) 0.05 'I I I I I q.l Figure 1. Solutions to Equation 6 for li+l/l; ::: [1.2 : I -ttl I II I r 0.5 1.0

PAGE 44

w co 6 5 +> 4 4-:r: le... 3 a z: o e... 2 1 o AV13RAGE..AREA (thousand square feet) o 1 2 3 4 5 6 7 I I I \ II I I I I I I I I / VI I I ./' I I I I i I It/ /1 i I i I I I i I I I I I I I I i I I 1/ I I i / Volume I I I I I I i IYi / I I I I \/ I I I I I I / V i If-i Area /! I / t I V/ i I I I 1 /1 I i I I V I I I VI I I I I I i I IIV'! i I I I I I r' I I I J l I I i I o 5 10 15 20 25 30 35 POND VOLUME (thousand cubic feet) I Figure 2. Depth vs. Average Area and Volume for Test Pond 8 Y I : I I. I I I I I I i I [ : I I I 40

PAGE 45

that the method utilizing equation 6 should give good results. The comparison between the experimental data and the theoretical prediction is shown in Figure 3. The agreement is, indeed, quite good. CALCULATION OF THE RESPONSE OF THE GROUNDWATER TABLE The reaction of the water table aquifer to vertical infiltration from a seepage pond is characterized by the growth of a groundwater mound. For analysis we consider a circular pond constructed in homogeneous and isotropic soil. The aquifer is underlaid by a horizontal impervious layer and initially has a horizontal free surface. When the wetting front of vertical infiltration reaches this free surface the mound begins to form. From this time the infiltratiori is assumed to continue at a constant rate. 39

PAGE 46

100 90 ,:t-80 +l ClJ ClJ 70 4U ..... 60 ..c ::s u -c 50 C rtI III ::s 0 40 ..s::: +l .j:::> 3 30 a a .....J I.L. z: 20 ...... 10 0 V Ot ::el'ved Je ptr, 0 / C31 cul a :.eJ nfl 0 .... 1 ---..... V / .-. i v tJ ( ( V 0.3/ '" ,.J / 0'( ;V/ V D r V 0 1'---, 1 9 10 11 12 13 14 15 2 4 5 6 7 8 3 10 9 8 7 .--.. 6 +l ClJ ClJ 4-........ 5 :t: I-a.. I.J.J 4 0 3 2 1 0 TIME (hours) Figure 3. Comparison of Experimental and Calculated Storage Routing I

PAGE 47

In a cylindrical coordinate system centered on the seepage pond as shown in Figure 4 the axisymmetric flow in the saturated mound is described by where nr ar + w n Ks = saturated hydraulic conductivity n = specific yield w = rate of vertical infiltration. (7) Both wand n are functions of radial distance, r. The value of specific yield is greatly reduced in the zone of infiltration under the pond. Based on the experimental evidence of Bodman and Coleman [1943] it seems reasonable to assign a specific yield to the'zone of infiltration with a value of 20% of the specific yield elsewhere. The rate of infiltration, w, is zero when r is greater than the pond radius, Ro' Noting that h = S + a, equation 7 can be rearranged: .' '" Ks (s + a) a2S + Ks (aSJ 2 + f.s + a) as + w at n arzn ar] rn l ar n (8) To make the equation dimensionless, introduce the following dimension less variables: Substitution of these variables along with some manipulation and omitting the primes we obtain: as a2 (2 ) n"IT arz-S /2 41 (9)

PAGE 48

s -,t h a r Figure 4. Groundwater Mound

PAGE 49

In this form, the non-linear terms are distinctly separated from the linear terms. They are treated differently in the numerical analysis. 2 = = 0 at r = 0 ar ar 2 (10) and S = 0 at r = 00 The initial condition is: S = 0 at t = 0 (11) To solve the problem by finite differences, we must set up an appropriate space time grid. The time increments are designated by ;IS and the space increments by j1s. The time derivative of S is approximated by a forward difference. (12 ) The non-l inear spatial derivatives are represented by central difference operators: a ar (S2/2) '" [S2,. J.+l 2S2 .. + S2. +11/2!J.r2 ,J J (13) (14 ) The linear space derivatives are the most influential terms in the equation. If the solution were to become numerically unstable, it would 43

PAGE 50

probably be because of these terms. Therefore, they are approximated by the mean of the finite difference representations on the (1 + l)th and the (i)th time rows. This;s the implicit method of solution developed by Crank and Nicolson [1947J which has good stability and convergence characteristics: '" ((Si+l,j+l Si+l,j_l)/(2Llr) + (S;,j+l S;,j_,)/(2Llr)] a2S [ 3rT '" 2 S;+l,j+l 2S;+1,j + S;+1,j_l)/Llr2 + (S. '+1 2S .. + S .. 1 )/Llr2 ) 1 ,.1 1 ,J 1 ,JWhen these finite difference approximations are substituted into equation 9, a tridiagonal system of equations is generated which can be solved by Gaussian elimination. (15 ) (16 ) The computer program developed from these finite difference operators was tested with several combinations of grid mesh ratio and grid size. 2 It was found to give stable results with a ratio of Llt/Llr 0.11. The mesh size chosen was Llr = 0.03. To simulate the boundary condition at infinity, the calculations were carried out to r = 25 at which point S was required to equal zero always. For comparison with the linearized analytical solution of Hantush [1967], the finite difference program was run with constant specific yield, n. The results of this run are shown by the dashed lines in 44

PAGE 51

Figure 5. They are practically identical to the results obtained by the method of Hantush. When n was allowed to vary as a function of r, the results were very much different, as shown by the solid lines in Figure 5. DISCUSSION AND CONCLUSIONS An analysis of the hydraulic operation of a seepage pond should be conducted in two phases. The first phase is concerned with the rate at which water seeps out of the pond by vertical infiltration through its bottom. The equation for vertical unsaturated flow developed by Green and Ampt [1911J can be used to describe the outflow in a storage routing process. The storage routing is simple enough to be done by hand when the soil is assumed to be homogeneous .. The Green and Ampt equation can also be adapted for use in soil where the hydraulic conductivity varies monotonically with depth [Bouwer, 1969J which is a very common occurrence. This calculation is about as simple as the storage routing procedure but the combination of variable ponded depth and variable hydraulic conductivity would be cumbersome enough to make the use of a computer or a programmable calculator desirable. The second phase of the analysis deals with the response of the water table to recharge from the seepage pond. The determination of adequate capacity of the aquifer requires use of the design charts for each specific case. The design chart itself is based on the numerical solution of the non-linear differential equation. The importance of the finite difference solution of the groundwater mound is that it makes possible the consideration of an axisymmetric variation of the specific yield. It would also be possible to include the effects of axisymmetric variation in other soil properties. 45

PAGE 52

+ Vl 1 14 P 0 3 I II ,p 0 3 j" .. ,. l,T-"II .. i-j--"I''I"'I,"lT--l .. .. I ............. ......... ........ :::: :: : .:.. =. :::'r::r: -. ....... '"T''' "l",r::I":r"T"r:r:1 < >:: :<::::-:' ::>< -:: :::::::::::::: '.: :: :::: :::: : :::: :: :: .. ::. III .. .... .. ...... ... 1 .......... ................ .. 1.1211, ;Vl,.,> ... !iP .. .. t ...... .. .... l.. vJ R . ... .j .... ./V.. .. ::j:. .... ", ........... .. :::+::: :::: :::: ::::,:::. p = -;: j' :::;'::.::<::" .::. :::.. .. : 1: V:: ::: :;::1:::: :::: :::: :::: :::::::. .... .......... I .. Ka .... I .. j./ ..... ............... .. ...... .. :::: :::: :::: :::: :::::::: s : .. ::::i ,':: :::: :::: ::.: ......... j ..................... .. .. ++ +--+-1 I 1 1 1 1 1 1 1 1, : 1 .... ( ,, .. tI,II, ...... .... ...... .. .. 1;, .. 1 ......... 1 = C3 ...... L.U ::I: ....... z ::::> en 0 ::E: Vl Vl L.U ....J Z 0 .... Vl Z UJ ::E: ...... ... /. :41"ll:::::n p= 1. 08 : : .!:: : ",A. J .. 1.:, I.. .. ". ..................... 1 .. 1 ...... .......... .. .j ............ .... ......... .. 1 ......... ......... 0... I ,.j.-' .... I .... t ... li ... ........ "I .... '.... .. .. 1 :: .... :: -i" ... \ .......... .. I .. 'j .. .. I .. ,. r! ... .". .. ...................... .. ......... .: .... :. : :::: To. :::: :::: ::J::: :". :::: :::: :::::.. :1,<: : ':l::: ::':tj : :::: : .. :,::,.':. :H:" re :::: ::.: :J::.f:: 1.04 ... : ;: ::::1.1': :::1:': j ,. ,I. M' .... I" ..... t .... I ..... i . Cpn$tant n ---1 02 .. Vbdable n .02 .05 .1 .2 .5 2 I I 5 10 DIMENSIONLESS TIME, tl Figure 5. Mound Height at rl = O,CQmparison of Constant .and Variable: ni

PAGE 53

It should be noted that the saturated hydraulic conductivities referred to in the two phases of the analysis are, very often, not the same. Infiltration is concerned with vertical flow while the groundwater mound is mainly influenced by horizontal flow. The difference between the saturated hydraulic conductivities in the two directions usually stems from horizontal layering of the soil. NOTATION a h H I L n p r r' S' t initial thickness of aquifer [LJ area [L 2 J area of infiltration [L2J volumetric fraction of fillable pore space thickness of aquifer [LJ depth of ponded water [LJ inflow to pond [L3J hydraulic conductivity of saturated soil [L/TJ hydraulic conductivity of the transmission zone in unsaturated flow [LIT] depth of penetration of the wetting front [LJ specific yield .dimensionless recharge intensity, radial distance [LJ dimensionless radial distance, r/Ro radius of seepage pond [LJ rise of groundwater mound above water table [LJ dimensionless mound height, S/a time [TJ 47

PAGE 54

t' v w r dimensionless time, stored volume in pond [L3 ] infiltration rate [LIT] 1/In + H [L] cap;ll QD-poten-{iiil-iifflelcfcaJ1acity[L] 48

PAGE 55

REFERENCES Bodman, G. B. and Coleman, E. A., IIMoisture and Energy Conditions During Downward Entry of Water into Soilsll, Soil Science Society of America Proceedings, Vol. 8, 1943, pp. 116-122. Bouwer, H., IIInfiltration of Water into Nonuniform Soil II Proceedings of the American Society of Civil Engineers, Irrigation and Drainage Division, Vol. 95, No. IR4, pp. 451-462, 1969. Crank, J. and Nicolson, P., IIA Practical Method for Numerical Evaluation of Solutions of Partial Differential Equations of the Heat Conduction Type II Proceedings Cambridge Philosophical Society, Vol. 43, 1947, pp. 50-67 Green, W. H. and Ampt, G. A., "Studies on Soil PhYSics. I. The Flow of Air and Water Through Soils", Journal Agricultural Science, Vol. 4, 1911 pp. 1-24. Hantush, M. S., I'Growth and Decay of Groundwater-Mounds in Response to Uniform Perco1ationll, Water Resources Research, Vol. 3, No.1, 1967 pp. 227-234. Morel-Seytoux, W. J. and Khanji, J., IIDerivation of an Equation of Infiltration, Water Resources Research, Vol. 10, No.4, 1974, pp. 795-800. Weaver, R. J. and Kuthy, R. A., Field Evaluation of a Recharge Basin, New York State Department of Transportation, Engineering Research and Development Bureau, Research Report 26, 1975. Whistler, F. D. and Bouwer, H., "Comparison of Methods for Calculating Vertical Drainage and Infiltration for Soils", Journal of Hydrology, Vol. 10, No.1, 1970, pp. 1-19. 49

PAGE 56

CHAPTER 4 THERMAL CONVECTION IN A CAVERNOUS AQUIFER INTRODUCTION ---------------------the author ana lyzed i nstabil ity-cri terla re-lated to onset of thermoha 1 i ne convection in an aquifer whose properties are similar to those of the Boulder Zone of the Floridan Aquifer. Such an aquifer is characterized by extremely large pore size and transmissivity leading to very intensive solute and heat dispersion as well as to invalidity of the laminar Darcy law even when flow velocities are extremely small. In I it was found that for moderate Reynolds numbers the doubly diffusive convection can be approximated by the singly diffusive convec-tion, for in such cases mechanical dispersion is larger than molecular solute and heat diffusion. The objective of this study is to analyze transport phenomena in the cavernous aquifer subjected to singly diffusive convection. BASIC EQUATIONS The analysis is related to a flow field similar to the Boulder Zone (the deep regions) of the Floridan Aquifer. We assume that density gradients are induced by a single component only, referred to as temper ature. In the cavernous strata, turbulent effects, as well as mechanical heat dispersion, are induced by even extremely slow fluid motions. In I the basic equations related to the calculation and approximated by the Boussinesq approach were presented-as follows: 50

PAGE 57

where au. 1 0 ax. -, + pgni + (1 + b)ui = 0 1 aT lL.. = _a_ (E. lL..) + u,. a a. aX. J a J Ui = velocity vector p = pressure p = fluid density K = permeability = porosity Jl = vi scos ity T = temperature a = coefficient of thermal expansion po,TO = density and temperature of reference Y = a coefficient defined by pC + p s C s(l (1) (2) (3) (4) (5) Through a brief literature survey presented in I it was suggested that the friction function, b, can be approximated by: b = 0.014 Re (6) In an isotropic medium the dispersion tensor, E;j' can be expressed as a sum of the isotropic molecular diffusivity and the second order 51

PAGE 58

symmetric mechanical dispersion tensor, as in E = KO + lJ lJ lJ where Here subscriptst and refer to transversal and longitudinal components respectively. (7) (8) An assumption of singly diffusive convection is justified for a mUlti-component system for moderate Reynolds numbers. A model suitable for the description of the mechanical dispersion tensor in such cases was suggested by Saffman [1959J. According to this model Ei -= \I E* --\I s + 2 (Re) 2(s + 1) (s + 3) $ (s + 1)2 Re 2(1 -s)(s + 2)(s + 3) ($) (9) (10) where s is a power coefficient describing the dependence between the velocity and the pressure drop (s varies between unity, for laminar flow, and half, for turbulent flow). The flow field model considered in this study while unperturbed conditions prevail, is described in Figure This is a saturated porous layer of infinite horizontal extent bounded by two impermeable planes located at bottom and top of the aquifer. Temperatures on bottom and top of the porous layer are To and respectively. Through the 52

PAGE 59

A '7 ---7/ // Uo TEMPERATURE PROFILE IMPERMEABLE Figure 1. Schematical description of unperturbed conditions 53

PAGE 60

porous layer the fluid flows uniformly in the longitudinal, x, direction. y and z are transversal and vertical coordinates respectively. The flow field variables can be nondimensionalized as follows: TI= (T Ta)/ T_ UI = Ud/E t E!. = E .jEt 1 J 1 J (p -p gz)K xu PI 0 + 0 (11) + b) where ti is a unit vector in the x direction. Substituting the dimensionless variables of (11) in (1)-(3) and omitting the primes we obtain -RTn. + (1 + s)u. = 0 aXi 1 1 3.T + _a_ (E rat ui ax. ax. ij ax. 1 1 J Here R is the Rayleigh number defined by The power s in (9) and (10) according to I can be calculated through (12) (13) .( 14) s = tn(Re) ( ) tn(1 + b) + tn(Re) 15 However, (15) can be utilized only when Re > 2.77 if b = 0.014 Re. 54

PAGE 61

For unperturbed conditions (12) and (13) yield u. = 0 1 8 = 0 T = -z 2 P = P -Rz /2 o LINEAR STABILITY ANALYSIS E.;(i;' j} = 0 lJ (16 ) (I7) (18) (19) (20) Stability criteria of the flow field are determined through the linear stability analysis. The flow field is subjected to small disturbances in the velocity (u, v, w), temperature (6), dispersion tensor (E .. ), friction function (8), friction power coefficient (s) and lJ pressure. Disturbances are very small and through the linear stability analysis second order terms are negligible. For the linear stability analysis the module of the velocity vector in the perturbed flow field is approximated by (21) Substituting (21) in (15), (9) and (10) we obtain expressions for the principal components of the dispersion tensor. By applying (8), all components of the dispersion tensor in the perturbed flow field are obtained. 55

PAGE 62

It is convenient to express the velocity components by utilizing a scalar function Q as follows (22) where (23) Introducing flow field perturbations in (12) and (13), neglecting second order terms and eliminating the pressure perturbation we obtain where 02 02 02 02 v =-2+-2+-2 ax ay az E* t c2 = Ud o (24) (25 ) (26) Assuming vanishing values of 6 and Q on top and bottom of the aquifer, these disturbances can be expanded by the following normal modes (6,Q) = (61,Q1) sin (nz)exp [i(a x + a.y) + (0 + iw)tJ (27) x yr where 61 and Ql' are constants, ax and ay are the wave number components. For point of stability (or = 0) substitution of (27) in (24) and (25) yields the following secular equation 56

PAGE 63

2 + 2 R I a 1T 2 ( 2 2 2 + 1T2 + a lax cIa C21T ) xa iYWI= a (28) where x = + 1 (29) (a/ay)2 + 1 The minimal value of R satisfying (28) is the critical Rayleigh number. Therefore, (28) yields the following criteria of point of insta-bil ity x = 1 a = 1T o a = a a = a x y w = a R = 4i o (30) Hence convection cells are two dimensional rolls whose axes are parallel to the unperturbed velocity vector. Superposition of the unperturbed velocity and the convection velocity leads to a helical flow field. Convection currents conducted in two dimensional rolls were also obtained for the ordinary Benard problem [Malkus and Veronis, 1958; Schluter et a1., 1965J as well as for free convection in a porous layer while mechanical dispersion effects are negligible [Straus, 1974J. However, in those cases only nonlinear stability analysis associated with stability analysis of the steady convection motion yield such a result,. whereas, in our study the linear stability analysis indicated that phenomenon. In our case, calculations concerning the anisotropy of the dispersion led to the conclusion that convection cells should be two dimensional rolls. 57

PAGE 64

Convection sets out in planes where the effective Rayleigh number attains maximal values, namely. where the coefficients of hydrodynamic dispersion attain minimal values. Inertial effects associated with the invalidity of the laminar Darcy law introduce the friction function, D. in the expression for R but do not affect the linear stability analysis and predictions. FINITE AMPLITUDE DISTURBANCES AND NONLINEAR STABILITY ANALYSIS The effect of the convection motion on transport processes through the aquifer can be predicted through the solution of the nonlinear equations of motion and heat transport related to supercritical conditions. The convection motion is two dimensional, therefore, the velocity components can be expressed by the stream function as follows v = 8z w 8y (31) Substituting the finite amplitude disturbances in (12). (13) and eliminating the pressure perturbation we obtain (32) (33) where Band H are the friction and the heat advection spectra, D and F are two parts of the heat dispersion spectrum, all of which are of the form = (6 1 1 (34) 58

PAGE 65

= a ( ae F(E .. ,e) = -,,-E .. -,,-) lJ aX. lJ aX. 1 J (35) (36) (37) Here sand Eij are finite amplitude disturbances in the friction functi on and the di spers i on tensor' r'espec ti vely. As long as the convection velocity is smaller than the unperturbed velocity the absolute value of the velocity in the flow field can be approximated by (38) where (39) (40) Applying (40), (6) and (15) we obtain s = [b/(1 + b)JA (41) = (1 + + b)J S (1 + + b)ReJ SA (42) Introducing (42) in (9), (10) and the dimensionless form of the dispersion tensor we obtain after minor approximations 59

PAGE 66

E2 EQ, E* 1 + t {(--1)(1 + A) { [1 5 (UOd)2 Et Et 1-2 5 E* 2s + 5 E* + E! (3 -5)J Q, I)} (1 A)} u u (5 + 2)(5 + 3) (E* -t t J (43) As long as the Rayleigh number is not very high, A is very small and can be approximated by the first term of (40), leading to the following approximations (44) (45) (46 ) where b E (31 = ( t )2 2(1 + b) ua 0 2 sl = Et 5 2(u od)2 A EtEt 1 1 1 E:t 0,1 = .. + r; (s + 2 -s + 1 s + 3) Et 2(u od)2 1 ( 47) E 2 ER-A t 0,2 = 2 (-1) (uod) Et 60

PAGE 67

If (32) and (33) are subject to the following homogeneous boundary conditions = 0 at z = 0,1 (48) then the system (32) and (33) can be solved by means of a set of truncated eigenfunctions expressing the finite amplitude disturbances. Double Fourier series expansions may conveniently be applied for such purposes as follows 00 A ,1, -I: 'l' 'f' -p,q=l p,q sin(pay) e 0p ,q cos(pay) q=l (49) (50) The calculation can be simplified by using the complex variable presentation of sin(pay) and cos(pay) leading to = -i 'l'p.,q e ipay p=_oo q=l 00 e = I: p=_oo 8 e ipay p,q provided that 'l' --p,q 1 0/ -P.q = 2' p,q 8 =8 =18 p,q -p,q 2 p,q 61 (51) (52) (53) (54)

PAGE 68

The solution of the system (32) and (33) through the Fourier series expansion means the determination of the Fourier series coefficients which can be done by power series expansion. Such a method was first used by Kuo [1961J for the analysis of the ordinary Benard convection. Palm et al. RUbiO[1975 j through approach for the analysis of free convection with no dispersion and turbulence in porous media. Rubin and Christensen [1975J suggested some guidelines for the utilization of such an approach when analyzing instabilities induced by salinity gradients in a saturated porous layer. The advantage of this method lies in its simplicity. According to (30) the inception of convection motion is of the marginal instability of exchange, namely, steady convection follows point of instability. Therefore, the first term in (33) vanishes when steady convection is attained. In such conditions the toefficients 'p,q and 8 p ,q of the Fourier series expansion can be expressed through a power series expansion as follows: N (, ,8 )= p,q p,q n=l (,(n) 8(n)) n p,q' p,q n where n is a small parameter defined by (55) (56) The Rayleigh number is also expanded by a finite power series as follows: S 2' R = R + R n J o os J= (57) 62

PAGE 69

where S = N/2 (58) Other flow field perturbations like v, w, a, andcEij as well as perturbation spectra H, B, 0, and F can be also expanded in double Fourier (n) (n) (n) and power series. Expressions for the coefflclents Hp,q' Bp,q' Dp,q and F(n) are presented in the Appendix. p,q If the analysis is conducted for N = 1, it reduces to the linear stability analysis, yielding critical values of wave number and Rayleigh number as given in (30). There is a pos i ti ve rel a ti onshi p between increases in Rayl ei gh numbers and wave numbers. However, taking the assumption that under supercritical conditions the wave number remai ns constant does not significantly affect predictions of transport phenomena for quite a wide range of Rayleigh numbers [Straus, 1974]. Such an assumption is not required by the method used here but considerably simplifies the analysis. Substituting the series expansions in (32) and (33) we obtain R e(n) + R e(n-2i) + + B(n) = 0 O P1f -p.,q os II m p q Tp q P q i =1 (59) (60) The functions and e are generated by superpositions of trigonometric functions. For n = 1 the only nonvanishing coefficients Therefore, only coefficients and with even values of lpl+q nonvanishing values. Moreover the values of the subscripts Ipi and q are smaller or equal to n. Coefficients with other subscripts vanish. According to (60). 63

PAGE 70

(61) For PFO, (59) and (60) yield 222 2 R [(p +q ) -4p ] + E p,q Ro i=1 (62) For p=q=l e(n) + E e(n-2i) + )_2(H(n+2) + D(n+2) + F(n+2/R (63) 1,1 i=ll,l 1,1 os 1,1 1,1 1,1 os Through (61), (62) and (63) all the coefficients and e(n) are p,q p,q detennined. According to (63) the coefficients and e(n+2) 1,1' 2,20,2 must pe determined simultaneously. However, by simple arrangements can be expressed explicitly. Such a procedure avoids any trial process. Calculations of the mean hori zontal temperature and the Nussel t number followed the determination of the series coefficients. These para meters are given by N T = -z (64) N () n-l (i) (n i) Nu=lqen +al''Ii"'!!l"L 'n=l q=l O,q 1,-p,q -p,S Calculation of Nu determines the convergence of the method and the termination of the series expansion. Through the calculation, presented in the next section, we took N = 6, 8, 10, 12, 14 and 16 according to the variation in the Nusselt number. If Nu varied by less than 2% as N was increased from N to N + 2 then the expansion was terminated. 64

PAGE 71

RESULTS AND DISCUSSION According to (44)-(47) convection effects are not determined only by the Rayleigh number but also by the Reynolds (Re) and Prandtl (Pr) numbers, as well as by the ratio between the porous layer thickness and the characteris-ti cfJo-re -( d/-dp whole a-naly-s i cab-leonly when mechanical dispersion is at least comparable with the molecular diffusion. Th;s relationship is determined by the magnitude of Re and Pro It seems that if Pr 1, which is reasonable for practical purposes of hydrology [Somerton, 1958J, then dispersion effects are significant even for Re 3. However, if the Prandtl number is smaller than unity (if the solid fraction is a good conductor) then dispersion effects become significant at higher Reynolds numbers. The effect of Re, Pr and d/dp on the convection phenomenon is introduced in the analysis through 6 1 and (and = + which determine effects of turbulence and dispersion induced by the convection motion. Figure 2 demonstrates changes in these coefficients due to variations in Re, Pr and d/dp when the porosity is 0.4. The effect of Pr vanishes for high Reynolds numbers as in such cases mechanical dispersion affects transport processes more than the molecular diffusion. When the Reynolds number increases, according to model (9) and (10), the mechanical dispersion tensor becomes more and more isotropic. This phenomenon leads to a reduction in the value of Turbulent effects become more and more significant when the Reynolds number increases, leading to a positive relationship between increases in the Reynolds number and the coefficient 65

PAGE 73

1\ a2 1\ al (31 01 "-J 10-1 3 10-2 3 10-3 3 10-4 3 10-5 3 10-6 3 -----------------===== Pr -'!!!!!I!!II---__ --;:;.;z;ft -0.5 -iil O did =103 -_ p --------.------=--________ -+ ________ __________ __ 3 10 30 100 Re Figure 2. Description of 8, ,a, and &2 vs. Re for various values of Pr and d/d p(cp=O.4).

PAGE 74

R/Ro:5..1 ---dId =20 p ----dId:::: 10 p R/Ro= 10 R/Ro=4 1.0 0.8 0.6 0.4 0.2 -1.0 -0.8 -0.6 -0.4 -0.2 0.0 T Z Figure 3. Mean horizontal temperature profiles for various values of R/R and did Pr=l. o p Re=3) 68

PAGE 75

At high Rayleigh numbers a thick region in the center of the porous layer should achieve a nearly isothermal state in the mean. However, as presented in Figure 3 such conditions lead to a positive temperature gradient or reversal of temperatures. Such a phenomenon was also identified in the ordinary Benard problem [Kuo, 1961; Veronis, 1966J. Veronis [1966J tried to explain the origin of such a strange phenomenon. However, better choice of wave numbers diminishes this effect. Figure 3 also demonstrates the creation of boundary layers on top and bottom of the aquifer leading to invalidity of the continuum approach even for moderate Rayleigh numbers if d/d p is not very large. Figure 4 presents variations of Nusselt number with Rayleigh number for various values of d/d p According to this presentation, the net effect of the mechanical dispersion induced by the convection motion leads to a reduction in the heat transport through the aquifer. Disper sion and turbulence induced by the convection motion act as a stabilizing mechanism in the flow field (ironical interpretation). However, this mechanism is more complicated than just a stabilization. In the calculation it was found that the phenomenon was associated with increased values of the higher modes of the series expansion, leading to a reduction in the convergence of the series expansions. Through the calculation, the maximal values of A and V were continu ously checked in order to examine the validity of approximations (44)(47) and to follow changes in the Reynolds and the Peclet numbers. It 69

PAGE 76

" o Figure 4. Description of Nu vs. R/Ro for various values of d/dp(=Oo4, Pr=l, Re=3). Nu 7 6 5 4 3 2 o 2 .. ."..:d/dp= 00 ---,d/dp= 20 -----d/d = 10 I P 4 6 8 10 R/RO

PAGE 77

was found that for the range R/Ro :: 10, >max was much sma 11 er than unity whi justi es utili on of (44)-(47) and yields only minor changes in the lds number. CONCLUS IONS The anisotropic character of the dispersion tensor leads to convection motions conducted in two dimensional rolls whose axes are paranel the ullpe-durbe-dve lodty vector. By choos i ng a new defi ni ti on for the Rayleigh number, rbulence and dispersion effects can be introduced in the linear stability analysis no substantial complications in the calculations and results. Finite amplitude analysis for homogeneous boundary conditions can be by Fourier series and power series expansions. The significance of mechanical dispersion and tlwbulence induced by the convection motion depends on PI", Re and mainly on d/dp An increase in Re diminishes the effect of PI" and reduces effects of dispersion and turbulence associated with the convection motion. An increase in d/d p reduces significantly effects of dispersion and turbulence due to the convection motion. For d/dp 102 these effects practically vanish. However. for smaller values of d/d p these effects lead to a reduction in Nuo The singly diffusive analysis presented in this article can be applied when density gradients are induced by a single component or IlIIhen boundary condit; ons and effecti ve di spers; on coeffi ci ents are identical for all components in a system. The latter canditionis approximately satisfied when mechanical dispersion effects are larger than any molecular diffusion in the aquifer. 71

PAGE 78

APPENDIX where sl=q-t Expressions for coefficients of heat advection, friction and heat dispersion spectra e!p-k!,S2e !p-k! ( ) 4 n-1 i 1 00 ( .) (. ) B n -1T L L J'-J p,q 8 1 "4 j=l k ,m=-oo If'k,t If' m,h t,h=l {th[4k(p-k-m)+(p-k-m)2+ 3 s (v=1,2,3,4;S,6,7)] v +[3k2hs -t2m(p-k-rn)][,(nki) (v=3,4,5;1,2,6,7)]} v p--m,sv ( ) 4 n-1 i -1 ( ) (.' ) F n -1T L L r If' J If' l-J !p! ,q16 i=2 j=l k,m--oo k,t m,h t,h=l {th[(&l + ,sv(V=S,6,7;1,2,3,4)] +m[2&1 k 2 -&2t2 ) [e f I,s/ v=3,4,S; 1 ,2,6,7)] (68) (69) (70) ,sv(v=1,3,S,7;2,4,6)]} (71) In (70) and (7l) the subscript Sv may obtain seven different values as follows: sl=q-t-h s6=h:..t-q s2=q+t+h s7=t-h-q '72 (72)

PAGE 79

NOTATION a wave number components of wave number critical wave number b friction function B friction spectrum Bp-;,q. hmts in-sari Ofl-SFOr B c1,c2 C C s d d p o D ,O(n) p,q P.q E'!': E lJ lJ q,Et F F F(n) p,q p,q 9 H H ,H (n) p,q P.q K JI.,. 1 n 1 coefficients defined in (26) specific heat of fluid specH; c heat of solid porous layer thickness characteristic pore size part of heat dispersion spectrum coefficients in series expansions for D heat dispersion tensors (mechanical and hydrodynamical respect i ve 1 y ) longitudinal heat dispersion coefficients transversal heat dispersion coefficients part of heat dispersion spectrum coefficients in series expansions for F gravity acceleration heat advection spectrum coefficients in series expansions for H permeabil i ty unit vector in the longitudinal direction unit vector in the vertical direction 73

PAGE 80

N total number of terms in the series expansion Nu Nusselt number p pressure Po pressure at the coordinates origin Pr Prand-tl number ( \lId R Rayleigh number defined in (14) Ro critical Rayleigh number Ros parameter defined in (58) Re Reynolds number s power coefficient s = N/2 t time T temperature To temperature at z = a u longitudinal velocity perturbation Ui velocity vector Uo unperturbed velocity U module of velocity vector v lateral velocity perturbation V module of velocity perturbation w vertical velocity perturbation x coordinate 1 x.y,z coordinate system a coefficient of volumetric thermal expansion a 1,a2,a3 coefficients defined in (47), (a3 = a1 + a2)' 8 perturbation in the friction function 74

PAGE 81

6 1 coefficient defined in (47) y parameter defined in (5) 0.. Kronecker1s delta lJ 6T difference in temperature between bottom and top of the porous layer Eij dispersion tensor perturbation n small parameter defined in (56) 8 temperature perturbation 81 constant defined in (27) A (n) ff"' t" f 8p,q8p,q8p,q coe lClen s ln serles expanslons or e K thermal diffusivity of saturated porous medium parameter defined in (40) v viscosity v kinematic viscosity perturbation in s coefficient defined in (47) P fluid density Po density at z = 0 P s solid density Gr parameter expressing amplification of small disturbances porosity x parameter defined in (29) stream funttion (n) ff" t" f coe lClen s ln serles expanslons or w parameter expressing oscillations 75

PAGE 82

Q scalar function defined in (22) Q1 constant defined in (27) 76

PAGE 83

REFERENCES Kuo, H. L., IISo 1 u tion of the Non-linear Equations of Cellular Convection and Heat Transportll, J. Fluid Mech., 1.Q, 1961, pp. 611-634. Malkus, W. V. R., and Veronis, G., IIFinite Amplitude Cellular Convection", J. Fluid Mech., 1, 1958, pp. 225-260. Palm E., Weber, J. E., and Kevernold, 0., liOn Steady Convection in a Porous Medium", J. Fluid Mech., 54(1), 1972, pp. 153-161. -lH-r.-=,----,,-lIftQrrn -1ti"tbfl=lertAnrr.a:rtl ys ; 5 of Ce 11 u 1 a r Co nvec t lOn n Porou s Med i a II Int. J. Heat & Mass Trans., ]&, 1975, pp. 1483-1486. Rubin, H., "Onset of Thermohaline Convection in a Cavernous Aquifer", to be published in Water Resourc. Res., 1976. Rubin, H. and Christensen, B. A., IIConvection Currents Associated With Hydrodynamic Dispersion in a Porous Mediumll, 16th Congress of IAHR, Sao Paul, Brazil, 1975. Saffman, P. G., "A Theory of Dispersion in a Porous Medium", Fluid Mech., 6(3), 1959, pp. 321-349. SchlUter, A., Lortz, D., and Busse, F., liOn the Stability of Steady Finite Ampl itude Convectiontl, J. Fluid Mech., 23(1), 1965, pp. 129144. -Somerton, H. W., "Some Thermal Characteristics of Porous Rocks", J. Petrolr Tech., Note 2008, 1958, pp. 61-65. Straus, J. M., "Large Amplitude Convection in Porous Media", Fluid Mech., 64(1), 1974, pp. 51-63. Veronis, G., "Large Amplitude Benard Convection", J. Fluid Mech., 26(1),1966, pp. 49-68. 77

PAGE 84

INTRODUCTION CHAPTER 5 SEMINUMERICAL APPROACH FOR THE MATHEMATICAL MODELING OF SINGLY DISPERSIVE CONVECTION IN GROUNDWATERS Considerable interest ;s now focused in the State of ElDrida, as well as in other locations in the United States, for the possible utilization of deep saline aquifers for waste disposal. Vernon (16) delineates the properties of the deep zones of the Floridan aquifer which makes this stratum available for such application. Henry and Kohout (2) mention the fact that in a thick system. like the Floridan aquifer, the effect of geothermal activity should be considered, too. One of these authors has postulated in previous articles (3,4) that geothermal activity induces groundwater circulation in the Floridan aquifer. Singly diffusive convection is the convection motion induced by a single dissolved component (e.g. temperature or salinity) in a fluid layer. The hydrodynamics of this phenomenon in a saturated porous media has been studied while, in most instances, assuming that the fluid is initially at rest (7,8). However, groundwaters are generally subject to hydraulic gradients leading to slow, effectively horizontal flow of the subsurface water. This movement through a formation, similar to the Boulder Zone in the deep saline region of the Floridan aquifer (1,4), leads to intensive mechanical dispersion of heat and soluted in the aquifer, as well as inertial effects, as demonstrated by invalidity of the laminar Darcy law. In such a system molecular diffusion effects are usually less significant than mechanical dispersion. Convection phenomenon under such conditions may therefore, be called dispersive convection. A system is subject to singly dispersive convectionif one of the following criteria ;s satisfied: a) density gradients are 78

PAGE 85

in a single component; or b) all molecular diffusivities of the dissolved components are much smaller than the mechanical dispersion, and all of them have the same boundary conditions. The objective of this art"icle is to present a rather simple method by which flow conditions in such an aquifer can simulated. BASIC EQUATIONS The basic equations applied for the analysis are the equations of __ cQnti nuitY3_ mot "ion ,_i:ul(l di sperstolJ, sublect to the_Boussi -+ \7U = 0 0 \7p + pgh + k (l+b)u = 0 aT -;(=) Y3f + u'\7T \7-E\7T .. The coefficient y appearing in Equation 3 is defined by: (1) (2 ) (3) . ( 4 ) In a brief literature survey, presented in a previous study(lO), it was suggested that the friction function, b, appearing in Equation 2 can be approximated by: . . ( 5 ) It is assumed that the fluid density depends linearly on the dissolved component, which is the temperature, as follows: ...... ....... (6) 79

PAGE 86

Uo .. .. TEMPERATURE .. .. PROFILE Figure 1. Schematical description of the aquifer with no convection motion. 80

PAGE 87

The flow field model considered prior to the inception of the convection motion, as presented in Figure 1, is a cavernous aquifer consisting of a porous layer of infinite horizontal It is bounded by two impermeable planes on which the temperature is constant. Through the porous layer the fluid flows uniformly in the longitudinal x direction. The transversal and vertical coordinates are y, and z, respectively. A moving coordinate system is applied with the velocity uo/y in the x direction, and consider (7) as characteristic length, temperature, dispersion coefficient, velocity, time, and pressure, respectively. In such a manner the following dimensionless basic equations are obtained (in further expressions all variables are dimensionless): + vu = 0 . (8) vp RTh + = 0 (9) vT = V.(r'VT) (10) Here the Rayleigh R -agll TKd number, R, and the variable B, are defined as follows: B = b[(U/u o ) -l]/(l+b) (11) In an isotropic medium the dimensionless dispersion tensor can be expressed by: . (12) As long as there is no convection motion Equations 8-10 yield: + 0 B = 0 T = Rz2/2 u = -z p = p 0 Exx = E Eyy = Ezz = 1 Exy = Exz = Eyz = 0 (13) JI, . . 81

PAGE 88

THE FLOW FIELD STABILITY The stability of flow conditions presented by Equations 13 can be determined by analyzing the growth of small disturbances in the aquifer. These are disturbances in the velocity, u, temperature, e, dispersion tensor, E and pressure,P. Second order terms depending on these disturbances are negligible. It is convenient to express the velocity vector through the scalar vari-able Q asf611ows: where t :: rtx'i,7X'i,7 . . . . (14) By substituting the small disturbances in the basic equations, and applying the boundary conditions 2 e,\71Q = 0 at z = 0,1 ........................ (15) (where :: we obtain: R'i,7iQ = + . . . (16) 2 ++ 2 'i,7d = 'i,7\7 + 'i,7 ................. (17) By expand
PAGE 89

E STEADY STATE CONVECTION In the prey; ous sect; on it was plAoved that convect; on cell s two dimens'iona"J rolls. It is convenient to express the velocity by the stream function . . ( 20 ) where the square brackets symbolize the box product. Substituting Equation 20 in the basic equations, and eliminating the pressure perturbatiQn, we_obtain: R[rt, 1. VT] + = S) .. + v 2 e -crt, =: + + F(, e). where is the finite amplitude disturbance in the dispersion tensor. Variables B, H. 0 and F are nonlinear terms defined as follows: S) = e) = -[1. ve] = e) = -v(E:-ve). (21) (22) (23) (24) (25) (26) As long as convection velocity is smaller than the unperturbed velocity the absolute value of the flow field velocity can be approximated by: (27) where 2 -+-+ V ::: UU If 1 terms depending on high orders of can be neglected. In such conditions. the expressions for B and the dispersion tensor can be by: 8 =: 8 1 V 2 83 (28) (29) (30)

PAGE 90

where 2 6 1 = b/[2(1+b)uoJ .. a 1 = (1/2u o)(aEt/au o ) a = 3 .. . . (31) 2 a2 = (E9,-l)/u o + (1/2u o)[a(E9,-Et)/au oJ . (32) The system of the differential equations 21 and 22_ ;s s_ubject to the foll ow; ng boundary cond it i on s : e = 0 at z = 0,1 ........................ (33) Such boundary conditions are simple according to the definitions presented by Orszag (9). In such cases, accurate simulation of incompressible flows can be obtained by spectral methods. Assuming that and e are periodic in the horizontal direction, these variables can be presented by sets of truncated eigenfunctions as suggested by Veronis in similar studies (17,18). = ; sin(pay)sin(qwz) .................. (34) p.q=l p,q e = e cos(pay)sin(qwz). .................. (35) p=o p,q q=l The calculation can be simplified by using the complex variable presentation of sin(pay) and cos (pay) leading to = -i E_ exp(ipay) ................ (36) D--oo 00 e = E p=-oo e [sin(qwz)] exp(ipay) ................. (37) p,q q=l provided that 84

PAGE 91

The convergence of the expression for the Nusselt number also determines the truncation of the Fourier series expansion (the value of N). NUMERICAL CALCULATIQNS Experiments concerning heat transfer characteristics of porous rocks (e.g. 5,6) showed that the phenomenon ofheatdjspersion_duetothe fluid movement in the stratum is very similar in nature to the characteristics of solute dispersion in porous media. We may adopt, then, models of mechanical dispersion, which are available in the literature, for the quantitative evaluation of the effect of the convection phenomenon on the intensity of transport processes in the aquifer. Saffman (11) suggested the following expressions for the coefficients of dispersion. et 1 + s + 2 \) Pr 2 ( s + 1 )(s + 3 ) ( :e ) . . . . (42 ) (s + 1)2 + 2(1-s)(s+2)(s+3) (:e). . . (43) where s is a power coefficient describing the dependence between the flow velocity and the pressure gradient. In a previous article (10) it was suggested that this coefficient can be calculated through the following expression s = In(Re)/{ln[(l+b)Re]} (44) Equation 41+ ;an h", "'d' ,":,; ;2.77. Through the calculation, it was assumea that PI" = 1, which seems to be a reasonable value for limestone and dolomite aquifers (13). The convergence of the numerical integration was very moderate. Several approaches were applied to speed convergence of the calculations. These were: (a) Each calculation was conducted in two steps, in the first step we took 86

PAGE 92

and. after ining resul for was introduced into the program; (b) Results lower values Nand/or R were used as initial quanti es hi parameters. As a criterion for steady calculations, to ta { I ( ., 04 < Ii 0 (45) In our calcula ion, it was that such a criterion is very conservative. a criterion of one hundred times less conservative, vari 1% were obtained in the selt number. Through the ng oscillations were detected in the value the 2 ations of the Nusselt number with Rayleigh number va 0(15 values These are the maximal values of filL! obtained when va value of the wave number (12). The wave number increased with values of ions in the wave number induced minor changes in Nu for a constant ei to Fi 2, mechanical dispersion, due to the convection motion, leads to a ion in the intensity of transport processes through the aquifer. The calculation 1 icated that in small values of d/d p the magnitude of the dis in This s es ion cients increases, when convection occurs, leading to a reduction the tempet'ature gradients on top and bottom of the on in temperature gradients affects transport processes more increase in ue of the dispersion coefficients. It was also t "IO'IrJer va 1 ues d/d p give rise to the higher modes of the Fourier the convergence of the calculation, For very ions a 87

PAGE 93

co co u 7 6 5 4 3 2 o 2 4 6 d/dp:: 00 d/dp:: 20 d/dp:: 10 8 10 R/RO Figure 2. Description of Nusselt number vs. R/Ro for various values of Pr=l; Re=3).

PAGE 94

small values of did the analysis fails to follow the physical phenomenon, p even at moderate Rayleigh numbers, as the thickness of the boundary layers developed on top and bottom of the aquifer is of the same order of magnitude as the characteristic pore size. In such cases the continuum approach applied through the analysis is not satisfied. Through the calculation, the maximal value of was continuously checked in order to examine the validity of the approximations presented in Equations 27-32. It was found that fOt' R,LR<"IO-, the maximal-value of .7-. was mtfchsmaller-thanunity. o --CONCLUSIONS Singly dispersive convection in aquifers can be analyzed by expanding the flow field disturbances in eigenfunctions. The anisotropic character of the mechanical dispersion determines the plane in convection motions are conducted. According to the analysis. convection cells are two dimensional rolls whose axes are parallel to the unperturbed velocity vector. Analysis of steady state convection can be conducted by transforming the equations of motion and continuity to a set of general first order differential equations that can be integrated through available subroutines. The effect of mechanical dispersion and turbulence, induced by the convection motion, depends mainly on the ratio between the porous layer thickness and the characteristic pore size; Prandtl and Reynolds numbers have less significant effects on the physical phenomenon. The analysis of the singly dispersive convection presented in this study can be appl ied when density grad"ients are induced by a single or when bOL!ndar'y condition_ and effective dispersion coefficients are identical for all components in a mul{icomponent system. 89

PAGE 95

APPENDIX I EXPRESSIONS FOR THE COEFFICIENTS OF THE SPECTRAL FUNCTIONS N-l N-I k I Hp,q = H_p,q = (wa/2) 'k,S{S(p-k}[Slp_kJ, W y (Y=1,2;3}J + W y (Y=2;1,3}J} .................... (46) .... ___ 2 __ _!! ... .... .-Dp,q =D."p:q= 222 2 + 2cxl(n /a }wy cx2(p-k) J['p_k W (Y=1,2;3}J} ............ (47) y where w1 = q-s Slp_kl ,wy(y=2;1 ,3} = Slp_kl 'W1 Slp_kl ,w2 Slp_kl 'W3 ............ (48) 2 2 N-1 N-Ikl N-Iml B = -B = (w a /4) L L L p,q -p,q k,m=l-N s=l h=l k,s m,h 2 2 2 2J ( -. ) J {sh[4k(p-k-m)+(p-k-m) + 3(w /a }wy y-1,2,3,4,S,6,7 + [3(a2/n2}k2hWy-s2m(p-k-m}J['p_k_m,Wy(y=3,4,5;1,2,6,7}J} ...... (49) where W 1 = q-s-h ,w5 = s+h-q W 2 = q+s+h W6 = h-s-q W3 = q-s+h W 7 = s-h-q 90 W 4 = q+s-h . (51)

PAGE 96

NOTATION a = wave number; a o = critical wave number; ax' a = y wave number components; b = friction function defined in Equation 5; B = friction spectrum defined in Equation 23 ; B -coefficients in the Fourier series q.q C s C w = specific heats of soil and water, respectively; d = porous layer thickness; d = characteristic pore size; p D = part of heat dispersion spectrum defined in Equation 25; D p,q = coefficients in the Fourier series expanded for D; eQ,' et = longitudinal and lateral dispersion coefficients. respectively; F H = dimensionless dispersion coefficients (dispersion coefficients divided by e t existing prior to convection conditions); E = dispersion tensor; F = part of heat dispersion spectrum defined in Equation 26; = coefficients in the Fourier series expanded for F; P.q g = gravity acceleration; H = heat advection spectrum defined in Equation 24; = coefficients in the Fourier series expanded for H; P.q I = unit matrix; K = permeability; + Q, = unit vector in the longitudinal direction; -+ n = unit vector in the vertical direction; N = truncation parameter; Nu = Nusselt number; p = pressure; 91

PAGE 97

Po = pressure at z = 0; Pr = Prandtl number (=V/K); R = Rayleigh number defined in Equation 11, Ro = critical Rayleigh number; Re = Reynolds number (=
PAGE 98

\! kinematic viscosity; P = fluid dens ity; Po = fluid dens ity at z = 0; P s = solid density 0", 0"2 = parameters expressing growth and oscillation of disturbances;


PAGE 99

REFERENCES 10 Burke, R. G., "Sun Whips Florida's Boulder Zonell, Oil and Gas Journal, VoL 65, No.7, 1967, pp. 126-127 2. Henry, H. R., and Kohout, F. A., "Circulation Pattern of Saline Ground vlater Affected by Geothermal Heating as Related to Waste DisposalI', Proceed;n s of the 1st S osium on Under round Waste Mana ement and Environmental Implications, +. -D. Cook ed American Association of Petroleum Geologists, Tulsa, Oklahoma, 1972, pp. 202-221. 3. Kohout, F. A., IIA Hypothesis Concerning Cyclic Flow of Salt Water Related to Geothermal Heating in the Floridan Aquifer, New York Academy of Science Transactions, Vol. 28, No.2, 1965, pp. 249-271. 4. Kohout, F. A., "Groundwater Flow and the Geo.thermal Regime of the Floridan Plateauli, Symposium Geological History of the Gulf of Mexico Caribbean Antil1en Basin, Gulf Coast Association of Geological Societies Transactions, Vol. 17, 1967, pp. 339-354. 5. Kunii, D and Smith, J. M., "Heat Transfer Characteristics of Porous Rockll, American Institute of Chemical Engineering Journal, Vol. 6, No.1. 1960, pp. 71-78. 6. Kunii, D., and Smith, J. M., IiHeat Transfer Characteristics of Porous Rocks: II. Thermal Conductivities of Unconsolidated Particles with Flowing Fluids", American Institute of Chemical Engineering Journal, Vol. 7, No.1, 1961, pp. 29-34. 7. Lapwood, E. R., "Convection of a Fluid in Porous Media", Proceedings of Cambridge Philosophical Society. Vol. 44, 1948, pp. 508-521. 8. Nield, D. A., "Onset of Thermohal ine Convection in a Porous Mediumll, Water Resources Research, Vol. 4, No.3, 1968, pp. 553-560. 9. Orszag, S. A., "Numerical Simulation of Incompressible Flows Within Simple Bounda ri es; Accuracy", Journa 1 of Fl u i d Meehan; cs, VA 1. 49, Pt. 1, 1971, pp.75-112. 10. Rubin, H., "Onset of Thermohaline Convection in a Cavernous Aquifer", Water Resources "'search, Vc 1, P, No.2, 1976, pp. 141147. 11. ':::"";;;", P. G., IIA Theory of u;spersion in a Porous Medium", Journal of Fluid Mechanics, Vol. 6, Pt. 3, 1959, pp. 321-349. 12. Schluter, A., Lortz, D., and Busse, F., liOn the Stability of Steady Finite Amplitude Convection", Journal of Fluid Mechanics. Vol. 23, Pt. 1, 1965, pp. 129-144. 13. Somerton, H. W., "Some Thermal Characteristics of Porous Rocks", Journal of Petroleum Technology, Vol. 10, Note 2008, 1958, pp. 61-65. 94

PAGE 100

14. Straus9 J. M., "Finite Amplitude Doubly Diffusive Convection", Journal of Fluid Mechanics, Vol. 56, Pt. 2, 1972, pp. 353-374. 15. Straus, J. M., "large Amplitude Convection in Porous Mediall, Journal of Fluid Mechanics. Vol. 64, Pt. 1,1974, pp. 51-63. 16. Vernon, R. 0., "The Beneficial Uses of Zones of High Transmissivities in the Florida Subsurface for Water Storage and Waste Disposalu Information Circular No. 70, Florida Bureau of Geology, 1970. 17. Veronis, G IIlarge Amplitude Benard Convection", Journal of Fluid Mechanics, Vol. 26, Pt. 1, 1966, pp. 49-68. 18. "Effect-Df Solute {)n Thermal-Convect-ion!',-Journal of Fluid Mechanics, Vol. 34, Pt. 2, 1968, pp. 315-336. 95

PAGE 101

INTRODUCTION CHAPTER 6 NUMERICAL SIMULATION OF SINGLY DISPERSIVE CONVECTION IN GROUNDWATERS Considerable interest is now focused in the State of Florida as well -as other locatiofis lnthe U: S. for the PQssible utilization of the deep and saline aquifers for waste disposal. Vernon [1970] counted the properties of the deep zones of the Floridan Aquifer which he believed makes this stratum capable of such applications. Henry and Kohout [1972J mentioned the fact that in a thick system like the Floridan Aquifer the effect of geothermal activity should be considered too. That statement was based on previous articles [Kohout, 1965; 1967J postulating that geothermal activity induces groundwater circulation in the Floridan Aquifer. Only sophisticated numerical approach may be used for the simulation of all possible phenomena that may occur in an aquifer like the Floridan Aquifer subject to geothermal activity. However, the applicability of such an approach should be determined according to ranges of variations of the parameters determining the physical phenomena. The purpose of this study is to determine such ranges for the possible of finite difference numerical simulation. ,:Ol\lc_ The basic equations utilized for the analysis were presented in a previous study 1976bJ. These are the followii!.;) equations of continuity, motion and dispersion subject to the Boussinesq approximation 96

PAGE 102

where 'lou = a vp + (l+b)u = a aT -+ = y 3t + uvT = v(E-VT) -+ u, velocity vector p, pressure p, fluid density K, permeabi 1 i ty cp, porosity ]1, viscosity -+ unit vert; ca 1 in n, vector T, temperature E, heat dispersion tensor y, coefficient defined by y = pCwct> + Ps Cs(l-ct pCwct> the downward direction (1) (2) (3) (4) In a brief literature survey presented in a previous study [Rubin, 1976aJ it was suggested that the friction function, b, would be approximated by b = 0.014 Re (5) It is assumed that the fluid density depends linearly on the temperature as follows p = p [l-a{T-T )J o 0 (6) 97

PAGE 103

where a 9 coefficient of thermal expansion po.T o density and temperature of reference. The flow field model considered in this study prior tc the inception of convection motion is described in Fig. 1. This is a saturated porous layer of infinite horizontal extent bounded by two impermeable planes on which the temperature is constant. Through the porous layer the fluid flows uniformly in the longitudinal. x. direction. The transversal and vertical coordinates are y and z, respectively. By considering as units of length. temperature, dispersion coefficient and velocity. respectively. we nondimensionalize the flow field variables as follows )t!= ()t-uot1/y)/d TI ::: (T-To)/lIT it! = lfd/et ui 0 = uod/e t U' = Ud/et (7) t' = tet/d 2 Ei ::: E/e t the vctriables of (7) in (1)-(3) and omitting primes we obtain vp + + (l+S)u = 0 aT -+ = + uvT = v(EVT) at (8) (9) 98

PAGE 104

UNIFORM VELOCITY PROFILE, Uo TEMPERATURE PROFILE T + liT o x IMPERMEABLE Figure 1. Schematical description of the aquifer with no convection motion. 99

PAGE 105

Here the Rayleigh number R and the coefficient S are defined respectively as fall ows s = b[(U/uo)-IJ/(I+b) (10) -In an isotropic medium the dispersion tensor can be expressed.by where = = ++ 2 E = I + (E 1) uu/U .R, I, unit matrix As long as there ;s no convection motion, (8) and (9) yield -7-U = 0 s = 0 T = z E = E = 1 yy zz THE FLOW FIELD STABILITY p = p -Rz2 /2 o EXY = Exz = Eyz = 0 ( 11) (12) (13 ) Criteria of instability of the flow field can be determined easily by analyzing the growth of small perturbations in the flow field. These are di sturbances in the velocity vector (u). temperature (8), di spers i on tensor and pressure. Second order terms depending on these disturbances are negligible. Such conditions yield that the variation W. in the absolute value of the velocity vector is (14) Therefore can be obtained through the following expression 100

PAGE 106

(15) It is convenient to express the velocity components by utilizing a scalar function as follows -+ -+ u = where -+ -+ o = n x v x v (16) --i tu t 1 n g the sma 11 pe rtu rba t fOns-Tn-l1FT3T-ci nd -;1 imina tin g the pressure disturbance we obtain (17) where (19) The disturbance e and the vertical component of the velocity perturbation vanish on top and bottom of the aquifer. Hence (17) and (18) are subject to the following boundary conditions at z = 0,1 (20) Applying (20) and substituting 8 by through (17) in (18) we obtain (21) 101

PAGE 107

We expand Q in the following normal mode (22) Substituting (22) in (21) for point of instability (01 = 0) and separating real and imaginary parts we obtain R = (a2 + + (23) Ra 2 2 x 2 [(E-l)a lu o + aEt/auoJ y(a ) x (24) where (25) The critical Rayleigh number is the minimal value of R satisfying (23). Therefore instability is demonstrated by x = 1 ax = 0 ay = a = 0 R 2 a o = = 0 (26) Hence the principle of exchange of stabilities holds. Convection cells are two dimensional rolls whose axes are parallel to the unperturbed velocity vector. Superposition of the unperturbed and the convection motions yields a helical flow field. The anisotropy of the hydrodynamic dispersion led to the characteristics of the point of instability presented in (26). NUMERICAL CALCULATION OF THE STEADY STATE CONVECTION In the previous section we proved that the convection cells are two dimensional rolls. Hence, the convection velocity components can be 102

PAGE 108

expressed by the stream function as follows -+ -+ U ::: -JlX\l\)! -+ -+ \)! = -[!/', n, 'i7r1J ( 27) where the square brackets symbolize the box product. Substituting (27) in (1)-(3) and eliminating the pressure we obtain (28) (29) If the model correlating the dispersivities with the velocity vector is known,then (28) and (29) can be solved while being subject to the specific boundary conditions of the aquifer. In order to save computer time quantities, it is convenient to expand the absolute value of the velocity vector in a power series as follows (30) where V is the absolute value of the convection velocity (V = uu). The variable 6 appearing in (28), according to (5), can be expressed through 6 (31) where 6 1 = b (32) 2(l+b)U/ The dispersion tensor can be expanded in power and Taylor's series as foll ows 103

PAGE 109

where (34) (35) From all the terms appearing in the right hand side of (33) only the first two are important for the calculation of transport phenomena and flow conditions. It should be mentioned that approximations (30) (35) are valid only for convection velocities slower than the flow velocity induced by the horizontal hydraulic gradient. According to (22) the boundary conditions determining the convection motions are aW W, ay' T = 0 at z = 0 aW W'ay O. T = 1 at z = 1 (36 ) -Zl = 0 1jJ'az' ay at y = O,L where L is the width of the half convection cell (L = n/a). This parameter is determined through the maximization of the Nusselt number. The system (28) and (29) was solved in this study by a combination of the SOR (Successive Over Relaxation) method for the solution of (28) and the ADI (Alternating Direction Implicit) method for the solution of (29). The FD (Finite Difference) approximation for (28) by applying central difference approximations yields the following SOR scheme for the m+l iteration 104

PAGE 110

,/,(m+1) = ,/,(m) + {'i,(m+l) (m+l) + 1jJ(m) + 'V 'V w'Vi-l,j+1)Ji,j-l i+l,j 1,J+1 1,J ; ,j i ,j -R(h/2)(T. J'+l T .. 1) + + 1" 1.J1,J 1,J 1,J-( 37) where w is the relaxation coefficient. As a first approximation for its optimal value we took its value for the Poisson's equation given by [Smith, 1969J W = 2 [4 + + (38) Usually Sl is very small. Therefore, (38) is applicable, leading to a good convergence of (37). The solution of (29) according to the ADI method was conducted in two successive schemes; the first yields the approximate temperature distribution T* from the given temperatures Tn (at time t n ) implicitly in the y direction through the following scheme 105

PAGE 111

. 1 1 1,J + 1.J+ 1.J1 .J 1.J k/2 V;,j 2h + W;,j 2h 1 1,J+ 1,J2h n n aE aE T +1 -T. 1 + (-E .. 1 ,J 1.J ay az 1.J 2h n n n .+T-I: 1 T'+ l .-2T .. +T. 1 + (E ) .. 1,]+ 1,J l,J-+ (E) 1 ,J l,J 1-,J yy 1,J h 2 zz ;,j h 2 Tn Tn Tn +T" .:r J 1-. 1 1-. 1 1 1 1 ( ) 1+.l,J+ 1,J+ 1+ ,J-1-,J-+ 2 E ., 2 yz 1,J 4h (39) The next step is the determination of the solution for T"+l by applying an implicit scheme in the z direction as follows n+1 n+1 n+1 T. H T"* 1T"* 1 T + 1 -T. 1 1,J 1,J l,J+ 1,]-1 ,J 1.J k/2 + Vi,j 2h + W;,j 2h aE aE 1 '" .-..J::L). 1,J 1,J-( ay + az ;,j 2h aE aEzz + _._) ay az i ,j n+1 n+1 T'+l .-T. 1 1 ,J 1-,J 2h .+T-I: 1 (E) l,J 1,J 1,J-(E) + 2 + .. yy 1,J h zz 1,J 1 1 1 () 1,J 1-,J+ 1 ,J1,J-+ 2 E ., 2 yz 1.J 4h 106 (40)

PAGE 112

Schemes similar to (39) and (40) were applied in previous studies concerning the ordinary Benard convection as well as singly diffusive convection in porous medium [Az;z and Hellum, 1967; Holst and Aziz, 1972; Cabell; and De Vahl Davis, 1971J. The value of the derivatives at the convection cells boundaries were calculated according to approximations with round off error of O(h2 ) as follows (-3r .. + 4r2 -r3 .)/2h 1,J ,J .J (3r +1 -4r + r 1 .) /2h m.J m,J m-,J ( 41) (-3r. 1 + 4r. 2 -r. 3)/2h 1, 1. 1, where r is a dummy variable. The truncation error in (37) is of 0(h2). The truncation error in (39) and (40) is of 0(h2 +k2). The velocity and dispersion tensor components appearing as coefficients in (39) and (40) were calculated for the time t +1' Therefore, these coeff; ci ents were obtained by n "2 linear extrapolation or interpolation before being applied in (39) and (40). For the initial conditions we assumed that the temperature is distributed according to Ti j = ih + 0.2(0.5 jh)sin(irrh) 107 (42)

PAGE 113

Such an assumption was recommended by Combarnous [1970J in a study concerning natural convection without hydrodynamic dispersion effects. The computation process was ended when maximum changes in and T in one time step was smaller than 105 However, we also followed changes in the Nusselt number, The Nusselt number was calculated along the row presenting the centerline of the cell through the following expression Nu -1 fL (-wT + E lL + E lL)dy [ zz az zy ay o (43) Through the calculations the time step was adjusted to assure convergence of the numerical computation. The mesh size was determined according to the variations in the Nusselt number. If the value of Nu was changed by less than 2% by changing M to M+2 (where M is the number of intervals in the vertical then it was assumed that the mesh size is sufficiently small. For R/Ro 6 the mesh size varied in the range 0.1 ; h < 0.05. In order to speed up the convergence, values of and T, obtained for negligible effects of the convection motion on the dispersion tensor were used as initial values for calculations considering these effects. NUMERICAL RESULTS AND DISCUSSION For the identification of the relationship between the velocity vector and the dispers;vities, we applied the model suggested by Saffman [1959J. According to this model e:t -v s + 2 (Re) 2{s+1){s+3) (44 ) 108

PAGE 114

0* _'-JI, (Re) v (45) where s is a power coefficient describing the dependence between the flow velocity and the pressure gradient. In a previous article [Rubin, 1976aJ it was suggested that s can be calculated through the following expression s = Q,n(Re) Q,n c( 1 +b ) Re J (46) The transversal and longitudinal dispersion coefficients are given by et 1 e"t -=-+v Pr v e Q, 1 e1 --= -+v Pr v where Pr is the Prandtl number. Through the calculations we assumed ( 47) Pr = 1.0 which seems to be a reasonable value for limestone and dolomite aquifers [Somerton, 1958J. In Fig. 2 we present the dependence between the Nusselt number and the convection cell width for various values of R/Ro and d/d p The effect of d/d p is significant when its value is comparatively small. It reduces the value of the Nusselt number as well as the width of the convection cell. The half width of the convection cell is the value of L associated with the maximal Nusselt number. This assumption was based on investigation concerning the ordinary Benard problem [SchUlter et al., 1965J as well as diffusive convection in porous media [Straus, 1974J, and cannot be proved through simple numerical experiments. However. through the calculation it was found that for wide ranges of L, the Nusse It number changed very 1 ittl e as shown in Fi g. 2. 109

PAGE 115

Fig. 2. Nu 4.0 --3.0 2.0 o. Nusselt numbers vs. half cell size for various values of R/R and did 0.4, Re = 3.0, Pr = 1.0) o p i jdJ.d_:Loo ____ : I p! :-p-: .-'cifd410 p. 0.2 0.4 0.6 110 .. : I 0.8 1.0 1.2 L

PAGE 116

At high Rayleigh numbers boundary layers are created at the boundaries of the convection cell. The simulation of processes occuring in these regions requires fine mesh size. Such requirements limit the applicability of the schemes (37)-(40) related to a constant mesh size for the whole gri d < Terms depending on the effect of the convection motion on the dispersion tensor have convectional nature. The magnitude of such terms mainly depends on t-he parameter d/d p' Therefore, simulation of processes related to small values of did require fine mesh size although in such p cases the Nusselt number is smaller than what expected in the diffusive convection. In Fig. 3 we present profiles of mean horizontal temperatures for various values of the Rayleigh number and the parameter d/dp In small values of d/dp mechanical dispersion due to the convection motion is significant. Fig. 3 demonstrates that such an effect is associated with a reduction in the values of the temperature gradients at the convection cell boundaries. In Fig. 4 is shown the dependence between the Nusselt number and the Rayleigh number for various values of d/d p According to this figure mechanical dispersion induced by the convection motion leads to a reduction in the effectiveness of transport processes through the aquifer. Such effects can be significant when the ratio between the aquifer thickness and the characteristic pore size is comparatively small. The validity of the calculations was checked by conventional methods as described in the previous section (changing time steps and intervals while following variations in the Nusselt number) as well as comparing 111

PAGE 117

Fig. 3. T z Profiles of mean horizontal temperatures for various values of did when R/R = 1, 3. ( = 0.4, Re =: 3.0, Pr = 1.0) p 0 112

PAGE 118

NU 4 3 2 1 1 2 3 4 5 6 R/R 0 Fig. 4. The Nusse1t number VS. R/R for various values of did (cp = 0.4, Re = 3.0, Pr = 1.0) 0 p 113

PAGE 119

the results with those obtained through Fourier series expansion. It was shown by Drs zag [1971J that spectral methods for numerical simulation of incompressible flows within simple boundaries lead to more accurate so1ution of the matnematical problem. In the case of thermal convection it is possible to utilize series expansions leading to analytical solution of the problem [Rubin, 1975; 1976bJ. By comparing our numerical results with the analytical ones we could present further evaluation of the applicability of the numerical schemes applied in this study. From such a comparison it was found that if mechanical dispersion effects are negligible the combination of the SDR and ADI can be applied up to R/Ro 6, The applicability of this procedure reduces if mechanical dispersion due to the convection motion is significant, CONCL US IONS A combination of the SOR and ADI methods can be applied for the analysis of singly dispersive convection in porous media. The procedure is practically applicable up to R/Ro 6. If mechanical dispersion due to the convection motion is significant the applicability of this method reduces, The requirement for coordinate intervals and time steps are determined by the formation of boundary layers at the convection cell boundaries as well as by the production of terms having the convection character, 114

PAGE 120

NOTATION a b k K L m M N Nu wave number critical wave number wave number components friction function defined in (5) specific heat of solid specific heat of water porous layer thickness characteristic pore size longitudinal dispersion coefficients (hydrodynamical and mechanical. respectively) transversal dispersion coefficients dimensionless longitudinal dispersion coefficient dimensionless transversal dispersion coefficient dispersion tensor coordinate interval unit matri x time step permeabi 1 i ty ha If cell wi dth unit vector in the longitudinal. x direction number of iteration number of intervals producing the cell height unit vector in the vertical, z direction number of intervals producing the cell width Nusselt number 115

PAGE 121

p R s t T T -+ u u v v w w -+ x x,y,z pressure pressure at the coordinates origin Prandtl number (=V/K) Rayl ei gh milTiber critical Rayleigh number Reynolds number power coefficient time temperature mean horizontal temperature temperature at z = a velocity vector velocity existing prior to the inception of the convection motion absolute value of the velocity vector velocity in y direction absolute value of the convection velocity velocity in z direction velocity perturbation in the x direction coordinates vector coordinates coefficient of volumetric thermal expansion a;(i=1 .. 4) coefficients defined in (35) 13 variable defined in (10) 131 coefficient defined in (32) y coefficient defined in (4) r dummy variable 116

PAGE 122

= E e \l \I p x w operator defined in (16) difference in temperature between bottom and top of the aquifer dispersion tensor resulting from the convection temperature perturbation viscosity kinematic viscosity fluid density density of reference solid density parameter of instability porosity parameter defined in (25) stream function relaxation coefficient function defined in (16) coefficient defined in (22) 117

PAGE 123

REFERENCES Aziz, K. and J. O. Hellum, IINumerical Solution of the Three Dimensional Equations of Motion for Laminar Natural Convectionll, The Phys. of Fluids, lQ(2) 1967, 314-324. Cabelli, A., and G:OeVahl Davis, IIA Numerical study of the Benard Cell, J. Fluid Mech 45(4), 1971, 805-829. Combarnous, M. A., IIConvection Naturelle et Convection Mixte dans Une Couche Poreuse Horizontalell, Revue Generale de Thermique, No. 108,1970, 1355-1375. Henry, H. R. and F. A. Kohout, "Circulation Pattern of Saline Groundwater Affected by Geothermal Heating as Related to Waste Disposal!!, Proceedings of 1st Symp. on Underground Waste Management and Environmental Implications, ed. by T. D. Cook, pp. 202-221, Amer. Assoc. Petrol Geologists, Tulsa, Oklahoma, 1972. Holst, P. H. and K. Aziz, IITransient Three Dimensional Natural Convection in Confined Porous Mediall, Int. J. Heat Mass Transfer, J2., 1972,73-90. Kohout, F. A., "A Hypothesis Concerning Cyclic Flow of Salt Water Related to Geothermal Heating in the Floridan Aquiferll, New York Acad. Sci. Trans., 1965, 249-271. Kohout, F. A., IIGroundwater Flow and the Geothermal Regime of the Floridan Plateau", Symp. Geological History of the Gulf of Mexico Caribbean Antillean Basin, Gulf Coast Assoc. Geol. Socs. Trans., 1I, 1967, 339-354. Orszag,S. A., "Numerical Simulation of Incompressible Flows Within Simple Boundar; es: Accuracyll, J. Fl ui d Mech., 49 (1), 1971, 75-112. Rubin, H., "On the Analysis of Cellular Convection in Porous Mediall, Int. J. Heat Mass Transfer, 1975, 1483-1486. Rubin, H., 1I0nset of Thermohaline Convection in a Cavernous Aquifer", Water Resources Res. 1976a, 141-147. Rubin. H., II Ther'ma1 Convection in a Cavernous Aquifer", to be published in Water Resources Research. 1976b. Saffman, P. G., "A Theory of Dispersion in a Porous Mediumll, J. Fluid Mech., 1959. 321-349. Schluter. A., D. Lortz, and F. Busse, liOn the Stability of Steady Finite Amp 1 itude Convect; on ". J. Fl u; d Mech., 23 (1). 1965, 129-144. Smith, G. D., Numerical Solution of Partial Differential Eguations. Oxford University Press, 1969. p. 149. Somerton, H. W., IISome Thermal Characteristics of Porous Rocks", J. Petrol. Tech., lQ, 1958, Note 2008, 61-65. 11Q

PAGE 124

Straus, J. M., "Large Amplitude Convection in Porous Medial!, J. Fluid Mech 64(1), 1974, 51-63. Vernon, R. 0 liThe Beneficial Uses of Zones of High Transmissivities in the Florida Subsurface for Water Storage and Waste Disposalll Florida Bur. Geol. Inf. Circ., No. 70, 1970. 119

PAGE 125

CHAPTER 7 SUMMARY AND CONCLUSIONS 1. The conceptua 1 des i gn framework of an urban dra i nage sys tem shou 1 d be basedO'nadequate attention to legal principles, and management programs. 2. An outlined diagram for urban drainage design is suggested. This diagram consists of three different phases: project feasibility, preliminary design and detailed design. The use of this diagram promotes more effective management through systematic and sequential consideration of design factors. 3. The adequacy of seepage pond to divert storm water to groundwater depends on two, not necessarily components: the ability of the pond to seep water and the growth and decay of groundwater mounds. Both of these aspects should be considered by the designers. 4. A simple meth6d based on the Green and Ampt equation is suggested for the analysis of the pond's ability to seep water. 5. The growth and decay of groundwater mounds can be calculated by a method developed by Hantush. However, this method is not always sufficiently accurate. Therefore, numerical solution of the nonlinear differential equation is frequently required. Design charts based on the numerical solution were developed in this study and can be used as an aid for pond designers. 6. Dispersion of pollutants in the Floridan aquifer, as well as possible invalidity of the laminar Darcy law, should be considered. 120

PAGE 126

7. Possible conditions for doubly or singly dispersive convection in the Floridan Aquifer should be considered. 8. Convective motions may significantly affect dispersion of contaminants in the aquifer. 9. A variety of methods may be used for the analysis of dispersive convection in the aquifer. The more efficient methods are those utilizing spectral expansion of flow field perturbations. These methods are, however, usually limited to simple boundary conditions. 121

PAGE 127

ERRATA SHEET Page 82 Lines 19, 20 should read: By expanding in the following normal mode, (18) Page 84 Line 18 should read: =-i [_ exp(ipay) ............ (36) p--oo q=l Page 86 Lines 19, 20 should read: Equation 44 can be applied for Re 2.77. Through the calculation. it was assumed that Pr = 1, which seems to be a reasonable value for limestone and Page 90 Line 16 should read: [w 81 k I' Wv ( v= 1 ,3,5,7 ; 2 ,4,6) J} . . . . (50) v p--m Page 92 Line 22 should read: 8 = temperature perturbatiun; Page 94 Lines 32, 33 should read: 11. Saffman. P. G., "A Theory of Dispersion in a Porous Medium", Journal of Fluid Mechanics, Vol. 6, Pt. 3, 1959, pp. 321-349. Page 96 Line 22 should read: BASIC EQUATIONS Page 98 Line 18 should read: Substituting the dimensionless variables of (7) in (1)-(3) and omitting