PB85-2153A1
UPCOMING OF A SALT-WATER/FREST-WATER INTERFACE BELOW
A PUMPING WELL
Oklahoma State University
Stillwater, OK
Jun 85
U.S. DEPARTMENT OF COMMERCE
National Technical Information Service
NTIS
-------
EPA/600/2-85/066
June 1985
UPCONING OF A SALT-WATER/FRESH-WATER
INTERFACE BELOW A PUMPING WELL
by
Jan Wagner
School of Chemical Engineering
Oklahoma State University
Stillwater, Oklahoma 74078
Douglas C. Kent
Department of Geology
Oklahoma State University
Stillwater, Oklahoma 74078
CR811142
Project Officer
Carl G. Enfield
Robert S. Kerr Environmental Research Laboratory
U. S. Environmental Protection Agency
Ada, Oklahoma 74820
ROBERT S. KERR ENVIRONMENTAL RESEARCH LABORATORY
OFFICE OF RESEARCH AND DEVELOPMENT
U.S. ENVIRONMENTAL PROTECTION AGENCY
ADA, OK 74820
-------
TECHNICAL REPORT DATA
(Please read Instructions on the reverse before completing)
1. REPORT NO.
EPA/600/2-85/066
3. RECIPIENT'S ACCESSION-NO.
PB3 q x i i? 3 •'. - /AS
4. TITLE AND SUBTITLE
UPCONING OF A SALT-WATER/FRESH-WATER INTERFACE
BELOW A PUMPING WELL
5. REPORT DATE
June 1985
6. PERFORMING ORGANIZATION CODE
7. AUTHOR(S)
8. PERFORMING ORGANIZATION REPORT NO.
Jan Wagner, Douglas Kent
9. PERFORMING ORGANIZATION NAME AND ADDRESS
Oklahoma State University
Stillwater, OK 74078
10. PROGRAM ELEMENT NO.
ABRD1A
11. CBKXKXKTCKXXNXXtO. UOOp. Mgr.
CR811142
12. SPONSORING AGENCY NAME. AND ADDRESS
Robert S. Kerr Environmental Research Laboratory
Office of Research and Development
U.S. Environmental Protection Agency
Ada, OK 74820
13. TY.IJE Of REPORT AND PERLPD COVERED
Final 9/83 - 2/85
14. SPONSORING AGENCY CODE
EPA/600/15
15. SUPPLEMENTARY NOTES
Carl G. Enfield, Project Officer
16. ABSTRACT
Analytical solutions for the upconing of an abrupt salt-water/fresh-water
interfere beneath a pumping well and for the concentration profile across a moving
interface are developed for two types of upconing problems. The first considers
the position of the interface and the salinity of the pumped water for a specified
pumping rate. The second type of problem addresses the pumping schedules to
prevent salinization of a well or to reach a predetermined salinity in the pumped
water.
- An interactive Fortran computer code has been developed to obtain solutions to
both types of problems. The user is provided with options to modify the definition
of a given problem, and, therefore, can gain some insight into the effects of
geometry and physical properties on the rate and extent of upconing and the
salinization of a well.
17.
KEY WORDS AND DOCUMENT ANALYSIS
DESCRIPTORS
b.lDENTIFIERS/OPEN ENDED TERMS C. COSATI Field/Croup
13. DISTRIBUTION STATEMENT
RELEASE TO PUBLIC
19. SECURITY CLASS (This Report)
UNCLASSIFIED
21. NO. OF PAGES
78
20. SECURITY CLASS (Thispage)
' UNCLASSIFIED
22. PRICE
EPA Form 2220-1 (9-73)
-------
DISCLAIMER
The information in this document has been funded wholly or in part
by the United States Environmental Protection Agency under assistance
agreement number CR-811142 to Oklahoma State University. It has been
subject to the Agency's peer and administrative review, and it has been
approved for publication as an EPA document.
-------
FOREWORD
EPA is charged by Congress to protect the Nation's land, air, and water
systems. . Under a mandate of national environmental laws focused on air and
water quality, solid waste management and the control of toxic substances,
pesticides, noise, and radiation, the Agency strives to formulate and imple-
ment actions which lead to a compatible balance between human activities and
the ability of natural systems to support and nurture life.
The Robert S. Kerr Environmental Research Laboratory is the Agency's
center of expertise for investigation of the soil and subsurface environment.
Personnel at the Laboratory are responsible for management of research pro-
grams to: (a) determine the fate, transport and transformation rates of
pollutants in the soil, the unsaturated zone and the saturated zones of the
subsurface environment; (b) define the processes to be used in characterizing
the soil and subsurface environment as a receptor of pollutants; (c) develop
techniques for predicting the effect of pollutants on ground water, soil and
indigenous organisms; and (d) define and demonstrate the applicability and
limitations of using natural processes, indigenous to the soil and subsurface
environment, for the protection of this resource.
This project was initiated to develop an interactive computer model which
could be utilized to predict the impact of density gradients on water flow.
The model should be useful in regard to salt water intrusion and salt water
injection problems. This model assumes idealized conditions of a homogeneous
isotropic media. Care should be utilized if there is significant heterogeneity
in the.aquifer.
Clinton W. Hall
Director
Robert S. Kerr Environmental
Research Laboratory
m
-------
ABSTRACT
Analytical solutions for the upconing of an abrupt salt-water/fresh-
water interface beneath a pumping well and for the concentration profile
across a moving interface are developed for two types of upconing prob-
lems. The first considers the position of the interface and the salinity
of the pumped water for a specified pumping rate. The second type of prob-
lem addresses the pumping schedules to prevent salinization of a well or to
reach a predetermined salinity in the pumped water.
An interactive Fortran computer code has been developed to obtain
solutions to both types of problems. The user .is provided with options
to modify the definition of a given problem, and, therefore, can gain some
insight into the effects of geometry and physical properties on the rate
and extent of upconing and the salinization of a well.
IV
-------
TABLE OF CONTENTS
Abstract 1
Table of Contents 2
List of Tables 3
List of Figures 4
Introduction 5
Section I - Mathematical Development 6
Upconing of an Abrupt Interface Dispersion 6
Superposition of Dispersion on the Upconing of an
Abrupt Interface 12
Concentration in Pumped Water 13
Section II - Applications 15
Case I - Estimation of Interface Elevations for a
Given Pumping Rate 15
Case II - Estimation of the Maximum Permissible
Pumping Rate to Prevent Salinization of the Well 16
Example Problem - Upconing Below a Coastal Collector Well 17
Section III - Program UPCONE 27
Basic Input Data 27
Edit Commands 31
References 37
Appendix A - Example Problems 38
Appendix B - Description of Program UPCONE 60
Appendix C - Listing of Program UPCONE 61
Appendix D - Listing of Function Subroutines 71
-------
LIST OF TABLES
Table 1 - Parameters for Semadar 1 - Test B 18
Table 2 - Salinity/Maximum Pumping Rate and Salinity/Time
Relationships for Semadar 1 (Xcr f 0.4d) '. 25
Table 3 - Edit Commands for UPCONE ; 36
-------
LIST OF FIGURES
Figure 1 - Upconing of an abrupt interface below a pumping well 7
Figure 2 - Predicted and observed interface elevations for Test B of
Semadar 1 20
Figure 3 - Predicted and observed upconing and recovery curves for
Test B of Semadar 1 21
Figure 4 - Concentration profile across initial transition zone 23
Figure 5 - Predicted and observed chloride concentrations in pumped
water for Test B on Semadar 1 24
VII
-------
INTRODUCTION
Relatively simple analytical models can often be used to solve ground-
water contamination problems, depending upon the complexity of the system
and the availability of field data. Analytical models can also be used to
gain some insight to the expected behavior of a complex system before pro-
gressing to the application of more sophisticated numerical models. In
general, relatively few input parameters are required to define a problem
using an analytical model and numerical results can be calculated in a few
seconds. Analytical models are well suited for interactive use, and in some
instances can be programmed on hand-held calculators.
This report presents an analytical solution for the upconing of an
abrupt salt-water/fresh-water interface below a pumping well. Dispersion
phenomena arising from the displacement of a moving interface, or a finite
transition zone between the invading and displaced fluids, can be superim-
posed on the analytical solution for the position of an abrupt interface.
An interactive Fortran computer code has been developed which enables the
user to modify input parameters and to control the computational sequence.
This interactive approach enables the user to gain insight into the effects
of geometry and physical properties on the rate and extent of upconing and
salinization of a well.
-------
SECTION I
Mathematical Development
McWhorter (1972) presented the equations which describe the flow in
saturated aquifers which are underlain by a zone of saline water and pointed
out the difficulties in obtaining solutions to these problems. The complex-
ity of the flow phenomenon has led many investigators to idealize the system
as a fresh-water zone separated from an underlying salt-water zone by a sharp
interface. In other words, the two fluids are assumed to be immiscible.
Schmorak and Mercado (1969) followed this approach and accounted for the mix-
ing of the two fluids by superimposing the effects of dispersion on the
transient solution for the position of an abrupt interface.
Upconing of an Abrupt Interface
The following discussion is based on the studies of Bear and Dagan as
reported by Schmorak and Mercado (1969). The basic assumptions underlying
the theoretical development are: (1) the porous medium is homogeneous and
nondeformable, (2) the two fluids are incompressible, immiscible, and sep-
arated by an abrupt interface (a geometric surface), and (3) the flow obeys
Darcy's law. The non-linear boundary condition along the interface between
the two fluids constitutes the major difficulty with the immiscible formula-
tion of the problem. Bear and Dagan used the method of small pertubations
to obtain an approximate solution for the position of the interface which
served as a tool for obtaining analytical solutions for cases involving
small deviations from an initially steady interface.
For the case of upconing beneath a pumping well partially penetrating
a relatively thick confined aquifer as shown in Figure 1, Schmorak and
2
-------
UJ
o
z
CO
Q
x=o
STATIC HEAD
DYNAMIC HEAD
..-IMPERVIOUS
FRESH WATER
DENSITY:/?,
ZONE OF
SUDDEN
RISE
VALIDITY
OF LINEAR
APPROXI-
MATION
ls»d
INITIAL INTERFACE POSITION
SALT WATER
DENSITY: ps
\; IMPERVIOUS
N
z~
o
UJ
_l
UJ
ZCr =XCr
RADIUS, r
00
Figure 1. Upconing of an abrupt interface below a pumping well.
-------
Mercado (1969) presented Bear and Dagan's solution for the position of the
interface as a function of time and radial distance from the pumping well
•.
as
Q
27T(Ap/p).K d
X
1
(1 + R2)15
fa.
1
T 2 + R2]15
where R and T are dimensionless distance and time parameters defined by
f •* X
K
R =4
K (2)
and
(Ap/p)K
Other notations are defined as follows (also refer to Figure 1) :
d distance from the bottom of the well to the initial
interface elevation (L)
K ,K horizontal and vertical permeabilities, respectively
X Z (L/t)
Q well pumping rate (L /t)
r radial distance from well axis (L)
t time elapsed since start of pumping (t)
X rise of the interface above its initial position (L)
Ap/p dimensionless density difference between the two
fluids, (pg - pf)/pf
9 porosity of the aquifer.
Application of the method of small pertubations restricts changes in
the interface elevation to relatively small values . In terms of the physi-
cal problem, this restriction implies d«lf. and d«l . Although the govern-
ing differential equations have been formulated for a confined aquifer, the
results can be applied to unconfined systems if the drawdown is negligible
compared to the saturated thickness of the fresh-water zone.
-------
The linear relationship, Equation 1, between the rise of the inter-
face and the pumping rate is limited to a certain "critical rise," X
cr
This limitation arises from linear approximation of the boundary conditions.
As the interface approaches this critical rise, the rate of rise increases.
Above the critical rise the interface reaches the pumping well with a sud-
den jump. Muskat (1946) defines the zone of accelerated rise for X/d > 0.48
and the critical rise within the limits of X/d ^ 0.60 to 0.75. Schmorak and
Mercado (1966) recommend application of the linear approximation for X/d < ^
0.5. Sahni (1972) investigated the zone of instability of the interface us-
ing both numerical and physical models and recommended design criteria for
skimming wells.
An abrupt interface such that (1) salinization of the pumping well
occurs only for X > X =fd where f is the fractional critical rise, and
cr
(2) Equation 1 is valid for 0 _< X _< X will be assumed in this report.
CL
Thus, the maximum permissible pumping rate which will ensure salt-free
water can be obtained from Equation 1.
For r = 0 and t -»• °°
X(°'M) = 2.d(ApVp)K (4)
X
and
(5)
For a decaying mound, an imaginary recharge well is superimposed at
* *
time t = t corresponding to the end of the pumping period, and for t > t
X(r,t)
where
Q
2ir(Ap/p)K d
X
1
K>'
+ .f fc
i
2 2\h
1+r) + R
-------
(Ap/p)K
Values of R and T are evaluated using Equations 2 and 3 respectively.
Dispersion
The upconing process as treated above assumes that the two fluids are
immiscible and that the interface between them is abrupt. Actually, the' in-
terface is diffuse and a transition zone exists between the two fluids in
which the concentration varies from the concentration in one fluid to the
concentration in the other fluid over a finite distance. This transition
zone is related to dispersion processes which alter the concentration pro-
file across the moving interface.
Bear and Todd (1960) approximated the concentration profile as a func-
tion of position, X; the "interface" position, X; the equivalent total dis-
tance the interface is displaced, |AX|, independent of direction; and the
dispersivity, D • The correlation is given by
e(X) = \ ERFC X"X_
2(D I AX I)
m1 '
where e(X) is a dimensionless, or relative, concentration defined as
C - Q
(8)
e(x) =
C = measured concentration at X
C = concentration of the invading fluid
s
C, = background concentration of the displaced fluid.
Now, Equation 9 is a normal distribution function with a mean X and a stan-
dard deviation
a = (2DJAXI)1* (10)
or
-------
p fy - 1 — 1
r IX > xj- i
~~ (2~)
oo
EXP
X
(x-x) 2]
2a2 J
dx
(11)
From the definition of the error function,
X-X 1
{X _> x} = - ERFC
(12)
a /2~J
The two-parameter distribution is completely defined once the mean, X, and
standard deviation, a, are known.
The mean of the distribution is assumed to be the rise of the inter-
face as determined from Equation 1, or
X = X| = X(r,t) (13)
The standard deviation is defined as
= |XU=0.841-XU.159
(14)
Now, the width of the transition zone is a function of the total distance
traveled, |AX|, (independent of direction) and the dispersivity as given
by Equation 10, or
20 = 2(2D
m
(15)
For an interface with an initial transition width, 2a , when raised by a
o
distance, AX,
The concentration distribution function then becomes
e (X) =4 ERFC
(16)
(17)
or
e(X) = -z ERFC
X-X
2a 2 + 4D I AX I
o m1 '
(18)
-------
Two important points should be noted concerning the preceding discus-
sion of dispersion. First, the "initial width" of the transition zone has
been defined as two standard deviations of the dimensionless concentration
distribution. This definition has been adopted for convenience and serves
to. define the standard deviation of the concentration distribution across
the initial transition zone. Secondly, the dispersion concept should be
limited to the zone below the critical depth. This point will be consid-
ered in more detail in the following paragraphs.
Superposition of Dispersion on the Upconing of an Abrupt Interface
The position of the interface as a function of time and radial distance
from the well is evaluated using Equation 1, which assumes an abrupt inter-
face between the two fluids. This elevation is assumed to correspond to
x| _n c» or the mean of the concentration distribution across the transition
zone. In other words,
X = X(r,t)
(l+T) +
(6)
assuming an abrupt interface. The effect of dispersion arising from the
displacement of the interface by a distance
|AX| = |x(r,t £ t*)| + |X(r,t > t*)| (20)
is superimposed to estimate the concentration distribution across the inter-
face using Equation 18, or
X - X
e(X) = -=• ERFC
2a 2 + 4D I AX
o m'
(18)
The only difficulties in the approach occur for e(X) = 0.0 and e(X) = 1.0.
Since
e(X) = 0 for X -»• »
8
-------
and
e(X) = 1.0 for X ^ 0 (22)
the transition zone would have an infinite width in theory. To overcome
this physical impossibility, the width of the transition zone is arbitar-
ily set at five standard deviations. This range includes approximately
99 percent of the area under the concentration distribution curve. Thus
e(X) = 0 for X = 2C + 2.5.0- (23)
and
e(X) = 1 for X = X - 2.5 o^ (24)
Note that these limits differ from those used to define the "initial
width" of the transition zone defined by Equation 14, or
(25)
Concentration in Pumped Water
The increase in concentration, or salinization, of pumped water is
probably due to the intrusion of invading fluid above the critical depth.
Data for two pumping tests on a coastal aquifer in the Ashqelon area of
Israel indicated that the increase in salinity of the pumped water was
approximately proportional to the average salinity above the critical depth
(Schmorak and Mercado, 1969).
Previous discussion has emphasized that the linear approximation for
the interface elevation is limited to elevations below the critical eleva-
tion and that the dispersion concept should be limited to the zone below
the critical depth. The complex mixing and flow phenomena above the criti-
cal depth, near the well screen, and within the well pipe are approximated
expirically using the approach followed by Schmorak and Mercado (1969).
-------
The average dimensionless concentration of the transition zone above
the critical rise, e(X > X ), is approximated as one-half the concentra-
cr
tion at the critical depth, or
e(X > X ) = 0.5 e(X ) (26)
cr cr
The concentration in the pumped water, e , is determined from dilution of
w
the average transition-zone concentration above the critical depth with
displaced fluid, or
e = I(X > X ) (27)
w T cr
where is an interception coefficient, or the fraction of transition zone
fluid in the total volume pumped.
10
-------
SECTION II
Applications
Two types of upconing problems are considered. The first involves the
description of the expected interface elevation and the salinity of the pump-
ed water as a function of time for a given pumping rate. The second problem
addresses the maximum rate at which a well can be pumped without exceeding a
specified salinity in the pumped water. Both types of problems are discussed
in the following paragraphs.
Case I - Estimation of Interface Elevations for a Given Pumping Rate
Case I problems are solved in a fairly straight-forward manner. Once
the physical properties of the aquifer and the initial conditions have been
specified, Equation 1 or Equation 6 is solved for the position of the abrupt
interface, or the mean of the transition zone, i.e.,
X(r,t) = X|£=0_5 = X (13)
In terms of elevations,
Z(r,t) = X(r,t) + ZQ (28)
where Z is the initial interface elevation.
o
The concentration profile across the transition zone is evaluated using
Equation 18. • The salinity of the pumped water is determined from the dilution
of the average transition zone salinity above the critical rise. From Equa-
tions 26 and 27
e = 0.5 <}) e(X ) (29)
w cr
and from the definition of dimensionless concentration (Equation 9)
C = e (C - C, ) + C, (30)
w w s b b
11
-------
where C is the concentration in the pumped water.
w
Case II - Estimation the -Maximum Permissible Pumping Rate to Prevent
Salinization of a Well
Case II problems present some difficulty as the pumping rate, Q, is
unknown. Thus, the elevation of the interface, X, and the total displace-
ment of the interface, |AX|, are also unknown. Equation 18 must be solved
for the maximum permissible rise in interface elevation such that
e(X =fd) < e (31)
cr — max
where
£max = 0 (32)
Assuming a constant, steady pumping rate the total displacement of the
interface will be equal to the rise of the interface or
| AX | = X
and Equation 18 can be written_as
e ' = -r ERFC
max 2.
X - X
1
f 2 _
20 + 4D X
o m
(33)
Equation 33 must be solved for X using trial-and-error procedures. The maxi-
mum permissible rise in the interface elevation,
X < X =fd (34)
max — cr
is then corrected for the concentration profile as
X* = X - (fd-X) (35)
max max
and
Z* = X* + Z (36)
max max o
12
-------
The maximum permissible steady-state pumping rate is then obtained using
Equation 5, or
Q" = 2ad(Ap/p)K X" (37)
max x max v '
where X depends only upon the critical rise and the dispersion pattern.
The time required to reach a predetermined salinity in the well can be
estimated by rewriting Equation 1 as
t(C ) = TTTT^T- f 1 I— 1 (38)
26d
(AP/P)KZ
i
1 - f2ir(Ap/P)K d
1 X
* >\
X \
max
/Q
1
Substituting Equation 37 into Equation 38 yields
t(C
26d
V (Ap/p)K_ -- ,
vmax y
- 1
(39)
which can be used to estimate the time required to reach a predetermined sa-
linity in the pumped water for. pumping rates, Q, greater than the maximum
*
steady-state pumping rate, Q
J f t- & xmax
An interactive computational code has been developed to calculate inter-
face elevations, and concentrations for both Case I and Case II problems using
the approach described above. The computer program is discussed in Section IV
of this report.
Example Problem - Upconing Below a Coastal Collector Well
The application of the analytical model will be demonstrated using the
field data for Test B on the coastal collector well Semadar 1 in the Ashqelon
area of Israel (Schmorak and Mercado, 1969). Test B consisted of pumping
Semadar 1 at a rate of 348 m /day for a period of 84 days. The upconing and
decay of the salt-water/fresh-water interface were monitored by measuring the
13
-------
TABLE 1
Parameters for Semadar 1 - Test B
3
Fresh-water density, p. 1.00 g/cm
Salt-water density, p 1.03 g/cm
s
Porosity, 9 0.33
Horizontal permeability, K 14.7 in/day
X
Vertical permeability, K 14.7 m/day
z
Initial interface elevation, Z -30.75 m MSL
o
Distance from bottom of well to initial
interface, d 15.5 m
Fractional critical rise, f 0.4
Chloride concentration in'salt water, C 22,000. ppm Cl
S
Chloride concentration in fresh water, C, 145. ppm Cf
Dispersivity, D 0.5 m
Initial width of transition zone, 2a 3.5 m
o
Interception coefficient,
-------
salinity profiles in four observation wells. Samples of the pumped water
were collected periodically and analyzed for chlorides. The properties of
the fluids and the aquifer as estimated by Schmorak and Mercado are summa-
rized in Table 1.
Program UPCONE was used to calculate the transient interface elevations
and chloride concentrations in the pumped water. The input data dialog and
printed results for this problem are presented in Appendix D.
The predicted and observed interface elevations after 16, 57, and 84
days of pumping are shown in Figure 2. The observed interface elevations
correspond to elevations for a 50 percent relative concentration of sea
water, or
C - C,
C -
C '
Cb
The predicted and observed interface elevations match fairly well with the
exception of the values for observation well T-2. Well T-2 is located 4.5
meters from the pumping well and penetrates to an elevation of 31.10 meters
below MSL. However, this well is screened to an elevation of only 29.02
meters below MSL, and Schmorak and Mercado (1969) indicate that saline water
was entrapped at the bottom of the well from a previous pumping test. The
predicted rise of the interface at well T-2 also approaches the critical
elevation for a fractional critical rise of 0.4. Thus, the observed inter-
face elevation could be in a zone of accelerated rise. Finally, the reported
concentration gradients were very steep in well T-2, and small errors in con-
centration measurements could lead to large errors in. estimating the position
of the interface.
The predicted upconing and recovery curves for each of the four obser-
vation wells are shown in Figure 3. With the exception of well T-2 the
15
-------
SEMADAR 1
T/2
-------
20
24
28
32
0)
i_
CD
I I I I i
I
T-2
r =4.5m
L CRITICAL ELEVATION
FOR Xcr/d = 0.4
H 1 1 1 1 I-
-i 1 1-
H 1
I 24
111
28
LU
co 32
z
LU
S 24
-LU
co 28
z
o
< 32
LU
_J
LU
24
28
32
T-3
r =12.4m
H 1 1 1 1 1
T-4
r =16.7m
H 1 1 h
H 1-
T-5
r =33.9m
50 100
TIME (days)
150
Figure 3. Predicted and observed upconing and recovery curves for Test B of
Semadar 1.
17
-------
predicted upcoming curves follow the observed interface elevations fairly
well. The recession curves for all four observation wells are also in fair
agreement with the field data.
Initial observed relative concentrations for wells T-3, T-4, arid T-5
are plotted on Figures 4a and 4b. The predicted concentration profile across
the transition zone using an initial transition zone width, 2o , of 3.5 meters
is also shown. This value represents an average of the widths of the tran-
sition zones.at the three observation wells.
The parameters listed in Table 1 were used to predict the concentration
of chlorides in the pumped water as a function of time. The results of the
simulation are summarized in Figure 5 and agree very well with the observed
values.
No effort has been made in this report to "optimize" model input para-
meters or to quantify the "goodness of fit" between observed and predicted
f
values of elevations, concentration profiles, or salinity of the'pumped water.
However, a qualitive comparison of the predicted and observed values as shown
in Figures 2, 3, and 5 indicate that the assumptions incorporated in the ana-
lytical model approach the field conditions.
Program UPCONE was also used to develop salinity/maximum pumping rate
relationships and salinity/time relationships for Semadar 1. These relation-
ships correspond to Case II types of problems. The corrected critical inter-
* *
face elevations, Z , and maximum steady-state pumping rates, Q , for sev-
max f f e> xmax
eral values of predetermined salinity in the pumped water are presented in
Table 2. The time required to reach the specified salinity in the pumped
water were also calculated for two optional pumping rates. These pumping rates
correspond to Test A and Test B of Semadar 1.
18
-------
26
o
28
(0
2
o 30
_l
UJ
CD
HI
_l
UJ
<0
b.
o
o
E
^^
-i
(0
2
$
o
_l
UJ
00
z
o
H
<
>
UJ
32
34
®T-3
AT-4
HT-S
1 2 5 10 20 40 60 80 90 95
RELATIVE CONCENTRATION (PERCENT)
.2 .4 .6 .8
RELATIVE CONCENTRATION
1.0
Figure 4. Concentration profile across initial transition zone,
(a) Relative concentration on probability scale.
(b) Relative concentration on arithmetic scale.
19
-------
800
2o
I- E 600
OC a
I- a
Z~
HI DC
O u 400
O >
UJ O
Q W 200
5^
o o
i I
20
40
TIME ( days )
60
80
Figure 5. Predicted and observed chloride concentrations in pumped water for Test B on Semadar 1.
-------
Relative
Concentration
of Salt Water
(Dimensionless)
0.001
.003
.005
.010
.020
.030
TABLE 2
Salinity/Maximum Pumping Rate and Salinity/Time
Relationships For Semadar 1 (X = 0.4d)
cr
Chloride
Concentration
In Pumped
Water
(ppm CL)
166.85
210.56
254.27
363.55
582.10
800.65
X
max
(meters)
2.85
3.73
4.30
5.36
7.20
9.491"
max
(m /day)
79.57
117.40
141.66
187.34
266.28
364.76f
Time
Test A
Q = 575 m3/day Q
(days)
3.73
5.95
7.58
11.21
20.01
40.25
Test B
= 348 m3/day
(days)
6.88
11.81
15.92
27.05
75.59
_
(t) X > X 6.20 meters
max cr
-------
The example problems using the data for Test B of Semadar 1 described
above are intended to support, in general, the validity of the theoretical
approach and to demonstrate the application of the analytical model to a
typical upconing problem. The reader interested in methods which might be
used to develop input parameters for the model are referred to the Schmorak
and Mercado (1969) discussion of the field investigation and interpretation
of the field data.
22
-------
SECTION III
Program UPCONE
Program UPCONE evaluates the position of an abrupt salt-water/fresh-
water interface beneath a pumping well as a function of time and radial
distance from the well. The program has also been written to (1) superim-
pose the effects of dispersion on the abrupt interface to estimate the con-
centration profile across the interface and the salinity of the pumped
water, or (2) estimate the maximum pumping rate or time required to reach a
specified salinity in the well. The program has been designed for interac-
tive use and requires data input under two modes of operation—"Basic Input
Data" and "Edit."
Basic Input Data
Basic input data are required to initiate a new problem using the UP-
CONE program. The data entries include the problem title, the physical prop-
erties of the aquifer and the two fluids, and the geometry of the system.
The user is prompted for the required data through a series of input commands
described below. Numeric data should be entered through the keyboard with
decimal points, and multiple data entries should be separated by a comma(s).
The first basic input command is
ENTER TITLE
7
Any valid keyboard characters can be used. The first 60 characters will be
retained for further problem identification.
23
-------
The next two input commands are used to define the units of all variables
used in the calculations and printouts of the results. Any consistent set of
units may be used. The two commands are
ENTER UNITS FOR LENGTH (2 CHARACTERS)
?
ENTER UNITS FOR TIME (2 CHARACTERS)
?
Any valid keyboard characters can be used. The first two characters will be
retained for identification of length and time units.
The next series of input commands are used to specify the physical prop-
erties of the fluids and the aquifer. Input data errors which may interrupt
the computational sequence are detected by the program and a command is issued
to reenter the data for the appropriate variable. The series of commands are
as follows:
ENTER FRESH-WATER AND SALT-WATER DENSITIES
? ?
The densities, or specific gravities, of the two fluids may be entered in any
units so long as the units are identical for both entries.
ENTER AQUIFER POROSITY
?
Enter the volume void fraction as a decimal value greater than zero and less
than one.
ENTER HORIZONTAL PERMEABILITY (L/t)
24
-------
ENTER VERTICAL PERMEABILITY (L/t)
?
Horizontal and vertical permeabilities must be entered with dimensions of
L/t in the units requested. Numerical values for both entries must be
greater than zero.
ENTER INITIAL INTERFACE ELEVATION (L)
7
The initial interface elevation may be either positive, zero, or negative,
depending upon the location of the reference elevation with respect to the
initial abrupt interface (or the elevation of the mean concentration of the
transition zone). The elevation must be entered in the units requested.
ENTER DISTANCE FROM BOTTOM OF WELL TO INITIAL INTERFACE (L)
?
The entry must be positive and in the units requested.
ENTER FRACTIONAL CRITICAL RISE
The fractional critical rise must be a decimal value greater than zero and
less than one.
The next basic input command is used to select the option of perform-
ing concentration calculations. The command is
CONCENTRATION CALCULATIONS? (Y/N)
If the user does not respond with Y, the problem title and parameters which
have been specified are listed. The critical rise and the maximum steady-state
25
-------
pumping rate which will maintain the interface at the critical elevation are
evaluated and the results are printed. The program then enters the "Edit"
mode.
If the response to the last basic input command is Y, the user will be
requested for additional basic input data required to carry out the concen-
tration calculations. The first command is
ENTER UNITS FOR CONCENTRATION (6 CHARACTERS)
?
Any valid keyboard characters can be used. The first 6 characters will be
used to specify concentrations for following data entries and printed results.
The following commands are:
ENTER SALT-WATER AND BACKGROUND CONCENTRATIONS (M/L3)
The concentration of any desired component in the invading fluid and the dis-
placed fluid are entered. If the user wishes to work in terms of dimension-
less concentrations, enter a value of 1.0 for the salt-water concentration and
a value of 0.0 for the background concentration.
ENTER DISPERSIVITY (L)
?
The dispersivity has dimensions of L and must be entered in the units re-
quested. Numerical values must be greater than zero.
ENTER INITIAL WIDTH OF TRANSITION ZONE (M)
26
-------
The width of the transition zone is defined as two standard deviations of
the concentration distribution function across the transition zone, as dis-
cussed in Section II of this report. For an initially abrupt interface
enter a value of zero.
ENTER INTERCEPTION COEFFICIENT
The interception coefficient is the fraction of transition zone water in
the total volume pumped. Enter a decimal fraction greater than zero and
less than one.
At this point the program will list the problem title and parameters
as they are currently specified along with the critical rise and the maxi-
mum permissible pumping rate for an abrupt interface. The program then
enters the "Edit" mode.
Edit Commands
Once the basic input data have been entered, the problem as currently
defined is listed and the program enters the "Edit" mode. The edit commands
are listed in Table 3- The request for information is
ENTER NEXT COMMAND
?
One of the responses from Table 1 should be given. If the response is in-
correct or improperly formated the statement
ERROR IN LAST COMMAND—REENTER
27
-------
is issued. Error messages for invalid numeric data will be issued as de-
scribed under Basic Input Data.
The request for information will be repeated until one of the responses
EL, PR, LI, NP or DN is entered.
EL will initiate the calculation of interface elevations
for a specified pumping rate. The user is given the
option of calculating the concentration profile across
the transition zone by responding to the following
command
CONCENTRATION CALCULATIONS ? (Y/N)
If the response is Y and the data required for concen-
tration calculations have not been entered previously,
the user is prompted for the required information using
"Basic Input" commands.
For the initial use of the EL edit command, the following
requests for data are issued:
ENTER PUMPING RATE (CU L/t) AND PERIOD (t)
? ?
• > •
Enter the well pumping rate and the pumping period in
the units requested. Both entries must be positive and
separated by a comma.
Two additional requests for data are used to define the
coordinates of the observation points. The first is
ENTER TFIRST, TLAST, DELTAT (t)
23
-------
Input units for the time variables must be in the units
requested. TFIRST must not be negative value. A zero
entry for DELTA! will result in interface elevations at
a single time. The second request is
ENTER RFIRST, RLAST, DELTAR (L)
? ? ?
. , . , .
The numerical values used to define the radial coordi-
nates of the observation points may be positive or nega-
tive. The results of the calculations will be printed
from RFIRST to RLAST.
The pumping rate and observation coordinate parameters
are listed and a request
CONTINUE ? (Y/N)
is issued. If the response is not Y, the program returns
to the edit mode. If the response is Y, the program pro-
ceeds to the computation of interface elevations at the
specified times and radial distances from the well and
prints the results. If concentration calculations were
requested, the concentration in the pumped water and the
concentration distribution across the transition zone are
also evaluated and printed for the specified times.
On subsequent use of the EL edit command, the pumping rate
and observation coordinate data as currently defined will
be listed and a request to continue will be issued.
29
-------
PR initiates the calculation of maximum permissible inter-
face elevations and pumping rates for a specified con-
centration in the pumped water. If the data required
for concentration have not been entered previously, the
user is prompted for the required data using "Basic Input
Data" dialog. The following request for information is
then issued:
ENTER MAXIMUM PERMISSIBLE CONCENTRATION IN PUMPED WATER (M/L3)
Enter the concentration in the units requested. The numer-
ical value must not be negative.
The program then lists the problem as currently defined
and evaluates the maximum permissible interface elevation
and pumping rate. The following request for information
is then issued:
ENTER OPTIONAL PUMPING RATE (L3/t)
If a value greater than the maximum pumping rate is entered,
the time to reach the specified concentration in the pumped
' water will be calculated, and the command reissued. If a
value less than the maximum pumping rate is entered the pro-
gram returns to the edit mode.
LI will list the problem as currently defined.
30
-------
NP will request a complete new problem definition using the
"Basic Input Data" dialog.
DN will terminate the program.
Although many tests for valid input data and properly formulated edit
commands have been embedded in the computer code, the user is encouraged to
correct "keyboard errors" before the data are transmitted. This practice
will serve to minimize the frustration of program termination as a result
of fatal errors during execution of the numerical computations.
31
-------
TABLE 3
«.•
Edit Commands for UPCONE
Command Variable Changed/Execution
CO Salt-water and background concentrations
CR Fractional critical rise
DI Dispersivity
DT Depth from bottom of well to initial interface
FD Fluid densities
1C Interception coefficient
KX Horizontal permeability
KZ Vertical permeability
• OB Time and radius coordinates
PO Porosity
QP Pumping rate and time
- RC Radius coordinates
TC Time coordinates
TW Initial width of transition zone
ZO Initial interface elevation
EL Interface elevation calculations
DN Terminate program
LI List problem definition
NP New problem
PR Pumping rate calculations
32
-------
REFERENCES
Abramowitz, M. and I. A. Stegun. 1966. Handbook of Mathematical Functions
with Formulas, Graphs, and Mathematical Tables. National Bureau of
Standards Applied Mathematics Series 55, U. S. Department of Commerce,
1046 pp.
Bear, J. and D. K. Todd. 1960. "The Transition Zone between Fresh and
Salt Waters in Coastal Aquifers." Contribution No. 29, Water Re-
sources Center, University of California, Berkeley, CA.
Carnahan, B., H. A. Luther and J. 0. Wilkes. 1969. Applied Numerical
Methods. John Wiley and Sons, New York, NY.
»
McWhorter, D. B. 1972. "Steady and Unsteady Flow of Fresh Water in
Saline Aquifers." Water Management Technical Report No. 20,
Engineering Research Center, Colorado State University, Fort
Collins, CO.
Schmorak, S. and A. Mercado. 1969. "Upconing of Fresh Water-Sea Water
Interface Below Pumping Wells, Field Study." Water Resources Re-
search, Vol. 5, No. 6, pp 1290-1311.
33
-------
APPENDIX A
Example Problems
The following pages contain the documentation of the upconing simulation
used to generate the predicted interface/time and concentration/time
relationships for the example problem discussed in Section II of this report.
34
-------
ENTER TITLE
?
SEMADAR 1 — TEST B
ENTER UNITS FOR LENGTH (2 CHARACTERS)
?
M
ENTER UNITS FOR TIME <2 CHARACTERS)
?
DY
ENTER FRESH-MATER AND SALT-WATER DENSITIES
?,?
1.00,1.03
ENTER AQUIFER POROSITY
?
0.33
ENTER HORIZONTAL PERMEABILITY (M /DY)
?
14.7
ENTER VERTICAL PERMEABILITY (M /DY)
?
14.7
ENTER INITIAL INTERFACE ELEVATION (M )
7
-30.75
ENTER DISTANCE FROM BOTTOM OF WELL TO INITIAL INTERFACE (M )
7
15.5
ENTER FRACTIONAL CRITICAL RISE
?
0.4
CONCENTRATION CALCULATIONS ? (Y/N)
N
35
-------
SEMADAR 1 — TEST B
DENSITY OP FRESH WATER ' 1.0000
DENSITY OF SALT WATER 1.0300
AQUIFER POROSITY .3300
HORIZONTAL PERMEABILITY 14.7000
VERTICAL PERMEABILITY (M /DY) 14.7000
INITIAL INTERFACE ELEVATION (M ) -30.7500
DISTANCE FROM BOTTOM OF WELL TO INTERFACE (M ) 15.5000
FRACTIONAL CRITICAL RISE . 4OOO
CRITICAL RISE (M ) 6.20OO
CRITICAL ELEVATION (M ) -24.5500
MAXIMUM STEADY-STATE PUMPING RATE (CU M XDY) 266.2818
ENTER NEXT COMMAND
7
EL
CONCENTRATION CALCULATIONS ? (Y/N)
ENTER PUMPING RATE
?,?
34B.,84.
ENTER TFIRST, TLAST, DELTAT
?,?,?
0.,57.,16.
ENTER RFIRST, RLAST, DELTAR (M )
?,?,?
0. ^40. ,5.
PUMPING RATE (CU M /DY)
PUMPING PERIOD (DY)
TFIRST
RFIRST
.0000 TLAST <= 57.0000 DELTAT =
.0000 RLAST - 40.0000 DELTAR •
NOTE! INTERFACE WILL RISE TO CRITICAL ELEVATION IN
CONTINUE ? (Y/N)
348.0000
84.0000
16.0000
5.0000
75.59 DY
36
-------
PUMPING RATE! 348.00 CD M /DY FOR 84.00 DY
INTERFACE ELEVATIONS (M )
R (M >
.00 5.00 10.00 15.00 20.00 25.00 30.00
T (DY)
.00 -30.75 -30.75 -30.75 -30.75 -30.75 -30.75 -30.75
16.00 -27.44 -27.75 -28.42 -29.09 -29.60 -29.95 -30.18
32.00 -26.05 -26.41 -27.23 -28.08 -28.78 -29.30 -29.67
48.00 -25.29 -25.66 -26.52 -27.45 -28.22 -28.82 -29.26
57.00 -24.99 -25.37 -26.25 -27.18 -27.98 -28.60 -29.OB
37
-------
INTERFACE ELEVATIONS (M ) (CONTINUED)
R (M )
35.00 40.OO
T (DY)
.00 -30.75 -30.75
16.00 -30.34 -30.45
32.00 -29.94 -30.13
48. 00 -29.60 -29.84
57.00 -29.43 -29.70
ENTER NEXT COMMAND
OB
ENTER TFIRST, TLAST, DELTAT (DY)
0.,160.,5.
ENTER RFIRST, RLAST, DELTAR (M )
4.5,0.' ,0.
ENTER NEXT COMMAND
EL
CONCENTRATION CALCULATIONS 7 (Y/N)
N
PUMPING RATE (CU M /DY) 348.0000
PUMPING PERIOD (DY) 84.OOOO
TFIRST • .OOOO TLAST = 16O.OOOO DELTAT = 5.0OOO
RFIRST • 4.5000 RLAST - .OOOO DELTAR = .OOOO
NOTEt INTERFACE WILL RISE TO CRITICAL ELEVATION IN 75.59 DY
CONTINUE ? (Y/N)
38
-------
PUMPING RATE: 348.00 CU M /DY FOR
INTERFACE ELEVATIONS (M )
B4.00 DY
R (M )
T (DY)
4.50
.00
5 . 00
1 0 . 00
15.00
20.00
25.00
30 . 00
35.00
40.00
45.00
50 . 00
55.OO
60.00
65.00
70.00
75.00
80.00
85.00
9O.OO
95.00
100.00
105.00
110.00
115.00
120.00
125.00
130.00
135.00
140.00
145.00
150.00
155.00
160.00
-30.75
-29.45
-28.52
-27.81
-27.27
-26.83
-26.47
-26. 18
-25.93
-25.71
-25.53
-25.36
-25.22
-25.09
-24.98
-24.89
-24.79
-25.00
-26. 13
-26.94
-27.55
-28.01
-28.37
-28.67
-28.91
-29. 11
-29.27
-29.41
-29.54
-29.64
-29.73
-29.81
-29.88
* NOTEi CRITICAL ELEVATION OP -24.55 M EXCEEDED AT R=0 AND T=
75.59 DY
ENTER NEXT COMMAND
7
RC
ENTER RFIRST, RLAST, DELTAR
-------
PUMPING RATE: 348.00 CU M /DY FOR 84.00 DY
INTERFACE ELEVATIONS (M )
R (M )
» 12.40
(DY) »
*
.00
5.OO
10.00
15.00
20.00
25.00
30.00
35.00
40.00
45.00
5O.OO
55.00
60.00
65. OO
70.00
75.00
80.00
85. OO
9O.OO
95.00
100. OO
105.00
110.00
115.00
120.00
125.00
130.00
135.00
140.00
145.00
150.00
155. OO
160.OO
-30.75
-29.99
-29.36
-28.85
-28.42
-28.06
-27.76
-27.50
-27.28
-27.08
-26.91
-26.76
-26.63
-26.51
-26.40
-26.30
-26.22
-26.30
-26.96
-27.49
-27.92
-28.28
-28.57
-28.82
-29.02
-29.20
-29.34
-29.47
-29.58
-29.68
-29.77
-29.84
-29.91
* NOTE: CRITICAL ELEVATION OF -24.55 M EXCEEDED AT R=0 AND T= 75.59 DY
ENTER NEXT COMMAND
RC
ENTER RFIRST, RLAST, DELTAR (M )
16.7,0.,0.
ENTER NEXT COMMAND
7
EL
CONCENTRATION CALCULATIONS 7 (Y/N)
N
PUMPING RATE (CU M /DY) 348.0000
PUMPING PERIOD (DY) 84.000O
TFIRST • .0000 TLAST •» 160.0000 DELTAT •> 5.0000
RFIRST - 16.7000 RLAST - .0000 DELTAR •> .0000
NOTE: INTERFACE WILL RISE TO CRITICAL ELEVATION IN 75.59 DY
CONTINUE' ?
Y
40
-------
PUMPING RATE: 34S.OO CD M /DY FOR
INTERFACE ELEVATIONS (M )
84.00 DY
* R (M )
»
T (DY) *
16.70
.00
5.00
1O.OO
15.00
20.00
25.00
30.00
35.00
40.00
45.OO
50.00
55.00
60.00
65. OO
70. OO
75.00
80.00
85.00
90.00
95.00
100.00
105.00
110.00
115.00
120.00
125.00
130.00
135.00
140.00
145.00
150.00
155.00
160.00
-3O.75
-30.23
-29.76
-29.36
-29.00
-28.70
-28.44
-28.21
-28.00
-27.83
-27.67
-27.53
-27.40
-27.29
-27. 19
-27.09
-27.01
-27.04
-27.48
-27.87
-28.20
-28.49
-28.73
-28.94
-29. 12
-29.27
-29.41
-29.52
-29.63
-29.72
-29.80
-29.87
-29.93
* NOTE: CRITICAL ELEVATION OF -24.55 M EXCEEDED AT R=0 AND T- 75.59 DY
ENTER NEXT COMMAND
RC
ENTER RFIRST, RLAST, DELTAR (M )
33.9
ENTER NEXT COMMAND
EL
CONCENTRATION CALCULATIONS ? (Y/N)
PUMPING RATE (CU M /DY)
PUMPING PERIOD
TFIRST • .0000 TLAST = 160.000O DELTAT «=
RFIRST - 33.9000 RLAST • .OOOO DELTAR =
NOTE: INTERFACE WILL RISE TO CRITICAL ELEVATION IN
CONTINUE ? (Y/N)
348.0000
84.0000
5.0000
. 0000
75.59 DY
41
-------
PUMPING RATE I 348.00 CU M /DY FOR 84.00 DY
INTERFACE ELEVATIONS (M )
* R (M )
* 33.90
T (DY) »
.00
5.00
1 0 . 00
15.00
20.00
25 . 00
30.00
35.00
40.00
45.00
50.00
55.00
60.00
65.00
70.00
75.00
80.00
85.00
90.00
95.00
100.00
105.00
110.00
115.00
120.00
125.00
130.0O
135.00
140.00
145.00
150.00
135.00
160.00
-30.75
-30.62
-30.48
-30.34
-30.20
-30.07
-29.94
-29.82
-29.70
-29.59
-29.49
-29.40
-29.31
-29.23
-29. 15
-29.08
-29.02
-28.98
-29.05
-29. 14
-29.23
-29.32
-29.41
-29.49
-29.58
-29.65
-29.72
-29.79
-29.85
-29.91
-29.96
-30.01
-30.05
* NOTEl CRITICAL ELEVATION OF -24.55 M EXCEEDED AT R=0 AND T= 75.59 DY
ENTER NEXT COMMAND
" ? .
OB
ENTER TFIRST, TLAST, DELTAT (DY)
?,?.?
0.,84.,5.
ENTER RFIRST, RLAST, DELTAR (M )
?,?,?
0.
ENTER NEXT COMMAND
?
EL
CONCENTRATION CALCULATIONS 7 (Y/N)
Y
ENTER UNITS FOR CONCENTRATION (6 CHARACTERS)
7
PPM CL
ENTER SALT-WATER AND BACKGROUND CONCENTRATIONS (PPM CD
?,?
22000.,145.
ENTER DIBPERSIVITY (M >
?
0.5
ENTER INITIAL WIDTH OF TRANSITION ZONE (M )
?
3.5
ENTER INTERCEPTION COEFFICIENT
?
0.08
42
-------
SEMADAR 1 — TEST B
DENSITY OF FRESH WATER 1.OOOO
DENSITY OF SALT WATER 1.0300
AQUIFER POROSITY .3300
HORIZONTAL PERMEABILITY 14.7000
VERTICAL PERMEABILITY (M /DY) 14.7000
INITIAL INTERFACE ELEVATION 15.30OO
FRACTIONAL CRITICAL RISE . 40OO
CRITICAL RISE (M ) 6.2000
CRITICAL ELEVATION (M > -24.550O
MAXIMUM STEADY-STATE PUMPING RATE (CU M /DY) 266.2818
CONCENTRATION IN SALT WATER (PPM CD
BACKGROUND CONCENTRATION (PPM CD
INITIAL WIDTH OF TRANSITION ZONE
-------
PUMPING RATEl 348.00 CU M /DY FOR 64.OO DY
INTERFACE ELEVATIONS
-------
T
(DY)
CONCENTRATION IN WELL AND PROFILES BENEATH WELL (PPM CD
ELEVATION FOR C/(E) (M )
C(WELL)
E(WELL)
145.6 2330.5 4516.0 6701.5 8887.0 11072.3
( .0)
.00 145.17
( .0000)
5.00 155.81
( .0003)
10.00 192.67
( .0022)
-26.4
( . 1) < .2)
-2B.5 -29.3
( .3) ( .4) ( .5)
-29.8 -30.3 -30.8
-24.0* -26.6 -27.5 -28.2 -28.8 -29.3
-22.4* -25.3 -26.3 -27.1 -27.7
-28. 3
15.00 244.28 -21.3* -24.4* -25.5 -26.3 -26.9 -27.6
( .0045)
20.00 297.22
( .007O)
25.00 345.51
< .0092)
30.00 387.60
( .0111)
35.00 423.69
( .0128)
40.00 454.52
( .0142)
-20.5* -23.7* -24.8 -25.6 -26.3 -27.0
-19.8* -23.1* -24.3* -25.1 -25.9 -26.5
-19.3* -22.6* -23.9* -24.7 -25.5
-18.8* -22.3* -23.5* -24.4* -25.2
-18.5* -22.0* -23.2* -24.1* -24.9
45.0O 480.92 -18.2* -21.7* -23.0* -23.9* -24.7
< .0154)
-26.2
-25.9
-25.6
-25.4
50.00 503.65
( .0164)
55.00 523.35
( .0173)
-17.9* -21.5* -22.7* -23.7* -24.5* -25.2
-17.7* -21.3* -22.6* -23.5* -24.3* -25.1
60.00 540.53 -17.4* -21.1* -22.4* -23.3* -24.2* -24.9
( .0181)
65.00 555.62
( .0188)
70.00 568.94
( .0194)
75.00 580.79
( .0199)
-17.3* -20.9* -22.2* -23.2* -24.0* -24.8
-17.1* -20.8* -22.1* -23.1* -23.9* -24.7
-17.0* -20.7* -22.0* -23.0* -23.8* -24.6
80.00 591.38 -16.8* -20.6* -21.9* -22.9* -23.7* -24.5*
< .0204)
84.00 599.06
< .0208)
-16.7* -20.5* -21.8* -22.8* -23.6* -24.4*
* NOTEl THE DISPERSION CONCEPT SHOULD BE .LIMITED TO THE ZONE BELOW
THE CRITICAL ELEVATION OF -24.5500 M
45
-------
CONCENTRATION IN WELL AND PROFILES BENEATH WELL (PPM CD
T
(DY)
.00
5.00
10.00
15.00
20.00
25.00
30.00
35.00
40.00
45.00
50.00
55.OO
6O.OO
65. OO
70.00
75.00
80.00
84.00
C(WELL)
E(WELL)
145. 17
( .0000)
155.81
( .0005)
192.67
( .0022)
244.28
( .0045)
297.22
< .0070)
345.51
( .0092)
3B7.60
( .0111)
423.69
< .0128)
454.52
( .0142)
48O.92
< .0154)
503.65
< .0164)
523.35
( .0173)
540.53
( .0181)
555.62
( .0188)
568.94
( .0194)
580.79
< .0199)
591.38
< . 0204 )
599.06
< . O2OB)
ELEVATION FOR C/(E> (M
11072.5 13258.0 15443.5 17629.0
('.5) < .6) ( .7) ( .8)
-30.8
-29.3
-28.3
-27.6
-27.0
-26.5
-26.2
-25.9
-25.6
-25.4
-25.2
-25. 1
-24.9
-24. B
-24.7
-24.6
-24. 3«
-24. 4»
-31.2
-29.9
-28.9
-28.2
-27.7
-27.2
-26.9
-26.6
-26.3
-26. 1
-26.0
-25.8
-25.7
-25.5
-25.4
-25.3
-25.2
-25.2
-31.7
-30.4
-29.5
-28.9
-28.4
-28.0
-27.6
-27.4
-27.1
-26.9
-26.8
-26.6
'-26.5
-26.4
-26.3
-26.2
-26. 1
-26.0
-32.2
-31. 1
-30.3
-29.7
-29.2
-28.8
-28.5
-28.2
-28.0
-27.8
-27.7
-27.5
-27.4
-27.3
-27.2
-27.1
-27.0
-27.0
)
19814.5
( .9)
-33.0
-32.0
-31.3
-30.8
-30.3
-30.0
-29.7
-29.5
-29.3
-29. 1
-29.0
-28.8
-28.7
-28.6
-28.5
-28.5
-28.4
-28.3
22000.0
(1.0
-35. 1
-34.6
-34.2
-33.8
-33.5
-33.3
-33. 1
-32.9
-32.8
-32.7
-32.5
-32.5
-32.4
-32.3
-32.2
-32.2
-32. 1
-32.1
« NOTE« THE DISPERSION CONCEPT SHOULD BE LIMITED TO THE ZONE BELOW
THE CRITICAL ELEVATION OF -24.55OO M
ENTER NEXT COMMAND
?
NP
46
-------
ENTER TITLE
7
SEMADAR 1 SALINITY/TIME RELATIONSHIPS
ENTER UNITS FOR LENGTH (2 CHARACTERS)
7
M
ENTER UNITS FOR TIME (2 CHARACTERS)
7
DY
ENTER FRESH-WATER AND SALT-WATER DENSITIES
?,?
1.00,1.03
ENTER AQUIFER POROSITY
7
0.33
ENTER HORIZONTAL PERMEABILITY (M /DY)
7
14.7
ENTER VERTICAL PERMEABILITY (M /DY)
7
14.7
ENTER INITIAL INTERFACE ELEVATION (M )
7
-30.75
ENTER DISTANCE FROM BOTTOM OF WELL TO INITIAL INTERFACE
7
15.5
ENTER FRACTIONAL CRITICAL RISE
7
0.4
CONCENTRATION CALCULATIONS ? (Y/N)
Y
ENTER UNITS FOR CONCENTRATION (6 CHARACTERS)
?
PPM CL
ENTER SALT-WATER AND BACKGROUND CONCENTRATIONS (PPM CL)
?,?
22000.,145.
ENTER DISPERSIVITY (M >
?
0.5
ENTER INITIAL WIDTH OF TRANSITION ZONE (M )
7
3.5
ENTER INTERCEPTION COEFFICIENT
7
0.08
47
-------
SEMADAR 1 SALINITY/TIME RELATIONSHIPS
DENSITY OF FRESH WATER 1.0000
DENSITY OF SALT WATER 1.0300
AQUIFER POROSITY .3300
HORIZONTAL PERMEABILITY (M /DY) 14.7000
VERTICAL PERMEABILITY (M /DY) 14.7000
INITIAL INTERFACE ELEVATION (M > -30.7500
DISTANCE FROM BOTTOM OF WELL TO INTERFACE (M ) 15.5000
FRACTIONAL CRITICAL RISE . 40OO
CRITICAL RISE (M ) 6.2000
CRITICAL ELEVATION (M ) -24.550O
MAXIMUM STEADY-STATE PUMPING RATE (CU M /DY) 266.2818
CONCENTRATION IN SALT WATER (PPM CD
BACKGROUND CONCENTRATION (PPM CD
INITIAL WIDTH OF TRANSITION ZONE (M )
DISPERSIVITY (M )
INTERCEPTION COEFFICIENT
22000.0000
145.0000
3.5000
. 5000
.0800
ENTER NEXT COMMAND
PR
ENTER MAXIMUM PERMISSIBLE CONCENTRATION IN PUMPED WATER (PPM CD
166.85
48
-------
SEMADAR 1 SALINITY/TIME RELATIONSHIPS
DENSITY OF FRESH WATER 1.0000
DENSITY OF SALT WATER 1.0300
AQUIFER POROSITY .3300
HORIZONTAL PERMEABILITY (M /DY) 14.7000
VERTICAL PERMEABILITY (M /DY) 14.7000
INITIAL INTERFACE ELEVATION
575.0000
3.7256
ENTER OPTIONAL PUMPING RATE (CU M /DY)
0.
ENTER NEXT COMMAND
PR
ENTER MAXIMUM PERMISSIBLE CONCENTRATION IN PUMPED WATER (PPM CD
7
210.56
-------
SEMADAR 1 SALINITY/TIME RELATIONSHIPS
DENSITY OF FRESH WATER 1.0000
DENSITY OF SALT WATER 1.0300
AQUIFER POROSITY .3700
HORIZONTAL PERMEABILITY 15.5OOO
FRACTIONAL CRITICAL RISE .4000
CRITICAL RISE (M ) 6.2000
CRITICAL ELEVATION (M ) -24.550O
CONCENTRATION IN SALT WATER (PPM CD 22000.0000
BACKGROUND CONCENTRATION (PPM CD 145.0000
INITIAL WIDTH OF TRANSITION ZONE (M ) 3.50OO
DISPERSIVITY (M ) .5000
INTERCEPTION COEFFICIENT .OBOO
MAXIMUM CONCENTRATION IN PUMPED WATER (PPM CD 210.5600
MAXIMUM RELATIVE CONCENTRATION .0030
MAXIMUM INTERFACE ELEVATION (M ) -2B.0165
MAXIMUM PERMISSIBLE PUMPING RATE (CU M /DY) 117.4OO9
ENTER OPTIONAL PUMPING RATE (CU M /DY)
7
348.
PUMPING RATE (CU M /DY)
TIME TO REACH CMAX (DY)
348.0000
11.810O
ENTER OPTIONAL PUMPING RATE
-------
SEMADAR 1 SALINITY/TIME RELATIONSHIPS
DENSITY OF FRESH WATER . 1.0000
DENSITY OF SALT WATER 1.030O
AQUIFER POROSITY .3300
HORIZONTAL PERMEABILITY (M /DY) 14.700O
VERTICAL PERMEABILITY (M /DY> 14.7000
INITIAL INTERFACE ELEVATION (M ) -30.7500
DISTANCE FROM BOTTOM OF WELL TO INTERFACE (M ) 15.5OOO
FRACTIONAL CRITICAL RISE .4000
CRITICAL RISE (M ) 6.2000
CRITICAL ELEVATION
-------
SEMADAR 1 SALINITY/TIME. RELATIONSHIPS
DENSITY OF FRESH WATER 1.0000
DENSITY OF SALT WATER 1.0300
AQUIFER POROSITY . .3300
HORIZONTAL PERMEABILITY 14.7000
INITIAL INTERFACE ELEVATION (M ) -30.750O
DISTANCE FROM BOTTOM OF WELL TO INTERFACE (M ) 15.5000
FRACTIONAL CRITICAL RISE .4000
CRITICAL RISE (M ) 6.2000
CRITICAL ELEVATION (M ) -24.5500
CONCENTRATION IN SALT WATER (PPM CD 22000.0000
BACKGROUND CONCENTRATION (PPM CD 145.000O
INITIAL WIDTH OF TRANSITION ZONE (M ) 3.500O
DISPERSIVITY (M > .5000
INTERCEPTION COEFFICIENT .080O
MAXIMUM CONCENTRATION IN PUMPED WATER (PPM CD 363.5500
MAXIMUM RELATIVE CONCENTRATION .010O
MAXIMUM INTERFACE ELEVATION (M ) -26.3880
MAXIMUM PERMISSIBLE PUMPING RATE (CU M /DY) 187.3442
ENTER OPTIONAL PUMPING RATE (CU M /DY)
7
348.
PUMPING RATE (CU M /DY) 348.0000
TIME TO REACH CMAX (DY) 27.0509
ENTER OPTIONAL PUMPING RATE (CU M /DY)
7
573.
PUMPING RATE (CU M /DY) 575.00OO
TIME TO REACH CMAX (DY) 11.21O7
ENTER OPTIONAL PUMPING RATE (CU M /DY>
•7
0.
ENTER NEXT COMMAND
7
PR
ENTER MAXIMUM PERMISSIBLE CONCENTRATION IN PUMPED WATER (PPM CD
?
582. 1
52
-------
SEMADAR 1 SALINITY/TIME RELATIONSHIPS
DENSITY OF FRESH WATER 1.OOOO
DENSITY OF SALT WATER 1.0300
AQUIFER POROSITY .330O
HORIZONTAL PERMEABILITY (M /DY> 14.700O
VERTICAL PERMEABILITY (M /DY) 14.7000
INITIAL INTERFACE ELEVATION (M ) -30.7500
DISTANCE FROM BOTTOM OF WELL TO INTERFACE (M ) 15.5000
FRACTIONAL CRITICAL RISE .4000
CRITICAL RISE (M ) 6.2000
CRITICAL ELEVATION (M ) -24.5500
CONCENTRATION IN SALT WATER (PPM CD 22000.0000
BACKGROUND CONCENTRATION (PPM CD 145.0000
INITIAL WIDTH OF TRANSITION ZONE (M ) 3.5000
DISPERSIVITY (M ) .5000
INTERCEPTION COEFFICIENT . OBOO
MAXIMUM CONCENTRATION IN PUMPED WATER (PPM CD 582.1000
MAXIMUM RELATIVE CONCENTRATION .0200
MAXIMUM INTERFACE ELEVATION (M ) -24.5511
MAXIMUM PERMISSIBLE PUMPINS RATE (CU M /DY) 266.2366
ENTER OPTIONAL PUMPING RATE (CU M /DY)
348.
PUMPING RATE (CU M /DY) 348.OOOO
TIME TO REACH CMAX (DY) 75.5346
ENTER OPTIONAL PUMPING RATE (CU M /DY)
575.
PUMPING RATE (CU M /DY) 575.OOOO
TIME TO REACH CMAX (DY) 20.0023
ENTER OPTIONAL PUMPING RATE (CU M /DY)
•?
0.
ENTER NEXT COMMAND
7
PR
ENTER MAXIMUM PERMISSIBLE CONCENTRATION IN PUMPED WATER (PPM CD
?
800.65
53
-------
SEMADAR 1 SALINITY/TIME RELATIONSHIPS
DENSITY OF FRESH WATER 1.0000
DENSITY OF SALT WATER 1.0300
AOUIFER POROSITY .3300
HORIZONTAL PERMEABILITY
-------
SEhADAR 1 SALINITY/TIME RELATIONSHIPS
DENSITY OF FRESH WATER 1.0000
DENSITY OF SALT WATER 1.03OO
AQUIFER POROSITY .3300
HORIZONTAL PERMEABILITY (M /DY) 14.7OOO
VERTICAL PERMEABILITY (M /DY) 14.7000
INITIAL INTERFACE ELEVATION (M ) -30.75OO
DISTANCE FROM BOTTOM OF WELL TO INTERFACE (M ) 15.5000
FRACTIONAL CRITICAL RISE .4000
CRITICAL RISE (M ) 6.2000
CRITICAL ELEVATION -24.5500
CONCENTRATION IN SALT WATER (PPM CD 22000.0000
BACKGROUND CONCENTRATION (PPM CD 145.0000
INITIAL WIDTH OF TRANSITION ZONE (M ) 3.5000
DISPERSIVITY (M > .5000
INTERCEPTION COEFFICIENT .0800
MAXIMUM CONCENTRATION IN PUMPED WATER (PPM CD 1019.2000
MAXIMUM RELATIVE CONCENTRATION .0400
MAXIMUM INTERFACE ELEVATION (M ) -13.1997»
MAXIMUM PERMISSIBLE PUMPING RATE (CU M /DY) 753.7639
• NOTE! THE DISPERSION CONCEPT SHOULD BE LIMITED TO THE ZONE BELOW
THE CRITICAL ELEVATION OF -24.55OO M
(MAXIMUM CONCENTRATIONS IN PUMPED WATER LESS THAN 582.1O PPM CD
ENTER OPTIONAL PUMPING RATE (CU M /DY)
•7
0.
ENTER NEXT COMMAND
?
DN
Stop - Program terminated.
C>
55
-------
APPENDIX B
Description of Program UPCONE
Program UPCONE has been written in an unextended Fortran computer code in
an effort to make the program transportable between computer systems. The
code consists of a main program and two function subroutines. The program has
been documented internally through the liberal use of comment statements.
The main program is divided into three sections. Section I provides for
the "Basic Input Data" as described in Section III of this report. The
numerical evaluation of interface elevations, concentrations, and maximum
pumping rates is accomplished in Section II of the main program which contains
the computational algorithms for Case I and Case II types of problems.
Section III provides for problem redefinition and control of execution under
the "Edit" mode.
Two function subroutines are required to calculate the concentration
distribution across the transition zone or to evaluate the maximum pumping
rate for a specified salinity of the pumped water. Function ERFC (Z) is a
rational approximation of the complimentary error function of the argument Z.
Function IERFC (X, Z, ZMIN, ZMAX) uses a fieguLa. ^aLbi. root finding
technique to find the argument, Z, for a specified value of the complimentary
error function, X. The last two arguments, ZMIN and ZMAX, define the lower
and upper limits of the initial search interval.
56
-------
APPENDIX C
Listing of Program UPCONE
57
-------
PROGRAM UPCONE
VERSION 1.0
JAN WAGNER
SCHOOL OF CHEMICAL ENGINEERING
OKLAHOMA STATE UNIVERSITY
STILLWATER, OK 74078
TELEPHONE (405) 624-5280
JANUARY, 1982
DIMENSION A(30),C(11),CW(51),E(11),EW(51),1C(20).KFLG(11).
1R(26),T(51 ),2(26,26)
REAL*4 KX.KZ
DATA KHAR1.KHAR2.NY/' ','*', 'Y'/
DATA IC/'FD' , 'PO' , 'KX', 'KZ' , 'ZO', 'DT' , 'CR' , 'CO', '01' , 'WT', '1C'
1'OB', 'TC', 'RC','OP','LI','EL','PR','NP','ON'/
READ DEVICE: NI=1
NI = 1
NO = 1
WRITE DEVICE: N0=1
MAXIMUM NUMBERS OF OBSERVATION POINTS FOR TIME AND RADIUS
HAVE BEEN SET AT 50 AND 25, RESPECTIVELY
DIMENSION CW(M),EW(M),R(N),T(N),Z(M,N)
— WHERE M=MAXTIM+1 AND N=MAXPTS+1
MAXTIM = 50
MAXPTS = 25
INITIALIZE PROGRAM FLOW PARAMETERS
1 IEDIT = 1
KCON = 1
KELE = 1
*..** SECTION I -- BASIC INPUT DATA
READ TITLE
WRITE(NO,5)
5 FORMAT('1',3X,'ENTER TITLE',/,3X,'?')
READ(NI,1O) (A(I),1=1,3O)
10 FORMAT(30A2)
DEFINE UNITS
WRITE(NO,15)
15 FORMAT(3X,'ENTER UNITS FOR LENGTH (2 CHARACTERS)'./,3X,'?')
READ(NI,2O) IL
20 FORMAT(A2)
WRITE(NO,25)
25 FORMAT(3X,'ENTER UNITS FOR TIME (2 CHARACTERS)',/,3X,'?')
READ(NI,20) ITU
FLUID DENSITIES
29 WRITE(NO,30)
30 FORMAT(3X,'ENTER FRESH-WATER AND SALT-WATER DENSITIES',/,3X,'?,?
REAO(NI,35) RHOF.RHOS
35 FORMAT(2F10.0)
40 IF(RHOS.GT.O.O) GO TO 55
WRITE(NO,45)
45 FORMAT(3X.'SALT-WATER DENSITY MUST BE GREATER THAN ZERO',
1' -- REENTER',/,3X.'?')
READ(NI.SO) RHOS
50 FORMAT(FIO.O)
GO TO 40
55 IF(RHOF.GT.0.0.AND.RHOF.LT.RHOS) GO TO 69
WRITE(NO.SO)
60 FORMAT(3X,'FRESH-WATER DENSITY MUST BE GREATER THAN ZERO',
1/.6X,'AND LESS THAN SALT-WATER DENSITY -- REENTER'./,3X,'?')
READ(NI,50) RHOF
GO TO 55
69 GO TO (70,700), IEDIT
UPCN001
UPCN002
UPCN003
UPCN004
UPCN005
UPCN006
UPCN007
UPCN008
UPCN009
UPCN010
UPCN011
UPCN012
UPCN013
UPCN014
UPCN015
UPCN016
UPCN017
UPCN018
UPCN019
UPCN020
UPCN021
UPCN022
UPCN023
UPCN024
UPCN025
UPCN026
UPCN027
UPCN028
UPCN029
UPCN030
UPCN031
UPCN032
UPCN033
UPCN034
UPCN035
UPCN036
UPCN037
UPCN038
UPCN039
UPCN040
UPCN041
UPCN042
UPCN043
UPCNO44
UPCN045
UPCN046
UPCN047
UPCN048
UPCN049
UPCN050
UPCN051
UPCN052
) UPCN053
UPCN054
UPCN055
UPCN056
UPCN057
UPCN058
UPCN059
UPCN060
UPCN061
UPCN062
UPCN063
UPCN064
UPCN065
UPCN066
UPCN067
UPCN068
UPCN069
UPCN070
58
-------
POROSITY
70 WRITE(NO,75)
75 FORMAT(3X,'ENTER AQUIFER POROSITY',/.3X,'?')
80 READ(NI,50) PO
IF(PO.GT.O.O.AND.PO.LT.1.0) GO TO 89
WRITE(NO,85)
85 FORMAT(3X,'POROSITY MUST BE GREATER THA.N ZERO AND LESS THAN',
1' ONE -- REENTER',/,3X,'?')
GO TO 80
89 GO TO (90,700), IEDIT
HORIZONTAL PERMEABILITY
90 WRITE(NO,95) IL.ITU
95 FORMAT(3X, 'ENTER HORIZONTAL PERMEABILITY ( ' ,A2, ' /',A2, ' ) ',
1/.3X,'?')
100 READ(NI,50) KX
IF(KX.GT.O.O) GO TO 110
WRITE(NO,105)
105 FORMAT(3X,'HORIZONTAL. PERMEABILITY MUST BE GREATER THAN',
1' ZERO -- REENTER',/,3X,'?')
GO TO 1OO
110 GO TO (111.700), IEDIT
VERTICAL PERMEABILITY
111 WRITE(NO,112) IL.ITU
112 FORMAT(3X,'ENTER VERTICAL PERMEABILITY (',1A2.'/',A2.') ',
1/.3X,'?') __
114 READ(NI,50) KZ
IF(KZ.GT.O.O) GO TO 119
WRITE(NO,115)
115 FORMAT(3X,'VERTICAL PERMEABILITY MUST BE GREATER THAN ZERO',
1' — REENTER',/,3X,'?')
GO TO 114
119 GO TO (120.70O), IEDIT
INITIAL INTERFACE ELEVATION
120 WRITE(NO,125) IL
125 FORMAT(3X,'ENTER INITIAL INTERFACE ELEVATION (',A2,')'./,3X,'?')
READ(NI.SO) ZO
GO TO (129,700), IEDIT
DISTANCE FROM BOTTOM OF WELL TO INITIAL INTERFACE
129 WRITE(NO,130) IL
130 FORMAT(3X,'ENTER DISTANCE FROM BOTTOM OF WELL TO INITIAL ',
1'INTERFACE ('.A2,') './.3X,'?')
135 READ(NI,50) D
IF(D.GT.O.O) GO TO 144
WRITE(NO,140)
140 FORMAT(3X.'DISTANCE MUST BE GREATER THAN ZERO -- REENTER'./,
13X,'?' )
GO TO 135
144 GO TO (145,700), IEDIT
FRACTIONAL CRITICAL RISE
145 WRITE(NO,150)
150 FORMAT(3X,'ENTER FRACTIONAL CRITICAL RISE'./ , 3X,'?')
155 READ(NI,50) THETA
IF(THETA.GT.O.O.AND.THETA.LT.1.0) GO TO 164
WRITE(NO,160)
160 FORMAT(3X,'FRACTION MUST BE GREATER THAN ZERO AND LESS THAN',
1' ONE — REENTER',/,3X,'?')
GO TO 155
164 XCR = THETA*D
ZCR = XCR + ZO
MAXIMUM STEADY-STATE PUMPING RATE
OMAXSS = 6.283185*((RHOS-RHOF)/RHOF)*KX*D*XCR
165 GO TO (169,700). IEDIT
DATA FOR CONCENTRATION CALCULATIONS
UPCN071
UPCN072
UPCN073
UPCN074
UPCN075
UPCN076
UPCN077
UPCN078
UPCN079
UPCN080
UPCN081
UPCN082
UPCN083
UPCN084
UPCN085
UPCN086
UPCN087
UPCN088
UPCNO89
UPCN090
UPCN091
UPCN092
UPCN093
UPCN094
UPCN095
UPCN096
UPCN097
UPCN098
UPCN099
UPCN100
UPCN101
UPCN102
UPCN103
UPCN104
UPCN105
UPCN106
UPCN107
UPCN108
UPCN109
UPCN110
UPCN111
UPCN112
UPCN113
UPCN114
UPCN115
UPCN116
UPCN117
UPCN118
UPCN119
UPCN120
UPCN121
UPCN122
UPCN123
UPCN124
UPCN125
UPCN126
UPCN127
UPCN128
UPCN129
UPCN130
UPCN131
UPCN132
UPCN133
UPCN134
UPCN135
UPCN136
UPCN137
UPCN138
UPCN139
UPCN140
59
-------
169 WRITE(NO,170)
170 FORMAT('0',2X,'CONCENTRATION CALCULATIONS ? (Y/N)')
READ(NI, 175) ICON •
175 FORMAT(A1)
IF(ICON.NE.NY) GO TO 275
176 KCON = 2
WRITE(NO,180)
180 FORMATOX,'ENTER UNITS FOR CONCENTRATION (6 CHARACTERS)',
1/.3X,'?')
READ(NI,185) IM1.IM2.IM3
185 FORMATOA2)
SALT-WATER AND BACKGROUND CONCENTRATIONS
188 IF(KCON.EQ.2) GO TO 189
WRITE(NO.180)
READ(NI.185) IM1.IM2.IM3
189 WRITE(NO,190) IM1.IM2.IM3
190 FORMATOX.'ENTER SALT-WATER AND BACKGROUND CONCENTRATIONS (',
13A2,') ',/,3X,'?,?')
READ(NI,35) CO.CB
195 IF(CO.GE.1.0) GO TO 205
WRITE(N0.200)
200 FORMATOX,'SALT-WATER CONCENTRATION MUST BE GREATER THAN',
1' OR EQUAL TO ONE -- REENTER',/,3X,'?')
READ(NI,50) CO __
GO TO 195
205 IF(CB.GE.O.O.AND.CB.LT.CO) GO TO 214
WRITE(NO,210)
210 FORMATOX,'BACKGROUND CONCENTRATION MUST BE GREATER THAN',
1' OR EQUAL TO ZERO',/,6X,'AND LESS THAN SALT-WATER',
2' CONCENTRATION -- REENTER',/,3X,'?')
READ(NI.SO) CB
GO TO 205
214 GO TO (215,700,215,215). IEDIT
DISPERSIVITY
215 WRITE(NO,220) IL
220 FORMATOX,'ENTER DISPERSIVITY ( ' , A2 , ' ) ',/,3X,'?')
225 READ(NI.SO) DM
IF(DM.GT.O.O) GO TO 234
WRITE(NO,230)
230 FORMATOX,'DISPERSIVITY MUST BE GREATER THAN ZERO -- REENTER',
1/.3X,'?')
GO TO 225
234 GO TO (235.700,235,235), IEDIT
INITIAL WIDTH OF TRANSITION ZONE
235 WRITE(NO,240) IL
240 FORMATOX,'ENTER INITIAL WIDTH OF TRANSITION ZONE C.A2,') ',
1/.3X,'?')
245 REAO(NI.SO) S02
IF(S02.GE.O.O) GO TO 254
WRITE(NO,250)
250 FORMATOX,'INITIAL WIDTH MUST BE GREATER THAN ZERO -- REENTER',
1/.3X,'?')
GO TO 245
254 GO TO (255,700,255,255), IEDIT
INTERCEPTION COEFFICIENT
255 WRITE(NO,260)
260 FORMATOX,'ENTER INTERCEPTION COEFFICIENT',/. 3X,'?'')
265 READ(NI,50) PHI
IF(PHI.GT.0.0.AND.PHI.LT.1.0) GO TO 274
WRITE(NO,270)
270 FORMATOX,'COEFFICIENT MUST BE GREATER THAN ZERO AND LESS THAN'
1' ONE -- REENTER',/,3X,'?')
GO TO 265
274 GO TO (275,700,275,510), IEDIT
UPCN141
UPCN142
UPCN143
UPCN144
UPCN145
UPCN146
UPCN147
UPCN148
UPCN149
UPCN150
UPCN151
UPCN152
UPCN153
UPCN154
UPCN155
UPCN156
UPCN157
UPCN158
UPCN159
UPCN160
UPCN161
UPCN162
UPCN163
UPCN164
UPCN165
UPCN166
UPCN167
UPCN168
UPCN169
UPCN170
UPCN171
UPCN172
UPCN173
UPCN174
UPCN175
UPCN176
UPCN177
UPCN178
UPCN179
UPCN180
UPCN181
UPCN182
UPCN183
UPCN184
UPCN185
UPCN186
UPCN187
UPCN188
UPCN189
UPCN190
UPCN191
UPCN192
UPCN193
UPCN194
UPCN195
UPCN196
UPCN197
UPCN198
UPCN199
UPCN200
UPCN201
UPCN202
UPCN203
UPCN204
UPCN205
UPCN206
UPCN207
UPCN208
UPCN209
UPCN210
60
-------
LIST PROBLEM DEFINITION
275 CONTINUE
WRITE(NO,280) (A(I),1 = 1,30), RHOF,RHOS,PO,IL,ITU,KX,IL,ITU , KZ
280 FORMAT('1',3X,30A2,//,
16X,'DENSITY OF FRESH WATER',25X,F10.4,/,
'DENSITY OF SALT WATER',26X,F10.4,//,
'AQUIFER POROSITY'.31X.F10.4,/,
'HORIZONTAL PERMEABILITY ('.A2.'/'.A2,') ',15X,F10.4,/,
56X,'VERTICAL PERMEABILITY (',A2.'/',A2.') ',17X,F10.4)
WRITE(NO,281)IL, ZO,IL,D,THETA,IL,XCR,IL,ZCR
281 FORMAT('0',
15X,'INITIAL INTERFACE ELEVATION (',A2,') ',14X,F10.4,/
26X,'DISTANCE FROM BOTTOM OF WELL TO INTERFACE (',A2,')
'FRACTIONAL CRITICAL RISE ',22X,F10.4,/,
'CRITICAL RISE C.A2,') ',28X,F10.4,/,
26X,
36X,
46X,
36X,
46X,
F10.4,/.
56X,'CRITICAL ELEVATION (',A2,') '.23X.F1O.4)
WRITE(NO,285) IL.ITU,QMAXSS
285 FORMAT('0',5X,'MAXIMUM STEADY-STATE PUMPING RATE (CU ',A2,'/',
1A2,') '.F10.4,/.'0'./,'0')
IF(ICON.NE.NY.OR.KCON.EQ.1) GO TO 700
WRITE(NO,290) IM1,IM2,IMS,CO,IM1,IM2,IM3,CB,IL,S02,IL,DM.PHI
290 FORMAT('0',5X,'CONCENTRATION IN SALT WATER C.3A2,') ',10X,F10.4./
16X,'BACKGROUND CONCENTRATION (',3A2,') ',13X,F10.4,/,
26X, 'INITIAL WIDTH OF TRANSITION ZONE ('A2,' ) ',9X,F10.4,//,
36X,'DISPERSIVITY C,A2,') ',29X,F10.4,/,
46X,'INTERCEPTION COEFFICIENT',23X,F10.4,/,'0',/.'0')
GO TO (700,700,301), IEDIT
*,„„* SECTION II -- NUMERICAL EVALUATION OF INTERFACE ELEVATIONS
CASE I PROBLEMS -- EVALUATE INTERFACE ELEVATIONS AND CONCENTRATION
30O CONTINUE
IEDIT = 3
PARAMETERS FOR INTERFACE ELEVATION CALCULATIONS
301 IF(KELE.EQ.2) GO TO 329
KELE = 2
PUMPING RATE AND PERIOD
302 WRITE(NO,305) IL,ITU,ITU
305 FORMAT('0',2X, 'ENTER PUMPING RATE (CU ' ,A2 , ' / ' .A2, ') AND'.
1' PERIOD (',A2,')',/,3X,'?,?')
READ(NI,35) Q.TPUMP
306 IF(Q.GT.O.O) GO TO 308
WRITE(NO,307)
307 FORMAT(3X,'PUMPING RATE MUST BE GREATER THAN ZERO -- REENTER',
1/.3X, '?')
READ(NI,50) 0
GO TO 306
308 IF(TPUMP.GT.O.O) GO TO 310
WRITE(NO,309)
309 FORMAT(3X,'PUMPING PERIOD MUST BE GREATER THAN ZERO -- REENTER'
1/.3X.'?')
READ(NI,50) TPUMP
GO TO 308
310 GO TO (700,700,312), IEDIT
COORDINATES OF OBSERVATION POINTS -- TIME AND RADIUS
311 IEDIT = 1
312 WRITE(N0.313) ITU
UPCN211
UPCN212
UPCN213
UPCN214
UPCN215
UPCN216
UPCN217
UPCN218
UPCN219
UPCN220
UPCN221
UPCN222
UPCN223
UPCN224
UPCN225
UPCN226
UPCN227
UPCN228
UPCN229
UPCN230
UPCN231
UPCN232
UPCN233
UPCN234
UPCN235
UPCN236
UPCN237
UPCN238
UPCN239
UPCN240
UPCN241
UPCN242
UPCN243
UPCN244
UPCN245
UPCN246
UPCN247
UPCN248
UPCN249
UPCN250
UPCN251
UPCN252
UPCN253
UPCN254
UPCN255
UPCN256
UPCN257
UPCN258
UPCN259
UPCN260
UPCN261
UPCN262
UPCN263
UPCN264
UPCN265
UPCN266
UPCN267
UPCN268
UPCN269
UPCN270
UPCN271
UPCN272
UPCN273
UPCN274
UPCN275
UPCN276
UPCN277
UPCN278
UPCN279
UPCN280
.61
-------
313 FORMAT(3X,'ENTER TFIRST, TLAST, DELTAT (',A2,') ',/.3X,'?,?,?')
READ(NI,315) TF.TL.DELT
315 FORMAT(3F10.0)
DELT = ABS(DELT)
316 IF(TF.GE.0.0.AND.DELT.LE.1.OE-06) GO TO 320
IF(TF.GE.O.O) GO TO 318
WRITE(NO,317)
317 FORMAT(3X,'TFIRST MUST NOT BE LESS THAN ZERO -- REENTER',
1/.3X,'?')
READ(NI,50) TF
GO TO 316
318 IF(TL.GE.O.O) GO TO 321
WRITE(NO,319)
319 FORMAT(3X,'TLAST MUST NOT BE LESS THAN ZERO -- REENTER',
1/.3X,'?')
READ(NI,50) TL
GO TO 318
320 TL = TF
321 GO TO (322,700.322), IEDIT
C
322 WRITE(NO,323) IL
323 FORMAT(3X,'ENTER RFIRST. RLAST, DELTAR (',A2,') ',/,3X,'?,?,?')
READ(NI.325) RF.RL.DELR
325 FORMAT(3F10.0)
DELR = ABS(DELR)
GO TO (700,700,329), IEDIT
C
329 WRITE(NO,330) IL.ITU,0,ITU,TPUMP,TF,TL,DELT,RF,RL,DELR
330 FORMAT('0',5X,'PUMPING RATE (CU ',A2,'/',A2,') ',23X,F10.4,/,
16X,'PUMPING PERIOD C.A2.') ',27X.F10.4,//,
26X,'TFIRST =',F10.4.3X.'TLAST =',F10.4,3X,'DELTAT =',F10.4,/.
36X,'RFIRST ='.F1O.4.3X,'RLAST =',F1O.4,3X,'DELTAR =',F10.4)
IF(O-LE.QMAXSS) GO TO 332
TCR » ((2.«PO*D)/(((RHOS-RHOF)/RHOF)*KZ))»((1./(1.-QMAXSS/Q))-1.)
WRITE(NO,331) TCR,ITU
331 FORMAT('0',2X,'NOTE: INTERFACE WILL RISE TO CRITICAL ELEVATION',
1' IN'.F10.2.A3)
332 WRITE(N0.333)
333 FORMATC'O',2X.'CONTINUE ? (Y/N)')
READ(NI,175) JFLOW
IF(JFLOW.NE.NY) GO TO 700
C
C RADIUS COORDINATES
335 CONTINUE
IR = 1
R(IR) = RF
IF(DELR.LE.1.OE-06) GO TO 345
DIF = RL - RF
IF(ABS(DIF).LE.1.OE-06) GO TO 345
IF(DIF.LE.O.O) DELR = -DELR
NPTS ' DIF/DELR
REM = DIF - DELR-FLOAT(NPTS)
TOL = 1.OE-06*ABS(DIF)
NPTS = NPTS + 1
IF(NPTS.LE.MAXPTS) GO TO 337
WRITE(NO,336) NPTS.MAXPTS
336 FORMAT(3X,I3,' RADIUS OBSERVATION POINTS EXCEED MAXIMUM Of'.14)
GO TO 700
337 CONTINUE
DO 340 IR=2,NPTS
R(IR) = R(IR-1) * DELR
34O CONTINUE
IR » NPTS
IF(A8S(REM).LT.TOL) GO TO 345
IR = IR + 1
R(IR) - RL
345 CONTINUE
TIME COORDINATES
IT = 1
UPCN281
UPCN282
UPCN283
UPCN284
UPCN285
UPCN286
UPCN287
UPCN288
UPCN289
UPCN290
UPCN291
UPCN292
UPCN293
UPCN294
UPCN295
UPCN296
UPCN297
UPCN298
UPCN299
UPCN300
UPCN301
UPCN302
UPCN303
UPCN304
UPCN305
UPCN306
UPCN307
UPCN308
UPCN309
UPCN310
UPCN311
UPCN312
UPCN313
UPCN314
UPCN315
UPCN316
UPCN317
UPCN318
UPCN319
UPCN320
UPCN321
UPCN322
UPCN323
UPCN324
UPCN325
UPCN326
UPCN327
UPCN328
UPCN329
UPCN330
UPCN331
UPCN332
UPCN333
UPCN334
UPCN335
UPCN336
UPCN337
UPCN338
UPCN339
UPCN340
UPCN341
UPCN342
UPCN343
UPCN344
UPCN345
UPCN346
UPCN347
UPCN348
UPCN349
UPCN350
62
-------
T(IT) = TF
IF(DELT.LE.1.OE-06) GO TO 355
DIP = TL - TF
IF(ABS(DIF).LE.1.OE-06) GO TO 355
IF(DIF.LE.O.O) DELT = -DELT
NPTS = D1F/DELT
REM = DIF - DELT*FLOAT(NPTS)
TOL = 1.OE-06«ABS(DIF)
NPTS = NPTS + 1
IF(NPTS.LE.MAXTIM) GO TO 347
WRITE(NO,346) NPTS.MAXTIM
346 FORMAT(3X,I3,' TIME OBSERVATION POINTS EXCEED MAXIMUM OF',14)
GO TO 700
347 CONTINUE
DO 350 IT=2,NPTS
T(IT) = T(IT-1) + DELT
350 CONTINUE
IT = NPTS
IF(ABS(REM).LT.TOL) GO TO 355
IT - IT + 1
T(IT) = TL
355 CONTINUE
TMAX = TL
IF(TF.GT.TL) TMAX=TF
COEF = 0/(6.2832»((RHOS-RHOF)/RHOF)»KX*D)
CONR = SORT(KZ/KX)/D
CONT = ((RHOS-RHOF)/RHOF)«KZ/(2.0*PO»D)
DO 365 1=1, IT
TAU = CONT»T(I)
TAU1 = CONT*(T(I)-TPUMP)
IF(T(I).LE.TPUMP) TAU1=0.0
XRO = COEF*(1 .0/(1 .0+TAU1 ) - '1 .0/( 1 .O+TAU) )
ZRO = XRO + ZO
DO 360 0=1. IR
RDIM = CONR*R(J)
Z(I,0) = COEF*(( 1
.0/SORT(( 1 .0+TAU1 )**2 •*• RDIM»*2))
360 CONTINUE
IF(ZRO.GT.ZCR) IPR=2
365 CONTINUE
370 IT = I
PRINT INTERFACE ELEVATIONS
WRITE ( NO , 375 ) 0 , I L , ITU , TPUMP , ITU . I L , I L
375 FORMAT( ' 1 ', 18X, 'PUMPING RATE : ' , F 12 . 2 , ' CU ' . A2 , ' / ' , A2 , ' FOR'
1F10. 2, A3, /, '0' , 18X, 'INTERFACE ELEVATIONS C,A2,') ',//,
23X, ' *' ,/,3X, ' * R ( ' ,A2, ' )' )
LIM1 = 1
LIM2 = 7
380 IF(LIM2.GT.IR) LIM2=IR
WRITE(NO,385) (R( L ) , L=LIM1 , LIM2 )
385 FORMAT(3X,' * '.7F9.2)
WRITE(NO,390) ITU
390 FORMAT(3X,'T C,A2,') » ' , / , 12X , ' « ' )
DO 400 1=1. IT
DO 393 L=LIM1 .LIM2
KFLG(L) = KHAR1
IF(Z(I,L) .GT.ZCR) KFLG(L)=KHAR2
393 CONTINUE
WRITE (NO, 395) T( I ) , ( Z( I , L) , KFLG( L) , L=LIM1 , LIM2)
395 FORMAT (5X.F8. 2, 1X , 7( F8 . 2, A1 ) )
40AX.GE.TCR)
1
405 FORMAT ( '0' ,2X,
* NOTE: CRITICAL ELEVATION OF',F8.2,A3,
1' EXCEEDED AT R=0 AND T-'.F8.2.A3)
IF(LIM2.EQ.IR) GO TO 415
LIM1 = LIM1 + 7
LIM2 = LIM2 + 7
UPCN351
UPCN352
UPCN353
UPCN354
UPCN355
UPCN356
UPCN357
UPCN358
UPCN359
UPCN360
UPCN361
UPCN362
UPCN363
UPCN364
UPCN365
UPCN366
UPCN367
UPCN368
UPCN369
UPCN370
UPCN371
UPCN372
UPCN373
UPCN374
UPCN375
UPCN376
UPCN377
UPCN378
UPCN379
UPCN380
UPCN381
UPCN382
UPCN383
UPCN384
UPCN385
UPCN386
UPCN387
UPCN388
UPCN389
UPCN390
UPCN391
UPCN392
UPCN393
UPCN394
UPCN395
UPCN396
UPCN397
UPCN398
UPCN399
UPCN400
UPCN401
UPCN402
UPCN403
UPCN404
UPCN405
UPCN406
UPCN407
UPCN408
UPCN409
UPCN410
UPCN411
UPCN412
UPCN413
UPCN414
UPCN415
UPCN416
UPCN417
UPCN418
UPCN419
UPCN420
63
-------
WRITE(NO,410) IL.IL
410 FORMAT('1',18X,'INTERFACE ELEVATIONS (',A2,') (CONTINUED)',//,
13X, ' " ,/,3X, ' T R (',A2, ' )' )
GO TO 380
415 CONTINUE
IF(ICON.NE.NY) GO TO 700
CONCENTRATION PROFILES
SO = S02/2.0
E(1) = 0.0
C(1) = CB
DO 420 K=2,11
E(K) = E(K-1) + 0.1
C(K) = E(K)*(CO-CB) + CB
420 CONTINUE
C
UP
CONT*TPUMP
XBAR1 = COEF*(1.0-1.0/(1.0+TAUP))
TAU1 = CONT*(T(I)-TPUMP)
XTOT = COEF«((1.0/(1.0+TAU1)) - (1.0/(1.0+TAU)))
XBAR2 = XBAR1 - XTOT
XBAR = XBAR1 + XBAR2
GO TO 440
435 XTOT = COEF*(1.0 - 1.0/(1.0+TAU))
XBAR = XTOT
440 CONTINUE
S1 = SORT(SO**2 + 2.0*DM*XBAR)
ARG » 10.0
IF(S1.GT.0.0) ARG = (XCR-XTOT)/(1.414214*51)
EZCR = 0.5*ERFC(ARG)
EW(I) " 0.5*EZCR»PHI
CW(I) = EW(I)*(CO-CB) + CB
Z(I. 1) = XTOT + 2.5*S1 * 20
Z(I. 11) = XTOT - 2.5*S1 + ZO
XLIM1 = 2.0
XLIM2 = 0.0
DO 445 K=2.5
CERF = 2.0*E(K)
CALL IERFC(CERF,ARG.XLIM1.XLIM2)
DIST » 1 ,41421«S1*ARG
2(1.K) = XTOT + DIST + ZO
L = 12-K
Z(I.L) = XTOT - DIST + ZO
445 CONTINUE
Z(I.S) = XTOT + ZO
446 CONTINUE
LIM1 = 1
LIM2 = 6
IN WELL AND PROFILES BENEATH WELL',
1' (',3A2,')'.//. 11X,'T C(WELL)'.14X,'ELEVATION FOR C/(E) ('.
2A2,') ' ,/.9X.'( ' ,A2. ' ) E(WELD')
WRITE(NO,430) (C(K),K=LIM1,LIM2),(E(K),K=LIM1,LIM2)
430 FORMAT(24X,6(F8.1,1X),/.24X,S(' (',F3 . 1 .')'))
431 DO 455 1 = 1 .IT
DO 449 K=LIM1,LIM2
KFLG(K) = KHAR1
IF(Z(I.K).GT.ZCR) KFLG(K) = KHAR2
449 CONTINUE
WRITE(NO,450) T(I).CW(I).(Z(I,K),KFLG(K),K=LIM1,LIM2),EW(I)
450 FORMAT(/,6X,F8.2,F9.2,1X.6(F8.1,A1),/.16X,('(',F6.4.')'))
455 CONTINUE
WRITE(NO,457) ZCR.IL
457 FORMAT('0',2X,'* NOTE: THE DISPERSION CONCEPT SHOULD BE LIMITED',
1' TO THE ZONE BELOW',/,11X,' THE CRITICAL ELEVATION OF',F12.4,A3)
IF(LIM2.EQ.11) GO TO 7OO
LIM1 = 6
LIM2 = 11
GO TO 447
UPCN421
UPCN422
UPCN423
UPCN424
UPCN425
UPCN426
UPCN427
UPCN428
UPCN429
UPCN430
UPCN431
UPCN432
UPCN433
UPCN434
UPCN435
UPCN436
UPCN437
UPCN438
UPCN439
UPCN440
UPCN441
UPCN442
UPCN443
UPCN444
UPCN445
UPCN446
UPCN447
UPCN448
UPCN449
UPCN450
UPCN451
UPCN452
UPCN453
UPCN454
UPCN455
UPCN456r
UPCN457
UPCN458
UPCN459
UPCN460
UPCN461
UPCN462
UPCN463
UPCN464
UPCN465
UPCN466
UPCN467
UPCN468
UPCN469
UPCN470
UPCN471
UPCN472
UPCN473
UPCN474
UPCN475
UPCN476
UPCN477
UPCN478
UPCN479
UPCN480
UPCN481
UPCN482
UPCN483
UPCN484
UPCN485
UPCN48G
UPCN487
UPCN488
UPCN489
UPCN490
64
-------
NCENTRATION IN PUMPED WATER
C
500 CONTINUE
IEDIT = 4
IF(KCON.NE.2) GO TO 176
510 WRITE(NO,515) IM1.IM2.IM3
515 FORMAT('0',2X,'ENTER MAXIMUM PERMISSIBLE CONCENTRATION IN PUMPED'.
1' WATER ( ' ,3A2.')',/,3X,'?')
516 READ(NI,50) CMAX
IF(CMAX.GE.CB.AND.CMAX.LT.CO) GO TO 518
WRITE(N0.517)
517 FORMAT(3X,'CONCENTRATION MUST BE GREATER THAN OR EQUAL TO CB',/.
13X,'AND LESS THAN CO -- REENTER',
2/.3X.'?')
GO TO 516
518 CONTINUE
SO » SO2/2.O
EMAX = (CMAX-CB)/(CO-CB)
EXCR = EMAX/(0.5*PHI)
ELIM = 0.0
IF(SO.GT.O.O) ELIM = 0.5*ERFC(XCR/(1.1414214*50))
CINIT = 0.5*PHI*ELIM«(CO-CB) + CB
IF(EXCR.LT.ELIM) GO TO 550
IF(EXCR.LE.O.O.OR.EXCR.GE.1.0) GO TO 520
CERF = 2.0*EXCR
XLIM1 = 3.0
XLIM2 = 0.0
CALL IERFC(CERF,ARG,XLIM1,XLIM2)
B = -4.0*ARG*ARG*DM - 2.0*XCR
CON = -2.0*ARG*ARG*SO*SO + XCR*XCR
GO TO 521
520 B = -12.5*DM - 2.0*XCR
CON = -6.25*SO*SO + XCR*XCR
521 CONTINUE
ROOT = B*B - 4.0*CON
IF(ROOT.LT.O.O) GO TO 65O
XBAR1 = (-B-(ROOT**0.5))/2.0
XBAR2 = (-B+(ROOT*«0.5))/2.0
XBAR = XBAR1
IF(EXCR.GT.O.S) XBAR=XBAR2
ZMAX = XBAR + ZO
JFLG = KHAR1
IF(ZMAX.GT.ZCR) JFLG=KHAR2
OMAX = 6.283185*((RHOS-RHOF)/RHOF)*KX»D«XBAR
WRITE(NO,280) (A(I),I=1,30).RHOF,RHOS,PO,IL.ITU,KX,IL,ITU,KZ
WRITE(NO.281) IL.ZO,IL,0.THETA,IL.XCR,IL,ZCR
WRITE(NO.290) IM1,IM2,IM3,CO,IM1,IM2,IMS,CB,IL,S02,IL.DM,PHI
WRITE(NO.525) IM1,IM2,IM3,CMAX,EMAX,IL,ZMAX,JFLG,IL,ITU,OMAX
525 FORMAT('0',5X,'MAXIMUM CONCENTRATION IN PUMPED WATER C.3A2,')'
11X.F10.4./,
26X,'MAXIMUM RELATIVE CONCENTRATION',17X,F10.4,/,
36X, 'MAXIMUM INTERFACE ELEVATION (',A2, ')',15X,F10.4,A 1,/,
46X,'MAXIMUM PERMISSIBLE PUMPING RATE (CU '.A2,'/',A2,')',4X,
5F10.4)
IF(ZMAX.LE.ZCR) GO TO 530
WRITE(NO,457) ZCR.IL
CLIM = 0.5«(0.5*PHI)*(CO-CB) + CB
WRITE(NO,527) CLIM,IM1,IM2,IMS
527 FORMAT(3X,'(MAXIMUM CONCENTRATIONS IN PUMPED WATER LESS THAN',
1F10.2,A3,2A2,') './.'O')
(CU ',A2,'/',A2,')'
530 WRITE(NO,535) IL.ITU
535 FORMAT('0',2X,'ENTER OPTIONAL PUMPING RATE
1/.3X,'?' )
READ(NI,50) OP
IF(QP.LE.QMAX) GO TO 600
TIME = ((2.«PO*D)/(((RHOS-RHOF)/RHOF)*KZ))*((1./(1.-QMAX/QP))-1.)
UPCN491
UPCN492
UPCN493
UPCN494
UPCN495
UPCN496
UPCN497
UPCN498
UPCN499
UPCN500
UPCN501
UPCN502
UPCN503
UPCN504
UPCN505
UPCN506
UPCN507
UPCN508
UPCN509
UPCN510
UPCN511
UPCN512
UPCN513
UPCN514
UPCN515
UPCN516
UPCN517
UPCN518
UPCN519
UPCN520
UPCN521
UPCN522
UPCN523
UPCN524
UPCN525
UPCN526
UPCN527
UPCN528
UPCN529
UPCN530
UPCN531
UPCN532
UPCN533
UPCN534
UPCN535
UPCN536
UPCN537
UPCN538
UPCN539
UPCN540
UPCN541
UPCN542
UPCN543
UPCN544
UPCN545
UPCN546
UPCN547
UPCN548
UPCN549
UPCN550
UPCN551
UPCN552
UPCN553
UPCN554
UPCN555
UPCN556
UPCN557
UPCN558
UPCN559
UPCN560
65
-------
WRITE(NO,545) IL,ITU,OP,ITU,TIME UPCN561
545 FORMAT('0',2X,'PUMPING RATE (CU ',A2,'/'.A2,')',6X,F10.4,/ UPCN562
1 , UPCN563
GO TO 530 UPCN564
C UPCN565
550 WRITE(NO,555) S02,IL,CINIT,IM1,IM2,IM3 UPCN5S6
555 FORMAT('0'.2X,'CMAX EXCEEDED FOR AN INITIAL TRANSITION '. UPCN5S7
1'ZONE WIDTH OF'.F10.4.A3,/. UPCN568
23X,'INITIAL CONCENTRATION IN PUMPED WATER IS',F10.4,A3,2A2) UPCN569
600 CONTINUE UPCN570
GO TO 700 UPCN571
C UPCN572
650 WRITE(NO,655) UPCN573
655 FORMAT(3X.'IMAGINARY ROOT OBTAINED FOR XMAX') UPCN574
GO TO 700 UPCN575
C UPCN576
C UPCN577
C *«**» SECTION III -- PROBLEM REDEFINITION AND CONTROL OF EXECUTION UPCN578
C UPCN579
C UPCN580
700 CONTINUE UPCN581
IEDIT = 2 UPCN582
XCR = THETA*D UPCN583
ZCR = XCR + ZO UPCN584
QMAXSS = 6.283185«((RHOS-RHOF)/RHOF)*KX*D*XCR UPCN585
WRITE(NO,705) UPCN586
705 FORMAT(//,3X,'ENTER NEXT COMMAND',/,3X.'?') UPCN587
710 READ(NI,715) NEXT UPCN588
715 FORMAT(A2) UPCN589
C UPCN590
DO 720 1=1,20 UPCN591
IF(NEXT.EO.IC(I)) GO TO 730 UPCN592
720 CONTINUE UPCN593
WRITE(NO,725) UPCN594
725 FORMAT(3X.'ERROR IN LAST COMMAND — REENTER',/,3X,'?') UPCN595
GO TO 710 UPCN596
730 GO TO (29.70,90,111.120,129,145,188,215.235,255, UPCN597
1311,312,322,302,275.300,5OO,1,10OO),I UPCN598
C UPCN599
1000 STOP UPCN600
END UPCN601
66
-------
APPENDIX D
Listing of Function Subroutines
67
-------
FUNCTION ERFC(Z) ERFC001
C RATIONAL APPROXIMATION OF THE COMPLIMENTARY ERROR FUNCTION ERFC002
C SEE SECTION 7.1 OF... ABRAMOWITZ AND STEGUN (1966) ERFC003
IF(ABS(Z).GT.7.5) GO TO 30 ERFC004
C THE FOLLOWING IDENTITIES ARE USED TO HANDLE NEGATIVE ARGUMENTS ERFC005
C ERFC( Z) = 1 - ERF(Z) ERFC006
C ERFC(-Z) = -ERFC(Z) = 1 + ERF(Z) ERFC007
C ERFC008
X = Z ERFC009
C NEGATIVE ARGUMENTS ERFC010
IF (Z.LT.0.0) X = -Z ERFC011
ERFC = 1.0/((1.0 «• 0.070523*X + 0.042282*(X*«2) ERFC012
1 + 0.009270*(X«»3) + 0.000152*(X«*4) ERFC013
2 + O.000276*(X*»5) + 0.000043*(X«*6))**16) ERFC014
C NOTE: 2-ERF(X) = ERFC(-X) = ERFC(Z) FOR Z<0 ERFC015
IF (Z.LT.0.0) ERFC = 2.0 - ERFC ERFC016
RETURN ERFC017
C ERFC018
C FOR Z>7.5, ERFC(Z)<2.32E-22 AND IS SET TO 0 ERFC019
30 ERFC « 0.0 ERFC020
C FOR Z < -7.5 ERFC(Z) IS SET TO 2.0 ERFC021
ERFC =2.0 ERFC022
RETURN ERFC023
END ERFC024
68
-------
C INVERSE COMPLIMENTARY ERROR FUNCTION IERF001
SUBROUTINE IERFC(CERF,X,XLH.XRH) IERF002
C INVERSE COMPLIMENTARY ERROR FUNCTION - IERF003
C A REGULA FALSI ROOT-FINDING TECHNIQUE IS USED TO IERF004
C LOCATE X FOR A SPECIFIED VALUE OF CERF. XLH AND IERF005
C XRH DEFINE THE INITIAL SEARCH INTERVAL. IERF006
C SEE CARNAHAN, LUTHER AND WILKES (1969) FOR A IERF007
C DISCUSSION OF THE METHOD. IERF008
FLH = CERF - ERFC(XLH) IERF009
FRH = CERF - ERFC(XRH) IERF010
10 XEST = (XLH*FRH - XRH*FLH)/(FRH - FLH) IERF011
FEST = CERF - ERFC(XEST) IERF012
IFUBS(FEST) .GT. 1 .OE-O4) GO TO 20 IERF013
X = XEST IERF014
RETURN IERF015
20 IF(FEST.GT.O.O.AND.FRH.GT.O.O) GO TO 30 IERF016
XLH = XEST IERF017
FLH = CERF - ERFC(XLH) IERF018
GO TO 10 IERF019
30 XRH = XEST IERF020
FRH = CERF - ERFC(XRH) IERF021
GO TO 10 IERF022
END IERF023
69
------- |