United States EPA-600/7-82-037!
Environmental Protection
Agency May 1982
dEPA Research and
Development
VERIFICATION AND TRANSFER OF
THERMAL POLLUTION MODEL
Volume VL User's Manual for
One-dimensional Numerical Model
Prepared for
Office of Water and Waste Management
EPA Regions 1-10
Prepared by
Industrial Environmental Research
Laboratory
Research Triangle Park NC 27711
-------
RESEARCH REPORTING SERIES
Research reports of the Office of Research and Development, U.S. Environmental
Protection Agency, have been grouped into nine series. These nine broad cate-
gories were established to facilitate further development and application of en-
vironmental technology. Elimination of traditional grouping was consciously
planned to foster technology transfer and a maximum interface in related fields.
The nine series are:
1. Environmental Health Effects Research
2. Environmental Protection Technology
3. Ecological Research
4. Environmental Monitoring
5. Socioeconomic Environmental Studies
6. Scientific and Technical Assessment Reports (STAR)
7. Interagency Energy-Environment Research and Development
8. "Special" Reports
9. Miscellaneous Reports
This report has been assigned to the INTERAGENCY ENERGY-ENVIRONMENT
RESEARCH AND DEVELOPMENT series. Reports in this series result from the
effort funded under the 17-agency Federal Energy/Environment Research and
Development Program. These studies relate to EPA's mission to protect the public
health and welfare from adverse effects of pollutants associated with energy sys-
tems. The goal of the Program is to assure the rapid development of domestic
energy supplies in an environmentally-compatible manner by providing the nec-
essary environmental data and control technology. Investigations include analy-
ses of the transport of energy-related pollutants and their health and ecological
effects; assessments of, and development of, control technologies for energy
systems; and integrated assessments of a wide range of energy-related environ-
mental issues.
EPA REVIEW NOTICE
This report has been reviewed by the participating Federal Agencies, and approved
for publication. Approval does not signify that the contents necessarily reflect
the views and policies of the Government, nor does mention of trade names or
commercial products constitute endorsement or recommendation for use.
This document is available to the public through the National Technical Informa-
tion Service, Springfield, Virginia 22161.
-------
EPA-600/7-82-037f
May 1982
VERIFICATION AND TRANSFER
OF THERMAL POLLUTION MODEL
VOLUME VI: USER'S MANUAL FOR ONE-DIMENSIONAL
NUMERICAL MODEL
By
Samuel S. Lee, Subrata Sengupta
and Emmanuel V. Nwadike
Department of Mechanical Engineering
University of Miami
Coral Gables, Florida 33124
NASA Contract No. NAS 10-9410
NASA Project Manager: Roy A. Bland
National Aeronautics and Space Administration
Kennedy Space Center
Kennedy Space Center, Florida 32899
EPA Interagency Agreement No. 78-DX-0166
EPA Project Officer: Theodore G. Brna
Industrial Environmental Research Laboratory
Office of Environmental Engineering and Technology
Research Triangle Park, North Carolina 27711
Prepared for:
U. S. Environmental Protection Agency
Office of Research and Development
Washington, D. C. 20460
-------
PREFACE
Emphasis continues to be placed on the use of digital computers in
solving nonlinear hydrodynamic and thermodynamic equations of fluid
flow. This publication of the thermal pollution group at the University
of Miami presents the solution of one such problem. This problem deals
with the use of a numerical one-dimensional model in predicting the tem-
perature profiles of a deep body of water. Although this model can be
applied to most lakes, a specific site (Lake Keowee, S. C.) application
has been chosen and described in detail. The programs are written in
fortran V and could be modified by the user. Some of these modifications
are suggested either in the text or in the specific programs.
A detailed derivation of the equations integrated has been left out;
however, to improve readability of the final equations, the meaning of
the terms and variables occurring in these equations are included.
This research was performed at the thermal pollution laboratory at
the University of Miami. Funding was provided by the National Aeronau-
tics and Space Administration (NASA-KSC) and the Environmental Pro-
tection Agency (EPA-RTP).
-------
ABSTRACT
A user's manual for a one-dimensional thermal model is described.
The model is essentially a set of partial differential equations which are
solved by finite difference methods using a high speed digital computer.
The main equations integrated are discussed. The programs are written
in fortran V and an example problem is discussed in detail.
iii
-------
CONTENTS
Forward
Preface }i
Abstract iii
Figures v
Symbols vi
Acknowledgments vii
1. Introduction 1
2. Recommendations 3
3. Program Description and Flow Chart 4
Description of program algorithm * ] 4
Background 4
Algorithm . 4
Main program and subroutines 6
4. Description of Program Symbols 9
Introduction 9
Description of main variables 10
5. Preparation of Input Data 12
6. Plotting Programs and Execution Elements 13
Plotting programs 13
Input data 13
References 17
Appendices 18
A. Example Problem ^8
Site description 18
Problem statement 18
Calculations of parameters and input data 18
B. Fortran Program Listings 27
iv
-------
FIGURES
Number Page
1 Ideal!zed deep body of water 5
2 Flow chart (calculation) 7
3 Flow chart (plots) 16
4 Lake Keowee 20
5 Sample output - Lake Keowee, 1971 2i|
6 Sample plots - measured average temperature profiles
(Stations 500-506) vs predicted temperature profiles,
Lake Keowee, 1971 25
7 Sample plot 26
-------
SYMBOLS
h
A(z)
Kz)
Q(z)
T
P
KZO
w*=
R
av
Vertical coordinate measured
upward from deepest point
of the lake. .As a subscript
it marks the vertical compo-
nents of a vector.
Depth of lake
Horizontal cross-sectional
area at height Z
Bottom-surf ace source of
mass per unit area
Bottom-surface source of
heat per unit area
Temperature (°C)
Density of water
Vertical velocity
Eddy diffusivity
Eddy diffusivity under
neutral condition
(T, } Friction velocity
empirical constant
Richardson number
Volumetric coefficient of
expansion of water
Surface shear stress
C
H?z)
A1
8
C C
71' .
n
6
2f
K
3B
Heat capacity
Heat source/unit volume
Average value of W*
Half of the annual variation W*
,C3/C4/C_ Phase angles
SoTar radiation incident on
the water surface
Average value of $
Half the annual variation of
Extinction coefficient °
Absorption coefficient
Volumetric discharge
Condenser temperature change
Discharge temperature
Surface heat flux
Surface heat exchange
coefficient
Equilibrium temperature
Average value of T_
Half the annual variation of T_
Surface temperature
Bottom surface heat flux
Lake surface radius
Area variation with depth
vi
-------
ACKNOWLEDGMENTS
This work was supported by a contract from the National Aeronau-
tics and Space Administration (NASA-KSC) and the Environmental Pro-
tection Agency (EPA-RTP).
The authors express their sincere gratitude for the technical and
managerial support of Mr. Roy A. Bland, the NASA-KSC project manager
of this contract, and the NASA-KSC remote sensing group. Special
thanks are also due to Dr. Theodore G. Brna, the EPA-RTP project
manager, for his guidance and support of the experiments, and to Mr.
S. B. Hager, Chief Engineer, Civil-Environmental Division, and Mr.
William J. McCabe, Assistant Design Engineer, both from the Duke Power
Company, Charlotte, North Carolina, and their data collection group for
data acquisition. The support of Mr. Charles H. Kaplan of EPA was
extremely helpful in the planning and reviewing of this project.
vii
-------
SECTION 1
INTRODUCTION
It is important that the thermal behavior of heated discharges and
their receiving basins be clearly understood.
A numerical model that can be used for predicting the seasonal thermo-
cline of a deep body of water is very useful in studying the environmental
impact of thermal discharges from power plants. This is not only required
for existing power plants but also for planned units. Thus, a predictive
capability is essential to the licensing procedure. Monitoring programs
cannot satisfy these needs, but from time to time, play a vital role in the
calibration and verification of mathematical models.
The one-dimensional, thermal numerical model, described in this
manual, features the effects of area change with depth, nonlinear inter-
action of wind-generated turbulence and buoyancy, absorption of radiative
heat flux below the surface, thermal discharges and the effects of vertical
convection caused by discharge. The main assumption in the formulation
of this model is horizontal homogeneity.
This model can be applied to most stratified deep bodies of water.
This stratification has a seasonal cycle and is an important natural
characteristic of a body of water. The body of water could be divided
into any number of slices. The temperature of each slice is predicted
by the model. The surface slice exchanges heat with the environment
of known climatic conditions while the bottom slice is assumed perfectly
insulated. Condenser cooling water is extracted from any one of the
slices and heated by the power plant. The discharge is injected into a
slice of the same temperature as the dischage.
The main function of the model is the prediction of the temperature
profiles in a deep body of water for any number of annual cycles. How-
ever, predictions connot be made on hourly basis - a feature usually
handled by a more sensitive three-dimensional model. This is the main
limitation of the model.
The procedure used in writing this manual is as follows:
Description and flow chart of the main program are given in Section
3, where the subroutines are also described. In the next section, a
list of the variables and dimensions are given. The next three sections
7
-------
show how a typical run is prepared, executed and plotted. An example
case is discussed in Appendix A, while Appendix B gives the fortran
source program listings.
-------
SECTION 2
RECOMMENDATIONS
The main disadvantage of a one-dimensional thermal model lies in
the fact that resolution is sacrificed for computational speed. Three di-
mensional models are bulky and time consuming but have much better
resolution, however, when long term simulations are necessary, a one-
dimensional mode! is recommended.
The model described here can be modified to include the single
effects of the various quantities involved in the surface heat transfer
phenomenon rather than using the equilibrium temperature concept.
This is particularly recommended for the user who is interested in
modeling the long term effects of one (for example, evaporation) of
the quantities involved in the surface heat transfer processes.
Furthermore, the model can be easily adapted to handle connected
multiple domains. This recommendation is discussed in the text.
-------
SECTION 3
PROGRAM DESCRIPTION AND FLOW CHART
DESCRIPTION OF PROGRAM ALGORITHM
Background
A view of an idealized deep body of water is shown in Figure 1.
This basin is divided into eleven slices. The inner nine slices are of
equal thickness, DZ, while the top and bottom slices are of thickness
DZ/2. The thickness, DZ, is determined from the depth of the basin
and the number of slices used. The temperature of each slice is as
shown in Figure 1; the horizontal lines correspond to the center of each
slice.
The condenser cooling water (CCW), if any, could be taken from
any slice. In Figure 1, the CCW is extracted from the center of Slice 2
which is at temperature T_. The discharge temperature, TD, is the sum
of T, and the increase in temperature through the condenser. TQ is in-
jectea into a slice of equal temperature or treated as a surface outfall if
T~ is greater than the highest temperature of the basin.
The basin also gains or loses heat from the surface as a result of
changing climatic conditions which are required as input data. These
could vary every time step, daily or monthly.
Algorithm
The problem is an initial value problem, so the values of dependent
variables are assumed known initially. The governing and associated
equations are discussed in the next section. The governing equation is
parabolic and mathematically represents a diffusion process with vertical
convection.
The values of the dependent variables at successive time steps are
obtained by using a forward-time Dufort-Frankel scheme.
The sequence in which calculations are performed is as follows:
(Refer to Summary of Variables - next section.)
1. The dependent variables, T, KZ/ W*, Ay, p, T£ and Kg, are ini-
tialized. The area of each slice is calculated ana then tne time step
-------
T12.= TSf- A12
" From
Condenser
(Cooling Water)
If TD = T 9
Tl , Al
To Condenser
(Cooling Water)
Figure 1. Idealized deep body of water
-------
is calculated. The heading of the beginning year is printed. The
values of the variables, K^, W*, AV, p, T- and K., are then calcu-
lated. The temperatures of the slices are finally calculated. If the
temperature profile is unstable, mixing of the unstable portion of
the profile is undertaken.
2. During the next time step, the temperatures are updated, and the
dependent variables are calculated again.
3. The values of the temperature T, eddy diffusivity KZ/ number of
days and surface heat transfer coefficient K- are printed every time
step, every day or normally at the end of each month. At the end
of the present year, the title of the new year is printed and compu-
tions continue as listed above. These steps are shown in a flow chart,
Figure 2. The results are stored on a magnetic tape and plotted
when necessary.
Description of Main and Subprograms
The fortran calculation programs consist of a main program (NASA)
and seven subroutines (YEARS, EQUIL1, STORE, CCW, SMOOTH, MIXIT
and AREAS).
1. MAIN: The main program handles the input data, calls the subrou-
tines and does the temperature calculations. Two alternatives are
given for handling the input data; these are either read through
cards or 5n-data files or through a block-data arrangement given
at the beginning of the main program. For users interested in the
block-data package, the following caution is necessary: Whenever
a data or set of data is changed, the main program must be recom-
piled!
2. YEARS: This subroutine prints the year heading. It is called at
the beginning of a new year.
3. EQUIL1: This subroutine reads the dewpoint temperature, wind
speed and solar radiation. It then computes the surface heat trans-
fer coefficient and the equilibrium temperature. Depending on how
the data has been averaged (e.g. days, months or years); it is
called as often as needed.
4. STORE: This subroutine stores the calculated data on magnetic
tape designated as Unit 8. The stored data could be read by the
plotting subroutine called READER. This subroutine and other plot
programs are described later.
5. CCW: This subroutine supplies the condenser cooling water data.
The data is also converted to the required units by this subroutine.
6. SMOOTH: This subroutine finds the largest value of the eddy dif-
-------
Increase
Time
No
Initialize the Dependent
Variables
Calculate the Surface
Areas of the Slices
i
Calculate Time Step
I
Print Start Year
Calculate K7, W*, A ,
p, TE. and Kg v
1
Integrate Governing
Equation for Temperature
I
If Temperature Profile is
Unstable, MIX
1
Increase
Time
Increase
Month 5
Year
No
New Year?
Yes
Print Year
Title
No
Time = 30 Days?
Yes
Print Time 5
Dependent
Variables
Is Time
Equal to
Desired
Time
STOP
Yes
Figure 2. Flow chart (calculation)
-------
fusivity and uses it to calculate the variable time step. It also
smoothens the calculated eddy diffusivity for unstable temperature
gradients. It is called every time step.
7. MIXIT: This subroutine looks for unstable temperature gradients
and mixes or stabilizes the temperatures. It is also called every
time step.
8. AREAS: This subroutine handles the surface areas of each slice
and converts the values to the required units. It is called only
once at the beginning of the computations.
9. INPUT: This is an in-data element containing all input data.
-------
SECTION H
DESCRIPTION OF PROGRAM SYMBOLS
Introduction
The programs have been written to calculate, as a function of depth,
thermal diffusivity and temperature profiles over complete annual cycles.
The equation integrated is
A(z)lt{pCPT) -f^1 -<>cATv) +QA'
The above equation requires two boundary conditions and one ititial con-
dition.
The initial condition is an input quantity supplied by the user and
equals the homothermal temperature of the basin. The boundary condi-
tions are:
1. At the surface;
KzI|Z=h=KS(TE-TS>
where Z = vertical coordinate measured from the deepest point
Tg = equilibrium temperature
TC = surface temperature
Kg = surface heat exchange coefficient
2. At the bottom;
Perfect insulation is assumed,
3T - n f31
TF|Z=0 - ° (3)
Calculations of the temperature profiles are made by numerical integration
of Equation (1). Calculations start with the homothermal conditions and
a forward explicit scheme is used.
Each time step, the surface temperature, TS = T12, is calculated
-------
and then the temperature of each slice is calculated. Solar radiation is
absorbed at the surface slice and the unabsorbed portion is transmitted
exponentially to the slices below.
The empirical relations involved in this manual are summarized below.
A full discussion is given in the final report. Lee et al. (1980).
Description of Main Variables
1. Density, p, fortran variable - ROW:
p = Aj + BjT + C^T2 (4)
where A. = density at 0°C
= 1.02943 gm/cc
B. = constant
= -0.00002
C- = constant
= -0.0000048
2. Eddy diffusivity, KZ/ fortran variable = XKZ
Kz = KZQ(1 + ^R.)-! (5)
and
\ - W*z 3Z
where R. = Richardson number
a' ~ 0.1, an empirical constant, fortran variable - SIGMA
g = acceleration due to gravity, fortran variable - G
W*= friction velocity, fortran variable - FRVEL
= (ts/p)
fortran variable for
B2(T - 4) + C2(T - 4)2 (6a)
a.,, AV
where A- = 0, volumetric coefficient of expansion at 4°C, fortran
variable - A1
B, = constant, fortran variable - A 2
z = 1.538 x 10"3
C- = constant, fortran variable - A3
= -2.037 x 10~7
ct can also be estimated by using Equation (4).
10
-------
where K_Q = eddy diffusivity under neutral condition (varies with
time), fortran variable - XKZO
KZO = A3 * B3 Sln(5t * C3} Ct is in days) (6b)
where A. = average value of KZQ, fortran variable - R9
B- = half annual variation of K_o/ fortran variable - R10
C- = phase angle, fortran variable - R8
3. Heat source, H, fortran variable - F6
H = n(1 - 6)A(Z)oexp(-n(Z - h)) (7)
where 3 = 0. 5, fraction of the solar radiation absorbed at the surface
n = 0. 75, solar radiation absorption coefficient
-------
SECTION 5
PREPARATION OF INPUT DATA
The input data is stored in an in-data file - INPUT. Alternatively,
it could be punched on cards. The input data is read in with an open
format. The main variables read are: dewpoint temperature, wind speed
and solar radiation. In some cases where the dewpoint temperature is
not available, the relative humidity, air temperature and a pschometric
chart are used to find the dewpoint temperature. If this involves a lot
of chart reading, subroutine EQUIL1 could be modified and the dewpoint
temperature calculated from a known equation supplied by the user.
if the latter case is used, then the input data base is enlarged to read
air temperature, relative humidity, wind speed and solar radiation. A
detailed input list of the constants is given in Appendix A.
12
-------
SECTION 6
PLOTTING PROGRAMS AND EXECUTION ELEMENTS
DESCRIPTION OF PROGRAMS
The fortran plotting routine consists of one main program (PLOTTER)
and one subroutine (READER).
PLOTTER: This program calls the calcomp fortran subroutines (re-
fer to a Calcomp plotting manual for details) and the subroutine (READER)
which reads the calculated results from a magnetic tape designated as
Unit 8. (See Item A. 4.) A flow chart is shown in Figure 3.
READER: Reads the calculated data stored on Unit 8 (magnetic
tape).
Execution Elements
Two execution elements are used, one for executing the calculated
results and the other for executing the plots.
DO-IT: This element compiles and prints the main program (NASA)
and then prepares an entry point table, maps the necessary programs
and subprograms, calls the in-data element containing the input data
and finally, executes the calculations. This is done as follows for a
UNIVAC 1100 computer at the University of Miami.
Only one magnetic tape is necessary.
1. @ ASG, AX FILE.
The 'FILE1 is assigned for the run.
2. 9 ASG, T 8., 16N, TAPENAME
A magnetic tape file named '8.' is being assigned. The tape is 9-
track, and the reel number is 'TAPENAME1. The calculated results
are stored on this tape.
3. @ PRT, S FILE.NASA
The main program is printed.
13
-------
1. 9 PACK FILE.
The 'FILE' is packed.
5. 0 PREP FILE.
The entry point table is prepared.
6. @ MAP, S
7. IN FILE.NASA
8. LIB FILE.
9. END
10. 0 XQT
11. 0 ADD FILE.INPUT
12. 0 FIN
PLOT-IT: Similar to DO-IT, but handles the plotting executions.
For a UNIVAC 1100 computer the following cards are necessary. Two
magnetic tapes are necessary.
1. 0 ASG, AX FILE.
2. 0 ASG, T 8., 16N, TAPENAME
3. ©ASG, T 11., 16, PLOTTAPE
A magnetic tape file named '11.' is being assigned. The tape is 7-
track, and the reel number is 'PLOTTAPE1. The plots are stored on
this tape.
4. 0 PRT, S FILE.PLOTTER
The plot program is printed.
5. 0 PACK FILE.
6. 0 PREP FILE.
7. 0 MAP, S
8. IN FILE.PLOTTER
9. LIB FILE.
10. END
-------
11. 0 XQT
12. @ ADD FILE.INPUT
13. 9 FIN
15
-------
INITIALIZE COUNTERS AND
ENTER STARTING YEAR
CALL SUBROUTINE 'READER'
TO READ MAGNETIC TAPE
READ MEASURED DATA IF
JRUN=1, OTHERWISE JRUN=2
PLOT TEMPERATURE PROFILE
FOR THE MONTH
PLOT STRATIFICATION CYCLE
FOR ALL YEARS
PRINT 'ALL PLOTS ARE NOW
COMPLETED NORMAL JOB EXIT'
NO
INCREASE
YEAR
PRINT 'PLOTS
FOR YEAR
COMPLETED
IS MONTH
DEC?
YES
STOP
YES
YEAR LESS
THAN 1980?
NO
Figure 3. Flow chart (plots)
16
-------
REFERENCES
Duke Power Company. Oconee Nuclear Station Environmental Summary
Report 1971-1976. Vol. 1. November 1977.
Sengupta S., Lee S. S. and E. V. Nwadike. A One-Dimensional Variable
Cross-Section Model for the Seasonal Thermocline. Proceedings of
the Second Conference on Waste Heat Management and Utilization.
p. 1X-A-3. December 1978.
Lee, S. S., Sengupta, S. and E. V. Nwadike. Verification of a One-
Dimensional Model for the Seasonal Thermocline at Lake Keowee.
NASA Contract NAS 10-9410. 1980.
17
-------
APPENDIX A
18
-------
APPENDIX A
EXAMPLE PROBLEM
The model described in this manual was verified using monthly-
averaged data supplied by Duke Power Company for Lake Keowee, South
Carolina. Accordingly, the data discussed below apply to Lake Keowee.
SITE DESCRIPTION
Lake Keowee is located 40 km west of Greenville, South Carolina.
It is the source of cooling water for Oconee Nuclear Station (ONS). It
was formed from 1968 through 1971 by damming the Little and Keowee
rivers. A connecting canal (maximum depth 30.5 m) joins the two main
arms of the lake. Flow out of the lake is through the Keowee Hydro
Station. Lake Keowee also exchanges water with Lake J ocas see-pumped
storage station. The three-unit ONS with a net capacity of 2580 Mwe
started operating in July 1973. ONS operated on annual gross thermal
capacity factors of 11, 28, 69 and 59% in the years 1973 through 1976,
respectively. From 1977 to 1979 the factors varied from 65 to 75%. A
map showing the geometry of the lake is given in Figure 4.
PROBLEM STATEMENT
Calculation of Parameters and Input Data
1. The fortran variable DM (I, J) is a two-dimensional array containing
the temperatures at the connecting channel between Lake Keowee and
the Jocassee-pumped storage station. The data is averaged monthly.
The units are in degrees Celcius (°C). I is the year counter and J
is the month counter. The inputs for the first year are punched on
the first card, the next year on the next card, and so on. Accord-
ingly, each card contains twelve inputs in open format (real floating
point numbers).
2. The following fortran variables/constants are also read in with open
format, five on one card.
IYEAR: starting year - 1971 (could be changed).
DZ: thickness of an inner slice (ft) - (maximum depth of lake)/(10.0)
XKZL: lower limit of the eddy diffusivity (ftVday) - corresponds to
19
-------
to
o
Figure 1. Lake Keowee
-------
the thermal diffusivity of solid water (15 ft2 /day).
H: maximum depth of lake, ft (150 ft).
G: acceleration due to gravity (ft/sec2).
PI: ir = 3.1415926.
A1: corresponds to A2 in Equation (6a); A1 = 0 °C .
A2: corresponds to B2 in Equation (6a); A 2 = 1.538 x 10~ 5 °C.
A3: corresponds to C2 in Equation (6a); A3 = -2.037 x 10~7 °C.
A4: corresponds to A1 in Equation (4); A4 = 1.02943 gm/cc °C.
A5: corresponds to B1 in Equation (4); AS = 0.00002 gm/cc °C.
A6: corresponds to C1 in Equation (4); A6 = -0.0000048 gm/cc °C2.
(NOTE: The units for A4 through A6 are automatically converted to
consistent units in the main program.)
TO: homothermal temperature of lake (initial condition); TO = 7.8 °C,
C : specific heat; C = 1. 8 BTU/lb °C.
SIGMA: see Equation (5); SIGMA =
-------
solar radiation. These can either be punched on cards or stored in
an in-data element. They are read every month. Each card contains
three members. For example: for January-March 1971 (Lake Keowee),
the data are
3.0, 6.69, 167.0
0., 9.3, 264.4
6.3, 9.28, 264.4
The first number on each line (each card) is the dewpoint tempera-
ture in °C. The second one is the wind speed in ft/sec. The third
quantity is the solar radiation in BTU/ft2day. If DATA1 = 1, a fourth
number must be included on each line (every card). This fourth
quantity Is the computed friction velocity for each month.
NOTE: The in-data element described above is called INPUT. (See
Fortran Source Program Listing, Appendix B.)
22
-------
Sample Output and Sample Plots
23
-------
- 1971 •
*
"OMH IS
3C 6.28 155.13
30 7.80. 7.80 7.di;. 7.HC 7. 80 7. 90 7.3,; 7.80 7.80 7.80 7.8C 7.8Q
BC1.7J 601.70 B01.70 831. 7G 8cl.7c 8T1.70 331.70 801.70 8C1.70 8ol.7n SC1.70 801.70 8C1.70
"OMH IS FES.
6C 12.19 165.52
bC 7.!9 7.19 7.3-» 7. .'9 7.39 7 . <9 7.39 7.39 7 . *9 7.39 7 . t 9 7.«9
/C1.S7 7C1.57 701.57 7G..57 7i-1.5J 701.57 731.57 73;. 57 701. ?7 fdl.57 7CI.=7 701.57 7C1.57
•GMH IS" fABCH '1971 -C . '* ' ' " 'C ^ 'C ' 3 '" -0
•iC .7.77 16M.C2
*C 7.59 J.6C 7.oC 7.t? 7.6b 7.77 7.93 8.16 H.Tt h.tb 8.5fc S.65
627. frC, 570. cO ^M.JJ :7t.6I 15C.r4 97.70 73.fi8 76.81 1 1M .c 5 221. .'0 -.16.75 59, .37 627. bu
•" »3 .0 «0 .0 .C .0---.G .0 .0 "
1971 '^ >u
t2C 7.71 J.71 7.73 7.fcO 7.92 «.!« 8.^5 9.35 1C. =2 11.22 1 1 . t J ll.fcfl
cuc.r.;_ J8.-4S ie.ts ia^s ni8.«s is^s is. is is.«5 is.us 10. oj 126.1-* «.97.si 6G3.ci
«u .C .0 .0 .3 .Q .0 .0 .3 .n n n
MOMh IS HAY 1971 'J *u •" •"
I5c 18.12 1H3.C8
15C ^-79 1.79 7.o3 7.93 8.12 8.50 9.21 lols? 1^.25 17. 65 18. 43 18. bl
625. Su 15. CO 15.00 15. QO 15.00 15. UC 1S.SG 15.00 15.00 15.00 19.22 34,3.56 625.90
10MH IS° JUNE •?»»! '° 'D '° '° 'C -°. •" -° -° . -0
loc i7.C7 21J.15
1?C 7>?? r^ I**I 7'99 6'15 8'<'7 9'11 JO. 31 12. SC 16. 3« 22.31 J2.89 23.1*.
698.61 15. CO 15.00 15. CO IE.CC ' 15. CC IS. CO IS. 00 15. OQ 15. CO '17.12 t4C.35 696.61
MONTH IS° JULY '^971 '° '° '° '° '° «° '° '° -° -0
210 27.15 2C1.71
210 8.H 8.11 8.2t 6.51 9.03 1C.C1 11.77 11.71 19.39 2 3 . 8 7 21.17 21 31
798.29 15. CO 15.00 15. CG 15. OC IS. 00 15.CO 15.00 1S.OQ IS.CO 32.60 622.68 798.29
•G «Q «o «o .c .0 .0 . .a .0 .0 n n
-ONIh li AUf.. 1971 -° •"
210 21.32 2C2.12
21C 8.18 8.18 8.6-4 9.C-. 9.77 11.15 13. 1C 16.93 22.11 21.2! 21. It 21. S«
d-*a.M2_ 15. co ^ is.oo is.rc ,.^-nc 15. GO is.r.c; ib.oo i5.ro js.ao -ia.^3 771.11 »9«.M^
•w •- •« .^ .^- .^ .0 .T .0 .P n n
-OMh li SEPT. 1971 >L >U «U
270 18.66 201.52
270 ». 98 8.98 9.2C 9.7" 13.73 12.12 11.07 IS. 91 22.67 22.98 23.07 23.12
972. J-, 1S.LO 15. OC !£.;„ 1E.CC 15.UG 1S.CC 15.30 IS. 'JO ;fc . M5 125. C9 926.57 972. J9
• c .C .c .0 ,c .0 ,n .n .0 .c o T
"ONIh IS OCT. 1971 'L *J 'a
3:,C 8. 31 167.98
9 ''tl.CO ''".00 '-5i.Cfl ""li.OC ll'lE.OO l3-t^.,C "'J^OO "'ti.il "i^.H "M.S* 19,,'L95 '7
-OMH IS NOV.
3JC 8.73 168.76
33C 1C. 12 1C. 12 1C. /I 11. 18 12.71 11.21 H.2J 11.23 11. 2J 11.23 11.23 11.23
971.11 I5.CO 15.00 15.QO 15.00 15.00 23.16 H3G.29 971.11 97i.H 97u.u 971.11 971.11
•C .0 .0 .0 .C .0 .0 .0 .0 .0 0
M3NTH IS DEC. 1971 '
36C 3.59 132. !5
36C 11.16 11.17 11. 5d 11.95 11.96 11.96 11.96 11.96 1). -J6 11.96 1 1 . V 6 11.96
9C1.10 171. S2 92.91 3G.79 18.07 716.19 931. 1C 901.10 9C1.1Q 901. 1C 931.10 901.10 901.10
Figure 5. Sample output - Lake Keowee, 1971
-------
TEMPERRTURE PROFILES FOR LRKE KEQWEE
(OEPTH IS MEfiSURiD FROM THE DEEPEST POINT OF THE LRKEJ
rSTHTIQNS 5QQ-50SJ
1971
o_
JRN
30 20.00
TEMP.ICJ
Q_
ES
MRY
_ IQ 20.00
a TEMP.CC)
o.
FEB.
uj U.QO 20.00
a TEMP.CC]
o
—10
*
o_
JUNE
10 20.00
a TEMP.fCJ
o,
Mo
a.
MRRCH
£] Ll.OO 20-00
a TEMP.(C)
o
i-iO
*°
o_
JULY
_ 10 20-00
a TEMP.(C)
RPRIL
a
"o
.00
TEMP
20.00
CO
TEMP
20-00
CCJ
o
,-•0
a.
SEP
) 2G.QO
• — w o t r* i
Cm « I L J
-
OCT
2Q.QO
[MP.tCJ
NOV .
30 20.00
TEMP.(C)
0
0.
DEC
TEMP
20-00
fCJ
Flqure 6. Sample plots - measured average temperature profiles (Stations
500-506) vs predicted temperature profiles. Lake Keowee, 1971
25
-------
STRR7
iFI-CRTIGN CYCLE FOR LRKE KEGWcc 1971-1979
N»
cn
o
o
Solid Lines (No Discharge)
Broken Lines (Discharge - Mid-layer Temperatures)
Equilibrium Temperatures
Surface TemperaturosA)|<)
Mid-layer Temperatures
36-00
TZ.OO
108.00 144.00 180-00
TIME IN DflYS *10l-
6.00 252.00 288-00 3Z4 .00
Figure 7. Sample plot
-------
APPENDIX B
FORTRAN PROGRAM LISTING
27
-------
DATA (MONTHS (J) , J- 1,12) /'JAN.', 'FEB. ', 'MARCH*, *
C' JULY' , "AUG. * , 'SEPT.* , * OC T .*, * NOV.*, *OEC. */
NASA SYM CREATED ON 12 A uG 80 AT I<4:l7:0b
1 C ONE DIMENSIONAL MODEL FOR THE SEASONAL THERMOCLINE
? C
3 C
4 DIMENSION T(20),AV(20),CB(2Q),Z(2a),A(20), XK 2(201, ROW (20> , TNI 20)
5 DIMENSION DM ( 20 ) , T 2 (20 ) , X TOO ( 1 0, 360)
6 DIMENSION L-EL T EM( 1 2 ) . QP( 12 )
7 CHARACTER*** MONTHSU2)
ft C
9 C
10 DATA (MONTHS (J) , J- 1,12) /'JAN.', 'FEB. ', 'MARCH*, * APRIL* , »MAY *," JUNE*,
11
12 C
13 (.
14 C IF YCJ NEED TO STORE RESJLTS ON MAGNETIC TAPE READ JRJl^l
15 C OTHERWISE JRUN-2.
16 C
17 READ 1 ,JRUN
18 READ 1 ,IYEAR,DZ,XKZL,H, G
19 READ 1 ,PI , Al ,A2,A3 ,A<*
2(1 READ 1 ,A5 ,A6 ,TO,CP .SIGMA
21 READ 1 ,R6 , W7 , R 8 ,R9 ,R1 0
22 1 FORMAT ( )
23 MMI=U
24 Z (1 ) -0.
2r> JIM-1
26 TDL-D.
27 DVE^C.
28 CALL AREAS (A)
29 J-l
30 JWrl
31 JJ-U
12 NDAYS-0
33 NDAYSl^O
34 TIMf-0.
35 TIME1-U*
36 TIME7-0.
37 TIHE3=Q.
3? T1MC4 -0.
39 TErTO
4H DO 20 1=1 ,12
m T(I)-TO
42 T2 (I )-TO
^:5 20 CONTINJC
M4 DO 22 1-2,11
M5 Z (I ):DZ/2. +( 1-2 )*DZ
-------
NJ
to
46
47
48
49
50
u
56
57
56
59
60
Ll
62
c3
64
C5
66
67
L8
69
7n
71
72
73
74
75
76
77
78
79
bO
el
£2
S3
hH
65
66
«9
vn
91
CONTINUL
Z (12)-H
DT r(0.4*DZ**2)/iaG
QP2-574.07363*(60.
CALL YLARS(SELTEM,
CALL CCW(QP.DELTEM
h-Q
OMEGA-?.*PI/360.
T (12 ) -10
T 12 -- T 0
JTOT -1
M i " 1
ROWU2 )-(A4+A5*T(1
RO*CP-ROW(12)*CP
CALL EUUIL1 (TN ,TE,
IF (MJ.EG.1 )DELTM2 =
FRVEL-*A6*(T112J )**2)*62.4
1-4 .)*A3*(T( 12)-H.)**2
CMA*AV ( I2)*f, *( (H-Z( 12)
4.*T (10) ) /( 1 .5
(9 )-
CP
Z,XKZU,XKZL, NDAY Sl,TM2,T12,T»OTltDZ)
DO V69 1=1 ,12
IF(XKZ(I).LT.XKZL)XKZ(I)-XKZL
IF4XKZII).LT.XKZO)X HZ (D-XKZO
CONTINUF
DO 91 L-2
F 1 =OT
, 11
d (l)*CP*A ( I) )
.*(XKZ(I)
IF (!YEAR.LL.1973)DLLrM2-U.Ci
IF (IYEAR.LL.1973)CP2-0.0
F 31-RO* (I )
F41=(HOU(I)*CP*QP2/(1.5*02))*DELTMi:*(T(I+l)-T
F4r(RO»i(I)*CP*OP(Jw)/(1.5*OZn *DEL 1EH ( JW ) .* t T (
)*OELTEM(JW)
-------
HHIi- IJHIJ. hbl
u r~ r £ * T
< I )Nlr( 1)1 Ihl
21* I- 1 Z6 OU 109 Chi
(2T) 1:?! 1
jnNIlMOD 626
( I )N1=( I) Zl v
31' I- I 6c'6 00
3WllrZ3WIi ITT
n + 3WIlr3WIl 0^1
(V'NillIXlwnV3 631
(z mar si 9i yzt
( TSM03**f )/( 31*1 SNO^+ (OT >N1- I I I ) Ni#* h ) - ( Z\ ) Nl cj T L?.\
91 01 0*) 9?I
»)) = {znNl t»I S?t
31 01 oy h?l
(Ot)Nl=OI3i I? I o
( Til MU 1131 0,3 t
l)r IS NOD 61 t
XX/10SH*r301-3l Bit
0+S'hJr XX ill
S£:*Dztfl3a 911
*0«- 2*6r r J gtl
0* Z/IN3CU* < el) NilrMi hll
( 2)1^( H Ml U I
3fmNOO 16 211
X+£JfZJ)r(I>Ni U I
*I-M1X( S'l 9* U JI UT T
•I=M1X(S*3TIIJI 6HI
i-iwxf zai*3Ti ti 11 ji sni
Or tWX( ?ai* 19M I) I) JI i.[)l
•T rWXtTl'llM T) 1) M SOI
•o-wx
-------
146
117
148
149
153
151
1 52
Ib3
l'->4
155
1 56
157
158
Ib9
100
161
16?
163
1 tl
165
166
167
168
169
170
171
17?
173
1 74
175
1 76
1 77
178
179
IdQ
Ifcl
It.?
U 3
Io4
Ih5
Ic6
Ifc7
It; 8
Ifc9
190
I'-'l
501
5f!2
66
VR8
313
7TQ
IF (NDAYS. LE. 360)T IME2-T I ML 2-360. LI
IF(NOAYS.CE.360)T IME-TIME-360.0
IF(NDAYS.GE.360)T lML4 = TlMt4-360.0
IF(NUAYS.GE.36G)T1ME5-TIML5-360.0
If (NOAYS.LE.360)JJ=0
IF(NDAYS.bE.360)J«rl
IF (I YEAR.GT.1979) bO TO 99
IF(NDAYS.6E.360)IYEAR-IYEAR+1
IF (NDAYS.bE.360)CALL CCW(UP,DELTEM,I YEAR,OT)
IF(NDAYS.GE.360)CALL YEARS(SFLT£M,QQPP,IYLAR)
IF(NDAYS.GE.360)JTOTrjTOT+1
IF(NDAYS.GE.360)J1M-JIM+1
IF(TIM£4.GE.1.0)GC TO 5U1
GO TO b02
MMIrMMI+1
XTUD »JIMtMhI)rTD
TIME4-TIME4-1.
CONTINUE
IF(NDAYS.Gt.36Q)NOAYSro
DO 66 Ir2 , 10
CbU )-(T (2)-T (1) )/7.5
Cbjll»=(Tll2l-TCllll/7.5
IMTIML1.GC.3U.) GO TO 98
GO TO 33
PRINT 988
FORMA] (IX
T1ME4-0.
MM1 -L
JJ-JJ+1
Jw -J« +1
NUAYS1 =TIHE3
M J"M J+l
DELTM2-OM (MJ) -T(5)
IF (MJ. bE. 12 )H J-l
CONTINUE
DO 7CiQ irl ,12
T (1 ) rTN (I )
IF (JRUN. EQ
, (CP < JWJ1 , JU J-l. 12)
,12F10.1)
,2) GO
(T ,AU
TO
,CB
111
tZ, A
XKZ.ROW.TN , DM. T2, MONTHS T2
-------
X
U1
(M
LJ •
>• a.
•-t JE
o:
o •
— • x
_! ^ u
o a './>
-/> — • — » »-<
X * CM Q
• (-«—••
O • —fSJUJ >-
Z X CJ-«UJ _J
I— « (M ••< » • X
-^ 1 1 >• ** z
10
LULU *•
i-< ~ a s o
o
o
H-I
—
<
h-
1
Q.
JE
O
CJ
u.
o
X ~> UO — • H- ^£ <>— CK
•>*7>-(/) »f»-Q< • »U-C\JO
O l^j
>— 33
s:
i-i js^a:i-«»-»"-i-ta:uia:a! JssQiu1^
o «* a: o cc ocr ecu. 0:0 o oo 01-10 i. u. octo ^: t-^
•^ a •— i
« CM a
CM
32
-------
Ul
noonn
m
1 1 t/> Z O
-toa Nf>- aaoa u i/>c
& -J CO •
*
* tn i> tn a o ui rv> •»• m
*• tn j> j-. •# t-
Q rn
(/I l/l
X J>
O X -<
O
33
rn
-4 .2 CD J>
J> 3D —4
r-^o m
r c r:
n*» -*
XJ5s.uO
i rn o i\>
«/>o o
• ct
n—• t-i
O-« 'JT On
2 f> O
tr. J>
• -H-H>
X -I
Hm
v >->
n i> ,^
r~ x>"
IT,
(/I M
Ki
O-sl
T)
-------
1
2
3
14
5
6
7
8
9
10
11
12
13
It
15
16
17
13
19
20
21
22
23
24
25
26
27
2<3
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
43
49
50
SI
5?
S3
b4
55
56
57
5P
59
bO
fcl
fc2
63
64
t5
66
6?
ecu
c
c
c
c
c
c
10
12
SYM CREATED ON 12 AUG SO AT 13:00:09
THIS SUBROUTINE CONTAINS THE CONDENSER
COOLING WATER. ASSUMES THAT COMPUTATIONS
START IN 1971.
SUBROUTINE CCW(QP,UELTEM,I YEA 3,DT)
DIMENSION QP(12) ,D£LTEM(12)
IFCIYEAR.GT.1979)GO TO 11
IYEA=IYEAR-197Q
ACOST=10.0
,5,6,7,8,9),IYEA
3,4
GO T0(l,l ,
DO 10 1=1 ,12
QP (I )=U.O
DELTEM tI)=C.Q
GO TO 11
DO 12 1=1,6
QP(I)=0.0 "
DELTEM (I )=•(,.0
QP (71 = 1890.2*ACOST
QP (8 )=1910.3*ACOST
QP (9 )=2170.7*ACOST
QP (10=2232.5*ACOST
QP (11>=217C.7*ACGS
QP (12)=3284.6*ACOS
DELTEM (7>=5.3
DELTEM(8)=4.b
DELTEM(9)=5.3
DELTEM (10 )=7.3
DELTEM (11 )=7.7
DELTEM(12 )=4.1
GO TO 11
QP (1
QP (2
QP (3
QP (4
T
T
=3069.3*ACOST
=3069.4*ACOST
=2976.9*ACCST
3*ACOST
= 2807
QP(5 =2164.6*ACOST
QP (6
QP(7 =5334.6*ACOST
OP (8 =4727.1*ACOST
QP (9 )=5961.4*ACOST
QP(10)=4953.4*ACOST
QP(ll)=420t.i*ACOST
QP (12)=5225.6*ACGST
DELTEM (1 )=H.2
DELTEM (2 )=7.4
173=8.4
) = a.o
)=2.7
)=6.0
DELTEM
DELTEM
DELTEM
OELTEM
DELTEM
DELTEM(8)=4.8
DELTEM (9)=5.8
DELTEM(1C)=3.
DELTEM(11)=7.
OELTEM (12 )=5.
GO TO il
(4
(5
(6
(7
(8
3P (1
QP (2
QP (3
QP (4
QP £5
QP (6
OP
QP
)=4612
)=3694
)=5456
)=S57G
5=6494
)=6574
(7 )=71Q4
(6 )=7513
QP(9 )=72C1
QP (1C)=6993
4*ACOST
9*ACOST
3*ACOST
6*ACOST
3*ACOST
2*ACOST
2*ACCST
l*ACQST
6*ACOST
.4VACOS7
-------
69
70
71
7?
7?
74
75
76
77
73
79
oO
ol
52
83
«4
bS
'36
67
68
89
va
91
9?
93
94
->5
96
97
99
99
ICQ
1C1
1C2
1C3
1C4
iqs
Iu6
1C7
1C3
1C9
110
111
1 12
113
114
115
116
117
118
119
1Z3
121
122
123
124
125
126
127
123
129
130
1 31
1 J2
1 33
134
135
136
137
1 38
129
QP 111 ) =74&7.1*ACOST
QP (12) =685C.9*ACGST
DELTEM (1 )=b.3
)=4.6
)-c.2
= e.3
J-c.6
}=8.e
(2
DELTEM (3
OELTEM <4
OELTEM
OELTEM
(
OELTEM (7 )=b.3
DELTEM (8 )-7.8
DELTEM (9 ) = 7.4
DELTEM
DELTEM (11 )
DELTEM (12 )
GO TO 11
QP (1 )=6069
OP (2 1=4440
QP 13 )r4P74
OP (4 )=4272
GP (5 J-3970
OP (6 >=S197
OP (7 >=583G
QP (8 )=724e
(10)=7.7
=b.5
=9.4
(4 }-
(5
(6
QP (10 1=563
QP (1 1) =58G
QP (12 ) =491
DELTEM (1jr
DELTEM(2) =
DELTEM
DELTEM
OELTEM
OELTEM
DELTEM(7
DELTEM (S
DELTEM (9
OELTEM (10 )
DELTEM (11 )
DELTEM <12 )
GO TO 11
QP (1 1=5045
QP (2 1=4985
OP (3 1=5113
QP (4 1=6013
QP (5 1=6302
QP (6 )=43B5
QP 17 1=5033
QP (S )=57CS
QP (9 1 r6964
QP(101=675
QP (1 1 1 =469
QP (121=535
OELTEM (1 1 =
DELTEM (2 1 =
DELTEM (31 =
DELTEM (4 1 =
DELTEM (5 }-
OELTEM(6)=
DELTEM (71 =
DELTEM (B ) =
DELTEM (9) =
DELTEM(101
OELTEM(11)
DELTCM(121
GO TO 11
QP (1 1=6176
QP(2 1=6444
QP(3 1=5195
QP (4 1=4811
OP (5 1=4984
OP (6 1=5659
QP (7 1=7053
.3*ACOST
.2#ACOST
.3*ACOST
.1*ACOST
.7*ACOST
.6*ACOST
.C*ACOST
.3*ACCST
.4*ACOST
7.6*ACOS1
9.2*ACQST
4.8*ACOS T
10.6
7.3
7.1
i.l
5.6
9.3
7.4
c.-5
b.O
=7.8
=6.7
=8.4
.3*ACOST
.2*ACOST
.5*ACOST
.6*ACCST
.4*ACOST
.3*ACOST
.6*ACOST
.9*ACOST
.C*ACOST
<+.7*ACOST
7.6*ACOS T
4.6*ACOS T
12.5
11.4
10.4
11. H
9.4
b.4
7.4
5.0
b.O
= 3.fa
-6.2
=7.9
,7*ACOST
.fc*ACOST
.7*ACQST
.8*ACOST
.2*ACOST
,9*ACOST
.b#ACOST
35
-------
140 . QP(8>=7914.9*ACOST
141 QP (9 )=6557.3*ACOST
142 QP (10 )::7407.4#ACOST
143 QP(11)=6065.1*ACOST
144 QP(12}=65a3.5*ACOST
145 OELTEM (1>-9.0
146 DELTEM (2 ) = 11.C
147 DELTEM(3>=13.2
148 OELTEM(4 } = 9.7
149 DELTEM(5)=10.1
150 DELTEM (6>=&.!
151 DELTEM (75 = 7.9
152 DELTEM(8)=7.5
153 DELTEM (9> = 7.6
154 DELTEM(10 )=6.2
155 DELTEM (11 >=S.4
156 OELTEM(12 1-7.2
157 GO TO 11
158 9 QP(1}=7207.7*ACOST
159 QP(2}=7319.9*ACOST
160. QP (3 )=7419.5*ACOST
161 QP(4)=7275.3*ACOST
162 QP(5)=4189.1*ACCST
163 QP(6>=5331 .2*ACOST
164 QP(7)=4733.3*ACGST
165 QP(8)=4733.3*ACOST
166 QP(9)r4733.2*ACOST
167 QP(10)=4733.3*ACOST
168 QP(11)=4733.3*ACOST
169 QP(12)=4733.3*ACOST
170 DELTEM(1)=10.3
171 DELTEM(2 )=10.4
172 DELTEM(3)=9.6
173 DELTEM 14 )=9.9
174 DELTEM I5) = c.2
175 DELTEM(61=7.1
176 DELTEM (7J-S.O
177 OELTEM (8) = 5.0
178 DELTEM(9)=5.C
179 OELTEM <1Q ) =5 . 0
180 DELTEM(11)=5.0
161 DELI EM (12 )=5.C
132 11 RETURN
163 END
36
-------
EQJIL1 SYM CREATED ON 11 JU N 80 AT 11:00:00
1 SUBROUTINE EQUIL1 C TN,TE,XK, TDEW,X TN,XTEf XXK,WIND,HSOL )
2 DIMENSION TN < 20 ) ,X TN( 20 J
3 2 RLADI5 ,1 )TDEW .WIND.HSOL
4 1 FORMATO
5 WINO-WlNO*U.<4b
6 HSOLrHSOL*3.6855
7 TM=(TN (12 I+TOEW)/2.0
8 FW-9.2+0.'*6*(W IND*>>2J
9 BETAro.35+0.015*TM+0.0012*( TM*42»
10 XK=4.5 +O.U5*TNJ12) *B£ T A*FW + 0. M7*F U
co
12 TE-TOEW+HSOL/XK
13 XTH=(XTN(121*TOEU) /2.0
14 XF»i-9.2*0.^6* (WIND**2 )
1^ XBETA=G.35+O.Q15*XTM+0.0012*(XTM**2I
16 XXK-4.5+0.05*XTN<12)*XBETA*XFW+0.q7*XFU
17 XXK-XXK*4 .232* I9./S.J
18 XTE^TDEW+HSOL/XXK
19 RETURN
?0 END
-------
1
2
3
H
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
23
29
7r*
*j •-*
31
32
33
34
35
36
37
38
30
40
41
42
43
44
45
46
47
43
49
50
51
52
53
54
55
56
57
53
59
cO
61
62
63
64
65
66
67
68
69
INPJT S
3.0,6.6
0.,9.3,
6.3, 5.2
7.5,8.7
17.2,7.
18.8,5.
20.,6.4
19.44,5
18.33,5
13.88,7
2.88,7.
5.5, a.3
1.67,6,
-2.22,9
1.11 ,9,
6.67,8.
11.11,7
13.13,7
IB. 7 7,6
22.22,6
18. 8,5.
11.5 ,7.
5.9,7.1
4.,6.8,
1.t7.22
-1.,7.3
1C. , 7.1
7.7,S.4
14.3,6.
20.25,3
22.2,5.
21.7,5.
^w•3 t o *
13.5 ,7.
7.2,8.1
3.2,5.6
S.2,5.8
C . , 5 . 8 ,
6.3,7.7
in.7,8.
17.2 ,6.
17. 8,6.
21.,5.2
2 1 . , 5 . 8
17.5 ,6.
10.2 ,5.
6.C , 7.2
3.8,0.9
3.0,6.3
3.5, 7.6
2.2, 9.6
7.2, 7.6
17.5 ,4.
19.3 ,5.
21.3 ,5.
21.0,5.
16.2,7.
1?.4 ,7.
7.9,6.9
i! • 0 » 7 »2
-1.0,7.
3.2, 8.5
3.9, 7.9
11.2 ,7.
14.G ,7.
1 ?.3 ,6.
19.8,5.
13.G ,6.
15.4 ,7.
YM CREATED ON 12 AJG 80 AT 13:01:13
5,167.0
264.4
3,264 .4
2,457.
5,480.5
65 ,478.
8,409.
.75 ,428.2
.77 ,329.
.02 ,261.3
53,247.7
,147.7
69,178.
.26 ,257.6
20,352.5
72,448.
.53,433.6
.95 ,564.3
.64 ,493.8
.07 ,453.5
47,386.3
17 ,298.1
3,220.9
148.
,162.7
,279.5
,348.5
4,449.3
83,449.5
.04 ,507.7
32,496.9
1 ,391 .6
3C3 ,338.4
1,341 .7
4,247.6
,154.
119m
226.9
,326.1
73 ,397.7
6,43t> .
98,559.3
,459.5
7,480.
74 ,339.2
7,302.5
,231. 1
,181.9
93 ,191.4
14 ,226.9
,326.1
,397.7
8,436.
82 ,559.3
1C,459.5
4 ,480.8
3,339.3
7,3C2.5
,231. 1
,181.9
4,2t79.8
,31G.9
,338.6
6 ,496 .9
3,448 .4
4 T4c3 .2
9,4S3 .3
65 ,48G.4
13,345.1
-------
70
71
72
73
74
75
76
77
78
79
80
81
82
S3
84
85
66
67
88
69
90
91
92
93
94
95
96
97
93
99
100
1C1
1C2
103
104
1C5
106
1C7
1C3
109
6*2 ,"7.2 1,2 8 7. 5
1.0, 7.27,237.5
-1.5 ,8.2,195 .0
-6.6 ,8.04 ,205.5
-2.78,8 .4 ,317.6
6.0, 1.7 ,328.5
JL0.2 ,7.6,427.3
15.4,6.2,473.
18. ,6. 7,543. 3
20.2,5.8,551.8
20.7,5.4,423.9
18.7,5.3,350.7
9.2, 7.2,286.6
7.0, 7.5 ,196.2
0.4 , 7.2 ,178.2
-2.8,7.9,227.
-5.0 ,6.8,3C8.
1.2, 7.6 ,408.
9.6, 7.6,429.
14. ,6.7 ,513.
19.4 ,4.7,593.
23. 3,5.7,568 .
2H. 8,5.1,461 .
15.5 ,5.7,365.
9.3,6.6 ,369.
9.0,5.8 ,232.
0.4, 7.3,191.
-3. 33,8 .6 ,206.
C.Q, 7.2 ,251.
b.O, 7.9,373.
9.2, 7.6 ,479.
14. ,6.7,513.
19.4 ,4.7,598.
23.8,5.7,56d.
2C. 8,5. 1,461 .
15.5 ,5.7,387.
9.3 , 6.6 ,369.
9.0, 5.8 ,232.
u.4 , 7.3 ,191.
u.4, 7.3 ,191.
39
-------
----- MIXIT SYM CREATED ON 12 AUG 80 AT 13:26:57
1 C
2 C
3 C THIS SUBROUTINE MIXES STABILIZES UNSTABLE
4 C TEMPERATURE PROFILES.
5 C
6 C
7 SUBROUTINE MIX IT ( T ,\, A )
8 DIMENSION TN(20),A12Q)
9 100 DO 10 1=1,11
10 IF(TN(I + 1 J.GE.TNUMGO TO 1
11 IF ( (TN II J-TNU+1 J ) .LT .0 .0)60 TO 1
12 T AV=UN(I + 1)+TN( I) )/2.
13 TN (1+1 )=TAi/
1^ TN(I)=TAV
15 I CONTINUE
16 10 CONTINUE
17 TMAX rAMAXl ITN (1) TTN(2} ,TN( 3) , TNj( HJ ,TN { 5)
18 C,TN(6 ) ,TN(7 ) ,TM3> ,TN
19 C(9) ,TN (10 ) ,TN (11 ) , TN( 12) )
20 IF(TN(12).LT.TMAX)CO TO lOu
21 3HO RETJRN
22 END
-------
u
5
6
7
8
9
10
11
1?
13
14
15
16
17
18
19
20
21
22
23
24
2S
26
27
28
29
30
31
3?
33
34
37
38
39
40
'•1
4?
43
44
45
PLOT
C
C
C
C
C
C
100
5
TER SYM CREATED ON 12 AJ6 80 AT 12:56:46
PARAMETER Nri/4fNN-12,NTIME-12,ND~110
DIMENSION IBUFUOQG)
DIMENSION TEMP<50) ,DEEP(50) , TEMPS (ND) ,
DIMENSION T(N) ,AV ,THF (ND» ,TSU (ND)
) , ZE D t ND i
E6I 50),E5«50),
ASURED DATA
NOT
BE MODIFIED
TZ, MONTHS, T2,QP,
CCP ,SIGMA,R3 ,R4 .r..,.
CXK20 ,TL,NDAYS,TN12
CNDAYS1 ,TIML,TIME2,
I COUNT-ICOCNT *1
IF (ICOUNT.bT.96)60
IF (J(xjN.EQ.2) GO TO
RE AD<5,8) (DEEP(INK
DO IS. KL-1 ,50
RLAD(5,fl) uEEP(KL)
CEC (KL » ,E7 (hL)
RE AD (b ,6 ) Afil ,BE1 ,C
DEEP (KH-
R6fR7.R8.R9. Rl J. uP 2 , FRE VEL , ROWCP f DT ,
,T12,Fl,F2fF3,F31,F41,F5,F6,F7,F8,TD,TD2,
TIME3,IYEAR,MJ,XK,TDD,J)
TO 333
2UO
J, TFMP( I NK ) , INK- l.NSTOP)
,El(KL),E2(KL),E3(KL),E4(KL)tE5(KL),
E1,DE1,£E1,FE1,LE1,HE1,OE1
E 1 (KL ]
E2 (KL
E3 (KL
E4 (KL
E b (KL
E6 (KL
E7 (KL
LJ U, LJ LJ UJ L^ UJ
1 1 1 1 1 1 II i t II II
-------
46 EDIKL)=OE1
47 IF(E3(KL).EQ.0.0)60 TO 16
48 IF (UEEPtKL I.LQ.t-l.llGO TO 16
<*9 TEMP(KL>-(EKKL)+E2(KL)+£3(KL>*E4
b? 15 CONTINUE
53 200 CONTINUE
54 16 NS10P=KL-1
55 IF (JRUN.EQ.2JGO TO 201
56 DO 222 JU-1,50
57 IF 1DEEPIKL >.EQ.(-l.nGO TO 223
58 READ(5,8)ALl,B£l,CEl,DEl,EEl,FEl,GEl,HElfOEl
59 IF (AE1 .EQ. (-1 ) ) 60 TO 223
60 2?2 CONTINUE
61 223 CONTINUE
62 201 CONTINUE
63 CONS2-1./0.3048
64 IFURUN.EQ.2)GO TO 202
65 00 9 INK-1 ,NSTOP
t6 DEEP(INK)=CONS2*DEEPri50.-DEEP( INK)
6ft DEEP (NSTOP+1)-0.0
69 DEEP (NSTOP+2)=2INN)/1.5
70 TEMP (NSTOP+1)rQ.O
jg 71 TEMP (NSTOP*2) -30.0 /I.5
7? 202 CONTINUE
73 6 FORMAT ( )
74 333 JOzJO+1
75 L-L + 1
/6 TSJ(L)^T<12)
77 X 2fjiXZL)43n.
7fi ZED(L)rX?D
79 . TEMPS (L ) = TEMP (II
bO TCQ(L)=TE
ei THF (L ) - (T (7) + T(6)) /2.
6? IBCDrMONTHS(JO)
fa 3 ' 7 (KN+1 )-0.0 "
b4 2 (NN+2 )-Z (NN 1/1.5
85 T (NN+1 UO.U
^6 T JNN+2 )-30./l.5
87 CALL AXIS(0.0,0.0,tHTEMP.(C»,-6,i.b,U.D,TI131,T(14il
fc8 CALL AXIS(0.0,0.0,9HDEPTH(FT),9,1 .5,90.0,7(131 ,Z(14))
89 CALL FLINE (T ,Z ,-NN,1,0,0)
'jD IF tlCOUNT. GT.Vt. ) GO TO 444
^l IF (JNUN.EQ.2)GO TO 203
9? CALL DASHL (TEMP ,DE L'P , NST OP , 1)
93 203 CONTINUE
94 444 CALL SYMBOL(1.0,0.5,0.14,I BCD,O.U,6)
•;5 CALL PLOT (2.25 ,0.0 ,-3)
-------
96
97
98
99
ICO
lul
10?
103
104
105
106
1C?
108
109
110
111
112
in
114
115
116
1 17
118
119
120
121
122
123
124
12r>
126
127
128
129
1 J-0
131
132
1 33
134
135
1 36
137
138
139
3
1
6
13
IF (JO.E0.4.OR.JO.Eu.8)GO TO 3
GO TO 1
CALL PLOT 1-9.0 ,-2.25,-3)
CONTINUE
CALL PLOT (-2.25 ,0. ti,-3)
CALL SYMBOL(-6.75,b.75,.14,41HTEMPERATLRE PROFILES FOR LAKE KEOUEE
C ,G.O ,41)
P1-JYEAR
MY^JYEAR
CALL NUMBER(999. ,999.,0.14,P1,O.U,U)
CALL SYMBOL<-6.75,6.5,0.1,54H(DEHTH IS MEASURED FROM THE DEEPEST P
COINT OF THE LAKE).0.0,54)
CALL PLOT IB.0,-9.25,-3)
PRINT 2,MY
FORMAT (IX, 'THE PLOTS FOR *,I 5,* ARE COMPLETE')
IF(M.EQ.9)GO TO 6
JYEARrJYEAR+1
GO TO 5
CALL PLOT (6.0 ,0.0,-3)
DO 13 1=1,96
DEEPSU )=ZLD(I)
DEEPS(97)=U.O
DEEPS 08 )r;j240.0/9 .0
TSUU09
TSUdlO
TEfc(lC9 rn.O
TEGdlCl
THF (109
-O.G
-3b./
l-U.O
:n4
THF(11U
TEMPS(97
TEMPS(98 J-35./5
-O.C
ZEDdlO r3240./'
CALL PLOTIC.0,2,
CALL AXIS (li.O ,0,
CALL AXIS(U.O,0
C(110 ) )
CALL FLINE (ZEU,TSJ,-108,1,2,11)
CALL FLINE (ZED,TEQ ,108,1,2,2)
CALL FLINE (ZED ,THF,-108,1,2,0)
If (JRUN.EQ.2)GO TO 204
CALL UASHL (DEEPS ,TEMPS,96, 1)
CONTINUE
0,-3)
0,12HTIME IN DAYS,-12,9.0,0.0,ZED<109 ) ,ZED(110))
0,16HTEMPERATURES (C),16,5.0,9L.,TSU(109),TSU
-------
140 t
mi c
142 C
143 C CHANGE TITLES TO SUIT NEEDS ( 4 LINES)
144 C
1**5 C
146 CALL SYMBOL (0.0,6. 11,0.1*4,
147 C46HSTR ATIF1CAT ION CYCLE FOR LAKE K£OkEE I 9 71 -1 979 , 0.0 ,46 )
148 CALL SYMBOL (0.0 ,0. 10,0 .10,87H 1971 1972 1973 197**
149 C 1975 1976 1977 1978 1979,0.0,871
150 WRITE(b,7)
1B1 7 FORMAT (IX ,'ALL PLOTS ARE NOW COMPLLTE«,//,» NOKMAL JOB EXIT*)
152 CALL PLOT (15.0 ,0.0,-3)
153 STOP
154 EfMb
-------
READER SYM CREATED ON 12 AU6 81) AT 13:21:45
1 C
2 C
3 C THIS SUBROUTINE READS THE MAGNETIC TAPE
4 C CONTAINING THE COMPUTED RESULTS.
5 C
6 C
7 SUBROUTINE READER < T,AV.CB i 2 tA.XKii, KOU, TN.DM.TZ, MONTHS ,T2, QP,
8 CCPlSlGMA,R3,RH,R5tR6,R71R6tR9tRllJ,CP2,FREVEL,ROWCP»DT,
9 CXHZO,TE,NDAYS.TN12,T12tFl,F2tF3TF3l,F41tF5fF6,F7,F6,TOfTD2,
10 CNDAYSl,TIME,TlME2,TIME3tIYEAR,HJfXK,TDD,JtNCASL,SF,EDEPT,\/OL)
11 DIMENSION T (20 ) ,AV (20) ,CB(20),Z(20),A( 20),XK2(20),
1? CROh(20)1TN(20),DM(20)tTZ(2U)|T2<2u>tUPU2)
13 CHARACTER*6 MONTHSC12)
14 1 CONTINUE
15 READ (b,ENU=l) IT( I Jl,IJ-1, 12),(A V(IJ) ,IJ= 1 ,12 ) ,
16 C(Cb(lJ)fIJzl,12).(^(IJ),IJ-l,12),(AUJ),IJ-l,12),
17 C(XKZ(IJ)flJ=Itl2lt(ROWlIJ)fIJ=1.12)vCTN(IJ)(IJ=l,lk)t
18 C(DM(lj),IJ = lt12)t(lZ(IJ)flJ-l,12>'t(MONlHS(IJ),lJ = l,12)t
19 CU2(IJ),IJ-1,12),
20 C(QP(IJ ) ,IJ = 1 ,12) ,
21 CCf,SIGHA,R3,R4.R5,ft6,R7,R8,R9,R1U,QP2,FREVLLtROWCP,DT,
i2 CXKZO,TE,NDAYS,fM2tTl2,FlfF2,F3,F31,F'4lfF5,F61F7fF8,TD,TU2,
23 CNIDAYS1,TIML,TIME2,TIME3,IYEAR,MJ,XK,TDD,J,NCASE,SF,EDEPT,VOL
21 RETURN
25 END
-------
SMOOTH SYM CREATED ON 12 AUG 80 AT 14:34:30
1 C
2 C THIS SUBROUTINE CORRECTS THE EDDY ulFFLSlVJTY
* L IF VARIABLE TIME STEP IS REQUIRED, »DU* SHOULD
4 C BE CHANGED 10 'DT* IN THL CALLING PROGRAM.
5 C
6 SUBROUTINE SHOOT HIXKZ,XKZU,XKZL,NuAYS1,TN12,TI 2,T,DT1 , DZ)
7 DIMENSION XKZ(20 ) , 1 (20)
8 DO 93 1=1,12
9 IF (XKZ (I ) . CT.XK2U) XKZ(I)=XKZU
10 IF (XKZ (I KLT.XKZL) XKZ(I)-XKZL
1] 93 CONTINUE
1? NEkrG
13 DO 96 1-2 ,12
14 IF (XKZ (I ). LQ.XK2L) NE*l-I
15 96 CONTINUE
16 IF tNEw .EO.O )GO TO 77
17 DO 55 1-1 ,NEW
18 XKZ(I)-XKZL
19 55 CONTINUE
20 77 CONTINUE
21 IF (NDAYSl.LE.6Q.OR.NDAYS1.&T.3QO) bO TO 29
22 IF (TN12.6E.T12) GO TO 19
23 IF UN12.LT.T12) GO TO 39
2<4 19 XHIN-AMIN1 (XKZ(l),XKZ(2),XKZ(3),XKZ(4),XKZ(5),XKZ(b),XKZ(7),XKZ
25 1 (B ),XKZ <9 > ,XKZ (1Q> ,XKZ(1 li ,XK7( 1^) )
26 DO 82 1=1 ,12
27 IF (XKZ (I ).LQ.XMINJ GO TO 61
28 82 CONTINUE
29 GO TO 29
JO fil I PIN-I
31 DO 7U 1=1 ,1MIN
32 XKZ (I ) =XKZ (IM1NI
33 7G CONT1.NJE
34 GO TO 29
3^ 39 XMAX=AMAXl(XKZ),XKZ(7),XKZ
36 l(h>,X*\Z«9),XKZ(10>»XKZ(ll),XKZ(12M
37 DO 62 1=1 .12
38 IF (XKZ ( I l.EO.XMAX) GO 10 61
39 62 CONTINUE
in GO TO 29
41 61 iMAXrl
42 DO 5G 1=IMAX,12
43 XKZ (II =XKZ UMAX)
44 *0 CUNT INjE
45 29 CONTINUE
46 200 XMAX-AHAX1(XKZ(1),XKZ(2),XKZ(3),XKZ(4J,XKZ(5),XKZ(6),XKZ(7),XKZ
47 1 (h ) ,XKZ (9 ) ,XKZ (10) ,XKZ tllll ,XKZ (12) )
48 DTI:(U.4*DZ**2J/XMAX
49 RLTJi
-------
STORE SYM CREATED ON 12 AJG faO AT 13:19:47 .
1 C
2 C
3 C THIS SUBROUTINE STORES THE COMPUTED RESULTS ON
4 C MAGNETIC TAPE.
5 C
6 C
7 SUBROUTINE STORE(T,AW,CB,I,A,XKZ,KUW,TN,DM,TZ,KONTHS ,T2,QP ,
8 CCP,SIbMA,R3,R4,R5,fl&,R7,R8,R9,R10t&P2,FREVEL,RUWCP,DT,
9 CXKiOtTEtNDAYStTN12,T121FlfF2.F3,F.$l|F,TN(20),DM(20),fZ(2u),T2(2G),
13 CQP(12)
11 CHARACTERS MONTHS(12)
15
16
17
18 C(DM(IJ)tiJ=lfi2).(lZ(IJ).IJ=l(12)|(MONrHS(IJ)tlJ=ltl2)
19
2H C(QP(IJ ) iTJ^l .12) .
21 CCPfSlGMAfR3,RH,RSfR6,R71Rfa,R9,RlUfwP2.FREVLLfRGWCPtDT,
22 CXK20,TE,NDAYS,TN12,Tl2,Fl,F2,F3,FJl,F41,F5,F6,f7,Fa,TD,TD2,
23 CNDAYSl,TlME,TlME2,lIME3,IYEAR,MJ,XK,TDD,J,lJCAStfSF,EDEPTtVOL
24 END FILE 8
2^ RETURN
26 END
COP (12 )
CHARACTERS MONTHS (12)
WRITE (8) a(IJ) »Io=l,12) , (Atf( IJ) ,IJr 1, 12) ,
C(X»sZaj!lIJ-l,12)f,(ROW(lJ)TlJ-l,i2),(TNlljT,ij-lIl2)l
C(DM(IJ),IJ-1,12) ,(lZ(IJIfIJ=lt 12) t(HON FHSCIJ)t U = li!2
C(T2(IJ ) ,1J-l,12) ,
-------
YEARS SYM CREATED ON 12 AUG SO AT 13:1U:C3
1 C
2 C THIS SUBROUTINE PRINTS THE YEAR TITLE.
3 C
H C
5 SUBROUTINE YEARS(S£LTEM,QQPP,I YEAR>
6 PRINT 99.IYEAR
7 99 FORM AT t59X ,17<•*•> ,/,59X,**',15X,•*'»/,59X,
8 C' ** ,2X , 'YEAR = • ,I4,2X ,»*',/, 5 9X , f **, 1 5X, « *•
9 C,/,59X ,17( •*•) )
10 RETURN
11 END
-------
TECHNICAL REPORT DATA
(Please rtod laaructtunt on the revene before completing)
t REPORT NO.
EPA-600/7-82-037f
2.
3. RECIPIENT'S ACCESSION NO.
4 TITLE AND SUBTITLE
Verification and Transfer of
Thermal Pollution Model; Volume VI. User's
Manual for One-dimensional Numerical Model
5 REPORT DATE
May 1982
6. PERFORMING ORGANIZATION CODE
7 AUTHORIS)
S.S.Lee, S.Sengupta, and E.V.Nwadike
B PERFORMING ORGANIZATION REPORT NO.
9 PERFORMING ORGANIZATION NAME AND ADDRESS
The University of Miami
Department of Mechanical Engineering
P.O. Box 248294
Coral Gables, Florida 33124
10. PROGRAM ELEMENT NO.
11. CONTRACT/GRANT NO.
EPA IAG-78-DX-0166*
12 SPONSORING AGENCV NAME AND ADDRESS
EPA, Office of Research and Development
Industrial Environmental Research Laboratory
Research Triangle Park, NC 27711
13. TYPE OF REPORT AND PERIOD COVERED
FinaIt 3/78-9/80
14. SPONSORING AGENCY CODE
EPA/600/13
15 SUPPLEMENTARY NOTES IERL-RTP project officer is Theodore G.Brna, Mail Drop
61, 919/541-2633. (*) IAG with NASA, Kennedy Space Center, FL 32899,
subcontracted to U. of Miami under NASA Contract NAS 10-9410.
16. ABSTRACT The six-volume report: describes the theory ofa tnree-aimen-
sional (3-D) mathematical thermal discharge model and a related one-
dimensional (1-D) model, includes model verification at two sites, and
provides a separate user's manual for each model. The 3-D model has two
forms: free surface and rigid lid. The former, verified at Anclote An-
chorage (FL), allows a free air/water interface and is suited for signi
ficant surface wave heights compared to mean water depth; e.g., estu-
aries and coastal regions. The latter, verified at Lake Keowee (SC), is
suited for small surface wave heights compared to depth (e.g., natural
or man-made inland lakes) because surface elevation has been removed as
a parameter. These models allow computation of time-dependent velocity
and temperature fields for given initial conditions and time-varying
boundary conditions. The free-surface model also provides surface
height variations with time. The 1-D model is considerably more econo-
mical to run but does not provide the detailed prediction of thermal
plume behavior of the 3-D models. The 1-D model assumes horizontal
homogeneity, but includes area-change and several surface-mechanism
effects.
17.
KEY WORDS AND DOCUMENT ANALYSIS
DESCRIPTORS
b.IDENTIFIERS/OPEN ENDED TERMS
c. COSATi I it-Id Croup
Pollution
Thermal Diffusivity
Mathematical Models
Estuaries
Lakes
Plumes
Pollution Control
Stationary Sources
13B
20M
12A
08H,08J
21B
13 DISTRIBUTION STATEMENT
Release to Public
19 SECURITY CLASS (Tras Report)
Unclassified
21 NO OF PAGES
56
30 SECURITY CLASS /Thispage;
Unclassified
22 PRICE
EPA Form 2220-1 (V73)
49
------- |