EPA-600/3-80-028b
                                              February 1980
      MODELING OF SIMULATED PHOTOCHEMICAL
         SMOG WITH KINETIC_MECHANISMS
    Volume 2.  CHEMK:   A Computer Modeling
         Scheme for Chemical Kinetics
                       by
                 G.  Z.  Whitten
                    H.  Hogo
      Systems Applications, Incorporated
              950 Northgate Drive
         San Rafael, California  94903
            Contract No. 68-02-2428
                Project Officer

                Marcia C. Dodge
  Atmospheric Chemistry and Physics Division
  Environmental Sciences Research Laboratory
Research Triangle Park, North Carolina   27711
  ENVIRONMENTAL SCIENCES RESEARCH LABORATORY
      OFFICE OF RESEARCH AND DEVELOPMENT
     U.S. ENVIRONMENTAL PROTECTION AGENCY
RESEARCH TRIANGLE PARK, NORTH CAROLINA   27711

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

-------
                                  ABSTRACT


     Mechanisms that describe the formation of photochemical  smog are devel-
oped using a computer modeling technique directed toward the  simulation of
data collected in two smog chambers:  an indoor chamber and a dual  outdoor
chamber.  The results of simulating 164 different experiemnts are presented
in Vol. I.  Individual compounds for which specific experiments were simu-
lated and mechanisms developed include the following:   formaldehyde, acet-
aldehyde, ethylene, propylene, butane, and toluene.  Experiments in both
chambers were simulated for all these compounds.   The  mechanisms reported
describe the decay of the precursor organic compound,  formation and decay of
secondary organics, conversion of nitrogen oxides, formation  of nitrates,
and the appearance and decay of ozone.  Special emphasis is given to the
chemistry of toluene.  Also included is a study of a generalized smog-based
or carbon-bond mechanism developed in a previous study.  Vol. II contains
the user's manual and coding for a chemical kinetics computer program, .CHEMK.

     This report was submitted to the U.S. Environmental Protection Agency
in fulfillment of Contract No. 68-02-2428 by Systems Applications,  Incorpor-
ated.  This report covers the period 23 August 1978 to 23 August 1979, and
work was completed as of 23 August 1979.
                                     iii

-------
                                 CONTENTS
Abstract	   iii
Figures	    vi
Tables	    vi

   1.  Introducti on	     1
   2.  CHEMK Input	     3
            Basic input for a new problem	     3
            Program options	    11

References	    22
Appendices

   A.  Program listing	   A-l
   B.  Input data—CHEMK format	   B-l
   C.  Input data—EPA format	   C-l
   D.  Output data	   D-l
   E.  Example use of variable photolysis and print options	   E-l

-------
                                 FIGURES

Number                                                             Page
  1     Data input flow chart  	   4
                                 TABLES
Number                                                             Page
  1     Input Card Description 	  16

-------
                                SECTION 1
                               INTRODUCTION


     Calculation of the time-dependent concentration profiles of a number
of reacting species in a complex reaction mechanism is central  to the
science of chemical kinetics.  Such calculations allow the validation of
kinetic schemes and reaction rate constants  through the comparison of
predicted concentrations with experimental data.  The need to rely on
approximation methods such as steady-state assumptions, quantum yields,
and a "rate-determining step" in making these calculations adds uncertainty
to predicted concentrations and could possibly conceal inaccuracies
in the kinetic mechanism if the experimental data were highly scattered.
These approximation methods are also time-consuming in that the user must
formulate them and evaluate their impact on  the predictions.

     Work by Whitten (1974) and Overton (1976) led to the development of
fast and efficient computer programs that make complex calculations common-
place.  These computer programs integrate the system of ordinary differen-
tial equations (DDEs) associated with a kinetic mechanism, given a set of
initial conditions.  The method used to integrate the system of ODEs was
developed by Gear (1971) and modified by Hindmarsh (1974).  Initial work
by Whitten (1974) led to the development of a computer program that could
handle 50 chemical species and 200 chemical  reactions.  That program,
CHEMK, used the modified Gear routines documented by Hindmarsh (1974).

     Since that time Spellmann and Hindmarsh (1975) have modified the Gear
routines to handle sparse matrices.  The Jacobian, or the matrix of partial
differential equations, associated with a given chemical kinetic mechanism
is usually sparse, meaning that it contains  more zero elements than non-
zero elements.  Thus, the use of a routine designed specifically for sparse
matrices would increase the efficiency of the computer program CHEMK in
terms of central processing (CP) time.  The work of Spellman and Hindmarsh

                                      1

-------
and our experience with the limitations of the first versions of CHEMK led
us to modify CHEMK.

     In this report, we present a new version of CHEMK that can handle
chemical kinetic mechanisms containing up to 89 species and up to 200
reactions.  The new version can also handle photolysis rate constants and
temperatures that vary with time, and reaction schemes with stoichiometric
coefficients.  The Gear routine used is that of Spellmann and Hindmarsh
referred to above.

     The version of CHEMK documented in this report requires approximately
44000 decimal words to load on a Control Data Corporation (CDC) 7600 com-
puter.  The efficiency of the new program has been increased about a factor
of 3 from the original versions of CHEMK.  Rather than present a detailed
description of the program and its various routines, the emphasis in this
report has been placed on providing a description of program inputs and
options available to the user.  A complete program listing is provided in
Appendix A.  Appendices B and C give the input data for an example run in
the CHEMK and the EPA format, respectively.  Appendix D gives the output
of that example run in EPA format, and Appendix E provides an example of
variable photolysis and print options.

-------
                                  SECTION 2

                                 CHEMK INPUT

     A brief description of each data input  card and  its argument list is
given below.  Both initial  program specification and  modified operation
are treated in this presentation.  The data  input flow is illustrated
schematically in Figure 1.

BASIC  INPUT  FOR A  NEW  PROBLEM

General Data  Card
     The arguments that can be entered on this card are:

     >  The last number of the reaction sequence.
     >  Choice of reaction format.
     >  Print option to print only instantaneous concentrations
        at specified intervals.
     >  Dilution factor.
     >  Number of time intervals in which the photolysis rate
        constants are changing.
     >  The time interval of each photolysis number.
     >  The number of photolysis reactions.
     >  The number of temperature values.
     >  The time interval of each temperature value.
     >  Option to plot the photolysis values as a function of time.

The last number of the reaction sequence should correspond to the number
on the last card in the reaction scheme, even though  it is not the highest
reaction number.  The next parameter (choice of reaction format) is  to
enable users of EPASIM (Overton, 1976) and users of the first versions of

-------
*«x
     i
 i  £

 is*
 < o w


IP!
i «•= =

    T
      i
                                                   t.
                                                   in
              °
                          -
                                *-
                               Is5
                                     a

                                     5
                                     u
                                     S
                 £  I
      g
                           i/i

                         512
                         U U
                         —• hU
                         IE G»
                         M W»


                         Is
                                    sili
                                             5
                                             3

-------
CHEMK to use the current version of CHEMK.  A more detailed discussion of
this option is presented in the section on kinetic mechanism cards.

     An option has been added to the CHEMK program that will control the
printing of the species' net rates of formation and the net rates of reac-
tion.  When this option is not used, only species concentrations are printed
at specified time intervals.

     To activate this option so that species net rates of formation and
net rates of reaction will be printed at the same specified time intervals,
the user enters any alphanumeric character anywhere in columns 11 through
14 of the General Data Card.

     The dilution factor represents a constant dilution rate for every
species such that:
                               dC.
where
          C. = concentration of species i
          D = dilution factor.

In many experiments all species are being diluted by a common first-order
factor.  If this factor is entered in the third field of the General Data
Card (the one that identifies the last reaction and on all subsequent
appearances of such a card), then this option will be activated.  If one
species (e.g., 02) is part of the dilution gas, a first order "flow out"
(e.g., 02> the only reactant with no products) reaction with a negative
rate constant equal to the dilution factor will cancel the dilution for
such a species.

-------
     The next three parameters (number of time intervals in which the
photolysis rate constants are changing, the time interval of each photoly-
sis number, and number of photolysis reactions) are used to input the con-
ditions of the time-varying photolysis rate constants.  If a blank or zero
is entered for the number of time intervals, then the time interval and
number of photolysis reactions should not be entered.  If the number of
time intervals is nonzero then the first card after the General Data Card
should contain the vector of reaction numbers that identify the reactions
with varying photolysis rate constants.  Up to seven numbers may be placed
on each card.

     The next card(s) will contain the initial photolysis value and the
values at the end of each time interval.  If, for example, the user inputs
5 as the number of time intervals and 60 (minutes) as the time interval,
then six photolysis values should be entered on the card (one for the
initial value and one for every 60 minutes thereafter).  Up to seven
photolysis values may be entered on one card.  The total  time in which the
photolysis rate constants are changing in the above example is 300 minutes.
If the simulation is carried for a longer period of time (e.g., 600 minutes)
then all photolysis rate constants are set to zero after 300 minutes.  The
maximum number of time intervals is 97 for each simulation.  If more values
are needed the MORE option may be used (see section on MORE option).

     If the number of temperature values is nonzero, the next set of cards
is the vector of temperature values (initial plus the end of each time inter-
val specified on the General Data Card).  The format for the temperature
values is similar to that for the photolysis values.  All temperature values
must be in units of degrees Kelvin (K).  All activation energies entered on
each reaction card must also be in units of degrees Kelvin.  The temperature
at each specified printout time is printed with the rest of the species
under the column TEMP.

-------
     If the user wishes to obtain a plot of the photolysis values as a
function of time, then a nonzero value is declared in columns 56 to 60.
A line printer plot is produced at the end of the species concentration
versus time profiles.  This option should be activated only when the plot
option is declared.

2.   Kinetic Mechanism Cards

     Each chemical reaction in the mechanism must be entered individually.
Entry requires (1) that the reaction be numbered, (2) that its components
be specified in alphanumeric form, and (3) that the rate constant at 298K be spec-
ified.  If the reaction is temperature-dependent, then the activation energy
(in degrees Kelvin) can also be specified.  A maximum of 200 reactions is
allowed by the program.

     The user has a choice of two different input formats for the reactions.
To accommodate users of EPASIM (Overton,  1976) CHEMK allows nonunity stoichi-
ometric coefficients for the products of the reactions.  The reaction input
format described in EPASIM and summarized in Table 1 is available for users
of EPASIM and the current version of CHEMK.  A maximum of three reactants and
three products (with stoichiometric coefficients) are allowed under the
EPASIM format.  To use the EPASIM format, the user must enter a nonzero value
in the second field of the General Data Card (see preceding subsection).

     Reaction identification numbers are not used under the EPASIM format
as described by Overton (1976).  Therefore, we introduced this parameter in
columns 17, 18, and 24 of the reaction card.  The first field represents the
reaction number from 1 through 99.  If the user has more than 99 reactions,
a modulus factor must be specified in column 24 to convert the reaction num-
ber to a number greater than 100.  For example, a reaction identification
number of 13 is given on the reaction card in columns 17 and 18 while column
24 is left blank.  A reaction identification number of 123 is given on the

-------
reaction card first by input of 23 in columns  17  and 18 and 1  in column 24.
For reaction 200,  columns 17 and 18 are left blank (representing zero)  and
a 2 is entered in  column 24.

     For users of earlier versions of CHEMK and new users a different reaction
input format is available.  A maximum of three reactants and four products
(without stoichiometric coefficients) is allowed  under this format.

     For the reacting species, any combinations of letters and numbers can
be used to identify the individual entries with the exception that the letter
M is exclusively reserved to designate the total  pressure and TEMP is reserved
for the temperature.

     To specify a "flow in" reaction, that is, the continuous addition of a
particular species, a reaction with only products can be used where the reac-
tion rate k is equal to the flow rate.  A reaction with a reactant but no
products implies a first-order "flow out" process.

     When preparing the data, one must ensure  that the rate constant is in
units consistent with the order of the reaction as specified in the reaction
mechanism.  Instead of input of a rate constant for each photolysis reaction,
a photolysis ratio is entered when variable photolysis is used (see the section
on the General Data Card).  This ratio is the  ratio of the actual photolysis
rate constant to the photolysis values entered on the cardfs) following
the General Data Card.  For example, if the user inputs photolysis values,
then each of the photolysis reactions will have a photolysis constant
determined at each instant  in time from the product of the  interpolated
photolysis value and the  ratio entered on the reaction card.  All photol-
 ysis  values must  be declared on the  card  after the  General  Data Card,
as  discussed earlier.   Therefore,  for  example, if the photolysis  values
represent.N02  photolysis  rate constants,  then  the N02 photolysis  reaction
identification  number must  be declared and  a  value  of  1.0 entered for
the rate  constant.   It  is also  wise  to keep in mind the  following points:
                                         8

-------
     >  H02 and HOO will  be treated as different species.
     >  2ND will be treated as a new species and not as two NOs
        unless the EPASIM format is used.
     >  The maximum number of different species is 89.

Title Card

     Each simulation is titled by entering the desired notation in columns 1
through 28 of this card.   Anything may be entered in the first four columns
of this card with the following exceptions:

     >  All blank:  This  card will terminate the program.
     >  CONT:  This is an option indicating that concentrations and
        species from a previous run will be used as initial con-
        ditions for a new run.
     >  MORE:  This choice indicates that new reactions or new photo-
        lysis and temperature values are to be added.
     >  PLOT:  This entry activates the plotting routine.

Output Specification Card

     The first two items on this card specify the number of reactants with
nonzero initial concentrations and the number of integration steps between
print cycles.  Normally,  the program suppresses the printing of the concen-
trations of all species as a function of the integration steps.  However,
the user may key the output to the frequency of the integration steps by
entering a negative number of steps between print-outs in the appropriate
slot on the data card.  Thus, if the user enters a value of -5, the pro-
gram will print the instantaneous concentrations every 5 iteration steps.
If the user only wants concentrations printed at specified intervals, then
this option should not be activated (i.e. leave this  slot blank).  To
obtain output at specific times, the user enters the  initial print time
and then either a fixed time increment or a fixed multiple to be used
between subsequent printings.  For example, a user may specify an initial

-------
print time of 10 minutes and set the time increment at 10 minutes.  In
this case, output will be printed at 10 minutes, 20 minutes, 30 min-
utes, and so on.  Or the user may specify an initial print time of 10 min-
utes and then set the print multiple at 2, in which case output is pro-
vided at 10 minutes, 20 minutes, 40 minutes, 80 minutes, and so on.

 Initial Reactant Identification Card(s)

      The vector containing the alphabetic designation of all reactants with
 a nonzero initial concentration is given on this card(s).

 Initial Reactant Concentration Card(s)

      The initial concentrations of all reactants, in appropriate units, are
 provided on this card in the same sequence as species are entered on the
 Initial Reactant Identification Card(s).

 Time Specification Card

      The simulation starting time and ending time are the first two arguments
 on this card.  Their time unit should be identical to that used to specify
 rate constants.  In addition, the absolute temperature and the error toler-
 ance for the integration step may be entered on this card.  Default values
                                    _2
 of these parameters are 298K and 10  , respectively.

      After each simulation the user may choose the following options:

      >  Start a new simulation with the same mechanism by first
         inserting a new title card.
      >  Terminate the program by inserting a blank card.
      >  Provide a new or changed reaction sequence and start a
         new calculation using the MORE option.
                                        10

-------
     -  Continue the calculation using the current results but
        with some parameter changed using the CONT option.
        Plot the computed output using the PLOT option.

All of these possibilities are discussed in the following section.

PROGRAM OPTIONS

More Than One Computation

     Rather than enter a terminating blank card, insert a new title card
followed by a new set of initial condition data.  The new data will provide
computational results using the previously specified reaction scheme.

To Continue with Previously Calculated Species

     If the new title card has CONT in its first four columns, the species
from the previous calculation and their concentrations will be used as initial
conditions for a problem in which either new reactions are entered or the
temperature or error tolerance is changed.  The next card must in this case
be the time specification card, having the starting and ending times, error
tolerance, and temperature.

New or Changed Reactions

     To introduce either new or changed reactions into the kinetic mechanism,
enter a card with the word MORE punched in columes 1 through 4 followed by
                                               i
the set of cards used for specifying the changes in the kinetic mechanism.
This set includes a first card similar to the General Data Card.  Any
parameter controlled on the original General Data Card can be changed
on this new version.  The difference between this card and the original
is that all blank fields will default to the parameters controlled on
the original card or to the present values changed by previous uses of
                                      11

-------
the MORE option.  For example the dilution factor may be the only change
desired.  To accomplish this a MORE card followed by a card with only the
new dilution factor in columns 21 through 30 would be used.  The next
card in this example might then be a CONT to continue the previous cal-
culation or a new title to start an entirely new calculation with the
same mechanism modified only by the new dilution factor.

     The most frequent utilization of this option is typically for addi-
tions or alterations to the original mechanism.  To alter any reaction a
new reaction card is entered with the same identification number as the
original.  The new card must contain the complete reaction in its new
form since it will replace the original.  New reactions must have identi-
fication numbers not previously used in the original mechanism.  In any
case the identification number of the last reaction card to be read in at
this time must appear on the new General Data Card in columns 1 to 3.

     A new variable photolysis vector, a new variable temperature vector,
or a new set of photolysis reactions controlled by the variable photolysis
vector can also be introduced using this option.  If a new photolysis reac-
tion is to be added to the original set that varies with the time dependent
photolysis vector, then the entire new set, including both the originals
and the new additions, must be read in.

Plot Option Card

     At the end of the computational cycle, a card with PLOT in its first
four columns must be entered.  This signals the program that the plot option
is to be used and plot data are  to follow.

Plot Title Card

     If a plot has been selected by entering the word PLOT on the above card,
then the user must supply a plot title card, which may have storage option

                                     12

-------
entry and an option to print the species  concentrations  at  certain  time steps.
The storage option entry determines  whether  the  plotted  output will  be stored
on a tape file for later access.  If the  output  is  to  be stored,  the user need
only insert any nonzero number in the option argument  space.

     The data is stored on a file titled  TAPE 7  in  unformatted form.  The
following parameters are written to  the file TAPE 7:

     >  NTIT—The plot title.
     >  NAME(L)--The Lth species being plotted.
     >  CLOW—The lowest concentration entered on the  plot
        specification card.
     >  CHIGH—The highest concentration  value entered on the
        plot specification card.
     >  TLOW—The lowest time value  entered  on the  plot
        specification card.
     >  THIGH--The highest time value entered on the plot
        specification card.
     >  NDAT--Number of observed data points.
     >  NPNT--Number of calculated data points.
     >  (TIME(J),J=1,NDAT)—The vector of time values  of the
        observed data.
     >  (DATA(J),J=1,NDAT)—The vector of observed
        concentration values.
     >  (SAVTIM(J),J=l,NPNT)--The vector  of  time values  of
        the calculated data.
     >  (SAVCON(L,J),J=1,NPNT)~The  vector of calculated concentrations.

The twelve parameters are written with a  single  WRITE  statement and should
be read in the same manner.

     If the user is interested in the actual concentrations of the  species
being plotted, then a nonzero value  entered  for  this option will  activate the
                                    13

-------
printing of the actual concentrations at the times closest to every one-
eightieth of the total simulation time.  Thus, if the total simulation time
is 400 minutes then calculated concentrations at 5-minute intervals are
printed.

Plot Specification Card

     The plotting routine can handle up to three species on one plot.  This
Is done by specifying the number of species to be plotted on the plot speci-
fication card.  The information entered on this card is:  The number of spe-
cies to be plotted, the concentration labels for the vertical axis (five labels
must be entered), the lowest and highest values of the concentrations, and
the lowest and highest values of the integration time.

Species  Information Card

     This card contains the name of the species to be plotted (the species
identification code), the number of data points to be plotted from experi-
mental data, and the plot output symbol.  If the number of experimental data
points is zero, the next card is not entered.  A species identification card
must be entered for each species to be plotted.

Experimental  Plot  Data Cards

        A unique feature of the PLOT program is its ability to overlay
experimental data over a computer concentration profile.  To do this, pairs
of time and concentration values must be entered for each species.  Four
pairs of these values appear on each card, and the number of entries corres-
ponds to the  value specified on the  previous data card.

End of Plot Data Card

     A blank card must be entered to signal the end of the plot data.  The
next card can start a new simulation.

                                     14

-------
End of Program Card

     A blank card entered instead of a new title or option card ends the
program.

     Table 1 provides an individual  description of each data card, its
argument list, format code, and column entries.

Calculation of the Round-Off Error

     The variable UROUND in block data ALPHA!  should be set to the round-
off error associated with each computer system.  Currently, UROUND is set
to the round-off error of 7.5 x 10    associated with a CDC 7600 computer
system.  To reset UROUND, the card (W.I6) in block data ALPHA1 should be
changed:

      DATA     UROUND/user's round-off error/     (card W.16)

     UROUND is calculated from the number of significant digits (N) used
for the mantissa of a floating point constant:

                           UROUND = 2~N

For example, on the UNIVAC 1110 computer, each word contains 36-BITs of
                                            -27
which 27 are used for the mantissa.  Thus, 2    is equal to approximately
7.5 x 10"9.  This is the value required for UROUND in CHEMK.
                                     15

-------
                   TABLE 1.   INPUT CARD DESCRIPTION
        Card
Column   Format
Description
General Data Card
Next card(s)
Next card(s)
  1-3     13      Last number of reaction sequence

  6-10    15      Choice of reaction input format:
                  Nonzero=EPASIM, default=CHEMK
                  format

 11-14    A4      Print option for net rates of
                  reaction and species rate of
                  formation:  default=printing is
                  suppressed

 15-20    —      Not read

 21-30   F10.0    Dilution factor
 31-35    15      Number of photolysis time units

 36-40    15      Time interval of each photolysis
                  time unit

                  Number of photolysis reactions

                  Number of temperature time units

                  Time interval of each temperature
                  time unit

                  Option to plot photolysis values
                  as a function of time

 61-80    —      Not read

  1-10   7110     If there are variable photolysis
 11-20            rate constants, the reaction
 21-30            identification vector identify-
                  ing number of each photolysis
                  reaction is read here
    •
 61-70
 71-80    —      Not read

  1-10   7F10.0   The initial photolysis and the
 11-20            value at the end of each time
 21-30            unit.  If photolysis occurs for
 31-40            N units of time, then there
                  should be N+l photolysis values
                     41-45    15

                     46-50    15

                     51-55    15


                     56-60
                     61-70

                     71-80
                  Not read

                  16

-------
                         Table 1  (Continued)
          Card
Column   Format
                 Description
Next card(s)
Kinetic Mechanism Card
(EPASIM format)
 1-10
11-20
21-30
31-40
61-70

71-80

 1-4,
 7-11,
13-16,
17-18
                             19-23
                             24
                             25-28


                             32-36


                             37-40


                             41-45
7F10.0   If variable temperature is
         used then the temperature
         values are read here.   If
         there are N units of time
         in which the temperature
         changes, then there should
         N+l temperature values
 A4
 A4
 A4
 12
        Not read

        The names  of the  reactants
        (up to three are  allowed)

        The identification  number
        of the reaction (number
        may be between 0  and 99; 0
        only if next field  is
        nonzero)

F5.0    The stoichiometric  coeffi-
        cient of the first  product;
        a blank implies a default
        value of 1
II      Modulus factor for  reaction
        identification numbers
        greater than 100  (values are
        0 or blank,  1, and  2);  if  a
        value of 1 is entered,  the
        reaction identification num-
        ber will be  between 100 and
        199 depending on  the value
        in columns 17 and 18
A4      The name of the  first
        product
F5.0    The stoichiometric  coeffi-
        cient of the second product

A4      The name of the  second
        product
F5.0    The stoichiometric  coeffi-
        cient of the third  product
                                      17

-------
                       Table 1  (Continued)
          Card
Column   Format
                          Description
Kinetic Mechanism Card
(CHEMK format)*
Title Card
Output Specification
Card
49-52     A4      The name of the third
                  product
55-64     F10.0   The rate constant at 298K
                  or, for photolysis reac-
                  tions, the ratio of the
                  photolysis value to reac-
                  tion number 1

66-72     F7.0    The activation energy in
                  degrees Kelvin (K)

          13      The reaction identifica-
                  tion number

          A4      The names of the reactants
                  (up to three allowed)


          A4      The names of the products
                  (up to four allowed)
 1- 3


 6- 9
11-14
16-19
21-24
26-29
31-34
36-39

41-50
                             51-60
 1-72
F10.0   The rate constant at 298K
        or the photolysis ratio
        for photolysis reactions

F10.0   The activation energy in
        degrees Kelvin (K)

18A4    The first four columns must
        not all be blank or contain
        one of the following:  MORE,
        CONT, PLOT
 1-3     13      Number of reactants with
                  initial nonzero concentra-
                  tions (maximum is 89)

21-30     110     Number of integration steps
                  between print.cycles;
                  default is 104 to suppress
                  printing keyed to step size.

31-40     F10.0   Specific time of initial
                  print (independent of step
                  size)
                                     18

-------
                           TABLE 1  (Continued)
         Card
Column   Format
                          Description
Initial  Reactant Card(s)
41-50     F10.0   Time interval
51-60     F10.0   Time multiple (note that
                  at least one of the above
                  two fields must be blank)

 1-4,    A4      Alphanumeric identification
11-14             of the reactants with ini-
21-24             tial nonzero concentrations
                  seven species per card
Initial  Reactant Con-
centration Card(s)
Time Specification Card
61-64

65-80

 1-10
11-20
21-30
61-70

71-80

 1-10


11-20


21-30
                  Not read

          F10.0   Initial concentrations  of
                  reactants.   Concentrations
                  must be given  in  the same
                  sequence  as  the  species are
                  given on  Initial  Reactant
                  Card(s)
                  Not read

          F10.0   Simulation starting  time.
                  Time must be in the  same
                  units as the rate constants
          F10.0   Simulation ending time.
                  Time must be in the  same
                  units as the rate constants
          F10.0   Temperature (constant
                  throughout the simulation
                  if variable temperature
                  profile is not used;
                  default=298K)
31-40     F10.0   Error tolerance; default=10

41-80     --      Not read
                                                                          -2
                                    19

-------
                           TABLE  1  (Continued)
         Card
Col limn	Format
Description
Options Card
Plot Option Cards:

  Plot Title
  Plot Specification
  Card
 1- 4     A4      Plot initiated if PLOT is
                  entered.   Additional  kin-
                  etic steps are processed
                  if MORE is entered.   A
                  new starting time is  used
                  if CONT is entered.   Pro-
                  gram terminates if blank
                  card is inserted.

 5-80     —      Not read
 1-12     3A4     Title of simulation
13-54     —      Not read

55-60             Option to print actual  con-
                  centrations of all  plotted
                  species
61-65     15      Data storage option.  If
                  nonzero, plot data is saved
                  on a file named TAPE?

66-80     —      Not read
 1-5     —      Not read

 6-7     12      Number of species to be
                  plotted; up to three species
                  may be placed on one plot

 8-15     —      Not read

16-19     A4      Concentration labels:
21-24             Five arguments must be
26-29             entered
31-34
36-39

41-50     F10.0   Lowest concentration value

51-60     F10.0   Highest concentration value

61-70     F10.0   Lowest time value

71-80     F10.0   Highest time value
                                    20

-------
                         TABLE  1  (Concluded)
         Card
Column   Format
Description
Species Identification
Card (one card for each
species)
 1-4    A4      Name of species to be
                  plotted

 6- 7     12      Number of data points to
                  be plotted from experimen-
                  tal data.  Value must be
                  between 1 and 80.  If set
                  at zero, no experimental data
                  are read



Experimental Plot Data



End of Plot Data

End of Program

8
9
10-80
1-10
11-20
21-30
31-40
61-70,
71-80
1- 4
5-80
1- 4
5-80
—
Al
—
F10.0
F10.0
F10.0

A4
—
A4
—
Not read

Symbol for plotted output
Not read
Time
Concentration
(A total of four sets
time and concentration
placed on each card)

Blank card
Not read
Blank card
Not read



of
are






*  Default reaction  format
                                    21

-------
                             REFERENCES
Gear, C. W. (1971), Numerical  Initial  Value Problems  In  Ordinary  Differen-
     tial Equations (Prentice-Hall,  Englewood  Cliffs,  New Jersey).

Hindmarsh, A. C. (1974), "GEAR:   Ordinary Differential Equation System
     Solver," Report UCID-30001,  Rev.  3, Lawrence Livermore  Laboratory,
     Livermore, California.

Overton, J. H. (1976), "Users  Guide  to EPASIM--A Chemical Kinetics  Simu-
     lation Program," TN-262-1643, Northrop Services,  Incorporated,
     Huntsville, Alabama.

Spellmann, 0. W., and A. C.  Hindmarsh  (1975),  "GEARS:  Solution of  Ordinary
     Differential Equations  Having a Sparse Jacobian  Matrix," Report UCID-
     30116, Lawrence Livermore Laboratory, Livermore,  California.

Whitten, G. Z. (1974), "Rate Constant  Evaluations Using  a New Computer
     Modeling Scheme," 167th National  Meeting, American  Chemical  Society.
                                   22

-------
  APPENDIX A
PROGRAM LISTING
       A-l

-------
                                                                         A    3
C****AAA******A*************AAA******************AAAAAAAAAAA***********A /\    4
C                                                                        A    5
C     CHEM< IS A FORTRAN PROGRAM WHICH, WHEN GIVEN A PREDETERMINED       A    6
C     KINETIC MECHANISM, COMPUTES THE CONCENTRATIONS OF THE VARIOUS      A    7
C     REACTANTS IN TIME. IT WAS WRITTEN BY GARY Z. WHITTEN OF SYSTEMS    A    8
C     APPLICATIONS INCORPORATED IN SAN RAFAEL, CA AND INCORPORATES       A    9
C     THE GEAR INTEGRATION PACKAGE OF A. HINDMARSH OF LAWRENCE LIVERMORE A   10
C     LABORATORIES IN LIVERMORE, CA.  THE PRINTER-PLOT ROUTINE WAS       A   11
C     PROVIDED BY DAVID C. WHITNEY OF SAI AND THE PROGRAM WAS CONVERTED  A   12
C     ANSI FORTRAN BY JIM MEYER OF SAI. THE  COMPUTER CODE THAT FOLLOWS  A   13
C     IS EXECUTABLE ON ANY IBM OR CDC MACHINE WITH A FORTRAN IV COMPILER A   14
C                                                                        A   15
C     TWO SEPARATE SETS OF DATA ARE NEEDED TO EXECUTE THE PROGRAM. ONE   A   16
C     SET CONTROLS THE COMPUTATIONAL PART OF THE PROGRAM AND IS          A   17
C     NECESSARY TO ESTABLISH THE INTEGRATION ROUTINE.  THE SECOND SET OF A   18
C     CONTROLS THE PLOTTER OUTPUT AND PROVIDES THE PARAMETERS NECESSARY  A   19
C     SPECIFY THE  OUTPUT FORMAT.  THE READER IS REFERED TO THE PROGRAM  A   20
C     FOR FURTHER INFORMATION.                                           A   21
C                                                                        A   22
C                                    GARY Z. WHITTEN                     A   23
C                                    SYSTEMS APPLICATIONS INCORPORATED   A   24
C                                    950 NORTHGATE DRIVE                 A   25
C                                    SAN RAFAEL, CALIFORNIA  94903       A   26
C                                    (415) 472-4011                      A   27
C                                                                        A   28
C**'********A*****A************AAAAA***'****A***AAAA**AAA**A*******'**l**AAA ^   £9
C                                                                        A   30
      COMMON /DAW NR,KR(200,7),A(200),S(200),ITITLE(7),TEM?,ERR,START, A   31
      lSTOPP,PC,NP,SIG(91)fIP(91),ITYPE(200),R(200),BK,SG,DILUT           A   32
      COMMON /NAMES/ SPECIS(91),REACT(91),NS                             A   33
      COMMON /FRPLOT/ NIT(3),SAVCON(90,80),SAVTIM(80),JGRID(89,40),NT    A   34
      COMMON /ALPHA/ IGO(4),IBLANK,MBLANK,OINTER                         A   35
      COMMON /APLOT/ JVERT(52,2),JBLANK,JSTAR,JPLUS,JBAR                 A   36
      COMMON /GEAR!/ T,GUESS,HMIN,HMAX,EPS1,UROUNO,NC,MF1,KFLAG1,JSTART A   37
      COMMON /GEAR2/ YMAX(IOO)                                           A   38
      COMMON /INOUT/ IN,IOUT,ITAPE                                       A   39
      COMMON /HEAT/ CV,Q(200),SC(200,7),ISC(200,3),ITEM3                 A   40
      COMMON /SPARSE/  IA(91),JA(1000)                                    A   41
      COMMON /STORE/ AST(35),IPL(7),TEMEND,NTEMP,TMI,NPHOT,PHI,IL,NFRST, A   42
      1IPH(30),QM(100),PM(100),PSTOP                                     A   43
      COMMON /PHOTR/ IPP,RDAT(80),RTIM(80),RR1,IN10                      A   44
      DIMENSION C(91),  IRS(200,7),  KRS(7),  SC1(7)                        A   45
      INTEGER  SPECIS,REACT,PHI,TMI,AST                                  A   46
C                                                                         A  47
C     ALPHAMERIC DATA  ARE PREASSIGNED  IN BLOCK DATA  ALPHA!                A   48
C                                                                         A  49
      DATA NPLOT/4HSAVE/,IBLK/1H /                                        A   50
      DATA IN2,IN3,IN5,IN6,IN7,IN8,IN9/7*1/                              A    51
      DATA XN4,XN6,XN9/3*1.0/                                            A   52
                                      A-2

-------
C                                                                        A   53
C     INITIAL PARAMETERS                                                 A   54
C                                                                        A   55
      IN=5                                                               A   56
      IOUT=6                                                             A   57
      ITAPE=7                                                            A   58
C                                                                        A   59
C     CLEAR PATTERN MATRIX AND SET THE FIRST ELEMENTS                    A   60
C     ALSO SET STOICHIOMETRIC COEFFICIENTS EQUAL TO 1                    A   61
C                                                                        A   62
      DO 5 0=1,200                                                       A   63
      DO 5 K=l,7                                                         A   64
      IRS(J,K)=IBLANK                                                    A   65
      SC(J,K)=1.                                                         A   66
    5 KR(J,K)=0                                                          A   67
      DO 10 J=l,200                                                      A   68
      A(J)=0.                                                            A   69
      R(J)=0.                                                            A   70
   10 S(J)=0.                                                            A   71
      NR=0                                                               A   72
      NS=0                                                               A   73
C                                                                        A   74
C     NX = LAST NUMBER OF THE REACTION SEQUENCE                          A   75
C     NPLOT = PLOT OPTION                                                A   76
C     DILUT= DILUTION FACTOR                                             A   77
C                                                                        A   78
      READ (IN.185) NX,NEPA,NPLIT,DILUT,NPHOT,PHI,IL,NTEMP,TMI,IPP       A   79
      GO TO 20                                                           A   80
   15 NS=NS-1                                                            A   81
      READ (IN.185) NX,IN2,IN3,XN4,IN5,IN6,IN7,IN8,IN9,IN10              A   82
      IF (IABS(IN2).NE.O) NEPA=IN2                                       A   83
      NPLIT=IN3                                                          A   84
      IF (ABS(XN4).NE.O) DILUT=XN4                                       A   85
      IF (IABS(IN5).NE.O) NPHOT=IN5                                      A   86
      IF (IABS(IN6).NE.O) PHI=IN6                                        A   87
      IF (IABS(IN7).NE.O) IL=IN7                                         A   88
      IF (IABS(IN8).NE.O) NTEMP=IN8                                      A   89
      IF (IABS(IN9).NE.O) TMI=IN9                                        A   90
      IF (lABS(INlO).NE.O) IPP=IN10                                      A   91
   20 IF (IL.LE.O.OR.IN7.LE.O) 60 TO 25                                  A   92
      READ (IN.200) (IPH(I),I=1,IL)                                      A   93
   25 IF (NPHOT.LE.O.OR.IN5.LE.O) 60 TO 30                               A   94
      NPHOT=NPHOT+1                                                      A   95
      READ (IN.190) (PM(I),I=1,NPHOT)                                    A   96
      PM(NPHOT+1)=2.0*PM(NPHOT)-PM(NPHOT-1)                              A   97
      PM(NPHOT+2)=3.0*PM(NPHOT)-2.0*PM(NPHOT-1)                          A   98
      PSTOP=FLOAT((NPHOT-1)*PHI)                                         A   99
   30 IF (NTEMP.LE.O.OR.IN8.LE.O) 60 TO 35                               A  100
      NTEMP=NTEM?+1                                                      A  101
      READ (IN,190) (QM(I),I=1,NTEMP)                                    A  102
                                     A-3

-------
      qM(NTEMs-H)*2.*QM(NTEMJ)-qM(NTef>-l)                                A  103
      QM(NTEMP+2)=3.0*QM(NTEMP)-2.0*qM(NTEMP-l)                           A  104
      TEMENO=FLOAT((NTEW>-1)*TMI)                                        A  105
   35 IF (NX.GT.O) WRITE (IOUT.130)                                       A  106
      IF (NX.LE.O) NS-NS+1                                               A  107
      IF (NX.LE.O) 60 TO 80                                              A  108
C                                                                        A  109
C     REACTION INPUT DATA                                                A  110
C                                                                        A  111
   40 IF (NEPA.LE.O) READ (IN.135) J,(IRS(J,I),I=1,7),A(J),S(J)          A  112
      IF (NEPA.LE.O) 60 TO 55                                            A  113
      READ (IN.220) (KRS(I),I=1,3),J,SC1(1),JJ,KRS(4),(SC1(LL-3),KRS(LL) A  114
     1,LL=5,6),RTE,ENERGY                                                A  115
      IF (JJ.GT.O) J=JJ*100+J                                            A  116
      DO 45 11-1,6                                                       A  117
   45 IRS(J,II)=KRS(II)                                                  A  118
      A(J)=RTE                                                           A  119
      S(J)=£NERGY                                                        A  120
      DO 50 11=4,6                                                       A  121
   50 SC(J,II)=SCl(II-3)                                                 A  122
   55 DO 60 1=1,7                                                        A  123
   60 IF (SC(J,I),LE.O.) SC(J,I)-1.                                      A  124
      DO 75 K=l,7                                                        A  125
      NL=K*5                                                             A  126
      NF=NL-4                                                            A  127
      DO 65 LK=NF,NL                                                     A  128
   65 AST(LK)=IBLK                                                       A  129
      IPL(K)=IBLK                                                        A  130
      IF (K.EQ.1.0R.K.EQ.4) 60 TO 70                                     A  131
      IF (IRS(J.K).NE.IBLANK) IPL(K-1)=JPLUS                             A  132
   70 IF (IRS(J,K).NE.IBLANK) CALL VALU (SC(J,K),K,NF,NL)                A  133
   75 CONTINUE                                                           A  134
      WRITE (IOUT.140) J,(AST(I),I=1,5),IRS(,J,1),IPL(1),(AST(I),I=6,10), A  135
     1IRS(J,2),IPL(2),(AST(I),1*11,15),IRS(J,3),IPL(3),(AST(I),1=16,20), A  136
     2IRS(J,4).IPL(4)f(AST(I),1*21,25),IRS(J,5),IPL(5),(AST(I),1=26,30), A  137
     3IRS(J,6),IPL(6),(AST(I),1=31,35),IRS(J,7),IPL(7),A(J),S(J)         A  138
      KR(J,1)=100                                                       A  139
      IF (J.GT.NR) NR
-------
      IF (ITITLE(l).EQ.IGO(l)) GO TO 15                                  A  153
      IF (ITITLE(l).EQ.160(2)) GO TO 115                                 A  154
      IF (ITITLE(l).EQ.160(3)) GO TO 85                                  A  155
      IF (ITITLE(l).EQ.IBLANK) STOP                                      A  156
      GO TO 90                                                           A  157
C                                                                        A  158
C     CALL PLOT                                                          A  159
C                                                                        A  160
   85 READ (IN,180) (NIT(I),I»1,3),IDT,KALCM>                            A  161
      NT=NT-2                                                            A  162
      CALL PLOT (NIT,NT,NS,SPECIS,SAVTIM,SAVCON,IDT,KALCM?,JGRID)         A  163
      GO TO 80                                                           A  164
C                                                                        A  165
C     OPTIONS CARD                                                       A  166
C                                                                        A  167
   90 READ (IN,155) N,NPRNT,TPRNT,TSTEP,TFACT                            A  168
C                                                                        A  169
C     TIME STEP SKIP OPTION                                              A  170
C                                                                        A  171
      IF (NPRNT*NPRNT.EQ.O) NPRNT=100000                                 A  172
C                                                                        A  173
C     TIME STEP LENGTH OPTION                                            A  174
C                                                                        A  175
      IF (TPRNT*TPRNT.EQ.O.) TPRNT=1.E10                                 A  176
C                                                                        A  177
C     CONCENTRATION OF SPECIES INITIALLY PRESENT                         A  178
C                                                                        A  179
      READ (IN.125) (REACT(I),I=1,N)                                     A  180
      READ (IN,190) (C(I),I=1,N)                                         A  181
C                                                                        A  182
C     STARTING AND ENDING INTEGRATION TIMES                              A  183
C                                                                        A  184
      READ (IN.175) START,STOPP.TEW,ERR                                 A  185
      IF (START*START.EQ.O.) START=0.                                    A  186
C                                                                        A  187
C     SPECIFICATION OF THE TEMPERATURE IF UNSPECIFIED IN INPUT           A  188
C                                                                        A  189
      IF (TEMP.LE.O.) TEMP=298.                                          A  190
C                                                                        A  191
C     SPECIFICATION OF THE ERROR BOUND IF UNSPECIFIED IN INPUT           A  192
C                                                                        A  193
      IF (ERR*ERR.EQ.O.) ERR=l.E-2                                       A  194
C                                                                        A  195
C     OUTPUT OF THE INITIAL CONDITIONS                                   A  196
C                                                                        A  197
      WRITE (IOUT.145) (ITITIE(I),1=1,7),(REACT(I),I=1,N)                A  198
   95 WRITE (IOUT,150) (C(I),I=1,N)                                      A  199
C                                                                        A  200
C     COMPUTE THE NET RATES OF REACTION                                  A  201
C                                                                        A  202
                                     A-5

-------
      IF (NTEMP.6T.O) TEM?=QM(1)                                         A  203
      IF (ITITLE(1).EQ.IGO(2).AND.NTEMP.GT.O) TEMP=TEMOLD                A  204
      CALL RATES (C,N)                                                   A  205
C                                                                        A  206
C     OUTPUT OF THE TEMPERATURE AND ERROR BOUND                          A  207
C                                                                        A  208
      WRITE (IOUT.160) TEM>,ERR                                          A  209
C                                                                        A  210
C     OUTPUT OF THE DILUTION FACTOR                                      A  211
C                                                                        A  212
      IF (DILUT*DILUT.EQ.O.) DILUTX).                                    A  213
      IF (DILUT.NE.O.) WRITE (IOUT.170) DILUT                            A  214
C                                                                        A  215
C     SET LIMITS FOR TIMED OUTPUTS                                       A  216
C                                                                        A  217
      IF (TSTEP.NE.O.) YMAX(1)=TSTEP                                     A  218
      IF (TSTEP.NE.O.) PC=TSTEP                                          A  219
      IF (TSTEP*TSTEP.EQ.O.) YMAX(1)=1.E10                               A  220
      IF (TSTEP*TSTEP.EQ.O.) PC*1.E10                                    A  221
      IF (TFACT.NE.O.) YMAX(2)»TFACT                                     A  222
      IF (TFACT*TFACT.EQ.O.) YMAX(2)=1.                                  A  223
C                                                                        A  224
C     RATE CONSTANTS                                                     A  225
C                                                                        A  226
      WRITE (IOUT,165) (R(I),M,NR)                                     A  227
C                                                                        A  228
C     WRITE INITIAL CONDITIONS OF THE CELL                               A  229
C                                                                        A  230
      IF (NPHOT.LE.O) 60 TO 105                                          A  231
      WRITE (IOUT.215) (IPH(I),I=1,IL)                                   A  232
      DO 100 1=1,IL                                                      A  233
      J=IPH(I)                                                           A  234
      YMAX(10*I)=R(J)                                                    A  235
  100 CONTINUE                                                           A  236
      WRITE (IOUT,225) (YMAX(10+I),M,IL)                               A  237
      NPT=NPHOT                                                          A  238
      WRITE (IOUT.205) PHI,(PM(I),I«1,NPT)                               A  239
  105 IF (NTEMP.LE.O) 60 TO 110                                          A  240
      NPT=NTEM?                                                          A  241
      WRITE (IOUT,210) TMI,(OM(I),M,NPT)                               A  242
  110 6UESS=1.E-10                                                       A  243
      T=START                                                            A  244
      IF (NPLOT.EQ.IBLANK) YMAX(3)=0.                                    A  245
      IF (NPLOT.NE.IBLANK) YMAX(3)=1.                                    A  246
      IF (NPLIT.NE.IBLANK) YMAX(4)=1.                                    A  247
      IF (NPLIT.EQ.IBLANK) VMAX(4)=0.                                    A  248
C                                                                        A  249
C     INITIALIZE PARAMETERS                                              A  250
C                                                                        A  251
      CALL YFIX (NS,TPRNT,C,NPRNT,INX)                                   A  252
                                      A-6

-------
      INX=4                                                              A  253
      TEMOLD=TEMP                                                        A  254
      GO TO 80                                                           A  255
C                                                                        A  256
C     CONTINUATION OF DATA                                               A  257
C                                                                        A  258
  115 READ (IN.175) START,STOPP,TEMP,ERR                                 A  259
      IF (TEMP.LE.O.) TEMP=298.                                          A  260
      IF (TEMOLD.NE.298..AND.TEMOLD.NE.O.) TENP=TEMOLD                   A  261
      IF (ERR*ERR.EQ.O.) ERR=l.E-2                                       A  262
      DO 120 1=1,NS                                                      A  263
  120 REACT(I)=SPECIS(I)                                                 A  264
C                                                                        A  265
C     INITIAL CONDITIONS                                                 A  266
C                                                                        A  267
      WRITE (IOUT.145) (ITITLE(I),I=1,7),(REACT(I),I=1,NS)               A  268
      N=NS                                                               A  269
      GO TO 95                                                           A  270
C                                                                        A  271
C                                                                        A  272
  125 FORMAT (7(A4,6X))                                                  A  273
  130 FORMAT (1H1.14H THE REACTIONS,86X,13HRATE CONSTANT.2X.15HACT. ENER A  274
     1GY (K))                                                            A  275
  135 FORMAT (I3,2X,7(A4,1X),2F10.0)                                     A  276
  140 FORMAT (/,2X,I3,2X,3(5A1,1X,A4,2X,A1),1H=,2X,4(5A1,1X,A4,2X,A1),1P A  277
     1E11.3.2X.E13.3)                                                    A  278
  145 FORMAT (1H1,30X,7A4,//,23H INITIAL CONCENTRATION ,//(10X,10(4X,A4, A  279
     14X)))                                                              A  280
  150 FORMAT (/(8X.1P10E12.3))                                           A  281
  155 FORMAT (I3,17X,I6,4X,3F10.0)                                       A  282
  160 FORMAT (/,34H THE TEMPERATURE OF THE CELL WAS =,1PE9.2,26H AND THE A  283
     1 ERROR TOLERANCE =,E9.2)                              <             A  284
  165 FORMAT (/,29H THE RATE CONSTANTS USED WERE,/,(/,8X,lP10E12.3))     A  285
  170 FORMAT (/,30H THE OVERALL DILUTION RATE WAS.1PE9.2)                A  286
  175 FORMAT (8F10.0)                                                    A  287
  180 FORMAT (3A4,43X,2I5)                                               A  288
  185 FORMAT (I3,2X,I5,A4,6X,F10.2,6I5)                                  A  289
  190 FORMAT (7F10.3)                                                    A  290
  195 FORMAT (7A4)                                                       A  291
  200 FORMAT (7110)                                                      A  292
  205 FORMAT (/44HOTHE N02 PHOTOLYSIS NUMBERS IN INTERVALS OF ,I3,15H TI A. 293
     1ME UNITS ARE,/,(1HO,1P10E13.3))                                    A  294
  210 FORMAT (/34HOTHE TEMPERATURES IN INTERVALS OF ,I3,15H TIME UNITS A A  295
     1RE/(1HO,1P10E13.3))                                                A  296
  215 FORMAT (/29HOTHE PHOTOLYSIS REACTIONS ARE,/(1HO,9I13))             A  297
  220 FORMAT (2(A4,2X)>A4,I2,F5.1,I1,A4,2X,2(F5.0,1X,A4,2X),F10.2,1X,F7. A  298
     12)                                                                 A  299
  225 FORMAT (57HOTHE PHOTOLYSIS RATIOS OF EACH OF THE REACTIONS ABOVE A A  300
     1RE/(1HO,4X,1P9E13.3))                                              A  301
      END                                                                A  302-
                                     A-7

-------
SUBROUTINE YFIX ( NO, TLAST,C,NQ, INDEX)
C
C
C THIS IS THE NORMAL OUTPUT ROUTINE
C

C




















C
C
C




















COMMON /SPARSE/ IA(91),JA(1000)
COMMON /DATA/ NR,KR(200,7),A(200),S(200),ITITLE(7),TEMP,ERR,START,
1STOPP,PC,NP,SIG(91) ,IP(91) ,ITYPE(200) ,R(200) ,BK,SG,DILUT
COMMON /NAMES/ SPECIS(91),REACT(91),NS
COMMON /FRPLOT/ NIT(3),SAVCON(90,80),SAVTIM(80),JGRID(89,40),NT
COMMON /6EAR1/ D,H,DUM(4),IDUM(4)
COMMON /6EAR2/ YMAX(IOO)
COMMON /INOUT/ INPP,IOUT,ITAPE
COMMON /HEAT/ CV,Q(200),SC(200,7),ISC(200,3),ITEMP
COMMON /STORE/ AST(35),IPL(7),TEMEND,NTEMP,TMI,NPHOT,PHI,IL,NFRST,
1IPH(30) ,QM(100) .PMdOO) .PSTOP
COMMON /PHOTR/ IPP,RDAT(80),RTIM(80),RR1,IN10
DIMENSION C(91), RT(200)
INTEGER SPECIS,REACT,TMI,PHI,AST
DATA IGO/4HCONT/
DATA NT1/1/
DATA START1/0.0/
N=NS-1
NS2*NS-2
NFRST=1

TIME INTERVALS

IF (ITITLE(l).NE.IGO) GO TO 30
IF (NT1.EQ.O) GO TO 30
TDC=(STOPP-START1)/80.
TD=START1+TDC
J=0
NTSV=NT
DO 25 1=1,80
5 J*J+1
IF (J.GT.NTSV) GO TO 20
IF (SAVTIM( J) .GE.TD) GO TO 10
GO TO 5
10 SAVTIM(I)=SAVTIM(J)
DO 15 K-1.NS
15 SAVCON(K,I)=SAVCON(K,J)
NT=I
GOTO 25
20 IF (I.EQ.l) NT=1
SAVTIM(I)=TD
25 TD=TD+TDC
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
B
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
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
A-8

-------
      60 TO 40                                                           8   51
   30 TDC=(STOPP-START)/80.                                               B   52
      TD=START                                                           B   53
      DO 35 1=1,80                                                       B   54
      TD=TD+TDC                                                          B   55
   35 SAVTIMdHD                                                       B   56
      START1=START                                                       B   57
      NT=1                                                               B   58
   40 NC=0                                                               B   59
C                                                                        B   60
C     INITIALIZE PARAMETERS                                              B   61
C                                                                        B   62
      NHS=0                                                              B   63
      MS=IABS(NQ)                                                        B   64
      TPRNT=TLAST+START                                                  B   65
      TSTEP=PC                                                           B   66
      TFACT=YMAX(2)                                                      B   67
      IF (TFACT.GT.l.) TSTEP=0.                                          B   68
      MT=50-((NS-1)/10)*3-(NR-1)/10                                      B   69
      NN=NS+1-MINO(NS/10,1)                                              B   70
      IF (YMAX(3).EQ.O.) NPLOT=0                                         B   71
      IF (YMAX(3).NE.O.) NPLOT=1                                         B   72
      IF (YMAX(4).EQ.O.) NOPT=0                                          B   73
      IF (YMAX(4).NE.O.) NOPT=1                                          B   74
      IF (NOPT.EQ.O) MT=50-((NS-1)/10)*2                                 B   75
      NTP=0                                                              B   76
      NT1=NPLOT                                                          B   77
      IF (NPLOT.EQ.O) NT=-1                                              B   78
      IN=INDEX                                                           B   79
      T=START                                                            B   80
      TNEXT-SWM.                                                     B   81
      CALL DIFFUN (NS2,TNEXT,C,RT)                                       B   82
      CALL SAVLIN (T.C.NS2)                                              B   83
      NFRST=2                                                            B   84
      IF (T-START) 60,60,50                                              B   85
   45 TNEXT=AMIN1(TPRNT,STOPP)                                           B   86
      IF (T.EQ.START) GO TO 55                                           B   87
   50 IF (NQ.LT.O.AND.IABS(NQ).NE.O) IN=3                                B   88
   55 CALL DRIVES (NS2,T,H,C,TNEXT,ERR,21,IN,IA,JA)                      B   89
      T=TNEXT                                                            B   90
      IF (T.GE.STOPP) TPRNT=STOPP                                        B   91
   60 IF (T.EQ.START) TIMNW=START                                        B   92
      IF (T.NE.START) TIMNW=TPRNT                                        B   93
      IF (T.EQ.START) GO TO 65                                           B   94
      IF (NQ.LT.O.AND.IABS(NQ).NE.O) TIMNW=T                             B   95
      NHS=NHS+1                                                          B   96
      IF (NHS.GE.MS.OR.T.GE.TPRNT) GO TO 65                              B   97
      GO TO 45                                                           B   98
   65 NHSO                                                              B   99
      CALL DIFFUN (NS2,TIMNW,C,RT)                                       B  100
                                     A-9

-------
      IF  (NOPT.LE.O) NTP=NTP+(N/10)+4                                    B  101
      IF  (NOPT.GT.O) lfTP*tfl?+(N/10)*3+(NR-l)/lCH-12                       B  102
      IF  (NTP.GT.MT.OR.T.EQ.START) WRITE (IOUT,120)                      B  103
      IF  (T.EQ.START) 60 TO 70                                           B  104
      IF  (NOPT.LE.O.AND.NTP.LE.MT) 60 TO 75                              B  105
   70  CONTINUE                                                           B  106
      IF  (NTP.6T.MT) NTP=0                                               B  107
      IF  (NN.6E.11) WRITE (IOUT,H5) (SPECIS(I),I=1,NN)                  B  108
      IF  (NN.LE.10) WRITE (IOUT,125) (SPECIS(I),I=1,NN)                  B  109
   75  CONTINUE                                                           B  110
      IF  (NS.6E.11) WRITE (IOUT,130) TIMNW,(C(I),I=1,10),H,(C(I),I*11,NS B  111
    12)JEMP,SG                                                         B  112
      IF  (NS.LE.10) WRITE (IOUT,130) TIMNW,(C(I),M,NS2),TEM?,SG,H      B  113
      IF  (NOPT.EQ.O) WRITE (IOUT.135)                                    B  114
      IF  (NOPT.EQ.O) GO TO 95                                            B  115
      RT(NS)=0.                                                          B  116
      RT(N)=0.                                                           B  117
      DO  80 1=1,NS2                                                      B  118
   80  RT(NS)=RT(NS)+RT(I)                                                B  119
      WRITE (IOUT.105)  (RT(I),I=1,NS)                                    B  120
      DO  90 M.NR                                                      B  121
      J=KR(I,1)                                                          B  122
      IF  (J.EQ.O)  RT(I)=0.                                               B  123
      IF  (J.EQ.O)  GO TO 90                                               B  124
      JT*ITYPE(I)                                                        B  125
      XT»1.                                                             B  126
      DO  85 L=1,JT                                                      B  127
      J=KR(I,L)                                                          B  128
      X>1.                                                             B  129
      IF  (J.6T.O.AND.J.LE.NS2)  XJ=C(J)**ISC(I,L)                         B  130
      IF  (ISC(ItL).EQ.-l) XJ=C(J)**SC(I,L)                               B  131
      IF  (J.EQ.NS) X>SG                                                 B  132
      XT=XT*XJ                                                          B  133
   85  CONTINUE                                                          B  134
      RT(I)-XT*R(I)                                                      B  135
   90  CONTINUE                                                          B  136
      WRITE (IOUT.110)  (RT(I),I=1,NR)                                    B  137
   95  IF  (TIMNW.EQ.START) 60 TO 55                                      B  138
      IF  (T.6E.STOPP)  RETURN                                            B  139
      IF  (IN.EQ.3.AND.T.LT.TPRNT)  GO TO 45                               B  140
      TPRNT=TPRNT*TFACT^TSTEP                                           B  141
      GO  TO 45                                                          B  142
C                                                                        B   143
C                                                                       B  144
  105  FORMAT (/,10H  NET RATES,2X,1P10E12.3,/,(12X,1P10E12.3))            B   145
  110  FORMAT (//,2X,22HTHE  REACTION RATES ARE,/,(1H .1P10E13.2))         B  146
  115  FORMAT (/,4X,5HTIME ,4X,10(4X,A4,4X),/,2X,8HINTERVAL,3X,10(4X,A4,4 B   147
    1X),/,(13X,10(4X,A4,4X)))                                           B  148
  120  FORMAT (1H1)                                                       B   149
  125  FORMAT (/,4X,5HTIME ,4X,10(4X,A4,4X),/)                            B  150
                                    A-10

-------
130 FORMAT (/,lPllE12.3/tllE12.3/,(12X,10E12.3))                        B   151
135 FORMAT (1H )                                                        B   152
    END                                                                B   153-
                                  A-ll

-------
SUBROUTINE RATES (C,N)
C
C
C THIS SUBROUTINE SETS THE INITIAL CONCENTRATIONS AND CALCULATES THE
C RATE CONSTANTS
C

C

























C
C
C
C
C






C
C
C



INTEGER SPECIS,REACT
COMMON /DATA/ NR,KR(200,7) ,A( 200) ,S(200),ITITLE(7), TEMP, ERR,START,
1STOPP,PC,NP,SIG(91),IP(91),ITYPE(200),R(200),BK5SG,OILUT
COMMON /NAMES/ SPECIS(91),REACT(91),NS
COMMON /HEAT/ CV,Q(200),SC(200,7),ISC(200,3),ITEMP
COMMON /ALPHA/ IGO(4),IBLANK,MBLANK>JINTER
DIMENSION C(91)
FCT=l./298.-l./TEMP
SIGG=0.
DO 5 1=1,91
5 SIG(I)=0.
DO 15 I=1,N
DO 10 J=1,NS
IF (SPECIS(J).EQ.REACT(I)) SIG(J)=C(I)
IF (SPECIS(O).EQ.REACT(D) GO TO 15
10 CONTINUE
SIGG=SIGG+C(I)
15 CONTINUE
N=NS
M=N-1
C(N)=0.
DO 20 1=1 ,M
C(N)=C(N)+SIG(I)
20 C(I)=SIG(I)
BK=0.

IF SIG(N) DOES NOT EQUAL ZERO, IT IMPLIES THAT THAT THE CONCENTRAT
SPECIES N HAS BEEN READ. OTHERWISE M IS THE SUM OF THE INITIAL
CONCENTRATIONS.

C(N)=C(N)+SIGG
IF (SIG(N).NE.O.) BK=SIG(N)-C(N)
IF (SIG(N).NE.O.) C(N)=SIG(N)
NP=0
DO 25 1=1, NR
IF (KR(I,1).EQ.O) GOTO 25

CALCULATE THE RATE CONSTANTS

IF (ABS(S(D).EQ.O.) R(I)=A(I)
IF (ABS(S(D).NE.O.) R(I)=A(I)*EXP(S(I)*FCT)
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
C
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
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
A-12

-------
25 CONTINUE                                                          C   51
   IF (NS.LE.9)  SPECIS(NS+1)=JINTER                                   C   52
   RETURN                                                            C   53
   END                                                               C   54-
                                 A-13

-------
      SUBROUTINE DIFFUN (N,T,X,XT)                                        D    1
C                                                                        D    2
C******AA*************************************************************** Q    3
C                                                                        04
C     THIS SUBROUTINE CALCULATES THE DERIVATIVE VECTOR OF THE ODEIS       D    5
C                                                                        D    6
   ***AAAAAAA*********************************************************** D    7
                                                                         D    8
      COMMON /DATA/ NR,KR(200,7),A(200),S(200),ITITLE(7),TEM5,ERR,START, D    9
     1STOPP,PC,NP,SI6(91),IP(91),ITYPE(200),R(200),BK,SG,DILUT           0   10
      CONON /NAMES/ SPECIS(91) ,REACT(91) ,NS                             0   11
      COWDN /HEAT/ CV,Q(200),SC(200,7),ISC(200,3),ITEW>                 D   12
      COWON /STORE/ AST(35),IPL(7) JEteND,NTE^TMI,NPHOT,PHI,IL,NFRST, D   13
     1IPH(30),QM(100),PM(100),PSTOP                                      D   14
      COMMON /PHOTR/ IPP,RDAT(80),RTIM(80)fRl,IN10                       0   15
      DIMENSION XT(N), X(N)                                              D   16
      INTEGER PHI.TMI                                                    D   17
      INTEGER SPECIS                                                     D   18
      P=BK                                                               D   19
      DO 5 I=1,N                                                         D   20
      XT(I)=-OILUT*X(I)                                                  D   21
    5 P=P+X(I)                                                           D   22
      SG=P                                                               D   23
      IF (NFRST.EQ.l) FCT=((3355.7046E-6)*TEMP-1.)/TEW>                  D   24
      IF (NFRST.EQ.l) GO TO 10                                           D   25
      IF (T.EQ.TOLD) GO TO 30                                            D   26
      IF (T.LE.l.) GOTO 30                                              D   27
   10 IF (NPHOT.LE.O) GO TO 30                                           D   28
      IF (NFRST.NE.l) GO TO 15                                           D   29
      PINT=FLOAT(PHI)                                                    D   30
      PINV=1./PINT                                                       D   31
   15 CONTINUE                                                           D   32
      IF (T.GT.PSTOP) GO TO 20                                           D   33
      IZ«IFIX(T*PINV)+1                                                  D   34
      Z=T*PINV-FLOAT(IZ-1)                                               D   35
      IF (T.LE.PINT) Rl=PM(l)+(0.5*PM(3)*(Z-l.h0.5*PM(l)*(Z-3.)-PM(2)*( D   36
     1Z-2.))*Z                                                           D   37
      IF (T.LE.PINT) GO TO 20                                            D   38
      Rl=PM(IZ)-K).25*Z*(5.*PM(IZ+l)-3.0*PM(IZ)-PM(IZ-l)-PM(IZ+2)+(PM(IZ- D   39
     11)-PM(IZ)-PM(IZ+1)+PM(IZ+2))*Z)                                    D   40
   20 IF  (Rl.LT.O.)  Rl*0.                                                D   41
      DO 25  IK=i,IL                                                      D   42
      IR=IPH(IK)                                                         D   43
      R(IR)=R1*A(IR)                                                     D   44
   25 IF  (R(IR).LT.O.) R(IR)=0.                                          D   45
   30 TNOW=TEf*>                                                          D   46
      IF  (NFRST.EQ.l) GO TO 35                                            D   47
      IF (T.EQ.TOLD) GO TO 50                                            D   48
      IF  (T.LE.l.) GO TO 50                                               D   49
   35 CONTINUE                                                           D   50
                                     A-14

-------
      IF (T.GT.TEMEND)  60 TO 50                                          D   51
      IF (NFRST.NE.l) GO TO 40                                           D   52
      TINT=FLOAT(TMI)                                                     D   53
      TINV-l./TINT                                                       D   54
   40 CONTINUE                                                           D   55
      IZ=IFIX(T*TINV)+1                                                  D   56
      Z=T*TINV-FLOAT(IZ-1)                                               D   57
      IF (T.LE.TINT)  TNOW=OM(l)-H(0.5*QM(3)*(Z-l.)-»-0.5*QM(l)*(2-3.)-QM(2)  D   58
     1*(Z-2.))*Z                                                         D   59
      IF (T.LE.TINT)  GO TO 45                                            D   60
      TNOW=QM(IZ)+0.25*Z*((5.*QM(IZ+l)-3.0*QM(IZ)-QM(IZ-l)-QM(IZ+2))+(QM D   61
     1(IZ-1)-QM(IZ)-QM(IZ+1)+QM(IZ+2))*Z)                                 0   62
   45 CONTINUE                                                           0   63
      IF (TNOW.NE.TEMP.AND.TNOW.GT.O.)  FCT=((3355.7046E-6)*TNOW-1.)/TNOW D   64
   50 DO 70 IR=1,NR                                                      D   65
      I-KR(IR,1)                                                         D   66
      IF (I.EQ.O.OR.R(IR).EQ.O.) GO TO 70                                D   67
C                                                                        D   68
C     CHECK FOR A ZEROTH ORDER REACTION                                  D   69
C                                                                        D   70
      RT=1.                                                              D   71
      JT=ITYPE(IR)                                                       D   72
      DO 55 L»1,JT                                                       D   73
      I"KR(IR,L)                                                         D   74
      IF (I.EQ.99.0R.I.EQ.O) GO TO 55                                    D   75
      IF (I.NE.NS) RT=RT*X(I)                                            D   76
      IF (I.EQ.NS) RT=RT*P                                               D   77
   55 CONTINUE                                                           D   78
      IF (ABS(S(IR)).NE.O..AND.T.GT.l.) R(IR)=A(IR)*EXP(S(IR)*FCT)       D   79
      RT=RT*R(IR)                                                        D   80
      DO 60 1=1,JT                                                       D   81
      I=KR(IR,L)                                                         D   82
      IF (I.EQ.O.OR.I.GT.N) GOTO 60                                     D   83
      XT(I)=XT(I)-RT                                                     D   84
   60 CONTINUE                                                           D   85
      DO 65 K=4,7                                                        D   86
      I=KR(IR,K)                                                         D   87
      IF (I.EQ.O) GO TO 70                                               D   88
C                                                                        D   89
C     I WILL BE NEGATIVE IF CLEAN HAS BEEN CALLED                        D   90
C                                                                        D   91
      IF (I.LT.O) GO TO 65                                               D   92
      XT(I)=XT(I)+SC(IR,K)*RT                                            D   93
   65 CONTINUE                                                           D   94
   70 CONTINUE                                                           D   95
      TEMP=TNOW                                                          D   96
      TOLD=T                                                             D   97
      RETURN                                                             D   98
      END                                                                D   99-
                                     A-15

-------
SUBROUTINE MATRX (C,KR)
C

C
C
C
C
C

MATRX CREATES A NR X 7 MATRIX OF THE REACTION SCHEME WHERE EACH
ROW REPRESENTS A REACTION. THE FIRST THREE COLUMNS ARE REACTANTS
AND THE LAST FOUR ARE THE PRODUCTS. THE ELEMENTS CORRESPOND
TO THE INDIVIDUAL SPECIES AND WILL BE USED AS SUBSCRIPTS.
C

C










C
C
C





C
C
C
















C
C

INTEGER SPECIS.REACT
COMMON /DATA/ NR,IR(200,7),A(200),S(200),ITITLE(7)JEW>fERR,START,
1STOPP,PC,NP,SI6(91) ,IP(91) ,ITYPE(200) ,R(200) ,BK,SG,DILUT
COMMON /NAMES/ SPECIS(91),REACT(91),NS
COMMON /GEAR1/ T,H,HMIN,HMAX,EPS1,UR(XJND,NC,MF1,KFLAG1,JSTART
COMMON /ALPHA/ 160(4), IBLANK.MBLANK.JINTER
COMMON /HEAT/ CV,Q(200),SC(200,7),ISC(200,3),ITEMP
DIMENSION C(91), KR(200,7)
NOLD=NS+1
DO 90 I=1,NR

SKIP REACTIONS ALREADY PROCESSED

IF (IR(I,2).EQ.NOLD.OR.IR(I,3).EQ.NOLD) IR(I,3)=99
IF (IR(I,2).EQ.NOLD) IR(I,2)*0
IF (IR(I,1).EQ.NOLD) IR(I,1)=99
IF (IR(I,1).EQ.NOLD) IR(I,3)*99
IF (IABS(IR(I,1)).NE.100) GO TO 90

IF LESS THAN THREE REACTANTS, FILL FIRST SLOTS.

IF (KR(I.l).NE.IBLANK) GO TO 5
IF (KR(I,3).NE.IBLANK) KR(I,1)=KR(I,3)
IF (KR(I,3).NE.IBLANK) SC(I,1)=SC(I,3)
SC(I,3)=1.
KR(I,3)=IBLANK
IF (KR(I.l).NE.IBLANK) GO TO 5
KR(I,1)=KR(I,2)
KR(I,2)=IBLANK
SC(I,1)=SC(I,2)
SC(I,2)=1.
5 IF (KR(I,2).NE.IBLANK) GO TO 10
IF (KR(I,3).NE.IBLANK) SC(I,2)=SC(I,3)
SC(I,3)=1.
IF (KR(I,3).NE.IBLANK) KR(I,2)=KR(I,3)
KR(I,3)=IBLANK
10 DO 15 K=4,7

GET RID OF M AS A PRODUCT
E
E

E
E
E
E
E
E

E
E
£
E
E
E
E
E
E
E
E
E
E
I
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
E
1
2

4
5
6
7
8
9
« f\
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
A-16

-------
C                                                                        E   51
      IF (KR(I,K).EQ.MBLANK) KR(I,K)=IBLANK                              E   52
   15 CONTINUE                                                           E   53
C                                                                        E   54
C     IF LESS THAN FOUR PRODUCTS, FILL FIRST SLOTS.                       E   55
C                                                                        E   56
      IF (KR(I,4).NE.IBLANK) GOTO 30                                    E   57
      DO 20 K»l,3                                                        E   58
      INDEX=8-K                                                          E   59
      IF (KR(I,INDEX).NE.IBLANK) GO TO 25                                E   60
   20 CONTINUE                                                           E   61
      INDEX=5                                                            E   62
   25 KR(I,4)=KR(I,INDEX)                                                E   63
      KR(I,INDEX)=IBLANK                                                 E   64
      SC(I,4)=SC(I,INDEX)                                                E   65
      SC(I,INDEX)=1.                                                     E   66
      IF (KR(I,1).NE.IBLANK.OR.KR(I,4).NE.IBLANK) 60 TO 30               E   67
      IR(I,1)=0                                                          E   68
      GO TO 90                                                           E   69
   30 CONTINUE                                                           E   70
      IF (KR(I,5).NE.IBLANK) GO TO 35                                    E   71
      IF (KR(I,7).NE.IBLANK) KR(I,5)=KR(I,7)                             E   72
      IF (KR(I,7).NE.IBLANK) SC(I,5)=SC(I,7)                             E   73
      SC(I,7)»1.                                                         E   74
      KR(I,7HBLANK                                                     E   75
      IF (KR(I,5).NE.IBLANK) GO TO 35                                    E   76
      KR(I,5)=KR(I,6)                                                    E   77
      KR(I,6)=IBLANK                                                     E   78
      SC(I,5)=SC(I,6)                                                    E   79
      SC(I,6)=1.                                                         E   80
   35 K=KR(I,6)                                                          E   81
      IF (K.NE.IBLANK) GO TO 40                                          E   82
      KR(I,6)=KR(I,7)                                                    E   83
      KR(I,7)=K                                                          E   84
      SC(I,6)=SC(I,7)                                                    E   85
      SC(I,7)=1.                                                         E   86
   40 DO 85 J=l,7                                                        E   87
      K=KR(I,J)                                                          E   88
C                                                                        E   89
C     PROCESS REACTANTS HERE                                             E   90
C                                                                        E   91
      IF (J.GT.3) GO TO 65                                               E   92
C                                                                        E   93
C     ALL M DEPENDENT REACTIONS ARE TO HAVE A 99 IN THE THIRD SLOT       E   94
C                                                                        E   95
      IF (K.NE.MBLANK) GO TO 60                                          E   96
      GO TO  (45,50,55),J                                                 E   97
   45 KR(I,1)=KR(I,2)                                                    E   98
      SC(I,1)=SC(I,2)                                                    E   99
   50 KR(I,2)=KR(I,3)                                                    E  100
                                     A-17

-------
      SC(I,2)»SC(I,3)                                                    B  101
      SC(I,3)-1.                                                         E  102
      KR(I,3)sMBLANK                                                     E  103
   55 IR(I,3)=99                                                         E  104
   60 K=KR(I,J)                                                          E  105
C                                                                        E  106
C     ZERO ORDER REACTIONS HAVE 99 IN FIRST SLOT                         E  107
C                                                                        E  108
      IF (KR(U).EQ.IBLANK) IR(I,1)=99                                  E  109
      IF (KR(I.l).EQ.IBLANK) K*99                                        E  110
C                                                                        E  111
C     ALL BLANKS ARE SET EQUAL TO ZERO                                   E  112
C                                                                        E  113
      IF (J.EQ.3.AND.K.EQ.MBLANK) K=99                                   E  114
   65 IF (K.EQ.MBLANK.OR.K.EQ.IBLANK) IR(I,J)*0                          E  115
      IF (K.EQ.IBLANK.OR.K.EQ.99) GO TO 85                               E  116
      IF (NS.NE.O) GO TO 70                                              E  117
      NS*1                                                               E  118
      GO TO 80                                                           E  119
   70 DO 75 L»1,NS                                                       E  120
      IF (K.NE.SPECIS(L)) GOTO 75                                       E  121
C                                                                        E  122
C     SLOT SET TO SPECIES NUMBER                                         E  123
C                                                                        E  124
      IR(I,J)=L                                                          E  125
      GO TO 85                                                           E  126
   75 CONTINUE                                                           E  127
C                                                                        E  128
C     IF NO SPECIES ARE FOUND, ADD ONE TO THE LIST                       E  129
C                                                                        E  130
      IF (SPECIS(NS).NE.ITEM>) NS=NS+1                                   E  131
   80 SPECIS(NS)=K                                                       E  132
      C(NS+2)-C(NS+l)                                                    E  133
      C(NS+1)=C(NS)                                                      E  134
      C(NS)=0.                                                           E  135
      IR(I,J)-NS                                                         E  136
   85 CONTINUE                                                           E  137
   90 CONTINUE                                                           E  138
      IF (SPECIS(NS).NE.ITENP) NS«NS+1                                   E   139
      SPECIS(NS)=ITEff                                                   E  140
      SPECIS(NS+1)=MBLANK                                                E   141
      NS=NS+1                                                            E  142
      DO 95  IK=1,NR                                                      E   143
      DO 95 MT-1,3                                                       E   144
      J=IFIX(SC(IK,MT)-KJROUND)                                            E  145
      ISC(IK,MT)=J                                                       E   146
   95 IF (SC(IK,WT)-FLOAT(J).GT.4.*UROUND)  ISC(IK,MT)=-1                 E  147
      DO 100  1=1,NR                                                      E   148
      IF (IR(I.l).EO.O)  GOTO 100                                        E  149
      ITYPE(I)»2                                                          E  150
                                     A-18

-------
C                                                                        E  151
C     CHECK FOR 99 CODE AND SUBSTITUTE M WHICH IS THE  LAST SPECIES        E  152
C                                                                        E  153
      IF (IR(I,1).EQ.99.AND.IR(I,3).EQ.99)  IR(I,1)=NS                     E  154
      IF (IR(I,1).EQ.NS.AND.IR(I,3).EQ.99)  IR(I,3)=0                     E  155
      IF (IR(I,2).EQ.O.AND.IR(I,3).EQ.99)  IR(I,2)=NS                     E  156
      IF (IR(I,2).EQ.NS.AND.IR(I,3).EQ.99)  IR(I,3)=0                     E  157
      IF (IR(I,2).EQ.O) ITYPEdH                                       E  158
      IF (IR(I,3).NE.O) ITYPE(I)-3                                       E  159
      IF (IR(I,3).EQ.99) IR(I,3)=NS                                      E  160
  100 CONTINUE                                                           E  161
      RETURN                                                             E  162
      END                                                                E  163-
                                    A-19

-------
      SUBROUTINE PLOT (NTIT,NPNT,NTOT,NA«,SAVTIM,SAVCON,IOT,KCP,JGRID)
  SUBROUTINE ******  PLOT
C
C
C THIS SUBROUTINE READS THE PLOT CARDS AND PLOTS THE RESULTS AS PART
C OF THE PRINTED OUTPUT ~ IT DOES NOT DRIVE A PLOTTER.
C
C SYMBOL DESCRIPTIONS -
C
C CGRID
C CHIGH
C CLOW
C CSPAN
C DATA
C IDT
C J
C JBLANK
C XONC
C JGRID
C JN
C JPLUS
C JSTAR
C JSYMB
C JVERT
C JX
           THE LENGTH OF THE VERTICAL AXIS,  PPM
           HIGHEST CONCENTRATION VALUE,  PPM
           LOWEST CONCENTRATION VALUE, PPM
           CONCENTRATION NORMALIZATION FACTOR
           CONCENTRATION DATA POINTS, PPM, UP TO 80
           IF NOT ZERO THEN PRINT RAW DATA USED FOR PLOTS
           DO-LOOP INDICES OR LOCAL POINTERS
           A HOLLERITH WORD OF FOUR BLANK CHARACTERS
           CONCENTRATION LABELS
           THE PLOTTING GRID
           HOLLERITH SYMBOL FOR NO DATA POINTS
           THE CHARACTER >+>
           THE CHARACTER >*>
           SYMBOL TO BE USED FOR PLOTTING SAVED POINTS
           VERTICAL LEGEND
           THE CHARACTER >X>
C K        DO-LOOP INDICES OR LOCAL POINTERS
C KCP      IF NOT EQUAL TO ZERO SAVE DATA FOR OFFLINE PLOTTING
C KCON     CONCENTRATION COORDINATE ON GRID
C KTIM     TIME COORDINATE ON GRID
C L        DO-LOOP INDICES OR LOCAL POINTERS
C M        DO-LOOP INDICES OR LOCAL POINTERS
C MAXCON   LIMIT ON NUMBER OF VERTICAL POINTS
C MAXPNT   MAXIMUM NUMBER OF SAVED TIME AND CONCENTRATION POINTS
C MAXTIM   LIMIT ON NUMBER OF HORIZONTAL POINTS
C N        DO-LOOP INDICES OR LOCAL POINTERS
C NAME     SPECIES NAMES, ONE PER SPECIES
C NDAT     NUMBER OF CONCENTRATION DATA POINTS
C NIN      THE FORTRAN INPUT UNIT (NORMALLY 5)
C NOUT     THE FORTRAN OUTPUT UNIT NUMBER (NORMALLY 6)
C NPLT     THE NUMBER OF SPECIES TO BE PLOTTED
C NPNT     NUMBER OF SAVED TIMES AND CONCENTRATIONS
C NTEST    SPECIES NAME FOR TESTING
C NTIT     USER-INPUT TITLE FOR PRINTOUT, 3 FOUR-CHARACTER WORDS
C NTOT     TOTAL NUMBER OF SPECIES
C NX       FLAG FOR CORRECT SPECIES TEST
C SAVCON   SPECIES CONCENTRATIONS, PPM, ONE PER SPECIES AT 80 TIMES
C SAVTIM   TIMES THAT CONCENTRATIONS ARE SAVED, MIN, UP TO 80 VALUES
C TGRID    THE LENGTH OF THE HORIZONTAL AXIS, MIN
C THIGH    HIGHEST TIME VALUE, MIN
C TIME     TIMES AT WHICH CONCENTRATIONS ARE  INPUT, MIN, UP TO 80
C TLOW     LOWEST TIME VALUE, MIN
F
F
F
F
F
F
F '
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
F
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
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
                                     A-20

-------
C TPRINT   TIMES FOR PRINTOUT ON HORIZONTAL AXIS, MIN                    F   51
C TSPAN    TIME NORMALIZATION FACTOR                                     F   52
C                                                                        F   53
C BEGINNING OF PROGRAM.                                                  F   54
C                                                                        F   55
C ENTRY POINT                                                            F   56
C                                                                        F   57
C SET DIMENSIONS OF INCOMING ARRAYS                                      F   58
C                                                                        F   59
      DIMENSION SAVCON(90,80), SAVTIM(SO), NTIT(3), NAME(91)             F   60
C                                                                        F   61
C SET DIMENSIONS OF LOCAL ARRAYS                                         F   62
C                                                                        F   63
      DIMENSION XONC(5), TIME(80), DATA(80), JN(3)                      F   64
      DIMENSION JGRID(89,40), TPRINT(9), KVERT(52)                       F   65
C                                                                        F   66
C DEFINE THE VERTICAL LABEL AND ESTABLISH ALPHAMERIC DATA                F   67
C                                                                        F   68
      COMMON /APLOT/ JVERT(52,2),JBLANK,JSTAR,JPLUS,JBAR                 F   69
      COMMON /PHOTR/ IPP,RDAT(80),RTIM(80),RR1,IN10                      F   70
C                                                                        F   71
C DEFINE MISCELLANEOUS DATA VALUES                                       F   72
C                                                                        F   73
      COMMON /INOUT/ NIN,NOUT,ITAPE                                      F   74
      DATA MAXTIM/89/,MAXCON/40/,MAXPNT/80/                              F   75
      DATA TGRID/88./,CGRID/40./,JX/lHX/,JNO/lH /                        F   76
      DATA KVERT/11*4H    ,4H P  ,4H H  ,4H 0  ,4H T  ,4H 0  ,4H L  ,4H  F   77
     1Y  ,4H S  ,4H I  ,4H S  ,4H    ,4H C  ,4H 0  ,4H N  ,4H S  ,4H T   F   78
     2,4H A  ,4H N  ,4H T  ,22*4H    /                                   F   79
C                                                                        F   80
C DEFINE VERTICAL AXIS VIA ASSIGNMENT STATEMENTS                         F   81
C                                                                        F   82
      IF (IDT*IDT.EQ.O) GO TO 10                                         F   83
      WRITE (NOUT.150)                                                   F   84
      WRITE (NOUT,145) (NAME(I),I=1,NTOT)                                F   85
      DO 5 J=1,NPNT                                                      F   86
    5 WRITE (NOUT.155) SAVTIM(J),(SAVCON(I,J),I=1,NTOT)                  F   87
   10 DO 15 J=l,40                                                       F   88
      JGRID(1,J)=JBAR                                                    F   89
   15 CONTINUE                                                           F   90
      JGRID(1,1)=JPLUS                                                   F   91
      JGRID(1,11)=JPLUS                                                  F   92
      JGRID(1,21)=JPLUS                                                  F   93
      JGRID(1,31)=JPLUS                                                  F   94
      JN(1)=JSTAR                                                        F   95
      JN(2)=OPLUS                                                        F   96
      JN(3)=JX                                                           F   97
      TSPAN=-10.                                                         F   98
C                                                                        F   99
C READ PLOT CONTROL CARD                                                 F  100
                                     A-21

-------
C                                                                        F  101
   20 READ (NIN,130) NPLT,XONC,CLOW,CHIGH,TLOW,THIGH                    F  102
C                                                                        F  103
C TEST FOR END PLOTTING                                                  F  104
C                                                                        F  105
      IF (NPLT.EQ.O) GO TO 105                                           F  106
      WRITE (NOUT.160) NTIT                                              F  107
C                                                                        F  108
C SET NORMALIZATION FACTORS AND VERTICAL LABELS                          F  109
C                                                                        F  110
      CSPAN*CGRID/(CHIGH-CLOW)                                           F  111
      TSPAN*TGRID/(THIGH-TLOW)                                           F  112
      JVERT(1,2)=JCONC(5)                                                F  113
      JVERT(11,2)=OCONC(4)                                               F  114
      JVERT(21,2)=JCONC(3)                                               F  115
      JVERT(31,2)=JCONC(2)                                               F  116
C                                                                        F  117
C SET HORIZONTAL TIME LABELS                                             F  118
C                                                                        F  119
      DO 25 J=l,9                                                        F  120
      TPRINT(J)=FLOAT(J-1)/8.*(THIGH-TLOW)+TLOW                       .   F  121
   25 CONTINUE                                                           F  122
C                                                                        F  123
C CLEAR GRID                                                             F  124
C                                                                        F  125
      DO 35 K-1.MAXCON                                                   F  126
      DO 30 J*2,MAXTIM                                                   F  127
      JGRID(J,K)*OBLANK                                                  F  128
   30 CONTINUE                                                           F  129
   35 CONTINUE                                                           F  130
      NX=0                                                               F  131
      DO 95 LK'l.NPLT                                                    F  132
      READ (NIN.135) NTEST.NDAT.JSYMB                                    F  133
C                                                                        F  134
C TEST NUMBER OF DATA POINTS AND READ DATA                               F  135
C                                                                        F  136
      IF  (NDAT.LE.O) GO  TO 45                                            F  137
      IF  (NDAT.LE.MAXPNT) GO TO 40                                       F  138
      WRITE (NOUT,185) MAXPNT                                            F  139
      GO TO 125                                                          F  140
C                                                                        F  141
C READ DATA POINTS                                                       F  142
C                                                                        F  143
   40 READ (NIN.140)  (TIME(J),DATA(J),J=1,NDAT)                          F  144
C                                                                        F  145
C TEST FOR CORRECT  SPECIES NAME                                          F  146
C                                                                        F  147
   45 DO 50 L-l.NTOT                                                     F  148
      IF  (NTEST.EQ.NAME(L)) GO TO 55                                     F  149
   50 CONTINUE                                                           F  150
                                     A-22

-------
      WRITE (NOUT,190) NTEST                                             F  151
      NX=NX+1                                                            F  152
      IF (NPLT.EQ.1.0R.NX.EQ.3) GO TO 20                                 F  153
      IF (NPLT.EQ.2.AND.NX.EQ.2) GO TO 20                                F  154
      GO TO 95                                                           F  155
C                                                                        F  156
C IF THERE ARE CALCULATED POINTS, GET THEIR COORDINATES                  F  157
C                                                                        F  158
   55 IF (NPNT.LE.O) GO TO 65                                            F  159
      DO 60 J=1,NPNT                                                     F  160
      KTIM=IFIX((SAVTIM(J)-TLOW)*TSPAN+1.5)                              F  161
      Y=(SAVCON(L,J)-CLOW)*CSPAN-0.5                                     F  162
      IF (Y.LT.0.0) GO TO 60                                             F  163
      KCON=IFIX(Y)                                                       F  164
      KCON=MAXCON-KCON                                                   F  165
C                                                                        F  166
C CHECK FOR BEING WITHIN GRID, THEN PLACE ON GRID                        F  167
C                                                                        F  168
      IF (KTIM.LT.2) GO TO 60                                            F  169
      IF (KCON.LT.l) GO TO 60                                            F  170
      IF (KTIM.GT.MAXTIM) GO TO 60                                       F  171
      IF (KCON.GT.MAXCON) GO TO 60                                       F  172
      JGRID(KTIM,KCON)=JSYMB                                             F  173
   60 CONTINUE                                                           F  174
C                                                                        F  175
C IF THERE ARE DATA POINTS, GET THEIR COORDINATES                        F  176
C                                                                        F  177
   65 IF (NDAT.LE.O) GO TO 75                                            F  178
      DO 70 J=1,NDAT                                                     F  179
      KTIM=IFIX((TIME(J)-TLOW)*TSPAN+1.5)                                F  180
      Y=(DATA(J)-CLOW)*CSPAN-0.5                                         F  181
      IF (Y.LT.0.0) GO TO 70                                             F  182
      KCON=IFIX(Y)                                                       F  183
      KCON=MAXCON-KCON                                                   F  184
C                                                                        F  185
C CHECK FOR BEING WITHIN GRID, THEN PLACE ON GRID                        F  186
C                                                                        F  187
      IF (KTIM.LT.2) GO TO 70                                            F  188
      IF (KCON.LT.l) GO TO 70                                            F  189
      IF (KTIM.GT.MAXTIM) GO TO 70                                       F  190
      IF (KCON.GT.MAXCON) GO TO 70                                       F  191
      IF (LK.EQ.l) JGRID(KTIM,KCONHSTAR                                F  192
      IF (LK.EQ.2) JGRID(KTIM,KCON)=JPLUS                                F  193
      IF (LK.EQ.3) JGRID(KTIM,KCON)=JX                                   F  194
   70 CONTINUE                                                           F  195
   75 IF (NDAT.GT.O) 60 TO 85                                            F  196
      WRITE (NOUT,195) NTEST,JNO.JSYMB                                   F  197
      GO TO 90                                                           F  198
   85 WRITE (NOUT.195) NTEST,JN(LK),JSYMB                                F  199
C                                                                        F  200
                                     A-23

-------
C SAVE DATA FOR OFFLINE PLOTS IF DESIRED                                F  201
C                                                                       F  202
   90 IF (KCP*KCP.EQ.O) GO TO 95                                        F  203
      WRITE (ITAPE)  NTIT,NAME(L),aOW,CHIGH,TLOW,THI6H,NDAT,NPNT,(TIME(J F  204
     1),>1,NDAT),(DATA(J),>1,NDAT),(SAVTIM(J),>1,NPMT),(SAVCON(L,J),J F  205
     2=1,NPNT)                                                          F  206
   95 CONTINUE                                                          F  207
      DO 100 K=1,MAXCON                                                 F  208
      WRITE (NOUT.165)  JVERT(K,1)$JVERT(K,2),(JGRID(J,K),J»1,MAXTIM)     F  209
  100 CONTINUE                                                          F  210
C                                                                       F  211
C PRINT THE HORIZONTAL  AXIS AND LABELS                                  F  212
C                                                                       F  213
      WRITE (NOUT.170)  JCONC(l)                                         F  214
      IF (THIGH.GT.80.) WRITE (NOUT.175) TPRINT                         F  215
      IF (THIGH.LE.80.) WRITE (NOUT.180) TPRINT                         F  216
      GO TO 20                                                          F  217
C                                                                       F  218
C END OF SUBROUTINE ~ RETURN TO CALLER                                 F  219
C                                                                       F  220
  105 IF (IPP.LE.O)  RETURN                                              F  221
      IF (IN10.LE.O) RETURN                                             F  222
      IF (TSPAN.LE.O.)  RETURN                                           F  223
      OVERT(1,2HH0.80                                                 F  224
      JVERT(11,2)=4H0.60                                                F  225
      JVERT(21,2)*4H0.40                                                F  226
      JVERT(31,2)=4H0.20                                                F  227
      JCONCUHHO.OO                                                   F  228
      CSPAN=CGRID/0.8                                                   F  229
      WRITE (NOUT.200)  NTIT                                             F  230
      DO 110  I=1,MAXCON                                                 F  231
      DO 110  J=2,MAXTIM                                                 F  232
  110 JGRID(J,I)=JBLANK                                                 F  233
      DO 115  I=1,NPNT                                                   F  234
      KTIM=IFIX((RTIM(I)-TLOW)*TSPAN+1.5)                                F   235
      Y=(RDAT(I)-aOW)*CSPAN-0.5                                        F  236
      IF (Y.LT.Q.O) GO TO  115                                            F   237
      KCON'IFIX(Y)                                                       F   238
      KCON=MAXCON-KCON                                                  F  239
      IF (KCON.LT.l) GO TO 115                                           F   240
      IF (KTIM.LT.2) GO TO 115                                           F  241
      IF (KTIM.GT.MAXTIM)  GO TO 115                                      F  242
      IF (KCON.GT.MAXCON)  GO TO 115                                      F  243
      JGRID(KTIM,KCON)=JSTAR                                             F  244
  115 CONTINUE                                                          F  245
      DO 120  K=1,MAXCON                                                  F  246
      WRITE  (NOUT,165) KVERT(K),JVERT(K,2),(JGRID(J,K),>1,MAXTIM)       F  247
  120 CONTINUE                                                          F  248
      WRITE  (NOUT.170) XONC(l)                                          F  249
      IF (TPRINT(9).GT.80.) WRITE  (NOUT,175) TPRINT                      F  250
                                     A-24

-------
      IF (TPRINT(9).LE.80.) WRITE (NOUT,180)  TPRINT                      F  251
      RETURN                                                             F  252
  125 CONTINUE                                                           F  253
      STOP                                                               F  254
C                                                                        F  255
C LIST OF FORMAT STATEMENTS                                              F  256
C                                                                        F  257
C                                                                        F  258
C                                                                        F  259
  130 FORMAT (5X,I2,8X,5(A4,1X),4F10.0)                                   F  260
  135 FORMAT (A4,1X,I2,1X,A1)                                            F  261
  140 FORMAT (8F10.0)                                                    F  262
  145 FORMAT (/,4X,5HTIME ,4X,10(4X,A4,4X),/,(13X,10(4X,A4,4X)))         F  263
  150 FORMAT (1H1)                                                       F  264
  155 FORMAT (1P11E12.4,/,(12X,10E12.4))                                 F  265
  160 FORMAT (1H1,//////////,62X,3A4,/,32X,7HSPECIES,2X,5HEXPT.,2X,4HSIM F  266
     1.)                                                                 F  267
  165 FORMAT (16X,2A4,89A1)                                              F  268
  170 FORMAT (20X,A4,1H+,8(11H	+))                              F  269
  175 FORMAT (I25,2X,8F11.0,/,62X,14HTIME (MINUTES),/)                   F  270
  180 FORMAT (I25,2X,8F11.3,/,62X,14HTIME (MINUTES),/)                   F  271
  185 FORMAT (33H PROGRAM CANNOT HANDLE MORE THAN ,I4,28H PLOT POINTS - F  272
     1 JOB ABORTED.)                                                     F  273
  190 FORMAT (33X,13HSPECIES NAME ,A4,26H NOT FOUND IN SPECIES LIST)     F  274
  195 FORMAT (33X,A4,6X,A1,5X,A1)                                        F  275
  200 FORMAT (1H1,//////////,62X,3M,//)                                 F  276
      END                                                                F  277-
                                     A-25

-------
   SUBROUTINE VALU (VAL.N.NF.NL)                                      G    1
   COWON /STORE/  AST(35),IPL(7),TEI€ND,NTEMP,TMI,NPHOT,PHI,IL,NFRST, G    2
  1IPH(30),OM(100),PM(100),PSTOP                                      G    3
   CCWON /GEAR1/  T,H,FWIN,H1AX,EPS1,UROUND,NC,MF1,KFLAG1,JSTART      G    4
   DIMENSION AC(5),  AD(5),  BC(5)                                      G    5
   INTEGER AC,AD,BC,AST                                              G    6
   DATA IBLK/1H /,IPER/1H./,IZRO/1HO/                                 G    7
   AVAL'ALOGIO(VAL)                                                   G    8
   IF (AVAL.EQ.O.) RETURN                                            G    9
   DO 5 1=1,5                                                        G   10
   ACUHBtK                                                        G   11
   AD(I)=IZRO                                                        G   12
 5 BCUHZRO                                                        G   13
   IAD=IFIX(AVAL+UROUND)                                              G   14
   IF (IAD) 30,10,10                                                 G   15
10 IREM=3-IAD                                                        G   16
   IAD=IAD+1                                                         G   17
   JREM=IFIX(VAL-HJROUND)                            .                  G   18
   REM-VAL-FLOAT(JREM)                                               G   19
   IF (REM.GT.UROUND.AND.IREM.GT.O)  GO TO 15                         G   20
   CALL CONVT (JREM,AC,5)                                            G   21
   GO TO 50                                                          G   22
15 CALL CONVT (JREM,AD,IAD)                                          G   23
   JREM=IFIX(REM*(10.**IREM)+0.1)                                     G   24
   CALL CONVT (JREM,BC,IREM)                                         G   25
   DO 20 J=1,IAD                                                     G   26
20 AC(J)=AD(J)                                                       G   27
   AC(IAD+1)=IPER                                                    G   28
   DO 25 K=1,IREM                                                    G   29
25 AC(K+IAD+1)=8C(K)                                                 G   30
   GO TO 40                                                          G   31
30 IAD=IABS(IFIX(AVAL-0.1))-1                                        G   32
   IVAL=IFIX(VAL*10000.+0.1)                        '                 G   33
   CALL CONVT (IVAL,AC,5)                                            G   34
   AC(1)=IPER                                                        G   35
   DO 35 1=1,IAD                                                     G   36
35 ACU+1HZRO                                                      G   37
40 IF (AC(S).NE.IZRO) GO TO 50                                       G   38
   DO 45 K-1,4                                                        G   39
   L=6-K                                                             G   40
45 AC(L)=AC(L-1)                                                      G   41
   AC(1}=IBLK                                                         G   42
   60 TO 40                                                           G   43
50 DO 55  I=NF,NL                                                      G   44
   K-I-NF+1                                                           G   45
55 AST(I)=AC(K)                                                       G   46
   RETURN  •                                                           G   47
   END                                                                G   48-
                                  A-26

-------
      SUBROUTINE CONVT (NUM.L.N)                                         H    1
C                                                                        H    2
C  SUBROUTINE CONVT CONVERTS INTEGERS TO ALPHANUMERICS                   H    3
C     FOR PRINTING                                                       H    4
C                                                                        H    5
C  ASSUMES VALUE OF INTEGER IS POSITIVE                                  H    6
C                                                                        H    7
      DIMENSION L(5), JDIGIT(IO)                                         H    8
C                                                                        H    9
      DATA JDIGIT/1HO,1H1,1H2,1H3,1H4,1H5,1H6,1H7,1H8,1H9/               H   10
      DATA JBLANK/1H /                                                   H   11
C                                                                        H   12
      NI=NUM                                                             H   13
      DO 5 M,N                                              -           H   14
      L(I)=JBLANK                                                        H   15
    5 CONTINUE                                                           H   16
C                                                                        H   17
      DO 10 K=1,N                                                        H   18
      I=N-K+1                                                            H   19
      NEXT=NI/10                                                         H   20
      NDX=(NI-NEXT*10)+1                                                 H   21
      l(I)=JDIGIT(NDX)                                                   H   22
      IF (NEXT.LE.O) GO TO 15                                            H   23
      NI=NEXT                                                            H   24
   10 CONTINUE                                                           H   25
C                                                                        H   26
   15 RETURN                                                             H   27
      END                                                                H   28-
                                     A-27

-------
   SUBROUTINE SAVLIN  (T,C,N)                                          I    1
   COWON/FRPLOT/  NIT(3),SAVCON(90,80),SA\rriM(80),JGRID(89,40),NT    I    2
   CCW3N /PHOTR/ IPP,RDAT(80),RTIM(80),RR1,IN10                      I    3
   DDCNSION C(91)                                                    I    4
   DATA NFRST/1/                                                      I    5
   IF (NFRST.EQ.l)  60 TO 5                                            16
   IF (T.EQ.TOLD) RETURN                                              I    7
 5 IF (NT.LT.O) RETURN                                                I    8
   NFRST*2                                                           I    9
   TOLD=T                                                            I   10
   IF (NT.6T.80) RETURN                                              I   11
   IF (T.LT.SAVTIM(NT)) RETURN                                        I   12
   SAVTIM(NT)-T                                                     -  I   13
   RTIM(NT)=T                                                        I   14
   DO 10 M,N                                                       I   15
10 SAVCON(I,NTK(I)                                                  I   16
   RDAT(NT)=RR1                                                       .1   17
   NT*NT+1                                                           I   18
   RETURN                                                            I   19
   END                                                               I   20-
                                  A-28

-------


c
c
c
c
c
c
c
c
c
c
c
c
c



c
c
c
c
c




c
c
c
c
c
c
c
c




c
c
c
c







SUBROUTINE SPLNA (N,X,Y,J,D,C,W)
DIMENSION X(10), Y(10), 0(2), C(30), W(30)

OVER THE INTERVAL X(I) TO X(I+1), THE INTERPOLATING
POLYNOMIAL
Y=Y(I)+A(I)*Z+B(I)*Z**2+E(I)*Z**3
WHERE Z=(X-X(I))/(X(W)-X(I))
IS USED. THE COEFFICIENTS A(I),B(I) AND E{I) ARE CONFUTED
BY SPLNA AND STORED IN LOCATIONS C(3*I-2),C(3*I-1) AND
C(3*I) RESPECTIVELY.
WHILE WORKING IN THE ITH INTERVALJHE VARIABLE Q WILL
REPRESENT Q=X(I+1) - X(I), AND Y(I) WILL REPRESENT
Y(H-1)-Y(I)


Q=X(2)-X(1)
YI=Y(2)-Y(1)
IF (J.EQ.2) GO TO 5

IF THE FIRST DERIVATIVE AT THE END POINTS IS GIVEN,
A(l) IS KNOWN, AND THE SECOND EQUATION BECOMES
MERELY B(1)+E(1)=YI - 0*D(1).

C(1)=0*D(D
C(2)=1.0
W(2)=YI-C(1)
GO TO 10

IF THE SECOND DERIVATIVE AT THE END POINTS IS GIVEN
B(l) IS KNOWN, THE SECOND EQUATION BECOMES
A(1)+E(1)=YI-0.5*0*Q*D(1). DURING THE SOLUTION OF
THE 3N-4 EQUATIONS,A1 WILL BE KEPT IN CELL C(2)
INSTEAD OF C(l) TO RETAIN THE TRIDIAGONAL FORM OF THE
COEFFICIENT MATRIX.

5 C(2)=0.0
W(2) =0.5*Q*Q*D(1)
10 M=N-2
IF (M.LE.O) GO TO 20

UPPER TRIANGULARIZATION OF THE TRIDIAGONAL SYSTEM OF
EQUATIONS FOR THE COEFFICIENT MATRIX FOLLOWS—

DO 15 I=1,M
AI=Q
Q=X(I+2)-X(I+l)
H=AI/Q
C(3*I)=-H/(2.0-C(3*I-1))
W(3*I)=(-YI-W(3*I-1))/(2.0-C(3*I-1))
C(3*I+1)=-H*H/(H-C(3*I))
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
J
• J
J
J
J
J
J
J
J
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
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
A-29

-------
r
i«
C
C
C
r
r
C
C
r
v»
r
L.
C
C
C
C
r
i*
W(3*I+1MYI-W(3*I))/(H-C(3*I))
YI»Y(I+2)-Y(I+l)
C(3*I+2)»1.0/(1.0-C(3*I+1))
15 W(3*I+2MYI-W(3*I+1))/(1.0-C(3*H-1))
E(N-l) IS DETERMINED DIRECTLY FROM THE LAST EQUATION
OBTAINED ABOVE, AND THE FIRST OR SECOND DERIVATIVE
VALUE GIVEN AT THE END POINT.
20 IF (J.EQ.l) 60 TO 25
C(3*N-3MQ*Q*0(2)/2.0-W(3*N-4))/(3.0-C(3*N-4))
GO TO 30
25 C(3*N-3MQ*D(2)-YI-W(3*N-4))/(2.0-C(3*N-4))
30 M-3*N-6
IF (M.LE.O) GO TO 40
BACK SOLUTION FOR ALL COEFFICENTS EXCEPT
A(l) AND B(l) FOLLOWS--
DO 35 II=1,M
T=M-I 1+3
35 C(I)-M(I)-C(I)*C(W)
40 IF (J.EQ.l) GO TO 45
IF THE SECOND DERIVATIVE IS GIVEN AT THE END POINTS,
A(l) CAN NOW BE COMPUTED FROM THE KNOWN VALUES OF
8(1) AND E(l). THEN A(l) AND B(l) ARE PUT INTO THEIR
PROPER PLACES IN THE C ARRAY.
C(1H(2)-Y(1)-W(2)-C(3)
C(2)-W(2)
RETURN
45 C(2)=W(2)-C(3)
RETURN
END
J
J
J
J
.— .1
J
J
J
... J
J
J
J
J
J
J
.— .1
J
J
J
J
J
J
... 1
J
J
J
J
... J
J
J
J
J
J
J
51
52
53
54
EC
99
56
57
58
CO
•jy
60
61
62
63
64
65
fifi
oo
67
68
cq
UJ
70
71
72
73
74
/H
75
76
77
78
7Q
fy
80
81
82
83
84
85-
A-30

-------
   SUBROUTINE  DRIVES  (N,TO,HO,YO,TOUT,EPS,MF,INDEX,IA.JA)             K    1
   DIMENSION IA(1), JA(1), YO(N)                                      K    2
   DIMENSION Y(91,6)                                                  K    3
   COMMON /6EAR1/ T,H,HMIN,WAX,EPSC,UROUND,NC,MFC,KFLAG,JSTART       K    4
   COMMON /GEAR2/ YMAX(100)/GEAR3/ERROR(100)/GEAR4/W1(90,3)           K    5
   COMMON /GEARS/ IW1(91,9)/GEAR6/W2(2000)/GEAR7/IW2(2000)            K    6
   COVMON /GEARS/ EPSJ,IPTI2,IPTI3,IPTI4,IPTR2,IPTR3,NGRP             K    7
   COMMON /GEAR9/ HUSED,NQUSED,NSTEP,NFE,NJE,NZA,NPL,NPU,NZL,NZU,NZRO K    8
   COMMON /INOUT/ INP,LOUT,ITAPE                                      K    9
   DATA  NMX/90/,LENW2/2000/,LENIW2/2000/                              K   10
   NGP=0                                                             K   11
   IF  (INDEX.EQ.4) GO TO  15                                           K   12
   IF  (INDEX.EQ.O) GO TO  30                                           K   13
   IF  (INDEX.EQ.2) GO TO  35                                           K   14
   IF  (INDEX.EQ.-l) GO  TO 40                                         K   15
   IF  (INDEX.EQ.3) GO TO  45                                           K   16
   IF  (INDEX.NE.l) GO TO  135                                         K   17
   IF  (EPS.LE.O.) GO  TO 120                                           K   18
   IF  (N.LE.O) GO TO  125                                              K   19
   IF  ((TO-TOUT)*HO.GE.O.) GO  TO  130                                  K   20
   MITER=MF-10*(MF/10)                                                K   21
   IF  ((MITER.NE.1).AND.(MITER.NE.2))  GO  TO 15                        K   22
   NP1=N+1                                                            K   23
   NZA=IA(NP1)-1                                                      K   24
   MAX=LENIW2/2                                                      K   25
   IPTI2=MAX+1                                                       K   26
   CALL  SORDER (N,IA,JA,IW1,IW1(1,5),MAX,IW2,IW2(IPTI2),IER)          K   27
   IPTI2=NZA+1                                                       K   28
   IF  (IPTI2+NZA-1.GT.LENIW2)  GOTO 145                              K   29
   DO  5  I=1,NP1                                                      K   30
 5 IW1(I,2)=IA(I)                                                     K   31
   DO  10 1=1,NZA                                                      K   32
10 IW2(I)=JA(I)                                                      K   33
   CALL  NSCORD (N,IW1(1,2),IW2,IW1(1,3),IW2(IPTI2),IW1,IW1(1,5),IW1(1 K   34
  1,8))                                                               K   35
   MAXPL=(LENIW2-NZA)/2                                              K   36
   IPTI3=IPTI2+MAXPL                                                  K   37
   MAXPU=LENIW2-IPTI3+1                                              K   38
   CALL NSSFAC (N,IW1(1,2),IW2,MAXPL,IW1(1,3),IW2(IPTI2),IW1(1,4),MAX K   39
  lPU,IWl(l,5),IW2(IPTI3),IWl(l,6),Y(lt6),IWl(l,9),Y,Y(l,2),Y(l,3),IW K   40
  21(l,7),IWl(l,8)fY(l,4),Y(l,5),IER)                                 K   41
   NPL«M(N,4)                                                      K   42
   NPU=IW1(N,6)                                                      K   43
   NZL=M(N+1,3)                                                     K   44
   NZU=IW1(N+1,5)                                                     K   45
   IPTR2=NZA+1                                                       K   46
   IPTR3=IPTR2+MAXO(NZA,NZL)                                         K   47
   IF  (IPTR3+MAXO(NZA,NZU)-1.6T.LENW2) GO TO 145                     K   48
15 DO 20 1=1,N                                                       K   49
   YMAX(I)=ABS(YO(I))                                                K   50
                                 A-31

-------
   IF (YMAX(I).EQ.O.)  YMAX(I)=1.E-10                                  K   51
20 Y(I,1)=YO(I)                                                       K   52
   NC»N                                                              K   53
   T=TO                                                              K   54
   H-HO                                                              K   55
   NZROO                                                            K   56
   TST=EPS*1.E-10                                                    K   57
   00 25 1=1,N                                                       K   58
25 IF (Y(I.l).GT.TST)  NZRO-NZROH                                     K   59
   NZRO=MAXO(NZRO,1)                                                  K   60
   NOLD=NZRO                                                         K   61
   HMIN=ABS(HO)                                                       K   62
   HMAX=ABS(TO-TOUT)*10.                                              K   63
   HWX=AMIN1(WAX,20.)                                               K   64
   EPSC=£PS                                                          K   65
   J^OMF                                                     ••      K   66
   JSTARTX)                                                          K   67
   NO=N                                                              K   68
   NMX1=NMX+1                                                        K   69
   EPSJ=SOJ^T(UROUND)                                                  K   70
   NHCUT=0                                                           K   71
   GO TO 50                                                          K   72
30 H1AX=ABS(TOUT-TOUTP)*10.                                           K   73
   FWAX*AMIN1(WAX,20.)                                               K   74
   GO TO 80                                                          K   75
35 HMAX=ABS(TOUT-TOUTP)*10.                                           K   76
   HMAX=AMIN1(HMAX,20.)                                               K   77
   IF ((T-TOUT)*H.GE.O.)  GO TO 150                                   K   78
   GO TO 85                                                          K   79
40 IF ((T-TOUT)*H.GE.O.)  GO TO 140                                   K   80
   JSTART*-!                                                         K   81
   NC*N                                                              K   82
   EPSC=EPS                                                          K   83
   ^FC=MF                                                            K   84
45 CONTINUE                                                          K   85
50 CALL STIFFS (Y,NO,IA,JA,W1,NMX,IW1,NMX1)                           K   86
   KGO=1-KFLAG                                                       K   87
   GO TO (55,95,110,100),KGO                                         K   88
55 CONTINUE                                                          K   89
   D=0.                                                              K   90
   NZRO=0                                                            K   91
   DO 70 1*1,NC                                                      K   92
   IF  (Y(I.l).GE.O.) GOTO 65                                        K   93
   NGP=NGP+1                                                         K   94
   DO 60 J=l,6                                                       K   95
   K=(J-1)*N+I                                                       K   96
60 Y(K,1)=0.                                                         K   97
65 CONTINUE                                                          K   98
   IF  (Y(I,1).GT.TST) NZRO=NZRO+1                                    K   99
   AYI=ABS(Y(I,1))                                                   K   100
                                  A-32

-------
    YMAX(I)=AMAX1(1.E-10,AYI)                                           K  101
 70 D=D+(AYI/YMAX(I))**2                                               K  102
    NZRO=MAXO(NZRO,1)                                                   K  103
    IF (NZRO.NE.NOLD)  JSTART=-1                                         K  104
    D=D*(UROUND/EPS)**2                                                K  105
    DO 75 II-l.NC                                                      K  106
 75 YO(II)=Y(II,1)                                                      K  107
    CALL SAVLIN (T.YO.N)                                               K  108
    IF (D.GT.FLOAT(N)) GO TO 115                                       K  109
    IF (INDEX.EQ.3)  GO TO 150                                          K  110
    IF (INDEX.EQ.2)  GO TO 85                                           K  111
 80 IF ((T-TOUT)*H.LT.O.) GO TO 45                                     K  112
    CALL INTERP (TOUT,Y,NO,YO)                                          K  113
    GO TO 160                                                          K  114
 85 IF (T.GE.TOUT)  GO TO 90                                            K  115
    IF (((T+H)-TOUT).LE.O.) GO  TO 45                                   K  116
    H=(TOUT-T)*(1.+4.*UROUND)                                           K  117
    JSTART=-1                                                          K  118
    GO TO 45                                                           K  119
 90 JSTART=-1                                                          K  120
    H=AMIN1(H,1.)                                                      K  121
    GO TO 150                                                          K  122
 95 CONTINUE                                                           K  123
100 IF (NHCUT.EQ.10) GO TO 105                                          K  124
    NHCUT=NHCUT+1                                                      K  125
    HMIN=.1*HMIN                                                       K  126
    H=.1*H                                                             K  127
    JSTART=-1                                                          K  128
    GO TO 45                                                           K  129
105 WRITE (LOUT.165)                                                   K  130
    IF (KGO.EQ.4) WRITE (LOUT,180) T                                   K  131
    STOP                                                               K  132
110 WRITE (LOUT.170) T,H                                               K  133
    STOP                                                               K  134
115 WRITE (LOUT.175) T                                                 K  135
    KFLAG=-2                                                           K  136
    STOP                                                               K  137
120 WRITE (LOUT,185)                                                   K  138
    STOP                                                               K  139
125 WRITE (LOUT,190)                                                   K  140
    STOP                                                               K  141
130 WRITE (LOUT,195)                                                   K  142
    STOP                                                               K  143
135 WRITE (LOUT,200) INDEX                                             K  144
    STOP                                                               K  145
140 WRITE (LOUT.205) T,TOUT,H                                          K  146
    STOP                                                               K  147
145 WRITE (LOUT.210)                                                   K  148
    STOP                                                               K  149
150 TOUT=T                                                             K  150
                                   A-33

-------
      DO 155 M,N                                                       K  151
  155 YO(I)=Y(I,1)                                                       K  152
      CALL SAVLIN (TOUT,YO,N)                                            K  153
  160 INDEX=KFLAG                                                        K  154
      TOUTP=TOUT                                                         K  155
      HO=HUSED                                                           K  156
      IF (KFLAG.NE.O) HO=H               .                                K  157
      RETURN                                                             K  158
C                                                                        K  159
C                                                                        K  160
  165 FORMAT (//44H PROBLEM APPEARS UNSOLVABLE WITH GIVEN INPUT//)       K  161
  170 FORMAT (//35H KFLAG « -2 FROM INTEGRATOR AT T » .E16.8.5H  H =,E16 K  162
     1.8/52H  THE REQUESTED ERROR IS SMALLER THAN CAN BE HANDLED//)      K  163
  175 FORMAT (//37H INTEGRATION HALTED BY DRIVER AT T = ,E16.8/56H  EPS  K  164
     1TOO SMALL TO BE ATTAINED FOR THE MACHINE PRECISION/)               K  165
  180 FORMAT (//35H KFLAG = -3 FROM INTEGRATOR AT T » .E16.8/45H  CORREC K  166
     1TOR CONVERGENCE COULD NOT BE ACHIEVED/)                            K  167
  185 FORMAT (//28H ILLEGAL INPUT.. EPS  .LE. O.//)                       K  168
  190 FORMAT (//25H ILLEGAL INPUT.. N .LE. O//)                          K  169
  195 FORMAT (//36H ILLEGAL INPUT.. (TO-TOUT)*H  .GE. O.//)               K  170
  200 FORMAT (//24H ILLEGAL INPUT.. INDEX =,I5//)                        K  171
  205 FORMAT (//44H INDEX = -1 ON INPUT WITH (T-TOUT)*H .GE. 0./4H T =,E K  172
     116.8,9H   TOUT =,E16.8,6H   H *,E16.8)                             K  173
  210 FORMAT (//42H INSUFFICIENT WORKING STORAGE IN IW2 OR W2//)         K  174
      END                                                                K  175-
                                     A-34

-------
   SUBROUTINE  STIFFS  (Y,NO,IA,JA,W1,NMX,IW1,NMX1)                      L     1
   COMMON  /GEAR1/ T,H,HMIN,HMAX,EPS,UROUND,N,MF,KFLAG,JSTART           L     2
   COMMON  /GEAR2/ YMAX(100)/GEAR3/ERROR(100)                           L     3
   COMMON  /GEARS/ W2(2000)/GEAR7/IW2(2000)                             L     4
   COMMON  /GEARS/ EPSO,IPTI2,IPTI3,IPTI4,IPTR2,IPTR3,NGRP              L     5
   COMMON  /GEAR9/ HUSED,NQUSED,NSTEP,NFE,NJE,IDUMMY(5),NZRO            L     6
   COMMON  /DATA/ NR,KR(200,7),A(200),S(200),ITITLE(7),TEM),ERR,START,  L     7
  1STOPP,PC,NP,SIG(91),IP(91),ITYPE(200),R(200),BK,SG,DILUT            L     8
   COMMON  /HEAT/ CV,Q(200),SC(200,7),ISC(200,3),ITEMP                  L     9
   DIMENSION Y(NO,1),  IA(1),  JA(1),  W1(NMX,1),  IWl(NMXl.l)             L    10
   DIMENSION EL(13),  TQ(4), RT(3)                                      L    11
   DATA EL(2)/1./,OLDLO/1./                                            L    12
   KFLAG=0                                                            L    13
   TOLD=T                                                              L    14
   IF (JSTART.GT.O) GO TO 50                                           L    15
   IF (JSTART.NE.O) GO TO 10                                           L    16
   CALL DIFFUN (N,T,Y,W1)                                              L    17
   DO 5 I-l.N                                                          L    18
 5 Y(I,2)=H*W1(I,1)                                                    L    19
   METH=MF/10                                                          L    20
   MITER=MF-10*METH                                                    L    21
   NQ=1                                                               L    22
   L=2                                                                L    23
   IDOUB=3                                                  .         L    24
   RMAX=1.E4                                                           L    25
   RC=0.                                                              L    26
   CRATE=1.                                                            L    27
   HOLD=H                                                              L    28
   MFOLD=MF                                                            L    29
   NSTEP=0                                                            L    30
   NSTEPJ=0                                                            L    31
   NFE=1                                                              L    32
   NJE=0                                                              L    33
   IRET=3                                                              L    34
   GO TO 15                                                            L    35
10 IF (MF.EQ.MFOLD) GO TO 25                                           L    36
   MEO=METH                                                            L    37
   MIO=MITER                                                           L    38
   METH=MF/10                                                          L    39
   MITER=MF-10*METH                                                    L    40
   MFOLD=MF                                                            L    41
   IF (MITER.NE.MIO)  IWEVAL=MITER                                      L    42
   IF (METH.EQ.MEO) GO TO 25                                           L    43
   IDOUB=L+1                                                           L    44
   IRET-1                                                              L    45
15 CALL COSET  (METH,NQ,EL,TQ,MAXDER)                                  L    46
   LMAX=MAXDER+1                                                       L    47
   RC=RC*EL(1)/OLDLO                                                   L    48
   OLDLO=EL(1)                                                        L    49
20 FN=FLOAT(NZRO)                                                      L    50
                                 A-35

-------
   EDN=FN*(TQ(1)*EPS)**2                                              L   51
   E=FN*(TQ(2)*EPS)**2                                                L   52
   EUP=FN*(TQ(3)*EPS)**2                                              L   53
   BND=FN*(TQ(4)*EPS)**2                                              L   54
   EPSOLD-EPS                                                         L   55
   NOLD=NZRO                                                          L   56
   GOTO (30,35,50),IRET                                              L   57
25 IF (EPS.EQ.EPSOLD.AND.NZRO.EQ.NOLD) GO TO 30                       L   58
   IRET=1                                                             L   59
   GO TO 20                                                           L   60
30 IF (H.EQ.HOLD) GO TO 50                                            L   61
   RH=H/HOLD                                                          L   62
   H=HOLD                                                             L   63
   IREDO=3                                                            L   64
   GO TO 40                                                           L   65
35 RH=AMAX1(RH,HMIN/ABS(H))                                           L   66
40 RH=AMIN1(RH,HMAX/ABS(H),RMAX)                                      L   67
   Rl-1.                                                              L   68
   DO 45 J=2,L                                                        L   69
   R1=R1*RH                                                           L   70
   DO 45 I-l.N                                                        L   71
45 Y(I,J)=Y(I,J)*R1                                                   L   72
   H=H*RH                                                             L   73
   RC=RC*RH    .                                                       L   74
   IDOUB-L+1                                                          L   75
   IF (IREDO.EQ.O) GO TO 290                                          L   76
50 IF (ABS(RC-1.).GT.0.3) IWEVAL41ITER                                L   77
   IF (NSTEP.GE.NSTEPJ*20) IWEVAL*MITER                  .             L   78
   T=T+H                                                              L   79
   DO 55 J1-1.NQ                                                      L   80
   DO 55 J2=J1,NQ                                                     L   81
   J»(NO*J1)-J2                                                       L   82
   DO 55 I-l.N                                                        L   83
55 Y(1,J)-Y(I,J)+Y(I,M)                                             L   84
60 DO 65 1=1,N                                                        L   85
65 ERROR(I)=0.                                                        L   86
   M=0                                                                L   87
   CALL DIFFUN (N,T,Y,W1(1,2))                                        L   88
   NFE=NFE+1                                                          L   89
   IF (IWEVAL.LE.O) GO TO 140                                         L   90
   IWEVAW)                                                           L   91
   RC=1.                                                              L   92
   NJE-NJE+1                                                          L   93
   NSTEPJ=NSTEP                                                       L   94
   CON=-H*EL(1)                                                       L   95
   ISV=M                                                              L   96
   LSV=L                                                              L   97
   NZ=IA(Nfl)-l                                                       L   98
   DO 70 I=1,NZ                                                       L   99
70 W2(I)=0.                                                           L  100
                                 A-36

-------
    DO 125 IR=1,NR                                                     L  101
    IF (KR(IR,1).EQ.O.OR.KR(IR,1).EQ.99)  GO TO 125                      L  102
    MT=ITYPE(IR)                                                        L  103
    DO 85 1=1,W                                                        L  104
    IF (KR(IR,I).EQ.N+2)  60 TO 85                                      L  105
    JX=I+l-I/3*3                                                        L  106
    KX=I+2-I/2*3                                                        L  107
    J=KR(IR,JX)                                                         L  108
    K=KR(IR,KX)                                                         L  109
    XJ=1.                                                              L  110
    IF (J.EQ.O) GO TO 75                                               L  111
    XJ=Y(J,1)                                                           L  112
    IF (J.EQ.N+2)  XJ=SG                                                L  113
 75 XK=1.                                                              L  114
    IF (K.EQ.O) GO TO 80                                               L  115
    XK=Y(K,1)                                                           L  116
    IF (K.EQ.N+2)  XK=SG                                                L  117
 80 RT(I)=R(IR)*XJ*XK                                                  L  118
 85 CONTINUE                                                           L  119
    DO 120 K=1,MT                                                      L  120
    I=KR(IR,K)                                                          L  121
    IF (I.EQ.N+2)  GO TO 120                                            L  122
    DO 100 L«1,MT                                                      L  123
    J=KR(IR,L)                                                          L  124
    IF (J.EQ.N+2)  GO TO 100                                            L  125
    M=IA(J)-1                                                           L  126
 90 M=M+1                                                              L  127
    IF (I-JA(M))  90,95,90                                              L  128
 95 W2(M)=W2(M)-RT(L)        "                                         L  129
100 CONTINUE                                                           L  130
    DO 115 L=4,7                                                        L  131
    J=KR(IR,L)                                                          L  132
    M=IA(I)-1                                                           L  133
    IF (J) 115,120,105                                                 L  134
105 M=M+1                                                              L  135
    IF (J-JA(M))  105,110,105                                           L  136
110 W2(M)=W2(M)+RT(K)*SC(IR,L)                                         L  137
115 CONTINUE                                                           L  138
120 CONTINUE                                                           L  139
125 CONTINUE                                                           L  140
    DO 135 J-l.N                                                        L  141
    KMIN=IA(J)                                                          L  142
    KMAX=IA(J+1)-1                                                     L  143
    DO 130 K=KMIN,KMAX                                                 L  144
    W2(K)=W2(K)*CON                                                    L  145
    IF (JA(K).EQ.J) W2(K)=W2(K)+1.-CON*DIIUT                           L  146
130 CONTINUE                                                           L  147
135 CONTINUE                                                           L  148
    CALL NSCORA (N,IA,JA,W2,IW1(1,2),W2(IPTR3),W2(IPTR2),IW1,IW1(1,7), L  149
   1IW1(1,8))                                                           L  150
                                  A-37

-------
    CALL NSNFAC (N,IW1(1,2),IW2,W2,IW1(1,3),IW2(IPTI2),IW1(1,4),W2(IPT L  151
   lR2),Wl(l,3),IWl(lt5),IW2(IPTI3),IWl(l,6),W2(IPTR3),Wl,IWl(l,7),M L  152
   2(1,8),IER)                                                         L  153
    M*ISV                                                              L  154
    L=LSV                                                              L  155
    IF (IER.NE.O)  60 TO 160                                            L  156
140 DO 145 I-l.N                                                       L  157
    IF (M.LE.O) 60 TO 145                                              L  158
    IF (-H*W1(I,2)*10..6T.Y(I,1))  60 TO 155                            L  159
145 W1(I,1)»H*W1(I,2)-(Y(I,2)+ERROR(I))                                L  160
    CALL NSBSLV (N,IWl,IWl,M(l,3),IW2(IPTI2),IWl(li4),W2(IPTR2),Wl(l L  161
   1,3),IW1(1,5),IW2(IPTI3),IW1(1,6),W2(IPTR3),W1(1,2),W1,W2)           L  162
    D*0.                                                               L  163
    DO 150 1=1,N                                                       L  164
    ERROR(I)=ERROR(I)+W1(I,2)                                          L  165
    D=D+(W1(I,2)/YMAX(I))**2                                           L  166
150 W1(I,1H(I,1)+EL(1)*ERROR(I)                                       L  167
    IF (M.NE.O) CRATE=AMAX1(.9*CRATE,D/D1)                             L  168
    IF ((D*AMIN1(1.,2.*CRATE)).LE.BND) 60 TO 175                       L  169
    Dl=0                                                               L  170
    M=MH                                                              L  171
    IF (M.EQ.3) 60 TO 155                                              L  172
    CALL DIFRJN (N,T,W1,W1(1,2))                                       L  173
    GO TO 140                                                          L  174
155 NFE=NFE+2                                                          L  175
    IF (IWEVAL.EQ.-l) GO TO 170                                        L  176
160T*TOLD                                                             L  177
    RMAX=2.                                                            L  178
    DO 165 J1=1,NQ                                                     L  179
    DO 165 J2=J1,NQ                  .                                  L  180
    J=(NQfJl)-J2                                                       L  181
    DO 165 I-1.N                                                       L  182
165 Y(I,J)«Y(I,J)-Y(I,J+1)                                             L  183
    IF (ABS(H).LE.HMIN*1.00001) GO TO 285                              L  184
    RH».25                                                  '          L  185
    IREDO=1                                                            L  186
    GO TO 35                                                           L  187
170 IWEVAL=MITER                                                       L  188
    GO TO 60                                                           L  189
175 IF  (MITER.NE.O)  IWEVAL=-1                                          L  190
    NFE=NFE4«                                                          L  191
    D«0.                                                               L  192
    DO 180 I=1,N                                                       L  193
180 WHERROR(I)/YMAX(I))**2                                          L  194
    IF  (D.GT.E) GO TO  195                                              L  195
    KFLAG=0                                                            L  196
    IREDOO                                                             L  197
    NSTEP=NSTEP+1                                                      L  198
    HUSED=H                                                             L  199
    NQUSED=NQ                                                          L  200
                                   A-38

-------
    DO 185 J=1,L                                                       L  201
    DO 185 1=1,N                                                       L  202
185 Y(I,J)=Y(I,J)+EL(J)*ERROR(I)                                       L  203
    IF (IDOUB.EQ.l) 60 TO 205                                          L  204
    IDOUB=IDOUB-1                                                      L  205
    IF (IDOUB.GT.l) GO TO 295                                          L  206
    IF (L.EQ.LMAX) GO TO 295                                           L  207
    DO 190 1=1,N                                                       L  208
190 Y(I,LMAX)=ERROR(I)                                                 L  209
    GO TO 295                                                          L  210
195 KFLAG=KFLAG-1                                                      L  211
    T=TOLD                                                             L  212
    DO 200 J1=1,NQ                                                     L  213
    DO 200 J2=J1,NQ                                                    L  214
    J=(NCH-J1)-J2                                                       L  215
    DO 200 I=1,N                                                       L  216
200 Y(I,J)=Y(I,J)-Y(I,J+1)                                             L  217
    RMAX=2.                                                            L  218
    IF (ABS(H).LE.HMIN*!.00001) GO TO 275                              L  219
    IF (KFLAG.LE.-3) GO TO 265                                         L  220
    IREDO=2                                                            L  221
    PR3=1.E+20                                                         L  222
    GO TO 215                                                          L  223
205 PR3=1.E+20                                                         L  224
    IF (L.EQ.LMAX) GO TO 215                                           L  225
    Dl=0.                                                              L  226
    DO 210 1=1,N                                                       L  227
210 D1=D1+((ERROR(I)-Y(I,LMAX))/YMAX(I))**2                            L  228
    ENQ3=.5/FLOAT(L+1)                                                 L  229
    PR3=((Dl/EUP)**ENQ3)*1.4+1.4E-6                                    L  230
215 ENQ2=.5/FLOAT(L)                                                   L  231
    PR2=((D/E)**ENQ2)*1.2+1.2E-6                                       L  232
    PR1=1.E+20                                                         L  233
    IF (NQ.EQ.l) GO TO 225                                             L  234
    D=0.                                                               L  235
    DO 220 1=1,N                                                       L  236
220 D=D+(Y(I,L)/YMAX(I))**2                                            L  237
    ENQ1=.5/FLOAT(NQ)                                                  L  238
    PRH(D/EDN)**ENQl)*1.3+1.3E-6                                     L  239
225 IF (PR2.LE.PR3) GO TO 230                                          L  240
    IF (PR3.LT.PR1) GO TO 240                                          L  241
    GO TO 235                                                          L  242
230 IF (PR2.GT.PR1) GO TO 235                                          L  243
    NEWQ=NQ                                                            L  244
    RH=1./PR2                                                          L  245
    GO TO 255                                                          L  246
235 NEWQ=NQ-1                                                          L  247
    RH=1./PR1                                                          L  248
    GO TO 255                                                          L  249
240 NEWQ=L                                                             L  250
                                   A-39

-------
    RH=1./PR3                                                          L  251
    IF (RH.LT.1.1)  GO TO 250                                           L  252
    DO 245 I-l.N                                                       L  253
245 Y(I,NEWO+1KRROR(I)*EL(L)/FLOAT(L)                                L  254
    GO TO 260                                                          L  255
250 IDOUB=10                                                           L  256
    GO TO 295                                                          L  257
255 IF ((KFLAG.EQ.O).AND.(RH.LT.1.1)) GO TO 250                        L  258
    IF (NEWQ.EQ.NQ) GO TO 35                                           L  259
260 NQ=NEWQ                                                            L  260
    L-NQfl                                                             L  261
    IRET=2                                                             L  262
    GO TO 15                                                           L  263
265 IF (KFLA6.EQ.-9) GO TO 280                                         L  264
    RH=10.**KFLAG                                                      L  265
    RH=AMAX1(HMIN/ABS(H),RH)                                           L  266
    H=H*RH                                                             L  267
    CALL DIFFUN (N,T,Y,W1)                                             L  268
    NFE=NFE+1                                                          L  269
    00 270 1=1,N                                                       L  270
270 Y(I,2)=H*W1(I,1)                                                   L  271
    IWEVAL*MITER                                                       L  272
    IDOUB=10                                                           L  273
    IF (NQ.EQ.l) GO TO 50                                              L  274
    NQ=1                                                               L  275
    L*2                                                                L  276
    IRET=3                                                             L  277
    GO TO 15                                                           L  278
275 KFLAG=-1                                                           L  279
    GO TO 295                                                          L  280
280 KFLAG=-2                                                           L  281
    GO TO 295                                                          L  282
285 KFLAG=-3                                                           L  283
    GO TO 295                                                          L  284
290 RMAX=100.                                                          L  285
295 HOLD=H                                                             L  286
    JSTART=NQ                                                          L  287
    RETURN                                                             L  288
    END                                                                L  289-
                                   A-40

-------
   SUBROUTINE COSET (METH,NQ,EL,TQ,MAXDER)                             M    1
   DIMENSION PERTST(12,2,3),  EL(13),  TQ(4)                             M    2
   DATA PERTST/1.,1.,2.,1...3158,.07407,.01391,.002182,.0002945,.0000 M    3
  13492,.000003692,.0000003524,1.,1.,.5,.1667,.04167,1.,1.,1.,1.,1.,1 M    4
  2.,1.,2.,12.,24.,37.89,53.33,70.08,87.97,106.9,126.7,147.4,168.8,19 M    5
  31.0,2.0,4.5,7.333,10.42,13.7,1.,!.,!.,1.,1.,1.,1.,12.0,24.0,37.89, M    6
  453.33,70.08,87.97,106.9,126.7,147.4,168.8,191.0,1.,3.0,6.0,9.167,1 M    7
  52.5,1.,!.,!.,!.,!.,!.,!.,I./                                       M    8
 5 MAXDER=5                                                           M    9
   GOTO (10,15,20,25,30),NQ                                          M   10
10 ELUH.O                                                          M   11
   GO TO 35                                                           M   12
15 EL(1)=6.6666666666667E-01                                          M   13
   EL(3)=3.3333333333333E-01                                          M   14
   GO TO 35                                                           M   15
20 EL(1)=5.4545454545455E-01                                          M   16
   EL(3HL(1)                                                        M   17
   EL(4)=9.0909090909091E-02                                          M   18
   GO TO 35                                                           M   19
25 EL(1)=0.48                                                         M   20
   EL(3)=0.7                                                          M   21
   EL(4)=0.2                                                          M   22
   EL(5)=0.02                                                         M   23
   GO TO 35                                                           M   24
30 EL(1)=4.3795620437956E-01                                          M   25
   EL(3)=8.2116788321168E-01                                          M   26
   EL(4)=3.1021897810219E-01                                          M   27
   EL(5)=5.4744525547445E-02                                          M   28
   EL(6)=3.6496350364964E-03                                          M   29
35 DO 40 K*l,3                                                        M   30
40 TQ(K)=PERTST(NQ,METH,K)                                            M   31
   TQ(4)=.5*TQ(2)/FLOAT(NO+2)                                         M   32
   RETURN                                                             M   33
   END                                                                M   34-
                                  A-41

-------
   SUBROUTINE NSCORA (N,IA,JA,A,IAP,JAWORK,AWORK,C,IR,ICT)             N    1
   INTEGER IA(1),JA(1),IAP(1),JAWORK(1),C(1),IR(1),ICT(1)              N    2
   REAL A(1),AWORK(1)                                                  N    3
   00 5 K»1,N                                                         N    4
   I(X<(K)                                                           N    5
 5 IR(ICK)=K                                                          N    6
   JMIN-1                                                             N    7
   DO 15 K=1,N                                                        N    8
   ICK-C(K)                                                           N    9
   JW\X=JMIN+IA(ICK+1)-IA(ICK)-1                                      N   10
   IF (JMIN.GT.JMAX) SO TO 15                                         N   11
   IAINK=IA(ICK)-1                                                    N   12
   DO 10 J=JMIN,JMAX                                                  N   13
   IAINK»IAINK+1                                                      N   14
   JAOUTJ=JA(IAINK)                                                   N   15
   JAOUTJ=IR(JAOUTJ)                                                  N   16
   JAWORK(J)=JAOUTJ                                                   N   17
10 AWORK(J)=A(IAINK)                                                  N   18
15 JMIKWMAX-H                                                        N   19
   DO 20 I-l.N                                                        N   20
20 ICT(I)«IAP(I)                                                      N   21
   JMIN-1                                                             N   22
   DO 30 I-l.N                                                        N   23
   ICK-C(I)                                                           N   24
   JMflX«JMIM*IA(IOC+l)-IA(ICK)-l                                      N   25
   IF (JMIN.GT.JMAX) GO TO 30                                         N   26
   DO 25 J*JMIN,JMAX                                                  N   27
   JAOUTJ=JAWORK(J)                                                   N   28
   ICTJ=ICT(JAOUTJ)                                                   N   29
   A(ICTJ)=AWORK(J)                                                   N   30
   ICT(JAOUTJ)=ICTJ+1                                                 N   31
25 CONTINUE                                                           N   32
30 JMIN=JMAX+1                                                        N   33
   RETURN                                                             N   34
   END                                                                N   35-
                                  A-42

-------
   SUBROUTINE NSNFAC (N,IA,JA,A,IL,JL,ISL,L,D,IU,JU,ISU,U,X,IRL,JRL,I 0    1
  1ER)                                                                 0    2
   INTEGER IA(1),JA(1),IL(1),JL(1),ISL(1)                              0    3
   INTEGER IU(1),JU(1),ISU(1),IRL(1),JRL(1)                            0    4
   REAL A(1),L(1),D(1),U(1),X(1)                                       0    5
   REAL LKI                                                           06
   IER=0                                                              0    7
   00 5 K=1,N                                                         0    8
   IRL(K)=IL(K)                                                        0    9
 5 JRL(K)=0                                                           0   10
   DO 90 K=1,N                                                        0   11
   X(K)=0.                                                            0   12
   11=0                                                               0   13
   IF (JRL(K).EQ.O)  GO TO 15                                          0   14
   I=JRL(K)                                                           0   15
10 I2=JRL(I)                                                          0   16
   JRL(I)=I1                                                          0   17
   11=1                                                               0   18
   X(I)=0.                                                            0   19
   1=12                                                               0   20
   IF (I.NE.O) GO TO 10                                               0   21
15 JMIN=ISU(K)                                                        0   22
   JMAX=JMIN+IU(K+1)-IU(K)-1                                          0   23
   IF (JMIN.GT.JMAX) GO TO 25                                         0   24
   DO 20 J=JMIN,JMAX                                                  0   25
   JUJ=JU(J)                                                          0   26
20 X(JUJ)=0.                                                          0   27
25 JMIN=IA(K)                                                         0   28
   JMAX=IA(K+1)-1                                                     0   29
   DO 30 J=JMIN,JMAX                                                  0   30
   JAJ=JA(J)                                                          0   31
30 X(JAJ)=A(J)                                                        0   32
   1=11                                                               0   33
   IF (I.EQ.O) GO TO 50                                               0   34
35 IRLI=IRL(I)                                                        0   35
   LKI=-X(I)                                                          0   36
   L(IRLI)=-LKI                                                       0   37
   JMIN=IU(I)                                                         0   38
   JMAX=IU(I+1)-1                                                     0   39
   IF (JMIN.GT.JMAX) GO TO 45                                         0   40
   ISUB=ISU(I)-1                                                      0   41
   DO 40 J=JMIN,JMAX                                                  0   42
   ISUB=ISUB+1                                                        0   43
   JUJ=JU(ISUB)                                                       0   44
40 X(JUJ)=X(JUJ)+LKI*U(J)                                             0   45
45 I=JRL(I)                                                           0   46
   IF (I.NE.O) GO TO 35                                               0   47
50 IF (X(K).EQ.O.) GO TO 95                                           0   48
   DK=1./X(K)                                                         0   49
   D(K)=DK                                     •                       0   50
                                  A-43

-------
   IF (K.EQ.N)  GO TO 90                                               0   51
   JMIN-IU(K)                                                          0   52
   JMAX=IU(K+1)-1                                                     0   53
   IF (JMIN.GT.JMAX) GO TO 60                                         0   54
   ISUB=ISU(K)-1                                                      0   55
   DO 55 J=JMIN,JMAX                                                  0   56
   ISUB=ISUB+1                                                        0   57
   JUJ-JU(ISUB)                                                       0   58
55 U(JH(JUO)*DK                                                     0   59
60 CONTINUE                                                           0   60
   1*11                                                               0   61
   IF (I.EQ.O) GO TO 85                                               0   62
65 IRL(I)-IRL(IH                                                    0   63
   I1=JRL(I)                                                          0   64
   IF (IRL(I).GE.IL(I+1)) GOTO 80                                    0   65
   ISLB=IRL(I)-IL(I)+ISL(I)                                           0   66
   J-JL(ISLB)                                                         0   67
70 IF (I.GT.JRL(J)) GO TO 75                                          0   68
   J=JRL(J)                                                           0   69
   GO TO 70                                                           0   70
75 JRL(I)=JRL(J)                                                      0   71
   JRL(J)=I        .                                                   0   72
80 1=11                                                               0   73
   IF (I.NE.O) GO TO 65                                               0   74
85 ISLK=ISL(K)                                                        0   75
   IF (IRL(K).GE.IL(K-H)) GO TO 90                                    0   76
   >JL(ISLK)                                                         0   77
   JRL(K)=JRL(J)                                                      0   78
   JRL(J)=K                                                           0   79
90 CONTINUE                                                           0   80
   RETURN                                                             0   81
95 IER=K                                                              0   82
   RETURN                                                             0   83
   END                                                                0   84-
                                  A-44

-------
   SUBROUTINE NSBSLV (N,R,C,IL,JL,ISL,L,D,IU,JU,ISU,U,X,B,Y)           P    1
   DIMENSION R(l),  IL(1),  JL(1),  IU(1),  JU(1),  C(l),  ISL(l),  ISU(l)    P    2
   DIMENSION L(l),  X(l),  B(l),  U(l),  Y(l),  D(l)                       P    3
   INTEGER R,RK,C,CK                                                  P    4
   REAL L                                                             P    5
   DO 5 K=1,N                                                         P    6
   RK=R(K)                                                            P    7
 5 Y(K)=B(RK)                                                         P    8
   DO 15 K=1,N                                                        P    9
   JMIN=IL(K)                                                         P   10
   JMAX=IL(K+1)-1                                                     P   11
   YK=-D(K)*Y(K)                                                       P   12
   Y(K)=-YK                                                           P   13
   IF (JMIN.GT.J1AX) GO TO 15                                         P   14
   ISLB=ISL(K)-1                                                       P   15
   DO 10 J=JMIN,JMAX                                                  P   16
   ISLB=ISLB+1                                                        P   17
   JLJ=JL(ISLB)                                                       P   18
10 Y(JLJ)=Y(JLJ)+YK*L(J)                                               P   19
15 CONTINUE                                                           P   20
   K=N                                                                P   21
   DO 30 1=1,N                                                        P   22
   SUM=-Y(K)                                                          P   23
   JMIN=IU(K)                                                         P   24
   JMAX=IU(K+1)-1                                                     P   25
   IF (JMIN.GT.JMAX) GO TO 25                                         P   26
   ISUB=ISU(K)-1                                                      P   27
   DO 20 J=JMIN,JMAX                                                  P   28
   ISUB=ISUB+1                                                        P   29
   JUJ=JU(ISUB)                                                       P   30
20 SUM=SUM+U(J)*Y(JUJ)                                                P   31
25 Y(K)=-SUM                                                          P   32
   CK=C(K)                                                            P   33
   X(CK)=-SUM                                                         P   34
30 K=K-1                                                              P   35
   RETURN                                                             P   36
   END                                                                P   37-
                                  A-45

-------
      SUBROUTINE YSMER (A,K,A1)                                          Q    1
      INTEGER A,A1(2)                                                    Q    2
      COWON /INOUT/ INP.LOUT.ITAPE                                      Q    3
      WRITE (LOUT.5) A,K,A1(1),A1(2)                                     Q    4
      RETURN                                                            Q    5
C                                                                       Q    6
C                                                                       Q    7
C                                                                       08
    5 FORMAT (1X,A10,I6,2A10)                                            Q    9
      END                                                               Q   10-
                                    A-46

-------
   SUBROUTINE INTERP (TOUT,Y,NO,YO)                                    R    1
   COVWON /GEAR1/ T,H,DUMMY(4),N,IDUMMY(2),JSTART                     R    2
   DIMENSION YO(NO), Y(NO,1)                                           R    3
   DO 5 I=1,N                                                         R    4
 5 YO(I)=Y(I,1)                                                        R    5
   L=JSTART+1                                                         R    6
   S=(TOUT-T)/H                                                        R    7
   Sl=l.                                                               R    8
   DO 15 J=2,L                                                        R    9
   S1=S1*S                                                            R   10
   DO 10 1=1,N                                                        R   11
10 YO(I)=YO(I)+S1*Y(I,J)                                              R   12
15 CONTINUE                                                           R   13
   RETURN                                                             R   14
   END                                                                R   15-
                                  A-47

-------
   SUBROUTINE SPARS (IA,JA,N)                                          S    1
   COMMON /DATA/  NRtKR(200,7),A(200),S(200),ITITLE(7),TEMP,ERR,START,  S    2
  1STOPP,PC,NP,SI6(91),IP(91),ITYPE(200),R(200),BK,S6,DILUT           S    3
   DIMENSION IA(N), JA(N)                                              S    4
   DO 5 I«1,M                                                         S    5
 5 IA(IH                                                            S    6
   KT=0                                                               S    7
   IA(N+1)=1                                                          S    8
   JA(1)=0                                                            S    9
   DO 70 IR»1,NR                                                       S   10
   IF (KR(IR,l).EQ.N+2)  GO TO 70                                      S   11
   IF (KR(IR,1).EQ.O.OR.KR(IR,1).EQ.99)  GO TO 70                       S   12
   MT=ITYPE(IR)                                                        S   13
   DO 65 K=1,MT                                                        S   14
   I*KR(IR,K)                                                         S   15
   IF (KR(IR,K).EQ.N+2)  GO TO 65                                      S   16
   DO 30 L=1,MT                                                        S   17
   J=KR(IR,L)                                                         S   18
   IF (KR(IR,L).EQ.N+2)  GO TO 30                                      S   19
   K1=IA(J)                                                           S   20
   K2-IA(J*1)-1                                                        S   21
   IF (K1.GT.K2)  GO TO 15                                             S   22
   DO 10 MHC1.K2                                                      S   23
10 IF (I.EQ.JA(M))  GO TO 30                                           S   24
15 DO 20 M=J,N                                                        S   25
20 IA(M*1)«IA(MH)+1                                                  S   26
   KT=KT+1                                                            S   27
   KD=KT-K2                                                           S   28
   K2*K2+1                                                            S   29
   DO 25 M«1,KD                                                       S   30
25 JA(KT+2-M)=JA(KT+l-M)                                              S   31
   JA(K2)=I                                                           S   32
30 CONTINUE                                                           S   33
   Kl-IA(I)                                                           S   34
   DO 60 1=4,7                                                        S   35
   K2-IA(W)-1                                                       S   36
   J=KR(IR,L)                                                         S   37
   IF (KR(IR,L).EQ.N+2)  GO TO 60                                      S   38
   IF (J) 60,65,35                                                    S   39
35 IF (K1.GT.K2) GO TO 45                                             S   40
   DO 40 M=K1,K2                                                      S   41
   IF (J.EQ.JA(M)) GO TO 60                                           S   42
40 CONTINUE                                                           S   43
45 DO 50 M=I,N                                                        S   44
50 IA(M*1)-IA(M*1)+1                                                  S   45
   KT=KT+1                                                            S   46
   KD=KT-K2                                                           S   47
   K2=K2+1                                                            S   48
   DO 55 M-l.KD                                                       S   49
55 JA(KT+2-M)=JA(KT+l-M)                                              S   50
                                  A-48

-------
   JA(K2)=J                                                           S   51
60 CONTINUE                                                           S   52
65 CONTINUE                                                           S   53
70 CONTINUE                                                           S   54
   DO 80 I-l.N                                                        S   55
   K1=IA(I)+1                                                         S   56
   K2=IA(I+1)-1                                                       S   57
   IF (K1.GT.K2) GO TO 80                                             S   58
   MT=K2-K1+1                                                         S   59
   DO 75 K=1,MT                                                       S   60
   DO 75 M=K1,K2                                                      S   61
   IF (JA(M).GT.JA(M-l)) GO TO 75                                     S   62
   J=JA(M-1)                                                          S   63
   JA(M-1)=JA(M)                                                      S   64
   JA(M)=J                                                            S   65
75 CONTINUE                                                           S   66
80 CONTINUE                                                           S   67
   M=N                                                                S   68
   DO 95 I=1,M                                                        S   69
   IF (lA(I-H).GT.IA(I)) GO TO 95                                     S   70
   NM=I+1                                                             S   71
   NN=N+1                                                             S   72
   KMIN=IA(NM)                                                        S   73
   KMAX=IA(NN)                                                        S   74
   DO 85 J=KMIN,KMAX                                                  S   75
   KM=KMAX+KMIN-J                                                     S   76
85 JA(KM)=JA(KM-1)                                                    S   77
   KNOW=IA(I)                                                         S   78
   JA(KNOW)=I                                                         S   79
   DO 90 LL=NM,NN                                                     S   80
90 IA(LL)=IA(LL)+1                                                    S   81
95 CONTINUE                                                           S   82
   RETURN                                                             S   83
   END                                                                S   84-
                                  A-49

-------
   SUBROUTINE SORDER (N,IA,JA,P,Q,MAX,V,L,IER)                        T    1
   INTEGER IA(1),JA(1),P(1),Q(1),V(1),L(1)                            T    2
   INTEGER S,SFS,PI,PJ,PK,VI,VJ,VK,QVK,DTHR,DMIN                      T    3
   IER*0                                                             T    4
   DO 5 S=1,MAX                                                      T    5
 5 L(S)-S*1                                                          T    6
   SFS=1                                                             T    7
   L(MAX)0                                                          T    8
   DO 10 K-l.N                                                       T    9
   P(K)=K                                                            T   10
   Q(K}*                                                            T   11
   V(K)=1                                                            T   12
10 L(K)=0                                                            T   13
   SFS=SFS+N                                                         T   14
   DO 50 K=1,N                                                       T   15
   JMIN-IA(K)                                                        T   16
   JMAX=IA(K+1)-1                                                    T   17
   IF (JMIN.GT.JMAX+1)  GO TO 145                                     T   18
   KDIAG=0                                                           T   19
   DO 45 J=JMIN,JMAX                                                 T   20
   VJ=JA(J)                                                          T   21
   IF (VJ.NE.K) GO TO 15  -                                          T   22
   KDIAG=1                                                           T   23
   GO TO 45                                                          T   24
15 LLK=K                                                             T   25
20 LK=LLK                                                            T   26
   LLX=L(LK)                                                         T   27
   IF (LLK.EQ.O) GO TO 25                                            T   28
   IF (V(LLK)-VJ) 20,30,25                                           T   29
25 LLK-SFS                                                           T   30
   IF (LLK.EQ.O) GO TO 150                                           T   31
   SFS<(SFS)                                                        T   32
   V(K)=V(KH                                                       T  33
   V(LLK)=VJ                                                         T   34
   L(LLKK(LK)                                                      T  35
   L(LK)=LLK                                                         T   36
30 UK=VJ                                                            T  37
35 LK=LLK                                                            T   38
   LLK=L(LK)                                                         T   39
   IF (LLK.EQ.O) GO TO 40                                            T  40
   IF (V(LLK)-K) 35,45,40                                            T  41
40 LLK=SFS                                                           T  42
   IF (LLK.EQ.O) GO TO 150                                            T   43
   SFS=L(SFS)                                                         T  44
   V(VJH(VJ)+1                                                      T   45
   V(LLK)=K                                                           T   46
   L(LLK)=L(LK)                                                       T   47
   L(LK)=LLK                                                          T   48
45 CONTINUE                                                           T   49
   IF (KDIAG.EQ.O) GO TO 160                                          T   50
                                  A-50

-------
 50 CONTINUE                                                           T   51
    J=0                                                                 T   52
    DTHR=0                                                             T   53
    DMIN=N                                                             T   54
    1=0                                                                 T   55
 55 1=1+1                                                               T   56
    IF (I.GT.N)  GO TO  140                                              T   57
    JMIN=MAXO(J+1,I)                                                    T   58
    IF (JMIN.GT.N)  GO  TO 70                                            T   59
 60 DO 65  J=JMIN,N                                                     T   60
    VI=P(J)                                                             T   61
    IF (V(VI).LE.DTHR)  GO TO 75                                         T   62
    IF (V(VI).LT.DMIN)  DMIN=V(VI)                                       T   63
 65 CONTINUE                                                           T   64
 70 DTHR=DMIN                                                          T   65
    DMIN=N                                                             T   66
    JMIN=I                                                             T   67
    GO TO  60                                                           T   68
 75 PJ=P(I)                                                             T   69
    P(J)=PJ                                                             T   70
    Q(PJ)=J                                                             T   71
    PI=VI                                                               T   72
    P(I)=PI                                                             T   73
    Q(PI)=I                                                             T   74
    LI=VI                                                               T   75
 80 LI=L(LI)                                                           T   76
    IF (LI.EQ.O)  GO TO 105                                             T   77
    VK=V(LI)                                                           T   78
    LLK=VK                                                             T   79
    LJ=VI                                                               T   80
 85 LJ=L(LJ)                                                           T   81
    IF (LJ.EQ.O)  GO TO 100                                             T   82
    VJ=V(LJ)                                                           T   83
    IF (VJ.EQ.VK) GO TO 85                                             T   84
 90 LK=LLK                                                             T   85
    LLK=L(LK)                                                          T   86
    IF (LLK.EQ.O) GO TO 95                                             T   87
    IF (V(LLK)-VJ)  90,85,95                                            T   88
 95 LLK=SFS                                                             T   89
    IF (LLK.EQ.O) GO TO 155                                            T   90
    SFS=L(SFS)                                                          T   91
    V(VK)=V(VK)+1                                                      T   92
    V(LLK)=VJ                                                          T   93
    L(LLK)=L(LK)                                                        T   94
    L(LK)=LLK                                                          T   95
    GO TO  85                                                           T   96
100 IF (V(VK).GT.V(VI)) GO TO 80                                       T   97
    1=1+1                                                               T   98
    QVK=Q(VK)                                                          T   99
    PI=P(I)                                                             T  100
                                  A-51

-------
    P(QVK)*PI                                                          T  101
    Q(PI)*QVK                                                          T  102
    P(I)=VK                                                            T  103
    Q(VK)=I                                                            T  104
    60 TO 80                                                           T  105
105 LI-VI                                                              T  106
110 IF (L(LI).EQ.O) GO TO 135                                          T  107
    LI-L(LI)                                                           T  108
    VK=V(LI)                                                           T  109
    LLK=VK                                                             T  110
    QVK=MINO(Q(VK),I)                                                  T  111
115 LK=LLK                                                             T  112
    LLK=L(LK)                                                          T  113
    IF (LLK.EQ.O) GO TO 120                                            T  114
    VJ=V(LLK)                                •                          T  115
    IF (Q(VO).GT.QVK) GO TO 115                                        T  116
    V(VK)=V(VK)-1                                                      T  117
    L(LKK(LLK)                                                       T  118
    L(LLK)=SFS                                                         T  119
    SFS=LLK                                                            T  120
    LLK=LK                                                             T  121
    GO TO 115                                                          T  122
120 IF (Q(VK).LE.I) GO TO 130                                          T  123
    IF (V(VK).LE.DTHR) GO TO 125                                       T  124
    IF ((DTHR.LT.V(VK)).AND.(V(VK).LT.DMIN)) DMIN=V(VK)                T  125
    GO TO 110                                                          T  126
125 J=MINO(Q(VK)-1,J)                                                  T  127
    DTHR=V(VK)                                                         T  128
    GO TO 110                                                          T  129
130 L(LK)=SFS                                                          T  130
    SFS=L(VK)                                                      '    T  131
    GO TO 110                                                          T  132
135 L(LI)=SFS                                                          T  133
    SFS=L(VI)                                                          T  134
    GO TO 55                                                           T  135
140 RETURN                                                             T  136
145 CALL YSMER  (3HROW,K,13H OF A IS NULL)                              T  137
    GO TO 165                                                          T  138
150 CALL YSMER  (3HROW,K,16H EXCEEDS STORAGE)                           T  139
    GO TO 165                                                          T  140
155 CALL YSMER  (6HVERTEX.VI.16H EXCEEDS STORAGE)                       T  141
    GO TO 165                                                          T  142
160 CALL YSMER  (6HCOLUMN,K,19H.. DIAGONAL MISSING)                     T   143
165 IER=1                                                              T  144
    RETURN                                                             T   145
    END                                                                T  146-
                                   A-52

-------
   SUBROUTINE NSSFAC (N,IA,JA,MAXPL,IL,JL,ISL,MAXPU,IU,JU,ISU,P,V,IRA U    1
  1,JRA,IRAC,IRL,JRL,IRU,JRU,IER)                                      U    2
   INTEGER IA(1),JA(1),IRA(1),JRA(1),IL(1),JL(1),ISL(1)                U    3
   INTEGER IU(1),JU(1),ISU(1),IRL(1),JRL(1),IRU(1),JRU(1)              U    4
   INTEGER P(1),V(1),IRAC(1)                                           U    5
   INTEGER VI,VJ,VK,PK,PPK,PI,CEND,REND                               U    6
   IER=0                                                              U    7
   DO 5 K=1,N                                                         U    8
   IRAC(K)=0                                                          U    9
   JRA(K)=0                                                           U   10
   JRL(K)=0                                                           U   11
 5 JRU(K)=0                                                           U   12
   DO 10 K=1,N                                                        U   13
   IAK=IA(K)                                                          U   14
   IF (IAK.GT.IA(K+1))  GO TO 200                                      U   15
   IF (JA(IAK).GT.K) GO TO 10                                         U   16
   JAIAK=JA(IAK)      •                                                U   17
   JRA(K)=IRAC(JAIAK)                                                 U   18
   IRAC(JAIAK)=K                                                      U   19
10 IRA(K)=IAK                                                         U   20
   JLPTR=0                                                            U   21
   IL(1)=1                                                            U   22
   JUPTR=0                                                            U   23
   IU(1)=1                                                            U   24
   DO 195 K=1,N                                                       U   25
   P(l)=l                                                             U   26
   V(1)=N+1                                                           U   27
   LSFS=2                                                             U   28
   VJ=IRAC(K)                                                         U   29
   IF (VJ.EQ.O) GO TO 30                                              U   30
15 PPK=1                                                              U   31
20 PK=PPK                                                             U   32
   PPK=P(PK)                                                          U   33
   IF (V(PPK)-VJ) 20,220,25                                           U   34
25 P(PK)=LSFS                                                         U   35
   V(LSFS)=VJ                                                         U   36
   P(LSFS)=PPK                                                        U   37
   LSFS=LSFS+1                                                        U   38
   VJ=JRA(VJ)                                                         U   39
   IF (VJ.NE.O) GO TO 15                                              U   40
30 LASTI=0                                                            U   41
   I=K                                                                U   42
35 I=JRU(I)                                                           U   43
   IF (I.EQ.O) GO TO 60                                               U   44
   PPK=1                                                              U   45
   JMIN=IRL(I)                                                        U   46
   JMAX=ISL(I)+IL(I+1)-IL(I)-1                                        U   47
   IF (LASTI.GT.I) GO TO 40                                           U   48
   LASTI=I                                                            U   49
   LASTID=JMAX-JMIN                                                   U   50
                                  A-53

-------
   IF (JL(JMIN).NE.K) LASTID=tASTICH-l                                 U   51
40 IF (JMIN.6T.JW) SO TO 35                                         U   52
   DO 55 J=JMIN, J1AX                                                  U   53
   VJ=OL(J)                                                           U   54
45 PK=PPK                                                             U   55
   PPK=P(PK)                                                          U   56
   IF (V(PPK)-VJ) 45,55,50                                            U   57
50 P(PK)=LSFS                                                         U   58
   V(LSFSHJ                                                         U   59
   P(LSFS)=PPK                                                        U   60
   PPK=LSFS                                                           U   61
   LSFS-LSFS+1                                                        U   62
55 CONTINUE                                                           U   63
   60 TO 35                                                           U   64
60 PI-P(l)                                                            U   65
   IF (V(PI).NE.K) 60 TO 225                                          U   66
   IF (LASTI.EQ.O) 60 TO 65                                           U   67
   IF (LASTID.NE.LSFS-3) 60 TO 65                                     U   68
   IRLL=IRL(LASTI)                                                    U   69
   ISL(K)»IRLL+1                                                      U   70
   IF (JL(IRLL).NE.K) ISL(K)-ISL(K)-l                                 U   71
   IL(K-H)«IL(K)-HASnD                                               U   72
   IRL(KHSL(K)                                                      U   73
   60 TO 80                                                           U   74
65 ISL(K)=JLPTR+1                                                   .  U   75
   PI-P(l)                                                            U   76
   PI-P(PI)                                                           U   77
   VI-V(PI)                                                           U   78
70 IF (VI.GT.N) 60 TO 75                                              U   79
   JLPTR=JLPTR+1                                                      U   80
   IF (JLPTR.6T.MAXPL) 60 TO 230                                      U   81
   JL(JLPTR)*VI                                                       U   82
   PI-P(PI)                                                           U   83
   VI-V(PI)                                                           U   84
   60 TO 70                                                           U   85
75 IRL(K)«ISL(K)                                                      U   86
   IL(K+1)-IL(K)+OLPTR-ISL(1C)+1                                       U   87
80 P(l)=l                                                             U   88
   V(1)=N+1                                                           U   89
   LSFS=2                                                             U   90
   JMIN*IRA(K)                                                        U   91
   JMAX»IA(K+1)-1                                                     U   92
   IF (JMIN.6T.JMAX) 60 TO  100                                         U   93
   DO 95 J=OMIN,JMAX                                                  U   94
   VJ=JA(J)                                                            U   95
   PPK=1                                                               U   96
85 PK=PPK                                                              U   97
   PPK=P(PK)                                                           U   98
   IF (V(PPK)-VJ) 85,205,90                                           U   99
90 P(PK)=LSFS                                                          U   100
                                  A-54

-------
    V(LSFS)=VJ                                                         U  101
    P(LSFS)=PPK                                                        U  102
 95 LSFS=LSFS+1                                                        U  103
100 LASTI=0                                                            U  104
    I-K                                                                U  105
105 I»JRL(I)                                                           U  106
    IF (I.EQ.O) 60 TO 130                                              U  107
    PPK=1                                                              U  108
    JMIN=IRU(I)                                                        U  109
    JMAX-ISU(I)+IU(I+1)-IU(I)-1                                        U  110
    IF (LASTI.GT.I) GO TO 110                                          U  111
    LASTI=I                                                            U  112
    LASTID=JMAX-JMIN                                                   U  113
    IF (JU(JMIN).NE.K) LASTID=LASTI[>1                                 U  114
110 IF (JMIN.GT.JMAX) GO TO 105                                        U  115
    DO 125 J=JMIN,JMAX                                                 U  116
    VJ=JU(J)                                                           U  117
115 PK=PPK                                                             U  118
    PPK=P(PK)                                                          U  119
    IF (V(PPK)-VJ) 115,125,120                                         U  120
120 P(PK)=LSFS                                                         U  121
    V(LSFS)=VJ                                                         U  122
    P(LSFS)=PPK                                                        U  123
    PPK=LSFS                                                           U  124
125 LSFS=LSFS+1                                                        U  125
    GO TO 105                                                          U  126
130 PI=P(1)                                                            U  127
    IF (V(PI).NE.K) GO TO 210                                          U  128
    IF (LASTI.EQ.O) GO TO 135                                          U  129
    IF (LASTID.NE.LSFS-3) GO TO 135                                    U  130
    IRUL=IRU(LASTI)                                                    U  131
    ISU(K)=IRUL+1                                                      U  132
    IF (JU(IRUL).NE.K) ISU(K)=ISU(K)-1                                 U  133
    IU(K+1)=IU(K)+LASTID                                               U  134
    IRU(K)=ISU(K)                                                      U  135
    GO TO 150                                                          U  136
135 ISU(K)=JUPTR+1                                                     U  137
    PI=P(1)                                                            U  138
    PI=P(PI)                                                           U  139
    VI=V(PI)                                                           U  140
140 IF (VI.GT.N) GO TO 145                                             U  141
    JUPTR=JUPTR+1                                                      U  142
    IF (JUPTR.GT.MAXPU) GO TO 215                                      U  143
    JU(JUPTR)=VI                                                       U  144
    PI=P(PI)                                                           U  145
    VI=V(PI)                                                           U  146
    GO TO 140                                                          U  147
145 IRU(K)=ISU(K)                                                      U  148
    IU(K+1)»IU(K)+JUPTR-ISU(KH                                       U  149
150 I=K                                                                U  150
                                   A-55

-------
155 Il-JRL(I)                                                          U  151
    CEND-ISL(I)+IL(I+1)-IL(I)                                          U  152
    IF (IRL(I).GE.CEND) 60 TO 160                                      U  153
    IRLI«IRL(I)                                                        U  154
    J-JL(IRLI)                                                         U  155
    JRL(I)=JRL(J)                                                      U  156
    ORL(J)-I                                                           U  157
160 1=11                                                               U  158
    IF (I.EQ.O) GOTO 165                                              U  159
    IRL(I)-IRL(IH                                                    U  160
    GO TO 155                                                          U  161
165 I *                                                                U  162
170 Il-JRU(I)                                                          U  163
    REND*ISU(I)+IU(I+1)-IU(I)                                          U  164
    IF (IRU(I).GE.REND) GO TO 175                                      U  165
    IRUI-IRU(I)                                                        U  166
    .WUflRUI)                                                         U  167
    JRU(I)=JRU(J)                                                      U  168
    JRU(J)=I                                                           U  169
175 1=11                                                               U  170
    IF (I.EQ.O) GO TO 180                                              U  171
    IRU(I)«IRU(IH                                                    U  172
    GO TO 170                                                          U  173
180 I=IRAC(K)                                                          U  174
    IF (I.EQ.O) GO TO 195                                              U  175
185 I1=JRA(I)                                                          U  176
    IRA(I)=IRA(IH                                                    U  177
    IF (IRA(I).GE.IA(I+1)) GO TO 190                                   U  178
    IRAI«IRA(I)                                                        U  179
    IF (JA(IRAI).GT.I) GO TO 190                                       U  180
    JAIRAI=OA(IRAI)                                                    U  181
    JRA(I)=IRAC(JAIRAI)                                                U  182
    IRAC(JAIRAI)=I                                                     U  183
190 1=11                                                               U  184
    IF (I.NE.O) GO TO 185                                              U  185
195 CONTINUE                                                           U  186
    ISL(N)=OLPTR                                                       U  187
    ISU(N)=JUPTR                                                       U  188
    RETURN                                                             U  189
200 CALL YSMER  (3HROW,K,13H OF A IS NULL)                              U  190
    GO TO 235                                                          U  191
205 CALL YSMER  (3HROW,K,20H HAS DUPLICATE ENTRY)                       U  192
    GO TO 235                                                          U  193
210 CALL YSMER  (3HROW,K,17H HAS A NULL PIVOT)                          U  194
    GO TO 235                                                          U  195
215 CALL YSMER  (3HROW,K,19H  EXCEEDS JU STORAGE)                        U  196
    GO TO 235                                                          U  197
220 CALL YSMER  (3HCOL,K,20H HAS DUPLICATE ENTRY)                       U  198
    GO TO 235                                                          U   199
225 CALL YSMER  (3HCOL,K,17H HAS A NULL PIVOT)                          U  200
                                   A-56

-------
    GO TO 235                                                          U  201
230 CALL YSMER (3HCOL,K,19H EXCEEDS JL STORAGE)                         U  202
235 IER=1                                                              U  203
    RETURN                                                             U  204
    END                     -                                           U  205-
                                  A-57

-------
   SUBROUTINE NSCORD (N,IA,JA,IAWORK,JAWORK,C,IR,ICT)                  V    1
   INTEGER IA(1),JA(1),IAWORK(1),JAWDRK(1),C(1),IR(1),ICT(1)           V    2
   DO 5 I«1,N                                                         V    3
 5 ICT(I)=0                                                           V    4
   IAWDRK(1)=1                                                        V    5
   DO 15 K-l.N                                                        V    6
   ICK-C(K)                                                           V    7
   JMIN=IAWORK(K)                                                      V    8
   JMAX=JMIN+IA(ICK+1)-IA(ICK)-1                                       V    9
   IAWORK(K+1)=JMAX+1                                                 V   10
   IF (JMIN.GT.JMAX) GO TO 15                                         V   11
   IAINK=IA(ICK)-1                                                    V   12
   DO 10 J»JMIN,JMAX                                                  V   13
   IAINK-IAIMM                                                      V   14
   JAOUTJ=OA(IAINK)                                                   V   15
   JAOUTJ=IR(JAOUTJ)                                                  V   16
   JAWDRK(J)=JAOUTJ                                                   V   17
10 ICT(JAOUTJHCT(JAOUTJH                                          V   18
15 CONTINUE                '                                           V   19
   IA(1H                                                            V   20
   DO 20 I=1,N                                                        V   21
   IA(M)=IA(IHCT(I)                                               V   22
20 ICT(I)=IA(I)                                                       V   23
   DO 30 I=1,N                                                        V   24
   JMIN=IAWORK(I)                                                      V   25
   JMAX=IAWORK(I+1)-1                                                 V   26
   !F (JMIN.GT.JMAX) GO TO 30                                         V   27
   DO 25 J=JMIN,JMAX                                                  V   28
   JAOUTJ=JAWORK(J)                                                   V   29
   ICTJ*ICT(JAOUTJ)                                                   V   30
   JA(ICTJ)«I                                                         V   31
25 ICT(JAOUTJ)=ICTJ+1                                                 V   32
30 CONTINUE                                                           V   33
   RETURN                                                             V   34
   END                                                                V   35-
                                  A-58

-------
 BLOCK DATA ALPHA!                                                  W    1
 COMMON /ALPHA/ IGO(4),IBLANK,MBLANK,JINTER                         W    2
 COMMON /APLOT/ JVERT(52,2),JBLANK,JSTAR,JPLUS,JBAR                 W    3
 COMMON /GEAR1/ T.GUESS,HMIN,HMAX,EPS1,UROUND,NC,MF1,KFLAG1,JSTART  W    4
 COMMON /HEAT/ CV,Q(200)>SC(200,7),ISC(200,3),ITEMP                 W    5
 COMMON /STORE/ AST(35),IPL(7),TEMEND,NTEMP,TMI,NPHOT,PHI,IL,NFRST, W    6
1IPH(30),QM(100),PM(100),PSTOP                                      W    7
 COMMON/PHOTR/ IPP,RDAT(80),RTIM(80),RR1,MO                      W    8
 DATA IGO(1)/4HMORE/,IGO(2)/4HCONT/,IGO(3)/4HPLOT/,IGO(4)/4H    /,! W    9
1BLANK/4H    /,MBLANK/4HM   /.JINTER/4HINTV/                        W   10
 DATA JVERT/12*4H    ,4H C  ,4H 0  ,4H N  ,4H C  ,4H E  ,4H N  ,4H  W   11
IT  ,4H R  ,4H A  ,4H T  ,4H I  ,4H 0  ,4H N  ,4H    ,4H P  ,4H P   W   12
2,4H M  ,75*4H    /                                                 W   13
 DATA JBLANK/4H    /,JSTAR/1H*/,JPLUS/1H+/,JBAR/1HI/                W   14
 DATA ITEMP/4HTEMP/,TEMEND/0.0/,PSTOP/0.0/,IN10/1/                  W   15
 DATA UROUND/7.5E-15/                                               W   16
 END                                                                W   17-
                                A-59

-------
       APPENDIX B
INPUT DATA--CHEMK FORMAT
           B-l

-------
             APPENDIX B
      INPUT DATA--CHEMK FORMAT
69 SAVE
1 R03
2 0
3 RO 03
4 R02 03
J R02 0
6 OH 03
7 1102 03
o on no2
9 Oil CO
10 no HO
1 1 no H03
12 H02 (103 ICO
is ro 1102
14 n02 1102
IS 0 PAR
16 Oil PAH
17 0 OLE
'fl 0 OLE
. ? Oil OLE
20 03 OLE
S 1 03 OLE
S2 0 ETII
23 0 nil
24 Oil ETII
23 03 ETII
26 RO AC03
27 NO nnn2
2(1 no RA02
24 RO HE02
TO RO ME02
31 HO ME02
32 03 HD02
S3 03 RA02
34 03 HE02
33 ou cum
36 OU CARD
37 CARD
30 aa
39 aa
-0 CARD
* PAH x
42 1(02 AC03
•S3 PAfl
•;•> 1102 AC03
->3 1102 ME02
46 MO CRIC
47 .102 CHIC
4.1 CARO CRIC
•»? no nciic
a> .102 ncnc
a i CAIUI note
13 >:iiu:
nil I:MIU
H4 Clllli
15 item;
oo IICIMI
97 III:IK;
ia tn:uc
vt uu ARO
60 UU AIIO
6 1 OU A1\O
. 000324
RO 0
03
H02
R03
RO
1102
on
I1N03
H02 CO2
N02 H02
R02 H02

R02 OU

KE02 Oil
HE02
HEU2 AC03 X
CAIUJ PAR
RA02
CARO CRIC
CARD ncnc x
HC02 1102 CO
CAIIO I'AR
RD02
CARO CRIC
R02 HE02 C02
R02 CARD CARO 1102
N03 CARD CARO 1103
H02 CARO neoa x
R02 CARD 1102
RRAT
CARD CARB H02
CARB CARB H02
CARB H02
11UU CO
AC03 X
00
1102 1102 CO
NE02 X U02 CO
CO

PAN
AC03 H02


R02 CARD
.103 CARD

!»02 CAtUJ P,\I\
.•O3
2. 700 E* 03
2.?oor»on
4.200L+04
O.OOOE-03
8.000E-03
6.000002
6.0IIOO02
1.200I1+04
2.400E-03
3.UOOE+03
1.200E+04
1 . 200F>04
4.000E+03
8 ••••£* 93
S.000E+02
B.«OOE+««
2.000E+92
S.OOOE^OO
9.SOOOU3
9.300E+03
I.OOOE*00
I.OOOH+02
I.OOOE»02
3.SOOE-01
I.OOOE*03
2.000C+03
2.QOOE-02
4.000E+03
4.000C+03
l.200t+04
O. 000l>03
2.000C*03
1 . 20Or.« 04
n.oo(in»oo
2.ouoi>o:i
ft .Tnoi'to'j
U. 4IIIII lli'l
V. (MIDI. Hl|
1 i &OOI1*U«
II. linn uu
4. U&1IO03
a.»o«K»oi
«.000t>03
I.60OE«03
1,SOOE»04



1400.*
2430.0

1000.0
1S23.0




-10600.0






























12300.0


















62 CARD W
                B-2

-------
 *•  HO   Y         RO   CAM PAH        3.«0«E+01
 04  HO   Y         103  AERO            1.S00E+OI
 ••  HO3  Y         CAM CAM            3.S00E+04
 M  03   Y         AERO                 6.009E-01
 47  OH   GLY       B02  OC   CO         1.000E+04
 •a  cc             Y    Y    Y          1.099E+92
 M  CLY            KE02 B02  CC         2.928E-01
•DM
 TO       SAVE      .••6324
  I  1102            HO   0              .3
 •7  CAM           OA                  .00064
 40  CAHB           CO                  .90932
 M  HO   HE02      R02  CAHB HE02 X    2779.
 30  HO   JJE02      H02  CAHB 802       923ft.
 ««  CLY            KE02 B02  CC        .»72
 Tl  RX             OB                  .»3
 72                 co                  .003
 73  03                                 .001
 TO  H20                                -.••6324
BC-237
 10                           69.       60.
•O        H02       OLE       PAR       AHO       CARS      B20
RX        CO        ETH                             .. .
 .377      .1*6      .09       7.38      .177      .1*12     24000.
 .994      .9        .87S
          36«.      3*3.2
PLOT
EC-237
      3         •••• «.2» «.49 «.6« •.86  9.       .8        •.        4«9.
03   23 0
•.         .992       IB.        .992      39.       .91       48.        .»42
M.        .11       78.        .IBS      90.       .264       198.       .338
 12».       .391       138.       .482       180.      .811       168.       .887
 1M.       .894       198.       .62       2i«.      .64       228.       .68
 24*.       .683      288.       .63       270.      .648      288.       .638
 3M.       .628      318.       .613      330.      .693      348.       .894
 36*.       .884
 HO  23 »
 • .         .377       IS.        .272      30.       .169      48.        .*83
 M.        . »4      78.        .024      90.       .«1B       108.       .913
 12*.       .911       138.       .01        ISO.      .097       168.       .007
 IB*.       .006       198.       .006      210.      .006       228.       .006
 24«.       .••«      288.       .006       270.      .008       288.       .008
 30*.       .993      318.       .008       330.      .098       348.       .90S
 360.       .008
 •02  28 2
 0.        .106       18.        .198      39.        .29       48.        .383
 M.       .368      78.        .38       99.        .324      108.       .297
 120.       .268      138.       .234      180.      .208       168.       .171
 180.       .142      198.       .117      210.      .096      22S.       .07
 240.       .074      288.       .068      270.      .063      285.       .08
 300.       .039      318.       .038      330.      .037      348.       .037
 369.       .036
                      BLANK CARD
                      BLAHK CARD
                                  B-3

-------
     APPENDIX C
INPUT DATA-ERA FORMAT
          c-i

-------
                           APPENDIX  C
                 INPUT  DATA-ERA  FORMAT
 n
KB
o
HO
no3
«os
01
003
01
08
NO
NO
NO2
NO
B03
0
OH
O
O
oa
03
03
0
0
OH
03
NO
NO
NO
NO
rt
NO
NO
03
00
03
OU
OH
CAM
CAM
FUR
.102
PAN
BOS
U02
NO
N02
CAM
no
NO2
CAM
CHIC
CRIC
CHIC
ncnc
ncnc
HUM:
 note
OH
OH
on
 CAM
no
NO
Noa
03
OH
CC
CLY
   SAVE
O3         3
03         4
O          3
03         4
03         7
N02        8
CO         *
NO        10
Noa       ii
N03   020 12
D02       13
U02       14
P.XR       IS
PAR       16
OLE       17
OLE       IB
OLE       l«
OLE       20
OLE       21
ETO       22
ETtl       23
ETO       24
CHI       23
ACO3      26
RI102      27
RA02      28
HE02      3»
          30
HC02      31
HE02      33
ROO2      33
RAOS      34
NEO2      33
CAM      36
CAM      37
          38
          34
          40
          41
X         43
AC03      43
 .000334
    HO
    03
    1(03
    H03
    HO
    BO2
    OB
    BH03
    003
3.  N02
3.  H02

    R02

    HE02
    HE03
    HE03
    CAM
    RA02
    CAM
    CAM
    HE03
    CAM
    MOS
    CAM
    H03
    KM
    N02
    R02
    HEOS
    N02
    NRAT
2.  CAM
3.  CAM
    CAM
    003
    AC03
    OA
3.  803
AC03
HE02
CRIC
CRIC
CRIC
ncnc
ncnc
new
                 CO

                 PAH
                 AC03
ARO
ARO
ARO
V  .
T
T
T
T
CLT
          43
          46
          47
          48
          49
          SB
          31
          32
          S3
          34
          as
          36
          S?
          SB
          S«
          60
          41
          62
          63
          44
          49
          66
          47
          68
          6«  S.  T
          rt
    NO2
    HOS

    N02
    N03

    CO

 3.  003

    HC02
    HC03
    CAM
    •02
    BOS
    OH

    M
    M3
 3.  CAM
    AERO
                            COS



                            OB

                            00

                            AC03
                            PAR

                            CRIC
                            ncnc
                            002
                            PAR

                            CRIC
                            HE03
                         3.  CAM
                         3.  CAM
                            CAM
                            X
                            CAM

                            B02
                            802
                            HOS
                            CO
                            X
                            CO
                             H02
    CAM
    CARD

    CARD
    CAM
    C02

    00
    BOS
 2. 002
a.  v
    CLY
    CLY

    CAM
    AIM
                CC
CO
COS
CO

X
v

PAH
                            00

                            06
                                               1.000E*00
                                               4.400E*06
                                               S.300OOI
                                               4.M6C-62
                                               1.340e«04
                                               r.r0oe*0i
                                 1480.
                                 «4«0.
                                                         isaa.




X



X
CO



COS
003
DOS
pp

DOS







B03








PAR
PAR

1.S60E-03 -10*00.
1.200E+04
l.SOOE+04
2.000E+41
2I700E+03
2.700E»03
4.2OOE+04
B.OOOE-03
B.OOOE-03
6.000E*02
6.000E*02
1.200E»04
2.400E-03
3.800E«03
1.200E+O4
1.200E+04
4.000E+03
1 . OOOE+03
a.000£*93
S.OOOE+02
3.000E«>«0
3. OOOE+02
S.OOOE«00
9.300E+03
l!oO0E»00
1. OOOE+02
l.OOOE+62
3.300E-01
2l OOOE+03
a.aooE-oa tssoe,
4.000E+03
4.000E+03
1.200E+04
a.OOOE+«3
2.000E+03
1.200E+O4
a.OOOE+03
2.00UE+03
                                               2.400E»02
3.400E»OS
4.230E»0S
a.SOOE+6!
6.AOOE+03
1.600E+63
I.S00E*«4
l.000E»0a
3.M6E+41
                      1.000E«04

                      3!023C-«t
                               C-2

-------
•DM
 74   t   SAVE       .88*324
Mtt              t      HO          O                  .8
CARB            38      OB          CO                ,MM4
CARB            41      CO                            .66632
HO    HEM      29      1102         CARB        PP    2778.
RO    KE03      31      B02         CARB        H02   9230.
CLY             76      HE02        102         CC    .972
RX              71      OH                            .63
                72      CO                             .663
03              73                                     .Ml
B20             74                                    -.•••324
EC-237
  !•                           60.       66.
HO        HOa       OLE       PAR       ARO       GARB      H30
RX        CO        tTB
 .377      .1*6      .99       7.38      .177      .1*13     24666.
 .••4      .9        .878
          36*.      363.2
 run
 CC-337
       3         •••• «.2« 6.4* •.«• «.M  ».       .S        t.        4M.
 O3    29  O
 • .         .M2      18.        .t«2      3«.       .•!       48.       .»42
 «0.        .11       78.        .108      90.       ,264      1*8.      .938
 126.       .391      138.       .463       IB*.      .811      168.      .887
 IBO.       .894      198.       .62       210.      .64       228.      .68
 24*.       .689      238.       .68       37».      .648      263.      .638
 3««.       .620      318.       .413      939.      .6t3      348.      .894
 3*6.       .884
 HO   23  H
 • .         .37?       18.        .272      36.        .169      48.       .683
 6«.        .64      78.        .424      90.        .fl!8       105.      •••13
 126.       .611       138.       .81        180.       .697       168.      .8*7
 180.       .0*6       198.       .666       210.       .006      229.       .8*6
 246.       .866      388.       .806      276.       .60S      283.       .8*8
 3*6.       .668      318.       .888       336.       .668      348.       .888
 366.       .888
 H02  29 2
 6.        .166      18.        .198      36.        .29       48.        .383
 66.       .368      78.        .98       98.        .384      108.       .297
  136.       .268      138.       .234      186.       .268      168.       .171
  186.       .142      198.       .117      21*.       .696      228.       .67
 246.       .674      288.       .668      270.       .663      288.       .68
 3*6.       .639      318.       .638      336.       .637      348.       .637
 368.       .836
                      BLAR CAM)
                            CARD
                                   C-3

-------
                 APPENDIX  D
                 OUTPUT DATA
(using EPA input format for the chemical  reactions)
                        D-l

-------
                           APPENDIX  D

                          OUTPUT DATA
      33   S  S
5     u u   u u          u
S     O 9   9 a          9
U     a is   9 tt          o
      » •»   9 a          9
§
t»
U


2
• • — M«->-999O-9OO9O9999999999O999
II     I      I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I  I
                                                        U  U U U
                                           x u
8 S ? =
o u S o < £ u
n tt el n a n a
rinfiontfifi ft ooococ.
o
I]
S « s
s § !
i i i f

- 8
* 3
y r
Hi

.§
2
n
n
e
9B

W U X U
n
« Cl S 0
2 1 S I

                                                       n ct  n
        fi  n    ti
                    2 = 2 2 * 2 « t  2  *  ?,  «  S 8 S S S S S S  S
                                0-2

-------
 e
 i
       o
                                          o  o
                                          i   i
                                 o  e
                                 i   i
            o  o>  e  o
            i   i   i   i
                           e  o
                           i  i
    •  tt  e  n  n  e
    e  e  e  e  e  e
    *  -*  +  *  *  *
    u  u  u  u  u
       e  o  e  o
       e  e  e  o
       -  e  e  «a
                      n
                        s  ;
                        n  n
                        e  e
o  n  »
e  e  o
n  n  ri
e  e  e
ct  -
e  e
N  n
e  o
                                                                              o  •
                                                                              +  *
                                                                              u  u
n  » B

22!
 sis
 e  a  ci
   u  u  u
   e  e  o
   o  e  e
   o  e  o

P.  —  —  ~
                                                      i  ?
                                                      a.  a.
                                                                              ri
                                                                              o  o
                                                                              u  u
   ri  n  ci
   222°
               x  S
                                             §  g
                                          .

                                      u  o
                    CI
                    o
                    u
                    fl CI    >- ».
                 = 00    J J
                 o = e  >• o 3
H  g  q
<  S  5
   u  o  u

   N  N
               n        ri
               o     ri  o
               U0OU
               
-------
?'••••••
 I I  I I I I
4 2 a 2 2 2 2

       *
     3  3
     3  I
siiii.
>• »i
• * » »
2 S  8 S 8
  : s s
              D-4

-------
 2

 >•
    •  «SSSe«SSS

    u'uuuuuuuwu
 ^
           a
           i
              S  §
       o     <  <
    e  u     u  u
              fl   b J
       Sa  o  e   o
       o  u  E   B
a  o
o  u
                 N


                 §
*  s
B  E  W  U

O
              2  g   d  S      S3
   —  B  —  <*  ~   »  ~   tt  m  
-------
   ;     2888 Xv 88

 ES     i lit |H
   I.     ••••«!••
   a     .; a - - ~ « -

   «•     c« « a ct a ~ M
   *     »*?*???

 M     U Hi 11



   8     22^:22:


 9!     !!!•!!!!

   I     22:2:2;


 M     Is IIHI


   ;M    : 2 2 2 2 2 :

 is!    iiiiiii
   • «    -< a a o • o •
   J *    ei -•<>•«•- a

   — "    « _ n • n •> -•
   aki IB    w u ill hi ia M u
   S I    ? S S 2 S S S
 •  £ &    ••«•••»



   * 5    i » » ? i * * *
 ^MU    feiUMUliiblMlil
        K8:898S9
        *Btl»H«5N



   2 -    ? ? • S ? ? ! •

 Mi:  1 I I I i I! i
   •>•»  ei — *•«!• — —
    • «l




   -;Sa»*---*-"
s  .EgS,^	

   S4 &J «4 7 *P^*Jfll^^JO
   !*• » Z!r C ^
-------
52,
                            fj O O O O O O
        O
        I
  u     III
  ->    O
C— -J    O
   -
                rj
                3
                O
                O
                a
                            OOCO3OO
          oo

          ii
                  -f n
                  so

                  ii

                  §§
                  1 1 N
                9 - «>
                  I I
                U

                I
             *
             O
                   nci
                   00
                   UU
                   NO
        oor    o
-------
        tiff
        ess.
                                *   — c -
                                              » C r-
                                              c to
                                                ' '.
                                        .          .
                                              ~
                                              »««»

                                              — «ri
            e-c:--z

            W M H* iL w. —
            »r. :is tut ;
            erac- — e -r ti
                                                                                c r: .-3   ci c j
                                                                                c r. s   c c f
                                                                                -»*>   I. e>»
                                                                                                    e c £ i. - i c
                                                                      ci :i I- £ II i
  fan   B a o
  ee   eee
« i  i     iii
UUid   VklU

SiJ   §«s
c * -   ti e n
                                                          I I  I  I I  I I
                                                         uuuuuuu
                                                         « e> eo

                                                      *—»
                                             eee
                                              ii
                                              -
                                     » c =-
                                     ens
                                     » v »   e B e
       te-o-»
                                             at.cn
                                             eeee
                                             r»ene
                                             ecs-r?
                                     eeee  eee
                                     o —«M
                                     oecn
                                             Oefr>
                                                  eeee  eee
                                                  i i  i *   t t  i
                                             Hg.  uuuu  uuu
                                           u < £  «N;iri  i»no
                                           — BU  — e«s>B  aee
                                           O«(-  o«-r.o  -v^
                           uuuuuuuu
                                                                                aaca  i»-r
                                                                                eeee
                                                                                iiit
                                                                                uuuu
                                                                                B —«t>
                                                          oeee
                                                          tiii
                                                          uuuu
                                                                                esoo
                                                                                n —«•—  — o — <
                                                                       slis
                                                          eooe
                                                          i i  t i
                                                          uuuu
                                                          FiSHC'
                                                          nee«»^
                                                       < — aesn — o 11 o-
                                                                       SB
                                                  ne
                                                  ee
                                                  4- +
                                                  uu
                                                  ec>
                                                  en
                                                  tin
                                                                      »«»*»««»
                                                                      eeeeeee
                                                                      i i  i i t  i i
                                                                      uuuuuuu
                                                                      « ti f n « ss »
                                                                                                     t— «» — —
                                                      oeeeeee
                                                       i i  i  i i  i i
                                                      M U U U M U (io>«— o —
                                                                                                    i  i i  i  i i  i
                                                                                                     uu uuu L;
                                                                                                     n n » — n c
                                                                                                     — o—
                                                                                                   e» — « N r « n u
                                                                                                      \»
                                                                                                   ecoeeeeo
                                                                                                    i i  i  i i  i i  i
                                                                                                   uuuuuuuu
                                                                                                   •ariGtic-etti-
                                                                                                   o —r»oe-o —a
                                                                                                   — *
                                                                                                   eeeeeeee
                                                                                                   i i  i i i  i i  i
                                                                                                 uuuuuuuuu
                                                                                                          e — rso
                                                                                                          
-------

           S3»
           OS
         f 1  1
         aau
                    1
                     eoaaoeo
                      i i i  i  i i  t
                     U U U U bl Cd U
                     r3rse>o o oci
                     f j o — o «i o a
                                               8e»
                                               oo
                                                       i  i  <
                                                      «uu
                                                      • 00
                                                      eot-
                                                      SNN
                                                               fl —
                                                                 I
                                                                           f I  I  I  f 1 I
                                                                           uuuuuuu
                                                                           .«e>oo — oo
                                                                           o o a o o e- o

                                                                           - * ei ci - cj —
  S J
  uu
 o r- ci    n N a
 ooo    ooo
 + 11     iii
 uuu    uuu
ecio

— ei —
                     OOOOOOO
                     I I I  I  I
                      ! — — ONCIO

                     N a » ei n - ei
                                u    o
                                _    U
                                       e-i
                                       UU
                                                      or-ci
                                                      000
                                                      nan
                                                      oeo
                                                               UU U
                                                               OOO
                                                                           n
                                                                           •» ^ s- o cira t
                                                                           i3O«or»r-o

                                                                           ••a ci ci si * c> —
n§
g=§=
   sS
         ooo
         uui

         823
         ooo
         iii
         a*W
         f. — a
         one
         -Q-
         NOTI
         ooo
         i  i  i
         uuu

n a o *
oooo
I I  I -f
uuuu

?§C10


r« — n ci
         — ci n ci
         0930
         n f- t- ci
         n t. -.• n
         rjr.no
        o ooo
         1  1 1  1
        -UUU
        oa-aa
        -an GO
        nci — o
        oano
        o ooo
        •TO CIO
        oeeio
        anoo
         UblU
         O>— f.
         030CI
                 -00*
                 ooo
         N»C1
         C1S-
         ON*
         eeo
         — ri o    — — v
                 enan
                 ooae

                 uuuu
— f -ci
e ooo

•' •* •'• u
n 3 NO

S-T3-

ci ci n ci
          Q
         0-00
                 ncir.0
                   i
         floor*
         ooeo
         i  i i i
         uuuu

         ?,??*-
        OOOQ
         till
        hlUUbl
        ^ rt C3 rt
        Nar-a-
uuuu
cin cie
NO nci
                 cm — a>
                 i  i  r
                 I
                              I  I  I I  I  f  I
                              UUUUUUU
                              n o ci o c = n
                              o u t- ci o e- a

                              o-a —nod —
                     •c r» * f o o n
                     o oooooo
                     I I  1  I I I I
                                     II
                            sis
000
I I  I
uuu
ooo
•f-CI
ooo

iii
NNCI
oo —
O 13 f
                                                      — O-&
                                             ooo
                                              I  I I
                                             uuu
                                             OO 3
                             ci N o r* N r- o
                             C O — •!• CJ — «•
                                                      NO*
                                                      OOO
                                                      I I  I
                                                      uuu
                                                      nno
                                                               ont-
                                                               uuu
                                                               »« o
                                                               o n tf>
                                                               — — cs
                                                               i  i
                   •  - — T  -. s •-. is a

                   ^o«j-'-cio-«>  2s
                                                          o N o « « o e>
                                                          ooooooo
                                                          i  i  i i i  i  i
                                                          uuuuuuu
                                                               ci f- a ci
                                                                           o - o o c e> o
                                                                            o 3 c a ci o
                                                                           r. a n to r> ci ti
                                                                           n«3cccc
                                                                           .
                                                                           uuuuuuuu
                                                                           o r. e> •* T ct ci a
                                                                           ci n o CT n •*« »
                                                                           sooeoooo
                                                                           i  i i i  i  t  i i
                                                                           u u u u u u u u
                                                                           Q —
                                                                          ci •» o n « •» n n
                                                                  oni3« cici «• o

                                                                '/> ei --•»-
                                                                         1
                                                                2 ei c t a r r- —
                                        D-9

-------
        fll»C
        ee —
        i i i
        ». —r.   *|»'«

        an*   —«««•
                          •eeeeeo
                          I I  I i  I i  i
                          UkJUUUUU
                          o » e c »
                                             »ee
                                             e ——
                                              i i i
                                                               ri • n « » n n
                                  B    e -
                                  -    us;
                                              •HE   enc
                                                               c :i» ? si is f
                                                               t» * — —!• — f
        8'eS   SiS
          II    II
rtctn


SSS    KKt-
                         •>eo
                         nS-c*-o-
                                              or-n   ne>n
                                              eae   eee
                                              + i i    itt
                                              UUU   UUU
                                              can   o — e-
                                              isi.0   nea
                                              eta»   atio
                                                                i i  i i i i i
                                                               uuuuuuu
        2§S   8SS
                                              ss:
                         •>Q
                         ee
        «B«S-   -O-
                            eeeee

                         uuuuwuu
                          • — * —
                                      	   UUU
                                      non   oeo
                                        >—   cio a
                                             BBC
                 SS
         poo   one
         :eo   -•»«•
        e»—   aeo
        oeo   eee
        tit   iii
        uuu   uuu
        no—   »s:i
         —f   I-
  »£§!
  	  nan
eeee  eea
i i i  *  iii
uuuu  uuu
e> ci t» ri  » n ci

n-ee  Sea

»ri —«j  -ti-e
       i i
        eeee  eeee
        iiii  i  i i  i
        uuuu  uuuu
        n — tt e  e c ci n
        VNOO  onNO
        BO — e  oocio
                  1
siis
            e
        i •* i i
        uuuu
         §t7»e>
         -cio
        O — !• —
               eeee
               i  i t i
               uuuu
                — OO
no
   pe
       man—
CO-J1C6
IO-3CCI
l-0t300»
CO-3V9C
               oooo

               uuuu
                goae
                t
       n —
                                              eee

                                              wuu
                         eeeeeee
                          i i  t t  i i i
                         uuuuuuu
                            ——-CIO
                                                               ooooooo
                                                               I  I I  I 1  I t
                                                               uuuuuuu
                         eeeeeee
                          i i  i i  i i i
                         uuuuuuu
                                                               i  i i  i i  i i
                                                               uuuuuuu
                                                               - »i. n o n o
                                U-«O
                         eeeeeeee
                         i i i i i i i *
                         uuuuuuuu
                                                               eeeeeeee
                                                               i  i i  i i  i i »
                                                               uuuuwuu
S
                          oeoeooo

                         HUUUUUUU
                                                                uuuuuuu
                                                                        c
                                                               e-oneoeov
                                                               -n — ««— co
                                                               eeeooeeo

                                                               UUUUUUUU
                                                                   — laoco
                                                                   ——tsee
                                                               oeeeeeee
                                                             — i i  i i  i i  i i
                         uuuiuuuuu
                            — fjo»n—  cs
                                     D-10

-------



















N
Cl
1
u











-OBCI
tfi





LJ

2
••
o ci



B
O
M
Cl

i
5
if

s
H



o
n








0
o






o







I
>
             O.O.C
D-ll

-------
          APPENDIX E

  EXAMPLE USE OF VARIABLE
PHOTOLYSIS AND PRINT OPTIONS
                E-l

-------
                          APPENDIX  E

             EXAMPLE  USE  OF VARIABLE
          PHOTOLYSIS AND PRINT  OPTIONS
 Ml                        14   *O   T             1
        I        |2        IT        18       I*       24       29
.•16      .13      .29      .39       .46       .8       .82
.93      .83      .9       .46       .3*       .29       .13
 01*
mn             i      no         o                i.
O               2      03,                          4400000.
.3   HO         3      H02                         23.8       I4SO.
hO2  03         4      H03                         .0474      2430.
H«J2  0          8      HO                          13000.
M3  HO         63.    H02                         20000.
»02  M3        7      H20S                        3910.      -861.
CSOS             8     H02        H03              6.03       10320.
k^OS             92.    BK03                        .02
hi)   R02       102.    BOHO                        .0000003
HOMO  UOHO      11      HO         BOS              .000018
noao            12     on         MO               . 10
OB   B02       13     IBI03                        16000.
OH   HO        14     BONO                        160OO.
002  HO        13     R02        OU               12000.
DOS  1102       16     BOOB                        7400.
UOOB            172.    OU*                         .00079
03             IB     010                         .002?
03             I*     0                           .069
OlD            2O     0                             4.4SE10 -97.4
010            212.    On                             6.8E09
OB   03        22     BOS                         76.9       1000.
03   802       23     OB                          2.03       SBO.
BCUO           24                                 .0046
UCaO           2S2.    U02                         .003
1KUO OB       26     U02                         14000.
HE02 HO       2O     HEO        H02              12000.
MEO            29     UCIIO       002              200000.
HCO   H02       3fl                                 18000.
HEO   «02       31     BCUO       BONO             4400.
HC02  B02       32                                 74OO.
REACTIVITY 1000 fPMC HETBAHE
  3                         60.       60.
HO        M2       C04
 .1        .1        1000.
          840.      303.      .01
fLOT
IUOCTIVITY 1*00 erne KETILXHE
      I       0.00 0.40  0.80 1.20 1.40 0.        1.6      0.        880.
03      0
           BLANK CAW)
           BLAUK CARD
                                E-2

-------
^,
i
u
s


n
0
u
e
n
»»»«l.  OS'S   —   fl
                                                            n«f.O~--fio*'^»f"!ST-*—   f»
                                                            —   —  —       •  ri  ci   fi   ei   N   ct   N   ri   M  :i  n  n   n
                                                                      E-3

-------

B
£
•.
S
0
f
»
I
S




•
^
•



i
» u
§ i
"
o

§u
1

-
e
»
n
f
•
M
Eg
5
1
o
n
e
i
S
C9
•
9
y
U rf

U ul
c S
!
0
I
ei
••
E
a
S
«•
n
i
•
•
;
i
ei
i
t>
8
u
IS
•

e
ui
e>
t»
0
ei
•
e
A
&
o
»

»

•
•V
«
u
«
N
V
«
I
»
«
8
u
e
ci
:
g
•
K
n
E
e
e
»
»•
f
S
*
u
o
«

»
e
+
u
e
e
<9
"
o



••

S
e
e
0
|
e*
I
•»
i
»
e
i
e
»
rt
I
«
CJ
e
ill
e
«
»

S
•^
u
s3
fl
O
«
5 S
4> 4>
2 u
S 0
ri e
— »

a N

» s
»
N
»
i
2 S
|
i
E
?
^ ^
u
g
Ik
e
S3
3 S
< W U
Sfc.
e e
= S
8
i
e
ca
8
J>
e
o

«
n
e
u
e
o
N
ei

*
e
^
1
r»
_
o
%

••

P
i

H
O
*
1
X
_l
<
>
1

s
2

3
S
E
u
e
e
e
o
i
s
0
i
0
o
1
hi
e
e
M
0
O
i
o
I
»
e
u
i
V
C)

o
u
0>
ti

e
i
u
e
S
rt
|
0
i
u
^J

e
u
§
ri

e

u
e
e
fe g
       n
       o
§n
e
ul
5
e






3
3
^*

S
e>
£
•j
S
a
3
V

e
2
u
5
u u u
o s e
* » «






u
T
•T





2
»
w)
i*»
?~

£
55
-J
P
d
_ *
•ri
• '
3
H
^
r™
u
p*
o
»
u
o
9
O
v




ER
5
:—

i^
rt
i
y
P
Cl
a
i .
u
e
0
^
~


o
i

a
«
T


                                      •o
                                      0)
               E-4

-------
             •      <8      >3      K      C-      •C
             O      O      O      O      O      O
              I        I       I       f       I      I
             w      u      u      u      u      u
ft           -      «•      t-      0>      fl      fl
o           »      »•      *      N      »      o
_:           «      ei      «      o      «      n
             o      o      «      a      a      e      r-
             —      oooooo

             U      U      hi      U      U      U      U
             O      -»««>      —      N
e           *      —      —      —      o      »      »
O           «      O      *      *      -a      *      o
             §or>    nn    en    *n    »n    »n    »n    »»
             oo>    oe    oo    oo    oe    oo    eo    oo
        4-     I +    14-     I *    14-     14-    14-     I4>    14-
O      u    uu    uu    uu    uu    uu    uu    uu    uu
e      o    fo    »o    no    »o    —o    no    —o    o —
O_    o    -to    -so    oo    — o    no    fie    fia    -    ei-    •*—    o—    »—    n-



             8     14^    14-     14-    (4>     14-    I4>     14-    14-
no,    u    uu    uu    uu    uu    uu    uu    uu    uu
o2    o    cio    ao    oo    o>o    oo    cia    -a    ao
EU    n    an    «n    «>n    «n    vn    an    en    nn
Bt-    o    ns    a>o    ao    eo    DO    oo    t-o    »o

      on    nn    en    cin    — n    no    —a    fin    an



             N—    »•—    t-o    i>-o*   >e»    w»    »o>    csa
             o—    a*-    a—    a—    oe    oa    oa    ee
              ii     it      ii     it      ii     it      ti     ii                ^~^
n           uu    uu    uu    uu    uu    uu    uu    uu               __
oo         fti    ao    no    Nn    e> r-    no    «o>    w                "-^
ci y         vri    «n    o a-    nn    —»-    «r-    «a    »»                    — —                C
                                                                                    i"
                                                                                    4->
             «e>    CO    aia    or-    at-    f-o    oo    nn                c
             oo    oo    oo    oo    oo    oo    oo    ao                Q

             uu    uu    uu    uu    uu    uu    uu    uu               O
                            on    no>    TO
                            as    -a         -ra    -oo
ou          (if    s> f i     aa    —a     e>a    — e>
EC          o «    * f i     o f i    — •»     no    fio

      oe     a—    ——     -»    n—     o'»    w-



        n     nn    on     tin    «n     fin    tin     -ft    -fi                 ?
        o     oo    oa     oa    oe     oo    30     oo    ee                 v*
        4     14-    14-     14-    t 4-     14-     14-     14-     14-                 GJ
        u     u u    uU     Uu    uu     uu    uu     uu    uiu                f&
  T    o     — o    •* o     no    90     oo    c^^     f 19    ^ t*
n=    a     no    ee     TO    oo     ao    ae>     Oe>    *«
ou    o     no    fie     fie    oo     QO    no     —e>    «a>                4J

                     »              -        -      *      »      »                 ^



              e>«    e*     »e>    an     5 ft    Sfi     a—    a—                 3
              oo    oe     oa    oo     oo    oo     oo    oe                O
              it     ii     it     it     ii      ii     ii      ii
  o          uu    uU     Uu    uU     uu    uu     Uu    f.i ^j
  a          n-r    nc>     «o    «\o     »•»    oo     on    fio

o =          oc^    on     on    oe     cl—    ao     ao    oo

      oo     fin    *»     a'fi    -N     -fi    —     r*n    aa     «o    on
z=  o       nci    se     sr    «-r     — ci    *N     o-r    rs —

      — o'     e« c i    » fl     * n    » *     — n    — N     — f I    — —
                                       E-5

-------
    t
                        u
                        e
                        n
                             n
                             i-
,  i
              s
              M
          S
          s
     I
    i

    to

e
1.
ss
1 +
uu

ss
uu
II
S3
I +

is
?!
uu
51
88
i ^>
uu
si
88
i +
uu
• a
    S3

sfe  si
         S3
     ss
     uu
                   oo
                   t +
                   uu
    ao

    »o
              ao

              -e
               as
               -e
                            uu

                            S3
    ss
n«
00
I. t
              oo
              I t
a*
oo

si

2   S-
                       uu

                       t»e<
                            uu
                                         TJ
                                          a>
    ss
    uu
             ss
             uu
          oo
          ii
                        *<•
     So
     ii
     ee
                   SS
                                          §
—n
oo
uu
oa
a»
88
uu
a v

on
eo
uu
—o
*•
(I«
ss
uu
S"
ei«
SS
MM
S?

88
uu
gs

                                          (A

 O  UU

 S  Sc~
 >=  ar>
10-3606
00-3610
              So

              uu
                  SS

                  uu
                  an
              s;
              • i
              uu
                            uu
                       e»o

                       aei
                            act
    uu
25  i?
    -w
         uu
         00
             uZ
             an
          o-

          iu

          S3
                  n-
     82
     i i
     uu
     oo

     §3
                            uu
                       oo  --
  PIFI
  eo

= uu
         net
         00
         uu
         s
-------
              o
              s
              e
              §
              S
              o
              §
              o
              o
              o
f.
                                                                         o>
                                                                         c
                                                                         0
                                                                         in

                                                                        &
— o
ii
                                                                        O
                                                            o
                                                            fl
                                                            N
O
un
                 o
                 ei
                             E-7

-------
* *
 *
                                     •o
                                     •o

                                     o

                                     o
                                     O.
                                     4J

                                     O
E-8

-------
                                  TECHNICAL REPORT DATA
                           (I'lcasc read Instructions on th.L reverse before completing)
1  REPORT NO.
   EPA-600/3-30-028b
                             2.
  TITLE AND SUBTITLE
  MODELING OF SIMULATED PHOTOCHEMICAL  SMOG  WITH KINETIC
  MECHANISMS   Volume 2.  CHEMK:  A Computer Modeling
  Scheme for Chemical Kinetics
                                                     5 REPORT DATE
                                                        February 1980	
                                                     6. PERFORMING ORGANIZATION CODE
                                                          3. RECIPIENT'S ACCESSION-NO.
7. AUTHOR(S)

  G.  Z.  Whitten and H. Hogo
                                                     8. PERFORMING ORGANIZATION REPORT NO.

                                                      EF78-107R
9. PERFORMING ORGANIZATION NAME AND ADDRESS
   Systems Applications,
   950  Northgate Drive
   San  rafael,  California
                    Incorporated
                      94903
                                                     10. PROGRAM ELEMENT NO.
                                                      1AA603   AC-054   (FY-79)
                                                     11. CONTRACT/GRANT NO.
Contract No.  68-02-2428
12. SPONSORING AGENCY NAME AND ADDRESS
   Environmental Sciences Research Laboratory-RTF, NC
   Office of Research and Development
   U.S.  Environmental Protection Agency
   Research Triangle Park, N.C. 27711
                                                     13. TYPE OF REPORT AND PERIOD COVERED
                                                       Final 7/78-9/79	
                                                     14. SPONSORING AGENCY CODE

                                                       EPA/600/09
15. SUPPLEMENTARY NOTES
16. ABSTRACT

        Mechanisms that describe  the  formation of photochemical  smog are developed
   using a computer modeling technique  directed toward the simulation of data collected
   in two smog chambers:  an indoor chamber and a dual outdoor chamber.   The results
   of simulating 164 different experiments are presented in Vol.  1.   Individual com-
   pounds for which specific experiments  were simulated and mechanisms developed
   include the following:  formaldehyde,  acetaldehyde, ethylene,  propylene, butane,
   and toluene.  Experiments in both  chambers were simulated  for  all these compounds.
   The mechanisms reported describe the decay of the precursor organic compound,
   formation and decay of secondary organic compounds, conversion of nitrogen oxides,
   formation of nitrates, and the appearance and decay of ozone.   Special emphasis is
   given to the chemistry of toluene.   Also included is a study  of a generalized
   smog-based or carbon-bond mechanism  developed in a previous study.   Volume 2
   contains the user's manual and coding  for a chemical kinetics  computer program,
   CHEMK.
17.
                               KEY WORDS AND DOCUMENT ANALYSIS
                  DESCRIPTORS
                                              b.lDENTIFIERS/OPEN ENDED TERMS
                                                                   c. COSATl Field/Group
   *
   *
   *
   *
   *
   *
Air pollution
Reaction kinetics
Photochemical reactions
Test chambers
Mathematical models
Computerized simulation
               13B
               07D
               07E
               14B
               12A
               09B
13. DISTRIBUTION STATEMENT

   RELEASE TO PUBLIC
                                        19. SECURITY CLASS (ThisReport!
                                         UNCLASSIFIED
            21. NO. OF PAGES
               113
                                                           JS (This page)
                                                                         22. PRICE
EPA Form 2220-1 (9-73)

-------