'Jp.iteo c 2"-;;s             -•r'rcn 2                Mew Jtirsev  "''lew York
Er~".';rc"r"»r ;a; PrfitSCtion     Z3 ~-3'1--r3! -'•.•;: 3           ?'., ir'o ; •; j
                       New /ork NV 10007        Virgin isiancs

FEDBAKS33   -  A COMPUTER PROGRAM FOR THE MODELLING OF FIRST ORDER
              CONSECUTIVE REACTIONS  WITH FEEDBACK UNDER A STEADY
              STATE MULTIDIMENSIONAL NATURAL AQUATIC SYSTEM
                 PROGRAM DOCUMENTATION AND USERS GUIDE
                                   George A. Nossa
                                   Chief, Environmental Systems Section
                                   Data Systems Branch
                                   Management Division
                                   U.S. Environmental  Protection Agency
                                   New York, New York  10007
                                   November 1978

-------
                                   -1-
                            TABLB OF CONTENTS

Abstract ------_______-_______________---  2
Introduction ------____----------__--------  2
Theory --------_____----___-_-_-_-------  2
Acknowledgements ----------------------------  6
References --------_-_---------___--------  6

Description of the Main Program and Accompanying subroutines ------  7
Program Flowchart  ---------------------------  8
Source Program Listing -------------------------16
Input Description  --------------------------- 3Q
Restrictions ------------------------------49

Application to the Delaware Estuary  ------------------50
          Description of Area of Study -----------------51
          Nitrogenous Point Source Loadings  --------------53
          Program Output ----'--------------------54
          Graphs of Computed vs Observed Values  ------------ 72

Program Execution at the COMNET facility ----------------83

Acknowledgements ---------------------------6,84

References to the Computer Program -------------------85

APPENDIX A
          Listing of Input data for Delaware Estuary Verification  -  -  - 86

-------
                                                 -2-.

           FEDBAK03 - A COMPUTER  PROGRAM FOR THE. .MODELLING OF FIRST ORDER CONSECUTIVE REACTIONS  WITH
                       FEEDBACK UNDER A'STEADY STATE' MULTIDIMENSIONAL NATBRAL AOTJATIC SYSTEM *


                                                 George A. Nossa
                                              Environmental Engineer
                                               Data Systems Branch
                                       U.S. Environmental Protection Agency
                                                 New York, N.Y.
 ABSTRACT

 The computer model  described  1s used to compute the
 steady-state distribution  of  water quality variables
 undergoing consecutive reactions with feedback and
 following first order kinetics.  The program has
 been developed 1n a general form but 1s specifically
 applicable to the reactions observed by nitrogenous
 species and the associated dissolved oxygen uptake
 1n the natural environment.

 The basis for this  model is the theory of conserva-
 tion of mass.  The  approach used to solve the equa-
 tions is.a finite difference  scheme developed by
 Thomann 14»5/, which has been shown to be a very
 effective tool 1n the field of water quality manage-
 ment.

 INTRODUCTION

 A computer model  has  been  developed to serve as  a
 useful  tool  1n the  prediction of water quality para-
 meters  which react  under first order kinetics  and  as
 a system of consecutive" reactions .where any parameter
 can react 1n .a feedback fashion. The problem setting
 assumes an aquatic  environment 1n which steady state
 conditions can be applied

 Lets consider a system of  five reactants:
  Y

  f
             KAI
"1
            *
            "13-
                               f
                     14
     SCHEMATIC
FIGURE 1
OF FIVE REACTANT SYSTEM
 In the above  system« the consecutive feedfoward and
 feedback reaction rates are presented.  All  the
 possible reaction loops for the first reactant have
 also been Included.  This particular system  will  be
 used in the theory development.

 THEORY

 The general estuarine advection/dispersion equation
may be written as:
    SN  =  E62N - U6N - KN + WN(X)
    fit   "fix?    «x
                                         where:
                                               N - concentration of constituent
                                               t = time, 1n tidal cycles
                                               E = tidally averaged dispersion coefficient
                                                   which includes the dispersive effects  of
                                                   tidal motion
                                               U = net advectlve velocity
                                               K = first order decay coefficient of constitu-
                                                   ent N
                                            W(x) - direct discharges of N
                                                       Assuming steady  state conditions: oN/dt
                                                       equation (1)  becomes:

                                                            0 = E d2N - U dN - KN + W(x) 	
                                                                                   0 and
                                                                                      (2)
                                         Direct solutions for the above equation  have been
                                         computed by O'ConnorO .2),  a computer  program has
                                         been documented by EPA,  Region II (3) which uses this
                                         technique to solve for water quality parameters.
                                         A second solution approach developed by
                                         solves the above differential  equation directly by
                                         replacing the derivatives  with finite-difference
                                         approximations.   This  approach is  used by computer
                                         program HAR03; documented  by EPA,  Region III6) to
                                         analyze systems  of consecutive reactions.  The
                                         Thomann solution technique will  also be the approach
                                         followed in program FEDBAK03.

                                         If we take the first reactant N-j on figune 1 and we
                                         incorporate all  the feedback loops, equation (2)
                                         becomes :

                                           0 = E d2N]   .  U dNi
                                                 --
                                                      where Kij are the appropriate first order reaction
                                                      rate constants.  The incorporation of the reaction
                                                      term KH allows for the first order decay of this
                                                      first component out of the system.

                                                      For the other remaining reactants equation (2) would
                                                      be written as:
                                           Q =        _     z _ K22N2+K12N2+K32N3+K42N4+
                                                          dx
                                                         dx
                                                                    K52N5+W2
                                                                                  <53N5+W3
                                                                                      (5)
* Reproduced  from:   Environmental  Modelling  and Simulation,  USEPA Office of
  Research and Development,  Washington,  D.C., July 1976  (EPA600/9-76-016)

-------
               . U dN4 .

                   dx~
                         K54N5+W4
                                                 .(6)
              . U

                         <45VW5
                                                 •(7)
The Thomann solution approach  divides the system
into completely mixed segments, as illustrated for a
one dimensional estuary in figure 2.
                       FIGURE 2
           SEGMENTED ONE DIMENSIONAL ESTUARY

where segments 1 and 'n1 each form an interface with
the boundaries.  Equation (3) for the "1" segment on
figure (2) can be written as:
           3x2"
                      dx
                                                 (8)
The derivatives on the equation above ean be replaced
by finite-difference approximations giving:
         dx
where: E'1fj
                                                 (9B)

                                                 (9C)-
          i i = subscripts used to denote the
                interface -between adjacent segments
                1 and- j
          I 4 a average length of segments 1 and
          1>J   j
         a., j and-*j j are weight factors to
                correct the concentrations from
                equation 9B
                                                 (10)
                                                 (ID
In order-for the ^final-concentrations to be positive,
it is required that:

                        	 02)
                                                         -3-

                                                           Substltuting  equations  (9A)  and  (9B)  Into  equation
                                                           (8) yields:
                                                                    ,1N2>1+V1K31 f1N3,
                                                                                                            (13)
                                                          Grouping terms in the above equation yields:
 Letting:

   »i.i-r-Qi-i,iai-i.rEi-i,i

   ai ,i=(5i ,i+l°i .i+T^i-l ,iB1-l
                                                                                                           (14)
                                                          The general  equation  for the  1th  segment  becomes:
                                                                                         ,TN4+V1K51 ,1H5. .  (18)

                                                          The use of this  finite  difference  approximation
                                                          scheme has a  numerical  dispersion, which  can  be
                                                          approximated  as  (5)
                                                              num
                      ..........................  (ISA)

Where £.,^15 the -numerical  dispersion.   This  is  par-
ticularly important to stream applications  where
advective'-velocities may be high  and  this effect  may
lead to distorted results.  -

Extension to multi -dimensional  analysis:

If we consider a grid of orthogonally straped
sections such as the one illustrated  in  figure 3:
                                                                                k4
                                                                              FIGURE 3-
                                                                    HYPOTHETICAtfTWO'-DIMENSIONAL SYSTEM

-------
Foll owing  the convention  that flows  entering  a  sec-
tion  is  negative and flow out of  a section  is posi-
tive,  a  mass  balance due  to  the transport and disper-
sion  of  material from section 1 to all  surrounding
sections k(s) 1s:(4»5,8)
            =°=
        dt
                (1=1,2. ..n)
                                              t1N5§i

                                                  (19)
This equation 1s  the equivalent  to  equation  (13)  for
the one  dimensional  case.   The generalization  of  the
advection  term is possible  since:
and
                                                  (20)
                                                  (21)
Using  equation (19),  if the  terms  containing  the de-
pendent  variable N-|tf  are grouped  on the  left hand
side and the  direct loads of this  component and the
terms  for the formation of NT  due  to other components
are placed on the right hand side, one obtains :fe;)
where
                                                  (22)
                                                 (22A)
For sections where flow  enters a section from the
boundary with  a  concentration C|j
                                              ... (23)
and the forcing  function  at the boundary is added to
the direct  loads at  that  section by:
                                                 (24)
For sections forming a boundary, where flow leaves
this section to an area-with a concentration CD:
                                             ••' (25)
and
The set of equations  for component _Ni for n number of
spacial sections  in the system described in figure 1
would be:
allNl ,l+a!2N1 ,
                                >4+. . .+alnN1 >n=
      N5,l

a21Nl ,l+a22Nl ,2+a23Nl ,3+
                                            ,n
                                               ,2
                                                                           ,2+a33Nl ,
                                                                                             ,n=wl ,3+
                                                             an!Nl ,l
                                                                             +a+- • -
                                                                        ,2an3
                                                                                                        ,n
                                                                                                         (29)
                                                                                                         (30)
                                                          In matrix notation equations (27) thru (30) can be
                                                          written as:

                                                             [A1](N1)-(W1MVK21](N2)+[VK31](N3)H-[VK41]

                                                                  (N4KVK51](N5) ........................ (31)

                                                          where:

                                                          [A]] 1s a square matrix of n order, containing the
                                                               a's as defined on equations (22A) and (22B),
                                                               note that the main diagonal has the reaction
                                                               rate constant K-p
                                                          (N-|),(N2),...(N5) are nxl  vectors of the reactant
                                                               over all sections
                                                          (W-]) 1s an nxl vector of the waste loads for react-
                                                               ant NT for all sections
                                                          [VK2l],[VK31],[VK4l] and [VK^leach of these is an
                                                               nxn diagonal matrix of the section volume and
                                                               the first order reaction coefficient at that
                                                               segment.

                                                          A similar analysis as above for the second reactant
                                                          N2 on figure 1 yields:

                                                             [A2](N2)=(W2)+[VK12](N1)+[VK32](N3)-i-[VK42]

                                                                  (N4KVK52](N5) ........................  (32)

                                                          where:
                                                       [A2] is an nxn matrix similar to [AiJ, but the
                                                            main diagonal contains the reaction rate
                                                            constant K22
                                                       (W2) is an nxl vector of direct waste loads for
                                                            component N2 over the n sections
                                                          For the other reactants N3, N4,
                                                          are -generated:
                                                                                          similar equations
                                                                   d:

                                                          [AS] (NaHwaHv*! 3] (NT )+[VK23] (N2)+[VK43]

                                                               (N4KVK53](N5) ........................  (33)

                                                          [A4](N4)=(W4)+[VK14](N1)H-[VK24](N2)+[VK34]

                                                               (N3KVK54](N5) ........................  (34)
                                                                  (N3)+[VK45](N4)  ........................ (35)

                                                          The above matrix equations  (31)  thru  (35) can be
                                                          written as a matrix  of Matrices: (5, 8)
         N5,2

-------
     -CVK21J-[VK31]-[VK4lJ-[VK5l!
             IA3  3-[VK43]-[VK35]

     -[VK24J-[VK34] [A4  J-[VK54]

K15] -CVKzsHVKssKVMs] [A5  J
                                     (N2


                                     <»3

                                     (N4

                                     (N5
     (N)
                (W).
(W2)

(W3)

(W4)
 (37)
re [A~] is the 5n x 5n matrix above and (R) and
 are 5n x 1 vectors.  The solution of the five
ctants over all the spadal sections are given by

(N)  =   [AT1 (W) ., ......................... (38)

II cation of the theory by the computer program

 program described follows a modular approach 1n
:h the user specifies to the main line program the
ions desired, and subroutines are called accord-
 y to perform specific tasks.

 steps to be accomplished can be summarized as:

 nput the physical characteristics of the system;
lamely, the geometry, temperature, hydro! ogic
:haracteristics, reaction schemes and correspond-
 ng reaction rates.
 alculate E' and °>'s for all the sections as
 escribed on equations (9C) and (10).  In order to
 andle the constraint stated on equation (12), the
 rogram tests the expression:
if such is the case
                         is recalculated as:
                                              C40)
 h places °^j well within the tolerable range.
 et up the system matrix [Af] by computing its
 lements  as given on equations (22A) -and {22B);-
 t should be noted that the difference between  the
 A-j] matrices in equations (31)  thru -{35 ) is  the
 ddition  of a separate main diagonal  term KtKjj'f.
 n order  to conserve space, this matrix is set  up "
 ithout this term and during the creation of  the
 atrix of matrices [A],  the appropriate VfKjj^
 erm is added.
 et up the matrix of matrices [AJ,  this 1s done by
 ombining an offline disk file containing all the
 iKj  terms for the system and -the [A^] matrix.-
 nput of  direct discharges into  the system and
 oundary  concentrations,  and from these compute
 he system s ource -vector 
-------
                                                    -6-
 If we assume these reactions to follow first order
 kinetics, the system can readily be solved using
 program FEDBAK03.  The computation of dissolved oxy-
 gen deficit can be accomplished two ways.  A deficit
 "species" can be defined (noted Ng above), the decay
 of which is the reaeration rate,  Ka.  The reaction
 schemes producing deficit are then defined, and the
 corresponding reaction rate would be the product of
 the stochiometric coefficient by the reaction rate
 of the reaction using up oxygen.  A second method
 to compute deficit concentrations as done for com-
 ponent NI in equations (22) thru (31), one obtains

    [B](D)1-3.43[VK23](N2)+1.14[VK34](H3) 	 (43)

 where [B] is a matrix similar to the [A-j] matrices
 of equations (31)'thru (35), except that the mafn
 diagonal term has the reaeration rate Ka instead of
 *11i (N2) and (N?) are nxl vectors of the steady-
 state concentration of these reactants.  (D)^ is an
 nxl vector of the deficit concentrations over all
 segments due to the oxydation of ammonia and nitrite.
 The solution to the deficit concentration over all
 space 1s given by:

    (D)i=3.43[VK23](N2)[B]-1+1.14[VK34](N3)[B]-l.(44)

 This method is used to compute deficit in program
 FEDBAK, by using the optional subroutine.

 The application to nitrification and dissolved oxygen
 deficit assumed first order kinetics for the bacter-
 ial reactions.  This should be confirmed by laboratory
 studies, or the nature of the system should be care-
 fully considered. This computer model  has been found
 to be very useful as a predictive tool  and in provi-
 ding insights to the behavior of nitrogen species  in
 the aquatic environment.  On new applications this
 will  ultimately depend on the applicatively of the
 underlying assumptions to the system of interest.

 ACKNOWLEDGEMENTS

 The author is grateful  to Steve Chapra  of the Great
 Lakes Environmental  Research  Laboratories and Richard
 W1nf1eld of the Manhattan College Department of Envi-
 ronmental  Engineering and Science for their review
 and comments on this paper.

 The data for test applications  of this model  were  made
 available by Dr.  Richard  Tortoriello of  the Delaware
 River Basin Commission. His  input to this  project  is
 gratefully acknowledged.

 REFERENCES

 (1) O'Connor, D.J.,  "Oxygen Balance of an Estuary".
 Jour.  San.  Enq.  Div. ASCE. Vol  86, May 1960,  pp 35-55

 (2) O'Connor, D.J.,  "The Temporal and Spacial Distri-
 bution of  Dissolved  Oxygen in Streams". Water Resour-
 ces Research,-Vol 3, No. 1, 1967, pp 65-79.

 (3) Chapra,  S.C.  and Gordimer S., ES001 A Steady Sta-
 te. One  Dimensional, Estuarine Hater Quality Model.
 USEPA,Region  II, New York, N.Y. September 1973.

 (4) Thomann, R.V., "Mathematical Model for Dissolved
Oxygen". Jour. San. Eng.-Div., ASCE, Vol 89, No. SA5,
October  1963, pp 1-30.

 (5) Thomann, R.V., System Analysis and water Quality
Management, Environmental  Science Division, New York,
1971.
(6) Chapra, S.C. and Nossa, G.A. HAR03 A Computer
Program for the Modelling of Water Quality Parameters
in*Steady State Multidimensional Natural Aquatic Sys-
tem., USEPA, Region II, New York,'N.Y. October 1974

(7) Stratton, F.E. and Me Carty, P.L., "Prediction of
Nitrification Effects on the Dissolved Oxygen Balance
of Streams" Eny. Sci. and Tech.. Vol. 1, No. 5, May
1967, pp 405^4Tcn

(8) O'Connor, D.J., Thomann, R.V.,  and D1 Toro, D.M.,
Dynamic Water Quality Forecasting and Management.
USEPA, Office of Research and Development,  Wash.,D.C.
Report No.  EPA-660/3-73-009, August 1973

-------
                 -7-
DESCRIPTION OF THE MAIN PROGRAM AND
      ACCOMPANYING SUBROUTINES

-------
                                                         -8-
                                                                   Mt-v/RE fc.

    R.£At>
CALL S«j




TO «»»
L
^x^
:D g«T{^.
,X1
C*t-<- iu«. DOOCf
To e.t»M«srtC BlF'tCtT

T*rAL "OCfi^lf /*>
BACU .3£cri£>A*


-------
                                    -9-
The Main Program

    FEDBAK03 has been programmed using a modular approach.  Wherever
possible, an attempt has been made to give each subroutine a unique function
The main program was designed to coordinate the various subroutines used  in
the model.  This inter-relationship among the various subroutines  is  illus-
trated in figure 6.

Subroutine DATA...

    The primary purpose of this subroutine is to input parameters  which do
not depend on the particular substance being modeled.  These include  system
parameters  (number of segments, indicators, etc.) and certain physical par-
ameters such as cross-sectional area, net advective  flow, etc., which are
only  input once when modeling a water body.  These are in contrast to par-
ameters such as loads and reaction rates which would vary with the constitu-
ent being modeled.

    The particular tasks accomplished by DATA are:

    1)  Read and write  system parameters
    2)  Initialize certain matrices to  zero
    3)  Read and write  interface parameters AREA, E  and Q and place them
        into the storage matrix ARRAY.
    4)  Read and write  the characteristic lengths  (LA)
    5)  Read chloride concentrations  (CHLOR), volume (VOL),  and temperature
         (TEMPER) at  each segment
    6)  Apply the  following  conversion  factors:

          Q(MGD) =  (CFS)  *  .6463 MGD/cfs

          VOL(MG) =  (106ft3)  *  7.48  gal/ft3

-------
                                   -10-
    At this point, a further word about the function of ARRAY  is  in order.
ARRAY is an N X 18 sized matrix which  is used to save space when  storing
the interface parameters AREA, E and Q.  Used in conjunction with IARAY,
ARRAY avoids the use of three N X N matrices when handling these  variables.
It does this as follows:

    IARAY(I,J)                      contains the section number of the
                                    section which forms the interface with
                                    section I.

    ARRAY(I,J)                      contains the area of interface
                                     (I,IARAY(I,J))

    ARRAY(I,J + 1)                  contains the E of interface
                                     (I,IARAY(I,J))

    ARRAY(I,J + 2)                  contains the Q of interface
                                     (I,IARAY(I,J))

    J                               is an integer variable which  is
                                    incremented by 3 from 1 to 16

    For computational purposes, every  interface must be accounted for in
terms of the section to which it applies.  This would involve duplicate
input - for instance, AREA  (I,J) is the same as AREA (J,I).  To simplify
input, only one of the values must be  input and subroutine DATA fills out
the remaining positions in ARRAY.

Subroutine FLOW...

    FLOW basically computes a balance of the flows into and out of a
section.  If there is more or less than a zero balance around the section,
the resultant value is printed as an excess flow.

    The purpose of this exercise is to both keep track of flows and to
insure that no mistakes were made when inputting the flow regime.  This can
be particularly useful when working in two or three dimensions.

-------
                                   -11-

5 ubrout in e S ETUP . . .

    The primary purpose of this subroutine  is  to  set  up  the  individual
flownet matrices  (A) ,  (described  in equations  31  through 35  in  the  theory
section) and output these on a direct  access device.   The tasks performed
by SETUP are:

    1)  Calculation of the bulk dispersion  coefficient,  E1

           According to the theory:

                                  E.  . A.  .
                        „•
                        E      =
           Where

                        E.  . = the dispersion  coefficient across  the  interface
                         !f J
                        A.  . = the cross-sectional area of the  interface
                         1 1 D
                        1.  . =  the average length of the adjoining
                         1'11    sections =  (1. + l.)/2


           The program calculates E '  (EPRIM) as follows :

                EPRIM( J) =AKRAY (I , JJ) * ARRAY (I , JJ+1) *417 . 166/ (LA(K, I) +LA (I ,K) )

           Where

                          ARRAY(I,JJ) = A.  .  (ft2)
                                         1 / 3
                          ARRAY (I, JJ+1) = E.  . (mi2/day)
                                           1 1 J
                          LA(I,K)  = 1.

                          LA(K,I)  = lk

           To conserve computer memory these lengths are stored as a  N X  6
           sized array and used in conjunction with array NLA to  identify
           the N,I interface.  This again avoids the use of an  N  X N  matrix
           to handle this variable.

                          417.1166 - conversion factor which is computed
                                     as follows :

           Original Units _ Conversion Factor _ New Units


           ft2*mile2          (5280 ft) 2 x 7.481 gal y 2 x MG      MGD
           day*ft*2            (mile) ^        ft3         10  gal

-------
                               -12-
2)   A value of alpha is calculated to adjust the advective term for
    sections of different length.  This value is checked against the
    positivity criterion (equation 12) and if it is not met, alpha is
    recalculated.

    Alpha(J) = 1.0-(EPRIM(J)/(2*ABS(ARRAY(I,JJ+2))))

    Which corresponds to

          o=l- E'/2Q

3)   The elements of the system matrix (A) are calculated and saved on
    mass storage

4)   The values for EPRIM and ALPHA are printed

5)   The section parameters:  volume, temperature and chloride con-
    centration are printed.

-------
                                    -13-
SUBROUTINE SOURCE

    This subroutine will read the reaction schemes and corresponding
reaction rate constants for each system segment and construct the  (VK)
system matrices.  It also reads all direct discharges and boundary  con-
ditions, to form the system waste vector  (W).  These steps are accomplished
in the following sequence:

    1.  The reaction coefficients are read and are corrected to the tem-
        perature in the section by the formula:


                         Ki - Ki 9(t-20)


        For the nitrification reactions, values of 9 are discussed  in
        reference  (2).  The corrected reaction rates are then printed.

    2.  The (VK) matrix is formed by multiplying by the corresponding
        section volumes and saved in mass storage.

    3.  The direct discharges of the various reactants into the system are
        read and used to construct the system source vector.  The contribu-
        tion to the source vector due to the boundary conditions are calcu-
        lated from the formulas given in equations (24) and  (26).

    4.  The system source vector is then printed out.

-------
                                        -14-
 Subroutine MATRIX

    This subroutine will:

    (1)   Combine the (VK)  matrices on mass storage with the (A)  matrices  to
         form the system matrix of matrices (A)  given in equation (37)
    (2)   Optionally print  the (A) matrix
    (3)   Read mile points corresponding to each of the segments

 Subroutine MATN

    This subroutine inverts the (A) matrix and multiplies it by  the source
    vector (W)  giving the solution vector (N).  This solution' vector represents
    the component concentrations at the various segments.  The printed output
    includes the mile points corresponding to the beginning of each segment.

 Subroutine SENS
    This subroutine will perform system sensitivity analysis on the reaction
    rate constants and/or direct discharges into the system.  The logical se-
    quence of this subroutine is:

 For Changes in direct discharges;

 1.)   Read and print revised loadings for section
 2.)   Modify original source vector to reflect the new loadings
 3.)   Multiply new waste vector by  the inverse of system matrix
 4.)   Print new solution vector
 5.)   Optionally,  the deficit analysis subroutine KBOD may be called

 For  Changes  in reaction rates;

 For changes of direct discharges and reaction rates, the first two steps
 above take place,  the procedure below is then followed:

 1.)   The  matrix of matrices (A)  is re-assembled using the  (A)and (VK)
      matrices on  mass storage files.
 2.)   The  revised  reaction constant(s) at all segments in the system are
      read and temperature corrected.          _
 3.)   The  above values are substituted in the (A) matrix. (Note:  Sub-
      sequent  calls  to this subroutine make the changes to the original (A)
      matrix and not from the one modified in this step.).  For cummulative
      effects  this new (A)  matrix should be copied to the mass storage file
      containing the old (A)  matrix._
4.)  Optionally print the revised  (A) matrix

This subroutine would then return  to the main line program which in turns
calls the MATN  and  optionally the  D.O.  deficit analysis routine.  Then if
desired it returns  back to this  routine for multiple runs of sensitivity
analysis.

-------
                                    -15-
Subroutine DODEF

   This subroutine will compute Dissolved Oxygen Deficit due to  selected
   reaction schemes.  For each section of the system this routine will de-
   termine the deficit due to each reactant and the resulting sum total
   deficit in the system.  The individual steps of this routine  are:

   (1)  Read the reaction coefficients at each section, and temperature cor-
        rect these by the expression:^

        KA(I)=KA(I)*1.02488(TEMPER(I)-20.)

   (2)  The temperature and reaeration rates are printed
   (3)  The saturation level of dissolved oxygen at each system  segment is
        computed from the expression.  ' '

        CS (!) = (!. -.000009*CHLOR(I))*(14.652-.41022*TEMPER(I)+-
               1.0079910*TEMPER(I)**2-.000077774*TEMPER(I)**3)

        Where:

          TEMPER(I) - is the temperature of the section
          CHLOR(I)  - is the chloride concentration (mg/1) of the section

   (4)  A system deficit matrix is formed by adding the  VKa  term to the
        (A) matrix on the mass storage file (N.B. the VK term is omitted
        for flexibility on the Flownet matrix (A) on mass storage)

   (5)  The above deficit matrix is  inverted and optionally printed.  The
        structure of the program logic requires a separate array to f&rm the
        system deficit matrix and compute its inverse.  In the repeated use
        of this subroutine during system sensitivity analysis, this inverse
        of the deficit matrix is saved for subsequent use.
   (6)  Background contributions of  deficit such as benthal demand or direct
        sources are read into the section, the deficit source vector is
        formed and multiplied by the inverse of the deficit matrix giving
        the system deficit response  due to direct loads.  These  results are
        then printed.
   (7)  Next, the number of reactions and the reaction paths which produce
        deficit as well as -the-"stochicmefcric coefficients are read.
   (8)  The (VKij)  matrix for the deficit reactions is formed by accessing
        the large (A) on mass storage.
   (9)  The deficit contribution to  each reaction scheme is computed as
        shown in equation (44)  for the oxydation of ammonia.
  (10)  The resulting deficit and the dissolved oxygen concentration is
        printed and accumulated for  the calculation of total deficit.

-------
                                    -16-
                          SOURCE PROGRAM LISTING
NOTE:  The source code and the sample test case can be obtained in a
       computer readable format for a nominal fee ($20.00 at present) by
       writing to the office originating this document.  This fee is waived
       for federal agencies or local governmental planning organizations.

-------
                                   -17-
   ******************************************************************

   THIS PROGRAM SOLVES A SYSTEM OF FIRST ORDER CONSECUTIVE REACTION
   WITH FEEDBACK IN A MULTI-DIMENSIONAL AQUATIC ENVIRONMENT

           *** FEDBAK03 ***

   THIS VERSION 2.0 CAN SOLVE UP TO A 60 SEGMENT SYSTEM, WHERE THE NO.
   OF COMPONENTS * NUMBER OF SEGMENTS IS EQUAL OR LESS THAN 180

   ******************************************************************
   INTEGER OUT
   DIMENSION TITLEC20)
   COMMON/BLK1/NR,N,NP,INPUT,OUT
   COMMON/OPT I ON/IPRINT,ISENS,KWDEF,KBOD,NRUNS
   DEFINE FILE 1(60,60,U,IREC1)
   DEFINE FILE 3(180,180,U,IREC3)
   ***************************************************************  *
                                                                    *
   ASSIGN DEVICE NUMBERS TO  CARD READER, LINE PRINTER ...          *
                                                                    *
   ******************************************************************
   INPUT=8
   OUT=6
   IDODSW=0
   READ*INPUT,200)TITLE
00 FORMAT (20A4)
   N= NO. OF SEGMENTS (SECTIONS)
   NR = NO. OF COMPONENTS (REACTANTS)
*
*
*
*
   READ(INPUT,9999)N,NR
9  FORMAT(2I2)
   WRITE (OUT,201) TITLE,N,NR
01 FORMAT(1H1,20X,20A4//28X,'NUMBER OF SECTIONS = »,13,10X,'NUMBER OF
  1 REACTION COMPONENTS = f,I2,//)
   NP=N*NR
   READ PROGRAM OPTIONS
   READUNPUT,1101)IPRINT,ISENS,K8QD,KWDEF
01 FORKAT(4I2)
   ******************************************************************
                                                                    *
   READ SYSTEM GEOMETRY — DISPERSIONS                              *
                                                                    *

   CALL DATA
   ******************************************************************
                                                                    *
  CHECK FLOW REGIMES                                                *
                                                                    *

   CALL FLOW
   ******************************************************************
MAIN001
MAIN002
MAIN003
MAING04
MAIN005
MAIN006
MAIN007
MAIN008
MAIN009
MAIN010
MAINOll
MAIN012
MAIN013
MAIN014
MAIN015
MAIN016
MAIN017
MAIN018
MAIN019
MAIN020
MAIN021
MAIN022
MAIN023
MAIN024
MAIN025
MAIN026
MAIN027
MAIN028
MAIN029
MAIN030
MAIN031
MAIN032
MAIN033
MAIN034
MAIN035
MAIN036
MAIN037
MAIN038
MAIN039
MAIN04Q
MAIN041
MAIN042
MAIN043
MAIN044
MAIN045
MAIN046
MAIN047
MAIN048
MAIN049
MAIN050
MAIN051
MAIN052
MAIN053
MAIN054
MAIN055
MAIN056
MAIN057
MAIN058
MAIN059

-------
                                    -18-
   SET UP SYSTEM MATRIX
                                               #
                                               *
                                               *
   CALL  SETUP
   ******************************************************************
                                                                     *
   READ AND PRINT  LOADS  +  REACTION  RATES —  COMPUTE  SOURCE  VECTOR    *
                                                                     *
   ******************************************************************
   CALL SOURCE
   ******************************************************************
                                                                     *
   NOW FORM LARGE  MATRIX OF  KV  TERMS  AND  (A)  MATRICES                *
                                                                     *
   ******************************************************************
   CALL MATRIX
    COMPUTE  INVERSE  OF  LARGE  MATRIX   *   SYSTEM  CONCENTRATION  VECTOR

    *********a
    CALL  MATN
    *********5

    TEST  FOR DISSOLVED  OXYGEN DEFICIT ANALYSIS
    IFCKBOD.GT.OJCALL  DODEF
40  IF( ISENSUOO, 100,801
    READ  NUMBER  OF  SENSITIVITY  RUNS  DESIRED
    ********************:
101  READ( INPUT,30DNRUNS
101  FORMATU2)
102  CALL  SENS
104  CALL  MATN
510  IF(KBOD.GT.O)CALL  DODEF
    IF(NRUNS.GT.O)GO TO 302
LOO  CALL  EXIT
    END
    SUBROUTINE  DATA
                                               *
                                               *
                                               *
                                               *
                                               *
    INPUTS DATA  WHICH
    CONSTITUENT  BEING
    HYDRODYNAMIC DATA
       SUBROUTINE DATA

WOULD NOT VARY WITH THE
MODELED, I.E. PRIMARILY
PARTICULAR
PHYSICAL AND
                                                                     *
                                                                     *
                                                                     *
                                                                     *
                                                                     *
                                                                     *
                                                                     *
>*********************************************************************
   I*******:
    INTEGER OUT
    REAL  LA(60,6)
    DIMENSION AREA(6),E(6),Q(6)
    COMMON/8LK1/NR,N,NP, INPUT,OUT
    COMMON/BLK3/VOL(60),TEMPER(60),CHLOR(60)
MAIN060
MAIN061
MAIN062
MAIN063
MAIN064
MAIN065
MAIN066
MAIN067
MAIN068
MAIN069
MAIN070
MAIN071
MAIN072
MAIN073
MAIN074
MAIN075
MAIN076
MAIN077
MAIN078
MAIN079
MAIN080
MAIN081
MAIN082
MAIN083
MAIN084
MAIN085
MAIN086
MAIN087
MAIN088
MAIN089
MAIN090
MAIN091
MAIN092
MAIN093
MAIN094
MAIN095
MAIN096
MAIN097
MAIN098
MAIN099
MAINIOO
MAIN101
DATA001
DATA002
DATA003
DATA004
DATA005
DATA006
DATA007
DATA008
DATA009
OATA010
DATA011
DATA012
DATA013
DATA014
DATA015
DATA016
DATA017

-------
                                    -19-
INTERFACE », 13, •-
                                                                  13.
   COMNON/COHOE/LA,SCALE(4),ARRAY(60,18)
   COMMON/OPT I ON/1 PR INT,ISENS,KWDEF,KBOD
   COMMON/UT1/IARAY(60,6)
   COMMON/BLK4/NLA(60,6)
31 FORMAT(SFIO.O)
CO FORMAT14F7.0)
01 FORMAT (///28X, 'ZERO SEGMENT NUMBER IN
  1« COMPUTATION DISCONTINUED')
00 FORMAT(6(I4,' -•,14,2X,F7.0,2X))
01 FORMAT(6(I3,F10.0))
00 FORMAT!              »1
  1 CHARAC. LENGTHS OF SEGMENTS (FT)
  IFACE  LENGTH    INTERFACE  LENGTH
  2E  LENGTH    INTERFACE  LENGTH')
11 FORMAT(3(3F6.0,I3))
02 FORMAT(//10X,'**SCALE FACTORS**'//10X,
  l'SCALE(l) =  ',F10.3,4X,'SCALE(2) =  »,F10.3,4X,
  2'SCALEt3) =  •,F10.3,4X,«SCALE(4) =  «,F10.3  )
00 FORMAT!»1«,2X,'INTERFACE',5X,'AREA',7X,'E',8X,«Q •,2(6X,'INTERFACE
  l^SX,'AREA*,9Xt>E>*7XttQ*)/,4X,'ROW-SEG',4X,•(FT**2)', 2X,
  2MMI**2/D)',2X,MCFS)',   6X, • ROW-SEG', 5X, • (FT**2 )',2X,
  3 M MI**2/D>»,2X,MCFS)',5X,» ROW-SEG', 5X, ' ( FT**2 ) ' ,2X, ' (MI**2/D)' ,
  4 2X,MCFS)»  )
100 FORMAT*  •   (•,I3,f-',I3,M',2X,F8.0,4X,F6.3,F8.1,2(5X,M',I3,
  1 '-'.IS,' )',2X^8.0,4X^6.3, F8.D)
                                      •///'  INTERFACE  LENGTH
                                         INTERFACE  LENGTH
                     INTER
                  INTERFAC
    READ AND WRITE THE SYSTEM PARAMETERS AND SCALE FACTORS
    READ
-------
                                     -20-
 02  ARRAY(I,JJ+2)=Q(J)*SCALE(3)
    ******** $*:). $ $ ** * ******:**:$*;**#**:£**:$: * $**$#**********************


    PROPAGATE AREAS,DISPERSIONS AND FLOWS FOR INTERFACES AT 
-------
                                   -21-
I**********$$********************************************************
   IF(KBOD.GT.O)READ(INPUT,31)(CHLOR(K),K=1,N)
   READtINPUT,31)   tTEMPER!I),VOL(I),I=1,N)
   IFITEMPER(2))14,13,14
.3  DO  14  1 = 1,N
   TEMPERU)=TEMPER{1)
L4  CONTINUE
   APPLY  CONVERSION  FACTORS  TO  FLOW  AND  VOLUME
   DO  12  1=1,N
   DC  23  J=U6
   JJ=(J-1)*3.+ 1
>3  ARRAYU,JJ+2) = ARRAY(I,JJ+2)*.6463
L2  VOL{I)=VOL(I)*7.48
   RETURN
   END
   SUBROUTINE  FLOW
                         SUBROUTINE  FLOW

   THIS SUBROUTINE CALCULATES  A FLOW BALANCE  AROUND  EACH  SECTION
   TO INSURE THAT THE FLOW REGIME HAS BEEN  INPUT  CORRECTLY
REAL LA(60,6)
INTEGER OUT
COMMON/8LK1/NR,N,NP,INPUT,OUT
COMMON/COHOE/LA,SCALE(4),ARRAY(60,18)
NMN=0
               SECTION ',13,' HAS AN EXCESS
                                            FLOW
14 FORMATi1
  !• CFS»)
02 FORMAT(*1«)
   DO 213 1=1,N
   QQ=ARRAYU,3)+ARRAY(I,6) + ARRAY(I,9)+ARRAYU,12)
  1+ARRAY(I,18)
   GQ=CQ/.6463
   IF4QQ-.001)216,37,37
16 IF
-------
                                    -22-

                                        ************^

   REAL LA(60,6)                                                     SETUP011
   INTEGER OUT                                                       SETUP012
   DIMENSION EPRIM(6),ALPHA(6)                                       SETUP013
   COMMON/BLK1/NR,N,NP,INPUT,OUT,AB(60,60)                           SETUP014
   COMMON/BLK3/VOL(60),T{60),CHLOR(60)                               SETUP015
   CCMMON/COHOE/LA,SCALE(4),ARRAY(60,18)                             SETUP016
   COMMON/OPTION/I PR INT,I SENS,KWDEF,KBOD                             SETUP017
   COMMON/UT1/IARAY(60,6)                                            SETUP018
   COMMON/BLK4/NLA(60,6)                                             SETUP019
23 FORMAT(6(I3,'-«,I3,F6.0,F6.3,2X))                                 SETUP020
24 FORMAT('i',49X, 'VALUES OF ALPHA  AND EPR IM» ,///,6( •  INTRFC EPRIM ALSETUP021
  1PHA  M,/,6(8X,'(MGD)',8X> }                                      SETUP022
42 FORMAT(4X,I4,9X,F12.2,16X,F5.2,23X,F10.0)                         SETUP023
43 FORMAT('1',3X,'SECTION »,8X,'VOLUME (MGAL)',8X,                   SETUP024
  1'TEMPERATURE  (CENT.)',8X,'CHLORIDE CONC  (MG/L)'   )                SETUP025
44 FCRMAT('l*,3X,iSECTION »,8X,'VOLUME (MGAL)«,8X,                   SETUP026
  1'TEMPERATURE  (CENT.)')                                            SETUP027
00 FORMAT(///,10X,'SEGMENT  (•,13,*-•,13,')  IS A BOUNDARY - BUT NO LENSETUP028
  IGTH WAS SPECIFIED')                                               SETUP029
01 FORMAT(//,10X,'PROGRAM TERMINATING DUE TO  ',13,'  BOUNDARY CONDITIOSETUP030
  INS WITH MISSING  LENGTHS'1                                         SETUP031
   WRITE  (OUT,24)                                                    SETUP032
   IEXIT=0                                                           SETUP033
   DO 28  1=1,N                                                       SETUP034
   DO 99  J=1,N                                                       SETUP035
99 AB(I,J)=0.0                                                       SETUP036

                                                                    *SETUP038
   CALCULATE EPRIM  AND  ALPHA FOR EACH INTERFACE                     *SETUP039
                                                                    *SETUP040

   DO 20  J=l,6                                                       SETUP042
   EPRIM(J)=0.0                                                      SETUP043
   ALPHA(J)=0.0                                                      SETUP044
   JJ=(J-1)*3+1                                                      SETUP045
   K=IARAY(I,J)                                                      SETUP046
   IF(K)20,20,10                                                     SETUP047

                                                                    *SETUP049
   GENERATE  (I,K) LENGTH ELEMENT                                    *SETUP050
                                                                    *SETUP051

10 XLEN1=0.0                                                         SETUP053
   DO 77  KK=1,6                                                     SETUP054
   IF(NLA(I,KK)-K)77,81,77                                           SETUP055
81 XLEN1=LA(I,KK>                                                    SETUP056
   GO TO  82                                                          SETUP057
77 CONTINUE                                                          SETUP058
!****«****************************************************************SETUP059
                                                                    *SETUP060
   GENERATE  (K,I) LENGTH ELEMENT                                    *SETUP061
                                                                    *SETUP062
82 XLEN2=0.0                                                         SETUP064
   DC 78 KK=1,6                                                      SETUP065
   IF(NLA(K,KK)-I)78,84,78                                           SETUP066
84 XLEN2=LA(K,KK)                                                    SETUP067

-------
                                    -23-

   GC TO 83                                                          SETUP068
T8 CONTINUE                                                          SETUP069
33 IFiXLENl+XLEN2)87,87,85                                           SETUP070
37 WRITE(OUT,100n,K                                                 SETUP071
   IEXIT=IEXIT+1                                                     SETUP072
   GO TO 17                                                          SETUP073
35 EPRIM(J)=ARRAYU,JJ)*ARRAY(I,JJ+l)*417.1166/sAB(ItI)+ALPHA(,EPRIM(2),ALPHA<2),    SETUP103
  2I,IARAYU,3),EPRIM(3),ALPHA(3>,I,IARAY(I,4),EPRIM(4),ALPHA(4),    SETUP104
  2I,IARAY(I,5),EPRIMt5),ALPHA{5),I,IARAYII,6),EPRIM<6},ALPHAC6)     SETUP105
   LINE=LINE-H                                                       SETUP106
   IFCFLOATCLINE/54)-FLOAT{LIN£)/54.J28,1028,1028                    SETUP107
28 WRITECOUT,23)                                                     SETUP108
   LINE=0                                                            SETUP109
28 CONTINUE                                                          SETUP110

                                                                    *SETUP112
   WRITE FLOWNET MATRIX ON MASS STORAGE DISK FILE 1                 *SETUP113
                                                                    *SETUP114
K********************************************************************SETUP115
   DC 30  1=1,N                                                      SETUP116
   WRITE(ltI)(AB(ItJ),J=l,N)                                         SETUP117
30 CONTINUE                                                          SETUP118

                                                                    *SETUP120
   WRITE SECTION PARAMETERS                                         *SETUP121
                                                                    *SETUP122
**********************************************************#**********SETUP123
   IF(KBOD.GT.O)WRITE(OUT,43)                                        SETUP124
   IF(KBOD.LE.O)WRITE(OUT,44)                                        SETUP125
   LINE=0                                                            SETUP126

-------
                                     -24-
                  WRITE(OUT,42)
                  WRITE
-------
    WRITE (VK)  COEFFICIENTS
    STORAGE  FILE 3
                                 -25-
                            FOR CORRESPONDING  REACTION  ON MASS
L01
 95
100
    DO  100  L=1,NS
    IP=NS*(J-1)+L
    READ(3'IP)(B(II),II=1,NP)
    JP=NS*(I-l)-»-L
    B(JP)=V(L)*K(L)
    IF(I-J)101. 95,101
    B(JP)=-B(JP)
    WRITE(3'IPMBUC),IC=1,NP)
    CONTINUE
    GO  TO 200
    READ  DIRECT COMPONENT DISCHARGES INTO  THE  SYSTEM
 89 DO 61 JK=1,NP
 61 W
-------
                                    -26-
GO TO 3
                                                                      SOURC105
GENERATE LA(K,K) LENGTH ELEMENT
                                                                     *SOURC107
                                                                     *SOURC108
                                                                     *SOURC109
 35
                                                 BOUNDARY CONDI
   XLEN1=0.0
   DO 77  KK=1,6
   IF(NLAU,KK)-L)77,88,77
   CONTINUE
   WRITEtOUT,102)
   CALL EXIT
   FCRMAT{//,10X,'PRQGRAM TERMINATING DUE TO  A
  INS WITH MISSING LENGTHS')
88 XLEN1=LA
-------
                                    -27-
    DO  60  L=1,NS                                                      MATRX019
    LL  =  NS *(J-1)+L                                                  MATRX020
»********************* jjc*^************************************** ****** *MATRX 021
                                                                     *MATRX022
   FORMATION OF THE LARGE MATRIX WHICH CONSISTS OF KV AND A MATRICES.*MATRX023
                                                                     *MATRX024
 60

 93
 20
L15

 70
 80

 89
 12
 39
JOO
 21
 20
B(LL)=B(LL)+A(L)
DO 93 JC=l,NP
AC(JJ,JC)=B(JC)
CONTINUE
CONTINUE
IFUPRINT)90,90,70
WRITEIR,(AC(IR,JC),JC=1,NP)
FORMAT(IX,'ROW NO.',13,IX,10F12.2/<12X,10F12. 2))
CONTINUE
MATRX026
MATRX027
MATRX028
MATRX029
MATRX030
MATRX031
MATRX032
MATRX033
MATRX034
MATRX035
MATRX036
MATRX037
MATRX038
MATRX039
                                                                     *MATRX041
    PTMIL IS AN ARRAY CONSISTING OF THE MILEAGE POINTS OF THE SECTION*MATRX042
                                                                     *MATRX043
 90
>00
READ  (INPUT, 200HPTMIL(I),I = 1,NS)
FORMAT  (8F10.0)
RETURN
END
SUBROUTINE MATN
                           SUBROUTINE MATN

    THIS SUBROUTINE INVERTS THE MATRIX AC AND MULTIPLIES
    IT BY THE FORCING FUNCTION WS(I)  GIVING SOLUTION VECTOR BU>
INTEGER OUT
COMMON/UT1/IPIVO(180),PIVOT(180)
CCMMON/UT2/INOEX<180,2)
COMMON/BLK1/N1,N2,N,INPUT,OUT,AC(180,180)
COMMON/OPTION/1 PR INT
COMMON/BLK4/PTMIL(60),B(180>
COMMON/BLK5/WS(180)
DO 900 1=1,N
DO 900 J=1,N
AC(ItJ)-AC(I»J)/1000000.
DO 20 J=1,N
IPIVO(J)=0
DO 21 1=1,2
INDEXfJ,I)=0
CONTINUE
DO 550 1=1,N

SEARCH FOR PIVOT ELEMENT
MATRX045
MATRX046
MATRX047
MATRX048
 MATN001
 MATN002
 MATN003
 MATN004
 MATN005
 MATN006
 MATN007
 MATN008
 MATN009
 MATN010
 MATN011
 MATN012
 MATN013
 MATN014
 MATN015
 MATN016
 MATN017
 MATN018
 MATN019
 MATN020
 MATN021
 MATN022
 MATN023
 MATN024
 MATN025
 MATN026
 MATN027
 MATN028
 MATN029

-------
                                    -28-
    AMAX=0.0
    DO  105  J=1,N
    IF  (IPIVO  (J)-l)
SO

BO
85
 DO
 05
50
00
60
 CO
 50
 150
>30
105
740
901
                    60,  105, 60
   DO  100 M=1,N
   IF  (IPIVG (M)-l)  80,  100, 740
   IF  (ABS(AMAX)-ABS(AC(J,M)))85,
   IROW=J
   ICOLU =M
   AMAX=ACU,M)
   CONTINUE
   CONTINUE
   IPIVO (ICOLU )=IPIVO (ICOLU
                                   100,  100
    INTERCHANGE ROWS TO PUT PIVOT ELEMENT ON DIAGONAL

    IF  UROW-ICOLU > 150,  260,  150
    DO  200 L=1,N
    SWAP=AC(IROW,L)
    AC(IROW,L)=AC(ICOLU,L)
    AC(ICOLU,L)=SWAP
    INDEX(I,1)=IRQW
    INDEX(I,2)=ICOLU
    PIVOT(I)=AC( ICOLU, ICOLU)

    DIVIDE PIVOT ROW BY PIVOT ELEMENT
    AC(ICOLU,ICOLU)=1.0
    DO 350 L=1,N
 50 AC ( ICOLU, L) = AC( ICOLU, L)/PIVOT( I)

    REDUCE NON-PIVOT ROWS
   DO 550 L1=1,N
   IF(L1-ICOLU ) 400, 550, 400
   S=AC
-------
                                    -29-
01 WRITE
88 FQRMATt'l',50X,'SYSTEM SOLUTION VECTOR',/, 25X,
  1'C OMPONENT   CONCENTRATION  (MG/L)',
  2/,«  SECTION',5X,»MILE PT',4X,8I8)
   DO 61 1=1, N2
61 WRITE
-------
                                    -30-
30 READUNPUT,31)LOAD,ISEG,IREACT
31 FORMAT(F10.0,212)
   IF
-------
                                    -31-
  2 'THETA',5X,13I8,/,4(13X,13I8,/M
   DC 115  1=1,NS
   READ(1»IJROW
   DC 29 J=1,N
   JJ=NS*(J-1)+ I
   READ(3«JJ)ROW2
   DO 60 L=1,NS
   LL = NS  *{J-1)+L
   FORMATION  OF THE  LARGE  MATRIX  WHICH CONSISTS OF  KV  AND ROWMATRICES*
                                                                     *

 60 ROW2(LL)=RCW2(LL)+ROWU)
    DO  93  JC=1,NP
 29 CONTINUE
 15 CONTINUE
    DO  25  I=1,NRATES
    READdNPUT,2>II,KK,THETA
  2 FORMAT(2I2,1X,F5.3)
    READ  REVISED DECAY  RATE  KIII.KK)  FOR  ALL  SEGMENTS
    IFdI.GT.O.AND.KK.GT.OJGO  TO  21
    WRITE(OUT,5)
  5  FOfcMAT(«   PROGRAM TERMINATED  DURING  SYSTEM  SENSITIVITY ANALYSIS
   1DUE TO ERRORS  IN REACTION  RATES  DATA INPUT'  )
    CALL EXIT
 21  READdNPUT,3)(K(L),L=l,NS)
  3  FORMAT(SFIO.O)
    TEMPERATURE  CORRECT  NEW  DECAY  CONSTANTS
    IF(THETA.LE.O)GO  TO  82
    DO  83   L=ltNS
 83  K(L)=K-20)
 82  WRITE(OUT,4)II,KK,THETA,(K
-------
                                    -32-
13

16
>5
ro
30
12
39
30
99
   NOTE   —   PROGRAM  ASSUMES  CHANGES  ARE  MADE FROM  *** ORIGINAL ***
             SYSTEM  MATRIX,  IF  CUMMULATIVE  EFFECTS  ARE DESIRED
             SAVE  THIS  NEW  SYSTEM  MATRIX  ON OFFLINE DISK  FILE 1
   IROW=IROW+1
   ICOL=ICOL-H
   CONTINUE
   CONTINUE
   NRUNS=NRUNS-1
   IF(IPRINT)90,90,70
   WRITEtOUT.80)
   FORMAT(1H1,///,43X,'REVISED SYSTEM MATRIX BEFORE SOLUTION*,///)
   WRITE(OUT,89)(I,1=1,NP)
   FORMAT(/,12X,10112,/)
   DO 39 IR=1,NP
   WRITE(OUT,12)IR,(AUR,JC),JC=1,NP)
   FORMAT(IX,•ROW NO.•,13,IX,10F12.2/<12X,10F12.2)>
   CONTINUE
   RETURN
   CALL EXIT
   END
   SUBROUTINE DODEF
* SENS155
* SENS156
* SENS157
* SENS158
  SENS159
  SENS160
  SENS161
  SENS162
  SENS163
  SENS164
  SENS165
  SENS166
  SENS167
  SENS168
  SENS169
  SENS170
  SENS171
  SENS172
  SENS173
  SENS174
  SENS175
  SENS176
 OODEF001
*DODEF002
   THIS SUBROUTINE WILL COMPUTE DISSOLVED OXYGEN  DEFICIT
   BY SELECTING COMPONENTS CAUSING DEFICIT WITHIN THE  SYSTEM
                                                                    *DODEF004
                                                                    *DODEF005
                                                                    *DODEF006
   THE OUTPUT WILL BE D.O. DEFICIT DUE TO EACH COMPONENT AND        *DODEF007
   TOTAL DEFICIT DUE TO ALL COMPONENTS IN EACH SECTION OF THE SYSTEM*DODEF008
10
   INTEGER OUT
   REAL KA
   DIMENSION STOCH{10»,ICOEF1(10),ICOEF2(101,BMAT(60,60),CS(60)
   DIMENSION VK(60),RES(60),PROD(60),DODSUM(60),SROW(60),WDODEF(60)
   COMMON/UT1/IPIVO(180),PIVOT(180)
   COMMON/UT2/INDEX(180,2),KA(60)
   COMMON/BLK1/NR,N,NP,INPUT,OUT
   COMMON/BLK3/V(60),TEMPERC60),CHLOR160)
   COMMON/BLK4/PTMIL(60),XK1(180)
   COMMON/OPTION/IPRINT,ISENS,KWDEF
   EQUIVALENCE 
-------
                                    -33-
                                                                    *DODEF038
   IF KA ARE SAME IN ALL SECTIONS, INPUT ONLY VALUE OF FIRST SECTION*OODEF039
                                                                    *DODEF040
******** ************#*#£*****$**$**$*$*$*****************************DODEF041
   IF(KA(2))2t2,3                                                    DODEF042
 2 DO 4  1=2,N                                                       DODEF043
   KA(I)=KA(1)                                                       DODEF044
 4 CONTINUE                                                          DODEF045
*********************************************************************QQDEF046
                                                                    *DODEF047
   TEMPERATURE CORRECT REAERATION RATES       .                      *DODEF048
*********************************************************************OODEF050
   DO 6  I=l,N                                                       DODEF051
   KA(I)=KA(I)*1.024**(TEMPER(I)-20.)                                DODEF052
 6 CONTINUE                                                          DODEF053
***********************************#*********************************OQDEF054
                                                                    *DODEF055
   PRINT REAERATION RATES                                           *DODEF056
                                                                    *DODEF049
                                                                    *DODEF057
*********************************************************************DQQEFO58
 3 WRITE(OUT,8)                                                      DOOEF059
 8 FORMATC1',40X,'DISSOLVED OXYGEN DEFICIT ANALYSIS',//,5X,         DODEF060
  l'SEGMENT',10X,'TEMPERATURE (CENT.»»,8X,'REAERATION RATE (I/DAY)' JDODEF061
   WRITE(OUT,9J{I,TEMPER(I),KA{I),I=1,N)                             DODEF062
 9 FORMAT(5X,I4,14X,F6.2,21X,F12.6)                                  DODEF063
*********************************************************************DQDEF064
                                                                    *DODEF065
   COMPUTE SATURATION LEVEL OF DISSOLVED OXYGEN IN EA. SECTION      *DODEF066
                                                                    *DQDEF067

   DO 7  1 = 1,N                                                       DODEF069
   CS(I)=(l.-.000009*CHLORtI))*<14.652--41022*TEMPERU)+             DODEF070
  1.0079910*TEMPER(IJ**2.-.000077774*TEMPER(I)**3.)                  DODEF071
 7 CONTINUE                                                          DODEF072
   COMPUTE DEFICIT  MATRIX                                           *DOOEF073
   DO 101 1=1,N                                                      DODEF074
   READ(1§ IMBMATdf J)v J-ltN)                                        DODEF075
   BMATCI,I)=BMAT(I,I)+V(I)*KA(I)                                    DODEF076
01 CONTINUE                                                          DODEF077

                                                                    *DODEF079
   PRINT DISS. OXYGEN DEF. MATRIX IF DESIRED                        *DODEF080
                                                                    *DODEF081
*********************************************************************DOOEP082
   IF(IPRINT)90,90,70                                                DODEF083
70 WRITECOUT,280)                                                    DODEF084
80 FORMATdHl,///,SOX,'DEFICIT MATRIX BEFORE INVERSION',///)         DODEF085
   WRITE(OUT,89)tI,I=l,N)                                            DODEF086
89 FORMAT(/,12X,10112,/)                                             DODEF087
   DO 39 IR=1,N                                                      DODEF088
   WRITE(OUT,162)IR,
-------
                                    -34-
 10 DC 900  1 = 1,N
   DC 900  J=l,N
 10 8MATU, J) = 8MATU,J)/1000000.
   DC 20 J=1,N
   IPIVO(J)=0
   DO 21 1=1,2
 >1 INDEX(J,I)=0
 !0 CONTINUE
   DO 550 1=1,N

   SEARCH FOR PIVOT ELEMENT
   AMAX=0.0
   00 105 J=1,N
   IF (IPIVO  (J)-l)
60, 105, 60
 iO DO 100 M=1,N
   IF (IPIVO  (M)-l)  80,  100, 740
 30 IF (ABS(AMAX)-ABS(BMAT(J.M)))85,  100,  100
 35 IROW=J
   ICOLU =M
   AHAX=BMATCJ,M)
 30 CONTINUE
 35 CONTINUE
   IPIVO UCOLU  )=IPIVO  (ICOLU  )+l

   INTERCHANGE ROWS  TO PUT  PIVOT  ELEMENT  ON DIAGONAL

   IF UROW-ICOLU  )  150,  260,  150
 50 DO 200 L=1,N
   SWAP=BMAT(IROW,L)
   BPAT(IROW,L)=BMAT(ICOLU,L)
 00 BMAT(ICOLU,L)=SWAP
 60 INDEX(I,1)=IRQW
   INDEX(I,2)=ICOLU
   PIVOT(I)=BMAT(ICOLU,ICOLU)

   DIVIDE PIVOT  ROW  BY PIVOT ELEMENT

   BMATlICOLU,ICOLUJ=1.0
   DO 350 L=1,N
 50 BMAT{ICOLU,L)=BMAT(ICOLU,L)/PIVOT(I)

   REDUCE NON-PIVOT  ROWS

   DC 550 L1=1,N
   IFCL1-ICOLU ) 400, 550,  400
 00 S=BMAT(L1,ICOLU)
   BMAT(L1,ICOLU)=0.0
   DO 450 L=1,N
 50 BMAT(L1,L)=BMAT(L1,L)-BMAT(ICOLU,L)*S
 50 CONTINUE

   INTERCHANGE COLUMNS

   DO 710 1=1,N
   L=N+1-I
   IF (INDEX(L,1)-INDEX(L,2))  630,  710, 630
i30 JROW=INDEX(L,1)
   JCOLU =INDEX(L,2)
 DQDEF097
 DODEF098
 DODEF099
 DODEF100
 DODEF101
 DODEF102
 DODEF103
 DODEF104
 DODEF105
*DODEF106
*DODEF107
*DODEF108
 DODEF109
 DODEF110
 DODEF111
 DODEF112
 DODEF113
 DODEF114
 DODEF115
 DODEF116
 DODEF117
 DODEF118
 DODEF119
 DODEF120
*DODEF121
*DODEF122
*DODEF123
 DODEF124
 DODEF125
 DODEF126
 DODEF127
 OODEF128
 DODEF129
 DODEF130
 DODEF131
*DODEF132
*DODEF133
*DODEF134
 DODEF135
 OODEF136
 DOOEF137
*DODEF138
*DODEF139
*OODEF140
 DODEF141
 DODEF142
 DODEF143
 DODEF144
 DODEF145
 DODEF146
 DODEF147
*DODEF148
*OQDEF149
*DODEF150
 DQDEF151
 DODEF152
 DODEF153
 DODEF154
 DODEF155

-------
                                    -35-
    00 705 M=1,N
    SWAP=BMAT(M,JROW)
    BNAT(M,JROW)=BMAT(M,JCOLU)
    B!«AT(M,JCOLU)=SWAP
705 CONTINUE
710 CONTINUE
740 CONTINUE
    DO 901 1=1,N
    DO 901 J=1,N
901 BMATU,J) = BMATU, JJ/1000000.
    IOODSW=999
    PRINT INVERSE OF DEFICIT MARTIX IF DESIRED
    IF(IPRINT)139,139,370
370 WRITE(OUT,281)
281 FORMAT(1H1,//,40X,«INVERSE OF DEFICIT MATRIX1,///}
    WRITE(OUT,89X1,1=1,N)
    DO 339 IR=1,N
    WRITE
-------
                                     -36-
   ICOL=(ICOEF1UI)-1)*N+1
   LAST=IROW+N-1
   ISTRT=1
   00 104 J=IROW,LAST
   READ < 3' J)(SROW(K),K=lfNP)
   VK(ISTRT) = ABS(SROW(ICOD)
   ISTRT=ISTRT-H
   ICQL=ICOL+1
.04 CONTINUE
   MULTIPLY   (VK)* *STOCH (II) *VK ( I )
    LOC=LOC+1
L55 CONTINUE
    MULTIPLY  CVK)*STOCH*(SOL) BY INVERSE OF DEFICIT MATRIX
    DO 106  1=1, N
    PRODI I)=0.0
    DO 106  J=1,N
    PROD ( I )=PRCD ( I KBMAT ( I , J )*RES i J )
106 CONTINUE
    PRINT RESULTING DEFICIT FOR REACTION
    WRITE (OUT, 12 ) ICOEF1 (II), ICOEF2( 1 1 ) ,STOCH( 1 1 J
 12 FORMAT!//, 20X, 'DEFICIT DUE TO REACTION* , 12, •  TO ',12,
   1 10X,fSTOCHIOMETRIC COEFFICIENT OF • ,F6.2, //,5X,
   1 'SECTION*, 5X,»MILE PTf,5X, 'DEFICIT CONCENTRATION (MG/L)',5X,
   1 fD.O. SATURATION {HG/L )', 5X, 'DISSOLVED OXYGEN (MG/L)r )
    DO 108  1=1, N
    DO=CS(IJ-PROD(I)
    WRITE(OUT,13)I,PTMIL! I ) ,PROD( I ),CS( I ),00
108 CONTINUE
 13 FORMAT(7X,I3,8X,F5.1,13X,F7.3,25X,F7.3,19X,F7.3)
    SUM REACTION DEFICIT FOR TOTAL DEFICIT COMPUTATION
    DC 107  1=1, N
    DODSUM(I)=DODSUM(I)+PRCDII )
107 CONTINUE
103 CONTINUE
    PRINT TOTAL DISSOLVED OXYGEN DEFICIT
                                                                      DODEF229
                                                                      DODEF230
                                                                      DODEF231
                                                                      DODEF232
                                                                      DODEF233
                                                                     'DODEF234
                                                                     *DODEF235
                                                                     *OODEF236
                                                                     *DODEF237
                                                                     ^DODEF238
                                                                      DODEF239
                                                                      DODEF240
                                                                      DODEF241
                                                                      DODEF242
                                                                      DODEF243
                                                                     'DODEF244
                                                                     *DODEF245
                                                                     *DODEF246
                                                                     *DODEF247
                                                                     *DQDEF248
                                                                      DODEF249
                                                                      DODEF250
                                                                      DODEF251
                                                                      DODEF252
                                                                      DODEF253
                                                                      DODEF254
                                                                      DODEF255
                                                                      DODEF256
                                                                      DODEF257
                                                                      DODEF258
                                                                     'DODEF259
                                                                     *DODEF260
                                                                     *DODEF261
                                                                     *DODEF262
                                                                     CDODEF263
                                                                      DODEF264
                                                                      DODEF265
                                                                      DODEF266
                                                                      DODEF267
                                                                     'DODEF268
                                                                     *DQDEF269
                                                                     *DODEF270
                                                                     *DOOEF271
   WRITE(OUT,14)(ICOEF1(I),ICOEF2(I),  I=1,NRDEF)
                                                                      DODEF273

-------
                                    -37-

    WRITE(OUT,15)                                                      DODEF274
 14  FORMATC1',40X,'TOTAL DISSOLVED OXYGEN  DEFICIT  FOR  REACTIONS't     DODEF275
   1 10(/,76X,I2,'  TO «,I2))                                           DODEF276
 15  FORMAT(//,5X,'SECTION*,5X,«MILE PTf,5X,                            DODEF277
   2'TOTAL DEFICIT  CONCENTRATION (MG/L)'t5X,'DISSOLVED  OXYGEN  (MG/L)•)DODEF278
    DO 109  1=1,N                                                      DODEF279
    00=CS(I)-DODSUM(I)                                                 DODEF280
    WRITE(OUT,13)I,PTMIL(I),DODSUM{I),DO                               DODEF281
109  CONTINUE                                                          DODEF282
    RETURN                                                            DODEF283
    END                                                               DODEF284

-------
         -38-
INPUT DESCRIPTION

-------
1.  TITLE CARD
                                          -39-
COLUMNS
1-80
VARIABLE
Title
FORMAT
20A4
SYSTEM SIZE CARD
COLUMNS
1- 2
3- 4
PROGRAM OPTIONS
variable.
COLUMNS
1- 2
3- 4
5- 6
7- 8
SYSTEM GEOMETRY
COLUMNS
1- 7
8-14
15-21
22-28
VARIABLE
N
NR
- If option
VARIABLE
IPRINT
ISENS
KBOD
KWDEF
CONTROL CARD
VARIABLE
SCALE (1)
SCALE (2)
SCALE O)
SCALE (4)
FORMAT
12
12
desired, enter
FORMAT
12
12
12
12

FORMAT
F7.0
F7.0
F7.0
F7.0
                                                       DESCRIPTION
                                                       Users description of run
                                                       DESCRIPTION

                                                       Number of system segments
                                                       Number of components
                                                              (reactants)

                                               a positive integer in corresponding
                                                       DESCRIPTION

                                                       Print system matrix and its
                                                       inverse
                                                       Perform sensitivity analysis
                                                       of system
                                                       Perform dissolved oxygen defi-
                                                       cit analysis
                                                       Control variable to read
                                                       direct (background)  sources of
                                                       deficit in the deficit anal-
                                                       ysis subroutine
                                                       DESCRIPTION

                                                       Area scale factor ft /ou*

                                                       Dispersion scale factor
                                                       mi /day/ou

                                                       Flow scale factor CFS/ou

                                                       Length scale factor ft/ou
     *note:  ou  stands  for  the  "original  units"  in which  the parameter  to  be  con-
            verted  by the scale  factor  is  expressed.  The purpose  of  the scale
            factors is to  allow the user to input parameters  in  units  which  are
            different  from those  specified  on the following pages.   For instance,
            according  to this  program length should be input  as  feet.  However,
            it  may be  more convenient to enter  it in miles and set  (SCALE(4)  to
            5280.  The program would  then internally convert  the length from
            miles to feet

-------
                                          -40-
5.' INTERFACE PARAMETER CARDS  (_2 per segment;  total number = 2*N)

                            Format    Description

                                      Cross sectional area of
                            F6.0      interface with, section
                                      Dispersion coefficient of
                            F6.0      first interface with section
                                      Flow- * across the first
                            F6.0      interface with section
                                      Section forming the first
                            13        interface with section
                                      Cross sectional area of
                            F6.0      second interface
                                      Dispersion coefficient of
                            F6.0      second interface
                                      Flow across the second
                            F6.0      interface
                                      Section forming the second
                            13        interface
                                      Cross sectional area of third
                            F6.0      interface
                                      Dispersion coefficient of
                            F6.0      third interface
                                      Flow across the third inter-
                            F6.0      face
                                      Section forming the third
                            13        interface

Note:  The second interface parameter card is identical to
       above, but for the fourth, fifth and sixth interfaces.
       If the segment in question has an interface which
       forms a boundary across which mass is transported
       (that is, a boundary condition) interface parameters
       must be input for it.  To do this, input the appropriate
       AREA, E and Q for the interface and input the segment number
       as the IARAY.

       It is only necessary to input the parameters for a
       particular interface once.  In other words after
       inputting the parameters of the interface of segment 1
       with segment 2, it is not necessary to input the
       parameters of the interface of segment 2 with segment 1.
INTERFACE
Columns
1- 6
7-12
13-18
19-21
22-27
28-33
34-39
40-42
43-48
49-54
55-60
61-63
PARAMETER (
Variable
AREA
E
Q
IARAY
AREA
E
Q
IARAY
AREA
E
Q
IARAY
                               Units
                               ft

                               mi /day
                               cfs
                               ft

                               mi /day

                               cfs
                               ft

                               mi /day
                               cfs
* flow out of a segment is positive;
  negative.
flow into a segment is

-------
                                         -41-
 6.
CHARACTERISTIC LENGTHS  (1 card per section;  1 length for each  interface;
therefore up to 6 lengths per card)
     Columns   Variable
1- 3

4-13

14-16

17-26

27-29

30-39

40-42

43-52

53-55

56-65

66-68

69-78

NLA(I,JI)

LA(I,JI)

NLA(I,JJ)

LA(I,JJ)

NLA(I,JK)

LA(I,JK)

NLA(I,JL)

LA(I,JL)

NLA(I,JM)

LA(I,JM)

NLA(I,JN)

LA(I,JN)


13

F10.0

13

F10.0

13

F10.0

13

F10.0

13

F10.0

13

F10.0
                        Format     Description
                                        First segment which forms an
                                        interface with segment I
                                        Length of segment I with
                                        respect to section JI*
                                        Second segment which forms an
                                        interface with segment I.
                                        Length of segment I with
                                        respect to segment JJ
                                        Third segment which forms
                                        an interface with segment I
                                        Length of segment I with
                                        respect to segment JK
                                        Fourth segment which forms
                                        an interface -with segment I
                                        Length of segment I with
                                        respect to segment JL
                                        Fifth segment which forms
                                        an interface with segment I
                                        Length of segment I with
                                        respect to segment JM
                                        Sixth segment which forms
                                        an interface with segment I
                                        Length of segment I with
                                        respect to segment JN
                                                                         Units
                                                                    feet
                                                                    feet
                                                                    feet
                                                                    feet
                                                                    feet
                                                                    feet
 * note:   the following figure illustrates  the  interpretation  of
          characteristic lengths  as  used in this program.
                        JN
                        JJ
                                    JK
                                JI
                                              JL
                                              JM
Figure 7:  Depiction of characteristic lengths

-------
                                        -42-
   Where
           @   - LA(I,JK) and LA(I,JI)

           (D   - LA(I,JJ), LA(I,JM) LA(I,JL) and LA(I,JN)

           As can be seen, any other than orthogonally shaped segments
   would present a problem as to the measurement of LA.  It is therefore
   reiterated that oddly shaped segmentation schemes should be avoided
   except when absolutely necessary to depict the natural geometry of
   the system.

           Note also that for sections which form boundaries across which
   mass is transported, the characteristic length of the segment with respect
   to the boundary must be inputted.  This length is read in as LA(I,I), where
   I is the boundary section.   (see equations 23 and 25).

7. SECTION CHLORIDE CONCENTRATIONS.  Optional.  This data is only necessary if
   KBOD>1  (see card type #3)

   Columns     Variable     Format
    1-10

   11-20
CHLOR(l)

CHLOR(2)
F10.0

F10.0
Description

Chloride concentration  (mg/1) at the
1st section
Chloride concentration at the second
section
   71-80       CHLOR(8)     F10.0      Chloride concentration (mg/1) at the
                                       8th section

   The card above is repeated as necessary until all chloride concentrations
   for the system have been entered.

8. TEMPERATURE/VOLUME CARDS

   Columns     Variable     Format     Description

    1-10       TEMPER(l)     F10.0      Water temperature (°C) of 1st section
   11-20       VOL(l)        F10.0      Volume (10  ft3) of 1st section
   21-30       TEMPER(2)     F10.0      Water temperature (°C) of 2nd section
   31-40       VOL(2)        F10.0      Volume (106 ft3) of 2nd section
   71-80
VOL(4)
F10.0
Volume of 4th section
   This card is repeated until the temperature and volume for all sections have
   been entered.   If the temperature is constant in the system, it is only
   necessary to enter this parameter for the 1st section, as the program will
   then test for the 2nd section temperature and since it is blank it will
   propagate the value of the 1st section to all the sections in the system.

-------
 9.
                                        -43-
                       REACTION RATE COEFFICIENTS
REACTION SCHEME CARD.

Consider the reactions:
                                  -f D
     These reactants are represented numerically in the program as follows:
     The reaction rates will have indices referring to these reactant numbers:
                              K
                          K14

     The program will read these reaction indices as follows:

               VARIABLE
     COLUMNS   NAME         FORMAT     DESCRIPTION
      1-2     I

      3-4     J

      6-10     THETA
                       FORMAT

                       12

                       12

                       F 5.3
           1st reaction index of reaction rate
           constant Kij
           2nd reaction index of reaction rate
           constant Kij
           temperature correction coefficient;  the
           program assumes the rate constants are
           given at 20°C and can be converted to
           any temperature by the equation:
           Kij(T°C) = Kij(20°C) * THETA**(TEMP-20).
           If this field is left blank, the program
           will not correct the reaction rates for
           temperature.
10.  REACTION RATE CONSTANTS (Units: I/day)
COLUMNS
1-10
11-20
21-30
VARIABLE
NAME
K(l)
K(2)
K(3)
FORMAT
F10.0
F10.0
F10.0
DESCRIPTION
reaction rate constant for 1st segment
reaction rate constant for 2nd segment
reaction rate constant for 3rd segment
     71-80
          K(8)
F10.0
reaction rate constant for 8th segment

-------
                                        -44-
     The card image above is repeated until the reaction rate constant Kij;
     where i, j are the  indices specified on card #9, have been read.  The
     reaction rates for a new reaction are then entered by preparing a new #9
     card with the appropriate indices for the new reaction, then reaction
     rate cards follow with the reaction rate constant for all segments.

     This procedure is repeated until the reaction rate coefficients have been
     read for all reaction schemes.   A blank card is inserted after the last
     set of reaction coefficients, which signals an end to this type of input.

11.  DIRECT DISCHARGES

               VARIABLE
     COLUMNS   NAME         FORMAT     DESCRIPTION

      1-10     LOAD         F10.0      Waste load into segment (Ib/day)
     11-12     ISEC         12         segment number where waste load is
                                       being discharged
     13-14     IREAC        12         component (number)  making up waste
                                       load

     These cards are repeated for all point sources into the system, the input
     is terminated by a blank card.

12.  BOUNDARY CONDITION

               VARIABLE
     COLUMNS   NAME         FORMAT     DESCRIPTION

      1-10     BC           F10.0      Boundary concentration (mg/1)
     11-12     ISYS         12         reaction component number
     13-14     ISECT        12         section number

     These cards are repeated for each reaction component at all boundary
     segments.  The last card will again be a blank card to terminate this type
     of input.

13.  MILEAGE POINT OF SEGMENTATION

               VARIABLE
     COLUMNS   NAME         FORMAT     DESCRIPTION

      1-10     PTMIL(4)      F10.0      Mile point designation for interface at
                                       end of first segment
     11-20     PTMIL(2)      F10.       mile point at end of second segment
     71-80     PTMIL(8)      F10.       mile point at end of eighth segment

-------
                                        -45-
     These cards are again repeated if necessary to input remaining mile points.
     The program will read the same number of mile points as the number of  seg-
     ments defined for the system on card image #2.

     This mile point is used only to locate the final concentrations in the system;
     it will not be used for any numerical computation such as the length of a
     segment.

                        DISSOLVED OXYGEN DEFICIT ANALYSIS  (KBODX)  ON OPTION CARD)

14.   SEGMENT REAERATION RATE

               VARIABLE
     COLUMNS   NAME         FORMAT     DESCRIPTION

      1-5     KA(1)         F5.3       Reaeration rate (I/day)  for first segment
      6-10     KA(2)         F5.3       Reaeration rate for second segment
     76-80     KA(16)     '   F5.3       Reaeration rate for sixteenth segment

     This card image is repeated if necessary to read reaeration rate constants
     for remaining segments.

15.  BACKGROUND DEFICIT DUE TO BACKGROUND SOURCES (KWDEF> 0 ON OPTION CARD)

                VARIABLE
     COLUMNS    NAME        FORMAT     DESCRIPTION

      1-10      WDODEF(l)    F10.0      Direct discharge of deficit (Ib/day) into
                                       the 1st segment
     11-20      WDODEF(2)    F10.0      Direct discharge of deficit (Ib/day) into
                                       the 2nd segment
     71-80      WDODEF(8)    F10.0      Direct discharge of deficit (Ib/day)  into
                                       the 8th segment

     The above card is repeated until these direct discharges have been defined
     for the whole system.

-------
                                        -46-
1-6.   SELECTING REACTION SCHEMES PRODUCING DEFICIT
     COLUMNS

      1- 5

     11-12

     13-14

     15-20
VARIABLE
NAME

NRDEF

ICOEFl(l)

ICOEF2(1)

STOCK(1)
FORMAT

15

12

12

F6.0
DESCRIPTION

Number of reactions producing deficit in
the system
1st index of reaction rate constant
                                     IJ
2nd index of reaction rate constant
Stochiometric coefficient used in converting
component concentration to equivalent
oxygen utilized
     The subscripting of the "variable names" in the above cards refers to the
     first reaction producing deficit, the same procedure is utilized for all
     other deficit producing reactions:
21-22
23-24
25-30
ICOEF1(2)
ICOEF2(2)
STOCK (2)
12
12
F6.0
                                       1st index of reaction rate constant K
                                       2nd index of reaction rate constant K.
                                       Stochiometric coefficient for second
                                       reaction
                                                                            IJ
                                                             'IJ
     71-72     ICOEFK7)    12

     73-74     ICOEF2(7)    12

     75-80     STOCH(7)     F6.0
                        1st index of reaction rate constant K
                                                             IJ
                        2nd index of reaction rate constant K

                        Stochiometric coefficient for seventh reaction
     NOTE:   The reaeration rates constants and the reaction schemes producing
     deficit are read only once in the program.  During consecutive runs; such
     as in  system sensitivity analysis, there is an internal switch in the pro-
     gram to bypass this input and also the computation of the system deficit
     matrix inverse after the initial use of the DODEF subroutine.

-------
                                        -47-
                         5YSTEM SENSITIVITY ANALYSIS

                          (ISENSX)  ON OPTION CARD)

17.  MAIN PROGRAM CONTROL VARIABLE

               VARIABLE
     COLUMNS   NAME          FORMAT     DESCRIPTION

      1-2      NRUNS         12         Number of sensitivity analysis runs to be
                                        performed

18.  SUBROUTINE CONTROL  VARIABLE

               VARIABLE
     COLUMNS   NAME          FORMAT     DESCRIPTION

      1-2      NRATES        12         Number of reaction rate schemes to be
                                        revised
      3-4      NLOADS        12         Number of loads in system source vector
                                        to be changed

19.  CHANGES TO WASTE LOADS  (NLOADS > 0)

               VARIABLE
     COLUMNS   NAME          FORMAT     DESCRIPTION

      1-10     LOAD          F10.0      New waste load  (Ib/day)  into  segment
     11-12     ISEG          12         Segment number  of  modified waste load
     13-14     IREACT        12         Component number of waste load

     This card image is  repeated until the number of new  waste loads specified
     on card #18' as NLOADS have  been  read.

           CHANGES TO REACTION RATE COEFFICIENTS (NRATES  > 0)

20.  REACTION RATE SCHEME TO BE  MODIFIED

               VARIABLE
     COLUMNS   NAME          FORMAT     DESCRIPTION

      1- 2     II            12         1st index of  reaction  rate Kn,  „,
                                                                   J.X f KK
      3- 4     KK            12         2nd index of  reaction  rate K
                                                                   X J. f JS-K.
      6-10     THETA         F5.3        Reaction  rate temperature correction
                                        coefficient

-------
21.  NEW REACTION RATE CONSTANTS
               VARIABLE
     COLUMNS   NAME
             FORMAT
                                      -48-
           DESCRIPTION
      1-10
     11-20
K(2)
F10.0      New reaction rate constant K
           for 1st segment.
F10.0      New reaction rate constant K
           segment
                                                                   11,KK (see above)
                                                                   11,KK for 2nd
     71-80
K(8)
F10.0      New reaction rate constant K
           segment.
                                                                   11,KK for 8th
     This card image is repeated until the reaction rate constant specified in
     card #20 as K       has been read in for all segments,
                  iX f KK                      	
     The card image #20 and the sequence of card image #20 is repeated until the
     number reaction schemes specified on card #19 (NRATES)  has been satisfied.

     A new sensitivity run (if NRUNS > I) would begin by reading a new card image
     #18 (see schematic diagram for various inputs).

-------
                                  -49-
                           KESTRICTIONS


1.)   FEDBAK03 as presently written is limited to a system of sixty sections

2.)   Each section can have a maximum of six interfaces.

3.)   Only one interface of a particular section may act as a boundary.

4.)   The maximum number of reactants is a variable, which when multiplied
     by the number of sections, it cannot exceed 180.

     Restrictions (1) and  (4) can be easily expanded by minor modifications
     to the computer program.  The application of this model to the nitro-
     gen cycle by assuming first order kinetics of the bacterial reactions
     is a major assumption which should be verified by laboratory studies.

-------
               -50-
APPLICATION TO THE DELAWARE ESTUARY

-------
                                   -51-
The area modelled in the Delaware Estuary is the 86 mile stretch of tidal
river between Trenton, New Jersey and Liston Point, Delaware.  Figure 8.
shows this area and the corresponding system segmentation which was
followed.

This area was extensively investigated during the 1960's by the Delaware
Estuary Comprehensive Study.  The purpose of this study was to determine
the water quality objectives for the estuary, the cost and benefits of
these objectives, and the treatment levels necessary to achieve these
objectives.  The survey data used in this test case was gathered by Hydro-
science Inc.     for the Delaware River Basin Commission.  The intent of this
section is to demonstrate the user of this model its application by re-verify-
ing the results published in the forementioned report.  Interested readers are
suggested to review the Hydroscience report for an in depth discussion on the
rational for the reaction coefficients, as well as the differences between
computed and observed values.

-------
  -52-
          DELAWARE ESTUARY
         COMPREHENSIVE  STUDY
            SECTIONS  FOR
        MATHEMATICAL MODEL
    U.  S. ENVIRONMENTAL PROTECTION AGENCY

   Region II          New York. New York
FIGURE

-------
                                                  TABLE 2

Source
Del River at Trenton
Trenton
Philadelphia N. E.
Camden
Philadelphia S. E.
Texaco
Philadelphia S. W.
Mobil
Repauno
C.D.C.A.
Chester
Sinclair
Sun Oil
Wilmington
Carneys Pt.
Chambers
POINT SOURCE
DRBC
River Mile
134.0
132.0
104.2
97.9
96.6
93.9
90.7
87.5
85.5
85.0
80.8
80.6
78.8
72.2
71.2
69.8
NITROGEN LOADINGS
Seg
No
1
1
9
13
14
15
16
17
17
17
18
18
19
21
21
22
Organic N
(#N/day)
9,000
740.
6,620.
2,280.
4,210.
v- 	
4,320.
.-• 	
1 — :-
1,000.
1,000.

	
8,250.
	
	
TO THE DELAWARE ESTUARY
Ammonia-N Nitrite N
(#N/day) (#N/day)
1,000. 	
7,020. 10.
8,930. 10.
1,870. 330.
9,230. 	
1,000. 	
13,780. 	
3,600. 	
18,000. 	
1,230. 40.
1,270. 	
2,000.
7,900. 	
• 5,100. 250.
	 	
	 	
Nitrate N Total N
(#N/day) (#/day)
16,000 26,000
40. 7,810.
110. 15,670.
700. 5,180.
13,440.
1,000.
	 18,100.
	 3,600. i
Ul
U)
	 18,000. '
50. 2,320.
2,270.
	 2,000.
	 7,099.
800. 14,400.
3.0,500. ) 30,500.
' >

In order to reflect nitrogen inputs from natural runoff and smaller tributaries, it was assumed that 1,0 mg/1
ammonia nitrogen were contributed along the length of the estuary with an associated flow of 20 CFS/miles,
This is equivalent to a uniform loading of 107 Ibs/day  mile"  as shown on the system waste vector of the computer
program.

-------
                   NUN'BCK OF SECTIONS =  30           NUMBER OF REACTION COMPONENTS =
**SC~ALE FACTORS**~~
      rr^nnrcrm	sraxrnnr^	rrc'c'c	rrsrirrTr^	tTvvv	STTSTIIA)  =—rcuo.oco

-------
	 1





























JJTERFACE
ROW-SEG
I- 1)
I- 0)
2- 3)
2- C)
._3- -4L
3- 0)
4- 5)
4- 0)
5- 6)
5- 0)
L 6- 7 )
6- 0)
7- 8)
7- 0)
8- 9)
8- 0)
9- 10)
9- 0)
10- 11)
1C- 0)
11- 12)
11- 0)
12- 13)
12- 0)
13- 14)
13- 0)
14- 15)
14- 0)
15- 16)
15- 0)
16- 17)
16- 0)
17- 18)
17- 0)
18- 19)
18- 0)
19- 20)

20- 21)
20- 0)
21- 22)
21- 0)
22- 23)
22- 0)
23- 24)
23- 0)
24- 25)
24- 0)
25- 26)
25- 0)
26- 27)
26- 0)
27- 28)
27- 0)
AREA

-------
INTERFACE
ROW-SEG
( 28- 29)
( 28- 0)
- L.2
( 30- 0)
AREA
(FT**2)
235000.
0.
0.
307400.
0.
E
6.900
0.0
7.000
0.0
7.200
0.0
0
CCFSI
5627.0
0.0
5697.0
0.0
5722.0
0.0
	 INTERFACE_
ROW-SEG
C 28- 27)
( 28- C)
( 29- 28)
( 29- 0)
{ 30- 20)
i 30- 0)
AREA
CFT**2»
214500.
0.
235000.
0.
254500.
0.
E Q
JMI**2/D> «CFS)
6.800 -5624.0
0.0 0.0
6.900 -5627.0
0.0 0.0
7.000 -5697.0
0.0 0.0
INTERFACE 	
ROW-SEG
28- 0)
28-
29-
29-
30-
30-
0)
0)
0)
0)
0)
AREA
(FT**2)
0.
0.
0.
0.
0.
0.
E
Q
(MI**2/D) tCFS)
0.0 Q.O
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
 I
U1

-------
        SECTION
   	SECTION
        SECTION
   	SECTION.
        SECTION
   . _.	SECTION
        SECTION
   	SECTION-
        SECTION
   	 SECTION
        SECTION
   	SECTION
        'SECTION
   	SECTION.,
        SECTION
   	SECTION.
        SECTION
   	 SECTION
        SECTION
   	SECTION_
        S'ECTION
   	SECTION.
        SECTION
   	SECTION.
        SECTION
   	S.ECTIJQN_
        SECTION
   	SECTION
        S'ECTION
  1 HAS  AN  EXCESS FLOW
  _2__H A S _AN_£XCESS FLOW
  3 HAS  AN  EXCESS FLOW
                          Oh 0.303000^ 03 CFS
                             0.126000E 03 CFS
                          OF 0.108000E 03 CFS
	4_HA$. AN. EXCESS. FLOW	Of _Q.999983£ 01_CFS_
  5 HAS AN  EXCESS  FLOW   OF 6.330000E 02 CFS
  6_HAS. AH. EXCESS^ FLOW	OF..D.9000QOE.D2 CF$._
  7 HAS AN  EXCESS  FLOW   OF -.130000E 03 CFS
  9 .H6_S_JN  EXCESS  FLOW   OF QT?sono3FL_PZ_rFS_
  9 HAS AN  EXCESS  FLOW   OF 0.319997E 02 CFS
..10_HAS__AN_EXCESS__E.LQ.H	D.F__Qt200000E_.D2_CFS_
 11 HAS AN  EXCESS  FLOW   OF 0.230000E 03 CFS
 12 HAS AN  EXCESS  FLOW   OF 0.459999E 02 CFS
                          OF O.A99992E 01 CFS
                         J)F 0.100002E 02 CFS
                          OF 0.649996E 02 CFS
                         _QF__0.7090POE_03 CFS_
                          OF 0.369000E 03 CFS
                         _QF__jO. 1130QOE 03 CFS_
                          OF 0.749999E 02 CFS
                         J)F 0.320001F 02 CFS
 13 HAS AN  EXCESS  FLOW
_14_HAS_AN_EXCESS  FLOW
 15 HAS AN  EXCESS  FLOW
 _16 HAS AN  exCESS_FigW_
 17 HAS AN  EXCESS  FLOW
 18 HASSAN  EXCESS  FLOW
 19 HAS AN  EXCESS  FLOW
 20 HAS AN  EXCESS  FLOW
 21 HAS AN  EXCESS  FLOW
_2 2 H A S_ A N  E X C E S S_£L QW_
                          OF 0.110001E 02 CFS
                          OF 0.3190COE 03 CFS
 23 HAS AN  EXCESS  FLOW .  OF 0.199981E 01 CFS
. 2 -V _H AS_ Aji . E X C E S S__F L OW	OF 0.20 0019E_0 l_C£S_
 25 HAS AN  EXCESS  FLOW   OF 0.299972E 01 CFS
_2_6_HAS AN  EXCESS  FLOW   OF 0.600020E 01 CFS
 27 HAS AN  EXCESS  FLOW
..2_8_HAS_AJ\LEXCES5_ELOiL
 29 HAS AN  EXCESS  FLOW
 30 HAS AN  EXCESS  FLOW
                          OF 0.999983E 01 CFS
                          OF 0.300010E 01 CFS
                          OF 0.699999E 02 CFS
                          OF 0.250000E 02 CFS
                                                                                                                                              01
                                                                                                                                              -j
10

9

e

-------
INTERFACE
2 -
3 -
4 -
5 -
6 -
7 -
8 -
9 -
10 -
11 -
12 -
13 -
15 -
16 -
17 -
18 -
19 -
20 -
21 -
22 -
23 -
24 -
25 -
! 26 -
: 27 -
28 -
29 -
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
LENGTH
21000.
20000.
20000.
20000.
20000.
20000.
10000.
10000.
10000.
10000.
10000.
10000.
10000.
10000.
20000.
20000.
20000.
20000.
20000.
20000.
10000.
10000.
10COO.
10000.
10000.
10000.
10000.
10000.
20000.
INTERFACE
1 - 2
2 -
3 -
4 -
5 -
6 -
7 -
8 -
9 -
10 -
11 -
12 -
13 -
14 -
15 -
16 -
17 -
18 -
19 -
20 -
21 -
22 -
23 -
24 -
25 -
26 -
27 -
28 -
29 -
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
LENGTH
21000.
20000.
20000.
20000.
20000.
20000.
10000.
10000.
10000.
10000.
10000.
10000.
10000.
10000.
20000.
20000.
20000.
20000.
20000.
20000.
10000.
10000.
10000.
10000.
10000.
10000.
10000.
10000.
20000.
INTERFACE LENGTH
1-0 0.
2 -
3 -
4 -
5 -
6 -
7 -
8 -
9 -
10 -
11 -
12 -
13 -
14 -
15 -
16 -
17 -
18 -
19 -
20 -
21 -
22 -
23 -
24 -
25 -
26 -
27 -
28 -
29 -
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
o;
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
INTERFACE LENGTH
1-0 0.
2 -
3 -
4 -
5 -
6 -
7 -
8 -
. 9 -
10 -
11 -
12 -
13 -
14 -
15 -
16 -
17 -
18 -
19 -
20 -
21 -
22 -
23 -
24 -
25 -
26 -
27 -
28 -
29 -
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
INTERFACE LENGTH
1-0 0.
3 -
4 -
5 -
6 -
7 -
8 -
9 -
10 -
11 -
12 -
13 -
14 -
15 -
16 -
17 -
18 -
19 -
20 -
21 -
22 -
23 -
24 -
25 -
26 -
27 -
28 -
29 -
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
INTERFACE LENGTH
1-0 0,
2 -
4 -
5 -
6 -
7 -
8 -
9 -
10 -
11 -
12 -
13 -
14 -
15 -
16 -
17 -
18 -
19 -
20 -
21 -
22 -
23 -
24 -
25 -
26 -
27 -
28 -
29 -
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0. i
0 w
- 00
0. i
0.
0.
30 -  30
20000.
30 -  29
20000.
                                          30 -
                                                    30 -
                                                        0.
                                                   30 -
0.
                                                                                                        30 -
                                                                                                            0.

-------
 INTRFC EPRIH  ALPHA
	(MGD)  _   _
  l-
          26. 0.993
INTRFC EPRIM  ALPHA
     	(_MGD)_
                        1-
                                64.  0.985
 INTRFC EPRIM ALPHA
	(HGP)	
  I-   0
                                                        0.  0.0
 INTRFC EPRIM ALPHA
	LMG_D)	
  1-  0    6. 0.0
 INTRFC EPRIM ALPHA
	IMGD)	
           o.
                                                                                          i-
                                                                                                       0.0
 INTRFC EPRIM ALPHA
	(_MGQ»	
           0. 0.6
                                                                                                                 1-
2-
3-
4-
5-
6-
7-
8-
9-
10-
11-
12-
13-
14-
15-
16-
17-
18-
19-
20-
21-
22-
23-
*- 24-
'- 25-

27-
28-
29-
30-
3 89.
4 103.
5 119.
6 142.
7 230.
8 931.
-9 1709.
10 2542.
11 3429.
12 4474.
13 5231.
14 6671.
15 5525.
16 4701.
17 5251.
18 5913.
19 6900.
20 8190.
2112129.
2219008.
2321064.
2*23230.
2524253.
2625540.

2830420.
2922545.
3018577.
3023080.
0.980
0.978
0.974
0.969
0.951
0.797
0.500
0.500
0.500
0.500
0.500
0.500
0.667
0.500
0.500
0.500
0.500
0.500
0.333
0.500
0.500
0.5CO
0.500
0.500

0.500
0.66/
0.500
0.500
2-
3-
4-
5-
6-
7-
8-
9-
10-
11-
12-
13-
14-
15-
16-
17-
18-
19-
20-
21-
22-
23-
24-
25-

27-
28-
29-
30-
1 64.
2 89.
3 103.
4 119.
5 142.
6 230.
7 931.
8 1709.
9 2542.
10 3429.
11 4474.
12 5231.
13 6671.
14 5525.
15 4701.
16 5251.
17 5913.
18 6980.
19 8190.
2012129.
2119008.
2221064.
2323230.
2424253.

2627034.
2730420.
2822545.
2918577.
0.985
0.980
0.978
0.974
0.969
0.951
0.797
0.500
0.500
0.500
0.500
0.500
0.500
0.667
0.500
0.500
0.500
0.500
0.500
0.333
0.500
0.500
0.500
0.500

0.500
0.500
0.667
0.500
2-
3-
4-
5-
6-
7-
8-
9-
10-
11-
12-
13-
14-
15-
16-
17-
18-
19-
20-
21-
22-
23-
24-
25-

27-
28-
29-
30-
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

0
0
0
0
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.

0.
0.
0.
0.
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0-
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
Or\
0.0
0.0
0.0
0.0
2-
3-
4-
5-
6-
• 7-
8-
9-
10-
11-
12-
13-
14-
15-
16-
17-
18-
19-
20-
21-
22-
23-
24-
25-

27-
28-
29-
30-
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

0
0
0
0
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.

0.
0.
0.
0.
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
o.b
0.0
0.0 .

0.0
0.0
0.0
0.0
2-
3-
4-
5-
6-
7-
8-
9-
10-
11-
12-
13-
14-
15-
16-
17-
18-
19-
20-
21-
22-
23-
24-
25-
7 A —
27-
28-
29-
30-
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

0
0
0
0
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.

0.
0.
0.
0.
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0 '
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
Cft
0.0
0.0
0.0
0.0
2-
3-
4-
5-
6-
7-
8-
9-
10-
11-
12-
13-
14-
15-
16-
17-
18-
19-
20-
21-
22-
23-
24-
25-
"7 A—
27-
28-
29-
30-
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

0
0
0
0
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.

0.
0.
0.
0.
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0 |
0.0 Ln
0.0 f '
0.0
On
0.0
0.0
0.0
0.0

-------
SEC r i or.1
1
2
3
•'«
5
6
7
8
9
10
11
I?.
• 13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
VCLUf-E (KGAL)
If.' 10. 16
2122. 72
3 4 4 0 . H 0
3979. 36
4757.28
5654.88
34C3.40
3769. 92
39K6 .84
4353.36
4712.40
4899.40
5191.12
6C21.4C
13912.80
15164.40
16336.32
17922. C8
20136.16
21931.36
113C9.76
1 1773.52
12701.04
13404.16
13838.00
14391.52
15363.92
1-6815.04
36t>22.07
42C37.59
TEMPERATURE (CENT.)
30. CO
3C.CC
28. CC
23. CC
28. CC
27. 5C
26.50
?6.CO
26.00
26. CC
26.00
27.00
27. CC
27.00
27. CC
26. 5C
26.50
26.50
26.50
26. CO
26.00
26. CO
26. CC
26. CO
26. CC
26. CC
26. CC
26. CO
26.00
26. CO
CHLORIDE CCNC IMG/LI
14.
14 .
15.
15.
15.
15.
16.
17.
19.
20.
21.
22.
25.
40.
51.
43.
ICO.
115.
235.
4CO.
650.
663.
1350.
1420.
1660.
1900.
.2180. |
2400. g
2860. i
3350.

-------
SECTION NUKBE
KEALIIGN THETA 1 2 ~' " 3 	 4~
14 15 16 17
27 28 29 30
1 TO 2 0.0 0.10CO 0.1000 O.ICOO 0.10CO
O.ICOO O.ICOO 0.1000 O.ICOO
O.ICOO C.1COO O.ICOO O.ICOC
2 TU 3 1.C80 C.2375 0.2375 0.2036 C.2036
C.0017 C.C017 O.C016 O.C825
	 " 0.1746 0.1746 O;i746"- 0.1746
4TC4 0.0 0.0 0.0 0.0 C.C
O.C500 0.0500 C.O " 0.0
C.C O.C5CO O.C500 O.C500
3 TC 4 0.0 0.3000 0.3COO 0.3000 '0.3000
O.C020 O.C020 C.COi:0 0.1100
C.3COO 0.3CCO ' C.3CCO" 0.3COC
2 TO 2 0.0 C.12CC 0.120C 0.12CO C.1200
C.12CC 'C.12CO 0.12CO 0.1200""
0.1200 0.1200 0.1200 0.1200
3 1C 3 U.O 0.3000 " 0.3COC' 0.3COO 0.30CO
0.3COO 0.3COO 0.3CCO 0.3COO
<;•
18
O.ICOO
C.1COO
0.2036
0.0825
O.C
C.C
""0.1000
0.1100
0.12CO
0. JUOO
C.3COO
R
6
19
0.1000
0.1CCO
" 0.1959
0.1S14
0.0
O.C
0.3000
0.3COO


C
0.
0
0.
0
0
0.
0.1200 0
TJTTzrro o~:
0. iUOO
0.3CCO
u
0.

r
20
.1CCC
1000
.1814
1746
.0
0
.3000
3COO
.12CO
120^0
.3000
3000


0
0.
0
0.
0
C.
0
0.
0

8
21
.1COO
1COO
.1746
1746
.0
Cf 	
.3COO
3COO
.1200
(T.12CO
0
0.
.3COO
3COO

9
22
0.1CCC 0
O.ICOO 0.
••'0:1746 0
0.1746 0.
0.0 0
O.C 0.
0.3CCO 0
0.3COC 0.
0.12CC 0
0~.T2CO 0.
"O.'JUC'O 0
0.3COO 0.

'10 11
23 24
.1000 0.1000
1COO 0.1000
.1746 0.1746
1746 0.1746
.0 0.0
0 0.0
.3000 C.3UCO
3CCO 0.3000
.1200 0.12CO
1200 0.12CC
.3000 0.30UO
3000 0.3COO

25
0.1CCO
0.10CO
0.1865"
0.1746
0.0
C.O
0,3000
0.3CCO
C.12CC
"0 . I2t C
0.3CCU
0.3CCO



0
C

0
C

0
C

13
26
0.1CCO
.1CCC
C.CC17
. 1746
O.C5CC
.0
C.CC2C
.3CCO
0.12CC
. 12CC
0. JCLU
.3CCO

-------

LCAC
1C165.0
405.0
405. 0
4C5.0
4C5.0
405.0
2C2.0
202. 0
6822.0
202.0
2C2.0
202-0
2482.0
4412.0
405.0
4725.0
1405.0
14C5.0
4C5.0
4C5.0
8452.0
202.0
202.0
2C2.0
202.0
202.0
2C2.0
202.0
405.0
405.0
8445.0
405.0
405.0
4C5.0
4C5.0
4C5.0
202.0
2C2.0
9132.0
202.0
202.0
202. 0
2C72.0
9432.0
1405.0
	 	 141H5.0
23235.0
3675.0
8305.0
. . .. .. .... 4C5.0
5302.0
	 ' ' " " ' " 2C2.0
202.0
2 C 2 . 0
202.0
- - ' 	 ------ 202.0
202.0
JI SCi(AȣES
ScCTTON
1
__ 2 -
3
4
5
6
7
•8
9
10
11
12
13
14
15
16
17
ie
19
20
21
22
23
25
26
27
28
29
30
1
2
3
4
5
6
7
e
9
10
11
12
13
14
15
16
71
18
19
20
21
22
23
24
25
26
27
ILB/OflY)
• 1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1 ,
1 0^
l "y
i
i
i
i
2
2
2
2
2
2
2
2
2
2
2
t
2
2
2
2
2
2
2
2
2
2 " 	
2
2 	
2
2
2

-------
  _2C2.0	28	
   405.0      29          2
   4C5.0      30          2
    1C". 0    "   I          3
  _1C .0	9	'   3
   33Cf.(f      13          3"
    -40.0       17          3
   25C.O       21          3
 1604C.O        1          4
  "TTCTS	T
   7CO.O       13
	^CT(T      IT
 313CO.O       21

-------
                  SYSTEM  eCUNLTARY  CONCENTRATIONS
SECTION              COMPONENT              BOUND.  CCNC.  (MG/L)

-------
CCNPCNtNT
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
SECTION
1
2
3
4
5
6
7
8
9
1C
11
12
13
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
3C
1
2
3
4
5
6
7
8
9
10
11
12
13
15
16 -
17
18
19
20
21
22
23
24
25
26
27
28
LOAD (Lfi/DAY)
1C 165.0
40'j.O
405.0
405.0
405.0
202. C
202.0
6822.0
2C2.C
202.0
202.0
2482.0
4412.0
405. 0
4725.0
1405.0
1405.0
40^.0
405.0
3452.0
202.0
202. C
202.0
202. C • • &
202.0 ui
202.0 '
202.0
405.0
405.0
844b.O
405.0
405.0
405.0
405.0
405.0
2C2.C
202.0
9132.0
202. C
2C2.0
202.0
20/2.0
9432.0
1405.0
14185.0
0.0
3675.0
83Cb.O
405. 0
5302.0
2C2.C
202.0
202.0
202.0
202,0
202.0
202.0

-------
2
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
4
4
A
A
4
A
A
A
A
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
30
1
2
3
5
6
7
a
9
10
11
12
13
1*
15
16 '
17
18
19
2C
21
22
23
24
25
26
27
28
29
30
1
2
3
4
5
6
7
ft
9
1C
11
12
13
14
15
16
17
ie
19
20
21
22
23
24
25
26
27
26
29
4C5.Q
10.0
O.C
. o.o
c.c
0.0
0.0
0.0
C.G
10.0
0.0
0.0
c.o
330.0
0.0
0.0
• ' o.o
40.0
0.0
0.0
0.0
250.0
0.0
C.C
O.C
0.0
0.0
c'.o ^
0.0 
0.0 '
0.0
16040. C
0.0 •
0.0
O.C
O.C
0.0
0.0
0.0
11C.C
0.0
23235.0
O.C "" "
, 7CO.O
0.0
O.C
O.C
50. C
O.C
0.0
3130C.O
o.o - - •
C.C
0 . C
O.C
~0.0
0.0
O.C
0.0

-------

-------
SECTION PILE PT
1
2
3
4
5
6
8
9
10
•11
12
13
14
15
16
1 7 "
18
19
20
21
22
23
24
25
26
27
28
29
30
125.51
125.72
121.93
118.14
1 14.35
110.56
108.67
1C6.78
104.89
103. CO
101.11
99.22
97.33
95.44
91.65
87.86
84.07
80.28
76.49
72.70
70.81
68.92
67.03
65.14
63.25
61.36
59.50
57.58
50 .00
1
0. 568
0.569
"13.573
0.592
0.608
0.615
0.659
C.721
2
0.0 2'
0.476
0.482
0.501
0.515"
0.522
0.562"
0.632
1.C07 C.979
1.021 0.970
1.C39
1.104
1.201
1.264
1.257
1.246
" "1.228
1.23C
1.217
1.202
1.185
1.113
1.047
0.975
0.894
0.803
0.705
0.604
"0.452
0.223
0.977
1.C45
" 1.U1
1.275
1 . 261
1.3C3
1.205
1.199
1.195
1.112
1.057"
'O.st!3
" 0.917"
0.848
0.693
" 07605"
0.516
"0:336"
0.190
3
0.077
0.155
0.206
0.251
7J72 87
0.305
	 07TZ5""
0.3^6
0.457
0.454
0.419
0.236 "
0.198
07T2§"~
0.114
0.240
0.346
0.519
C.563
075537"
0.532
" '0.501
0.476
0.395
O.H47
0.298
0.223
0.109
4
C.915
C.940
I.CC6
1.135
1.3C3
1.495
1.727 '
1.939
2.339 .
2.826
3.376
3.253
2.546 / ' ' '
2.705
2.473
2.518
2.896 •
3.554
4.410
5.C98
5.424
5.350
5.230
5.012
4.684
4.247 " cri
3.7C7 f
3.1C6
2.251
1.C82

-------
SEGMENT
1
2
3
4
5
6
7
8
9
1C
' 11
12
13
14
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
TEfPER/STURE (CENT.)
3C.CC
3C.OO
28.00
28.00
28. CO
27. 5C
26.50
26.00
26.00
26. CO
26. CO
27.00
27. UU
27.00
27. CO
26.50
26.50
26.50
26.50
26.00
26. CO
26.00
26.00
26.00
"• 26.00
• 26.00
26.00
26.00
26.00
26.00
REAeRATICN RATE (I/DAY)
0.226175
C .228175
C.2176C5
C.2176C5
0.217605
0.215040
C.21CCCO
C.2C7525
0.207525
C. 207525'
C.2C7S25
C.2125C5
C.2125C5
0.2125C5
" C.212SCS
0.21COOO
C.21CC06
C.21CCCO
C.21CCCO
C.2C7525
C. 207525
C.2C7525
C.2C7525
0.207525
0.2C7525 . . 1
0.207525 _2
C. 207525 I
0.207525
C. 207525
0.2C7525
2 TO  3
secTiCN
i
2
3
4
5
6
7
8
9
10
11
" ' " "" 12 ""
13
14
15
16"
17
"18 ~~
19
20
21
22
23
MILE PT DEFICIT CONCENTRATION
129.5 0.275
125.7
121.9
iie.i ~
1 14 .4
11C. 6
ice. 7
106.8
1C4.9
" 1C 3 . 0
1C1.1
S9.2
97.3
S b . 4
91.6
67.9
64 . 1
80.3
76.5
12. 1
7C.8
68.9
67.0
0.579
O.B05
1.C23
1.210
1.335
1.465
1.5B8
1.916
2.078
2.068
1.922
1 .396
1 .044
C.724
C.657
1.159
1.645
2.380
2,6^1
2.596
2.516
2.417
(MO/LI o.o. SATURATION c'/U
7.436
7.436
7.722
7.722
7.722
7. 796
7.944
8. 020
e.020
e.c2o
8.020
7.869
7.869
1 . dbB
7.fi67
/.942
7.938
7.937
7.929
/ . *i^l
7.974
7.97J
7.924
DISSOLVED OXYGEN (MG/L)
7.161
6.857
6.917
6.699
6.512
6.461
6.480
6.432 	 "
6. 104
b.942 	
5.952
"" "5.948 ' - 	
6.473
6.e^ 	
7.143
7 . 2 « 6
6.780
6.292 " " '
5.549
b.3 l~i • "• 	
5.379
5.457
5.507

-------
25
26
2 1
28
29
30
63. '3 	
6 I./,
59 . 5
57.6
53. S
50. 0
*! . 1 1 -i
1.913
"1.685
1 .450
1.C86
C.535
7,'JOl
7.864
7.848
7.815
7.779
•^ . i a a
5.971
6.178
6.398
6.728
7.244
            'DEFTim DIE  TC'TTETCTromnrD—4~
                         "TTCUHTTJpirrR IC^CUEFFICIENT OF    flTT
SECYKN
1
2
.3
4
5
6
7
a
9
10
11
12
13
	 14
15
16
17
" 	 18
19
20
21
22
23
25
26
27
28
29
129.5
125.7
121.9
118.1
110.6
108.7
106.8
104.9
1C J . 0
101. 1
S9.2
97. 3
95.4
91 .6
67.9
84. 1
" 80/3
76.5
" 72.7 ""
7C. 8
68.9
67.0
65.1
63.3
61.4
59.5
53.8
0.019
C.066~
0.128
C.2C1
0.277
0.344
C.4C3
C.449
0.516
C.538
0.343
0.251
C. 164
C.124
0.178
0.302
0.550
0.691'
0.725
0.726
0.715
0.689
C.648
0.594
0.529
C.45S
0. J47
TIGN 
-------
2 TC 3
3 TC 4

SECTION
1
2
3
4
5
6
7
8
S
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
MILE PT
129.5
125.7
121.9
118.1
114.4
110.6
ice. 7
1C6.8
1C4.9
1C3.0
101.1
99.2
97.3
9b.4
91 .6
87.9
84.1
8C.3
76.5
72.7
7C.8
6 8. -9
67.0
65.1
63.3
61.4
59.5
57.6
53.8
50.0
TOTAL DEFICIT CONCENTRATION I^G/L)
0.295
C.645
C.933
1.224
1.467
1.679
1.868
2.037
2 . 4 3 i
2.632
2.6C6
2.4CO
1.739
1.295
o.B&y
0.780
1.337
1.947
2.S30
3.312
3.320
3.242
3.132
2.972
2.761
2.507
2.214
I. 90S
1.433
0.707
DISSOLVED OXYGEN (MG/LJ
7.142
6.791
6.789
6.498
6.235
6.117
6.077
5.983
5. sea
5.388
5.413
5.470
6.130
6.573
6.978
7.162
6.601
5.990
4.999
4.681
4.654
4.732 ' ^j
4.792 7
4.947
5.140
5.377
5.649
5.940
6.382
7.073

-------
       2D
       00
       40
                                    0
                                                                  e
       30
     10
        10
       oo
          130
                           ?)
                           o
                                                      |  0
130      120
110
100
90
80
70
60
50
                                               e
                                               o
                                               n
                                                  O-JUy 30, 664
                                                  O-AUG K). 1964
                                                  • -AUC H. I9M
                                                  OWPC DATA
                                                  ALL LWS
                                                  FLDW= JOOO CFS
                                                  TEMP -=26 -ZT C
                                              19
                                              I
                           Q
                           O
                                                    a
                                                    e
                                            s    r^
                                             j_
                  110-       CO      90       80
                           DISTARCE-MILES.
                                   70
                                   60
                                    50
                                                                                  i
                                                                                  ^j
                                                                                  NJ
FIGURE 9:  Observed vs Computed  (Dashed Line)  Nutrient Profiles, August 1964

-------
     20
   o>
   6 1.0



   'f

     00
     5.0
       J30
      30
   to
   i
      1.0'
        130
         •   8
         El  0
                  i-
                                  Q
                                  •
120
no
100
   90
80
70
60
50
                                                     I
                                                     8
                               T
               o

               9
                        •

                        •

                        O
                                             5
                                3
                               f
                                                O
                                                o
                                                               o
                                                               9
                                                               o
                                                               a
                                           •
                                           o
 120
 110
100
    90      80
DISTANCE-MLES
         70
          60
         50
                                                                                           I
                                                                                          ~j
                                                                                          ui
FIGURE  10:  Observed vs Computed (Dashed Line)  Nutrient Profiles, August  1964

-------
                                     -74-
Th e following output pages are a second method of computing deficit profiles
for the Delaware Estuary, where the computation of dissolved oxygen deficit
is accomplished by defining a deficit "species" as part of the system.

In this case it is defined as the NS component.  Since ammonia and nitrite
(N2 and N3 correspondingly) are the species producing deficit, a reaction
scheme for these two species is defined.  The reaction rates for each of
these is the product of the stochiometric coefficient for the oxydation
reaction times the reaction rate constant of the species that consumes
oxygen.

It should be noted that this method computes total deficit and not the con-
tribution of each component.  For a large system, this technique significantly
increases the amount of computer time and main memory required, since it in-
creases the order of the various matrices by the number of segments in the
system.

The pages that follow begin with the "System Reaction Rate Constants", since
the preceding pages would be identical to the example shown previously.

-------


REACTION THETA
__1 TO

2 TO
4 TO

2 TO
3 TO

5 TO
3 TO

2 TO
CO
» 3 TO

2 0.0
S E C T I
1
14
27
0,1000
2
15
28
0.1000
0.1000 0.1000
0.1000 0.1000
3 1.080
4 0.0

5 1.080
5 0.0

5 0.0
4 0.0

2 0.0
3 0.0

0.2375
0.0017
0.1746
0.0
0.0500
0.0
0.8146
0.0058
0.5987
0.3420
0.0020
0.3420
0.2282
0.2125
0.2075
0.30CO
0.0020
0.3000
0.1200
0.1200
0.1200
0.3000
0.3000
0.3000
0.2375
0.0017
0.1746
0.0
0.0500
0.0500
0.8146
0.0058
0.5987
0.3420
0.0020
0.3420
0.2282
0.2125
0.2075
0.3000
0.0020
0.3000
0.1200
0.1200
0.1200
0.3000
0.3000
0.3000
0 N
3
16
29
0.1000
N U M
4
17
30
0.1000
0.1000 0.1000
0.1000 0.1000
0.2036
0.0016
0.1746
0.0
0.0
0.0500
0.6984
0.0056
0.5987
0.3420
0.0020
0.3420
0.2176
0.2100
0.2075
•0.3000
0.0020
0.3000
0.1200
0.1200
0.1200
0.3000
0.3000
0.3000
0.2036
0.0825
0.1746
0.0
0.0
0.0500
0.6984
0.2828
0.5987
0.3420
0.1250
0.3420
0.2176
0.2100
0.2075
0.3000
0.1100
0.3000
0.1200
0.1200
0.1200
0.3000
0.3000
0.3000
B E
5
IS
0.1000
0.1000
0.2036
0.0825
0.0
0.0
0.6984
0.2828
0.3420
0.1250
0.2176
0.2100
0.3000
0.1100
0.1200
0.1200
0.3000
0.3000
R
6
19
0.1000
0.1000
0.1959
0.1814
0..0
0.0
0.6720
C.6222
0.3420
0.3420
0.2150
0.2100
0.3000
0.3000
0.1200
0.1200
0.3000
0.3000

7
20
0.1000
0.1000
0.1814
0.1746
0.0
0.0
0.6222
0.5987
0.3420
0.3420
0.2100
0.2075
0.3000
0.3000
0.1200
0.1200
0.3000
0.3000

S
21
0.1000
0.1000
0.1746
0.1746
0.0
0.0
0.5987
0.5987
0.3420
0.3420
0.2075
0.2075
0.3000
0.3000
0.1200
0.1200
0.3000
0.3000

9
22
0.1000
0.1000
0.1746
0.1746
0.0
0.0
0.5987
0.5987
0.3420
0.3420
0.2075
0.2075
0.3000
0.3000
0.1200
0.1200
0.3000
0.3000

10
23
0.1000
0.1000
0.1746
0.1746
0.0
0.0
0.5987
0.5987
0.3420
0.3420
0.2075
0.2075
0.3000
0.3000
0.1200
0.1200
0.3000
0.3000

11 12 13
24 25 26
0.1000 0.1000 0.1000
0.1000 0.1000 0.1000
0.1746 0.1885 0.0017
0.1746 0.1746 0.1746
0.0 0.0 0.0500
0.0 0.0 0.0
0.5987 0.6466 0.0058
0.5987 0.5987 0.5987
0.3420 0.3420 0.0020
0.3420 0.3420 0.3420
0.2075 0.2125 0.2125
0.2075 0.2075 0.2075
0.3000 0.3000 0.0020
0.3000 0.3000 0.3000
0.1200 0.1200 0.1200
0.1200 0.1200 0.1200
l
0.3000 0.3000 0.3000 pJ
0.3000 0.3000 0.3000 I
10
 9
 d
 7
 6

-------
DIRECT   DISCHARGES    (LB/DAY)
LOAD SECTION
inifti.n i












tn










12
10
9
8
7
6
5
4
3
405.0
405.0
405.0
405.0
405.0
202.0
202.0
6822.0
202.0
202.0
202.0 '
2482.0
4412.0
405.0
4725.0
1405.0
1405.0
405.0
405.0
8452.0
202.0
202.0
202.0
202.0
202.0
202.0
202.0
405.0
405.0
8445.0
405.0
405.0
405.0
405.0
405.0
202.0
202.0
9132.0
202.0
202.0
202.0
2072.0
9432.0
1405.0
14185.0
23235.0
3675.0
8305.0
405.0
5302.0
202.0
202.0
202.0
202.0
202.0
202.0
2
3
4
5
6
7
8
' 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
71
18
19
20
21
22
23
24
25
26
27
COMPONENT
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
I
1
1
1
1
1
1
1
i ;
i
i
i
i
i
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2

-------
405.0
AO'S.O
10.0
10.0
330.0
40.0
250.0
16040.0
110.0
700.0
50.0
31300.0
29
30
1
9
13
17
21
I
9
13
17
' 21
2
2
3
3 ... _
3
3
4
4
4
4
4
10
 9_
 e
 7
 6
 5

-------
                                                                                       DUMD A BJC_C-
                                                            SECTION	COMPONENT	BOUND.  CONC.  1MG/L)
                                                                                                                                                                                                           I
                                                                                                                                                                                                          -J
                                                                                                                                                                                                         -£D
                                                                                                                                                                                                           I
12
I I
10
 9
 a
 7

-------














tn
t










12
1 1
0
9
B
7
e
5
4
)
COMPONENT
1
1
1
1
1
1
1
I
1
1
1
1
1
1
1
1
I
I
1
1
1
1
1
1
1
1
1
1
1
1
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
SECTION
1
2
3
4
5
6
7
a
9
10
11
12
13
14 '
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
LOAD (LB/OAY)
10165.0
405.0
405.0
405.0
405.0
405.0
202.0
202.0
6822.0
202.0
202.0
202.0
2482.0
4412.0
405.0
4725.0
1405.0
1405.0
405.0
405.0
8452.0
202.0
202.0
202.0
202.0
202.0
202.0
202.0
405.0
405.0
8445.0
405.0
405.0
405.0
405.0
405.0
202.0
202.0
9132.0
202.0
202.0
202.0
2072.0
9432.0
1405.0
14185.0
0.0
3675.0
8305.0
405.0
5302.0
202.0
202.0
202.0
202.0
202.0
202.0
202.0 .








*




l
ID
1
















-------














ft










12
1
0
3
a
7
6
3
4
J
2
3
3
1
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
4
A
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
4
zq
3O
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
A05.0
4O5.O
i-ri ft
0.0
0.0
0.0
0.0
0.0
0.0
0.0
10.0
0.0
0.0
0.0
330.0
0.0
0.0
0.0
40.0
0.0
0.0
0.0
250.0
0.0
0.0
0.0
0.0
0.0 |
0.0 oo
o.o V
0.0
0.0
16040.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
110.0
0.0
23235.0
0.0
700.0
0.0
0.0
0.0
50.0
0.0
0.0
0.0
31300.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0

-------
f
5
	 5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
5
2 5
2 5
3O
1
2
3
4
5
6
7
8
9
10
11
12
13
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30

o.o
o.o 	
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0 '
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0 (i:
o.o -i-
0.0 '
0.0
11
10
 9
a

-------
SECTION
1
2
3
4
5
6
7
8
9
10
11 '
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
= 27
29
30
fiIL£_£JL_
129.51
125.72
121.93
118.14
114.35
110.56
108.67
106.78
104.89
103.00
101.11
99.22
97.33
95.44
91.65
87.86
8^.07
80.28
76.49
72.70
70.81
68.92
67.03
65.14
63.25
61.36
59.50
57.58
53.79
50.00
C O M P O
I
0.568
0.569
0.573
0.592
0.608
0.615
0.659
0.721
1.007
1.021
1.039
1.104
1.202
1.264
1.257
1.246
1.228
1.230
1.217
1.202
1.185
1.113
1.047
0.975
0.894
0.803
0.705
0.604
0.452
0.223
N E N T
2
0.472
Q.476
0.482
_Q.*.501 	
0.515
0.522
0.562
_0 , 632 	
0.979
_0*97Q
0.977
1.045
1.161
1.275
1.261
1.303
1.205
1.199
1.195
1.112
1.057
0.983
0.917
0.843
0.774
0.693
0.605
0.518
0.386
0.190
CON
3
0.077
0,155
0.206
0.287
0.305
0.326
0.348
0.423
_0J,457_
0.454
0.419
0.286
0.198
0.125
O.J14
0.240
0.346
0.519
0.563
0.553
0.532
0.507
0.476
0.439
0.395
0.347
0.298
0.223
0.109
CENT
4
0.915
_Q.340
1.006
1.303~
1.495
1.727
1.939
2.339
2,8?6
3.376
3.253
2.946
—i±W5—
2.473
_J.518
2.896
3.554
4.410
5.098
5.424
5.350
5.230
5.012
4.684
4.247
3.707
3.106
2.251
1.082
RATION 
-------
                                 -83-
Execution at the WCC facility

The two versions of the program have been stored as core-image
modules on a partitioned dataset at the EPA Washington Computer
Center Facility  (COMNET).  Users of this facility may execute the
programs by using the JCL stream listed below.

//EPAIII   JOB  ('AAAA1,MIII,,'FEDBAK03',B,,20),TIME=(5,30)
//JOBLIB DD UNIT=3330-1,VOL=SER=USER66,DISP=OLD/
//  DSNAME=CN.EPAGAN.R2TW.WATER.MODELS
//SI  EXEC PGM=PPPPPPPP,REGION=300K
//GO.FT01F001 DD UNIT=SYSDA,DSNAME=+FILE1,DISP=(NEW DELETE),
//  SPACE=(720,(720,5)),DCB=(RECFM=FS,BLKSIZE=720,LRECL=720)
//GO.FT03F001 DD UNIT=SYSDA,DSNAME=+FILE3,OISP*(NEW,DELETE),
//  SPACE=(1280,(400,5)),DCB=(RECFM=FS,BLKSIZE=1280,LRECL=1280)
//GO.FT06F001 DD SYSOUT=A
//GO.FT07F001 DD SYSOUT=B
//GO.FT09F001 DD SYSOUT=A,DCB=BLKSIZE=120
//GO.FT08F001 DD *
       ** DATA CARDS GO HERE **
/*
       WHERE   III = USER'S INITIALS

               AAAA = ACCOUNT NUMBER

               PPPPPPPP=FEDBAK3A  FOR 120  SEGMENT MODEL

               PPPPPPPP=FEDBAK3B  FOR 180  SEGMENT MODEL

-------
                                     -84-
ACKNOWLEDGEMENTS

Some material presented in this report is based on extensions of the HAR03(4)
computer model, originally documented by Steve Chapra.

The typing was very professionally handled by Ms. Gloria Kasporwitz.

-------
                                   -85-
                  REFERENCFS FOR THE COMPUTER PROGRAM


(1)   Committee on Sanitary Engineering Research,  Solubility of Atmos-
     pheric Oxygen in Water, J.  Sanit. Eng.  Div., Am. Soc.  Civil Engrs.,
     SA 4,  41 (1960)

(2)   Nitrification of the Delaware Estuary,  prepared for the Delaware
     River  Basin Commission by Hydroscience, Inc., Westwood, N.J., June
     1969.

(3)   Thomann, R.V., Systems Analysis and Water Quality Management,
     Environmental Science Division, New York, 1971.

(4)   Chapra, S.C.  and Nossa G.A.  HAF.g3_ A.Computer Program for the Modelling
     of Water Quality Parameters  in .a. Steady State Multi dimensional
     Aquatic System., USEPA Region II, New York,  N.Y., October 1974.

-------
                              -86-
                       APPENDIX A
LISTING OF INPUT DATA FOR VERIFICATION OF DELAWARE  ESTUARY

-------
                                     -87-


 *** NITRIFICATION MODEL FOR DELAWARE ESTUARY *  1964 CONDITIONS *** W/ DEFICIT
3004
     1
 1000.  1.0     1.0  1000.
  6.5   .40 -3000.  1 15.8  .40  3303.   2

  21.4  .40  3429.  3

  24.6  .40  3537.  4

  28.5  .40  3547.  5

  34.1  .40  3580.  6

  41.4  .40  3670.  7

  49.6  .90  35AO.  8

  51,2  1.6  3565.  9

  55.4  2.2  3597. 10

  60.9  2.7  3617. 11

  65.0  3.3  3847. 12

  66.0  3.8  3893. 13

  72.7  4.4  3898. 14

 88.3   4.5  3908. 15

 98.0   4.6  3973. 16

104.9   4.8  4682. 17

113.4   5.0  5051. 18

126.3   5.3  5164. 19

142.8   5.5  5239. 20

150.4   5.8  5271. 21

151.9   6.0  5282. 22

162.9   6.2  5601. 23

176.8   6.3  5603. 24

181.7   6.4  5605. 25

188.4   6.5  5608. 26

196.4   6.6  5614. 27

214.5   6.8  5624. 28

235.0   6.9  5627. 29

-------
                                     -88-
254.5
7.0  5697. 30
307.4   7.2  5722. 30
1 21.
3 20.
4 20.
5 2C.
6 20.
7 20,
8 10.
9 10.
10 10.
11 1C.
12 10.
13 10.
14 1C.
15 1C.
16 20.
17 20.
ia 20.
19 20.
20 20.
21 20.
22 10.
23 10.
24 10.
25 10.
26 10.
27 10.
28 10.
29 10.
30 20.
30 20.
14.
19.
100.
1660.
30.
28.
26.
27.
26.5
26.
26.
26.
0102
.10
.10
.10
.10
0203 1.08
.11
.11
.05
.11
0^04
2 21.
1 20.
2 20.
3 20.
4 20.
5 20.
6 10.
7 10.
8 10.
9 10.
10 10.
11 10.
12 10.
13 10.
14 20.
15 20.
16 20.
17 20.
18 20.
19 20.
20 10.
21 10.
22 10.
23 10.
24 10,
25 10.
26 10.
27 10.
28 20.
29 20.
14.
20.
115.
1900.
242.
636.
533.
694.
2184.
1512.
1850.
4896.

.10
.10
.10
.10

.11
.11
.05
.11

15.
21.
235.
2180.
3 .
27.5
26.
27.
26.5
26.
26.
26.

.10
.10
.10
.10

.11
.11
.11
.11

15.
22.
400.
2400.
364.
756.
582.
805.
2396.
1574.
1924.
5620.

.10
.10
.10
.10

.11
.11
.11
.11

15.
25.
650.
2860.
28.
26.5
26.
27.
26.5
26.
26.


.10
.10
.10
.10

.11
.001
.11
.11

15.
40.
663.
3350.
460.
455.
630.
1860.
2692.
1698.
2054.


.10
.10
.10
.10

.11
.001
.11
.11

16.
51.
1350.

28.
26.
27.
26.5
26.
26.
26.


.10
.10
.10


.11
.001
.11


17.
43.
1420.

532.
504.
655.
2030
2932
1792
2248


.10
.10
.10


.11
.001
.11


    0.0
      0.0
0.0
0.0
                                            0.0
0.0
0.0
0.0

-------
-89-
0.0
0.0
0.0
3304
.30
.30
.110
.30
0202
.12
.12
.12
.12
0303
.30
.30
.30
.30
1C165.
405.
405.
405.
405.
405.
202.
202.
6822.
202.
202.
202.
2482.
4412.
405.
4725.
1405.
1405.
405.
405.
8452.
202.
202.
202.
202.
202.
202.
2G2.
405.
405.
8445.
405.
405.
405.
405.
405.
202.
202.
9132.
202.
0.0 0.0 o.O
o.o o.o o.o
0.0 0.0 .05

.30 .30 .30
.30 .30 .30
-110 .30 .30
.30 .30 .30

-12 .12 .12
.12 .12 .12
-12 .12 .12
.12 .12 .12

.30 .30 .30
.30 .30 .30
.30 .30 .30
.30 .30 .30
1 1
2 1
3 I
4 1
5 1
6 1
7 1
8 1
9 1
10 1
11 1
12 1
13 1
14 1
15 1
16 1
17 1
18 1
19 1
20 1
21 1
22 1
23 1
24 1
25 1
26 1
27 1
28 1
29 1
30 1
1 2
2 2
3 2
4 2
5 2
6 2
7 2
8 2
9 2
10 2
.05 .05 .05 0.0
0.0 0.0 0.0 0.0
.05 .05

.30 .30 .30 .30
.002 .002 .002 .002
.30 .30 .30 .30
.30 .30

.12 .12 .12 .12
•12 .12 .12 .12
.12 .12 .12 .12
-12 .12

.30 .30 .30 .30
.30 .30 .30 .30
.30 .30 .30 .30
.30 .30









































-------
                                      -90-
202.
202.
2072.
9432.
1405.
14185.
23235.
3675.
8305.
405.
5302.
202.
202.
202.
202.
202.
202.
202.
405.
405.
10.
10.
330.
40.
250.
16040,
110.
700.
50.
3130C,
129.51
104.89
84.07
63.25
.18
11 2
12 2
13 2
14 2
15 2
16 2
71 2
18 2
19 2
20 2
21 2
22 2
23 2
24 2
25 2
26 2
27 2
28 2
29 2
30 2
1 3
9 3
13 3
17 3
21 3
1 4
9 4
13 4
17 4
21 4
125.72
103.0
80.28
61.36

129.51
104.89
84.07
63.25
.18
125.72
103.0
80.28
61.36

121.93
101.11
76.49
59.50

118.14
99.22
72.70
57.58

114.35
97.33
70.81
53.79

110.56
95.44
68.92
50.0

108.67
91.65
67.03


106.78
87.86
65.14


    2      23   3.43 3 4 1.14






/*

-------