r/EPA
               United States
               Environmental Protection
               Agency
                     EPA-600/7-82-03 7a

                     May 1982
Research and
Development
               VERIFICATDNAND
               THERMAL POLLUTION MODEL
               Volume L Verificati
              TRANSFER OF
             on of
               Three -dimensional Free-surf ace Model
               Prepared for
              Office of Water and Waste Management

              EPA REGDNS 1-10
               Prepared  lj>y

               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-037a

                                             May 1982
           VERIFICATION AND TRANSFER
           OF THERMAL POLLUTION MODEL
  VOLUME I:  VERIFICATION OF THREE-DIMENSIONAL
    FREE-SURFACE MODEL
                       By
        Samual  S. Lee, Subrata Sengupta,
           S. Y. Tuann and C. R. Lee
      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
     This report is one in a series on the University of Miami's thermal
pollution models.  Much of the background,  formulations,  solutions tech-
niques and applications of these  models has  been summarized in a three-
volume report by  Lee and Sengupta (1978).  These mathematical models
were developed  by the Thermal Pollution  Group at the University of Miami
and were funded by the National Aeronautics and  Space Administration.
The  primary  aim was to have a package of mathematical models which has
general application in predicting  the thermal distribution of once-through,
power  plant heated discharge to  the aquatic ecosystem,  The joint  effort
was planned  so  that the calibration  and verification of these models de-
pend on  simultaneous remote sensing and ground-truth data acquisition
support.  The concept is the development of adequately calibrated  and
verified  models  for direct prediction of thermal dispersion by  the user
communities.   The intended user communities include the utility companies
and the  regulatory agencies at the  federal and state level.

     The purpose of the  present effort is to further verify these models
at widely different sites using minimal calibration, and  then to provide
the program  code and user's manual to the Environmental Protection
Agency  (EPA) for  future users.

     Two sites chosen were Anclote Anchorage on the west  coast of
Florida and Lake Keowee  in South Carolina.   The  free-surface model was
applied to Anclote  Anchorage and the rigid-lid to  Lake Keowee.  Two
data acquisition trips, one in the summer and the other in  the winter,
were carried out at each  site. The acquisition was a collaborated  effort
jointly by the UM, EPA,  NASA-KSC (Kennedy Space Center)  and the
corresponding utility company personnels.

     The two-year project consisted of two phases.  During Phase I,
the individual model was  modified and calibrated to fit  the corresponding
site; then the model was  verified against the remote sensing isotherms
and  in-situ measured velocity  and temperature data  under both summer
and  winter conditions.  In  general, the computed  isotherms were compar-
able with isotherms based on remote sensing.   In  Phase II, source pro-
grams  of both models were documented and  transferred to the EPA, and
the user's guides were prepared to familiarize the potential users.   The
results of the two-year effort are summarized in a set  of three final re-
ports.
                                   ii

-------
                              ABSTRACT
     To assess the environmental impact of waste heat disposal  by power
plant operation into  natural water bodies, mathematical models are essen-
tial, especially for predictive studies.  The Thermal Pollution Group at
the University of Miami has developed a package of models for this pur-
pose.  A joint effort with EPA, NASA, Duke  Power Company and  Florida
Power Company was  conducted to verify these models with remote sensed
IR data and in-situ  measurements.

     The free-surface  model, presented in this volume, is for tidal
estuaries and  coastal regions where ambient tidal  forces play an important
role in the dispersal of heated water.  The model is time  dependent, three
dimensional and can  handle irregular bottom  topography.  The vertical
stretching coordinate is adopted  for better treatment of kinematic condition
at the water surface.  The results include surface elevation, velocity and
temperature.

     The mode! has  been verified at  the Anclote Anchorage site of Florida
Power Company.  Two data bases at  four tidal  stages for winter and summer
conditions were used to verify the model. Differences between  measured and
predicted  temperatures are on an average of  less than 1°C.

-------
                             CONTENTS
Preface  	     "
Abstract	     l»
Figures	     v
Tables	     !x
Symbols   	*	     *
Acknowledgments   	     Xl

     1.  1 ntroduction  	     1
              Background  	     1
              University of Miami  models  	     2
              Description of Anclote Anchorage  	     *
     2.  Conclusions	     *>
     3.  Recommendations	     8
     H.  Mathematical Formulation  and Model Description  	     9
              General  background on  free-surface (tidal)  models .     9
              Governing equations and boundary conditions ......     12
              Uncoupled system	     17
              Computational  grid	     18
              Finite difference equations	     20
              Solution procedures  	     25
              Stability   	     27
     5.  Application to Anclote  Anchorage	     29
              I ntroduction	     29
              Choice of domain  and grid system	     30
              Summary of data	     30
              Calculation of  input	     32
              Results  	     38

References  	•	     43
                                    iv

-------
                               FIGURES


Number

   1   Anclote Anchorage location in the state of Florida	   52

   2   Map of Anclote Anchorage	   53

   3   Definition sketch of  o-coordinate	   54

   4   Grid arrangement in the horizontal projection  	   55

   5   Four cells in a vertical column with velocities  shown at
       definition point and  scalar variables at the center of cell ..   56

   6   Notations and variables used in calculations	   57
                                                                   58
   7   Grid work for the Anclote Anchorage  	

   8   Location  of stations  for in-situ measurement, June 1978 ...   59

   9   Velocity  from  in-situ measurement at 1710-1903,
       June 19,   1978  	   60

  10   Velocity  from  in-situ measurement at 0648-0812,
       June 20,   1978	   61

  11   Velocity  from  in-situ measurement at 1125-1245,
       June 20,   1978  	   62

  12   Velocity  from  in-situ measurement at 1450-1605,
       June 20,   1978  	   63

  13   Daytime  flight lines  on June 19 and 20,  1978  	   64

  14   Surface temperature in deg C from  in-situ measurement
       at 1710-1903,  June  19, 1978  	   65

  15   Surface temperature in deg C from  in-situ measurement
       at 0648-0812,  June  20, 1978  	   66

  16   Surface temperature in deg C from  in-situ measurement
       at 1125-1245,  June  20, 1978  	   67
                                    v

-------
                               FIGURES
Number                                                         Page

  17   Surface temperature in deg C from in-situ measurement
       at  1450-1605, June 20, 1978  	    68

  18   Location of stations for  in-situ measurement, January  1979     69

  19   Velocity from in-situ measurement at 1020-1340, January
       30,  1979	    70

  20   Velocity from in-situ measurement at 1440-1800, January
       30,  1979		    71

  21   Velocity from in-situ measurement at 1430-1640, February
       1,  1979	    72

  22   Surface temperature in deg C from in-situ measurement
       at  1020-1340, January 30, 1979	    73

  23   Surface temperature in deg C from in-situ measurement
       at  1440-1640, January 30, 1979	    74

  24   Surface temperature in deg C from in-situ measurement
       at  1430-1640, February  1, 1979	    75

  25   Surface temperature from in-situ measurement  at 1020-
       1340,  January 30, 1979   	    76
     *
  26   Surface temperature from in-situ measurement  at 1440-
       1800,  January 30, 1979   	    77

  27   Surface temperature from in-situ measurement  at 1430-
       1640,  February 1, 1979  	    78

  28   Semidiurnal  tide for June 19-20,  1978 at south end of
       AncJote Key  	    79

  29   Semidiurnal  tide for January  30-Fefaruary 1, 1979 at
       south end of Anclote  Key  	    80

  30   Surface velocity by modeling at 1030, June 20, 1978  	    81

  31   UW velocity by modeling at 1030,  June  20,  1978  	    82

  32   VW velocity by modeling at 1030,  June  20,  1978  	    83

  33   Temperature from IR  at  1130, June 20, 1978   	    84
                                   VI

-------
                              FIGURES

Number
  34   Surface temperature by modeling at 1030, June 20, 1978  ..   85
  35   Surface Velocity by modeling at 1430, June  20, 1978  	   86
  36   UW velocity by modeling at  1430,  June 20, 1978  	   87
  37   VW velocity by modeling at 1430,  June 20, 1978  	   88
  38   Temperature from IR at  1500,  June 20,  1978	   89
  39   Surface Temperature by modeling at 1430, June 20,  1978  ..   90
  40   Surface velocity by modeling at 1730, June  20, 1978  	   91
  41   UW velocity by modeling at 1730,  June 20, 1978  	   92
  42   VW velocity by modeling at 1730,  June 20, 1978  	   93
  43   Temperature from IR at  1730,  June 20,  1978 	   94
  44   Surface temperature by modeling at 1730, June 20, 1978  ..   95
  45   Surface velocity by modeling at 2030, June  20, 1978  	   96
  46   UW velocity by modeling at  2030,  June 20, 1978  	   97
  47   VW velocity by modeling at  2030,  June 20, 1978  	   98
  48   Temperature from IR at  2000,  June 20,  1978 	   99
  49   Surface temperature by modeling at 2030, June 20, 1978  ..  100
  50   Surface velocity by modeling at 1100, January  30, 1979  ...  101
  51   UW velocity by modeling at  1100,  January 30, 1979  	...  102
  52   VW velocity by modeling at  1100,  January 30, 1979  	  103
  53   Temperature from IR at  1130,  January 30, 1979  	  104
  54   Surface temperature by modeling at 1100, January 30, 1979  105
  55   Surface velocity by modeling at 1600, January  30, 1979   ..  106
  56   UW velocity by modeling at  1600,  January 30, 1979  	  107

                                   vii

-------
                              FIGURES






Number                                                        Page
	  ••••••.                                                        ••••M*...!..



  57   VW velocity by modeling at 1600, January  30, 1979  	  108




  58   Temperature from IR at  1700, January 30, 1979  	  109




  59   Surface temperature by modeling at 1600, January  30,  1979  110




  60   Temperature from IR at  1600, February  1, 1979  	  111
                                   vtl!

-------
                              TABLES





Number                                                       Page



   1    Climatic Data for Summer Run at Ancfote Anchorage  	   46



   2    Climatic Data for Winter Run  at Anclote Anchorage  	   49
                                   IX

-------
                              SYMBOLS
Ht
B
H
K
K
L
P
t
Vertical_eddy viscosity,
cm2 sec
HorlzontaJ eddy djffusivity,
cm2 sec
Vertical_eddy diffusivity,
cm2 sec             ,
Coriolis factor,  sec
Relative humidity  in fraction
of unit
Acceleration  of gravity,
cm sec
Local  water depth with re-
spect to  mean sea level, cm
Total  water depth, cm
Node  index  in the direction
of the x-axis
Node  index  in the direction
of the y-axis
Node  index  in the direction
of the z-axis
Surface heat exchange co-
efficient,  BTU ft"Z day
deg
Reference length, cm
Pressure, dynes cm
Time,  sec
                                        ave
                                       u
T     Water temperature, deg C
T.,    Air temperature, deg F
      Average of air and dew-
      point temperatures,  deg F
      Dewpoint temperature,
      deg  F
      Equilibrium temperature,
      deg  F
      Ambient surface tempera-
      ture, deg F
      Component of water  velo-  ,
      city  along x-axis, cm sec
      Wind speed,  mph
      Component of water  velo-  ,
      city  along y-axis, cm sec
      Component of water  velo-«
      city  along z-axis, cm sec
n     Displacement of the  free
      surface with  respect to the
      mean water level, cm _ -
p     Water density, gm cm
fl     Nondimensional vertical
      fluid velocity
cr     Nondimensional vertical
      coordinate
                                       w

-------
                        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.
Albert W. Morneault from Florida Power Company (FPC), Tarpon Springs,
and his  data collection  group for data acquisition.  The support of Mr.
Charles  H.  Kaplan of EPA was extremely  helpful  in the planning and  re-
viewing  of this project.
                                  XI

-------
                              SECTION  1

                            INTRODUCTION
BACKGROUND

     The  problem of disposing of the waste heat produced as a result of
generation of electrical energy, whether by fossil-fuel or by nuclear-fuel,
is a  dominant consideration in making power production  compatible with
ecological concerns.  For every unit of energy converted to electricity,
approximately two units are rejected in the form of waste heat.  The ul-
timate  heat  sink  for the heat removed  by the condenser-cooling-water
system is the earth's atmosphere.  The cooling water taken  from  natural
or man-made water bodies is  circulated once through the condenser,  and
the heated water is discharged back to the same bodies, which usually
are lakes, rivers, estuaries,  or coastal waters.   Eventually,  the heat is
transferred to the atmosphere through evaporation, radiation and conduc-
tion  over relatively large areas at the air-water interface.

     The  use of  natural water bodies as an intermediate means  for dispo-
sal of waste heat must  take into account the effect upon the environment
of the  circulation and  temperature rises  produced in the receiving water.
The  rate  of oxygen consumption by aquatic species increases with rising
water temperature; however,  the  ability of water to hold dissolved oxygen
decreases with rising water temperature.   There are possibilities of  im-
pairment of biological functions of fish and of breaking  important links
in the  food  chain.  The lethal effects of thermal pollution are sometimes
obvious;  the sublethal  effects on  hydrobiological systems and waste  assi-
milation capacities are not easy to foresee  unless interactive hydrothermal,
chemical  and biological  studies are conducted in an integrated fashion.

     Accurate understanding  of hydrothermal  behavior of the receiving
water bodies is an important  factor in  a power plant system for the  fol-
lowing  reasons:

I.  To  provide a priori information about the nature and extent of ther-
    mal impact on the  aquatic life forms.

2.  To analyse the circulation pattern  of the receiving water body so
    that  recirculation between intake and outlet and consequent decrease
    in  cooling efficiency can  be minimized.

3.  To assess the thermal impact  on the aquatic life forms existing in the
    receiving ecosystem so  that post-operational remedies can be done to

                                   1

-------
    reduce the hazards.

     It is therefore apparent  that not only environmental but planning
and designing interests also are at  stake.

     The above-mentioned objectives can only be met by having  large
data sets over the entire discharge flowfield  under various extreme hydro-
logical and meteorological conditions.  Measurement for temperature and
velocity made over the affected domain could be used  to develop maps  for
velocity and  temperature distributions.  These in-situ measurements can
serve for diagnostic and  monitoring purposes under limited circumstances;
however,  they are not relevant for predictive objective.  The physical
modeling is useful for hydrodynamicai behavior studies when  properly
verified; it is generally expensive and time-consuming.  Frequently,
physical  dimensions  necessitate distortion in  the model, making exact
dynamic  similitude impossible.   Mathematical modeling is tractable and
predictive; it has become a powerful means in simulating complex environ-
mental flows.

     In order to establish the rationale of model  formulation,  the physical
mechanisms underlying the heat dispersion from a heated discharge need
to be outlined.  The following mechanisms govern the heat dispersal.

I.   Entrainment of ambient fluid into the  thermal discharge.

2.  Buoyant spreading of discharged heated  effluent.

3.  Diffusion by ambient turbulence.

4.  Interaction with ambient currents.

5.  Heat loss to  the atmosphere through air-water interface.

The first four mechanisms  redistribute heat  and  momentum in the  re-
ceiving water body.   The last mechanism eventually transfers heat to  the
atmosphere.   It has been customary, therefore, to make assumptions and
approximations which  enable the model solvable.   For  example, the ambient
turbulence is considered by assuming the eddy viscosity and  diffusivity
dependent on the mean  velocity field.  There assumptions necessitate
careful calibration of  models to assure reliability.  Large data bases are
needed for proper calibration,  especially  for three-dimensional,  time-de-
pendent models.   Remote sensing is the  only available method of obtaining
large synoptic data  bases.  Sengupta et al.  (1975) has discussed the need
for remotely sensed data for adequate development of time-dependent,
hydrothermal models.

UNIVERSITY OF MIAMJ MODELS

      Critical evaluation of mathematical  models used for thermal  pollution
analysis has been made by Dunn et al.  (1975).   They  compared  the per-

-------
formance of various models in  predicting a standard data base.  A gene-
ra!  conclusion that can be made from their analysis  is that, though some
models may perform well under certain conditions, a generalized three-
dimensional mode! which accounts for wind, current, tide,  bottom  topo-
graphy and diverse meteorological conditions is yet  to be developed.
One of the first  three-dimensional models was by Waldrop and Farmer
(1973,  1974 a, b).  This model was essentially a free-surface formula-
tion.   One of the first three-dimensional models which adequately  accounts
for bottom topography and comprehensive meteorological conditions was a
rigid-lid model developed by Sengupta and Lick  (1974 a, 1976).  They
used a vertical stretching to convert the variable depth, a concept cus-
tomarily adopted  by numerical  weather forecasting peoples.

     The  thermal pollution  research  team at the University  of Miami has
for the past several years  been developing a  package of three-dimensional
mathematical  models which could have general application to problems  of
power  plant heated discharge  to the aquatic ecosystem.  The primary
motivation behind the  effort was to develop a series of models with mini-
ma! restrictive assumptions, enabling applications to diverse basin confi-
gurations and to various  driving forces of ambient flow.  The effort is
closely integrated with simultaneous  remote sensing  and ground-truth
data acquisition support.  Our aim is to develop adequately calibrated
and verified  models for production purpose;  that is, for direct applica-
tion by the user communities.   The  user communities are the power in-
dustries and the regulatory agencies like the  Environmental Protection
Agency and the Nuclear Regulatory  Commission.

     For the time-dependent,  free-surface mode! which is the main concern
of this volume, the governing  equations are the incompressible Navier-
Stokes equations, conservation of  mass,  energy and an  equation of state.
The conventional hydrostatic and Boussinesq approximations are made.
The mean velocity  field closure is  assumed in defining eddy viscosity
coefficient.   The boundary conditions at water-land  interfaces are no-
slip, no-normal velocity and adiabatic conditions.  At  the air-water inter-
face, wind stress and  heat transfer  coefficients are  specified.  At open-
to-sea  boundaries,  conditions  are specified for surface elevation and tem-
perature.   Otherwise,  normal  derivatives of temperature and velocity are
equal to  zero.  Initial  condition  is assumed to be equilibrium; that is
"cool start."  Conditions at intake and outlet  are completely specified  in
space and  time.

     The  features of the UM's  free-surface model can  be summarized as
follows:

1.   Gravity, earth's rotation,  nonlinear inertia and  bottom  friction terms
    are included  in the hydromechanic part.

2.   The driving forces include wind, tide, river outflow and power plant
    intake and discharge.

3.   Convection, diffusion and  heat transport at water surface are

                                    3

-------
    included in  the thermodynamic part.

H.  The a system of coordinate, devised  by a meteorologist for avoiding
    difficulties at the  free surface, Js  incorporated.

5.  It  predicts spatial and temporal variation of surface elevation due
    to  tide surge and  wind set-up.

6.  It  has means for graphically representing velocity and temperature
    fields.

     This free-surface model is general enough to be applied with
minor modifications to  a  large variety of  sites,  such as lakes,  rivers,
estuaries,  tidal  inlets  and coastlines.   This model has been successfully
applied to Lake Okeechobee for hydrodynamic study, Biscayne Bay for
dispersion study and Hutchinson Island for hydrotherma! study.  This
model can simulate the receiving water body in the far-field and de-
tailed features of thermal plumes and mixing zones in the near-field.

DESCRIPTION OF ANCLOTE ANCHORAGE

     The UM's thermal pollution team has been  studying the Anclote
Anchorage on the west central coast of Florida near the town of Tarpon
Springs  (Figure  1) since  the summer of 1978.  The Anclote power plant,
operated by the  Florida  Power Corporation, has two 515 MW, oil-fired
electrical generating units.  The once-through  cooling water for the
two units is to be drawn from  the  Anclote River  by six  pumps deliver-
ing a total of 930,000  gpm.  After a temperature rise of 6.1°C, this
water is diluted with a flow of 1,060,000  gpm at ambient temperature to
reduce  the temperature rise to  2. 8°C at  the outlet to Anclote Anchorage.
This mixing is done in a  1250 meter-long  man-made channel leading to
an outlet.

     The Anclote Anchorage (Figure 2) consists of shallow coastal water
separated from the Gulf of Mexico  by a series of barrier island parallel
to the  coast line. The Anchorage  has  a  relatively unrestricted exchange
of water with  the Gulf through  natural channels  to the north and to the
south of the Anclote Keys.  Depth ranges from 0. 3 to 3. 6 m, with a
mean of 1.8 m.  Shallow  regions of less than 0. 6  m comprise approximately
5 km in  length and 6  km in average width.   Currents in the Anchorage
are tidal- and wind-driven, with the tide entering from the south  stronger
than that from the north.

     Prior to power plant construction, the Marine Science Institute of
the University of South  Florida (USF)  was contracted by the Florida
Power  Corporation (FPC) to investigate the possible environmental im-
pact of the plant operation.  The USF's Anclote Environmental Project
beginning in 1970 was  comprehensive in nature.   Much of their efforts
was to obtain  a detailed  picture of the Anclote  environment prior to
alternation of the environment  by power  plant construction and operation.

-------
Valuable background data  were collected in  a  variety of areas,  including
general physical characteristics of the Anclote River and Anchorage,
suspended sediments, turbidity, temperature, salinity, water quality,
seagrasses, benthic invertebrates and fishes.   All  information was ap-
plied to preconstruction plant  design  in order to minimize environmental
impact.  This application resulted in major changes in intake and dis-
charge channel dredging,  outfall design and thermal dishcarge  character-
istics.   The same  team also monitored the undergoing changes of the
environment during and after  the dredging  operations for the intake  and
discharge channels.

     Unit I of the Anclote Power Plant commenced operation  during  the
fall of  1974, after  which a post-operation ecological monitoring program
for evaluating the nature and  degree of thermal impact was  carried out
and maintained by  the USF.   At the time of this study, the planned
Unit 2 was still  pending permission.

     Since the nature of thermal impact at  Anclote Anchorage is not
clear,  it is decided  that a joint effort by UM, NASA-KSC, EPA and FPC
to use a three-dimensional mode! with support of remotely sensed data
and in-situ data for calibration and verification may be appropriate to
study  this thermal  impact  at Anclote.  It is for this purpose that the
existing  free-surface model was adapted to  the Anclote site, as we  have
done in its application to Biscayne Bay.  As the original program was
developed for well-mixed shallow coastal waters, except  for  some modi-
fication to accommodate the tidally influenced  Anclote River  flow,  the
program  was ready for calibration.  We shall account for this application
in a later section.

-------
                              SECTION  2

                            CONCLUSIONS
     A numerical simulation of the hydrotherma! characteristics of a  well-
mixed, shallow, coastal water body at Anclote Anchorage under  the  effect
of waste heat disposal by a power plant in the summer and winter situa-
tions  is presented.  The model takes into consideration the effects of
Coriolis force, wind, tide, bottom  topography, power  plant intake and
discharge, river outflow and  surface heat transfer.  Results  obtained
with the model have been verified with  in-situ measurements  and 1R data.
Reasonable agreement  has been obtained.  From the experiences and re-
sults  of these simulation runs,  the conclusions may be summarized as
follows:

1.  Inputs for the description of open boundaries, discharge, wind  and
    heat transfer must be found experimentally from field data for ac-
    curate hydrotherma I predictions.

2.  The  shape of the  thermal plume was dominated by the stage of the
    tidal cycles, as clearly exhibited in  the plots  of IR and calculated
    isotherms.

3.  Tide plays a main role in the  Anchorage as a driving force;  its  in-
    fluence from the south  is stronger  than that from the north.  A
    dividing line is observed.  On the  northern side of this  line, the
    water flows north, while  on the south  side it  flows south.   At this
    ridge line,  the water has minimal transverse motion.   The location
    of this  line varies with time as the tide from  both ends is not in phase.

4.  It is important to impose  correct boundary conditions,  especially the
    correct tidal functions on the south and north boundaries, for ob-
    taining a good prediction of the thermal plume.

5.  Wind  does appear to be an important external force affecting surface
    currents.   It is found that the surface currents are close to the
    wind  direction when wind speed is  in excess of 15 mph.

6.  An estimation of the Rossby Number revealed that the nonlinear
    inertia  terms can  be safely neglected for  Anclote  site.  This is  be-
    cause the ratio of inertia force to the  Coriolis force  is small  compared
    with unity.

7.  The recirculation  of cooling water should  be prevented by proper
    design  of intake and discharge location, since the heated water

-------
    recirculating back to the intake will cause a reduction in the effi-
    ciency of the  power plant.

     Numerically the model behaves very well for both summer and winter
simulation runs.   The model is  able  to include very shallow water depths.
A problem associated with the analysis is computing the response of the
thermal plume to very strong winds such as hurricane, since the formula
used in this study to estimate the vertical eddy viscosity coefficient does
not involve wind speed.  A  brief investigation was made regarding the
effects of hurricane-force winds.  It was found that currents become un-
realistically large  in such  cases  unless the vertical eddy viscosity coeffi-
cient was increased with wind speed.  This general problem of computing
the response of shallow coastal  water to very strong winds requires addi-
tional research.

-------
                              SECTION  3

                          RECOMMENDATIONS
     One of the difficulties encountered in model verification with the
remote sensing temperature field was  the lack of measured current data
of comparable accuracy.  We believe that, at the verification state,  the
accuracy of velocity calculation needs to be assessed in  order  to clear
ambiguities about the limitations of  comparison of calculated surface iso-
thermal maps with the IR imagery mosaics.

     In the initial phase of mixing,  the plume shape is governed by the
volume of discharge, the geometry of outlet and the initial temperature
difference between  the  discharge and receiving waters.  Consequently,
for predicting the details of the near-field thermal plume a finer grid  is
needed near  the outlet.   Thus, a combination of grid  structures would be
desired for the computation with a course grid for far-field and a  fine
grid for  the  near-field.   This,  however, would result in higher computa-
tional costs.

     For any non-reactive, dissolved, chemical constituent, the governing
transport equation  is   completely analogous  to the temperature equation
of the present model.   Thus, without much endeavor, the model can be
extended to include the dispersal of these constituents.

     For situations of deep and stratified water,  the bouyancy effects
are important factors in shaping the plume.  Therefore, a coupled system
of momentum and energy equations  should  be considered as the basis of
formulation.

-------
                             SECTION H

      MATHEMATICAL FORMULATION  AND MODEL DESCRIPTION


GENERAL BACKGROUND ON FREE-SURFACE  (TIDAL) MODELS

     Recent concern about the ecological  future of our estuaries and
coastal waters has generated a need for practical and reliable methods of
predicting  the environmental impact of  the widespread use of natural  sea-
shores as industrial  zones.   The tendency of setting up  once-through,
seawater-cooled, fossil- or nuclear-fueled power generation plants is in-
creasingly  strong, due to the demand for electrical power at competitive
prices and prohibitive costs  of cooling  towers and man-made cooling lakes.
Thus, it is of special importance that mathematical models be developed
for simulating and predicting thermal pollution.  The estuaries,  tidal
inlets and  coastal waters which  serve as  the receiving body for  waste
heat are usually of complex configuration and topography.   The flow  is
driven by  tide,  wind, run off and buoyancy force.  Thus, a complete
hydrothermal model for coastal waters must, in addition to  solving the
three-dimensional equations of mass, momentum and energy, include the
salinity  equation to determine closely the local density.   A  model of this
type is  still far beyond feasibility.  Therefore, assumptions and approxi-
mations  must be made to render the system closed and tractable.

     Hindwood and Wallis  (1975)  have compiled a bibliography of 141.
papers concerning computer  models  for tidal hydrodynamics.   In  general,
the output from the  hydrodynamic/hydraulic model is recorded on magnetic
tape. This tape is then entered into the thermal dispersion model.  The
dispersion  model operates through successive solution of the finite differ-
ence equations  for the change of temperature (with time  and space) due
to diffusion (eddy mixing),  advection  (velocity transport),  and  heat
flux through the water surface.  Therefore, it is appropriate to  say  that
the hydrodynamic/hydraulic  model serves as the backbone to the thermal
pollution model of the well-mixed water.  Commonly, the  hydrodynamic
models are two-dimensional and  based on  the vertically-integrated equations
of motion and continuity for an  incompressible fluid.  Thus,  the  vertical
structure of  the circulation was not considered.   This procedure turned
out to be sufficient for the investigation  of tidal processes.

     Of the two-dimensional  models, we should briefly mention the Leen-
dertse (1967) model and Reid and  Bodine  (1969) model for their comprehen-
siveness and popularity.  In Leendertse's model,  the unknowns are verti-
cally-averaged velocities and water level; the bottom friction is in terms
of Chezy coefficient  which is to be calibrated.  Reid and Bodine  used

-------
vertically integrated transports and water level as unknowns and Darcy
friction coefficient.  Both models used the space-staggered  system which
gives the simplest scheme consistent with the control  volume approach to
deriving difference equations.  However, in  the Leendertse model, the
finite difference equations are expressed  in an  alternating direction im-
plicit (ADI) form, with two successive time-level operations being executed
during each time cycle; while in the Reid and Bodine model a time-split
explicit (leap-frog) method is used for  marching forward, the time step
is  limited by numerical stability  requirement  (cf. Platzman,  1958), less
than the value min(Ax, Ay)/ 2gD    ), in which D     is the maximum
depth.                          max             max

     In the context of two-dimensional tidaUy-driven  flows, we should
include the finite element  method models which  were developed in the last
few years.  This is a  method combining finite element for spacial discre-
tization and finite difference  for temporal discretization.   The variables
in  discrete element is  approximated by  simple polynomials whose  coeffi-
cients  are expressed  in terms of nodal  values of the variables and their
derivatives.   Correctly formulated, the physical conservation principles
are satisfied and, in  theory,  the element shape and interpolation function
are quite free as long  as  certain compactibility  conditions are satisfied.
The freedom of  using  irregular grid to fit complicated geometry  is the
primary  advantage  of  this method.

     At the present stage, FEM  (finite  element method)  models for tidal
hydrodynamics are all two-dimensional.   Wang's CAFE (Circulation Analy-
sis by Finite  Element)  model  and Brebbia's shallow water model  (in his
Finite  Element Hydrodynamic  Problem Orientated Language (FEHPOL)
package)  are  both productive.   CAFE  (Wang,  1978)  has a linear  triangu-
lar element for all variables,  vertically-integrated transports and water
level,  and has split-time method for time integration,   Patridge and
Brebbia  (1976) use a six-node, quadratic triangular element for all vari-
ables;  in this case,  they are vertically-averaged velocities and free-
surface elevation.  The 4th order Runge-Kutta (explicit) and Trapezoidal
Rule (implicit) methods are used for time integration.  In general, the
flexibility of the grid  layout,  the consistency of FEM  formulation  and the
easiness in taking into account the spatially  variable  properties are the
advantages.   However,  the time integration scheme  is very much problem-
oriented, and more fruitful research could lead to a better  time  integra-
tion method.  It is worthwhile to point  out that, in the  FDE (finite
difference method)  for two-dimensional  hydrodynamics, the  explicit me-
thod had advantage over  the implicit method in terms of computer time.

     In the classical approach of two-dimensional computations  of tides
and storm surges in a shallow sea, it is assumed that the velocities are
uniform  over  the vertical; however, the vertically uniform velocities  are
of little  evidence in the case of those propagation and transport processes
which  are essentially  vertically structured.  The application of a three-
dimensional model requires extensive computation on a large number  of
grid points if relatively good  resolution is needed.  This precludes the

                                    10

-------
use of the implicit method which requires that matrix solutions of a row
or column of variables  be found.   In fact, Leendertse et ai. (1973) deve-
loped a three-dimensional model which really uses the explicit leap-frog
method instead of the implicit methods such as the ADI method Leendertse
used in his two-dimensional  model.

     Leendertse's three-dimensional model is a vertically-integrated, multi-
layered model, in which the usual  assumptions of hydrostatic pressure,
imcompressibility and small density  variation are made.  The interfaces
are assumed as  fixed horizontal planes, while the bottom layer has its
lower face fit to bottom topography and the top  layer has its  upper face
representing the free surface.  The top layer has a time-variable height.
The other layers may be intersected by the bottom  and have a height
which  is dependent on  the bathymetry.  The  description of the finite-
difference equations from the differential  equations  is accomplished in
two steps:  First,  the  equations for the layer are derived  by vertically
integrating  the  variables over the  layer thickness, and subsequently,
finite difference approximations for the layer  equation are developed.
Thus,  the vertically-integrated momentum, heat, and salt balance equa-
tions can be presented for each layer.  From the continuity equation,
the time derivative of the free surface and  the vertical velocity at the
interfaces can be  derived.   The horizontal pressure gradient in each
layer  is approximated from hydrostatic equation with layer-averaged
density, which  is tied  to salt and heat through the equation of state.

     In Leendertse's model (Leendertse and Liu, 1975),  the balance of
momentum fluxes between the local and convective accelerations,  pressure
gradient, Coriolis force,  lateral diffusion, and Snterfacial stresses is
accounted for within each layer; likewise, the balance of heat fluxes  be-
tween  the time rate,  convection, lateral diffusion and cross-layer diffu-
sion  is  for each layer;  so is  the salt  balance.  The same lateral diffusion
coefficient is used  for  both heat and  salt.   Even postulations are intro-
duced  to express  the vertical momentum,  salt and heat exchange under
vertically stable or unstable stratifications;  these layer-averaged equa-
tions clearly point  out  the foremost problem of his model.  There are
many  parameters to be determined  and an enormous amount  of supporting
field data is needed.

     Leendertse et al.  (1973)  discuss  the  numerical  finite difference
solution scheme  in  some detail.  The explicit leap-frog method is used  to
avoid  the difficulties encountered by  ADI.   The  spatial grid structure  is
cell-like with u, v and w components at the center of the corresponding
normal faces and pressure (p), density (p),  salinity  (S) and  temperature
(T) at the center of the  cell.  For  programming  reasons,  the  bottom  must
be approximated in steps of layer thickness, causing  some numerical  pro-
blems at the jumps.  But the most troublesome one is the specification  of
boundary conditions for the seaward boundary of the model, since field
measurements do not provide enough  data to describe that  boundary
completely.   Numerical  experiments  involving the simulations of hydrody-
namic behavior of both Chesapeake  Bay and San Francisco  Bay are given

                                   11

-------
with a graphical  representation of velocity, salinity and temperature  results.

     A specific three-dimensional model has been developed by Sundermann
(1974) and applied  to the North  Sea.   In  the model,  the usual assumption
of a hydrostatic  pressure and the Boussinesq approximation for the turbu-
lent Reynolds stresses  are  made.  However, no horizontal  turbulent ex-
change of momentum is retained.  To fit the vertical velocity profiles of
a wind-generated surface current  and a compensating countercurrent near
the bottom, as observed in nature, it was found necessary to include a
bottom boundary layer  that can  be modeled by having a vertical  eddy
viscosity  dependent on the depth.  The spatial discretization is carried
out by means of a  cubic grid net, while the explicit method of Crank-
Nicolson  has been used for approximating the vertical diffusion term.

     The  a system  of coordinates, devised by  N.  A. Phillips (1957),  was
used  for  numerical forecasting.  Its advantages are that the kinematic
boundary conditions at free-surface and bottom are made  simple and  a
vertical stretching  is used  to avoid the difficulty  of using regular grid
net for  irregular bathymetry.  Briefly, the (x, y,  a, t)  system  is used
to replace (x, y, z, t), and  the free surface and  bottom  are trans-
formed into a = 0 and a = -1  surfaces.  The modified vertical velocity,
Q = dcr/dt,  is zero at both  these surfaces, while the actual vertical
velocity is w = dn/dt at the free surface and w =  0  at the bottom.

     Freeman et  al. (1972)  introduced u-transformation into their three-
dimensional free-surface model for wind-generated  circulation in a closed
region.   Sengupta  and Lick (1974) incorporated ^-transformation  in their
rigid-lid model for  wind-driven  circulation in  lakes.  However, the trans-
formation was used only for the advantage of vertical stretching  so that
a constant vertical grid size, A a,  can be used throughout the  domain.
This  rigid-lid approximation effectively eliminates  surface gravity wave
by  imposing a zero-actual  vertical velocity at the  water surface.   An
additional equation, a  Poisson equation for  surface pressure which con-
tains the rigid-lid  condition,  can be  derived.  Therefore,  at each time
step the Poisson equation  has to be solved  by iterative method.   However,
for large problems, particularly when normal derivative boundary condi-
tions are used,  this can be a time consuming part of the  calculation  pro-
cedure.   As mentioned, this  model is  not for tidally-driven coastal flow.

     The present free-surface model  is similar to  that of  Freeman et  al.
in using  o-coordinate  system.  However,  the horizontal turbulent exchange
of momentum terms is  neglected  while the vertical  turbulent exchange term
is approximated  by using the Dufort-Frankel differencing  to ease the con-
straint on numerical stability due  to  neglect of the horizontal  diffusion
terms.  This model was developed by Carter  (1977) for his study on wind-
driven flow in Lake Okeechobee.  Later,  Sengupta et al.  (1978) applied
it to  tidal flow in Biscayne Bay.

GOVERNING EQUATIONS AND BOUNDARY CONDITIONS

      Flow of water in  estuaries  and coastal areas  is  predominantly hori-

                                     12

-------
zontai.  Vertical velocities do occur and are important, as they character-
ize the  vertical  circulation.   However, the vertical acceleration is extremely
small  in comparison  with  the  gravitational acceleration, so that the vertical
acceleration is neglected  and the vertical equation of motion is then re-
placed by the hydrostatic assumption.  The  effect of density  variations
on the inertial and diffusion  terms  in the horizontal  equations of motion
is  neglected.  Density  variation is  retained in the lateral pressure gradi-
ent terms; that  is,  the effect of bouyancy is accounted for by allowing
density variations in the horizontal pressure gradients which  influence
the fluid motion through  the horizontal momentum equations.  Consequently,
the continuity equation replaces the mass conservation equation.   In the
estuaries, the density is influenced by the salinity and temperature; in
this model,  the  density will  be taken dependent  only on  the temperature.

     For turbulent closure,  the Boussinesq approximation for  the turbu-
lent Reynolds stresses  is made, and lateral dispersion of momentum is
further assumed unimportant when  compared to the rest  of  the terms in
the horizontal equations of motion.   The vertical  eddy viscosity, the
vertical and  horizontal  eddy diffusivity of temperature are assumed to  be
constants,  although the horizontal eddy  diffusivity  has orders of magni-
tude larger than the vertical eddy  diffusivity, being due to the  much
larger horizontal characteristic length,  L,  in comparison  with  the vertical
characteristic length,  H.

     With these  conditions, we can  write the set of governing equations
expressing the conservation of mass, mementum,  and energy in incompres-
sible  flow:
              9t  ' Tx '  T^ '  ~3z~   iy  '  p  3x"   "V3^ = °            (1)


              3v  , 3uv ,  3vv ,  3vw  ,  r.. .  1  3p   „  32v _ n
                                                        " U            (2)
                               -H + oe = 0
                               3z   p§   u                             (3)
                            3x    3y   3z                               (4)


            3T L 3uT  L  3vT t  3wT   „ ,32T  ,  32T,    „ 32T   „
                                P = P(T)                               (6)


                                    13

-------
where x,  y,  z = Cartesian coordinates positive eastward,  northward and
                 upward, respectively

      u,  v,  w = Respective  components of velocity in x, y,  z direction

      t        = Time

      p       = Pressure
      T

      f
       B,
       B.
= Density

= Temperature

= Coriolis parameter

= Gravitational acceleration

- Vertical eddy viscosity

= Horizontal eddy diffusivity

- Vertical eddy diffusivity
     The first three equations represent the equations of motion.   Equa-
tion (I)  is the equation of continuity, and Equation  (5) represents the
energy equation.  The equation  of state. Equation (6), expresses  the re-
lation  between the density and the temperature.  For a more complete
representation of an estuarine ecosystem, an equation expressing  the
balance of the salts dissolved  in the water should  be included  together
with balance equations of dissolved substances  that are important  to  the
analysis.   All these equations  have the same form  as Equation  (5).   If
salinity is included  then the density should be  related to  both salinity
and temperature.

     Boundary conditions for the above equations  must be specified:

     u = v = w = 0                  at z = -h(x, y)                     (7)

                              \
      p = 0
      3u
      92
_*   av   Ty
A  '  82 = pAT,
      3T
            K
                 S  at z = n(x, y, t)
      aif
                  (Te  -  V
                 y~ - w = 0
                                                       (8a)
                                                       (8b)
                                                       (8c)
                                                       (8d)

-------
where n      = Free-surface elevation above mean water level

       T ,  T  = Surface wind stresses
        x   y

       C      = Specific heat at constant pressure

      K      =: Surface heat transfer coefficient

       T      = Equilibrium  temperature

       T      = Surface water temperature

Also, boundary conditions at closed and open  boundaries must be known.

     One major  difficulty in the treatment of the  free-surface mode! is
the fulfillment of kinematic condition  at the free surface.  The approach
used  in the model is to follow a vertical stretching transformation sug-
gested by  Phillips  (1957) and used successively in lake circulation  studies
by  Freeman et al.  (1972).   Using  this o-transformation,  the free surface
becomes a  fixed, flat surface and  the variable depth bottom becomes a
fixed, flat bottom of constant depth.  Thus, this  method allows easy
adaptation  to  various bottom topographies, and a constant vertical  grid
size can be used throughout the domain.

     With respect to the mean water  level  (MWL),  z is n(x, y, t) at the
free surface and -h(x,  y) at the  bottom.  The a-transformation of the
vertical coordinate for  the free-surface model is obtained by introducing

                    _    z - TI(X. y, t)          z - n ..
                  a ~ h(x,  y) + n(x, y,  t)   H(x, y, t)               (9)


where H = h + n is the total depth measured  from local water surface.
Figure  2 shows  the (x, y,  a) coordinate system.   Note that the  value
of o decreases monotonically from  zero at the free surface to minus unity
at the bottom.  The modified vertical velocity,  ft,  is clearly zero at both
free surface and bottom.  From  Equation (9),  the actual vertical velocity,
w,  is related  to fl by


                w = (h + n) n + (0 + 1)|3 + ff(u|| + vf|)             (10)


     Differential transformation relationships are required to convert
equations in (x, y, z)  to those  in (x, y a).  The first derivatives can
be written  as


                 3F       _  3F          if n— +  -fLlhlE
                (}      = (}       " H1 3x    3x  3a
                                    15

-------
                                        l   3H .  3TK3F
                           r JF,      _ 1  3F
                           I3zjx,y,t   H  30                          (lie)


where F Is any dependent variable and H  = h + n is the total depth and
is independent of a.  The continuity equation can be written as
By integrating Equation (12)  from  the bottom to the surface and observing
the fact  that fi vanishes at either face,  we obtain the  first useful form of
the continuity equation.
     The  horizontal  momentum  equations are
                           3Hvu   ^u +  (1 +  } |u
             3t     3x      3y     " 3a-          3a


The energy equation  can  be written as

                                                              B
     The coupling of momentum and energy equations may be  retained via
the density which is assumed to be a  function of temperature  only.   An
empirical formula

                                  16

-------
                p  = 1.029431 - 0.000020T  -  0.000005T2                (

may  be used for sea water of salinity  38 parts per thousand.   The hydro
static equation can be integrated  to obtain a diagnostic equation  for p.
     The continuity equation,  Equation (12),  can also be integrated from
the free surface  to a to yield  the modified vertical velocity n  at a plane.
Thus,  a second useful  form of Equation  (12)  yields
- £ In - Ira fMH + 3Hv
~ H 3t   H\  (  3x-     3y
           -'0
                                                                      (19)


The actual vertical velocity w is recovered through relationship Equation
(10).  In solving  parabolic-type partial-differential  Equations  (14)-(16)
and differential Equation (13), we choose the initial conditions to  be zero
elevation from MWL and  zero  velocities.  Given the initial  temperature
field, the boundary conditions are u = v = w = 0 on all solid  surfaces,
and at a = 0, the wind stresses (t • , T )  are  exerted;  i.e.

                   A           A
                   2y3_u       _y §z  _
                   pH 3cr ~ Tx' pH 3a   Ty at ° ~ U                  (20a,b)

The adiabatic conditions are assumed  on all solid  surfaces, and at a =  0,
the heat flux is set proportional to the difference between surface tempera-
ture T  and  equilibrium temperature T ; i.e.

                   PC B   «T   ,,
                    ...P v  |£  = K  (T   - T )  at a  = 0                   f2;n
                     H    3cr     s  e    s                             t.^1-'

     The specification of velocities and temperatures at the open boundar-
ies where the tide enters  the  basin is more difficult.   In  the  present study,
the tides outside the basin are given, and the difference  of surface eleva-
tions across  the open boundaries is used to  determine the normal  velocities
there.  The  temperature at the open  boundaries  is set equal to ambient
temperature.  The velocities and  temperature at  the discharge is specified
according to plant operation.

UNCOUPLED SYSTEM

     The foregoing system may be uncoupled so  that  the  hydrodynamic
model for tidal-driven flow is  separated from  the  thermal  dispersion  model.
The former system consists of Equations (13), (14),  (15)  and  (19),  while
the latter contains Equation (16)  only. This decoupling amounts to  an

                                   17

-------
assumption that the density variations due to the temperature rises
caused by  waste heat  dispersion are negligibly small,  so that the water
is  of constant density and  the energy equation can  be uncoupled from
the continuity and momentum equations.   The lateral pressure gradient
is  thus  directly related to the gradient of the free surface.

     In  the present problem of shallow water flow strongly influenced by
tide, the thermal  dispersion is mainly caused by velocity transport and
the buoyancy effect is negligible with comparison to the tidal effect.
Numerically,  this  decoupling has two  obvious advantages.   Firstly, the
time step of  the coupled model is the same as that of  the  hydrodynamic
model.   This time step is limited by the criterion for  computational
stability by the explicit method; namely,  At  must satisfy the condition
cAt/As«  1,  where As represents one of  the space intervals in the three-
dimensional grid  and c denotes the maximum  characteristic speed.   In
the free-surface model,  c is identified with the  wave  speed and As with
horizontal  spacing; thus,
                            At  <     A*,  Ay)
                            "• *. .            '
                                        max
However, the stability analysis of linearized energy equation alone yields
the thermal time step to be
                         1 Ax  '  Ay  *  2Bh(Ax^ H

where u and v are the maximum particle velocities in  x and y direction
respectively.   For the present  problem, these time steps are 50 sec and
100 sec respectively;  therefore, in using At  =15 sec, we can store the
hydrodynamic results  every 20  steps to match At— =  300 sec used in the
detached calculation of temperature field.  Secondly,  the decoupling
allows one to  try  for  flow solution before it is used for temperature cal-
culations, and the appropriate  flow solution may be used for many tem-
perature solutions of  various initial and boundary conditions of tempera-
ture.

     The effect of stratification Is known to arrest the thermal dispersion.
For a  shallow basin having an  insignificant river  discharge  in comparison
with tidal flow, which is the case of  concern,  the basin is nearly well-
mixed; the density variation is small.   Thus, uncoupling is physically
sound and numerically beneficial.

COMPUTATIONAL GRID

     To represent the equations in finite difference form, a horizontal
staggered computing grid system is  used.   Its  plan version, shown in


                                     18

-------
Figure  4,  indicates the arrangement of field variables in the x, y plane.
Figure  5  shows a vertical fluid  column subdivided into four layers;  each
has  a constant  nondimensional thickness,  A a = 0.25.   The  u,  v and w
velocities are shown at their definition point respectively.   Having divided
the  region into cells by a series of grid  points which are spaced at dis-
tances of Ax,  Ay, and A a,  the time variable is differenced into increment
of At such that t = nAt, where n denotes the current time step.

     The water depth, h. .,  is  specified  at  full-grid  point  where both
horizontal indices, i and f/  are  integers;  however, since the total depth,
H - h + n, is  needed at half-grid point,  where both indices i and j have
half-integer values,  (i + $,  j +  i), the specified water depth  at full-grid
points is  averaged.  The  following notations are  used to indicate water
depth at  the u-,  v- and n-points.
                                                                       (23)
     A  notation for n at full-grid point is needed for calculating the
nonlinear inertia terms of the momentum equations:

               _      .,  n  ,   n     ,  n       j.ni
               Ei,j  = Hni,i  + Vl.j +  Vl.i-1  *  Yj-l3                 (24)


where n  indicates  the present-time surface elevation above the mean
water level.   With  Equations (22)  to  (24) and the calculated water level,
n..,  at present time  step, the present total water depths at u-, v- and
n-points are respectively given by the following  expressions.
                                                                       (25)
                                                                       (26)
                            H. .  = hn
                                    . .   -    .
                                    i.J     i.J                           (27)

In the momentum equations, values of u and v are required at  half-grid
points where values of these variables  are not  defined and thus, not
stored.  In these cases,  values are obtained by linear interpolation be-
tween  the values stored for that variable at the two neighboring  points.

                                     19

-------
Therefore, the following  average notations are adopted.


It is clear that notations  (28) and  (29)  give velocities at half-grid point
while notations (30)  and (31) give velocities at full-grid point.  No
superscript n  is needed for these notations.  These and other variables
are shown in  Figure 6,  where u-, v-  and n-points are marked.

FINITE DIFFERENCE EQUATIONS

     In the equations presented here, the space  subscripts are i, j,  k
and  are all integers, unless  otherwise noted.  The superscript n refers
to the time level; n is  present time, at  which  the predictive equations,
(13) through  (16),  are  used to  advance the fields of u, v, n  and T  to
new time  level n + 1.  In this study,  first derivatives with respect to
time are always represented  by  central time difference, i.e.,  for any
variable F,

                                   -  F?:1,
                  V3t'i,j,k        2At         v^" '                     (32)


where 0(At2) refers to  the order  of truncation error associated with this
differencing.  The  central time scheme is thus of second order accuracy,
since the neglected terms are terms multiplied by (At2)  or a higher order
of At.  This indicates that the finite difference approximation can be
made to approach the differential  by taking sufficiently  small At.

     The second derivatives  with  respect to a in Equations  (14) to (16)
are written  in the DuFort-Frankel format,  i.e.
                 n      F         F     -_ pn-l    pn
                      -  i.i.k-1    1.1.k    i.i.k    l,j
             3a2i,j,k                  (A a)2                             (33)

 This expression for  second order crderivatives when used with the central


                                   20

-------
time difference is known to avoid the inherent instability caused by using
the common form of second order difference;  i.e.
     The continuity equation,  (13),  is the predictive equation  for surface
elevation  n/ and thus applies at n-point.  For convenience, we introduce,
at level k of (i,j) -column, a variable
     Fn =
      k              Ax                          Ay

where the superscript n on the  right-hand-side  variables has been dropped.
Then the finite difference form of the continuity equation can be  written as

                           n+1     n-1
Using Simpson's  rule, the above equation gives
                 = "U1 + (2At)T CF1 *  4F2 ^  2F3  + 4F4 + F5]           (34)


It is to be noted that Fg = 0 since both  u and v vanish at the solid  bottom.

     Another  form  of continuity equation,  (19), can be used to calculate
the modified vertical velocity n at each level of the column.  That is


                               ...  ,- .     n+1    n-1
                     o      - "fk"1^Ag  n, .-  - n. •
                      i i k  -- H       -^ _ 1>J
                      lf]lJS      Hi,          2At
                                                                      (35)
where the trapezoidal  formula  has  been used and  H.  . is  from Equation (27);
i.e.                                               ''*

                       H. .  =H*  = h?. + n?,
                        i,]     i,]     i.]    i>]
                                     21

-------
In the computation, always the present  time  ft  only  is needed; therefore,
no superscript n is required for ft and,  thus,  one array is needed to
store J2.  On  the other hand, three arrays are necessary to store ri  at
the past (n-1), present  (n)  and advance (n+1) stages.

     The actual vertical  velocity itself,  w, is not  needed in the calcula-
tion.   However, when  it is desired, it can be obtained from Equation (10),
which  is rewritten  here in finite difference form.
              wi i k = Hi
               l.J.'K    l
                                              n+1    n-1
                                                 2At
                                           2Ax
                         IJ'                2Ay
                        n+1   n-1
                           2At       i.j.k      2Ax
                                  2Ay                                 (36)
     The predictive equations for velocities are given by  Equation (1H]
and  Equation (15);  in  finite difference form, these equations are


 Ui,j ui,j,k " Ux,1 ui,j.k _ -   vi,j,k + vi-l,j.k + vi-
           2At           ~   i,j                   4

                      .                n-+l     n-1
            "              a       " U     " "
                      Uj~              (A a)2                  " Xi,j,k

                                                                       (37)

 7n+1vn+1  - V11'1  n~1         U      4- U        +            +
 Ad—Llik	LJ—iJLiik_ = -fv>  . ^J^ —klli-k	i+l,j-l,k	i+1'3'k
           2At
                       v   i.i.k-l   Vi,3,k   vi,j,k
            Ay       V. .                  (Aa)2                  "i.j.k
                       i.J
                                                                       (38)

                                   22

-------
where, unless otherwise noted, ail the superscripts are n  and thus are
dropped  for convenience.  The nonlinear inertia terms X. . .  and Y. .  .
                                                      i,J,rvI,J,lN
are
         Hi.iuf.i.kuf,j,k
       ~
                         Ax
                               Ay
                                        n-KL    n+1      n-1    n-1
                                       "1,1  *  Vl,1  "  \,]  " Vl,i
                                                    4At

                                                                      (39a)

  i,j,k "                  Ay

        . (hm.i ^Vi.i^^i.i.kVi.i.k- (hi.i +Ei,i)ui7-n.i.kvf+i.j,k
       +                            Ax
                                                n+1    n-1    n-1
                                               "1,1-1" "l.j  " ni.i-l
                                                   4At

                                                                     (39b)
Unless otherwise noted, ail the superscripts  are n  and thus dropped,
since the inertia terms are calculated at  present time level  n.   In Equa-
tions (39a,b),  r\n   and n are presently obtained variables; since at each
time level n and B are calculated first by Equations (34)  and (35)  respec-
tively, u and v are followed by  Equations  (37) and (38).

     When the  calculations for velocities  are  to be  performed  using-the
predicative equations,  (37) and  (38), the variables u. . .  and v. .  .
which appeared respectively on the right-hand-side of ^Equations1'(37) and
(38) should be  rearranged, so that a form for velocities similar  to Equa-
tion (34)  appears.  In this form, all quantities on  the right-hand-side are
obtained based on  specified h, calculated u and v  at present and previous
time levels, newly  computed Q and n at three time  levels.   All the variables


                                  23

-------
and notations that  appeared in Equations (37)  through (39) have been
defined and are shown in  Figure 4.   However,  in Figure  4, the common
superscript n and  subscript k  have  been omitted for brevity.  Figure 4
also shows clearly  the horizontal extent  of variables involved in the cal-
culation of n, u and v in  the  cross-hatched region, namely n. ./ u. . .
and v. . .  (for  k = 1 to 4).                                 I/J    I/J'K
      i/J/ K
     Similar to the treatment of the  momentum  equations,  the  energy
equation, Equation  (16),  can be written In finite difference form as
                          KLprH-1    Hn-lTn-l
                        i.j   i,j,k "  i,j   1,1. k  _
                                 2At

         T11        9Tn    4. Tn       Tn        9T11    -I- T11
              .1. k "   i,j,k *  i-l,j.k  ,   ,]+l,k   ^Ai.1.k'    i.j-l.k,
                 n      _ Tn+l   _  Tn-l    Tn
                           1.
                                (A
-------
for clearity.  Equation (40),  similar^tp Equations  (37)  and (38), can
lead to a  predictive equation for T.  . . .

     In the model, since the  density is considered as a constant, the
equation of state, Equations  (17), is of no  use, and the system  is un-
coupled.  With that,  the temperature field can be solved separately pro-
vided  that all results from the hydrodynamic model  are available.  That
is,  the spatial and temporal variation of temperature is solved after the
spatial and  temporal  variation of velocity  fields is known.  Alternatively,
Equation  (40) can predict Tn '  after obtaining n    , 0, un   and v" f
at each time cycle.   The former involves  two separated programs while
the latter is a coupled program  only.   In coding the model above, both
possibilities are taken into account by using flag  statements  (see Users
Manual for Free-Surface Model,  1980).

SOLUTION PROCEDURES

     Clearly,  Equations (37)  and (38) are for u and v at interior points;
that is, they  are  not on  the  boundaries.  Figure  4  implies why  these
formulae are not for  the  boundary points.  Therefore, Figure 2  indicates
that the normal velocity along the boundaries must  be  either specified or
calculated by some other means.  For solid  lateral boundaries parallel to
x- or  y-axis, the normal velocities are specified to  be zero all the time.
At the river mouth and discharge outlet, the normal components are given
by the known flowrate and the average water  depth.   At the seaward
boundary, the imposing tide  is specified along a parallel line at half-grid
size away from the boundary.  The  difference of  water elevations at
half-grid  points across the open boundary is used to calculate the normal
velocity there.

     At the bottom,  the no-slip condition requires that u = v = 0 at all
times;  hence, there is  no need for doing Equations  (37)  and  (38) at K = 5.
At the free surface or k = 1, the wind-produced  shear stresses  are pre-
scribed as a function of space and  time, and can  be written  as
These  conditions demand a modification of predictive formula for un+1  and
v    at k = 1.  That is, the DuFort-Frankel format for second derivative,
Equation (33), is  replaced by (indices i and j are omitted and values of
k are indicated)
                 ,  .,
                 K=l
                                   25

-------
and a similar one  for v-momentum  equation.   In this scheme, the free-
surface,wind drag has  been incorporated to drive the current through
the u    and v    at the  free surface.

     The boundary conditions on vertical velocity w  are  the kinetic con-
dition


               w = ^              at a = 0

and the rigid bottom condition

               w - 0               at a = -1

These conditions expressed in terms of modified vertical  velocity fl are

               fl - 0               atc = 0ork =  1

               Q = 0               at a - -1 or k = 5

As stated earlier, these conditions have been incorporated into  the  first
integral form of the continuity equation; namely,  the predictive equation
for surface elevation. Equation (13) or Equation  (14) with F  = 0.   Conse-
quently, the finite difference formular for Q, Equation (35), is  for  k = 2,
3 and 4, while Qk=s1 and  8k=5 are set  to zero at all  times.

     The boundary conditions of energy equation  are adiabatic on solid
boundaries and known  temperature or  zero normal derivative of tempera-
ture on open boundaries.   It can  be seen that if  one neglects the  lateral
thermal diffusion terms, which are small in comparison  with convection
terms anyway, then no adiabatic condition is  needed on  the solid lateral
wall;  this is due to zero  normal  velocities on  the  solid  lateral  wall.   How-
ever, condition is still  needed on  the open boundaries.   At the bottom,
the heat flux is zero;  while at the surface, the heat flux is proportional
to (T  - T  ) .These conditions are required when  Equation  (40) is used
to calculate T     at k  = 1 and k = 4 respectively.

     As in most other  hydrodynamic models for transient problems,  the
computation routine works step-by-step  in time.   This  means  that the
computation proceeds through a  sequence of time  steps,  each  advancing
the entire  flow configuration  through a small, but finite, increment of
time,  At.   The results of the present  and the last steps act as a basis
for the calculation to proceed to the next one, whereby  the initial  con-
ditions can  develop,  within the limitations imposed by the boundary con-
ditions, into the subsequent  flow  configurations.  That is, provided that
the values of dependent variables  are  known  initially,  the values at  sub-
sequent times are obtained by using the explicit scheme.  The leap-frog
finite difference formulae. Equations (34), (36),   (37) and  (40), predict
surface elevation,  n, and  two horizontal velocities,  u and v, and water
temperature, T, at time level n+1. The vertical velocity w does not
appear in these equations and thus, is left out until needed for picturing
                                    26

-------
the flow configuration.  Each time cycle itself contains the following  recur-
sive processes.

1.  The countinuity equation in the form of Equation  (34)  is used  for n
    through  the mass conservation  of the fluid  column at (i,j).   At the
    same time, the modified vertical velocity fi    is also calculated at
    k =  2, 3 and  4,  since it is also based on the continuity equation but
    in a different form; namely, Equation  (35).

2.  The nonlinear inertia terms  that appear in the horizontal momentum
    equations, (37)  and (38),  are considered as driving terms;  that is,
    they are calculated from formulae (39a,b) by using known  results
    at the present time  (n). This  calculation usually takes a significant
    portion of the total computing time;  thus, it is  advised to drop this
    calculation whenever justified.   For  general tidal  flows, the inertia
    effect  is negligibly small in  comparison with the Coriolis force.  There-
    fore, if the Rossby  Number, which is the ratio of inertial force to
    Coriolis force, is very  close to zero, the calculation of X and Y  terms
    can be skipped  entirely.

3.  The horizontal momentum equatipns,  Equations (37)  and (38), are
    used to calculate u     and  vn   .  Here,  the specified  boundary  con-
    ditions,  such  as prescribed normal  velocities, surface wind  stresses
    and  specified  tide, come into play.

4.  The nonlinear terms, R, in  the energy equation,  Equation  (40),  is
    calculated by  formula Equation  (41).   Unlike the  nonlinear terms,
    X and Y of the momentum equations, R is to be included and calcu-
    lated under normal circumstances.

5.  Temperature at  advance step, Tn   ,  is calculated by Equation  (40),
    and  the  adiabatic condition, known  discharge temperature,  given
    ambient temperature and surface heat  transfer rate play a  part in
    determining the temperature field.

STABILITY

     The leap-frog  (explicit) method has a limit on  the size of the time
step.  Exceeding  this limit makes the computation unstable.  According  to
Platzman (1963), the maximum time  step for an inviscid  linear system, i.e.
the system with viscous and nonlinear inertia effects  neglected, the maxi-
mum allowable time step is


                        ..       r,f 2
                        Atmax = "I5   +

where  H     - maximum depth  in the problem.  If the Coriolis effect  is
also neglected, then


                                     27

-------
                             At
                               max
     The leap-frog method with the DuFort-Frankel  scheme for vertical
diffusion terms has been  adopted in the  present model.  As pointed out
by Forsythe  and Wasow (I960), this format is unconditionally stable  for
a pure diffusion  model, a system in which the viscous diffusion is the
only mechanism  responsible for transport.  Although,  the  present hydro-
thermal model is  considerably more complex than a pure diffusion model,
the use of the DuFort-Frankel format is  an important consideration for
a shallow water system because the vertical diffusion criterion  tends to
become relatively more restrictive as the vertical dimension becomes
smaller.
                                     28

-------
                             SECTION 5

              APPLICATIONS TO ANCLOTE ANCHORAGE
INTRODUCTION

     The present model has been successfully applied in  thermal disper-
sion studies at the Anclote Anchorage site.  The  site, located on the
west central coast of Florida near the town of Tarpon Springs, is a shal-
low channel between the mainland and Anclote Key which separates the
channel from the Gulf of Mexico, as shown in Figure  1.  Anclote Anchorage
is of interest since it is the  receiving water body for the Anclote Power
Plant cooling water discharge.   The channel is  relatively shallow  with
depths ranging from 0.3 to 3.6 m.  Shallow regions of less than  0.6 meter
comprise about 35% of  the Anclote Anchorage area, which is approximately
5 km  in length and 6 km in average width.  The  principal driving  mecha-
nism for current circulation is  tidal flux at the  north and south entrances
of the channel.  The tide is predominantly semidiurnal with  a mean range
of two feet.  Earlier measurements of temperature and salinity indicated
the currents flow  in and out through  both entrances.  However,  the water
exchange appears to be stronger in the south than  in the north, or the
currents generally flow in  the  north  direction during flood tide and flow
in the south direction  during ebb tide.

     The Anclote Power Plant operated by the Florida Power Corporation
has two  515 MW,  oil-fired, electrical generating units. Cooling water  is
drawn from the Anclote River through a man-made channel.   The six
pumps delivering a total of 1,990,000  gpm (125.6  m3/sec) are designed to
raise the water temperature of 2. 8°C above  the ambient water temperature.
The heated water  is discharged back into  the Anclote basin  through the
discharge channel with dredged submarine extension.  The designed
total flowrate is approximately  53 times the long-term average flowrate
of the Anclote  River.   At present, only Unit 1  is operative,  while  Unit
2 is awaiting permission.  The present flowrate, therefore, is 995,000
gpm (62.7  m3/sec).

     The model as  applied  to the Anclote Anchorage shows its capacity
for considering the effects of geometry and  bathymetry,  spatio-temporal
variation of the free surface, various boundary conditions, including tides
of different phase and  range,  surface heat transfer based on equilibrium
temperature concept, and changing meteorological conditions.  In addition,
turbulence has been considered by using the eddy  transport approxima-
tions.
                                   29

-------
CHOICE OF DOMAIN AND  GRID SYSTEM

     In this study, the Anclote Anchorage is designed  to contain the area
of water bounded  by the Anclote  Keys and sand barriers on the west,
mainland on the east,  and two imaginary E-W lines  drawn from the northern
tip of the barriers and the south end of Anclote  Key to the mainland.
These  two  imaginary boundaries are  considered to be so far that  the
thermal plume will not reach them even  at low tida! stage.  Its location
and schematization are shown in  Figure  7.  This  water  is open to tides
from the Gulf of Mexico at both ends.   The Anclote River,  where the
intake  of cooling water is  located, is also included  since the recirculation
of discharge flow  is of concern to the power company.

     The adopted  horizontal  grid  layout and its index are also shown in
Figure 7.   The  grid work Is allowed to  orient away from  north-south,
east-west  system, but in general, the y-axis of the grid  system aligns
with south-north.  Thus,  the  subscript i increases eastwardly while j
increases  northwardly.  The z-axis is chosen  upward from  mean water
surface while  the  subscript k  increases  downward from  the water surface.

     The selection of the grid size is governed by  several  constraints.
If the  grid  size is too large,  the  approximation of the channel in the
system will  be inaccurate, and at certain size,  the  computation will  become
meaningless.  However,  decreasing  the  grid size  will  lead to a considerable
increase in  the  computer time  since the  computations  must then be made
on  more points.   In addition,  the time step will decrease  because the
dispersive  properties of the computational method are related to  the  ratio
of the  time  step to the spatial grid  size.  After some numerical experiments,
the model on a  16x 14 x  5 grid with grid size Ax  =  Ay =  416. 75 m  and
A a = 0. 25 is considered  a good compromise between the resolution desired
for the region near the  discharge and the limitation of  the  computer
(UNIVAC  1108)  at  the University of Miami.  Also, care has been  taken
to have the intake and outlet at grid points.  The  velocities at these
grid points are  specified such that  the  flowrate and direction can be
easily  represented.

SUMMARY OF DATA

     Field  measurements of current velocities and water temperature have
been made for model calibration as well  as verification  purposes  for  the
hydrothermal  prediction.   Two field  trips were carried  out in the summer
and winter  respectively.   The procedures and  data results will  be discussed
briefly below.

June 1978 Data  Acquisition

     On June 19 and 20, a team was sent out to  the  field for measurement
of current velocity and temperature.  In coordination,  the  infrared  (IR)
scanner data was  obtained by  flights over the  channel.   The current
measurement and temperature  readings were obtained at 9 points; 4
                                   30

-------
points were located  at the north end, four others at the south end, and
one in the middle of the channel, as shown in  Figure 8.  At each point,
the measurements were done at different depths, namely  the water  surface
and some depths below the water surface.

     The in-situ measurements were carried out by  personnel on  three
boats.  The Barnes PRT-5 Radiometer was used to measure  the surface
temperature, and a  thermister (#9)  with floating mechanism  was used to
measure the temperature of the water surface.   Temperature profiles are
measured by lowering another thermister (#4) into  the water in increments
of 2 to 3 feet.  Current velocities are measured with a Bendix Model No.
665 w/readout current meter which  reads the magnitude and direction of
the current.  Figures 9 through  12 show the current velocities at 4 tidal
stages, which give a picture of the current patterns in  the  channel,
expecially at the boundaries.   One  may notice that the current is stronger
at the south boundary than at the north boundary.  Both boundaries can
exchange water with the Gulf of Mexico.

     The flights were made at 609.6 meter  or 2000  feet altitude.  The
black body of 1R  scanning window was  set at different ranges for each
flight set.  For the first flight,  the black  body range was 74° to 98°F,
the second was  75°  to 99°F, the third was 78°  to 102°F  and  the fourth
was 80° to 104°F.  All have a  full black body range of 24°F; it provides
a satisfactory resolution of 4°F for  each of the six colors between white
and black.  A  finer resolution was also obtained which reduced the  color
band  to 0.74°C  (1.33°F) temperature spread.   At 2000 feet  altitude, the
scanning width  is 941.8  meters or 3090  feet.  With  this scanning  width,
the whole stretch of the channel was  covered by ten east-west flights
as labeled in Figure  13.

     The in-situ measurement of water surface  temperature at the time
when  the airborne IR data was undertaken provides a calibration  of IR
temperature.  Figures 14-17 show the measured temperature  at  their
points at four tidal  stages of this mission.   Since the channel is  quite
shallow and only the surface temperature is of  primary concern in this
study,  only the surface  temperature is  presented.

January 1979 Data Acquisition

     Field measurements  by boat and  by IR scanning were schedualed
on January  30 and 31.  On previous experience, the current measurement
and temperature reading were  carried at different points from the pre-
vious  mission.   Figure 18 shows  the location  of all points where the data
were collected by three boats.  The route  for boat #1 was along  points
Ji' .§'  1' 1 and Jj' and for boat #2/  the  route was along  points 8,  7, jj,
W,  and JI  MeaTiwhile, boat #3 was working  on  the  region near" the coast
between  outlet and intake.  The water temperature and current velocity
were measured by the same instruments used in the previous data acqui-
sition.  The measurements were taken at the surface and successive
depths of 3 ft (1  m)  intervals.  The flights were coordinated at the same
                                   31

-------
time as for IR scanning, while the boats were collecting the data.   These
ten flight lines, shown in  Figure 13, were made at 609.6 meters or 2000
feet altitude.  The black body of IR scanning range was set at  44-80°F
so each of the  six color bands between  white and black would represent
a 6°F interval.

     On January 30,  the weather permitted both  morning and afternoon
data collection.  However, on January  31,  the conditions was so bad that
the mission had to be postponed to February 1 afternoon.  Figures 1 9
through 21 show the  current velocities at three tidal stages; namely,
flood tide, ebb tide and high tide respectively.   Figures 22 through 24
show the  surface temperature distribution  at corresponding  stages.-  To
further show the measured temperature fields,  isotherms were interpolated
from these in-situ data and  presented in Figure  25 through 27 correspond-
ing to each tidal stage respectively.

     In these two data acquisitions, although a plan of synchronized
measurements of current at tide changes was carried out in order  to
provide current data  at slack, flood and ebb, the effort was not so
successful due to technical difficulties  in obtaining reading and  other un-
expected circumstances.  Therefore, the current data could only be of
use as reference.

CALCULATION OF INPUTS

     The  important input parameters and some  specification of boundary
conditions, such as intake and discharge velocities, discharge temperature,
tidal condition,  river flow, surface heat transfer, and wind stress, will
be presented in this  section.

1.  Time step,  DT

         In order to  determine the time step,  DT,  the  stability  criterion
    has to be followed:

                            Dx    _     41760
    About  1/3 of this  value is reasonably safe  to use.

    Here, we use  DT  = 15 sec.

2.  Vertical eddy  viscosity,  A

         The vertical  eddy viscosity 5s  estimated by means of the  formula
    where  H is the local depth and C is an empirical constant.   This
    type of formula for horizontal  turbulent diffusion was originally
                                    32

-------
    suggested by Richardson  (1926).  Since that time, the 4/3 power has
    been substantiated by a considerable amount of empirical evidence and
    theoretical analysis.   The so-called constant C may vary with the in-
    tensity of the  turbulence  and is not well  established.  In this study,
    H is used as maximum depth and C is 0.002,  so  that we  have
                    Av = .002x(360)    = 6 cm2 /sec

    For shallow, well-mixed tidal water, about three times  the calculated
    value was found suitable.  Here, we use  A  = 20 cm2 /sec.

3.   Horizontal eddy diffusivity, B.

        The horizontal eddy diffusivity is calculated by the same formula
    as  mentioned above.  However,  the  maximum depth  is replaced by the
    maximum  length of the domain, which is  6 km in this study.   So we
    have

                Bh = .002x(600,000)It/3  = 100,000 cm2 /sec

4.   Vertical eddy diffusivity, B

        In this study, the turbulent Prandtl Number is assumed as  1.
    Thus, the vertical eddy  diffusivity is equal to the vertical eddy  vis-
    cosity or B  = A   =20 cm2 /sec.

5.   Surface heat transfer coefficient, K

        The procedures for  K  calculation are as follows:

    a.  T . = T  - (14.55 + 0.11VT )(1-f) - [(2.5 + 0. 007T  )(1-f)]3
         Q    o                  a                       3

       where Tj = dewpoint temperature, °F

              T  = air temperature, °F
               3

                f = relative  humidity in fraction  of  unit

    b.   6  = 0.255 - 0.0085Tave + 0. 000204Tave2

       where Tgve =  fl"s  +  Td) /2,  and B is  an  intermediate step

              T    = ambient surface temperature,  °F

    c.  f(u) = 70 + 0. 7u2

       where u =  wind speed, mph

    d.  KS  =  15.7 +  (p + 0. 26)f(u)

       where KS is the  surface heat transfer coefficient in BTU/ft2 day.


                                   33

-------
        The values of T , f,  u, T  needed for  this calculation are read
        from the climato&gical data3.

6.   Wind drag coefficient, C ,

         The wind  stresses on the free surface  are introduced into the
    mode! as


                   v 3ui     _            v 3 v i    _
                  pH 80'0=0  = Tx  "      pH 3~0'a=0 ~ Ty

    The  subscripts x and y indicate the shear stress acting in the x and
    y direction respectively.   The relation of these stresses to the wind
    speed at a certain height is very difficult to determine theoretically
    and its  value is usually based on semi-empirical formulae.   The well
    known form of  the  relationship between shear stress T and wind
    speed U  (usually measured at a height of ten meters] is

                               * = PaCdU2

    where p   is the air density and  C , is  a dimensionless drag coefficient.
    In this  study,  the drag coefficient  formulae obtained by Wu (1969)  is
    adopted, and its  value was given by

                     C . = 0.00125 U2   for  U <  1  m/sec
                      d
                        = 0.0005/CT   for  1 <  U < 15 m/sec

                        = 0.0026      for   U >  15 m/sec

7.   Intake and discharge velocities

         The  intake and discharge velocities are calculated according to
    the discharge flow rate from power plant data,  the grid size and the
    average depth  at the intake and discharge outlet.   The procedures
    are shown as follows:

    a.  Flowrate =  955,000 gpm  (from power plant physical data)
                 =  62.8 m3/sec
                 *  62. 8 x 10  cm3/sec

    b.  Both  intake and  discharge channels  are at 45° from N,
        therefore,

        31.4 x 10  cm3/sec is crossing  the Ax and Ay  at the  point of
        intake and  discharge.

    c.  The average  depth at intake and discharge outlet is approxi-
        mately 41 or  122 cm,  and the width is Ax = Ay  = 41760 cm;
        so the cross-sectional area  is 41760 x 122 cm2.


                                   34

-------
   d.  The average velocities:
               U    = V
                ave     ave
                               31.4 x  10
                          41760 x 122

e.  The velocity profiles are assumed as  shown.

         1

         2
                                       = 6.163 cm/sec
                                                     7.0 cm/sec
    Intake:   V3(14, 4,  k)  = 7 -  3 x

             U3(15, 3,  k)  = V3(14, 4,  k)

    Discharge:  V3(14, 8,  k) = 7 - 3 x
        k   3

            4

            5

    f.  To  allow for channel  storage during tide change we assume the
       intake and discharge velocities  to be sinusoidal, i.e.

                                                (EST -  7.625)]

                                               for k = 1,  2,  3,  4

                                               y^ (EST -  7.5)]

                   U3(15, 8, k) = -V3(14, 8, k) for k = 1,  2,  3,  4

       where 7. 625 and 7. 5 are taken  to be the values of the phase shift
       which  takes  into account the time to travel  from the south end of
       Anclote Key to the concerned point.

8.   Tidal condition  on June 19,  1978

    Simulated diurnal tide is shown in Figure 28,  where

    a.  Period = 12.5 hr

    b.  Stage  = short term average  sea level - MSL = 48 cm

    c.  Amplitude = £ short term average tide  range = 65 cm

    d.  Time shift = 7.125 hr
       i.e. at 7.125 am,  June 19,  1978, the tide at the south end of
       Anclote Key was zero.

    e.  W - E  lapse  = 0.014 hr/DX

       Wave propagation  speed  C =J 2gh  -J 2x980x360 = 850 cm /sec
       (H  = 360 cm is the maximum depth of the Anchorage.)
                                    35

-------
       The time  needed to travel one grid distance is

                          iff?- SO sec- 0.0,4 hr

       We use 0.014 hr per  DX  for phase shift  in W  -  E direction and
       the imposing tide at the  south entrance is

            n  =48 + 65 sin[~I- (EST - 7.125 - 0.014(1  -  1)]
              S               1 £•• 5
       1 = grid no. in W - E direction.

    f.  S - N  lapse = 0.15 hr

       Distance from south entrance to north entrance is  about  543000 cm.
       Time for wave to travel this distance is  54300°  =0.18 hr.
       We take  0.15 hr as phase difference betweln°the south and the
       north  boundaries; there,  the imposing tide at the  north entrance
       is

       n   =48+65 sin!—— (EST -  7.125  -  0.15 - 0.014(1  - 1)]
         n              1 it o
9.   Tidal condition  on January 30, 1979

        Simulated diurnal tide is shown in Figure 29.  The  calculation
    procedures are  of the same  as summer tidal conditions.

    a.  Period =  12.0 hr

    b.  Stage  = 36. 6 cm

    c.  Amplitude = i short-term average  tide range =  42. 7 cm

    d.  Time shift =  10 hr
       i.e. at 10 am, January  30, 1979, the tide at the south end of
       Anclote Key was zero.

    e.  W - E  lapse = 0.015 hr/DX

       This value  is slightly higher than the summer case since  the
       maximum  water depth in  the winter is less than the maximum water
       depth  in  the summer.  The imposing tide at the south entrance  is

            ns  = 36. 6 + 42.7 sin[j| {EST  - 10  - 0.015(1 - 1)]

    f.  S  -  N  lapse = 0.2 hr/DX

       This is the time for wave to travel from  south entrance to north
       entrance.   So the imposing tide at the north  entrance is

            TI  = 36.6 + 42.7 sin[|| (EST  - 10  - 0.2 -  0.015(1 - 1)]
              n                   12

                                   36

-------
10.   Anclote River flowrate and temperature

     a.  The distance traveled from South Anclote Key to Tarpon  Springs
        is 20 DX.  We  estimate a time lapse of 0.5 hr to account  for the
        retardation  due to buffering  effect of river storage and Anclote
        River's natural outflow.

     b.  The average current is estimated to be 20 cm/sec,
        therefore, we  take

               U3(16, 1,  k)  = 20 cosly— (EST  -  7.625)]

               V3(15, 1,  k)  =-20  coslj^L (EST - 7.625)]

               for k =  1,  2, 3, 4

     c.  The surface elevation at Tarpon Springs  is to be  calculated.

     d.  To be in accordance with given velocities at Tarpon Springs,
         the temperature there  is also assigned  and its value has  a
         24 hr period instead of 12.5 hr.   This  temperature on June
         19,  1978 is

               T3(15,  1, k) =26.9 + 0.5 sin[|| (EST -  12)]

         while  on  January 30, 1979 the temperature was assumed as

               T3(15,  1, k) =11.9+0.5 sin[|| (EST -  12)]

         where the 12  hr  shift  is to make the peak temperature occur
         at 1800.  Thus, the water in and  out at  Tarpon Springs has a
         temperature ranging from 26. 4 (before  dawn) to 27. 4 (late after-
         noon) in the summer and ranging  from  11.4 (before dawn) to
         12.4  (late afternoon)  in the  winter.

 11.  Discharge temperature

          On June 19-20,  1978, the recorded discharge temperature at
     daytime was  in the range of 29.3 -  30. 3°C, while on January 30-
     February  1,  1979,  the temperature range was 16. 4 - 15,2°C.  To
     account for the further drop of  discharge  temperature due to cooler
     ambient temperature at nighttime, we  assumed a sinusoidal variation
     of discharge temperature with diurnal period.

     a.  Discharge temperature is estimated for  June 1978

               T3(14,   8, k) =29.4+0.4 sin[|| (EST -  12)]

         Therefore, the highest discharge  temperature of  30. 3°C  happens
         at 6 pm and the  lowest 29. 3°C at  6 am.
                                    37

-------
     b.  Discharge temperature is  estimated for January 1979

               T3(14, 8,  k)  = 15. 8 + 0.6 sin[~ (EST  -  12)]

         The highest discharge temperature is  16. 4°C  at 6 pm and the
         lowest discharge temperature is 15. 2°C at 6 am.

12.  The Gulf  temperature

          The  Gulf water outside the Anclote Anchorage as well as the
     atmosphere is sink to the heat disposal from  the  power plant;
     therefore, the boundary conditions on temperature  at the north and
     south entrance are not considered  as adiabatic as in normal case of
     far-field thermal pollution problem.  Instead, we  specify the outside-
     Anchorage ambient temperatures.  Again,  they are  24 hr periodic
     and their values should  be  in accordance  with the measured  tempera-
     ture in the same neighborhood.  Here in compliance with measured
     data, we  use
                 Tab ~ 27*° + °'2 sintl  (EST ~

     for the summer simulation during June 19-20 I978 and  use


                 Tab = 11*8 + °'4 sin[15  (EST ~ 12)]

     for the winter simulation during January 30-February  1,  1979.

RESULTS

     Due to the fact that the archival data for Anclote site is  inadequate,
the following inputs are used for the computer runs:  data from NOAA
tide table,  solar radiation,  wind and  the power plant operating conditions.
Tables  1 and 2 show  the operating conditions of  Anclote Power Plant
during  June 19-20, 1978 and January 30-February 1, 1979 respectively.
The isotherms  obtained from the IR data from each flight are interpolated
by hand from mosaic digicolor films and then  plotted by the computer.
These isotherm plots  are presented in this section for easy comparison
with the predicted isotherms.  The average deviation of the calculated
temperatures from  IR temperatures is also indicated.   This  deviation is
calculated  by simply averaging the temperature differences  between the
measured and calculated temperatures at each point of the domain.  The
average deviation is given  by

               S2  = S[TB(i, j,  D - TlR(i, j)]2/Z (i,  j)
                   i'j                          i,j
where TB  is the calculated temperature, T1R  is the  IR temperature and
I(i, j)  is the number of surface half-grid points in  the domain.
U
     The T1R and  TB isotherms  are  then  compared to assess the accuracy
of the model in predicting the dispersion of the waste heat.  It  is recog-


                                     38

-------
nized that the principal factors  affecting the flow  pattern and the shape
of the thermal plume are the tide  driving through the north and  south
boundaries, wind effects on the water surface, bottom topography, heat
transfer through the air-water interface, and  the  intensity of eddy vis-
cosity and diffusivity.   The effects of each  of these factors on the flow
pattern have  been discussed in  detai! by Lee  et al.  (I978a, b).  These
effects  are important in understanding the numerical behavior of  the model,
and  play  a leading role in proving the capabilities of the model.  In this
report, however, only  the results of the verification runs are presented.
The results of the summer and winter verifications are discussed  in the
next  two  sections.

Summer Results

     The  figures from  the summer  simulation show the hydrothermal dis-
persion of waste heat under the conditions of  June 19-20, 1978, as de-
scribed in the previous section.  The simulation run started at 0400 EST,
June 18,  with initial  conditions of zero  velocity field, equilibrium  water
level and  constant ambient temperature  (so-called  "cool start").   The tide
and  heated discharge are then imposed.   After 20 hours of "warm up,"
the thermal plume can  be seen to  develop.   Experience indicates that the
initial condition is  not  important, as  its  effects die out in the  first one
or two  hours.  However, it  is convenient to begin a simulation through a
"cool start."

     Figure 30 shows the surface  flow pattern at 1030 EST,  June 20,
which corresponds to the high tide at south Anclote Key.   The current
velocity is relatively small and the mainstream flows in the  north  direc-
tion.   A recirculation can be seen to occur  as a result of the flow of
water from the Anchorage into the Anclote River.   The net effect is the
flow of part of the discharge into the river  entrance where the intake
structure is located.

     Figure 31 shows the resultant velocity  at high tide.   The result of
the velocity components, u and  w, is plotted on the vertical cross sections
at J = 4,  8 and 12.  It is believed that  the  currents at these  sections are
most affected  by the plan-form configuration,  bathymetry and  other fea-
tures such as the river mouth,  discharge and tidal boundary etc. The
u-w  velocity  profiles show clearly the effects of bottom friction in retard-
ing the flow.   Figure 32 shows the v-w  resultant velocity on the  vertical
cross sections at I = 4, 8 and 12.  In these two figures,  the vertical
velocity component w has been exaggerated  to make the vertical circula-
tion  of the current detect!ble.  One may notice that,  since  the velocity
vector at  the  top level  (a =  0}  is plotted right on the  water surface,
joining  the tails of these velocity vectors would show the free-surface
profile.

     Figure 33 shows the isotherms plotted from the IR data.  The com-
puter-plotted  points were obtained  through  visual  interpolation from the
mosaic IR image, as discussed earlier.   Due  to the difficulties  of  visual


                                   39

-------
interpolations and large grid  size, the computer-plotted isotherms  un-
avoidably differ somehow from the real isotherms in the digicolor film.
This figure shows that the heated water is recirculated back to  the intake
through  the  Anclote River.  This will cause an  increase in the input
temperature of the condenser cooling  water, thus  reducing the efficiency
of the power plant.

     Figure  32 shows the computed isotherms  at the same  tidal stage.
Comparing this figure with Figure 33, the good agreement of these
temperature  fields can  be easily seen.  The average deviation of the
simulated temperature from 1R temperature is  only  0. 359°C.

     Figures 30 through 34 form the first set of results of the
summer simulation.  This set  of plotted results, which  includes the
surface velocity, u-w  velocity,  v-w velocity,  IR temperature and the
calculated temperature, coincides with the high  tide at  south Anclote
Keys.  Figures  35 through 39 is the second set of results of the same
simulation run, at 1430 EST,  June 20, or at a time corresponding  to
maximum ebb tide at south Anclote Keys. Therefore,  this set shows
the current and temperature  fields of the Anchorage four  hours after
the first set of results was taken.  It is  also  seen that the  Anclote
River empties into the  Anclote basin while the basin drains  northward
and  southward.  There are strong ebb  flows occurring at  both the north
and  south end of the Anclote Keys barriers.  The dividing  line  lies along
the j = 10 grid line.  On the northern side  of this line, the water flows
north, while on the southern side,  it flows  south.   At  this  ridge  line,
the water has minimal transverse motion.  The location of  this Tine varies
with time as the  tides from both ends are not in phase.  In reality,  most
of the flow  in the Anclote Anchorage  comes  primarily from the south  end;
therefore,  the ridge line is observed  to  exist  at a  region  close to  the
northern end of the Anchorage.  Also,  this line is observed to move
southward  during ebb  tide and northward during  flood tide.  The flow
pattern corresponds with the observed  flowfield at a similar tidal stage.

     Figure  36 shows a strong  u-velocity component at the river mouth,
while Figure 37 shows  a strong v-ve!ocity component near the south
boundary.   Figures 38 and 39 are two corresponding isotherm plots at
this  tidal stage.   It can be seen that the thermal  plume moves according
to the stage of the tidal cycles, as would be expected.  The 27. 50°C
isotherm  from the simulated results covers a larger area than  that of the
IR plot.   Generally  speaking, the  simulated  results are in good agreement
with the IR data.

     Figures 40 through 45 is the third  set of results  of the simulation
run,  but at  1730 EST,  June 20, which corresponds to  low tide at  the
south Anclote Keys.  The current field is seen  to be generally quite  small
except at Anclote River region.  The river  is still undergoing its  outflow-
ing process and continues to  empty its  tidal storage into the basin, while
the open sea flows into the Anchorage through  both open  boundaries.
It has to be mentioned here that the  lower part of the  27. 50°C isotherm
from the IR  data is  not very  visible in  the IR digicolor film.  As  a result,


                                   40

-------
the contour shown in Figure  43 is estimated by  extrapolation  from  the
IR data.  The calculated isotherm at this  temperature,  as  shown in
Figure  44, has a peculiar shape at the same region.  This may be due
to the fact that this region, north of the  dredged ship channel leading
to Tarpon Springs, is relatively  shallow and becomes especially so during
low tidal  stage.  The sum effect is that the heat transport, due to con-
vection, becomes small.   This shallow region can be easily seen from
Figure  42 along the I =  8 cross section at a location about J = 4, 5 and
6.  If this isotherm is ignored, then the remaining isotherms,  shown in
Figures 43 and  44, have the  same sort of tongue-shape profiles, show-
ing some  degrees  of good agreement.

     Figures 45 through  49 is the last set of results of the summer
simulation.  The time is  2030 EST at the same day and  corresponds to
flood  tide at the south end of the Anchorage.   The incoming tide from the
south end drives northward into Anclote basin.   The dividing  fine dis-
cussed  earlier, does not appear,  as the northbound  tidal  current changes
course into  the northeast direction.  Both Figures 48 and  49 show that
the thermal  plume has been pushed back,  and the isotherms become more
compact.  These features could be explained by considering that the
changing  tide retards the convective transport of the thermal plume.
The 27. 50 and 28. 25°C  isotherms  in Figure 49 are seen to be  pushed
toward  the northeast direction rather than being pushed toward  the east
coast, as shown by the  corresponding IR  isotherms. Figure 48.  This
could be  because the north and south tidal conditions are  no  longer con-
sistent  with the real tidal conditions  after such  a long  time of simulation
(64.5 hours).  At this stage,  it  seems that the  northern tide could in-
duce  a  stronger current to push  the plume toward the  southeast direction.
If correct tidal data had  been provided, the results might have improved.

Winter Results

     The  winter simulation run started at  2200 EST, January  29, 1979.
The meteorological data and input conditions are as discussed in the pre-
vious section.  Here, we present  only two sets of results corresponding
to the successive  flood tidal and maximum  ebb tidal stages at  the Anchorage
on January  30, 1979.  These tidal stages were recorded during the  first
and second data acquisition missions  which took  place on January 30,  1979,
at about 1100 EST and 1600 EST  respectively.    It should  be pointed out
that,  due to sudden change in weather condition,  the second  and third
missions were about 48 hours apart.   In between these  missions, the sea
had become  very rough under very stormy conditions.  The effect of the
stormy weather can be clearly seen by comparing the two IR temperature
fields, Figures 58 and 60; the former was  taken at about  1600  EST,
January  30, while the latter was  taken at  about  1600 EST, Feburary 1.
There is therefore a 48 hour  time lag.  Figure  58 shows the thermal
plume at ebb tide  of a typical  winter day.   This is a great contrast to
Figure 60 which shows the plume during high  tide but  after the in-
fluence  of the stormy weather.  It is to^ be noted that the  isotherms
shown in  Figure 60 are  not of the same values as those shown  in Figure
58.   In  fact, there was a 4°C  drop in temperature, and a  2. 3°C  drop in


                                  41

-------
the discharged water temperature.  Particularly, the wind condition  had
gone through a very substantial change during  this 48 hour time gap.
This wind speed,  shown in Figure 60, was still  50% stronger than  that
shown  in Figure 58.  Therefore, it was no surprise to find that the
thermal plume, shown in Figure  60, was so compact and,shrunken.
The windy condition  and the  relatively low ambient temperature were
the apparent causes of this much dwarfed plume.   Observing the com-
pactness of the neighboring isotherms and judging the resolution power
of the  present grid system, one is tempted to suspect the ability of
modeling under this kind of extremely severe meteorological condition.
Experience shows that it is difficult to obtain such a compact thermal
plume through simulations  of  such adverse natural conditions.   Consider-
ing the roughness  of the grid system, this is actually not a surprise.
Therefore, for the winter  simulation only two sets of results, namely,
1100 EST and 1600 EST on January 30,  are  presented.

     Figures 50 through 54 is the first  set of results  of 1100 EST,
January 30,  1979, which corresponds  to the flood  tidal  stage at south
Anclote Key.   Figure 50 shows the surface-flow pattern from calculated
results.  The incoming tide from  the  south entrance drives directly  into
the Anclote River.   There is  a dividing line  across the Anchorage at
this tidal stage.  This  ridge  line is located along  the J = 9 grid line,
at which the incoming tides from  both north  and south meet.   This can
also be seen from the in-situ  measurement, as shown in Figure 19.  The
discharge at the outlet is  pushed head-on by this  flow.  As a result,
the thermal plume is squeezed by this flow.  The isotherm plots in
Figures 53 and  54  show this effect very clearly.  Although Figures  53
and 54 show similar tendencies,  the calculated isotherms seem not  to  be
in good agreement  with the IR isotherms;  particularly so for the 15. 40°C
isotherms.   However,  the  calculated plume areas are close to those ob-
tained  from the IR data.   Both temperature fields show the dispersion
of the  thermal  plume along the shore.  The recirculation of heated water
is clearly  indicated in both figures, showing that a recirculation actually
occurs at  this  tidal stage under the meteorological influence.

     Figures 55 through 59 show the  second  set of results of the  winter
simulation runs.  This  set of plotted  results  shows the current and
temperature fields  at 1600  EST,  January 30,  1979,  during ebb  tidal con-
ditions at the Anchorage.  Figure 55  shows  a strong current driving from
northeast  to southwest along  the channel.  No ridge line is observable,
and the discharged water  is not entrained by this main current; in  fact,
most of the discharged water  flows directly  through the river  mouth  into
the intake  of the power plant.  Thus, serious recirculation happens at
this time; this is also shown  in Figure 58 of the IR isotherms and in
Figure 59 of the calculated isotherms.  The calculated current pattern
shown  in Figure 55 is in good agreement with that of  the in-situ mea-
surements shown in Figure 20.  A similar agreement can be seen between
the calculated thermal plume shown in Figure 59 with that of the IR  data
shown  in Figure 58.  It can also  be seen that all  the calculated  isotherms
cover a slightly larger area than the  corresponding IR isotherms;  this
difference is shown to be  insignificant by a low deviation of 0.65°C  above
the IR measured temperature.


                                   42

-------
                            REFERENCES


Carter,  C. V.  The Hydrothermal Characteristics of Shallow Lakes.  Ph.D.
     Dissertation, University of Miami, Florida, 1977.

Dunn, W. E., Policastro, A.  J. and R.  A. Paddock.   Surface Thermal
     Plumes:   Evaluation of Mathematical Models for the Near and Com-
     plete Field.  Part One and Two,  Energy and Environmental Systems
     Division, Great Lakes Project,  Argonne National Laboratory, May
     1975.

Forsythe, G. E. and W.  R.  Wasow.   Finite-Difference Method for Partial
     Differential Equations.   J. Wiley  and Sons,  1960.

Freeman,  N.  G.,  Hale, A. M. and M.  B. Danard.  A Modified Sigma
     Equations Approach to  the Numerical Modeling of Great Lakes
     Hydrodynamics.  J. Geophysical  Research,  77(6), 1050-1060,  1972.

Hinewood, J.  B. and  I. G. Wallis.  Classification of Models of Tidal
     Waters.   J. Hyd. Div., ASCE, Vol. 101,  No. HY10,  1975.

Leendertse, J. J.  Aspects of a Computational Mode! for  Long-Period
     Water-Wave Propagation.  RM-5924-Pr, The Rand Corp.,  1967.

Leendertse, J. J., Alexander, R. C.  and S.  K. Liu.   A Three-Dimen-
     sional Model for  Estuaries and  Coastal Seas:  Volume I,  Principles
     of  Computation.  R-1417-OWRR,  The Rand Corp., December 1973.

Leendertse, J. J. and S. K. Liu.  Modeling of Three-Dimensional  Flows
     in  Estuaries.  Proc.  2nd Annual  Symposium of the Waterways, Har-
     bors and Coastal Engineering Division of ASCE, September 1975 or
     Symposium on Modeling Techniques - 75.

Lee, S. S. and S. Sengupta.  Three-Dimensional Thermal Pollution Models.
     NASA CR-154624, 1978(a).
     Volume  1 - Review of Mathematical  Formulations
     Volume  II - Rigid-Lid Models
     Volume  III - Free-Surface Models

Lee, S. S. and S. Sengupta.  Demonstration of Three-Dimensional Ther-
     mal Pollution Models at  Anclote Rivers, Florida, and  Lake  Keowee,
     South Carolina.  NASA  Project,  Mid-Term Report, Department of
     Mechanical Engineering, University of Miami,  Florida, December
     1978(b).
                                  43

-------
                            REFERENCES
Lee, S. S., Sengupta, S., Tuann, S.  Y.,  and C. R. Lee.  User's
     Manual for Three-Dimensional Free-Surface Model, NAS 10-9410.
     Final Report, August 1980.

Patridge, P. W.  and C. A. Brebbia.   Quadratic  Finite Elements in
     Shallow Water Problems.   Proc.  ASCE, Vol.  102, No.  HY9, Septem-
     ber,  1976.
     \
Phillips, N. A.  A Coordinate System Having Some Special  Advantages
     for Numerical Forecasting.  J. Meteorology.  14:184-185,  1957.

Platzman,  G. W.  A Numerical Computation of the Surge of 26 June 1954
     on Lake Michigan.  Geophysics.   6:407-438,  1958.

Platzman,  G. W.  The Dynamic Prediction of Wind Tides on Lake  Erie.
     American Meteoroligical Society, Vol. 4,  No. 26, September 1963.

Reid, R, 0. and B. R. Bodine.  Numerical Model for Storm Surges in
     Galveston Bay.  Proc. ASCE, Vol.  94, No.  WW1,  February 1968.

Richardson, L. F.  Atmospheric Diffusion Shown on  A Distance-
     Neighbor Graph.  Proc. Roy. Soc., 110:709, 1926.

Sengupta,  S.  and W.  J. Lick.  A  Numerical Model for Wind-Driven
     Circulation and Temperature Fields  in Lakes and Ponds.  FTAS/TR-
     74-98, Case Western  Reserve University, 1974(a).

Sengupta,  S., Veziroglu,  T.  N., Lee, S. S.  and N.  L. Weinberg.  The
     Application of Remote Sensing to Detecting  Thermal Pollution.  NAS
     10-8470,  Mid-Term Report, May 1975.

Sengupta,  S.  A Three-Dimensional Model for Closed Basins.  ASME Paper
     76-WA-HT21, 1976.

Sengupta,  S.,  Lee, S. S. and  H.  P.  Miller.   Three-Dimensional Numerical
     Investigations of Tide and Wind-Induced  Transport Processes in
     Biscayne Bay.  SEA  Grant Technical  Bulletin No. 39, University of
     Miami, 1978.

Sundermann, J.  A Three-DimensSona! Model  of a Homogeneous Estuary.
     Proc.  14th Coastal Engineering  Conference,  Vol. Ill, June 1974.

Waldrop, W. and R. C. Farmer.  Three-Dimensional  Flow and  Sediment
     Transport at River Mouths.   Coastal Studies Institute, Louisiana
     State  University, Baton Rouge, Technical Report No.  150, Septem-
     ber 1973.

                                44

-------
                            REFERENCES
Waldrop, W. and R.  C. Farmer.  T hree-Dimensional  Computation of
     Buoyant Plumes.  J. Geophysical  Research, Vol. 79, No.  9, 1974(a).

Waldrop, W. and R.  C. Farmer.  Thermal Plumes for industrial Cooling
     Water.   Proceedings of the 1971 Heat Transfer  and Fluid Mechanical
     Institute,  Davis,  L. R.  and R.  E. Wilson  (ed.), Stanford, California,
     Stanford University Press, June 1974(b).

Wang,  J. D.   Real-Time  Flow in Unstratified Shallow Water.  Proc. ASCE,
     Vol.  104,  No. WW1, February 1978.

Wu, J.  Wind Stress and Surface Roughness at Air-Sea Interface.  J.
     Geophysical Research.   74(2) : 444-455,  1969.
                                  45

-------
TABLE 1 CLIMATIC DATA FOR SUMMER RUN AT ANCLOTE ANCHORAGE
Date
6/18
•



















6/19


EST
04
05
06
07
08
09
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
01
02
03
(1)
Tair(°C)
23.3
23.3
23.3
23.3
25.0
26.7
28.3
30.0
31.1
32.2
31.6
31.1
30.6
30.0
29.4
28.8
28.3
27.2
26.1
25.0
23.9
22.2
21.7
21.1
(2)
Humidity
.90
.90
.90
.90
.84
.77
.70
.65
.60
.56
.57
.58
.59
.61
.64
.67
.70
.73
.77
.75
.74
.73
.76
.80
(3)
Wind
Speed
357.6
268.2
312.9
312.9
312.9
447.0
536.4
536.4
536.4
536.4
491.7
536.4
581.1
581.1
625.8
581.1
447.0
402.3
312.9
402.3
312.9
223.5
223.5
312.9
(4)
Wind
Direction
50..
50.
50.
50.
50.
70.
80.
80.
70.
90.
80.
80.
90.
80.
80.
80.
80.
80.
80.
SO.
80.
80.
60.
40.
. (5)
Solar
Radiation
0.0
0.0
0.0
0.05
0.40
0.75
1.05
1.40
1.60
1.70
1.60
1.50
1.30
1.10
0.'70
0.30
0.05
0.0
0.0
0.0
0.0
0.0
0.0
0.0
(6)
T
surf
26.4
26.4
26.5
26.5
26.7
26.8
26.9
27.0
27.0
27.1
27.1
27.2
27.2
27.2
27.1
27.0
26.9
26.9
26.8
26.8
26.8
26.8
26.8
26.8

-------
TABLE 1  (Continued)
Date
6/19




















EST
04
05
06
07
08
09
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
(1)
Tair(°C)
21.1
21.7
22.2
22.8
23.9
25.0
26.7
27.8
29.4
30.6
30.0
29.4
28.8
28.3
27.2
26.1
25.6
25.0
25.0
25.0
24.4
(2)
Humidity
.84
.82
.80
.79
.73
.67
.60
.56
.53
.50
.53
.57
.61
.65
.69
.74
.76
.78
.79
.79
.79
(3)
Wind
Speed
312.9
312.9
312.9
312.9
357.6
402.3
447.0
491.7
536.4
447.0
312.9
268.2
447.0
536.4
581.1
581.1
357.6
402.3
357.6
312.9
357.6
(4)
Wind
Direction
50.
50.
70.
70.
70.
80.
80.
90.
80.
80.
90.
80.
350.
350.
360.
350.
60.
120.
110.
110.
110.
(5)
Solar
Radiation
0.0
0.0
0.0
0.10
0.25
0.40
0.60
1.25
1.25
0.60
0.80
1.20
0.70
0.50
0.15
0.15
0.05
0.0
0.0
0.0
0.0
(6)
Tsurf
26.8
26.9
26.9
27.0
27.0
27.0
27.0
27.0
27.0
27.0
26.9
26.9
26.9
26.9
26.9
26. 8
26.8
26.7
26.7
26.6
26.5

-------
TABLE 1 (Continued!
Date
6/20




-














EST
01
02
03
04
05
06
07
03
09
10
11
12
13
14
15
16
17
18
19
20
(1)
Tair(°C)
24.4
23.9
24.4
25.0
25.0
25.0
25.0
26.1
27.2
28.8
28.9
28.4
30.0
30.6
30.6
30.0
30.0
29.4
28.8
28.0
(2)
Humidity
.79
.80
.81
.82
.82
.82
.82
.78
.74
.70
.66
.62
.59
.58
.57
.57
.56
.55
.53
.50
(3)
Wind
Speed
312.9
268.2
357.6
357.6
402.3
357.6
357.6
357.6
491.7
581.1
447.0
268.2
357.6
402.3
402.3
268.2
402.3
312.9
223.5
223.5
(4)
Wind
Direction
100.
100.
110..
110.
110.
110.
90. ,
90.
110.
110.
110.
110.
90.
90.
110.
100.
100.
110.
90.
90.
(5)
Solar
Radiation
0.0
0.0
0.0
0.0
0.0
0.0
0.05
0.20
0.60
0.30
0.40
0.40
0.60
0. 50
0.55
0.40
0.30
0.20
0.15
0.05
(6)
Tsurf
26.4
26.5
26.6
26.6
26.7
26.7
26.7
26.8
26.8
27.0
27.0
27.0
27.0
27.1
27.0
27.1
27.2
27.1
27.0
27.0

-------
TABLE 2  CLIMATIC DATA FOR WISTEB, RUN AT ANCLOTE ANCHORAGE
Date
1/29


1/30



















EST
22
23
24
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
(1)
Tair(°c>
8.8
8.8
8.8
8.3
7.7
7.2
6.7
6.1
5.5
5.0
6.7
8.3
10.0
11.6
13.8
15.5
17.2
18.3
19.4
18.8
18.3
17.2
16.1
(2)
Humidity
.98
.93
.89
.86
.87
.88
.89
.89
.89
.89
.75
.63
.54
.50
..47
.43
.39
.35
.31
.38
.45
.52
.59
(3)
Wind
Speed
350
350
350
350
350
325
300
280
270
250
280
310
350
375
400
400
400
375
350
310
280
250
230
(4)
Wind
Direction
320
320
320
330
340
350
360
10
20
40
60
80
100
140
180
220
220
220
210
200
170
140
190
(5)
Solar
Radiation
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.05
0.25
0.45
0.75
0.95
1.15
1.30
1.15
1.00
0.80
0.40
0.10
0.0
(6)
T
surf
11.8
11.8
11.8
11.7
11.7
11.7
11.7
11.7
11.7
11.7
11.7
11.8
11.9
12.0
12.1
12.2
12.4
12.2
12.1
12.1
12.0
12.0
11.9

-------
TABLE 2  (Continued)
Date




1/31






1/31











EST
21
22
23
24
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
(1)
Tair(
15.5
15.0
15.0
15.0
14.4
14.4
13.8
13.3
12.7
12.2
11.6
12.7
14.4
16.1
17.2
16.6
16.1
15.0
13.8
12.7
11.6
10.5
9.4
(2)
..Humidity
.66
.72
.70
.68
.66
.66
.65
.65
.77
.89
1.00
.97
.95
.93
.86
.79
.72
.74
.76
.77
.77
.77
.77
(3)
Wind
Speed
210
200
210
230
250
250
270
290
320
360
400
420
500
600
700
800
800
850
900
900
700
600
600
(4)
Wind
Direction
250
310
340
10
40
50
60
70
90
110
140
180
220
250
280
310
330
330
330
330
330
330
330
(5)
Solar
Radiation
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.05
0.15
0.25
0.35
0.45
0.55
0.50
0.50
0.40
0.30
0.15
0.05
(6)
T
surf
11.9
11.8
11.8
11.8
11.7
11.7
11.7
11.7
11.7
11.7
11.8
11.8
11.8
11.9
12.0
12.1
12.2
12.1
11.9
11.7
11.5
11.3
11.0
           50

-------
TABLE 2 (Continued)
Date





2/1

















EST
20
21
22
23
24
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
(1)
Tair(°c)
8.8
8.3
8.3
7.7
7.7
7.2
6.7
6.1
5.5
4.4
3.8
3.3
3.8
5.0
6.1
7.2
8.3
9.4
10.0
10.5
11.1
10.0
8.3
(2)
Humidity
.79
.81
.83
.80
.77
.74
.74
.73
.73
.73
.73
.73
.68
.58
.51
.45
.39
.33
.30
.32
.34
.40
.46
(3)
Wind
Speed
600
600
550
500
500
500
550
550
550
500
500
500
550
600
650
650
650
600
600
600
600
600
600
(4)
Wind
Direction
330
330
340
340
340
340
340
350
350
350
360
360
360
360
360
350
350
350
340
340
330
330
330
(5)
Solar
Radiation
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0,0
0.0
0.0
0.05
0,15
0.25
0.40
0.55
0.70
0.85
0.80
0.70
0.55
0.40
(6)
T
surf
10.7
10.4
10.0
10.0
10.0
10.0
10.0
9.9
9.9
9.8
9.8
9.7
9.8
9.9
10.0
10.1
10.2
10.2
10.1
10.1
10.0
10.0
10.0
        51

-------
                                           Jacksonville
                                                Daytona Beach
          Tallahassee
                            Gainesville
             Anclote
             Anchorage
Anclote R.
Tarpon Springs

Tampa
                  Anclote
                   Key
                                     Lake
                                     Okeechobe
                                        Ft.
                                        Myers  Ft
                                             Lauderdale
                                                     Melbourne
  GULF OF MEXICO
                       .JWest
                         Palm
                         Beach
                                                        *) Miami
Figure 1.  Anclote Anchorage location in the state of Florida
                                52

-------
U7
U)
                                                   Q>

                                                                                            km
                                                           rage

-------
_ z-n(x,y,t)
  h(x,y) + n(x,y ,
                           h(x,y)
        Fiaure  3.   Definition  sketch of a-coordinate
                              54

-------
                               All Point Have
                               Same i
i=
                                                                   All Points
                                                                   Have Same j
      Figure  4.   Grid arranqement in the horizontal projection
                •  (full-grid) depth  (h)  point;  — u velocity point;

                 i v  velocity point;  +  (half-grid) w,ft,p,T,p and  n point
                                  55

-------
         v(i,j,5)=0
                                                  ,J, 5)=0
Fiqure 5.  Four cells in a vertical  column with velocities
           shown at definition point, and scalar variables
           at the center of cell  f»l
                             56


-------
                                           (i)
Half-grid
point
                                                             Full-grid  point
 .  Fiqure  6.   Notations (left) and variables (right) used in calculation
               of.T (#), n  (-»),  u (>*) and v f^)"  within  fi,i)-block
               (cross-hatched)
                                      57

-------
Ul
00
                              Figure 7.  Grid work for the Anclote Anchorage

-------
Ul
to
                     Figure 8.   Location of stations  for in-situ measurement, June 1978

-------
en
o
                                                                                       Discharge Rates 994,000gpm

                                                                                       Wind Velocity: 12.5 mph

                                                                                                   (155°)
                        Figure 9.  Velocity  from in-sttu measurement at  1710-1903, June  19,  1978

-------
          4    56    7   8   9   10    II
                                                            Discharge  —	•> l m
                                                                      0 10  20cm/sec
                                                               Discharge Rate:  994,000 gpra
                                                               Wind Velocity: 8.0 mph(90°)
Figure 10.  Velocity  from in-situ  measurement at 06H8-0812, June  20,  1978

-------
en
ro
                                                                                     Discharge Race: 994,000gpm
                                                                                     Wind Velocity: 8.0 mph(110}
                     Figure  11.   Velocity from  in-situ  measurement at 1125-1245,  June 20, 1978

-------
O>
to
                                                                                     Discharge Rate: 994,000 gpm
                                                                                     Wind Velocity: 7.5
                     Figure 12.  Velocity  from in-situ measurement at 1150-1605, June 20.  1978

-------
Ol
              N
                     2000 feet altitude
                     3090 flight line overlap
                        5 feet resolution
                                                   FL  1
                                                   FL  2
                                                   FL  3
                                                   FL  4
                                                   FL 5
                                                   FL 6
                                                   PL 7
                                                   FL 8
 t
  ANCLOTE

ANCHORAGE

     1 km
                                                   FL 9
                                                   FL 10
                     Figure  13.  Daytime flight lines on June  19 and 20,  1978

-------
U1
                                                                               IDischarge Temp: 30.3°C
                                                                               Wind Velocity: 12.5 mph(355^
                   Fiqure  11.   Surface temperature in degree C from in-situ measurement at
                                1710-1903, June 19, 1978

-------
en
01

                   u\
                              26.9
                                      26.9
                                      27.J
27.2
        26.8
           t
                                               27.1
                                                   26.9
                                                                26.7
                                                                                            Tidal Cycle

                                                                                 'Intake

                                                                                     j Discharge Temp: 29.3°C
                                                                                 V-S~| wltld Velocity: 8.0 mph(90°)
                                                          10
                     Fiaure  15.   Surface temperature in degree  C from in-situ measurement  at

                                   0648-0812, June 20, 1978

-------
13	-,
II
 3 1	
           27.3
                   27.0
                   27.3
27.7
        27.1
27.6
                                 27. •
                  26.9
                                    9   tO   II   ,2 ~~I3~  ,-Jl4   IS
                                              jDischarge Temp:  29.6°C
                                               Wind Velocity: 8.0 nmh(110°
Figure 16.  Surface temperature in degree C from in-situ measurement at
              1125-1215,  June  20, 1978

-------
01
CO
                                                                              [Discharge
                                                                                 |  Discharge Temp: 29.6 C
                                                                                   Hind Velocity: 7.5 mph
                                                                                              (100°)
                            34    *    l>    7   89   10   11
                   Figure 17.  Surface temperature in degree C  from in-situ measurement at
                                 11(50-1605, June 20, 1978

-------
                                  1 K
                                                 5 »
en
10
                                 lx
                                            14
                                             H

                                            5x
IK
                                                    13
    6x
                                                    NA-5
                                                          n

                                                         12
                                                         9«
                                                            17
                                                             Ix
                                                   9   10   I)   12   13
                    Figure  18.   Location of stations for in-situ measurement, January  1979

-------
                                                                 Discharge: 994,000 gpro
                                                                 Wind Velocity: 4.0 mph(320°)
Figure  19.   Velocity from  in-situ measurement at  1020-1340, January  30,  1979

-------
                                                                  Discharge: 994,000 gpm
                                                                  Wind Velocity: 6.0 mph
Fiqure  20.   Velocity from in-situ measurement at  U40-1800, January  30,  1979

-------
                                                                  Discharge  	> 2 •
                                                                           0   10  20cm/sec
                                                               J.	.         Tidal Cycle



                                                                    Discharge: 994,000 gpm
                                                                         Velocity: 14.0 mph
                        6    7   S    9   10   II   12  13  r->\*  IS
Figure 21.   Velocity from  in-situ measurement  at  1430-1640, February 1,  1979

-------
                  12.1
                                  12.1
                                                   12.3
             12.6
                               12.3
                                        11.8
                                                 12.0
                                                13.6
                      14.0*
                                                     14.7
                                                               I
                                            13.1
                                   11.9
                                        12.2
                                            13.8
                                                         	r_

             11.8
                           11.8
12.0
                                   12.0
                                    11.9
              12.2
                               11.9
                                                             ^Discharge
                                                         	4-N--I
                                                                                  ikm
                                                                 !       Tidal Cycle

                                                                 I  Discharge Temp: 15.3°C
                                                              v 	j  Wind Velocity: 4.0 roph
                                                               N.I               (320°C)
                 56    789   1O
Figure 22.  Surface  temperature  in degree C from in-sltu measurement at
              1020-1310, January  30,  1979

-------
                       7   t   9  10   II   |2  13  j~>\4  IS
Figure 23.  Surface temperature in degree C from in-situ measurement at
            1 WO-1800, January  30,  1979

-------
en
                    Figure  21.  Surface temperature in  degree C from in-situ measurement at
                               1130-1610,  February 1, 1979

-------
Figure 25.   Surface temperature from in-situ measurement at
            1020-1310, January 30, 1979

-------
Figure 26.  Surface temperature from in-situ measurement at
            1 MO-1800,  January 30,  1979

-------
00
                          Fi«;jure 27.  Surface temperature from In-situ  measurement at
                                      U30-16UO,  February 1, 1979

-------
                -Tide at south end of Anclote Key
                                                         2 TT
                -Simulated tide for calculation  n= 48 +  65 x sin l~  (H>))
                                          -Average Level
Figure 28.  Semidiurnal tide  for June  19-20,  1978 at south end of Anclote Key

-------
CO
o
                                            Tide at south end of Anclote Key
                                                  simulated tide for calculation n=36.6+42.7xsin[~- ^t-^
              Figure  29.  Semidiurnal tide for January 30-February 1,  1980 at south end  of Anclote" Key

-------
                   TIHECJUNE 2Q.197S)>
                   HINO SPEEO(CH/SEC)«
                   MIND OIRECTiaN(OEG/N)»
                   BIR TEMPERBTUREtDEO-Oi
                   OISCHflRGE TEHPtOEG-Oi
                   OI3CH FLOHRflTECCUM/SEC) i
                   LENGTH SCflLEClCMs  X CHJi
                   VELOCITY SCflLEtCM/SEC)>
10. S
450.0
110.
28.8
29. S
62.7
41013.
52.49
k 4
I 4
% h
. 4
' '
. t
1
• * t f f » , . .
• »/*«»..,
t >>«,,,,.
>»»«...* ^
>'<«*»..—.
'»•«».» 4
'^<««444
                                                          \

                                                          \
                                                              »   «•
                   HIGH TIDE
Figure  30.   Surface velocity, Anclote Anchorage  by modeiina
                                 81

-------
                TIHEtJUNE 20.1978)1
                HIND  SPEEDtCM/SEC)i
                MIND  DIRECTIONOEG/N)i
                HIR TEHPERRTUREtOEG-CJi
                OISCHflRGE TEHP(OEG-CJ>
                OI3CH FLOHRBTEtCUM/SEC)i
                LENGTH SCRLEtlCflr X CM J»
                VELOCITY SCflLErCM/SECJ«
10.S
450.0
UQ.
28.8
29.5
62.7
41019.
52.43
                                                             J=  8
                                                             J=  -4
                HIGH TIDE
Figure  31.   UW velocity, Anclote Anchorage by modeling
                                82

-------
                 TIHEUUNE 2Q.1978)>
                 WIND SPEEOlCH/SECJi
                 MIND DIRECTIONOEG/N)i
                 flIR TEHPERflTUREtOEG-C)i
                 OISCHflRCE TEHPCOEO-CJi
                 OI3CH FLQHRflTE(CU»/3ECJi
                 LENGTH SCRLEtlCtts  X  CfM i
                 VELOCITY SCRLECCM/SEC)«
10.5
450.0
110.
Z8.S
29.5
62.7
41019.
52.49
                                                              1= 12
                                                              1= 8
                                                              1=  4
                 HIGH TIDE
Figure  32.   VW velocity, Anclote Anchorage by  modeling
                              83

-------
  TIHEtJUNE 20.1978J*
  HIND SPEEOtCM/SEC)i
  HIND OIRECTIONCOEO/N)«
  RIR TEMPERRTUREtDEC-CJi
  DISCHARGE TEf1P(D£0-C)i
  OISCH FLQMRflTEfCUM/SEC)»
  LEMOTH SCRLE(1CM= X Cf1)»
  VELOCITY SCflLEfCM/SECJi
11.5
420-0
110.
29.1
29.6
6Z.7
41019,
52.49
   HIGH TIDE
Figure 33.  Temperature  from IR
                  84

-------
                  TIJ1EUUNE 20.1978)1
                  HINO SPEED!CM/SECJ>
                  HIND OIRECTI8N(OEG/N)t
                  BIR TEWERflTURElOEG-CJi
                  DISCHflRQE TEHP(OEG-C)»
                  DI3CH FL0HRaTEtCUH/SEC)i
                  LENGTH SCRLEUCfls X  CHli
                  VELOCITY SCHLEtCM/SEC)»
10. S
450.0
110.
28.8
29.5
62.7
41019.
52.49
                   OEVIHTI9N FR8M IR TEMP«
                   HIGH  TIDE
0.359
Figure  34.   Surface temperature, Anclote Anchorage by modeling
                                  85

-------
                  TIME(JUNE 20.1978)1
                  HIND SPEED(CM/SEC)i
                  HIND OIRECTI8N(DEG/N)»
                  flIR TEHPERF1TURE£ DEO-C )i
                  DISCHARGE TEHP(OEG-C)i
                  OISCH FLQHRBTECCU«/3EC)i
                  LENGTH 3CflLE(lCM= X Ctt)>
                  VELOCITY SCHLEtCn/SECJi
14.5
400.0
110.
30.6
30.0
62.7
41019.
52.49








*
-







» - . / ^^^.----.N
•»«/X^^--*-^
' » J ///---**•
• M / / ' 	
' / / / / / " -" • 	 "~ N
4 / / / / ' ^ ^ " ' N
* * *
\ i, ',',', ','.:::
. . ,




» > %
» ^^
-

*

\ \
N \
                 /   '
                  EBB  TIDE
Figure  35.   Surface velocity, Anclote Anchorage by  modelinq
                                 86

-------
              TIHECJUNE 20.1978Ji
              MIND SPEECHCtl/SEC) i.
              HIND DIRECTIBWOEG/NJi
              HIR TEJIPERflTURECOEOrOi
              OI3CHHROE TEf1P{OEG-CJ«
              D13CH FLOHRflTECCUM/SEC)»
              LENGTH  SCflUEUCtls-^CfllJ
              VELQCITY SCflLE(CM/SEC J a
14.5
400.0
110.
30.6
30.0
62.7
41019.
52.49
                                                            Ja 12
                                                            Js 6
                                                            J= 4
               EBS   TIDE
Fiqure 36.  UW velocity,  Anclote  Anchorage  by modeling
                              87

-------
                TIMECJUNE 20.1978Jt
                HIND SPEEOCCM/SEC)i
                MIND OlRECTIQN(OEG/N)i
                HIR TEflPERflTURECOEG-CJi
                OISCHRRGE TEHPCOEG-CJ*
                OI3CH FL3HRflTECUm/3EC)>
                CENCTH SCflLEncn=  X  C«JJ
                VELOCITY SCflLECCM/SEC)»
14.5
400.0
110.
30.6
30.0
62.7
41019.
52.49
                                                             1= 12
                                                             1= 8
                                                             1= 4
                 EBB   TIDE
Figure  37.   VVV velocity, Anclote Anchorage by  modeling
                               88

-------
 TIHEUUNE 20.1978J>
 HIND SPEED(Cfl/SECn
 HINO OIRECTI3N{OEG/N)i
 RIR TEflPERflTUREtOEG-CJi
 OISCHflRGE TEflPfOEO-CJi
 OISCH FLQHRRTEICUH/SECJi
 LENGTH SCflLEriCHs X CflJi
 VEL3CITY SCflLECCM/SEC)i
15.0
400.0
110.
30.6
29.S
62.7
41019.
52.49
 EBB  TIDE
Figure 38.  Temperature from  !R
                 89

-------
               7IMEUUNE  20,1978)1
               HIND SPEEOtCtl/SECJi
               HIND OIRECTIONCOEG/HJi
               RIR TEflPERflTURECQEG-CJj
               OISCHRRGE  TEHP{DEO-C)i
               OISCH FLQHRRTE(CUM/SEC)»
               LENGTH SCflLEdCtls X Ct1)»
               VEtaCITTf SCflLEtCII/SEOi
14. S
400.0
no.
30.S
30.0
6Z.7
41019.
52.49
                                                                so
               OEVIflTION FRQfl IR  TEHP t
               ES8  TIDE
0.361
Figure  39.   Surface temperature,  Anclote Anchorage by modeling
                              90

-------
                 TIt1E(JUNE 20.1978)1
                 HIND SPEED CCt1/3ECJ»
                 HIND OIRECTIQNCOEG/N)!
                 RIR TEHPERflTUREtOEG-CJi
                 OI3CHRRCE TEHPCOEG-CJi
                 DI3CH FL8HRRTEtCUH/3EC)i
                 LENGTH SCRLE(1C«= X  CM) i
                 VELOCITY SCRLEt CM/SEC )i
17.5
310. Q
110.
29.4
30.2
62.7  ;
41019.
52.49
                                                    ..   »    *

                 »»»»
                                                      \
                                                               \
                                                           N\
                 L8H  TIDE
Figure 40.   Surface velocity, Anclote Anchorage by   Modeling
                                  91

-------
 TIMEfJUNE 2Q.1978}»
 HIND SPEEDtCn/3EC)i
 HIND OIRECTIQNdJEG/NJi
 fl!R TEMPERRTURECQEG-Cli
 DISCHARGE TEMP(DEG-CJt
 OI3CH FLQHRRTEtCUH/SECli
 LENGTH SCRLEUCris  X  CM It
 VELQCITY SCflLEtCM/SEC)»
                                             17.5
                                             310.0
                                             no.
                                             23.4
                                             30. Z
                                             SZ.7
                                             41019.
                                             52.49
                                              J= 1Z
                                              J= 8
                                              Js 4
  LGH  TIDE
Figure
UW velocity, Anclote Anchorage  by modeling
                 92

-------
              TIMEUUNE 20.1978)1
              HIND SPEED!CM/SEC)«
              HIND OIRECTIQN(DEG/N3i
              RIR TEMPERflTURE(DEC-C)i
              OISCHfiRGE TEMPJDEO-CJi
              OISCH FL3HRaTE(CUH/3ECJ>
              LENGTH SCflLE(lCh= X  Cf1)i
              VELOCITY SCHLEICH/SEC)«
17.5
310.0
110.
29.4
30.2
S2.7
41019.
52.49
                                                           1= 12
                                                           1= 8
                                                           1= 4
              LOU  TIDE
Figure  42.   VW velocity, Anclote Anchorage by modeling
                              93

-------
TIME!JUNE 20.19783.
HIND SPEEDCCn/SEOi
WIND DIRECTIQN(DE&/NJ>
HIR T£f1PERRTUREtOEC-C)j
OlSCHflRGE TEMPtOEO-CJi
DI3CH FL9HRHTEtrUH/3EC)»
LENGTH SCRLEUCtts X  CHJ»
VELOCITY SCflLE(Cf1/SEC3>
17.S
310- 0
110.
29.4
30.2
62.7
41019.
52.49
 L9H  TIDE
  Figure  43.   Temperature from IR

-------
                TIflEtJUNE 20.1378)1
                HINO  SPEEOCCM/SEC)*
                HIND  OIRECTIQNCOEG/NJi
                HIR TEMPERRTUREtOEG-C)»
                OISCHflROE TEMP(DEG-C)i
                OI3CH FL8HRflTE(C'Jtt/SECJ>
                LENGTH SCRLEI1CH= X CH)>
                VELOCITY SCflLECCJ1/SEC)«
17.5
310.0
110.
29.4
30.2
62.7
41019.
52.49
                                                                 so
                OEVIRTIQN  FROM  IR TEMP*
                LOW  TIDE
0.538
Figure 44.  Surface temperature, Anclote Anchorage by modeling
                              95

-------




.



»
-


TIHECJUNE ZQ,1978)i 20.5
HIND SPEED (Cfl/SECJi 225. 0
HIND QIRECTI9N(DEG/M)» 90-
RIR TEHPERRTUREtOEO-C)i 28.0
OI3CHRROE TEf1PtOEG-C)i 30.1
OI3CH FL8HRRTEfCUH/3EC)J 62.7
LENGTH SCRLEdCtls X CHJ« 41019.
VELQCITY SCRLEt CM/SEC J» 52.49

, , % ,, .. — ^ * * *

•• « i t t s 
f f , 	 .
lilt//,*..
/* t
////,,..
lit//,*...


* * 4

* • *
» * .
» *
b

i
. \
\_- ,
\
N \
          FLQO TIDE
Flqure 45.  Surface velocity, Anclote Anchorage by modeling
                         96

-------
         TIflECJUNE 20.1978)1
         HINO SPEEOtCM/SEC)«
         HIND OIRECTiamOEG/NJi
         flIR TEMPERRTURE(OEO-C)«
         OI3CHRRSE TEMPI OEG-Oi
         OI3CH FL8HRFITEtCUf1/3EC)»
         LENGTH SCRLEtlCMs  X  CM) i
         VELQCITY SCHLEtCM/SEOi
20.5
225.0
90.
28.0
30.1
62.7
41019.
52.49
                                                       = a
                                                      J= 4
         FLQD TIDE
Figure 46.  UW velocity,  Anclote Anchorage by modeling
                          97

-------
          TIHECJUNE 20.1978)1
          HIND SPEEDCCn/SEC)i
          HIND OIRECTIQN(OE&/N)i
          RIR TEllPERRTURECDEO-CJi
          OISCHflROE TEMP(OEG-C)J
          OISCH FLQHRRTEtCUM/3EC)i
          LENGTH SCRLEUCf1= X  CflJi
          VELOCITY SCRLECC«/SEC)«
20.5
225.0
30.
28.0
30.1
62.7
41019.
52.49
                                                       1= 12
                                                       Is 8
                                                       la 4
           FLOO  TIDE
Figure 47.  VW velocity,  Anclote Anchorage by modeling
                           98

-------
 TIHECJUNE 20.1978)»
 WIND SPEEOCCtt/SECJi
 HIND OIRECTIQNCDEG/N)!
 RIR TEMPERflTUREIDEC-CJi
 OISCHflRGE T£f1P(DEG-C)i
 OISCH FLQHRflTEfCUM/SECJ»
 LENGTH scm.Eticri= x cnj«
 VELOCITY SCflLECCM/SEC)i
20.0
225.0
90.
28.0
29.3
62.7
41019.
52.49
  FLQD TIDE
Figure 48.  Temperature  from IR
                  99

-------
                 T1HEUUNE 20.1978]i
                 HIND  3PEEOCCM/SEC) i
                 HINO  DIRECTION!OEG/N)i
                 flIR TEHPERRTUREfDEG-CJi
                 OISCHflRGE TEMP(OEO-C)»
                 OI3CH FL3HRRTEtCU«/3EC)i
                 LENGTH SCRLEllCns X CMJi
                 VELOCITY SCflLECCM/SECli
20. S
225.Q
90.
28.0
30.1
62.7
41019.
52.49
                 OEVIflTION FRQM IR  TEMP>
                 FLQO TIDE
 0.361
Figure  49.   Surface temperature, Anclote  Anchorage by modeling
                                TOO

-------
                TIMEfJBNURRY  30.1970) :
                HIMO SPEED(CN/SEOi
                NINO OIRECTlQNfCEC/NJ:
                flIR TEMPERHTUREf DEC-CM
                01SCHRRCE TEflPCDEG-C)!
                DISCH FLQHRflTEfCUM/SECJ»
                LENGTH 5CRLEriCn=  X  CMJ>
                VELOCITY SCflLEfCH/SECJ»
11.0
4QQ.Q
320.
12.5
15.6
62.7
41013.
52.49
                \   \  \

                \   \   \

                «   \   \
             !  /
                «   /   s

                f   /   /

                /  /  /

                /  /  /
    \
                                                       \
            \
                FLQO TIDE
Figure 50.   Surface velocity, Anclote Anchorage  by modeling
                                 101

-------
          TIflEf JflNUflRY 30.1979)>
          HIND SPEEOrCM/SEC)>
          HIND 01RECnON(OEG/N)»
          HIR TE«PERflTUREfOEO-C)«
          OlSCHflRGE 7EflPCDEG-C)>
          01SCH FLOHRRTEfCUM/SECJ«
          LENGTH SCRLEflCff= X  CH)»
          VELOCITY SCRLEfCn/SECJ«
11.0
400 .0
320.
12.5
15- S
G2.7
41013.
52.49
                                                       J= 12
                                                        J= 8
                                                        J= 4
           FLQO TIOE
Figure 51.  UW velocity,  Anclote Anchorage by modelina
                          102

-------
              T'lEtJflNUBRY 30.1979)>
              HIND SPEEOfCM/SECJi
              HIND DIRECTIQNfOEO/N).
              BIR TEflPERfiTUREfOEO-CJs
              OISCHflROE TEHP(OEG-CJ«
              OISCH FLOHRfUEf CUM/SEC J«
              LENGTH SCBLEflCNr X CHJ»
              VELOCITY SCRLEfCM/SECJ«
11.0
4QQ.O
320-
12-5
15.6
62.7
41Q13,
52.49
                                                          1= 12
                                                           1= 8
                                                           1= 4
              FLOD TIDE
Ficure  52.   VW velocity, Anclote Anchoraqe by  modeling
                          103

-------
 Tlft£{JRNUflRY  30.1979}!
 HiNC  SPEEO(Cn/SEC)J
 HIND  OIRECT1QN(OEG/NJJ
 flIR TEHPERflTUREfOEG-C)!
 QlSCHflRCE  TEflPtOEG-CJ:
 OiSCH FLQHRRTEfCUM/SECH
 LENGTH SCRLEHCHr X CM J s
 VELOCITY SCfiLEfCfl/SECJs
It .5
400.0
320.
12.7
IS.5
62.7
41019.
52.49
 PLOD TIDE
Figure  53.   Temperature  from 1R
                  104

-------
             11HECJBNURRY 30.1379)>
             H1NO SPEEDfCM/SEC)»
             H1NO DIRECT1QN(OEO/N)»
             R1R TEf1PERRTURE(CEC-C)»
             01SCHHRCE TEMPCDEO-C)>
             OISCH FLQHRflTE(CUM/SEC)i
             LENGTH SCflLEftCf1= X CM)
             VELOCITY SCRLEfCM/SEC).
11.0
400.0
320.
12.5
15-6
G2.7
41019.
52.49
              OEYIRT1QN FROM IR TEMP*
              FLQO TIDE
 0.742
Figure  54.   Surface temperature, Anclote Anchorage  by modeling
                              105

-------
         TlftEfJflNUflRY 30.1979)1
         HIND SPEED fCft/3EC)»
         HINO DIREC71QNtDEO/NJ>
         ftlR TErtP£RRtUREfOEO-C)«
         01SCHRRGE TEf1PCOEG-C}«
         01SCH  FLQHRRTEfCUM/SECJ»
         LENGTH SCRlEflCfls X CM)»
         VELOCITY 5CfllEfCft/SECJt
16.0
400.0
360.
14.0
1G.3
52.7
41019-
52.49
           /  /  /
        1  It//'-'   '
Figure  55.  Surface velocity, Anclote Anchorage by modelina
                       106

-------
           T1HEURNURRY 30.19791*
           HIND SPEEOfCM/SECJ»
           WIND DlRECTIQNfOEO/NJ»
           R1R TEHPERflTUREfOEC-CJ>
           DI5CHRROE T£HP(OEO-CJi
           OISCH FLQHRRTEfCUM/SECJ»
           LENGTH SCRLEflCNr  X  CMJ»
           VELOCITY SCRLEfCH/SECJ»
16.0
400.0
360.
14.0
16.3
62.7
41019.
52.49
                                                        Ja  12
                                                        J= 4
           EBB  TIDE
Figure  56.  UW velocity, Anclote Anchorage by  modeling
                           107

-------
          nflEfJflNURRf 30,1379)>
          HIND SPEEDfCfl/SECJi
          H1NO DIRECTIQN(OEQ/N)»
          flIR 7EflPERfnUREfOEC-C)»
          DISCHARGE TEMP(OEO-C)>
          01SCH FLQ«RflTEfCU«/SECJ«
          LENOTH SCflLEflCf1= X CM)i
          VELOCITY 5Cm.EfCH/SECJ«
16.0
400.0
360.
14.0
16.3
62.7
41013.
52.49
                                                       1= 12
                                                        r  4
          EBB  710E
Figure  57.  VW velocity, Anclote Anchorage by modelinq
                          108

-------
 TIMEURNURRY 30.1979) i
 HIND SPEEOfCM/SEC )>
 HIND DIRECTIQNfOEO/NJ:
 RIR TEflPERRTUREfOEG-C):
 01SCHRRGE  TEf1PtDED-CJ»
 DISCH FLOHRRTEfCUM/SECJs
 LENGTH SCRLEflCM= X CMJ:
 VELOCITY SCRLEfCfl/SECJs
17.0
400.0
360.
18.0
16.3
62.7
41019.
52.43
 £83  TIDE
Figure 58.  Temperature  from 1R
                  109

-------
              MIND SPEEQfC.fi/SECJi
              HIND OIRECTION(QEO/N)«
              HIR TEf1PERfnUREfDEC-C)J
              OISCHBROE TEMPfOEO-Oi
              01SCH FLQHRflTEfCUf1/SECJ«
              LENGTH SCflLEflCMs X  CH)«
              VELOCITY SCflLEfCM/SEC It
16.0
400.0
360.
14.0
IG.3
82.7
41013.
52.43
               OEVlflUQN PROM IR
               £88   TIDE
0.350
Figure  59.   Surface temperature, Anclote Anchorage by  modeling
                               110

-------
TlflEClST   FEB.  .1979Ji
HINO SPEEOtCtl/SECU
WIND OIRECTI9NCOEG/N}:
flIR TEflPERHTURE(OEG-C)«
DISCHflRGE TEHPtOEG-CJi
DISCH FLQHRRTECCUM/SEC) :
UENGTH SCfllEClCMr X CM) :
VELOCITY 5CRLEICM/SEC )J
16.0
600.
330.0
11-1
U.O
62.7
41019.
52.49
HIGH TIDE
Figure 60.  Temperature from !R
                  111

-------
TECHNICAL REPORT DATA
(Please read Inunctions on the reverie before completing)
1. REPORT NO. 2.
EPA-600/7-82-037a
4.T.TLEANDSUBT.TLE Verification and Transfer of
Thermal Pollution Model; Volume I. Verifica-
tion of Three-dimensional Free-surface Model
7.AUTHOR(s> s.S.Lee, S.Sengupta, S.Y.Tuann, and
C.R.Lee
9. PERFORMING ORGANIZATION NAME AND ADDRESS
The University of Miami
Department of Mechanical Engineering
P.O. Box 248294
Coral Gables, Florida 33124
12. SPONSORING AGENCY NAME AND ADDRESS
EPA, Office of Research and Development
Industrial Environmental Research Laboratory
Research Triangle Park, NC 27711
3. RECIPIENT'S ACCESSION- NO.
5 REPORT DATE
May 1982
6. PERFORMING ORGANIZATION CODE
8. PERFORMING ORGANIZATION REPORT C
10. PROGRAM ELEMENT NO.
11. CONTRACT/GRANT NO.
EPA IAG-78-DX-0166*
13. TYPE OF REPORT AND PERIOD COVERS
Final: 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-2683. (*) IAG with NASA, Kennedy Space Center, FL 32899,
subcontracted to U. of Miami under NASA Contract NAS 10-9410.
e.ABSTRACT  The six-volume report:  describes the theory of a 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 tw
 forms: free surface and rigid lid. The former,  verified at Anclote An-
 chorage  (FL),  allows a free air/water interface and is suited for sign
 fleant surface wave heights compared to  mean water depth; e.g., estu-
 aries and coastal  regions. The  latter, verified at Lake Keowee (SC), i
 suited for small surface wave heights compared  to depth (e.g., natural
 or man-made inland lakes) because surface elevation has been removed a
 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
a. DESCRIPTORS
Pollution
Thermal Diffusivity
Mathematical Models
Estuaries
Lakes
Plumes
13. DISTRIBUTION STATEMENT
Release to Public
b. IDENTIFIERS/OPEN ENDED TERMS
Pollution Control
Stationary Sources
19. SECURITY CLASS (This Report)
Unclassified
20. SECURITY CLASS (This page)
Unclassified
c. COSATI Field/Group
13B
20M
12A
08H,08J
21B
21. NO. OF PAGES
123
22. PRICE
EPA Form 2220-1 (9-73)
112

-------