United States
Environmental Protection
Agency
Environmental Research
Laboratory
Corvallis OR 97330
Research and Development
Mathematical
Cooling Tower
EPA-600/7-78-102
June T978
Energy/Environment

-------
                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.
 This document is available to the public through the National Technical Informa-
 tion Service, Springfield, Virginia  22161.
            For S.l. by thy Supenntcndcn, ol bocumrm,. L'.S. Govcrnm.nl Printing Ottce
                   Wi.hmiton. D.C :o«o: Slock No OS5-001 -01080-0

-------
                                             EPA-600/7-78-102
                                             June 1978
 MATHEMATICAL MODEL FOR MULTIPLE COOLING TOWER PLUMES


                          by
                    Frank H. Y. Wu
                   Robert C. Y. Koh

W. M. Keck Laboratory of Hydraulics and Water Resources
           California Institute of Technology
               Pasadena, California 91125
           Grant Number EPA (5) R-803989-01-1

                    Project Officer:
                      I
                   Mostafa A.  Shirazi
         Ecosystems Modeling and Analysis Branch
       Corvallis Environmental Research Laboratory
                 Corvallis, Oregon 97330
       CORVALLIS ENVIRONMENTAL RESEARCH LABORATORY
           OFFICE OF RESEARCH AND DEVELOPMENT
          U.S.  ENVIRONMENTAL PROTECTION AGENCY
                 CORVALLIS, OREGON 97330

-------
                           DISCLAIMER

This report has been reviewed by the Corvallis Environmental Research
Laboratory, U.S. Environmental Protection Agency and approved for
publication.  Approval does not signify that the contents necessarily
reflect the views and policies of the U.S. Environmental Protection
Agency, nor does mention of trade names or commercial products con-
stitute endorsement or recommendation for use.
                               ii

-------
                                FOREWORD



     Effective regulatory and enforcement actions by the Environmental


Protection Agency would be virtually impossible without sound scientific


data on pollutants and their impact on environmental stability and


human health.  Responsibility for building this data base has been


assigned to EPA's Office of Research and Development and its 15 major


field installations, one of which is the Corvallis Environmental Research


Laboratory (CERL).


     The primary mission of the Corvallis Laboratory is research on the


effects of environmental pollutants on terrestrial, freshwater, and
                               t
marine ecosystems; the behavior, effects and control of pollutants in


lake systems; and the development of predictive models on the movement


of pollutants in  the biosphere.


     This report presents a procedure for estimating the nature and


extent of dispersions of moist plumes that are emitted from cooling


towers.  The estimates are \iseful in assessing environmental impact of


these systems.
                                          A. F. Bartsch
                                          Director, CERL
                                111

-------
                             ABSTRACT
     A mathematical model is developed resulting in a computer program
for the prediction of the behavior of plumes from multiple cooling
towers with multiple cells.  A general integral method based on the
conservation of mass, momentum, energy (heat), and moisture fluxes
(before and after a plume merging), were employed in the prediction
scheme.  The effects of ambient stratifications of temperature, mois-
ture, and wind are incorporated in the model.

     An axisymmetric round plume is assumed to be emitted from each
individual cell before interference with neighboring plumes.  A finite
length slot plume in the central part and two half round plumes at both
ends of the merged plume were used to approximate the plume after merg-
ing. The entrainment and drag functions are calculated based on the modi-
fied merged plume shape.

     The computer output provides the predicted plume properties such
as excess plume temperature, humidity and liquid phase moisture (water
droplet), plume trajectory, width, and dilution at the merging locations
and the beginning and ending points of the visible part of the plumes.
Detailed printout and contour plots of excess temperature and moisture
distribution can also be obtained if desired.

     Based on comparison with laboratory data this model gives good
predictions for the case of dry plumes (no moisture involved). It
should be noted that several empirical coefficients are as yet not
accurately known.  Verification of this model for the wet plume (such
as for prototype cooling tower plumes) and the determination of the
values for these empirical coefficients to be used in prototype appli-
cations must await detailed comparison with field data.

     This report was submitted in fulfillment of Grant Number EPA (5)
R-803989-01-1 by California Institute of Technology under the sponsor-
ship of the U.S. Environmental Protection Agency.  This report covers
a period from 12/22/75 to 7/31/77, and work was completed as of 7/31/77.
                              iv

-------
                         TABLE OF CONTENTS







Chapter                                                              Page




   1.     INTRODUCTION                                                 1




   2.     THEORETIC MODEL                                              4




          2.1  Formulation                                             5




          2.2  Entrainment                                            13




          2.3  Merging Process                                        17




   3.     COMPUTER PROGRAM                                            28




   4.     RESULTS, COMPARISONS AND DISCUSSIONS                        31




   5.     CONCLUSIONS AND RECOMMENDATIONS                             57




   REFERENCES                                                         59
                             t



   APPENDIX A:  COMPUTER PROGRAM




   APPENDIX B:  EXPLANATIONS OF THE IMPORTANT SYMBOLS IN THE


                PROGRAM MTP




   APPENDIX C:  LISTING OF PROGRAM

-------
                             LIST OF FIGURES

Figure                                                              Page

  2.1       Top View of the Coordinate System Associated with         6
            the Cooling Tower Configuration

  2.2       Side View of the Coordinate System Associated with        6
            the Cooling Tower Configuration

  2.3       Merged Plume Shape                                       ^

  2.4       Definition Sketch                                        20
                                                                     21

  2.5       Modified Merged Plume Shape                              24

  2.6       Correction of Ay for the Merged Plume                    26

  3.1       Flow Chart of the Computer Program                       30

  A.I       Excess Temperature Distribution for the Case of          32
            4 Towers in Linear Array and 0° to Wind Direction
                               I
  4.2       Excess Humidity Distribution for the Case of 4           33
            Towers in Linear Array and 0° to Wind Direction

  4.3       Excess Liquid Phase Moisture Distribution for            34
            the Case of 4 Towers in Linear Array and 0° to
            Wind Direction

  4.4       Excess Temperature Distribution for the Case of          35
            4 Towers in Linear Array and 90° to Wind Direction

  4.5       Excess Humidity Distribution for the Case of 4           36
            Towers in Linear Array and 90° to Wind Direction

  4.6       Excess Liquid Phase Moisture Distribution for the        37
            Case of 4 Towers in Linear Array and 90° to
            Wind Direction

  4r. 7       Excess Temperature Distribution for the Case of          38
            4 Towers in Linear Array and 135° to Wind
            Direction

  4.8       Excess Humidity Distribution for the Case of 4 Towers    39
            in Linear Array and 135° to Wind Direction

  4.9       Excess Liquid Phase Distribution for the Case of 4       40
            Towers in Linear Array and 135° to Wind Direction


                                 vii

-------
                       LIST OF FIGURES (Continued)

Figure                                                              Pag

  4.10      Excess Temperature Distribution for the Case of          41
            5 Towers in Round Array

  4.11      Excess Humidity Distribution for the Case of 5           42
            Towers in Round Array

  4.12      Excess Liquid Phase Moisture Distribution for            43
            the Case of 5 Towers in Round Array

  4.13      Comparisons of Plume Trajectory, Width and Dilution      45
            Between the Present Theory and Fan's (1967)
            Experiments for F = 20 and K=8

  4.14      Comparisons of Plume Trajectory, Width and Dilution      46
            Between the Present Theory and Fan's (1967)
            Experiments for F = 40 and K = 8

  4.15      Comparisons of Plume Trajectory, Width and Dilution      47
            Between the Present Theory and Fan's (1967)
            Experiments for F = 80 and K = 16

  4.16      Comparisons of Plume Excess Temperature Between          48
            the Present Theory and Chan et al.'s (1974)
            Experiments for F = 4 and K=1.02

  4.17      Comparisons of Plume Trajectories Between the            50
            Present Theory and TVA (1968, TVA-11, single tower)
            Field Data in a Stable Atmospheric Condition
            O6/3z = 0.590K/100 m)

  4.18      Comparisons of Plume Trajectories Between the            51
            Present Theory and TVA (1968, TVA-14, two towers)
            Field Data in an Atmosphere with an Inversion
            (39/32 = 1.42°K/100 m)
                                viii

-------
                             LIST OF TABLES

Table                                                               Page

 A.I        Input Data Cards for Example Cases (Different            52
            Wind Directions to a Line Array (3 cases) and
            Round Array (one case) of Towers)

 4.2        Input Data Cards for Three Cases of Fan's (1967)         53
            Experiments with F = 20, K=8; F = 40, K=8; and
            F = 80, K= 16

 4.3        Input Data Cards for Chan et al.'s (1974)                54
            Experiment with F = 4, K = 1.02

 4.4        Input Data Cards for TVA (1968, TVA-11,  single           55
            tower) Field Data  in a Stable Atmospheric
            Condition (39/3z= 0.59°K/100 m)

 4.5        Input Data Cards for TVA (1968, TVA-14,  two              56
            towers) Field Data in an Atmosphere with an
            Inversion O8/3z=1.42°K/100 m)
                                 ix

-------
                             LIST OF SYMBOLS
 a              Entrainment coefficient associated with jet type entrainment

 a.             Entrainment coefficient associated with the effect of
               buoyancy

 a.             Entrainment coefficient associated with thermal type
               entrainment

 a,             Entrainment coefficient associated with ambient turbulence

 a              Finite length of the slot plume in the central part of
               the merged plume

 A              Finite length of the slot plume in the central part of
               the modified merged plume

b             Plume radius

b ,,b -       Plume radii at the ends of the merged plume

b             Plume width of the slot plume in the central part of the
              merged plume

B1,B2         Plume radii at the ends of the modified merged plume

BXZ           Plume half height of the cross-section normal to s-axis

BY            Plume half width

 C              Tracer concentration

 C.             Drag coefficient
 Q

 C              Specific heat at constant pressure
 pa

 e              Saturation vapor pressure
 O
 E              Volume rate of entrainment

 E^             Volume rate of entrainment associated with the effects
               of momentum and buoyancy for a round plume

 E.             Volume rate of entrainment associated with jet type
  8            entrainment for slot plume

 E2             Volume rate of entrainment associated with thermal type
               entrainment

-------
                       LIST OF SYMBOLS (Continued)




E_            Volume rate of entrainment associated with ambient

              turbulence



F             Sum of excess vapor and liquid phase moisture fluxes

              (H+W)

                                          r I f\

F2            Densimetric Froude number (m   /y28)



FrT            Local densimetric Froude number (U 2/[(p -p )bg/p ])
  L                                             p     a  p     o

Fr
  LC          Critical local densimetric Froude number



g             Acceleration of gravity



G             T-L W/C
                  v   pa


H             Excess vapor phase (humidity) moisture flux



HT            Total height of the merged plume



L             Latent heat of evaporation or condensation
 v                             , *


m             Kinematic momentum flux (f.u 2dA)
                                      ^ A p


M             Excess kinematic momentum flux



p.            Dry air absolute pressure



p             Total absolute pressure (p, + e )



P             Plume periphery



q             Specific vapor phase moisture (humidity)



Q             Plume volume flux



s             Plume trajectory



t             Temperature



t             Steam point temperature
 S


T             Excess temperature flux



U             Velocity



U '           Measure of ambient turbulent fluctuation
 a


U 1,U 2       Plume velocities at both ends of merged plumes


                                 xi

-------
                       LIST OF SYMBOLS (Continued)



U             Plume velocity at central part of merged plume
 s


W             Excess liquid phase moisture flux



W,            Plume width
 a


WD            Total width of the merged plume



x             x-coordinate (parallel to wind direction)



y             y-coordinate (normal to x-axis)



z             z-coordinate (vertical coordinate)



a             Entrainment coefficient for pure jet  (F -»• °°)



a             Entrainment coefficient for pure plume



3             Buoyancy flux Cf.gU (p -p )/p dA)
                              A  P  a  p   3.


F             Adiabatic lapse rate



6             Angle between the tangent of s and x-axis



P             Density



a             Liquid phase moisture



4"             Angle between the centerline of the (inclined) plume

              cross-section and the horizontal line parallel to y-axis



y             Volume flux ff.U dA)
                          ^ A p


(  )          in ambient
    3.


(  )          at tower exit



(  )          in plume
                                xii

-------
                         ACKNOWLEDGEMENTS







     The writers would like to express their gratitude to Professors




Norman H. Brooks and John F. Kennedy for providing valuable comments




and suggestions during the investigation.  They would also like to




thank Drs. Mostafa Shirazi and Larry Winiarski for several stimulating




discussions on the project.




     The work reported herein was supported by EPA Grant Number




(5) R-803989-01-1.
                                xiii

-------
                             CHAPTER 1



                            INTRODUCTION




     The use of multiple cooling  towers  is a common means of disposal of


waste heat  from large power plants.  A better understanding of the be-


havior and  interaction of plumes  from multiple towers would be useful


not only to cooling  tower design  and operation, but also in the assess-


ment of their environmental impact.  Due to the relatively close


proximity of neighboring tower exits, the individual plumes from multiple


cooling towers rapidly interfere  with one another, thus changing the


overall plume shape  and its mixing characteristics.  In addition, the


ambient stability  (temperature profile), humidity, wind velocity and


wind direction to  the tower array also influence the plume behavior.
                               »

In this report a mathematical model resulting 'in a computer program is


developed which will provide the  framework to allow reasonable predictions


to be made of the  characteristics of plumes from multiple cooling towers


under various ambient conditions.


     Throughout this report, no distinction will be made between a


plume and a jet.   Neither will the distinction be made between the multi-


ple plumes from the  several cells of a tower and the multiple plumes


from several individual towers.


     The plume from  a cooling tower is a buoyant jet which has been


studied in the past  by numerous investigators, such as Morton, et al.


(1956), Morton (1957), Slawson and Csanady (1967,1971), Koh and


Brooks (1975).  Several existing  single  tower plume models, including .


dry and wet plumes,  were developed by Briggs (1969), Hoult et al. (1969),


Abraham (1970), Fox  (1970), Csanady (1971), Wigley and Slawson (1971,1972),

-------
Hirst (1971), Hanna (1972), Weil (1974), Wigley (1975a,b), and




Schatzmann (1977).



     A model for multiple plumes was first developed by Koh and Fan




(1970) in analysing the analogous problem of disposal of waste heat into




the ocean by multiple port diffusers.  The integral method was used and




the individual round buoyant jets were approximated by a two-dimensional




slot jet after interference.  A transition region during merging was




considered.  Although a discontinuous centerline temperature resulted at




the point of merging, the overall predictions of the plume properties




were quite satisfactory.  Jirka and Harleman (1974) approached the




multiple port diffuser problem by replacing it with an equivalent slot




jet having the same mass and momentum fluxes as the multi-port discharge.




As expected this method generally over-estimates dilution except for




those instances when the plumes  are  very close to each other initially.




Briggs (1974) added an enhancement factor to his single tower equation




for plume rise by considering the effects of number of towers and tower




spacing.  Meyer et al.  (1974) modified Briggs' equation and used a




"peak factor" to develop a model which can give a fairly good prediction




of visible plume length.  However, the prediction  of the plume




trajectory is less accurate.  Davis  (1975) developed a mathematical




model for calculating plume rise and dilution from multiple cell




mechanical draft cooling towers with the wind normal to the tower array.




The entrainment function used includes the effects of plume interference




and the  changing entrainment surfaces during merging.  The model also




provides calculation techniques  for various modes of plume development.




The plume properties remain continuous when the calculations proceed

-------
smoothly from one zone to another.  The values of the coefficients in




his entrainment function still need to be determined from suitable




laboratory and field experiments.




     Data from multiple towers are very scarce.  Chan et al. (1974)




made laboratory simulations of plumes from multiple towers using water




as the fluid medium.  Two sets of normalized excess temperature distri-




butions were presented.  These are useful for model comparisons.  The




studies by Carpenter et al. (1968) and Slawson et al. (1975) at TVA




gave field data on plume trajectory and some plume properties.  More




complete data including plume trajectory, width, dilution, and ambient




conditions (i.e. temperature, humidity and wind velocity profiles, wind




direction to tower array, tower configuration, etc.) however are required




for proper verification.




     The model developed in this report is based on a general integral




method applied to the conservation equations for mass, momentum, energy,




and moisture fluxes.  An axisymmetric round plume is assumed initially




for each tower exit.  As the plumes merge, combinations of round and




slot plumes are employed to simulate the shape of the resulting merged




plumes.  The merging criteria, merging processes, changes of plume




shape and entrainment functions are a part of the model and are discussed




in Chapter 2.  Some results of model predictions and comparisons with




laboratory and field data are presented and discussed in Chapter A.  A




computer program has been written to perform the calculations and is




included in Chapter 3 and Appendices A and C.

-------
                             CHAPTER 2





                         THEORETICAL MODEL






     The present mathemtical model is developed for the prediction of




plume properties from multiple cooling towers.  These include plume




temperature, moisture (vapor and liquid phases), excess temperature,




excess moisture, velocity, width, dilution, trajectory and visible




plume length.  For ease of application to practical situations, this




model is capable of handling rather arbitrary vertical profiles of




ambient temperature, humidity, and velocity; arbitrary but steady wind




direction to the tower array; and randomly arranged tower configurations.




     The assumptions made in developing this model are as follows:




     1.  The flow is fully turbulent.  Molecular transport can be




neglected in comparison with turbulent transport so that there is no




Reynolds number dependence.




     2.  Longitudinal turbulent transport is small compared with




longitudinal advective transport.




     3.  Pressure is hydrostatic throughout the flow field.




     4.  The cross-plume profiles are similar for plume velocity,




temperature, density, humidity and liquid phase moisture.




     5.  The Boussinesq assumption is valid.  This implies that the




variations of fluid density throughout the flow field are small compared




with the reference density chosen.  The variations in density are only




considered in the buoyancy term.




     Using these assumptions, a general integral model for multiple




cooling tower plumes based on the conservation of mass, momentum, energy




and moisture fluxes along the plume trajectory is developed.  By providing




                                 4

-------
the ambient conditions and  the empirical equations for entrainment and




drag, the conservation equations are  integrated stepwise for the center




line properties along the plume trajectory.  Before interference the




plumes are assumed to be individual,  axisymmetric, round buoyant jets.




During the merging process, a combination of a finite slot jet in the




central part and  two half round jets  at both ends is assumed to be the




cross-sectional shape of the merged plume as an approximation but only




for the calculation of entrainment and drag.  Finally, the completely




merged plume gradually tends to become round in cross section again,




whereupon the individual axisymmetric analysis is reapplied.  The for-




mulations of the  basic plume conservation equations, entrainment function,




merging criterion and merging process are presented in the following




sections.
                             *





2.1  Formulation




     The coordinate system  chosen with a typical cooling tower configura-




tion is shown in  Figures 2.1 and 2.2.  The x-axis is parallel to the




steady ambient wind direction,  s  is the coordinate along the plume




path and 8 is the angle between the tangent to  s  and the x-axis.  The




individual plumes from the  cooling tower cells are presumed to be dis-




charged vertically into a stratified atmosphere, and bent over due to




the effect of ambient wind.  The plume properties are defined as velocity




U, density p, temperature t, specific humidity q, and liquid phase




moisture a.   Here, the specific humidity (vapor phase moisture) q and




liquid phase moisture a are defined as the ratio of mass of vapor (or




liquid) phase moisture to the total mass of the mixture in a unit




volume.  The subscripts o, p, and a are used for the values at tower




                                5

-------
                                Y
Figure 2.1  Top View of the Coordinate System Associated with the
            Cooling Tower Configuration
                  U0COS
e
                                Up= U0COSC7+U
 Figure 2.2  Side View of  the Coordinate  System Associated with the Cooling
             Tower Configuration
                                 6

-------
exit, in the plume, and in the ambient atmosphere, respectively.  The



plume volume flux Q, kinematic momentum flux M, temperature deficiency



flux T, vapor phase moisture (or humidity) deficiency flux H, and liquid



phase moisture deficiency flux W are defined as follows:



                 Q = JA Up dA                               (2.1)
                 M = /  U   dA                              (2.2)
                      A  p
                 T = /A (tp - ta) Up dA                    (2'3)
                 H = JA (q  ~ qa) Up dA                     (2.4)
                 W = JA (Op -' oa).Up dA                     (2.5)





Any 'rainout1 of liquid droplets in the ambient atmosphere will be



neglected so that o  is equal to zero.
                   Si


     The conservation of mass equation is:
where E is the volume rate of entrainment of ambient fluid.  The function



used for E is an empirical expression including the effects of plume



geometry, local mean velocity, buoyancy and ambient turbulence.  The



detailed form of E is presented in Section 2.2.



     The conservation equation for horizontal momentum flux is:






     IF < JA Pp  Up2 ^  - cos 6} = pa UaE + % p& U/ sin2  6 Cd Wd | sin e|



                                                            (2.7)

-------
where the first term on the right hand side of equation  (2.7)  is due  to

the momentum of the entrained ambient fluid and the second term is due

to the horizontal drag force (Fan, 1967) on the plume; Cd is a drag

coefficient to be determined empirically and W, is the width of the

plume; the absolute value of sin6 is used in order to account  for both

the ascending and descending parts of the plume.

     The conservation equation for vertical momentum flux is:


                 dl" {'A "p UP2 dA ' sine} =


     «/A (pa-pp> dA - «/A Vp1* + I Pa "a  ' S^ 9 CdWd COs9

                                                             (2.8)


where the first term on the right hand side of equation  (2.8)  is the

buoyancy force due to the density difference between the plume and

ambient fluid; the second term is the effect of the weight of  liquid

droplets suspended in the plume; the third term is the vertical drag

force on the plume with the negative and positive signs  corresponding

to 0$6$TT/2 and -ir/2*6<0, respectively.

     Koh and Fan  (1970) have shown that the governing equation for the

flux of a tracer quantity C is as follows:

                 d                         dc
                 dT/A  
-------
and the liberation of  latent heat due  to  condensation of water vapor



inside the plume  (Csanady,  1971), the  conservation equation for energy



(or heat) flux  is derived as follows:
A
3- /A (t  - t ) U  dA = -
as 'A   p   a   p
dz
              + r
                                        U  dA + 3- /
                                      A  p         J
                                    3    .
                                    ds JA
                                                            pa



                                                            (2.10)
                                                               a  U  dA
                                                                P  P
where L  is the latent heat of evaporation  (or condensation) i.e. ,




Lv(t) =  [597.31 - 0.57x(t°C)] x 4.1868 Jg'1 for t*0°C and



L  = [677 + 0.622 x  (t°C)] x 4.1868 for t<0°C, and C   is the specific
 if                                                  p 3.


heat of  air at constant pressure.




     The conservation equation for  (vapor and liquid phase) moisture




flux is:                      .
d

T~
ds
                                                dq
                   ) U  dA + 3-  f . 0  U  dA - -- ^ sin6  f. U  dA
                                 J                         J
                .
                App
                                                           .
                                                          Ap
                                                             (2.11)
     We assume a top-hat distribution of plume properties across the



plume.  Therefore, U  , T , p , q  , and o  are constants inside the
plume.  In particular, U  is defined as follows:
U  = U  cos9 + U
 P    a
                                                             (2.12)
where U is the 'net' plume velocity relative to the ambient wind.



     Applying the relationship between temperature and density for a constant



pressure process (pressure variations in the plume cross section  is



neglected) we can write:

-------
                p     t         p   p     t  -t
                a    p      ,  a— p    p
                —  =  —c-  ;  and  	"=-  =  —^7
                PP    
-------
                 ap = 0          for qp <. qsp  (dry plume)




                                                             (2.21)





                 qp = qsp  (tp)   f°r S * qsp  (wet plume)





where q   is the plume saturation specific humidity.
       sp



     Equation (2.21) implies that when the plume is unsaturated, no




liquid phase moisture due  to condensation need be considered and when




the plume is saturated, the plume humidity is equal to the plume satura-




tion specific humidity.  To calculate the saturated specific humidity,




thermodynamic equilibrium  between liquid and vapor is assumed.  The




Clausius-Clapeyron equation (or tabulated values) can then be used to




calculate the saturated humidity.  The equation used in  the model to




calculate the saturated specific humidity q  is (Linsley et  al.,
                                           S
                             t


1975):





                  0.622e (t)             0.622 e(t)
       /   \            ^       —.

     Vt>p' = P.. - 0.378  e (t) " P,(z) + e  (t) - 0.378  e  (t)
                t          S       Q       S             S



                                                             (2.22)






where t and p are absolute temperature and pressure, respectively; P




is the total pressure which is the sum of the dry air pressure  Pd(z)




and the saturated vapor pressure e (t).
                                  S



     Since the variation of pressure is small  (for instance, AP = 1% of P




corresponding to Az = 100M), the absolute pressure at sea level (i.e.,




Pd = 1013.25 mb) is used in equation (2.22), which becomes:




                            0.622 eg(t)


                 qs(t) = 1013.25 + 0.622 e (t)               (2.23)
                                          S


     An approximate expression of e  (t) was  developed by Richards  (1971)
                                   S


as



                                11

-------
     e  (t) = 1013.25 x exp(13.3185 t -1.9760 t 2-0.6445  t  3-0.1299  t  A)
      S                           V         V           v          v


where




     
pa
            ' V'o - 'ao - c-
                             pa



    F  =H  +W  =fD2U  (q
     o    o    o   4  o   o  Mo
            = Qo
    z  = 0
     o
                                                             (2.24)
                                             pa

                            L,
                                 12

-------
where the subscripts o  and  ao  are  associated with  the values at  tower


exit  (or initial values)  and the ambient  at the  same level, and


Ly(t) =  [597.31 - 0.57  x  to(°C)] x 4.1868 Jg'1,  and C a =  1.005  Jg~loK~




2.2  Entrainment


     The entrainment of ambient fluid  into the plume is a  function  of


plume geometry, local mean  velocity, buoyancy, and ambient  turbulence.


The entrainment function  first proposed by Morton  et al.  (1956)  is:
                       E,   =  a  2rrb  U                           (2.25)
                        Ir          p
where a is entrainment coefficient  determined  from  experiments;  b  is  the


round jet radius; and U   is  the'jet centerline velocity.


     Based on the integral conservation equations of mass, momentum,


energy, and mechanical energy, and  assuming  similar profiles,  Fox  (1970)


and Hirst (1971) derived  an  entrainment function for round jets  which


includes the effect of buoyancy to  the entrainment.  It reads  as follows:


                             a2
                 Elr = (al + Fr~ Sin9) 2lTb Up                (2'26)
                               Li


where a, and a2 are entrainment coefficients,  and Fr.  is  the local

                                             2
densimetric Froude number defined as Fr  = U  /[(p_ -  p ) b g/P0] =
                                       LI     p     a   p       o
  2
U  A[(t  - t ) t bg/t  t  ],  and g is gravitational  acceleration.  Based
 p     p    a   o    pa

on experimental results, a,  was determined to  be 0.057 for pure  round


jet with Gaussian profile distribution (i.e.,  Fr, -»• ~).   Hirst (1971)


suggested the value of a~ =  0.97.   This appears to  be  too large  when  his


results are compared with experiments.  A better estimate of the value of


                                  13

-------
a» can be made from the work of List and Imberger  (1973), who,  based on
dimensional analysis and experimental data, derived a  similar expression
for E..  for round jets  (Koh and Brooks, 1975)
                 Elr =  (0.057 + ^P) 2*b Up                (2.27)

       9    c / o
where F  = m   /y26 with m being  the kinematic momentum flux
       r\
- /. U   dA, y the volume flux =.J. U  dA, and 3  the buoyancy flux
   A 6 p   p
       *    a
     Comparing equations  (2.26) and  (2.27)  in a quiescent  ambient
                                 2
(i.e. sine = 1) and calculating F  by using the Gaussian  similarity
                                                   2
profiles for U  and p  , one finds a2 = 0.083 Fr /F  =  0.4775.   Hence
equation (2.26) becomes

                 E,  =  (0.057 + °;4775   sine) 2irb  U          (2.28)
                  Ir             FrL               p

     The entrainment coefficient a in equation  (2.25)  has  the  extreme
values for pure jet  (i.e., FrL -* ») a  = 0.057 and pure plume  (i.e.,
FrT •»• 0) a  = 0.082 in  a  quiescent ambient.  A critical value  of FrT
  L       p                                                         L
may be determined from

                  0.082 = 0.057 + ^
 which gives
                  FrLC - 19'X
      For Gaussian similarity profiles of plume properties it will be
 assumed that
                                                              (2.29)
      Elr = (0.057 -1-        sine) 2*b U  for FrL > 19.1
            0.082 • 27rb U                for Fr. < 19.1
                         H                     Lt
                                 14

-------
which implies that Fr  = 19.1 was considered a small number below which
                     LI


the entrainment of a buoyant jet is similar to that of a pure plume.



     Experimental results for two-dimensional slot jet are not sufficiently



comprehensive to obtain a similar entrainment expression [i.e., equation



(2.29)] (Koh and Brooks, 1975).  Therefore, based on the experimentally



determined entrainment coefficient, the following form will be used for



the slot jet with Gaussian profile distribution, viz.,





                 Els = 0.14  • 2A • U                        (2.30)





where A is the length of the slot jet.



     The entrainment functions embodied in equations (2.29) and  (2.30)



are based on Gaussian profiles of plume properties.  In this study,



top-hat similarity profiles are' assumed.  The difference in the resulting



entrainment functions [equations (2.29) and (2.30)] due to this is a



factor of /2~.  Therefore, equations (2.29) and (2.30) can be rewritten



as follows:





     E,  = (0.0806 + °16753 sine) 2irb U   for FrT > 19.1
      ir              rr.              p        L*

                        L                                   (2.31)


         = 0.1160   2ub U                 for FrT $ 19.1

                         P                      L






     Els = 0.198 • 2AUp                                     (2.32)





     As the plume bends over towards the direction of the ambient wind,



the plume velocity is about equal to the wind velocity.  Then the entrain-



ment should be nearly as if the plume were a two-dimensional thermal  in



a stagnant atmosphere.  This entrainment is proportional to the jet



periphery and the velocity of the thermal.  Abraham  (1970) proposed the


                                 15

-------
following form






                 E, = a» P U  sine cos6                      (2.33)
                  £    J    Si





where P is the jet periphery, cos6 is arbitrarily chosen to diminish the




thermal type of entrainment closed to the initial stage of the vertical




jet, and a., is the entrainment coefficient for a line thermal.  For




large Reynolds number, the experimentally determined value of a- is




0.5 (Richards, 1963).  But a better value suggested by Koh and Chang (1973)




from their plume measurements and numerical model is 0.3536, which will




be used in this present model.  Thus




                 E0 = 0.3536 P U  sine cos6                  (2.34)
                  2             a




     Another type of entrainment is associated with ambient turbulence,




and expressed as





                 E3 = a4 P IT                                (2.35)






where a, and U' are the entrainment coefficient and a measure of turbulent
       f      3



velocity fluctuations.  Based on dimensional analysis, Briggs (1969)




found that U' is associated with eddy energy dissipation in the inertial
            cl



subrange and gave an estimate of a^ = 1.  In practice, the root-mean-




square value of the ambient wind velocity fluctuation may be used to




approximate IK, which is equal to a few percent of the mean wind velocity




under normal atmospheric conditions.




     Finally, we may combine equations  (2.31), (2.32), (2.34) and (2.35)




to construct a complete entrainment function as






                 E = P {a|u|+0.3536 Ujsin9 |cos6+1.0 ITJ     (2.36)






                                 16

-------
where  for a round jet,




                  P  =  2ub
                  a  =  0.0806  + °'p'"-J|sine|  for  FrT  >  19.1
                                 f i _               \_i
                      0.1160                 for  Fr.  $  19.1
                                                  Li
and for a  slot jet,
                 P =  2A


                 a =  0.198


Here U is the net velocity  in  the plume relative  to  the ambient velocity;


the absolute value is used  here  to account for both  the ascending and


descending parts of the plume.


     In the literature on buoyant jets, various investigators employing


the integral approach have  devised differing entrainment functions.  A


recent survey for the round jet  can be found in Wright (1977).  The


entrainment function  expressed in equation (2.36) and incorporated in


the present model is  but one possible expression.  Should a different


form be shown to be superior in  the future, the model can readily be


modified.




2.3  Merging Process


     The individual plumes  from  the multiple cells of a cooling tower


typically merge within a relatively short distance from the exits.


Before the plumes merge, equations for individual round buoyant jets


are applied in this model to calculate the plume behavior.  When several


individual plumes are merged, the resulting plume cross-section is no


longer round, but rather tends to be elliptical In shape.  In this model,



                                 17

-------
this merged plume is approximated by a slot jet in the central part and




two half round jets at the two ends of the merged plume as shown by the




solid lines in Figure 2.3.  The nonuniform size of the plume is due to




the effect of the wind direction with respect to the tower configuration.




In general, it is necessary to consider all types of plume merging




including all the possible combinations between individual round plumes




and modified merged plumes as shown in Figure 2.4.  The basic merging




criterion considered here is that the plume cross-sections are in con-




tact with each other.  An additional criterion is incorporated for the




merging between two individual round plumes: the area of the trapezoid




should be equal to the sum of the areas of the two half round plumes




as circled by the dashed lines in Figure 2.3.  When the plumes satisfy




these merging criteria, they are merged.  The fluxes of the merged plumes




are summed to maintain the conservation of fluxes.  Moreover, the new




shape and the new centroids of the merged plumes are determined,




and the integration of the equations is continued.  Upon merging, the




entrainment and drag functions are altered due to the change in plume




shape.  The merged plume shape is characterized by the radii Bl and B2 of




the two half round plumes, length of the slot jet A, and the angle <|»




(shown in Figure 2.3) between the centerline of the (inclined) plume




cross-section and the horizontal line parallel to the y-axis.  As the




plumes merge a new set of Bl, B2, A and  the plumes are




classified  into two categories:  one will be called horizontal for which




the total width (WD) of the new merged plume is larger than the total




                                18

-------
Figure 2.3  Merged Plume Shape
                                19

-------
   	L
                 (a)
Figure 2.4  Definition Sketch
                                 20

-------
             (g)
(h)
Figure 2.4  Definition Sketch
                               21

-------
height (HT), otherwise, it will be called vertical.  Both  types  are

illustrated by Figures 2.4(a), (c) and  (e); and  (b), (d) and  (f),

respectively.  These two categories are not different in substance

and the distinction is made primarily for coding convenience  in  the

computer program.  Bl and B2 are chosen to be the radii of the left  end

(or lowest end) and the right end (or highest end) of the  horizontal (or

vertical) plumes, respectively.  The left end of the plume is the end

closer to the x-axis.  For the cases of merging among individual round

plumes as shown by Figure 2.A(a) and (b), $ is the angle between the

horizontal line and the line connecting the two centers of the merged

ending plumes.  For the cases of merging between merged plumes,   is the

average of angles $. and $2  of the merging plumes 1 and 2, respectively,

as shown in Figures 2.4(c) and (d).  For the cases of individual round

plume 2 joining the merged plume 1 as shown in Figures 2.4(e),  (f),  (g)

and (h), $ is assumed to maintain the original value of the merged plume

2 since the resulting merged plume is envisioned to be dominated by  the

merged plume 2.  After Bl, B2 and  are determined, A can  be  calculated

by the following equations.


For horizontal plumes:

                 A = (WD-Bl-B2)/cos<|>     for cos * 0

                 A = HT-B1-B2           for cos    for sin^, V 0

                 A = WD-B1-B2           for sin<|>  = 0          (2.37)


When  the new shape of  the merged plume  is determined, then the  calcula-

tion  can be  performed  forward one integration  step  for  the round jets at
                                  22

-------
the ends and the slot jet  (with unit fluxes found by dividing the total




fluxes of the slot jet by  the finite length A) in the central part.




This will result in new values for the radii of the round jets at the




ends of the merged plumes, b , and b 2 and the half width and length of




the central part slot jet  b  and a.  Because of the different entrain-
                           s



ment rates for round and slot jets, the calculated plume cross-section




determined by b , , b ~ » b  and a n^Y not be smooth enough to represent a




realistic shape.  The discontinuities occurring at the junctions of the




round and slot jets are demonstrated by the dashed line curve in Figure




2.5.  In order to eliminate the discontinuity and to obtain a modified




smooth plume cross-section described by Bl, B2, and A, the following set




of equations  is  proposed:
                 0.5* (brl2 Url  + br22 Ur2) + 2bsaUs=






                [O.Sir (Bl2 + B22) + A (Bl + B2)] • U        (2.38)
                 a + brl + br2 = A + Bl + B2                (2.39)






                 B1/B2 = brl/br2                            (2.40)






where U , , U 2, U  and u are the plume velocities corresponding to the




half round jets with radii brl and br2, the slot jet with half width bg,




and the overall merged plume defined by Bl, B2 and A, respectively.




     Equation (2.38) describes the redistribution of the volume flux




from the calculated merged plume to the proposed modified plume. Equa-




tion (2.39) maintains the same plume length between the calculated and




modified plumes.  Equation (2.40) keeps the same ratios of the radii of




the two half round plumes between calculated and modified plumes.



                                23

-------

Figure 2.5  Modified Merged Plume Shape
                                 24

-------
     After the modified plume  shape  is  determined,  the  half  width,  BY

and the half height,  BXZ of  the merged  plume,  indicated in Figure  2.3,

can be determined by  the following equations  for  the  purpose of  checking

plume merging at the  next step:


      BY = 0.5x  (A  •  cos<{> +  Bl +  B2)     for cos<{>  H  0

         = Bl                            for cos*  =  0  and  Bl  £ B2

         = B2                            for cos4>  =  0  and  B2  > Bl

     BXZ = 0.5x  (A  •  sin<|> +  Bl +  B2)     for sin  *  0

         = Bl                            for sin<(>  =  0  and  Bl  i B2

         = B2                            for sin<)>  =  0  and  B2  > Bl

                                                             (2.41)


     Due to the uneven change of  Bl, B2  and A  for each  integration step,

the y-coordinate of the plume centroid also needs to  be readjusted.

The amount of adjustment Ay  noted in Figure 2.6 is

        T
      = [
          A+B1+B2     A             A   +B1   +B2
   Ay=    -1    •]   J  x  -    + B2j+1 -  1+1   J*1   -1  i
      = 0.5x|cos<)>|x   (Bl. - B2.)x -^— + B2   - I
                     L   J     1     A.J       J+1

where j and j+1 refer to the  calculation steps.

     With the modified merged plume cross-sectional  shape,  the  entrain-

ment and drag force  can be  determined  and the conservation  equations

integrated.  During  the calculation, barring  further merging, the length

of the slot jet A generally will be reduced and the  radii of  the two

ending round plumes will be increased.   Finally, when A diminishes to

zero, the shape of the merged plume cross-section  becomes practically
                                25

-------
                                                   Aj+l +B2j + l
                                                                     j + l
                                                                        Aj+l   H
                                                   0.5X(Aj+Blj-B2j)X—	  —
                                                   0.5X(Aj+Blj+B2j)-B2j

                                                  = 0.5X(Aj+Blj-B2j)
Figure 2.6  Correction of Ay for the Merged Plume
                                 26

-------
round.  At that point, a round plume is again adopted to carry through




the final stage of plume calculation.
                                 27

-------
                             CHAPTER 3






                         COMPUTER PROGRAM






     The governing equations for predicting the dynamic behavior of




multiple cooling tower plumes were presented in Chapter 2.  No analytical




solution can be obtained due to their complexity as well as the arbitrary




ambient conditions in the governing equations.  Therefore, a computer




program written in Fortran IV language was developed to solve these




equations.  A standard fourth order Runge-Kutta method was employed in




the solution.  The inputs to the program include tower exit conditions,




ambient conditions, tower configuration, entrainment and drag coeffi-




cients, and some control parameters.  The basic routine of the computer




program begins with inputting data, setting initial conditions (sub-




routine SETIC) and calculating the first plume (from the tower with the




smallest value of x) by setting an indicator IND(I) = 1 for i




individual round plume (subroutines RUNGS and DERIVR).  As the calcula-




tion continues, the subroutines CHKNWP, ALIGN and PLMERG are called to




check  for the appearance of  any new plumes, to align the existing plumes




at approximately the same x-coordinate, and to check for the merging




among  the existing plumes, respectively, along the direction of the




plume  trajectory.  If new plumes appear (whenever x exceeds the x-




location of downstream tower exits), the results for such new plumes are




calculated stepwise until the stage is reached to necessitate the




checking of  the merging criterion.  If the plumes merge, the indicator




of the i   and the j   plumes are changed to IND(I) = 2 and IND(J) = 0




 (J>I).  In the subroutine RESETI, the fluxes of the merged plumes are




added  together, and the initial conditions for the merged plumes are



                                28

-------
reset.  The subroutines DERIVR and DERIVS are called to calculate the

plume half widths and velocities of the round and slot jets in order

to determine the shape of the modified merged plume.  Then, subroutine

DERIVE is used to calculate the dynamic properties of the modified merged

plume.  The calculation stops when the integration step number is equal

to the desired (input) step number.  The outputs include the input

information, and the calculated plume properties such as temperature,

excess temperature, moisture, excess moisture, half width and trajectory

at visible, merging and final stages of the plumes.  The detailed list-

ings and examples of the input and output are presented in Appendix A.

The general structure of the computer program is described by the flow

chart shown in Figure 3.1.  Some important variables in the text and
                              «
program are compiled and listed in Appendix B.  The complete computer


program is presented in Appendix C.
                               29

-------
                             Input
            Call BLANK to initialize some storages
            Calculate gradients of ambient profiles
                              I
             Call SETIC to set initial conditions    |
                              I
        Call  RUNGS and DERIVR to calculate the. results
                  of individual round plumes
        Call CHKNWP to check if  any new plume appears
   Call  ALIGN to align the existing plumes at  approximately
                     the same  x-coordinate
           Call  PLMERG to check if  any plumes merged
                               I
Call RESETI (including call RUNGS, DERIVR, DERIVS and DERIVE)
to reset the initial conditions and calculate
the shape and results of the modified merged plume
1


{ yes
Printout minimum amount of results
1

•* Further output needed?
I yes
Call OUTPUT for detailed printout and contou
1

r plots

                          »   END
Figure 3.1  Flow Chart of the Computer Program
                      30

-------
                                 CHAPTER 4






                   RESULTS, COMPARISONS AND DISCUSSIONS






     In this chapter, four example cases are presented.  The results of




the present theory are also compared with the laboratory results of Fan




(1967), Chan et al.  (1974), and field data from TVA (Carpenter et al., 1968).




     In the example  runs, a line array of four cooling towers and a round




array of five cooling towers are considered.  For the cases with the line




array, three wind directions (i.e. 0°, 90°, and 135° with respect to the




tower array) are chosen.  The input data cards are shown in Table 4.1,




which include the number of towers, the desired number of calculation




steps, control parameters, tower, configuration, ambient levels, temperature




profile, humidity profile, wind velocity profile, tower exit conditions,




coefficients of contour plots and heading of plots.  Normally, the outputs




consist only of the  input information, the results at the merging points,




and those at the beginning and ending points of the visible phases of plumes.




However, detailed printouts and/or contour plots can also be provided by the




program upon request.  The contour plots of excess temperature, humidity




and liquid phase moisture for these examples are shown in Figures 4.1




through 4.12.  The plots represent the distribution of the highest values




projected onto the X-Z plane.  Detailed explanations of the input and




output parameters are presented in Appendix A.




     Three sets of data from Fan (1967) for a single jet, one set from Chan




et al. (1974) for six towers and two sets from Carpenter et al. (1968) for




a single tower and multiple towers are chosen for comparison with the model.
                                 31

-------
                 EXCS TEMP  H TWRS 1 flRRflYS 0  OEG TO  WIND
          I2h
         10
          8
          0
             1 I  I  I  i  I  I  I  I  I  I  I  I  I  I I I  I  I  I  I
                                 TTTTTI
i  i  i  i
i  i  i  i  i  i  i i i
                                                        i  i  i  i  i  r
               o
             6
             x/o
                                 8
10
12
Figure 4.1  Excess Temperature Distribution for the Case of 4 Towers  in
           Linear Array and 0° to  Wind Direction
                            32

-------
                  EXCS HUMI  4 TWRS 1 flRRflYS 0  DEO  TO WIND
              I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I
           12
           10
           8
        o
        \
        rvi
         _r. I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I I  I  I  I  I  I  I  I  I

                0       2       4       6       8       10       12


                                       X/0
Figure 4.2   Excess Humidity Distribution for the  Case of 4 Towers  in

            Linear Array and 0°  to Wind Direction

-------
                  EXCS LIQM U TWRS 1 flRRflYS  0 OEG TO WIND
               I  I  I  I  I  I I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I I I  l_
          -2
               I  I  I  I  I  I  I  i  I  i i i i I  i  i  i  i  i  i  i  i  i  i  i  i  i  i
                                        6
                                        X/D
8
10      12
Figure 4.3  Excess Liquid Phase Moisture Distribution for the Case of
           4 Towers in Linear Array and 0° to Wind Direction
                              34

-------
                   EXCS  TEMP 4 TWRS  1 flRRflYS 90  OG TO  WIND
               i  i  I  I  I  I  I  l  l  l  l  I  l  I  I  I  I I I I  l I  l  l  I  I  I  i
                                                                12
Figure 4.4  Excess Temperature Distribution for  the Case of 4 Towers  in
           Linear Array and  90° to Wind Direction
                             35

-------
             I2h
             10 h
             8
          a
          \
          M4
                    EXCS HUM I  »A THRS 1 flRRRYS  90 DG TO HIND
                    I  I  I I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  l/T
             0
             2<  I  I  I  I  I I i l  I  l  l  I  I  I  I  I  I  l  I  I  l  I  I  I  I  I  l  I
                  0       2       4       6       8       10      12
                                          X/0
Figure  4.5  Excess Humidity Distribution for  the Case of 4 Towers  in
            Linear Array and 90° to Wind Direction
                                36

-------
                    EXCS LIQM 14  THRS  1  RRRRYS 90 DG TO  WIND

                  I  I  I  I I  I  I I  I I  I  I I I  I  I  I  I  I  I  I  I  I  I  I  I  I
            12
             10
             8
          a
          \
          rvi
             0
           -2
                'l I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I
                  0
6

X/D
8
10
Figure 4.6  Excess Liquid Phssa Moisture Distribution for  the Case of A

            Towers in Linear Array  and 90° to Wind Direction
                               37

-------
                  EXCS  TEMP 4  TWRS  1  RRRflY 135  DEG TO WIND
               I  I  I I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  M  I  I  I  I I
                   I  i  I  i  i  i  i i  i  i  i
                                                               12
Figure  4.7  Excess Temperature Distribution  for the Case of  4 Towers in
           Linear Array and 135° to Wind  Direction
                             38

-------
                  EXCS HUMI 4  TWRS  1  flRRRY  135 OEG TO WIND
            12
            10
            8
         o
         X.
         rvj
            0
           -2
               I  I  I  I  I  I  I  I  I  I  I  I
    II I  I  I
                 I  I  I  I  I  I  I  I  |  |  i  I  I  I  I  I  I  I  I I  I  I  I  I  I  I  I
                 02468

                                        X/0
10       12
Figure 4.8   Excess Humidity  Distribution for  the Case of 4  Towers in

            Linear Array and 135° to Wind Direction
                              39

-------
                 EXCS LIQM 4  THRS 1  flRRflY  135 DEG TO  WIND
              I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I I I I I I  I I  I  I  I  II
        o
        fvl
          -2
               I  I  I  I I I I  I  I  1  I  I  I  I  I  I  I  I  I_L I  1  I  I  I  I  I  I
                                        6       8
                                       X/0
10      12
Figure 4.9   Excess Liquid Phase Moisture Distribution for the Case of 4
            Towers in Linear  Array and 135°  to Wind Direction
                              40

-------
                        EXCS TEMP 5 TWRS IN  ROUND  flRRflY

                      I  I  I  I  I  I  I  I I I  I  I  I  I I
         o
         \
         ivl
                      	I  I I  I  I  I  I  I  I  I  I
                  i  i ......
0
                                                                12
Figure A.10  Excess Temperature Distribution for  the Case of 5 Towers

             in Round Array
                              41

-------
                        EXCS HUMI 5  TWRS IN  ROUND  RRRflY
                ' ' '  '	I  I  I  I i i i l  l  l  i  »  t  i  i  i  i
                 0246       8       10       12
                                        X/0
Figure 4.11  Excess Humidity Distribution for the Case of 5 Towers in Round
            Array
                               42

-------
                         EXCS LIQM 5  TWRS  IN ROUND flRRflY
            12
            10
             8
                I  I  I  I  I  I  I  I  I  I  I  I  I  I I I I I  T I  I  I  I  I  I  I
          o
          x
          rvi
            -2
                  0
I  I  I  I  I  |  |  |  |  I I I I  I  I  I  I  I  I  I  I  I  I  I  I
    24       6       8       10       12
                   X/0
Figure 4.12  Excess Liquid Phase Moisture Distribution  for the Case of
             5  Towers in Round Array
                                43

-------
The input data cards for these cases are shown in Tables 4.2 through 4.5.



The ambient conditions for the data from Fan and Chan et al. are



uniform.  The velocity and temperature at the tower exit and ambient



are chosen to satisfy the given values of densimetric Froude number F



and velocity ratio K.  The predicted results of dilution, plume trajectory



and width for Fan's cases are shown in Figures 4.13, 4.14 and 4.15, and



are compared with his experimental results.  The comparisons are



generally good.  The predicted excess temperature distribution



for the six tower case from Chan et al. is shown in Figure 4.16, together



with the experimental results.  In this case the six towers are in one



line array and the ambient flow to it is normal.  The contour plot is for



the distribution of values in the central x-z plane.  It seems that the



present model tends to overpredict the excess temperature.  The main



reason might be because of the neglect of the effect of the mixing in



the plume wake zone in the present model.  However it should also be



noted that the experimental results of Chan et al. may have been influ-



enced by the blockage of the flow by the model towers due to the finite



width of the experimental flume.



     The field results from TVA at the Paradise power plant include two



cases.  One is for a single tower in a stable ambient (TVA-11, potential


                       3fl                     3A
temperature  gradient  •£ = 0.59 °k/100m, 0 < ~ * 1.0 °k/m).  The other
                       0Z                     9Z


is for  two towers in an ambient with a temperature inversion (TVA-14,


96
"fa* 1-42  °k/100m > 1.0 °k/m).  Since only the average temperature



gradient and average wind velocities at a few levels were



available only rough estimates of ambient temperature and wind



velocity profiles were constructed based on the limited data.

-------
        200
         150
              D
              100
                                                 50
                                                     = 50
                                     •TRACE OF JET BOUNDARY
                                       (0.51cm ORIFICE, PHOTO N0.7650.l)~|

                                 _• _ PRESENT THEORY

                                  o  EXPERIMENT (FAN,1967)
                                 	  FAN'S THEORY
           I   1   I  I

                                       I   i  i
                                                                 50
                                                      y
                                                      D
                                                                 100
                             F=20   k=8
                                  60 (d)
 100
 80
 60
 40
 20
  0
         • PRESENT THEORY
         o EXPERIMENT (FAN, 1967)
    0
20
40
                                         8°
                                      100
120
140
Figure 4.13  Comparisons of Plume Trajectory, Width and Dilution Between the
Figure          P    Theory and Fan's (1967) Experiments for F = 20 and K = 8
                                45

-------
250           200
    I   I  I   i   |   '   r
                                     x
                                     D
                              150
                     100
                          50
                                                             =50
                                                = 100
                                 = 150
                              TRACES OF JET BOUNDARIES (0.51cm ORIFICE)
     _•_ PRESENT THEORY    ~-*~ STATIONARY JET EXPERIMENT (PHOTO NO. 7646.4)
     O   EXPERIMENT (FAN ISST)**'^ TOWED JET EXPERIMENT (PHOTO NO. 7656.4)
     	 FAN'S THEORY
     I  I  I   I  I   I   i  I   I   I   I   i   I  I  I   i   I   I   i  I   i   i   i
                                                   50
                                                                          100
100
 80
 60
 40
 20
                             F=40  k=8  60(g)
         • PRESENT THEORY
         o EXPERIMENT (B\N, 1967)
            20
40
                                •
                                o
60        80
    S/D
                                                   100
120
140
  Figure 4.14  Comparisons  of  Plume Trajectory, Width and  Dilution Between the
               Present Theory  and Fan's (1967) Experiments for F - 40 and K - 8
                                 46

-------
  250
                                       x

                                       D
200
                               150
                             100
                          50
                    I   i   i   I
        - PRESENT THEORY

       O  EXPERIMENT (FAN,1967)

          FAN'S THEORY
                                                                  = 50
                                                   = 100
                       = 200
                    = 150   TRACES OF JET BOUND ARIES (0.51cm ORIFICE)
                      /—/STATIONARY JET EXPERIMENT(PHOTON0.7646.li
                      X..ATOWED JET EXPERIMENT (PHOTO NO. 7657.4)
                                I
                                                           50
                                                                            100
                               F=80  k = !6   60(j)
  100
   80
   60
 Q

QO
  40
   20
   0
           • PRESENT THEORY

           O EXPERIMENT (FAN, 1967)
                       •
                       o
                                 •
                                 o
                    •
                    o
              20
       40
60S/D   8°
100
120
140
  Figure 4.15   Comparisons of Plume Trajectory, Width and Dilution Between the
                Present Theory and Fan's  (1967) Experiments for F = 80 and k = 16
                                  47

-------
*«
00
                      0
                    -2
                                  1    I   T
I    I   I    I   I    I   I    I    I   I    T
                                 	EXPERIMENT

                                 	PRESENT THEORY
                               PLANE
                               MEASURED
                              I    I   I    I    I   I    I   I    I    I   1    I   I
                              0
                                       8
               Figure  4.16  Comparisons of  Plume Excess Temperature Between the Present Theory and Chan et al.'s

                           (1974) Experiments  for F = 4 and K =  1.02

-------
Three ambient relative humidity profiles (100%, 70%, 0%) associated with




the relative temperature profiles were generated and tested.  The exit




plume humidities were assumed 100% (saturated) except for one dry plume




case which is 0% for ambient and exit humidities.  The input data cards




are presented in Tables 4.4 and 4.5.  The predicted results and the compari-




sons are shown in Figures 4.17 and 4.18.  From the variations of plume




trajectory for different ambient humidity conditions, the effect of the




ambient humidity can be seen to be quite significant.  Similar conclusions




could be drawn for the effect of ambient temperature and wind velocity.




The present model overpredicts the plume trajectories for TVA's cases.




This could be due to the incomplete information of the ambient conditions




and the neglect of drift in the tower initial conditions.  Adequate




ambient and source conditions are mandatory for proper model validation.
                                49

-------
 01
 o
    600
500
    400
     300
    200
     100
              o  EXPERIMENT(TVA-II,I968)
            PRESENT THEORY
                   A
                   •
                   V
                   D
                  EXIT HUMI
                    100%
                    100%
                    100%
                      0%
AMB HUMI
  100%
   70%
    0%
    0%
                        200
                             400
                  600
               x(m)
                                                                       800
1000
Figure 4.17
Comparisons of Plume Trajectories Between the Present  Theory and TVA  (1968, TVA-11, single tower)
Field Data in a Stable Atmospheric Condition (39/3z-0.59°K/100 m)

-------
       500
       400
        300
        200
        100
         0
   o  EXPERIMENT (TVA-14,1968)
PRESENT THEORY   EXITHUMI  AMB HUMI
       A            100%
       •            100%
                    100%
                     v
                     G
 100%
 70%
  0%
                      0%
  0%
           0
             200
400
                                                       x(m)
600
800
1000
Figure ^.18  Comparisons of  Plume Trajectories Between the Present Theory and TVA (1968, TVA,  two towers)
            Field Data in an Atmosphere with an Inversion (36/^z=l.42°K/100 m)

-------
Table 4.1  Input Data  Cards  for  Example Cases (Different Wind Directions to a Line
           Array  (3  cases) and Round Array  (one case) of Towers)
"CASE  (1)   * TOWERS IN ONE ARRAY   0 DEGREES TO  WIND
     4  120   30  30  11
— -0-" 003 10 011111-
.00000 11. 45000 2?. 90000 34.35000
• .00000 .oonoo .00000 .00000
.00000 100.00000 200.00000
10.50000 10.30000 10.10000
70.00000 70.00000 70.00000
4.10000 5.22504 5.76969 4.10000 5.22504
5.76969 4.10000 5.22504 5.769f>9
	 9.44880 10.26800 31.90000 .02821 .00000
5.00000 5.00000 n 0 0 0
EXCS TEMP 4 TWRS 1 ARRAYS 0 OFG TO WIND 2/0
tXCS HUMI 4 TWRS 1 ARRAYS 0 DEG TO WIND Z/D
EXCS LIQM 4 TWPS 1 ARRAYS 0 DEG TO WIND Z/D
"CASE (2)" 4 TOWERS IN ONE ARRAY ' 90 DEGREES TO WIND
4 120 30 30 11
	 0 0 0 3 10 0 1 111 1
.00000 .00000 .00000 .00000
.00000 11.45000 2?. 90000 34.35000
.00000 loo.ooooo 200.00000
lO.bOOOO 10.30000 10.10000 	 	 ~ "" '"
70.00000 70.00000 70.00000
4.10000 5.22504 5.76969 4.10000 5.22504
5.76969 4.10000 5.22504 5.76969
9.44880 10.26800 31.90000 .02821 .00000
b. 00000 5.00000 0000
~F"X<~S TEMP 4 TWRS 1 ARRAYS 90 OG TO WIND 7/0
EXCS HUMI 4 TWRS 1 ARPAYS 90 OG TO WIND Z/D
EXCS LIDM 4 TWRS 1 ARPAYS 90 DG TO WIND Z/0
CASE (3) 4 TOWERS IN ONE ARRAY 135 DEGREFS TO WIND
4 120 30 30 11
	 0 0 " 0 3 10 0 •• — 1 	 1 111
.00000 8.10000 If). 19000 24.29000
'-" "".00000 fl. 10000 1ft. 19000 24.29000
.00000 100.00000 200.00000
lO.bOOOO 10.30000 10.10000
70.00000 70.00000 70.00000
' 4V10000 5.22504 5.76969 4.10000 5.22504
5.76969 4.10000 5.22504 5.76969
-9.44880 10.26800 31.900CO .02821 .00000
5.00000 5.00000 0000
EXCS TTMP 4 TWRS 1 ARRAY 135 DEG TO WINDZ/D
tXCS *JMI 4 TWRS 1 ARPAY 135 DEG TO WINr,Z/D
"EXCS'LIQM 4 TWRS 1 AROAY 135 DEG TO WINDZ/D
CASE (4) 5 TOWERS IN ROUND ARRAY
- - 5 100 30 30 11 	
0 0 0 3 10 0 0 1 1 1 1
	 .00000 5.00000 10.00000 15.00000 20.00000
10.00000 .00000 20.00000 2.00000 12.00000
	 '.ooooo 100.00000 200. ooooo - - — •-•-
10.50000 10.30000 10.10000
	 70.00000 70.00000 70.00000
4.10000 5.22504 c. 76969 4.10000 5.22504
5.76969 4.10000 5.22504 5.76969 4.10000
9.44880 10.26800 31.90000 .02821 .00000
	 5.00000 5.00000 0 0 00
EXCS TEMP 5 TWRS IN ROUND ARRAY Z/D
EXCS HUMI 5 TWRS IN ROUND ARRAY Z/D
EXCS LIOM 5 TWRS IN ROUND ARRAY Z/D






5.76969 4.10000



X/D
X/D
X/D








5.76969 4.10000



X/D
X/D
X/D








5.76969 4.10000



X/D
X/D
X/D








5.76969 4.10000
5.22504 5.76969


X/D
X/D
X/D
                                                                            5.22504
                                                                            5.22504
                                                                            5.22504
                                                                            5.22504
                                       52

-------
Table 4.2  Input  Data Cards for Three Cases of Fan's (1967)  Experiments with
           F = 20,  K =  8; F = 40, K = 8; and F = 80, K -  16
  CASE (1)  60(0)  F=20  K=8

   ~ 1 350   1   1   1
     000210
      .0
    "~ Y00000~-50.00000
    31.87000  31.87000
      .00000    .00000
      .13700    .13700
      'V00760" -1.10000 "20.00000     700000.00000
  CASEJ2]	60
-------
    Table 4.3  Input Data Cards for Chan et al.'s (1974) Experiment
               with F = 4, K = 1.02
Case (1) 6 Towers in One Array 90 Degrees to Wind F=4, K=1.02
 6 100   21  17  11
 00    120
0
1  1
0  1



25





6
EXCS
.00000
.00000
.00000
.00000
.00000
.37621
.37621
.05690
.50000
.31250
TEMP 6
.
,
900.
25.
.
.
.
.
.
5.
TWRS
00000
06501
00000
00000
00000
37621
37621
38374
50000
06250
.00000
.13003



.37621
.37621
30.00000
1.00000
0 0
1 ARRAY 90 DEC
.00000
.19504



.37621
.37621
.00000
3.00000
0 0
TO WIND Z/D
.00000
.26006



.37621

.00000



.00000
.32507



.37621




X/D
                                                  .37621   .37621   .37621
                                                     X/D A3=0.3536 CD=1.5
                                 54

-------
Table 4.4  Input Data Cards for  TVA (1968, TVA-11, single tower) Field Data in
	a-Stable-Atmospheric Oondlt:ion-<^e/9z-fr;»9J>K/1^0 m)

	CASE.  (1)   rOO* EXIT HUMI~8T100« "AMB'HUMI
      1  300
               1
               8
                          0
           1
  000
~~  .00000                              .    - -
	.00000
    rOOOOO~I00700000"150.00000  200.00000  250.00000 300.00000 350.00000 400.00000
  20.00000   20.16000   20.24000   20.32000   20.40000  20.48000  20.56000  20.64000
100.00000  100.00000  100.00000  lOO.OOOnn  100.00000 100.00000 100.00000 100.00000
  4.10000    4.10000    4.10000    4.70000    5.30000   5.90000   6.50000   7.10000
"7.90000   17.10000  139.00000     .67883     .00000
 ;ASE
                       HUMI 5-
      1  300   1   1   1
  -•- — Q--- o~  0   8~~~1 ~"
-------
         Tnput Data Cards  for- TVA  (1968, TVA-1&, two towers) .Elfild_Dat«. Iti an
         Atmosphere with an  Inversion  (39/3z-1.42°K/100 m)
CASE" (1) "10095 EXIT HUMI & 100% AHB HUMI~
  ~2 300"~ 1  "1"~"1                                            "            ""
   00061010000
   -.00000  29.06000      '	                          "
    •OOOOO  54.65000
  —;ooooo  50.00000 100.00000 iso.ooooo 200.00000 250.00000
   7.00000   7.43000   7.86000   8.29000   8.72000   9.15000
  LOO.OOOOO 100.00000 100.00000 100.00000 100.00000 100.00000
   8.60000   8.80000   9.10000  10.50000  11.60000 J2.60000    8.60000    8.80000
  "^9.10000  10.50000  11.60000  12.60000   ~""
   7.90000  20.20000 149.00000    .73951    .00000
               EXIT HUMI"
               HUMI
2
0
•
•
*
7.
70.
8.
9.
7.
300 1
0 0
T)0000 "
ooooo
ooooo*
ooooo
ooooo
60000
10000~~
90000
1 1
6 1
29.06000
54.65000
50.00000
7.22000
70.00000
8.80000
10.50000
20.20000
0

100.
7.
70.
9.
11.
149.
1 0

OOOOO ]
44000
OOOOO
10000
60000 ~
OOOOO
0 0
[50.DOOOO
7.67000
70.00000
10.50000
12.60000
.73951
0
200
7
70
11

.00000
.89000
.00000
.60000
.00000


250.00000
8.11000
70.00000
12.60006


8.60000
                                                                          8.80000
CASE (3)  100* EXIT HUMI &   0% A"B
 2"300 —
 0    0
  .00000
  .00000
" .00000
 7.00000
~.ooooo
 8.60000
 9.10000
 7.90000
           1   1"  "1
           061
            29.06000
            54.65000
            50.00000
             7.22000
              .00000
             8.80000
            "10.50000
100.00000
  7.44000
—.ooooo
  9.10000
"11.60000
150.00000 200.00000 250.00000
            20.20000 149.00000
  7.67000
   iOOOOO
 10.5000ft
 12.60000
   .73951
 7.89000
  .00000
11.60000

  .00000
 8.11000
  .00000
12.60000
8.60000   8.80000
            "0* EXIT HUMI I   0% AMB "HUMI
2 JOO ]
0 0 (
.00000
.ooooo
.uoooo
7.00000
.00000
8.60000
9.10000
L 1 1
) 6 1
29.06000
54.65000
'50.00000
7.22000
.00000
8.80000
10.50000
0 1
100.00000
7.44000
.00000 ~
9.10000
11.60000
0 0
150^00000"
7.67000
.OOnno
10.500QO
12.60000
0
200.00000
7.89000
.00000
11.60000
250.00000
8.11000
.00000
12.60000
   7.90000  20.20000  149.00000
             .00000
             .00000
                                                                8.60000    8.80000
                                     56

-------
                             CHAPTER 5




                  CONCLUSIONS AND RECOMMENDATIONS




     In this study, a mathematical model and corresponding computer



program have been developed for the prediction of plume behavior from



multiple cooling towers.  Some comparisons between the predictions



based on the present model and the measured results from laboratories



and the field are made in order to  test the model.  The following



conclusions and recommendations are made based on this study.



     (1)  The model is developed for arbitrary vertical profiles of



ambient temperature, humidity, wind velocity, and arbitrary tower



arrangement.  The velocity defect for the downstream towers due to the



effect of the upstream towers and plumes can be included by specifying



different ambient velocity profiles for each plume.  A general expression



for the velocity defect of the downstream towers might be developed in



the future.


     (2)  The temperature range for which this model is valid is -50°C



to  140°C, because of the accuracy associated with the calculation of



the saturation humidity.


     (3)  A set of suggested values of entrainment and drag coefficients



have been incorporated in the computer program.  Because of the rapid



merging and usually rapid bending over of the plumes, the coefficients



a~  a  and C. are the most important ones.  Better estimates of their
 3s      a


values are needed such as by further experimental study or field program.



     (4)  The merging criteria and processes (defined in this model by



equations (2.38) (2.39) and (2.40)) could also be improved when further





                               57

-------
research results on plume interaction are available.




     (5)  The blockage and recirculation effects in the wake zone of the




towers and plumes have not been incorporated in the present model.




Future effort could be made to include these effects.




     (6)  Based on comparisons between model and laboratory results




(Fan, and Chan et al.), good predictions of dry plume behavior can be




obtained.  In order to verify the model for actual cooling tower plumes,




a more complete set of experimental data (including the plume width,




trajectory, dilution and detailed ambient profiles of temperature,




humidity and wind velocity)are required.  Therefore, a complete set of




field measurement are strongly recommended for validation of the model.
                                 58

-------
                            REFERENCES
Abraham, G.  (1970),  "The Flow of Round Buoyant Jets  Issuing Vertically
     into Ambient Fluid Flowing in  a Horizontal Direction," Advan. Water
     Poll. Res. Proc.  Int.  Conf. Water Poll. Res.  5th  Paper 111-15.

Briggs, G. A.  (1969),  "Plume Rise," AEC  Critical  Review  Series  Report.
     Report  number TID-25075.

Briggs, G. A.  (1974),  "Plume Rise from Multiple Sources,"  Cooling  Tower
     Environment-74.   CONF-74032, ERDA Symposium  Series.   National
     Technical  Exchange Service.  U.S. Dept. of Commerce,  Springhill, VA.,
     pp. 161-179.

Carpenter, S. R., F. W. Thomas, and R. E.  Gartrell (1968), "Full-Scale
     Study of Plume  Rise at Large Electric Generating  Stations," TVA,
     Muscle  Shoals.  Alabama.
Chan,  T. L., S.-T. Hsu, J.-T. Lin,  K.-H. Hsu,  N.-S.  Huang, S.  C. Jain,  C.  E.
     Tsai, T.  E.  Croley II, H.  Fordyce and J.  F.  Kennedy,  "Plume Recirculation
     and Interference  in Mechanical Draft Cooling Towers," Iowa Inst. Hyd.  Res.
     Rep. No.  160, 41  pp.

Csanady, G.  T.  (1971), "Bent-Over Vapor  Plume," J. Appl. Meteor.,  10,
     34-42.

Davis, L. R.  (1975), "Analysis of Multiple Cell Mechanical Draft Cooling
     Towers," Environ. Prot. Agency Rep.  Office of Research and Develop-
     ment, Ecological  Research Series, EPA-66013-75-039.

Fan, L. N. (1967), "Turbulent Buoyant Jets into Stratified or Flowing
     Ambient Fluid," California Institute  of Technology, W. M. Keck
     Laboratory of Hydraulics and Water  Resources, Rep. No. KH-R-15.

Fox, D. G. (1970), "Forced  Plume in a Stratified  Fluid," J. Geophys.
     Res., 75  (33), 6818-35.

Hanna, S. R.  (1972), "Rise  and Condensation of Large Cooling Tower Plumes,"
     J. Appl. Meteor., 11,  793-799.

Hirst, E. A.  (1971), "Analysis of Round, Turbulent, Buoyant Jets Dis-
     charged to Flowing Stratified  Ambients," Oak Ridge National
     Laboratory,  Report Number ORNL-4685.

Hoult, D. P., J. A. Fay, and L.  J.  Forney  (1969),  "A Theory of Plume
     Rise Compared with Field Observations," J. Air Pollut. Contr.
     Assoc.  19(9), 585-90.

Jirka, G. and D.R.F.  Harleman (1974),  "The Mechanics of Submerged Multi-
     port Diffusers  for Buoyant Discharges in Shallow Water," MIT Ralph
     M. Parsons Lab.  for Water Resources and Hydraulics,  Report Number 169.


                                 59

-------
Koh, R.C.Y. and N. H. Brooks  (1975), "Fluid Mechanics of Waste-Water
     Disposal in the Ocean," Ann. Rev. Fluid Mech., 7:187-211.

Koh, R.C.Y. and Y. C. Chang  (1973), "Mathematical Model for Barged Ocean
     Disposal of Wastes," Environ. Prot. Agency, Office of Research and
     Development, Environ. Prot. Tech. Series, EPA-660/2-73-029.

Koh, R.C.Y. and L. N. Fan (1970), "Mathematical Models for the Predic-
     tion of Temperature Distributions Resulting from the Discharge of
     Heated Water into Large Bodies of Water," Environ. Prot. Agency
     Rep. 16130 DWO 10/70, 219 pp.  (Also Tetra Tech, Inc., Rep. TC-170).

Linsley, Jr., R. K., M.  A. Kohler, and J.L.H. Paulhus (1975), Hydrology for
     Engineers , McGraw-Hill Book Company, Inc., New York, NY.

List, E. J. and J. Imberger  (1973), "Turbulent Entrainment in Buoyant
     Jets and Plumes," J. Hydraul. Div., Proc. ASCE, 99:1461-74.

Meyer, J. H., T. W. Eagles, L. C. Kohlenstein, J. A. Kagan, and W. D.
     Stanbro (1974), "Mechanical Draft Cooling Tower Visible Plume
     Behavior:  Measurements, Models, Predictions," Cooling Tower
     Environment-74.  CONF-74032, ERDA Symposium Series.  National
     Technical Exchange Service, U.S. Dept. of Commerce, Springhill, VA,
     pp. 307-352.

Morton, B. R. (1957), "Buoyant Plumes in a Moist Atmosphere," J. Fluid
     Mech., 2, 127-144.

Morton, B. R., G. I. Taylor and J. S. Turner (1956), "Turbulent Gravita-
     tional Convection from Maintained and Instanteous Sources," Proc.
     Roy. Soc. London, Ser. A 234:1-23.

Richards, J. M. (1963),  "Experiments on the Motion of Isolated Cylin-
     drical Thermals through Unstratified Surroundings," Intern. J. of
     Air and Water Pollut., 7, pp. 17-34.

Richards, J. M. (1971),  "Simple Expression for the Saturation Vapor
     Pressure of Water in the Range -50° to 140°C," Brit. J. Appl. Phys.,
     4, L15-L18.

Schatzmann, M. (1977), "A Mathematical Model for the Prediction of Plume
     Rise in Stratified Flows," Sonderforschungsbereich 80, University
     of Karlsruhe, W. Germany.

Slawson, P. R. and G. T. Csanady  (1967), "On the Mean Path of Buoyant,
     Bent-Over Chimney Plumes," J. Fluid Mech., 28, 311-322.

Slawson, P. R. and G. T. Csanady  (1971), "The Effect of Atmospheric
     Conditions on Plume Rise," J. Fluid Mech., 47, 33-39.
                                60

-------
Slawson, P. R., J. H. Coleman, and J. W. Frey (1975), "Some Observations
     on Cooling-Tower Plume Behavior at the Paradise Steam Plant,"
     Cooling Tower Environment—1974 (ERDA Symp. Series: CONF-740302),
     147-160.

Weil, J. C.  (1974), "The Rise of Moist, Buoyant Plumes," J. Appl. Meteor.,
     13, 435-443.

Wigley, T. M. L.  (1974), Comments on "A Simple but Accurate Formula for
     the Saturation Vapor Pressure over Liquid Water," J. Appl. Meteor.,
     13, 608.

Wigley, T. M. L.  (1975), "A Numerical Analysis of the Effect of Condensation
     en Plume Rise," J. Appl. Meteor., 14, 1105-1109.

Wigley, T. M. L.  (1975), "Condensation in Jets, Industrial Plumes and Cooling
     Tower Plumes," J. Appl. Meteor., 14, 78-86.

Wigley, T. M. L.  and P. R. Salwson  (1971), "On the Condensation of Buoyant,
     Moist Bent-Over Plumes," J. Appl. Meteor., 10, 253-259.

Wigley, T. M. L.  and P. R. Salwson  (1972), "A Comparison of Wet and Dry
     Bent-Over Plumes," J. Appl. Meteor., 11, 335-340.

Wright, S. J.  (1977), "Effects of Ambient Crossflows and Density Strati-
     fication of  the Characteristic Behavior of Round, Turbulent Buoyant
     Jets,"  California Institute of Technology, W. M. Keck Laboratory of
     Hydraulics  and Water Resources, Report No. KH-R-36, 254 pp.
                               61

-------
                            APPENDIX A


                         COMPUTER PROGRAM


     The computer program based on the model and listed in Appendix C

was tested on an IBM 370/158.  The detailed input and output information

are listed in this Appendix.  The input ambient wind velocity profiles

(AU(NP,MG)) for each tower are designed to allow consideration of the

velocity defect in the wake of upstream towers (i.e., the tower array

parallel to the ambient wind direction).  In addition, some suggested

input values are listed below for reference:

                 NS = 300
                               »
                 NX = 40

                 NY = AO

                 NCONT =11

                 IX(3) = IX(6) = IX(ll) =0

In this Appendix, the input sequence as well as the input and output

variables are tabulated, explained and related to the symbols used in the

text of this report.
                                A-l
                                 63

-------
                                  INPUT SEQUENCE
Symbol

NP,NS,NX,NY,NCONT
(IX(I),1=1,11)
(CX(I),I=1,NP)
(CY(I),I=1,NP)
(AZ(I),I=1,MG)
(AT(I),I=1,MG)
(AH(I),I=1,MG)
((AU(I,J),J=1,MG),I=1,NP)
DI(l),UO(l)f TO(1),HO(1),WO(1)
(DI(I),1=1,NP)
(JO(I),I=1,NP)
(TO(I),I=1,NP)
(HO(I),I=1,NP)
(WO(I),I=1,NP)
A1,A2,A3,A4,CD,TURBF
DX,DZ,XO,ZO
WIDTH, !IITE,MORE,NOMAP,ICENT
NOTICK
(HEDN(!.),k=l,10),(LABY(L),L=l,5),(LABX(M),M=l,5)
(HEDN(k),k=l,10),(LABY(L),L=l,5),(LABX(M),M=l,5)
(HEDN(k),k=l,10),(LABY(L),L=l,5),(LABX(M),M=l,5)
Parameter   Format   Subroutine








IX(1)=0+
IX(1)=1*
IX(1)=1*
IX(1)=1*
IX(1)=1*
IXdH*
IX(2)=1*
IX(3)=1*
IX(11)=1*
*
IX(8)=1*
IX(9)=1*
IX(10)=1*
514
1114
8F10.5
8F10.5
8F10.5
8F10.5
8F10.5
8F10.5
8F10.5
8F10.5
8F10.5
8F10.5
8F10.5
8F10.5
8F10.5
4F10.5
2F10.5
414
20A4
20A4
20A4
CTPS
MTP
MTP
MTP
MTP
MTP
MTP
MTP
MTP
MTP
MTP
MTP
MTP
MTP
MTP
OUTPUT
OUTPUT
OUTPUT
OUTPUT
OUTPUT
+ Skip card if IX(1)=1
* Skip card if the corresponding IX(I)=0, (1=1,2,3,11,8,9,10)
                                         A-2
                                          64

-------
In Program

NP

NS

NX, NY
[CONTA(NX,NY)]
NCONT
[CONTB(NCONT)]
EXPLANATION OF THE INPUT SYMBOLS

In Text    Remarks

           Total number of towers

           Desired number of calculation steps

           Horizontal and vertical grid sizes, respec-
           tively (for contour plot)

           Desired contour levels for plotting
IX(2)
IX(3)
IX (4)

IX(5)

IX(6)
IX(7)


IX(8)



IX(9)
IX(ll)
                             IX(1)=0
           IX(2)=0

           IX(2)=1


           IX(3)=0

           IX(3)=1

           IX(3)=-1

           MG


           INRPR


           IPNT=0

           LC=0

           LC=1

           IX(8)=0

           IX(8)=1


           IX(9)=0

           IX(9)=1


           IX(10)=0

           IX(10)=1


           IX(11)=0
                                 A- 3
                                   65
 Same  exit  conditions  of  all  the
 plumes  (Input  one  card only)
 Different  exit conditions  of  the
 plumes

 No  input card  of entrainment  and
 drag  coefficients  is  needed
 Input the  desired  entrainment  and
 drag  coefficients

 No  input card  of DX,  DZ, XO  and ZO
 is  needed
 Input the  desired  values of  DX, DZ,
 XO  and  ZO
 No  plot needed

 Number  of  vertical levels  for
 ambient conditions
 Interval of detailed printout  for
 plume 1
 For contour plot (always use  IPNT=0)

 For cluster array  (round array) of
 towers
 For line or random array of towers

 No  contour map plotted for plume
 excess temperature
 Contour map plotted for plume
 excess temperature

 No contour map plotted for plume
 excess humidity
 Contour map plotted for plume  excess
 humidity
 No contour map plotted for plume
 excess liquid moisture
 Contour map plotted for plume  excess
 liquid moisture

No input card of WIDTH, HITE,  MORE,
NOMAP, ICENT, and NOTICK is needed
 Input the desired values of WIDTH,
HITE, MORE, NOMAP,  ICE1TT and NOTICK

-------
In Program

CX(NP)

CY(NP)

AZ(MG)

AT(MG)


AH(MG)



AU(N'P.MG)


DI(NP)

UO(NP)

TO(NP)

HO(NP)


WO(NP)
DX.DZ



XO,ZO




WIDTH,HITE


MORE


NOMAP


ICENT


NOTICK
                 In Text
                   x

                   y
                   z
                   ta
                   Dc
                   U
A1,A2,A3
and A4
CD
TURBF
al»a2»a3
and a i«
Cd
U'
Remarks
x-coordinates of towers (in m)
y-coordinates of towers (in m)

Ambient levels (in m)

Ambient temperature profile corresponding to
AZ(I) (in °C)         "
Ambient relative humidty profile correspond-
ing to AZ(I) (in percentage of the saturation
humidity)
Ambient wind velocity profile corresponding
to AZ(I) for each tower (in m/sec)

Diameter of each tower (in m)

Exit velocity of each plume (in m/sec)
Exit temperature of each plume (in °C)

Exit specific humidity of each plume  (in
kg/kg)
Exit liquid phase moisture of each plume
(in kg/kg)
Entrainment coefficients (Default: Al=0.0806,
A2=0.4775, A3=0.3536, A4=0.)

Drag coefficient (Default: CD=1.5)

Intensity of ambient turbulent fluctuations
(in percentage, decimal; Default: TURBF=6.)

Increments of grid size in x and z directions,
respectively (Normalized by the diameter of
the first tower; Default:DX=0.5, DZ=0.5)

Location of the center of the top of  the
first tower in the grid (Normalized by  the
diameter of the first  tower; Default:XO=1,
Z0=2.)

Width and height of contour map, respectively
(Inches; Default: 8",  8")

MORE=1    Do not finish off the map (Default:
          'MORE=0)
NOMAP=1   Do not force grid to be square
          (Default:NOMAP=0)

ICENT=1   Do not center the title  (Default:
          ICENT=0)

NOCITK=1  Do not draw  tick marks  (Default:
          NOTICK-0)
                                A-4
                                  66

-------
In Program       In Text    Remarks

HEDN(IO)                    40 characters to plot as title on top
                            (Default: Blank)

LABY(5)                     20 characters to plot as label on vertical
                            left  (Default: Blank)
LABX(5)                     20 characters to plot as label on horizontal
                            bottom  (Default: Blank)
                                  A-5
                                    67

-------
                 EXPLANATIONS OF THE OUTPUT SYMBOLS
                       In Text
                       Cd« Ua
In Program

EXIT,COEF,TWLC,AMBL,
INPR,IPNT,LC,ETPL,
EHPL.EMPL.COPL

A1,A2,A3,A4,CD,
TURBF

NTHP

CX.CY
[(CX(I).CY(I)]
DIA,VELO,TEMP,HUMI,
LPMO [DI(I),UO(I),
HO(I)WO(I)]


NSTEP*
X.Y.Z                  x,y,z
[PX(NP,NS),PY(NP,NS),
PZ(NP.NS)]

PTEMP,PHUMI,PLQDMOIST  tp,qp,ap
[PT(NP,NS),PH(NP,NS),
PW(NP.NS)]
Remarks
Correspond to IX(i) to IX(ll)
respectively
                         ,32,33 ,a«t ,   Same as input symbols
EXCEST, EXCESH
[PET(NP,NS),PEH
(NP.NS)]

PCROSEC
[BXZ(I)]

SLOTLEN[A(I)]
PANGLE[PCOS(NP,NS)
PSIN(NP.NS)]

DILUTN
[PDIL(NP.NS)]

PVELO [PU(NP,NS)]
                       At, Aq


                       HT/2

                         A

                         e


                       Q/QO

                        u
                                      The Nth plume
                                      Same as input symbols
Tower diameter, exit values of
plume velocity, temperature,
specific humidity and liquid phase
moisture, respectively
The Ntn step of calculation referred
to each tower
The horizontal, lateral and vertical
coordinates of plume center

Plume temperature, specific humidity
and liquid phase moisture respectively

Excess plume temperature and specific
humidity

Half height of plume cross-section

Finite length of slot jet of the
merged plume
The angle between the tangent of
plume trajectory and the horizontal
line
Plume dilution
Plume velocity
  NSTEP refers to the number of calculation step of plume 1 when each
  plume first appeared
                                  A-6
                                   68

-------
                               APPENDIX B

         EXPLANATION OF THE IMPORTANT SYMBOLS IN THE PROGRAM MTP
In Program
PX(NP,NS)
PY(NP.NS)
PZ(NP,NS)
PS(NP,NS)
PCOS(NP.NS)
PSIN(NP.NS)
PQ(NP.NS)
PMX(NP.NS)
PMZ(NP.NS)
PF(NP.NS)
PG(NP,NS)
PH(NP,NS)
PW(NP.NS)
PT(NP,NS)
PET(NP.NS)
PEH(NP,NS)
PAN
   sin
   2-BY
    TI
    g
    al
    a2
    a3
    aA
Excess temperature
Excess humidity
    *
Net plume velocity
Average plume width
    *
Plume dilution
Normalized PET.PEH & PW (or PEW)
Contour levels
Indication of the status of each plume
Plume which has not been started
Single plume
Merged plume
All plumes are merged
Beginning step number for each plume
Merged plume pair
Merged plume step numbers associated
   with MP(I)
    *
Beginning step number for visible plume
Ending step number for visible plume
3.1415926
    *
    *
    *
    *
                                B-l
                                 69

-------
In Program
In Text
Remarks
CD
TURBF
UC
B
TP
HP
WP
ET
EH
MG
CFRL
ALV
CPA
ANG
ITHP
ENTRAN
IQ
IL
IK
KE(I)
    UT
 Fr
   LC
 -pa
    6

    E
            <3

             =3
            >3
YP(I)
YR1(I),YR2(I),
& YS(I)

YR1P(I),YR2P(I),
& YSP(I)
                   YCD-Q
                   Y(2)=Mcos8
                   Y(3)=Msin6
                   Y(4)=G
                   Y(5)=F
                   Y(6)=x
                   Y(7)=z
   it
   *
   *
Plume width
   *
   *
   *
Excess temperature
Excess humidity
Elevation level
   *
   *
   *
   *
ith plume
   *
All the plumes have not been
   completely merged
All the plumes have merged
All the plumes have merged and
   become a round plume again
Number of merged pairs
Plume step number
Ending step number of each plume
                       *
                       *
                    Derivatives of Y(I)  with respect to s
                    Y(I)  associated with the two half
                       round plumes and  the central slot
                       plume for the merged plume
                    YP(I) associated with the two half
                       round plumes and  the central slot
                       plume for the merged plume
* Refer to "List of Symbols"
                                B-2
                                 70

-------
    APPENDIX C
LISTING OF PROGRAM
         C-l
           71

-------
FORTRAN IV G LEVEL 20.7 VS
                       PAIN
DATE
6/07/77
l
-------
 FORTRAN IV G LEVEL 20.7 VS

  0001
                                         MTP
                                                           DATE =
                                                                     6/07/77
                                                                                  14:32
  0002
  0003
  0004
                   A
               ----- 3
                   4
                   5
                   6
                    SUBROUTINE  MTP( PX , PZ ,PO, PKX, FMZ, PG, PF, PCOS, PS IN, PENT, PU.PS.PB.P A,
                   lPI*P^IjP-t*£E±UP-i(*£Atlj£DU.,PY,PiL,CCNTA,CCNT6,CONTC,NP, NS ,NX ,NY ,
                   2NCONT.N4)
                    EXTERNAL DERIVE, DERI VS.DERIVF
                    DOUBLE PRECISION AQ, BO, CO, DO
                    COKMCK ySTCR£l/lND130),INDT(30),NOVlKl30),IS<30)
                   1       /STORE2/KP(30),MS(30),IIP(30)
                                                                         ,AU13Q.3Q1 ,
                                    AUGC30.30)
                           -^TflR£4/AL30),B_ll30),B2l30),BXZl30),BYI30)fPMCOS(30),
                                    PMSINC30),DW<30) ,NBV< 3d) ,NE V(30) tNCV(30)
                           -^OLNSIl/PAl.GRA.Al,A2,A3,A4.CD,TURaF,UC,B,TP,ET,HP,EH,WP,
                                    MG,MGl,CFRLtALV,CPA,ANG,IQ,ITHP,ENTRAN
                                      fNi itTn
                   8
                           /CCNST3/AQ,BQ,CC,DQ
                           V£CmTiJU/WIOTH,HITE,MOR£,NJMAP,ICENT,NOTICKtHEDNI10),
                                    LABY(5),LABX(5)
 0005
                    DIMENSION  PXCNP,NS),PZ«NP,NS),PQ(NP,NS),PMX(NP,NS) ,PMZ(NP,NS),
                   2PENT(NP,NS),PHINP,NS),PWINP,NSJ ,PEH ( NP.NS ) ,CX<301 ,CY( 30 ) , KE<30) ,
                   4Y(6),YP(8) ,YS(81,YSPC8I,YR1(8),YR1PI8),YR2(8),YR2P{8»,PDIHNP,NS
               .    -5CONTAlNX,AlYJJ.CONTfllNaiNT).a;NTClN4J,JX112J,MGPA2(30),MGSTll30J.
                   6 PC (NP , NS ) t LCHK ( 301 , AHP ( 30)
             X _ IMTIAII7F
 0006
 0007
                                                                                  __
                    CALL BLANKtPX,PZ,PGtPMX,PMZ,PB,PU,PFfPGPCOS,PSIN,PENT,PT,PS,PA7
                   lP£T«P-£H,PAN,PYfPOILtP-H,PW,PCfNP,NS,NX.NY,NCONT,N4,CONTA,CONTB,
                   2CONTCfMGPA2«MGSTl,DX,DZ,XOtZC,IPNT,KE,LCHKtAHP)
                   -4NPUT -CCNTAOL -P
                 	GO TO -35   	
                 34  READC5.1I  (DI( I),1=1,NP)
                                        C-3
                                          73

-------
FORTRAN IV G LEVEL 20.7  VS
                                        MTP
                                             DATE -
6/07/77
14:32
 0029
_QO30	
 0031
 0032
 0033
                   READ<5,1J  (TO(I),I*1,NP)
      READIS.ll  (HOU!,1*1.NP!
   35 JF1IXI2) ^LE. -0)  GO 10 32
C     INPUT  COEFFICIENTS  OF  ENTRAINMENT AND DRAG ;TURBF»TURB.  FLUCT.
£     TURBF=X(IN-DECIMAL) CF AMBIENT JUNO VELOCITY
      READ(5.1)A1,A2,A3»A4«CD,TURBF
 0035
-0036  _
 0037
-OO38
 0039
                  1  FORMATJ8F10.51
                -2  FORMATJ1214J —
                32  t,RlTE<£,555)
                5£5  FORMAT11H1J
                    MRITE(6.3)
                                     np TUP WAR i ABIES«j	
 0041
-0042

-0043-
 0044
-00.45.
 0046
-0047
                    WRITEC6.4)
                 -4  FORMAT-IIX^ lt£NGTJUJ!,  _T£MP^C, -M01STUREJKG/XG,   VELO:M/SEC,  -ANG1
                  1:DEG«)
                                  ^IStNXtNYtNCPNT	_.    .       		
                   WRITE(6t49)
                   if.a\ TCI (•t
                   WRITE(6.48I  IIX(J).J'1,11)
        	  _50 FORMATiiJUtN£U_CF-XO*ERS«i-,44,f—NO.-OF-CALCULAT10N STEPS-*.14.
                  1*  GRID  SIZE  «X,Yf=',I4,-  X',I4,«  CCNTOUR LEVELS**.I*//)
-0043	   49 fORMATiUC.lCONTKOL-P-ARAMETERS'J-  	    	
 0049           36 FORMAT*1  EXIT  COEF  TMLC  AHBL  INPR  IPNT   LC   ETPL  EHPL  '
	I'PMPt  rnpi • i	      	
 0050
-0051-
 0052
-0053
 0054
                48  FORMAT UX, 11114, 2X)//)
                  5  FORMATUX, 'COEFFICIENTS OF ENTRAINNENT  AND  DRAG* I
                  —WRITE U».AJAl-.A2«A3.*«,CDfTURBF   -----                       --
                  6  FCRMAT(1X,'A1»' ,F8.5,5X,'A2-',F8.5,5X,1A3«i,F8.5,5X,'A4-t,F8.5t5
 0055
-0056
 0057
-0058
 0059
                    hRITE«6,7»
                  7  FORMATiUC.iAMBJiliT _P-RCFJ-L£S'_I 	
                    WRITEI6.8)
                  B  FORMAT-i3X,LH£JGHT-!_,7J(,«J-EMP^12X, 'tlUMJOlTy • ,9X, 'VELOC ITV J
                    00 556 I-1,MG
 0061
-O062-
 0063
-0X344-
 0065
                    T»1.-373.16/TK
                    ES-0. 622*1013. 25*EST
                    AH(1I>AHP( II*0.01*FSH
                    pn o im\tv.rr
 0067
-0068
 0069
-0070
 0071
 nmy
                  9 WR1TE(6,10)AZ(I),AT(II,AH(I),AHP(I),AUI1,I)
                •OO-EORMAI41JUJ
                    HRITE(6.11I
                Jl -FORMATi/J	
                    HRITE(6.12)
                 i? PHRMATIi«t«Tnupp
 0073
-007*

-0075 -
 0076
 nn77
                    WRITEI6,18)
                   1*HUMI',9X,*LPMO*I
                ___ DO 47 J-.UNP -------  .-....-..
                 17 WRITE(6,13)I,CX(I).CY(I),DI(II,UO(I),TO(I),HO(II.MO(II
                                        C-4
                                          74

-------
FORTRAK IV G LEVEL 20.7 VS
                                        FTP
                                                          DATE
                                                                     6/07/77
                                                                                14:32
 -007d _. .. __.
              (
  0079
  0080
  0081
  0032
 -0083 	
  0084
  0085
  0086
  0087
  0086
 -DO 8 9	
  0090
                   1F12.8)
                   -WRITE(6, 242)	
                    CALCULATION OF TEMPERATURE,HUMIDITY
                    MG1=MG-1
                    DO 14 1=1,MG1
                    11=1+1
                    DV=AZ(II)-AZ(I)
                                                        AND VELOCITY GRADIENTS
                 14
                   AhG(I>«(AH(II)-AH(I))/DV
                   DO 47 J=1,NP
                   00 47  1=1, FG1
                   11=1+1
                   DV=AZ( II I-AZ(I)
  0092
  0093
  0094
-COS5
  0096
  0097

  0098
  OC99
-0100
                   DO 148 I»1,NP
                   -D02=01(IJ/2.
                   BXZ(I»=DC2
                   BYU)=002
                   B1(I)=C02
                   WRITE<6,15)
                15 FORMATUX, 'P-LUME     1  APPEARS AT NSTEP
             C     SET INITIAL CCNDITICNS,
                                                               -IS//)
                   ALV=FALV(TP»
                   CALL _S£LLf
                   INTEGRATION BY RUNGE-KUTTA  METHOD
0101
0102
0103
0104
-QIC 5 ....
1=0 ...
S=0.
DS=DI(ITHP)/20.
fin rn «;««
 0106

 0107

 01C8
-0109-
 0110
 -0111
 0112
 -0113
 0114
-0115-
 0116
                86 IF(NIP1 .GT. NP) GO TO  19
                  -CHECK IF -ANY-NE* PUJKE  APPEARS                                 __
                   CALL CHKNk.Pt PX.PZ,PQ,PMX,PMZ,PG,PF.PCCS.PS IN, PENT, PU, PS,P 8,P A,PT,
                —lPET.PH.PEh.PWtPAN,PY,POU_.PC,CX,NP,NS,DERlVR,CY,KE)
                19 IF(NII-1)16.105,105
                    ITHP=1
              ...  — ALV-FALVtTPJ
                    DS=0.1*PBt 1,
              _59fl-XALL RUNGS IS.
                    B1(1)=B
                                       yP-, L.DERI VR)
                   BXZ(1)*B
 -0118   —

-0119	

 -0120	


 0121
                   SOLUTIONS FOR SINGLE PLUME
                 -CALL SCLUTN(PXtPY.Pi.PQ.PMX,PMZ,PG,PF,PU,PENT,PS,PA,PB,PCOS,PSIN
                  lPT,PET,PHtPEH,PW,PAN,PDIL,PC,CY,Y,YP,S,ITHP,IKT,NP,NS)
                   r.n Tn -ju	.             	
                   CHECK IF ANY NEK PLUPE MERGING CCCURS
               105 CALL ALIGN(PJ(.PY.PX.PiJ,PKXtPMZ,PG,PF,PU.PENT.PS,PA,PB.PC.PCOSt
                  lPSIN,PT,PET,PH,PEH,PW,PAN,PDIL,CY,Y,YP,YS,VSP,YRl,YRlPfYR2,YR2P,
                  2DERIVE.DERIVS.DERIVR.S.LC.NP,NS,KE)
                   CALL PLMERG(PX,PY,PZ,PT,PET,Ph,PEH,PH,PA,PB,PAN,POIL,PCCS,PSIN,
                                        C-5
                                          75

-------
FORTRAN IV G LEVEL 20.7 VS
                                        MTP
DATE =
6/07/77
14:32:
 0122
.-0123	
 0124
 0125
 0126
                16 00 23 1=1, NI
                   1THP=1
             C
             C
                   IKT= IK+NCV IM I I-ISI I )»1
                   IFIINOm ^LT.U GO TO 23
                   DS«0.1*P3( I.IKT-1)
                   RESET INITIAL CONDITIONS AND CALCULATE PLUME PROPERTIES FOR  SINGLE]
                   AND MERGED PLUMES
 0128
                  lCERIVE.DEPIVS,DERlVRfPX,PZ,PQ,PMX, PMZ,PG,PF,PCOS,PSINtPENT,PUtPB,
                  2PA.PT.PET.PEH.PUtLCl
                   SOLUTICNS OF  PLUMES  INCLUDING  SINGLE  6 MERGED  PLUMES
                   CALL SCLUTNlPX,PY,PZ,Pg.P«X,PMZ,PGf PF, PU.PENT.PS, PA, PB.PCOS.PSIN,
                  lPT,PET,PH,PEHtPH,PAN,PDILtPCiCYtY,YPfSfIt IKT.NP.NSJ
                23— CCWIJNUE
0130
0131
0132
0133
0134
nj 3*
0136
0137
0138
4139
0140
0142
0143
0144
0145
0146
0147
0148
O149
0150
0151
0152
01*?
0154
-0155
0156
0157
0158
0159
0160
-0161
0162
4163
0164
0166
0167
0163
0169
0170
01 71
IL-0
78 DO 77 1=1, NP
IFUIK*NOVIKmj .LT. NSI GO TO 77
-INDUJ-0 	 -
77 CONTINUE
nr 7«; T=i,Kp
IFIIND(I) .GT. Oi GO TO 79
J5 CONTINUE 	 -
GC TO 80
75 JK=JK*1 - 	 _
IF(IO-2)86, 16,397
nr. nr HI 1 = 1, Np
IF(KE( I) .LE. 0) GO TO 82
	 -4FILCHK(I)^.LI-. -IJaO-I^J 81
KEU )=KEl I )•*•!
K=HGPA2(! )
PX( I,J)=PX(K,LI
PBCI !J)=PB
	 PC(I«JJ=PC1K,LJ 	
FCOS(I,J)=PCOS(K,LJ
PET(I,JI»PET(K,L)
Pf-H I ,JI«PEH(K,IJ
PH(ItJ)*PM(K,L)
GO TO 91
62 KEUI-NS
MRITE(6f808)
808 FODMAT( ix , 'BfSULTS OF THE VISIBLE 1
HRITEC6.225)
DO 807 I-1,NP
DO 803 J=1,KEN
IF IPHl 1 1 J) .GT.O.IGO TO 802
IF(NCV(II.EO.O)GO TO 803
NCVI IJ-0
CH«2. - ......
                                         C-6

                                           76

-------
FORTRAN

 0172
  0173
  0174
  0175
  0176
          IV  G  LEVEL  20.7 VS
                                         PTP
                                                           DATE =
                                                                     6/07/77
                                                                                  14:32:
  0178
  4179
 -0195
 0196
 -0197
 0198
—0199
 0200
                   f.r U( 11J J

             -B09 «RMAT(lX.2(I5,lXJ,3(i=8.2.1XJ,2(F7.2.1X),3(F8.5,lX),6(F7.2,lx),/l

               B02 JFINCVUJ..NE.OJGO TO B03
                   NBVU)*J
                          -J.	
                   CH=1
                            ,''',,
                         I,J|,PW(I,j),pBtIfJlfPCCI,J),PAU,J),PANU,J),PDlLU,J),
                 — -2PUI 1 1 J )
0180 241 FORMAT|1X,2( I5.1X) ,3(F8.2,1X) ,2 IF7.2,1X) 3IFB 5 IX) fe(F7 7 x 1
^?18J 803 CONTINUE ' * ' '-«fl*ll
0182
O183
0184
-0185 -
0186
.0187 .
0188
-0189
0193
-0191 —
0192
-0193
0154
IF(CH-1.) 804,805,807
-804 -WRITE 4 6, 233) I
233 FORMATUX, 'NO VISIBLE PART OF PLUME FOR PLUME*, I4/)
	 	 GC TJ3-807 	 	
8C5 HRITE(t.243)I,J
	 241 FflRMATJlXT IPI lIMft , J^ 1 CTli | VISIRIP AT T^JC Tea u I M« T t nm <-Tr«. ...
8C7 CONTINUE TIRMlNATJON ST£P',IW.
242 FORMAT*//)
224 FORMAT (IX, -RESULTS AT THE LAST STEP OF CALCULATION')
22£ FORMAT (IX.'NSTEP NTHP X Z Y PTEMP EXCEST
                  -!,•  PHUMI ---- EXCESH JU.QOMOIST-PAWHFWO PCROSEC SLOTLEN   PANGLE    .
                  2'DILUTK  PVELO'I
                   NTHP=i
                   wRITE(6,226INS,NTHP,PX(l,NS),PZll,NSI,PYIl,NS»,PT(l,NS),PET(l,NS)
                  2PDIL(1,NS),PU(1,NS»
              -226 FORMAT41X.2(J5*lXJ,3
-------
FORTRAN

 0001
IV G LEVEL 20.7 VS
                                        BLANK
                                                  CATE
                                                                    6/07/77
14:32:
           SUBROUTINE BLANK IPXtPZ,PQfPMX,PMZ,PB,PUtPFtPGtPCOSfPS IN,PENT,
          ]PTrPS,PAtPPTtPPH,PAktPYfPnil tPH.PH.PC.NP.NS.NX. NY.NCQNT.N4.CQNTA.
0002
0003

OOU4
0005
0006
0007
OOC8
.0009 .
0010
0012
0013 	
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
00^
0034
•0035
0036
0037
0038
0039



2CONTBtCONTC»HGPA2,MGSTl,DX,rt,XO,ZC,IPNTtKE,LCHK,AHP)
.DOUBLE PRECISION Ab, BO, CO, DO
CGPMON /STCRE1/IND(30),INCT(30>.NOVIM30),IS(30)
J. /STORE2/MP13QJ.MSC30I .I1PC30)
2 /STORE3/AH(30),ATI30»,A2<30I,AHG<30I,ATG(30),AU(30,30I,
A Aiirr(-»nf -*o^
3 /STORE4/A(30) ,B1 ( 30) ,62 ( 30) ,BXZ(30) , BY (30 J ,PMCOSt30» ,
4 - -J>MSJNL30J«OM(30J,NBV(30ltNEV(30) tNCVOO)
5 /CONSTi/PAI ,GRA,A1,A2,A3 , A4 .CD.TURBF ,UC, B,TP,ET, HP, EH, HP,
J> JlGtM£l*CFRL*ALV«CPAfANG,IQ,ITHP«ENTRAN
7 /CGNST2/IL.NIJ,IK.NI,NIP1
9 /CONST3/£Q,BO,CQ.CQ
IPB(NP.NS) ,PU(NP^NS»'pF(NP*NS/tPG0.
P01 I • Jt "0
PMX(I,J)*0.
PF-0.
PC(I,Ji«0.
PCOS(1,J)*0.
PUll ,JI-0.
P£T(I,J»«0.
DPMI 1 1 j) «n_


C-8
78

-------
FORTRAN IV G LEVEL 20.7 VS
                                        BLANK
                                                          CATE
                                                                    6/07/77
                                                                                14:32
  0040
 -OO41
  0042
  0043
  0044
  0045
  0046
-0047-
  0048
 .£049-
  0050
  0051
  0052
  £053-
  0054
  0055
  0056
 -0057-
  0058
                   PTU.JI-0.
                   PA(I,J)»0.
                   -PHCI,JJ*0.
                   PY(I,J)=0.
                   PANU,J)=0. .._
                   PDIL(IfJ)-C.
                   DO 2 1*1,NX
                 - £0 2 J=1,NY
                 2 CCNTA(I,JI=0.
                   DO 3 J=1,NCONT
                 3 CONTB(I)=0.
                 	DO S Ial,M&	
                 5 CONTC(I)*0.
                  -00-6 1«=U30
                   DC 6 J-1,30
                 6 AU(ItJ)=0.
                 	DO_7_l=i*2
  0060
  0061
  0062
  £063
  0064
_G£65.
  0066
  OC67-
  0068
  0069
  0070
--OOZ1-
  0072
  0073
  0074
  -OC75
  0076
_007_7_
  0078
 -OOJ^-
  0080
                   KE( I)=0
                   AT(!)=0.
                   Ah(I)=0.
                     •fi(j >=n.
                 7 AhG(IJ-0.
                   AC=0.000009153132
                   80=0.0002112502
                  -00=0.003660244
                   00=0.009494118
                   PAI-3.14159265
                  --GRA =9. 3066	
                   Al=0.0806
                 	A2=0.6753	
                   A3=0.3536
                   CD=1.5
                	ItRBFsO.
                   DX=0.5
 0082
                    X0=l.
                    20=?.
 0084
                    IPNT=0
                   _CfBL=A2/l O.116-AU.
0086
nnai
0088
000°
0090
nnol
0092
flC93
0094
-0095 	
IL>0
MTIsQ
IK = 1
	 NIjl>2
ITHP=1
	 IKO(1)S1 - - 	 	 ----
IS(U = 1
 0096
                   END
                                        C-9
                                         79

-------
FORTRAN IV G LEVEL 20.7 VS
                      FALV
DATE
6/07/77
 0001

 0002
 0003
 OOOt
 0005
 0006
-0007	
  FUNCTION FALV(TC)
 . JO DETERMlNE_LAT£M-h£AT. ALV	
  IF(TC .LT. 0.) GO TO 1
  FALV*(597.31-0.57*TC)**.1868
  GO TO 2
1 FALV-1677.01*0.622*TCJ*4.1868
2 RETURN
                                        C-10
                                          80

-------
FORTRAN  IV G LEVEL 20.7 VS
                      SETIC
                                                           DATE  -
                                                   6/07/77
 0001
 0003
_0004.
 0005
 0006
 0007
 .0008
 O009
.-D010-
 0011
-0012
  SUBROUTINE SETIC(Y,I,CX)
—COMMON -yr£aMSJl/J>AIJGRA,Al.A2,A3,A4,CD,TUftBF,UC,B,TP,ET,HP,£H,WP,
 1               MG,MGltCFRL,ALV,CPA,ANGtIO,ITHP,ENTftAN
 2       y~STaRE3/AH(30)tATi30).A£(30),At-G(301,ATGl30)tAU(30,301 ,
 *               AUG(3C,30)
 3       /STOR£5/DI(30),U013C),TO(30),HO(30),WO(30)
  OIKENSICN Y(8),CXI30)

  YI2)-0.                          "  "                          ~
        iii*nniT >
                       1I-WO(IJ»ALV/CPA)
  Y(6)=CX(I)
  RETURN
 -END  ,-
                                         C-ll
                                            81

-------
FORTRAN IV G LEVEL 20.7 VS
                            RUNGS
                                                          DATE -
                                                        6/07/77
1*532
 0001
 0002.
 0003
        SUBROUTI ME  RUNGS(X,H,H,Y,V PR I ME.INDEX. DERIV)
	    __DJMENSJOK JOBJ .YPRJME18) ,ZI8),hlI8)jW2(8)tW3C8).M4(8)
  C      RUNGS-RUNGE-KUTTA SOLUTION OF SET OF FIRST ORDER O.D.E. FORTRAN  I
  C      DIMENSIONS  *UST -BE SET FCR EACH PROGRAM
  C      X    INDEPENDENT VARIABLE
  C      H    INCREMENT DELTA  X.MAY BE CHANGED IN VALUE
  C      N    NUMBER  OF ECUATICNS
__£	Y  _P£P£NPFNT VARJARIF R1DCK	OKE DIMENSIONAL ARRAY        	
  C      YPRIME   DERIVATIVE  CLOCK    ONE DIMENSIONAL ARRAY
 £      THE  PROGRAMMER  KUST  SUPPLY INITIAL VALUES OF Yd) TO YCN)
  C      INDEX    IS  A  VARIABLE WHICH SHOULD BE SET TO ZERO BEFORE EACH
 _C      INITIAL  ENTRY TD THE SUBROUTINE. I.E.. TO SOLVE A DIFFERENT
  C      SET  OF EQUATIONS OR  TO START fclTH NEW INITIAL CONDITIONS.
_C   __LH£_PECGKAMK£R J4USI_WLJJ£ -A -SUBROUTINE CALLED .DERIVE WHICH     	
  c      COMPUTES THE  DERIVATIVES AND STORES Tt-EM
 _C      Th£  .ARGUMENT^1_LST -IS SUBROLTINE DERIVE(X.N.Y,YPRIME 1
        IF(INDEX) 5.5,1
   	1-DO 2-J=l.N	      	          -                    	
0005
OOP^
OOC7
0008
0009
0010
0011
001?
0013
0015
0016
0017
00 IS
0019
O020
0021
0022 — -
0023
n<)2^
0025
-0026
0027
O028
wl(I)=H»YPRIME(I)
y 71 TlsV I I I*lhi1 ( I I*.SI
A=X*0. 5*H
-CALL OERIV.UU*UZ,YPaiME)
DO 3 I=1,N
3 Z(I>=Y(I)*.5*W2(I)
CALL CERIV(A,N.Z, YPRIME)
-CO 4 1*1, N
V.3
-------
FORTRAN IV G LEVEL 20.7 VS
                                       OERIVR
                                                         DATE
                                                   6/07/77
                                                              15:27:
 0001
 0002
 0003
 000*

 0005
 0006
 0007
 0008
 0009
 0010
 0011
 0012
 0013
 0014
 0015
 0016
 0017
 0018
 0019
 0020
 0021
 0022
 0023
 0024
 0025

 0026
 0027
 0028
 0029
 0030
 0031
 0032
 0033
 0034
 0035
 0036
 0037
 0038
 0039
 0040
 0041
 0042
 0043
 0044
 0045
 0046
 0047
   SUBROUTINE OERIVR(S,N,Y,YP)
   DOUBLE PPECISICK *Q,BC,CC,CQ
   COMMON /STOREl/IKCt30),INDT<30»,NOVIM30»,IS<30l
  1       /STOftE3/AK30),AT(30),*Z(30},AHG(30),ATG<30),AJ(3C, 30J t
  A               ALG(30,30)
  3       /CONSTl/PMtGR*,Al,A2,A3,A4fCO,TURBF,UC,B,TP,ETt-IP,EH,*Pt
  4               KG, MGl tCFPL t ALV.CPA, ANG, 1C, ITHP.ENTRAN
  5       /CCNST3/ACtBQtCQtDQ
  6       /STORE5/CI(30I,UO(30),TO(30),HO(30),WO(30)
   DIMENSION Y(8),YP<8)
   DETERMINE AMBIENT CONDITIONS
   TPK«TP+273.16
   TOK»TO( ITHP)*273.16
   DO 88 1-2, MG
   IFIY<7) .GT. AZ(I)I GC TO  88
   11=1-1
   DZZ=Y(7)-«ZtII)
   TA=AT( II )*ATG(II )*CZZ
   HA=AH
-------
FORTRAN IV G LEVEL 20.7 VS
                        OERIVR
DATE
6/07/77
15:27:
 0048

 0049
 0050

 0051
 0052
 0053
 0054
 0055

 0056
 0057
 005S

 0059

 0060
 0061
 0062
 0063
 0064
 0065
 0066
 0067
 0068
 0069
 0070

 0071
 0072
 0073
 0074
 0075
 0076
 0077
 0078
 0079
    HUMIDITY AND TEMPERATURE
    CALL ITEP«TPK,hP,EST,Cl,C2l
    TO DETERMINE THE PLUME LIQUID PHASE MOISTURE
    WP«C3-HP
    IFCfcP .GT. O.I GC TO 79
    DRY PLUME
    MP«0.
    HP-C3
    TP»TA*Y4
    TPK-TP»273.16
    GO TO 178
    MET PLUME
 79 TP-TPK-273.16
178 ET-TP-TA
    TO DETERMINE ACIABATIC LAPSE RATE
    GAMA*FGAMA(TAK«t-At
    DETERMINE PLUME EKTRAINPENT
    RT-TPK/TAK
    PER-2.*PAI*B
    IFIPT .EC. 1.)  CC TO 99
    FRL*UC*UC*TPK/(GRA*TCK*ABS(RT-1.)»B)
    IFIFRL .GT. CFRL) GO TO 9
    A12*0.116
    GO TO 10
  9 A12«A1*A2*APSIN/FRL
    GO TO 10
 99 A12-A1
 10 ENTPAN«PER*IA12**BS«U)*»3*UA»APSIN*PCOS*A4*UP|
    EQUATICNS CF 'CCNSERVATICN OF VOLUME, MCM. ,ENERGE AND MOIST.  FLUXES
    YP(1»«ENTRAN*RT
    YP(2)-(UA*ENTRAN*PMC»/FSIN)*RT
    YP|3I-«RT-1.-HP|«GRA*LSL/UC-SIGN*FMC*PCOS*RT
    YP|4»«-fTG»GAMA»*9Y
    YPC5I— HG*SY
    YP(6)*PCOS
    YP(7)-PSIK
    RETURN
    END
                                        C-1A
                                          84

-------
FORTRAN IV G LEVEL 20.7 VS
                                       DErUVS
                                                         DATE -
                                                                   6/07/77
                                                               15:213
 0001
 0002
 0003
 0004
 0005

 0006
 0007
 0008
 0009
 0010
 0011
 0012
 0013
 0014
 0015
 0016
 0017
 0018
 0019
 0020
 0021
 0022
 0023
 0024
 0025

 0026
 0027
 0028
 0029
 0030
 0031
 0032
 0033
 0034
 0035
 0036
 0037
 0038
 0039
 0040
 0041
 0042
 0043

 0044
 0045
             C
             C
   SUBROUTINE DERIVS«StN,V,YP»
   DOUBLE PRECISION AQ.BCrCCtOO
   COMMON /STOREl/INCf30).INDT(30),NOVIM30),IS(30)
  1       /STORE3/AH(30l,AT(30«,«Z(3Cl,AHG(30l,ATG1301,AUJi0.30J,
  A               ALG(30,30)
  3       'CCNSTl/PAI,GR*,Al,A2,A3,A4,CD,TUReFtUC,B,TP,ET,rll>, BH.WF,
  *               PGtMGltCFRL,ALV,CPA,ANG,ICtITHP,ENTRA*
  5       /CCNST3/ACtBCtCQ,CQ
  6       /STORE5/DI(30>,UC(30«,TO{20»,HO(30»,WO(30)
   DIMENSION YIB).YP(B)
   TPK=TP*273.16
   DETERMINE AMBIENT CONDITIONS
   TOK«TO(ITHP)»273.16
   DO 68 I=2tMG
   IF(Y(7» .GT. AZ(IM GC TO 88
88
90
65
    YJ7)-AZ< Hi
TA=ATUI»**TG(II»*CZZ
HA=AHCII)*AHGim*DZZ
UA=AU(ITHPt II)*AUG(ITHP,II )*CZZ
TG=6TG(II»
HG*AhG(II)
GO TO 90
CONTINUE
DZl»Y«7)-AZCD*UA*UA*PSIN»PSIN*C.5
   Y4=Y(4)/LSU
   V5«Y(5I/LSU
   C1*ALV/CPA
   C3=Y5*HA
   TO ASSUME THE PLUME  IS SATURATED  AND CALCULATE THE SATURATED PLUPE
   HUMIDITY AND TEMPERATURE
   CALL ITER(TPKtKPtEST,CltC2)
   TO CETERMNE THE PLUME LIQUID  PHASE MOISTURE
   WP-C3-HP
   JFCKP -GT. O.t GO TO  79
                                       C-15
                                         85

-------
FORTRAN IV G LEVEL 20.7 VS
                        OERIVS
DATE
6/07/77
 0046
 00*7
 00*6
 00*9
 0050

 0051
 0052
 0053

 005*

 0055
 0056
 0057

 0058
 0059
 0060
 0061
 0062
 0063
 006*
 0065
 0066
    DRY PLUME
    MP-0.
    HP«C3
    TP-TAtV*
    TPK«TP*273.16
    GO TO  178
    WET PLUME
 79 TP-TPK-273.16
178 ET-TP-TA
    EH-HP-HA
    TO DETERMINE  ACIA6ATIC  LAPSE  RATE
    GAMA*FGAPAITAKtPM
    DETERMINE PLUME ENTRAPMENT
    TPK-TP+273.16
    RT-TPK/TAK
    ENTRAN»2.*J0.198»ABS(UMA3*UA*APSIN*PCOS»A**UP»
    EQUATIONS OF  CONSERVATION OF  VOLUKE,MOM..ENERGE  AND HOIS 1. FLUXES
    YP(1)*ENTRAN*RT
    YP(2»»-(TG«GAMA)*SV
    YP<5l«-HG»SV
    YP(6>«PCCS
    YP(7)-PS1N
    RETURN
    END
                                       C-16
                                          86

-------
FORTRAN IV G LEVEL 20.7 VS
                                       DERIVE
                                                         DATE
                                                   6/07/77
                                                                               15:27:
 0001
 0002
 0003
 000*
 0005

 0006
 0007
 0008
 0009
 0010
 0011
 0012
 0013
 0014
 0015
 0016
 0017
 0018
 0019
 0020
 0021
 0022
 0023
 0024
 0025

 0026
 0027
 0028
 0029
 0030
 0031
 0032
 0033
 0034
 0035
 0036
 0037
 0038
 0039

 0040
 0041
 0042
 0043
 0044
 0045
   SUBROUTINE DEPIVE(S,N,Y,YPI
   DOUBLE PRECISION AQ.BCtCOtCO
   COMMON /STORE1/INC130),INDTJ30I,NOVIM30I,ISC30I
  1       /STORE3/AH30«.AT(30»,AZ(30),AHG(30>,ATG(30),AU(30»30),
  A               AUC00.3C)
  2       /STCRE4/A(30),B1(30),B2(3C),BXZ(30),BY(30I,P*CJS( 30»,
  3               PrSIN<30)tDU(30),NeV(30),NEv(30l,NCV(JC»
  4       /CONSTl/PAI,GR*,Al.A2,A3,A4,CDfTUReFtUC,B,TP,ET,rt*.EH.WF,
  5               MGtMGl,CFRL,ALV,CPA,ANG,IC,ITHP,ENTRAN
  6       /CCKST3/AC,eQ,CO,CO
  7       /STORE5/OI<30ltUC<30)fTO<30l,HOI30),WC(30l
   DIMENSION Y(8I»YP<8»
   TPK-TP+273.16
   DETERMINE AMBIENT CONDITIONS
   TOK»TO(ITHP|*273.16
   DO 88 1*2,MG
   IF(V(7) .GT. AZdM GC TO 88
   11=1-1
   DZZ = Y<7|-AZUI)
   TA=ATUI MATGUI »*DZZ
   HA=AH(II)«AHG(II)*DZZ
   UA=AU(ITHP,II)+AUG(ITt-P,IIJ*CZZ
   TG=ATG(II»
   HG=AHG(II)
   GO TO 90
88 CONTINUE
   DZ1=Y(7)-*ZIMGII
   TA=AT(MG1)+ATG(KG1)*021
   HA=AH(MG1*»AHG(KG1»*CZ1
   UA=AU(ITHP,MG1)*«UG(ITHF.MG1I*DZ1
   TG=ATG(HG1)
   HG=AHG(MG1»
90 TAK»TA*273.16
   UP=CA*TURBF
   DETERMINE MOMENTUM,TRAJECTORY ANC VELOCITY OF PLUME
   PM=SORT(Y(2)*Y(2)«Y<3)**(3))
   PCOS»YI2)/FM
   PS1N»Y«3I/PM
   IFCPCOS .NE. O.I €0 TC 86
   ANG=90.
   GO TO 65
66 ANG=ATAMPSIN/PCCS)*180./PAI
85 SIGN=1.
   IFtPSIN .LT. 0.) SIGK»-1.
   APSIN=ABS(PSIN)
   SY=PSIN*Y(1»
   UC«PM/YC11
   U=UC-UA*PCCS
   USU=UC*(A«ITHP»*(ei(ITHP)*B2(ITHPII*0.5*PAI*
-------
FORTRAN IV C LEVEL 20.7 VS
                          DERIVE
DATE
6/07/77
15:27.
 0046

 0047
 0048

 0049
 0050
 0051
 0052
 0053

 0054
 0055
 0056

 0057

 0058
 0059
 0060
 0061
 0062
 0063
 0064
 0065
 0066
 0067
 0068
 0069
 0070
 0071
 0072
 0073
 0074
 0075
 0076

 0077
 0078

 0079
 0080
 0081
 0082
 0083
 0084
 0085
 0086
 0087
C     TO ASSUME THE PLUNE IS SATURATED AND CALCULATE THE SATURATED PLUPI
C     HUMIDITY AND TEMPERATURE
      CALL ITEDITPKfHP,EST,CltC2t
C     TO DETERMINE THE PLUME LIQUID PHASE MOISTURE
      WP-C3-HP
      IFChP .GT. 0.) GC TO 79
C     DRY PLUME
      WP'O.
      HP«C3
      TP-TA+Y4
      TPK-TP+273.16
      GO tO 178
C     MET PLUME
   79 TP-TPK-273.16
  178 ET-TP-TA
      EH=HP-HA
C     TO DETERMINE ADIABATIC LAPSE PATE
      GAMA*FGAPAITAK,I-A)
C     DETERMINE PLUME ENTRAINPENT
      TPK«TP+273.16
      RT-TPK/TAK
      PER=2.*A(ITHP»
      A12-0.196
      IFCPT .EG. l.J GO TO S9
      FRC-UC*UC»TPK/IUA*ENTRAN*PHC*«PSIN)*RT
      YP(3»-«RT-l.-WP»*GRA*tSL/UC-SIGN»PMC*PCOS*RT
      YPC4J —(TG*GAMA)*SY
      YPt5l— HG*SV
      YP(6)«PCOS
      YP«7)»PSIK
      RETURN
      END
                                        C-18
                                          88

-------
FORTRAN IV G LEVEL  20.7  VS
                                        PHI
                                                          DATE -
                                                                    6/07/77
                                                                                14:32
 0001
 -O002	
 0003
                   SUBROUTINE  PHI(I,JR,H,CY,PZ,NP,NSiIKQ,IKP)
                  —CCMMCN /5TflR£l/JWflL3.0JJJNDT(30J,NOVlK(30),IS(30)
                   1       /STORE4/A<30),61(30), 62(301,6X2(30),BY(30J.PMCOS130J,
                  -2.               J>MSINI30J,OM(30itNBV(30J.NEV(30) .NCVI30)
                   DIMENSION  PZINP,NSI,CY<30)
                   TO CALCULATE _THE  ANGLE OF THE INCLINATION OF  MERGED PLUME «PHI)
0004
_0fl£5_
0006
-0007 .
0008
.0009
0010
0011
0012
-0013 -
IJR*INDT(JR)
1P-1NDT1MJ
IFdJR .NE. IM) GO TO 1
JFIIJR -.££• J. _jAND _1M EO
IFCJR .EG. M) GC TO 6
Iff I IB -PC. ? AND IM FO
2 DY*ABS(CY(JRJ-CY(MI)
DS=SQRTIOY»DY*DZ*DZJ
	 	 FMCOSU)=BY/OS_ __

1) GO TO
2) GO TO
1

 -0015	
 0016
_OO17	
 0018
 -0019	
 0020
 .0021	
 0022
                   PMSIN(I)=DZ/DS
             _____ GO TO -5        __ . ___
                 1 IFIIJR .NE.  1) GC  TO  6
                 7 pjur.nst n»pw nsiMi
                   PMSIN( I)=PMSIN«M)
              ----- GC TO 5
                 e PMCOSC I)=PMCOS(JR)
                   GO TO 5
 0024
 0025-
 0026
 0027—
                   PMCOSi I)=SQRT(0.5*I1.*COSS)J
              ---- PMSINUJ
                 £ RETURN
                __ END
                                       C-19
                                         89

-------
FORTRAN IV G LEVEL 20.7 VS
                     ITER
DATE
                                                                   6/07/77
14:3.
 0001
.-0002	
 0003
 0004
 0005
 0006  ._
 OOC7
_DOOB	
 0009
 D010  	_
 0011
 D012
 0013
 SUBROUTINE ITEP tTPK.HP.EST ,C1,C2»
 T-l.-TS
 £ST-1013.25*EXPU13.3185-C1. 976* (0.6445*0. 1299*T)*T>*T)*T I
 ES«0.622*EST
 HP*ES/(1013.25*ES1
 FTPK=C2-C1*HP*273.16-TPK
.-ifiABS-iFipju— ^.L.n.r.n KFUJKN.
 DST«(13.3185-(3. 952* (1. 9335+0. 5196*T I*TI*TI*TS/TPK
 FTPKD*-1.*C1*>1P*OST*JHP-1.J
 TPK*TPK-FTPK/FTPKD
 GO TO 1
 END
                                          -G-20	
                                             90

-------
FORTRAN IV G LEVEL 20.7 VS
                                       FGAHA
                                                         DATE
                                                                   6/07/77
 0001

 0002
 0003
 0004
 0005
 0006

 0007
 0008
 0009

 0010
 .0011.
 0012
  FUNCTION FGAMA(TAK.HA)
— 10 J)£TERHJNE_AOJAaAUC J.APSE RATE
  T=1.-373.16/TAK
  EST-EXPi{ 13. 3185-11. 976* (0. 6445+0. 1229*T|*T )*T ) *T )
  £5-0.622*1013. 25*EST
  HAS=0.99*ES/ (1013.25+ES)
  IFIHA .GE. HAS) GO TC 1
  FGAMA=O.OOS76
  GO TO 2
1 RESP=EST/TAK
  SATURATED AOIABATIC LAPSE FATE
  FGAMA=O.OC976*(1.*5^20.*RESP)/(1.*839COOO.*RESP/TAK»
2-A£TURN ____  _____
  END
                                         C-21
                                           91

-------
FORTRAN IV G LEVEL 20.7 VS
                        SOLUTN
DATE -
6/07/77
14:32:
 0001

 0002



-0003 —
    SUBROUTINE  SOLUTNJPX,PY,PZ,PGtPMX,FMZ,PGtPFtPU,PENT,PStPA,PB,PCOSi
	lPSJN.PJ,P.ET..P.tUPEa.P.*fPAN.PDlL,PCtCYjYiYPtStl.JtNPtNSI
    CCHMON /STORE4/A(30)fBl(30)tB2(3Cl,BXZ(30)fBY(30ltPMCOS(30)i
   1               FMSINI30)(OMI30J,NBV<30),NEV(30)tNCV(30)
   2       /CONST1/PAIfGRA,Al,A2.A3tA4,COtTURBFfUCtBtTP.ETtHP.EH.MPt
   3               MG.MGl.CFRL.ALV.CPA.ANGfIQ.ITHP.ENTRAN
   4       /STOREI/IN0(30)11NOTI30)iNOVIK (30)tIS(30)
                                                      . PMXI
                  IPPZ(NPtNS),PG(NP,NS)fPF(NP,NS),PU(NP,NS),PENT(NPtNS),PS(NP,NS),
                  2PAINP.NSJ.PBJNP.NS).PCOS(NP.NS),PS1NlNP,NS)tPT(NP.NS),PET(NP.NS)f
                  3PH(NP,NS),PEH(NP,NS),PH(NP,NS),PAN(NPfNS),PD1L(NP,NS),CYC30),Yt8),
0004
0005
OOOb
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
PX(ItJ)*Y(6>
P7I TTnn = Vl 71
pY'
PUU tJ)*UC
PENTU,J4s£NTAAN 	 ..- 	
PS(I,J)-S
PA (I tJ)«A( I)
PC(ItJl=EXZ( I)
PPII ,.ll=R
PCOSd ,J)*YP(6I
J»SIN(I,JJ»YP47J^ 	
PT(ItJ)-TP
-PEHI,J)»ET 	
PH(I,J)*HP
ppm 1 1 ii =FW
PHII t J)«UP
PAN! It J)ZANC
PDIL(I,J)*Y(l>/Pti(Itl>
-RETURN . 	 	 -
END
                                         C-22
                                           92

-------
FORTRAN IV G  LEVEL  20.7 VS
                                        CHKNWP
                                                           DATE  -
                                                                     6/07/77
                                                                                 14:32:
 0001

 0002
 0003
    SUBROUTINE CHKNWP(PX,PZ,Pg,PMX,PMZ,PG,PF,PCOS,PSIN,PENT.PU.PS,PB,
	!PA,PT,P£T,PiL.P.£h,Pjl,F.Ak,PY,PClL,PC,CX,NP,NS,DERlVR,CY,KE)
    EXTERNAL DERIVR
    COMMON /STOP.E1/IND(30),INDT(30J,NOVIK(30),IS(30)
   1       /STORE4/A(301 ,81(30), 62(30), BXZ(30),BY(30).PMCOS(30),
   2               PMSIN(30),DW(30),NBVI30),NEVI30),NCV(30)
   3       /CONST1/PAI,GRA,A1,A2,A3,A4,CD,TURBF,UC,B,TP,ET,HP,EH,HP,
   4*               M/i u/~ i r* c ni  A i tt  ^ n *  * ±< r<  * **  • *.. n  ^ ».*«*..
 0004
            -C  -
   5        /CONST2/IL,NII,IK,NI,NIP1
   6        /STQR£5/DH30).UO(30),TO{30),HC(30I,WO(30)
   DIMENSION PX(NP,NS),PZ(NP,NS),PO(NP,NSI,PMX(NP,NSI,PMZ(NP,NS),
   1PB(NP,NS),PU(NP,NS),PF(NP,NSJ,PG(NP,NS),PCOS(NP,NS),PSIN(NP,NS)
   2PENT(NP,NS»,PA(NP,NS»,PT(NF,NS),CX(30),Ym,YP(8),PS(NP,NSJ,
  _3PHJNPjNSX,R»LLNP_,aSJj.E£LiJtf-,NSJjB£HiNPjNSJ ,PAN (NP , NS ) , PY (NP ,NS )1
   4PCIL(NP,NS),CY(301,K E(30),PC INP,NSI
   CHECK  JF  ANY NEk -PLUfE APPEARS
0005
0006 	
OOC7
onna
0009
0010 	
0011
D012
0013
0015
0016 __.
0017
0018
0019
r
0020
0021 	
0022
0023
0024
no?s
0026
Q0?7
IKK=IK*NOVIM1)-IS(1)
NIPP=NIP1
DP 77 -I = NIPP,NP
OS=DI(J)/20.
- 4F4ABS(CX4NIPP_)^CX(JJJ -GE. -CSJ -GO
40 FORMATdX, 'PLUME1, 2X,I3,'« APPEARS
N ! I =1 ._ j
IS( J)=IKK
-NOVIK(JJ=NOVJK(1J
MP1=J*1
IND«J)=1 	 . 	 -
IKS-1
TP-TO(J)
CALL SETIC(Y,J,CXI
S=0.
JTHP*-i
3C CALL RUNGS(S,DS,7,Y,YP,L,DERIVR)
.9 .
Jfl
AT
-


                                                           NSTEP='f I 5//)
             w     SCLUTICNS  FOR UNMERGED SINGLE ROUND PLUME
 0028	CAU._SQLUTN(P-X.PJL,Pi-.PQ,PMX,FM2,PG,PF, PU.PENT.PS, PA, PB, PCOS.PSIN,
                  !PT,PET,P»-,PEH,PWfPAN.POIL,PC,CY,V,YP.S,J,IKS,NP,NS)
                            jHf  TRTTPBlrN FHC STrPJ>-IMG-ThE_CALCULATJQN OF THE MEM	
0029
0030
O031
0032
4033
0034
0035
0036
0037
0038
C
	 	 	
	 23
22
	 19
ISSUED SINGLE PLUME
Jf-CP2IJ.JKS> »r.F- P*«l, IKKJJ^CD JD-23
IFdKS .GE. NS) GO TC 22
DS«B*0.1
_ALV=FAl V(TP) - - --
GO TO 30
KE(J)*1KS 	 	
CONTINUE
AETURN — -
END
                                         C-23
                                           93

-------
FORTRAN IV 6 LEVEL 20.7 VS
                      ALIGN
                                        DATE
6/07/77
14:32:
 0001              SUBROUTINE ALIGNJPX.PY.PZ,POfPKX.PMZ,PG.PF,PU,PENT,PStPA,PB.PC,
 	lPCOS,PSl/iL.P:i-.P£T*P.hUf£H.PHtPAN,PCJL.CY,YtYP,YS. YSPtYRl,VRlP,YR2,
                  2YR2PtDERIVE,DERIVS,DERIVR,S.LC,NF.NS.KE)
 0002              EXTERNAL OERIVEtOERIVS.DERIVR
 0003              COMMON /STORE1/INDI301,INDTI30I.NOVIK(301.ISC 301
                  J.       ycCNST_l/PAI,GRAtAl.A2fA3tA4tCD,TURBF,UC.B,TP,ET,HP,EH.WP,
                  2               MG,MGl,CFRLtALV,CPA,ANG,IO,ITHP,ENTRAN
 	3	aiC\ST2/ IL .NII .IX.NJ^AJPl _  	
 0004
4       /STORE4/A(30)fBl(301,B2(30)tBXZ(30>rBY(30),PMCQS(30),
5               PM51N(30),DH(30J.NBV(30JtNEV(30> ,NCV(30i
 DIMENSION  PX(NP,NS).PYINP.NS)tPZ(NP.NS).PQ(NP.NSItPMX(NP,NS),
aPM2(NPtNSJ«PG(NP,NSJ-tPFlNP,NS),PU
0\15
4)016
0017
-0018
O019
On?n
0021
_0022 	
0023
-0024
0025
p^MAx-pxu'ii!'1
DC 1 I=ltNI
1MzIK«MnilIKI I I-Itf 1 1
IFCPXMAX .LT. PX(I,If)) PXMAX=PX(I,IM)
	 -J -CCNTINbE 	 	
00 2 1=1, NI
JMINDUJ —LT-._H GO TO 2
IM=IK+NOVIK«I)-IS(I )
DS"=P8(I,IM)*0.1
JFUPJIMAX-PJU1.JMJ--OSJ ^.LT- O.J GO TO 2
IM=IM+1
IFUNOUJ^.GT. 41 -GO JO 8
IFIIM .GT. KEII1) GO TO 8
Rl 1 T)sPR( I TIM1
B2(I»=PBII tIM)
J>CU«IM)=PBU.1NJ 	 -
GO TO 9
R ITHP=[
 0027
  IKT=IM
  CALL A£SET_U_UJXlL,DS^S,NPtNS,CY.Y,YP,YS,YSP.YRl,YRlP,YR2.YR2P,   -
1DERIVE,OEPIVS,OERIVR,PX,PZ,PC,PMX,FMZ,PG,PF, PCOS,PSIN,PENT,PU,PBf
 0029
 CALL  SOLUTN(PX,PY,PZ,PJ,PMX,PMZ,PG,PF.PU,PENT,PS,PA»PB,PCOS,PSIN,
0030
0031
0032
0033
0034
9 NOVIK1 I)=NOVIK(I)*1
GC TO 11
2 CONTINUE
PETURN
END
                                          C-24
                                            94

-------
FORTRAN IV G LEVEL  20.7  VS
                      FLHERG
                                                          DATE
                                                  6/07/77
                                                                                14:32:
 0001
 0002
  SUBROUTINE  PLMERG«PX,PY, PZ,PT,PET,PH, PEH,PW,PA,PB,PAN,PDI L , PCOS ,
 lPSJN.CY,NP.,NS*Kfc.MGPA2JMGSTl.PC»LChK,PU)
  COMMON /STCRE1/INDI30) , I NOTC30) ,NOVIM30) ,15(30)
 1       /STORE2/MP(301,MS(3C1,IIP130»
 2       /STCRE4/A(30),BU30),e2<30l,BXZ«30l,BV(30),PMCOS»30),
 3               PHSIN130) ,OW(30) ,NBV t 30) tNEV ( 30) ,NCV(30)
 4       /CONST1/PA1,GRA,A1,A2,A3.A4,CO,TURBF,UC,B.TP,ET,HP,EH,HP,
 5 ________ M
 OO03
 6        /CCNST2/IL,NI I ,IK,NI,NIP1
 -DIMENSION PJUNP,NSJ,PZ(NP,NS),PB{NP.NS),PCOS7
0008
0009
0010
0011
0012
Opi i
0014
0015 	
0016
-0017 	
0018
001°
0020
0021
O022
0023
nnyL.
0025
0026
0027
0020
0029
OfMft
0031
^)032
0033
0034
0035
QlMfr
0037
0038
0039
0040
—0041
5LCHKC30) ,PU(NP,NS)
-C JIEARRAJUGE JHE-JTJIAJECTORIES OF -J HE EXISTING PLUMES
NJ*0
	 00 101-1=1, NI_ 	 	 __
IFUND(I)-1)101,102,102
10? I^r-= IK *Nnu TM I )-I S1I >
NJ=NJ*1
	 PXI(NJ)=PX
IRE(NJ)=I
101 CONT INUE
C CHECK PLUME MERGING
— IL=0
NI1 = NJ-1
nn i n3 i «] ,N| i
J 1*1*1
J3C 103 J=J1,KJ
PYD*ABSIPYI(I)-PYHJ)I
. P_iO=ABSlPiIiI J--=P-24U J J
PYZ=SORT
-------
FORTRAN IV G LEVEL 20.7  VS
                                        FLMERG
                                           DATE
6/07/77
                                                                                14:32
 0042
 -00.43	
 0044
 0045
 0046
 0047
-0048.
 0049
 0050
 0051
 0032
 0053
                    II-IRE(I)
               111
               113
               103

               205
     INDTtJJI-INDCJJJ
     JNDTdlJ'INOdl)
     RESET  IND(III*2.INDfJJ)«0 FOR JJ>II MHEN PLUMES II AND JJ MERGED
     JND(JJ)*0
     IIPHDI t 1
LS*MSd)
HT-IIP1LP1
NT>IIP(LS)
MGSTHLS)-MT
IFIKPIICI .r.T _ MTI r.p JO 116
LCHMLS)»1
116 WRITE(6.110)MT,LP.PX(LP,MTJ,PZ(LP»MT),PYILP,MT) ,PT(
1«>ET(LP,^T) (PHiLP.MT) fPFHItP.MT) tPyUP,MT)tPp((,P|MT>
2PAILP,MT),PAN(LPfMT).PDILILPtMT),PU(LP.MT)
i in FORMAT (ix.?it<:riyi T^IF«-?T I * I r ?'F7-7. lX)»3IF6i-$r 1^1
117 INDX*2
LP-LS
MT»NT . . .
GO TO 116


LP.MTI,
.PC(LP.MT),
.6IF7.2.1XJJ
                                          C-26
                                             96

-------
urv • n>H
0091
D092 .
0093
0094
0095
0096
OCS7
009S .
0099
0100

0101
0102
0103
0104
0105
0106
n i» o ICVCL iu«/ va PLHERG
119 FORMATC/I "
- - 	 ICS-COWTJNUE 	
00 777 1=2, NI
IF(INO(I)-1J 777,779,779
777 CONTINUE
- - MI=0
ITHP«1
779 DO 400 _JsL?rMP
IF(INO(I))16,400,16
-. 	 -40C CONTINUE
C ALL THE PLUMES ARE MERGED KHEN INOJ
	 _j ^Q ( i J *3
10=3
N 1=1
ITHP=l
	 16 -RETURN 	 _^ .
ENC
DATE * 6/07/77








~~

II=C FOR KKNPfrl






                                        14:32:
C-27
  97

-------
FORTRAN IV G LEVEL 20.7 VS
                                       RESETI
                                                          DATE
                                                                    6/07/77
19:42:
 0001
                   SUBROUTINE RESfTI 4 I, IKT ,OS.S,NP,NS,CV, Y, YP, VS, VSP, VR1. VitiP ,VR2,
                  1VR2P,DERIVE,OEJUVStOEPlVR,PX,PZ.PQtPM)
 0002
 0003
                  2PU,PB,PA,PT,PET,PEH,Ph,LC)
                  -EXTERNAL -OE«IV4^,4J€«IVS»0€«IVR  --
                   COMMON /STCREl/INCt30),lNDT430l,NOVIK430)tlS(30)
                          /STCRE4/A430»,BIC30),B2I30I,BXZ<30I,BYI30),PMCJSI 301,
                  4
                 _5
                          /CONSTl/PAI,GRAtAl,A2,43,A4,CD,TUReF,UC,B,TP,ET,f,IRP(30).COdl33J,Ylft),
 0005
 OO06
 0007
-0008	
 0009
 -0010	
 0011
 O012	
 0013
-001-4	
 0015

 0016
 O017
 0018
-OOI*	
 0020
 O021  ——
 0022
 0023	
 0024
-0025	
 0026
 0027	
 0028
 -0029	
 0030
-0031	
 0032
 -0033	
 0034
 0035	
 0036
                  3PA(NP,NS>tPCaS(NP,NSI,PSIN»PPX
-------
FORTRAN IV G LEVEL 20.7 VS
                                       RESETI
                                                         DATE
                                                 6/07/77
                                                                               19:42:
0042
-0043
0045
0046
0047
0048
-0049
0050
0051
0052
0053
0054
-0055
0056
0057
0058
O059
0060
-0061
0062
0063
0064
0065
0066
0068
O069
00/0
0071
0072
0074
0075
4)076
OOJ7
-9O7-8
0079
-0080
0081
-0082
0083
Ann*.
0085
fin at.
Y(5)*YI5)»PF(JJ,IKPI
APT = APT»PTUJ,IKP)
IJK«!JK+1 -
IRP< 1JK|=JJ
CYMANC=CY-fiY(JJ»
	 D 7UAkl<- - 07/11 luni.nw......
	 -400
404
401
244
	 -28
26
- — 	 34
	 	 	 gf ft
en
859
C
575
572
	 	

PZMINC=PZ(JJ,IKP)-BXZJJJ)
IFUNOTANC)-6C TO 26
JRZ-JJ
IF(PZMIN .LE. PZKINC) 60 TC 34
P7MT NT P 7 MI MT
CONTINUE -
IF(MM-IJK) 58,59,59
GO TO 60
APT=APT/FLCAT(IJK)
LC=C M?ANS CLUSTER OF TCWERS; LC-1 MEANS LIKE TOWER AitKA Y
CY( II=0.5*«CYMAK*CYMIN)
IKO«IK*NCVIK(M)-IS(M|
DCYAI=CYKAK-CYKIN
HJPZAI-PZMAN-PZ^ll. 	
IFiCCYAI .LT. CPZAI) GC TO 402
i fit rVMAU^ »•«"«"• «T f P 7HAN— rf H INI 1 CO TD 40?
IF(CY(M) .EO. CY(JR)) GC TO 407
 0088
TO  CALCULATE  THE  ANGLE  CF  THE  INCLINATION OF MERGED PLUXE  «PHI)

-------
FORTRAN IV G LEVEL 20.7 VS
                                       RESETI
                                                         DATE
6/07/77
0096
0097 	 	
0098
0099
0100
0101
0102
Al A**
0104
0105
0106
0107
fti fto
0109
0110
0111
0112
0113
/ii n £
0115
0116
0117
0118
0119
•0120
0/21
0122
0123
0124 •--
0125
fil 3t*
0127
O128
0129
0130
0131
ni T"
0133
0135
ni ift
0137
Q13B
GO TO 573
GO TO 573
402 JR=JRZ
M=MZ
1KQ=IK+NOVIKIMI-IS«M>
IKP=IK*NCVIK(JR|-IS( JP|
C
PZY=ABS(CYIJR)-CY(MI 1
TO CALCULATE THE -ANGLE -CF THE INCLINATION OF MERGED PLU-4E (PHI)
CALL PHI 
                   -VB1U*-BBB	
 0144
-0145
 0146
 O147
 0148
                   YRlf2I*BBB*PU(M,IKQ)*PCCS-PZ
-------
FORTRAN  IV  G  LEVEL  20.7 VS             RESETI

              C      RESET I.C. FOR THE HALF ROUND PLUME
                                                          DATE
                                                                    6/07/77
                                                                     19:42:
-0150-
 0151
 0152
 0153
 0154
 0155
-015*-
 0157
 0158
 0159
 0160
-04*1-
 0162
 €163
 0164

 0165
-01**-
 0167
-01*8
 0169
-0170
 0171
 0172
 O173	
 0174
 -0175 —-
 0176
-0177	
 0178
 0179
 0180
 0181
-01-82-
 0183
        VP2m5»-0.5*(VPl(;*»VR2l5Hi/A(I»
        YS(7|*Y<7>
        CALCULATE -NEW -HAL^-W!CT«- AND -VELO, -OF THE -HALF ftCUNO
        L»0
        AL-V-AtVl ---   ---------------
        CALL RUNCS -AKO ^ELO. OF
        L = 0
     -- ALV-FALVUPT^        ----
        CALL RUNGS(S,DS,5tYS,VSF,ltOERIVSJ
                    -rO^F S t VS , VST , L
                                                              SLCT JET
        BS=P
 0185
    	C-
             c
-01-86	
 0187
 0188	
        S=S-DS
       — CALCULATE -MERGEC -PLUME
        THE PREVIOUS STEP
                                               ^THE -KOOIFI EO PLUME SHAPE dETe^MINED C
                    CALL  PUNGS(S.DS,7,Y,YP,LiOERIVEI
 0189
-0190
 0191
 0192
 0193
 0194
 0195
-6196-
        DETERMINE 81,82 AND A CF THE MODIFIED PLUME SHAPE FCK OLiULATINC
       -EWTfAINKENT^NO-CRAG-FOP-NEXT
        UU-UC
        C61=BR2/eRl

           *U 5*07963* (eRl*BRl*URl»BR2*BR2*UR2»»2.*BS*M*US»/JU
      ~APllil5707963*«U*CBl*CEl)-CB4*CB4
                                         C-31
                                          101

-------
FORTRAN IV G LEVEL 20.7 VS

 0197              CP1=-CB3
RESETI
DATE
6/07/77
19X42:
0199
0200
0201
0202
.n?n^
0204
0205
0206
0207
0208
Ji?AQ
0210
0211
0212
0213
0214
0215
0216
0217
0218
0219
u??n
0221
0222
0223
02.24
0225
-A?? A
0227
O228
0229
0230
flyn
0232
0234
./loac
0236
n ? T 7 .
0238
— n?ao —
62P = B2(U
-AP=MII
C THE DETERMINED B1,B2 *NC A
-BlU»=<-BPl*SQRTALC
fi V 7 < IkmAAlAVl 1 A I i I
GO TO 416
BY(I)=AMAX1(B1(I)
GO TO 31
Cr AI rni ATC r T ur i c
29 L-0
AL V* FALV ( PT ( I IK1
I| 1

,B2(I))
ROUND -HAH> E - 	 ... _ .
» •
CALL RUNGS(StDS,1,Y,YF,LtOERIVR)
B1(II=B
BXZtll-B
31 RETURN
	 ciun



  C-32
   102

-------
FORTRAN IV G  LEVEL  20.7 VS
                                        OUTPUT
                                                          DATE
                                                                    6/07/77
 0001
 O002
  19:42:


IX. I
-0003-
                   COMMON VSTORE1/INC<30»,INDT<30»-,4IOVIK<30»,IS<30*
                   1        'CCNTUU/WICTH,HITE,MOPE,NCMAP,ICENT,KOTICK,HEONUOJ•
                   -2                LABY<5»,IABX<5>
                   3        /STORE5/DI<30»,UOC30l,TO(30l,HOC30»tWOI30)
 0005
-0006-
 0007
 0008
 0009
-001O-
 0011
                   lPmNP,NSl,PEH,PY-
                42 IFCFEWC .GT. O.I GO TC 47
               	tX
-------
CORTRAN
DATE «
6/07/77
19:42:
0041
nnt?
0043
0044
0045
0046
0047
nn&n
0049
0050
0051
0052
0053
nn*i±
0055
0056
0057
0058
0059
flf\4*f\
0061
OO6Z '
0063
006.4 -
0065
ftnx-A
0067
UUOO
0069


11
31
a
	 32
— •-•*
a?
37
JU
7
	 5
1ft
17
--15

PZ(ItJ»«PZU,Jl/OI(U
PC(I,J)=PC(I,J)/OI(1)
PSII,J)'PSU,J)/OI<1)
PW(I,J)-PMI,J)/PEWO
PET«I,J)=PETII,J»/PETC
IF( IX(I» .LE. 0» GO TC 5
CALL ^P02(CONTA,OXYOZ,PX,PZ,PC,PET,PCOS,PSIN,NP,NS,KE.XJ ,ZO,NX.N'Vi
WRITE(6,31)
FORMATi IX, 'EXCESS TEMPERATURE PLOTS//*
GO TO 37
WRITE16.22*
-FOR^ATIlXi^-exCESS'-SPECIHC-ftUMIOITV PLOT',//)
GO TO 37
HRITE(t,33)
RE API 5, 30 1 (hEOMK),K*l,10)t(LABY(L) ,L»1 ,5 ) , CLABXC M),M«1 ,3)
CALL CCNTU(CCNTA,CCNTE,CONTC,NX,KY,NCONB,IFM)
IF(I-10)1C,17,15
IF( IX( I))5,5,9
END

                                        C-34
                                         104

-------
FORTRAN  IV  G  LEVEL 20.7 VS
                                        FMAX
                                                          DATE -
                                                   6/07/77
                                                                                 14:32
 0001
 .4002 —
 0003
 0004
 0005
 0006
 0007
-flOQB _
 0009
 0010
   FUNCTION PKAX«PT,NP,f*S,KE)
	DIMENSION PUNP...NSL.KEOO)
   PMAX«PTU,1)
   DO 1 I*1,NP
   K=KE«I)
   DC 1 J=1,K
   IF(PT(I,J) .GT. PMAX) PMAX-PT(I,JJ
_l_CONTOJUU£	._  	
   RETURN
   -END
                                       C-35
                                        105

-------
FORTRAN IV G LEVEL 20.7 VS
                            GRD2
DATE
6/07/77
14:32-
 0001

 0002
 0003

 0004
 0005
.0006
 0007
 0008
 0009
 0010
 0011
_0012_
 0013
 0014
 0015
 0016
 0017
-0018-
 0019
 0020
 0021
 0022
 0023
 .0024.
 0025
 0026
 0027
 -0028
 0029
-0030-
 0021
 -0032
 0033
 0034
 0035
-0036.
 0037
 0038
 0039
 -0040-
 0041
-004-2-
 0043
 0044
 0045
 O046
 0047
 0048.
 0049
 4050
 0051
 0052
 0053
 -0054.
        SUBROUTINE GRD2JA,DX,OZ,PX,PZ.PB,PG,PCOS,PS1N,NP.NS.KE,XO,ZO,10,
    	1JD)	     _....
        COMMON /STORE1/INO(3C).INOT(3C) tNOVIKOOl tISI30l
        DIMENSION JUlD,JD)tRM(5JiC(5i,lC(4J,PX(NPfNS)tPZINPtNS)fPG(NPtNS)
       1PB/nX» 1.0005 ______
        PB2C-PB{N,K*U*PCOS(K,K*1)
        PB1S=PBIN,K)*PSIN(N,K)
        DZZ-PZ«N,K*14-P2tN,X)
        DXX=PX(N,K+1»-PX(N,K)
        DPBS-PB2S-PB1S
        OS«SO.RHOZZ*XZ
        CB=PB(N,K»1)-PB(N,KI
        OG=PGJN,KJ-PmN^K*J,J ____
        IF(PSIN(N,K*1) .LT. C.I GO TO 101
    105 IL*IC(3)+1
	GO TOJ.12
    104 IL=ICI3>+1
	GO JO -113
    103 lL
        C(3»-PZtN,K+l)*PB2C-(PX(N,K+l)-PB2S)*RM<3)
        IFiDXPB .EQ. 0.) GO TO 20
        C(4)«PZ(N,K»1»-PB2C-
-------
 FORTRAN IV G LEVEL 2C.7 VS
                                        GR02
                                                          GATE
                                                                    6/07/77
  0055
 .0056
  0057
  0058
  OC59
  0060
  0061
 -00£2-
  0063
  O064
  0065
  O066
  0067
 -OO43.
  0069
  O070
  0071
  O072

-0013-
  0074
  OC75
  0076
 -0077-
  0078
-0079-
  0080
  0081
  0082
  OC83
  0084
      1FCI  ILL  1)  GO TO  12
    _X=OX*f LCATL1-1)	
      IFII  .GT.  ICI3)>  GO TO 7
      GO  TO  8
    7 -ZU=RMI2)*X»CI2)
    8  JU=ZU/DZ*1.0001
 	IFU—GT. -1CL2JJ_GO_LO-10	
      IFII  .GT.  Kill)  GO TO 115
      2L-RMI3)*X»C13)
      GO  TO  11
 115  ZL=RM11J*X4C11J
      GO  TO  11
 	1C _ZL=AW4) «J(*£X4)
   11  IFUZL-OZ*FLOATIJU-1)J  .GT. O.J  GO TO 12
      1FIJU  .LT.  JL)  GO  TO 12
     10XX,DZZ,I,K, ID.JDI
 0086
 0087
 0088
 O089
 0090
 nnci
      GO TO 99
  101  1FUC<2)-IC<4J)  121,122,123
  122  IL=IC(2)*1
      GC TO 124
  122  IL=IC(4I+1
    _GC JDJ.25 _________
  121  IL=1CI2I+1
  125  KMJ4)=(DZZ-OPBCJV10XX*OPBSJ
      C«4)=PZ(N,K+1I-PB2C-(PX(N,K+11*PB2S)*RMC4J
  12«  OXOP«OXX-OP8S-    —
      IFIDXOP  .EC.  0.) GO  TO 21
      RJU3 J *1 DZ Z *JP-BC jy-DXDK- __ ______  .
      C ( 3 ) =P Z I N t K* 1 » * PB2C -(PX(N,K+1)-PB2S)*P.M(3)
  21  RM(2J=-PCOSiN,K*l)/J>SlN{N,K*H
      C(2)=PZUtK+ll-PX(N,K+l)*RM(2)
  . .— JR-ICI3I ---  --------
   •   IFUL ,GT.  IR) GO  TO  99
      DO _L2fl_LaiL»Jfi ________  ______
 0092
 -0093
 CO 9 4
 4)095
 0096
-CC5J-
 0098
 0099
 0100
 -0101.
 0102
-0103.
 0104
      IFII  .GT.  IDJ  GC  TO  128
    . JFI1—U.. -1J-GO-TO 1?H
      X=DX*FLOAT(I-1)
    -JFtl  ,GT.  1C1J.JJJJO  10 130
      IFCI  .GT.  IC12H  GC  TO 132
     GO TO  133
-J.32 -ZU=RMI1)*J
     GO TO  133
 130-21
 132 JO-ZU/DZ4-1.0001
 0106
 0107-
 01C8
...0109 .
 	GC TO  136	
 135 ZL=RMI2)*X*CI2)
 136 IF(IZL-OZ*FLOATUU-1I) .GT. 0.) GO TO  128
     JL=ZL/DZ*1.9999
                                          C-37
                                           107

-------
FORTRAN IV G LEVEL 20.7 VS
                      GR02
                                                          CATE  -
                                                  6/07/77
14:32:
 0110
 CALL PC2*in
CH=0.
172 MRITEC6.178J
176 FORMATI1H1I
DO 662 J=1.JO
L=JD*1-J
bt2 bRJT£16.3) < Al I , LJ*J *J P» J" »
3 FORMAT(lXt20(F5.2tlX)l
IF1CH -.EC- 0«J GO TO 173
IP-IP+20
GC TO 171
17' DO 310 N»l t^P . . ,
CO 310 I«l tNS
PX(NtI)=PX
-------
FORTRAN IV G  LEVEL  20.7  VS
                                        FD2
                                                          GATE
                                                                    6/07/77
                                                                                14:32
 0001

 0002

 0003
 .0004.
 0005
-OOO6-
 0007
 0008
 0009
 0010
 0011
 nni 7
                    SUBROUTINE  PD2(A,PX,PZ,PB,PG,PCOS,PSIN,RM,C,NP, NS.N,JL,JU,08,OG,
                  _1£Z-,JL. DXX ,U£Z-, J^K, la, -JO J
                    DIMENSION AUD,JDI,PX«NP,NS»,FZCNP,NS),PB»NP,NS),PG
-------
FORTRAN IV G LEVEL 20.7 VS
                       CONTL
DATE
8/03/77
10:1*:
 0001
 OOC2
 0003

 0004
 0005
 0006
 0007
 0008
 0009
 0010
 0011
 0012
 0013
 0014
 0015
 0016
 0017
 0018
 0019
 0020
 0021
 0022
 0023
 0024
 0025
 0026
 0027
 0028
 0029
 0030
 0031
 0032
 0033
 0034
 0035
 0036
 0037
 0038
 0039
 0040
 0041
 0042
 0043
 0044
 0045
 0046
 0047
 0048
 0049
 0050
 0051
 0052
 0053
 0054
 0055
   SUBROUTINE CONTUISET
   YMAX«LIM
   YNAX-YKAX
   YNIN-YKIN
   XNAX~XNAX
   XNIN-XNIN
   IF(KEV(2I.EQ.1I GOTO 10
   DX*XMAX-XMIN
   DY-YMAX-YMIN
   IFIDX.GT.OYI GOTO 11
   DUM'XSIZE
   XSIZE-YSIZE
   YSIZE-DUM
11 YSIZE«ANIN1IYSIZE.26.)
   OX-AHIM C XSI ZE/OXtYSI ZE/ CY I
   XSIZE-DX*XMIN)
                                         C-40
                                          110

-------
FORTRAN IV G LEVEL 20.7 VS
                                       CONTU
                                                         DATE
                                                     8/03/77
10:14
 0056
 OC57
 0058
 0059
 006C
 0061
 0062
 0063
 0064
 0065
 0066
 0067
 0068
 0069
 0070
 0071
 0072
 0073
 0074
 0075
 0076
 0077
 0078
 0079
 0080
 0081
 0082
 0083
 0084
 0085
 0086
 0087
 0088
 0089
 0090
 0091
 0092
 0093
 0094
 0095
 0096
 0097
 0098
 0099
 0100
 0101
 0102
 0103
 0104
 0105
 0106
 0107
 0108
 0109
 0110
 0111
 0112
 0113
 0114
 0115
     VD-YSIZE/X
-------
FORTRAN IV G LEVEL 20.7 VS
                          TOPO
                                      DATE
8/03/77
10X14S
 0056
 0057
 0058
 0059
 0060
 0061
 0062
 0063
 0064
 0065
 0066
 0067
 0068
 0069
 0070
 0071
 0072
 0073
 0074
 0075
 0076
 0077
 0078
 0079
 0080
 0081
 0082
 0083
 0084
 0085
 0086
 0087
 0088
 0089
 0090
 0091
 0092
 00 S3
 0094
 0095
 0096
 0097
 0098
 0099
 0100
 0101
 0102
 0103
 0104
 0105
 0106
 0107
 0108
 0109
 0110
 0111
RS-RBIJ)
XS=X (J)
YS-Y CM
RM-RR
XM»XX
YM-YY
 IF 
-------
FORTRAN IV G LEVEL 20.7 VS
                          TOPO
DATE
8/03/77
 0112
 0113
 0115
 0116
 0117
 0118
 0119
 0120
 0121
 0122
 0123
 0124
 0125
 0126
 0127
 0128
 0129
 0130
 0131
 0132
 0133
10116 IFUBS (XPA-XPB)-.00115003,5004,5004
 5003 IFUBS (YPA-YPB)-.OOIUOC,5004, 5004
 5004 CALL   PLOTZUPAtVPA.XPBfYPB)
  100 RC-CNTRINC+1)
      NC»NC*1
      60 TQ 80
  103 XPA • XCL
      YPA « YCL
      GO TO 99
  106 0-(RC-RM)/tRL-RM»
      XPA«XCH-0*»XCH-XCLI
      YPA»YCH-QMYCH-YCL»
      GO TO 99
  110  GO TO L,(112,118)
  112  ASSIGN 118 TO L
      RR   -RB(J-l)
      XX   -X (J-l)
      YY   »Y (K)
      GO TO 37
  118 CONTINUE
       RETURN
       END
                                            C-45
                                             115

-------
FORTRAN IV 6 LEVEL 20.7 VS             PLOT2              DATE -     8/03/77      1

 0001              SUBROUTINE PLOTZ1X1.Vl*X2fY2I
 0002              INTEGER EC
 0003              LOGICAL*! LXS
 0004              DIMENSION X(2>,Y(2),LXS<1)
 0005              COMMON/CONTUU/SKIUPI26),
                  *              XMINtXMAXt VMINtYMAXtXOFFtYOFFtXSIZEtYSIZEtNCOMT
 0006              XUI-XSIZE *X1+XOFF
 OOC7              X(2)*XSIZE *X2+XOFF
 0008              Y(1)»YSIZE *Y1*YOFF
 0009              Y(2)-YSIZE *V2*YOFF
 0010              CALL SYSPLT(X(UtVll),3)
 0011              CALL SYSf>LT(X(2),Y(2lt2l
 0012              RETURN
 0013              ENTRY FINDMTCLXS.HT,NBTI
 0014              MT-NBT
 0015              00 44 J-l.NBT
 0016              IF(EC(LXS(MTI.* *).EQ.O) RETURN
 0017           44 MT-HT-1
 0018              RETURN
 0019               END
                                          C-46
                                           116

-------
FDYN
         SUPPORT FOP FORTRAN DYNAMIC BIT AND  FULL WORD  ARRAYS.
                                PAG;
  LOG  OBJECT CODE  , ADOR1 ADDR2  STKT   SOURCE  STATEMENT
OCOOOO

OCOOOO 47FO FOOC
000004 06
000005 C7C5E3ESC5C3
OOOOOB 00
COOOOC 90E4 DOOC
000010 184F

000012 1821
000014 5830 2000

 000018 5800 3004
 00001C 8800 0002
 000020 5000 3000
000024 4510 4028
000028 OAOA
OC002A 5010 3004
00002E 4130 3000
000032 1B13
000034 8A10 0002
00003B 4110 1001
00028
00004
00000
00002
00001
 00003C 5010 3008
 000040 98E4 OOOC
 000044 92FF OOOC
 000048 41FO 0000
 00004C 07FE

 00004E 47FO FOOC
 000052 06
 OC0053 C609C5C5E4D7
 000059 00
 00005A 90E4 DOOC
 OC005E 184F

 000060 1821
 000062 5830 2000

 OC0066 9801 3000
  OC006A 4111 0000
  00006E OAOA
OOOOF
00000
00001
00002
00003
00004



OOOOC


ooooc

00000

00000
00000
00004
00002
00000

00028

00004
00000

00002
00001
ooooa


ooooc
ooooc
00000


ooooc


ooooc

0004E

00000
00000
00000


00000


2 X15
3 XO
4 XI
5 X2
6 X3
7 X4
8
9 GETVEC
10
11*
12*
13*
14*
15
16
17
IB
19
20
21
22
23
24*
25*
26
27
28
29
30
31
32
33
34*
35*
36*
37*
38 FREEUP
39*FREEUP
40*
41*
42*
43
44
45
46
47
48
49
50
51*
52*
53
EOU IS
EOU 0
EOU 1
EQU 2
EQU 3
EQU 4
ENTRY FREEUP
CSECT
SAVE (14, 4),, «
6 12IOil5>
DC AL1I6)
DC CL61 GETVEC1
STM 14,4,12(131
LR X4.X15
USING GETVEC, X4
LR X2tXl
L X3,0(,X2)
USING ANCHOR, X3
L XOtOIMEN
SLA X0,2
ST XO, LENGTH
GETMAIN R.LV-COI
BAL 1,**4
SVC 10
ST XltLENGTH*4
LA X3,0(.X3)
SR XI, X3
SRA XI. 2
LA Xl.l(tXl)
ST XI. INDEX
DROP X3
RETURN (14.4I.T.RC-0
LM 14.4.12(13)
HVI 12(13), X'FF1
LA 15.0(0,0)
BR 14
SAVE (14,4),,*
B 12(0.15)
DC ALK6)
CC CL6' FREEUP*
STH 14,4,12(131
LR X4.X15
USING FREEUP, X4
LR X2,X1
L X3.0(.X2l
USING ANCHOR. X3
LM XO, XI, LENGTH
CROP X3
FREEMAIN R.LV-IO) ,A-( 1 )
LA 1,0(11
SVC 10
RETURN I14.4),T,RC-0














•X4


•X3












• X3











• X4


*X3

«X3




                  ASM 0200 10.13 08/03/77
                                                                                                                       11
                                                                                                                       12
                                                                                                                       13
                                                                                                                       14
                                                                                                                       15
                                                                                                                       16
                                                                                                                       23
FREE A DYNAMIC ARRAY.
SAVE REGISTERS.
          BRANCH AROUND ID

          IDENTIFIER

          SAVE REGISTERS
SET UP A BASE REGISTER.

1X21- AIPARANETER LIST I.
1X31- AIARRAY ANCHOR!.

1X0)- DIMENSION OF ARRAY.
(XO)- LENGTH OF ARRAY IN BYTES.

(XII- A(ARRAV).
          INDICATE GETMAIN
          ISSUE GETMAIN SVC
1X31- AIANCHORIDI.
(XII- AIARRAYIDI - A( ANCHOR (1).
CONVEPT TO ARRAY INDEX MI TH
RESPECT TO THE ANCHOR I
ARRAY III - ANCHORIINDEX).
PLANT INDEX IN ARRAY ANCHOR.
RETURN TO CALLER.
          RESTORE THE REGISTERS
          SET RETURN INDICATION
          LOAD RETURN CODE
          RETURN
SAVE REGISTERS.
          BRANCH AROUND ID

          IDENTIFIER
          SAVE REGISTERS
SET UP A BASE REGISTER.

(X2I- AIPARANETER LISTI.
(X3>> MARRAV ANCHOR*.

(XO,X1)» FREEMAIN PARAMETERS.

RELEASE THE CORE.
     CLEAR THE HIGH ORDER BYTE
ISSUE FREEMAIN SVC
RETURN TO CALLER.
 00860000
 00880000
 00900000

 01180000
       55
       56
       57

       59
       60
       61

       62
 68800001
 69000001
       6 i
*      68
       69
       70
       71

 00260000
 00640000
 00700000
 00800000
      124
 00860000
 00880000
 00900000

 01180000
      125
      126
      127

      129
      130
      131
      132
 03130018
 03140018

-------
            FDVN
                     SUPPORT FOR FORTRAN DYNAMIC  BIT AND FULLKORO ARRAYS.
                                                                                          PAGE
              LOG  OBJECT CODE

            000070 9BE4 DOOC
            000074 92FF OOOC
            OC007B 41FO 0000
            00007C 07FE
            000000
            CCOOOO
            000004
            000008
            CCOOOO
            00000B
AODR1 ADDR2  STNT   SOURCE STATEMENT
OOOOC
000 OC
00000
                                                                                                              ASM 0200 10.13 08/03/77
      00000
                                        00008
54»
55»
56 »
57»
58 ANCHOR
59
60 DIMEN
61
62 LENCTt
63 I NCR
64 INDEX
69
LN
XVI
LA
BR
DSECT
CS
OS
ORG
CS
OS
EQU
END
14,4,121131
12(13), X*
15,010,01
14

F
F
ANCHOR
Zf
F
1NCR

FF«










          RESTCRE THE REGISTERS
          SET RETURN INDICATION
          LOAD RETURN CODE
          RETURN

UNUSED MORO.
BITS/ROM OR • OF FULLWORDS.

LENGTH OF THE ARRAY, BYTES.
BYTES BETWEEN ROMS IN A BIT
ARRAY.
INDEX FROM START OF ANCHOR TO
START OF ARRAY FOR FULLHORD
ARRAYS.
 00260000
 00640000
 00700000
 OOBOOOOO
      465
      469
      470
      479
      476
*     480
      481
*     482
*     483
      484
      486
M?
H« *•
00 00

-------
                                   TECHNICAL REPORT DATA
                            (Please read Instructions on the reverse before completing)
1. REPORT NO.

   EPA-600/7-78-102
2.
                              3. RECIPIENT'S ACCESSION NO.
4. TITLE AND SUBTITLE
 "MATHEMATICAL MODEL FOR  MULTIPLE COOLING TOWER  PLUMES"
                                                              REPORT DATE

                                                              June 1978
                                                            6. PERFORMING ORGANIZATION CODE
 . AUTHOR(S)
  Frank  H.  Y.  Wu & Robert  C.  Y.  Koh
                              8. PERFORMING ORGANIZATION REPORT NO.
9. PERFORMING ORGANIZATION NAME AND ADDRESS
  W.  M.  Keck Laboratory  of Hydraulics & Water  Resources
  California Institute of Technology
  1201  E.  California  Blvd.
  Pasadena, CA  91125
                              10. PROGRAM ELEMENT NO.

                               EHA540
                              11. CONTRACT/GRANT NO.
                               EPA (5)  R803989-01-1
12. SPONSORING AGENCY NAME AND ADDRESS
                                                            13. TYPE OF REPORT AND PERIOD COVERED
  U.S.  Environmental  Protection Agency    Corvallis,  OR
  Corvallis Environmental  Research Center
  200 SW 35th Street
  Corvallis, OR  97330
                              14
                                . SPONSORING AGENCY CODE


                               FPA/600/02
 15. SUPPLEMENTARY NOTES
 16. ABSTRACT

   A mathematical model  is  developed for the prediction  of plume properties such  as
   excess plume temperature,  humidity and liquid phase moisture (water droplet),  plume
   trajectory, width, and dilution at the merging  locations and the beginning and ending
   points of the visible part of the plumes.  Detailed printout and contour plots of
   excess temperature and moisture distribution can  also be obtained if desired.

   Based on  comparison with laboratory data this model gives good predictions for the
   case of dry plumes (no moisture involved).  It  should be noted that several empirical
   coefficients are as yet  not accurately known.   Verification of this model for  the wet
   plume  (such as for prototype cooling tower plumes) and the determination of the value:
   for these empirical coefficients to be used in  prototype applications must await
   detailed  comparison with field data.
                                 KEY WORDS AND DOCUMENT ANALYSIS
                   DESCRIPTORS
   Cooling Towers
   Plumes
   Mechanical  Draft
   Merging
   Thermal Pollution
   Cooling Water
 18. DISTRIBUTION STATEMENT


   Release to Public

 EPA  Form 2220-1 (R.v. 4-77)PREV.OUS EO.T.ON .s OBSOLETE
                                               b.lDENTIFIERS/OPEN ENDED TERMS
                    Mathematical Model
                                                19. SECURITY CLASS (This Report/
                  20. SECURITY CLASS (Thii page)
                     Unclassified
13B
                                             21. NO. OF PAGES

                                              134
                                            22. PRICE
                 119
                            U S GOVERNMENT POINTING OFFICE: 1978—797-35*'195 REGION 10

-------