-------
The work upon which this publication is based was performed
pursuant to Contract No. EHSD 71-22 with the National Air
Pollution Control Administration, Environmental Health Service,
Public Health Service, Department of Health, Education and
Welfare.

-------
                                CONTENTS
SECTION     	     PAGE
  I         INTRODUCTION                                             1
            A.    Perspectives on Modeling Needs                     1
            B.    Overview of Previous Contributions                 1
            C.    Advances Made in the Present Work                  5
  II        KINETIC IMPROVEMENTS IN THE PHOTOCHEMICAL DESCRIPTION    7
            A.    Additional Oxidation Chains                        7
            B.    Influence of Initial Composition                  14
  III       MATHEMATICAL AND METEOROLOGICAL IMPROVEMENTS            21
            A.    Application of Fade" Approximants to Smog
                  Chamber Calculations                      .        22
            B.    Application of Pade1 Approximants to the
                  Atmospheric Model Equation                        26
            C.    Meteorological Reformulation of Model             30
  IV        TESTS OF THE MODEL AGAINST LOS ANGELES BASIN DATA       38
            A.    Objectives of the Tests                           38
            B.    Review of the Previous Validation Study           38
            C.    Sensitivity Studies on 1969 Trajectories
                  with the Expanded Model                           40
            D.    History Analysis of the 1030 El Monte
                  Trajectory                                        45
  V         CONCLUDING REMARKS                                      51
            REFERENCES                                              57

-------
ii

-------
                              ILLUSTRATIONS
NO.
1
2
3
4
5
6
7
8

Schematic of Diffusion Model for Air Pollution Simulation
Expanded Kinetic Mechanism
Computed Curves Compared With Experimental Points
Hydrocarbon Decay With Two Reactivities
Ozone Buildup With Two Reactivities
Influence of NO :HC-Ratio on Oxidant Production
X
Estimated Trajectory of Air Mass Arriving at El Monte,
1030 hours, 29 Sept. 1969
Estimated Trajectory of Air Mass Arriving at El Monte,
1130 hours, 29 Sept. 1969
PAGE
4
9
12
15
16
17
33
34
9     Estimated Trajectory of Air Mass Arriving at El Monte,
      1230 hours, 29 Sept. 1969                                     35

10    Estimated Trajectory of Air Mass Arriving at El Monte,
      1330 hours, 29 Sept. 1969                                     36

11    High Reactivity Hydrocarbon Concentration at Huntington
      Park on Type 2 Days, 1968.  f = 1/4; r = 1/2                  39

12    Ozone Concentration at Huntington Park on Type 2 Days,
      1968.  f = 1/4; r = 1/2                                       41

13    Reactive Hydrocarbon History Along the 1030 Trajectory.
      f = 1/4; r = 1/2 (Ground level concentrations)                46

14    Nitric Oxide History Along the 1030 Trajectory,  f = 1/4;
      r = 1/2 (Ground level concentrations)                         47

15    Nitrogen Dioxide History Along the 1030 Trajectory.
      f = 1/4; r ° 1/2 (Ground level concentrations)                47
                                                                     iii

-------
ILLUSTRATIONS  (Cont.)
NO.                                                                 PAGE
16     Oxides  of  Nitrogen (NO + N02>  History  Along  the 1030
       Trajectory.   f  =  1/4;  r = 1/2  (Ground  level  concentrations)    48

17     Ozone History Along the 1030 Trajectory.   f  =  1/4;  r =  1/2
       (Ground level concentrations)                                  50
iv

-------
I.    INTRODUCTION

A.    PERSPECTIVES ON MODELING NEEDS
      Control of emission sources and of future land use offer possibili-
ties for managing air quality.  Under statutory mandate, regional imple-
mentation plans must be submitted to show how these controls will be
applied.  Mathematical models that give estimates of air quality in terms
of emissions in a particular airshed should play a key role in formulating
the plans.

      The work described in this report has resulted in a diffusion model
with photochemistry that is applicable to pollutants that undergo chemical
changes.  Although our further development of an earlier model  is still
confined to the Los Angeles basin data, the method is applicable to the
broad class of problems.  Most of the pollutants covered by criteria
documents participate in some forms of chemical reaction on time scales
comparable with air residence times in urban areas.  To account for this
we have retained both the meteorology effects of dispersion and advection
coupled with photochemical effects of primary pollutants producing oxidant
compounds.

B.    OVERVIEW OF PREVIOUS CONTRIBUTIONS
                              234
      Since literature reviews ' '  of air pollution modeling in general
are available, our remarks here are confined to the previous work that
was specifically directed toward reactive pollutants.  The discussion
centers on the relationship between the details emphasized by the model
and the levels of precision in either input or output data.  At the very
least, any model must include the physics and chemistry of the governing
phenomena,  but complexities which far exceed reasonable expectations of
specific validation are probably unjustified.  In abstract language this
says that we can fit any observed data by using enough adjustable param-
eters in the model, but that we leave the realm of physics if we do not
subject most of the parameters to independent experimental verification.

-------
      The models described by Calvert  and Wayne, et al.,   stress chemical
kinetic mechanism.  In the former, the photochemical smog  mechanism is
reduced to 17 steps and in the latter to 31 steps.  In the larger mecha-
nism care is taken to describe plausible reaction for the  C-H../NO /air
system, but it is ultimately applied to the atmospheric case where dozens
of smog-forming hydrocarbons are involved.  To preserve a  consistent level
of detail, the number of reactions must be increased by a  large factor;
however, this is avoided by adjusting the rate constants to obtain observed
results.  Homogeneously stirred air parcels following wind trajectories
are assumed in order to retain first-order, ordinary differential equations.

      Another way of restricting the mathematical complexity to first-
order, ordinary differential equations is given in Friedlander and
Seinfeld's paper.   This method reduces the chemical system to seven
functional steps to be consistent with the level of observational detail
                                           A
usually available for polluted atmospheres.   Unlike the models discussed
above, this one accounts for finite rate atmospheric diffusion.  Retention
of the mathematical simplicity is bought at the price of imposing the
Lagrangian similarity hypothesis.  This means that the vertical profiles
of concentration in the polluted layer belong to simple parametric fami-
lies.  Careful choice of the profile form functions may well allow this
method to be a fast and powerful computational tool.

                                               8 9
      Three-dimensional, time-dependent methods '  have been proposed
recently, but results for reactive atmospheres have not been reported at
this writing.  Simplified chemistry must be employed in each of these
approaches because of the emphasis on details of advection and diffusion.
The body of meteorological data for most air basins falls far short of
the input requirements for any transport formulation of this complexity.
Here again it is difficult to avoid the problem of allowing too many
*
 Even if more detail were included, there are no measured data on rate
 constants for most of the additional reactions.  Consequently, the model
 is provided with additional degrees of freedom "wherever the kinetic
 data are open to estimation".

-------
unspecified parameters to obscure the physically based portions of the
calculation.  One of the new methods employs an Eulerian coordinate frame
which introduces the chronic problem of artificial (or numerical)
diffusion.  Small time and space steps are the only means suggested
heretofore to remedy this problem.  The other method uses a discretized
representation of species concentrations in a Lagrangian frame.  For
multicomponent systems, the latter type of approach has previously been
limited because of the large quantity of rapid access memory storage
needed for each species.  Both of these problems have been recognized
for some time with three-dimensional, time-dependent calculations.  More
powerful computers and novel algorithms may well hold their solutions.

      The method developed at General Research Corporation is more modest
in its chemical detail than the models first described above, and it is
far simpler in its transport formulation than the three-dimensional,
time-dependent approaches.  An attempt has been made to keep the mathe-
matical descriptions near the level of detail and precision of observa-
tional data.  This has resulted in a lumped-parameter, chemical kinetic
mechanism and a Lagrangian air parcel with finite-rate vertical diffusion.
The chief advantages of this over the approach of Ref. 7 are complete
flexibility in specifying vertical diffusion and in accounting for
distributed emission sources.

      Our approach avoids the artificial numerical diffusion problem by
the choice of natural (or so-called "intrinsic") coordinates.  We align
the coordinate system with the natural airflow at each instant so that
there is practically no advection across any grid line.  This permits us
to follow air masses, to include vertical dispersion, and to account
for chemical changes simultaneously.  Figure 1 schematically illustrates
these points.  The large memory requirements do not appear in our approach
because we trace the history point-by-point and output the results sequen-
tially.  Balancing the needs for realistic chemical and meteorological
descriptions against efficiency of calculation, the GRC model has achieved

-------
                         SUNLIGHT IS GIVEN AS
                         A FUNCTION OF TIME
00
c^

-------
      It is useful to digress here to examine the relative magnitude of
the horizontal diffusion terms which are neglected.  Carbon monoxide will
be chosen because it undergoes negligible atmospheric reaction over time
scales of interest.  We shall form a ratio taken from the CO-balance on
the moving air parcel.   The ratio of CO-loss horizontally to CO-gain
vertically gives a measure of how serious our assumption is.  The air
                          2
parcel is taken to be 9 km  in horizontal area because this approximates
the spatial resolution of the source grid.  The vertical dimension is
taken to be 100 meters, typical of a high pollution day.  In order to
force a high value for the ratio, the time of 1000 hours PST is selected
because the numerator of the ratio is high (the air having just been
loaded up by morning traffic) and the denominator is low (peak traffic
has passed).  Horizontal loss is the turbulent diffusion coefficient
(-V300 m /min) times the average crosswind gradient (^3.8 x 10    mass
                                                               f  9
fractions per meter) times the area of all four sides (1.2 x 10  m ).
The horizontal gradient is maximized by using a concentration decrement
equal to the highest 10 A.M. reading in Huntington Park in 1969, which
was 19 ppm.  The characteristic distance over which the concentration
drops from 19 ppm to zero is taken as the order of 30 miles.  This gives
                     3
a loss of about 0.1 m /min.  Using a grid square in the central basin
(containing both freeway and surface street contributions), we get a
                  3                           2
gain of about 20 m /min averaged over the 9-km  area.  Consequently, the
error committed is of the order of 0.5% in the CO-mass balance.  We
consider this to be a high estimate for the reasons outlined above.
C.    ADVANCES MADE IN THE PRESENT WORK
      Following the philosophy expressed above, we have advanced our
techniques to the stage where they can be adapted to planning and imple-
mentation of air pollution abatement measures.  As described in the next
section, the photochemical description has been extended to account for
varying ratios of oxides of nitrogen to hydrocarbon.  Some hypotheses of
carbon monoxide interactions have been tested.  These kinetic improvements
add several reaction steps to our previously adopted mechanism.  With no

-------
changes in numerical methods, the costs in computing time would have
grown roughly in direct proportion.

      Numerical integration improvements have more than offset the
increases in chemical complexity.  Another section describes these along
with refinements in the meteorological realism of the model.  The moving
element shown in Fig. 1 replaces our previous formulation which was
limited either to a horizontally uniform slab or a two-dimensional
Eulerian coordinate system.  The increased computing efficiency applies
to both the homogeneous (smog chamber) and the diffusion (atmospheric)
runs because we generalized it to apply to partial differential equations,

      In the final sections of the paper, validation tests are presented
for Los Angeles basin data and the significance of the results are
discussed in terms of future applications of the technique.

-------
II.   KINETIC IMPROVEMENTS IN THE PHOTOCHEMICAL DESCRIPTION

A.    ADDITIONAL OXIDATION CHAINS
      It has been known for some time that hydrocarbon consumption rates
observed in chamber experiments cannot be explained by ozone and oxygen
atom attack alone.  This is discussed in the review by Altshuller and
Bufalini   where it is noted that Schuck and Doyle   termed the disparity
an "excess rate."  Our previous modeling work treated this by increasing
both the 0_ and 0-atom rates of reaction with hydrocarbon.

      To use known rate constants to a maximum extent in the present work,
we have added a hydrocarbon oxidation chain to reflect the attack of  RO
                                                          12 13
upon the hydrocarbon.  Because of its suspected dominance,  '   the
hydroxyl radical (OH) was assumed to be the only RO reacting with the
hydrocarbon.

      Following the inorganic cycle formed by

            hv + NO,, -»• NO + 0
                                                                     d>*
            NO + 03 -> N02 + 02                                       (2)

we have allowed each hydrocarbon oxidation to generate multiple  R0_
radicals as expressed by the b-factors in the steps
            0 + HC -v b1(R02)                                         (3)
 The two reaction steps bracketed by (1) are combined by assuming 0-atom
 stationarity internal to the logic of the computer program.  Hence
 reaction (1) nets the system hv + N02 -»• NO + 0».

-------
            OH + HC -* b2(R02)                                         (4)
            03 + HC -»• b3(R02)                                         (5)

Subsequent conversion of  NO  to  NO,,  was assumed to occur via

            R02 + NO •* N02 + d(OH)                                    (6)

where  d  expresses that  fraction of the conversion responsible  for
returning hydroxyl to the system.

      Note that reactions (4) and (6)  form something of a closed loop
and that a stability requirement like  b_d < 1  may be needed  to prevent
R02~runaway.  This is not a precise requirement because of the variety
of external factors that  influence the flow rate of free radicals around
the loop.  Some of the  R0_  is formed by the other two reactions, and
there are radical removal processes in the chain termination steps we
discuss in the following  paragraph.

      Both  R0_  and  OH  are held in  check by removal steps that
terminate the chains.  We continue to  use some lumped parameter  reactions
such as

            R02 + N02 -> c(PAN)                                        (7)

and some elementary reactions like

            OH + NO + M -> HONO + M                                    (8)

            OH + N02 + M  -»• HN03 + M            '                       (9)

This rounds out the mechanism we have  called Mark XII.  It is  diagrammed
on Fig. 2.
 PAN is used to  denote  peroxyacetylnitrate.

-------


he
i
i ^
NO * 1 *
1
-
*






1
j
^-t— , '
HC ' 1— ^
1
1 	 J






'



>2+M


;
, i

J












3


r

^


r~
i

i ^~~

NO2

|










^U

V





V

P AKI









^
y

(C








^

,

l
! I
^ i
i i

• HI








r

[























R(
(HC

|














3
V

t







i
1
|






2
,)

1







k


Figure 2.  Expanded Kinetic Mechanism

-------
      In summary,  the main  changes  from the previous scheme  are:
addition of radical  species  RO   (treated as OH), the provision for a
multiplicity  of hydrocarbons (not shown in equations), and the elimination
of ozone production  from an RO-  reaction.  The multiplicity of hydro-
carbons is done by adding parallel  reactions to  (3), (4), and (5).  The
omission of the ozone production  reaction is justified on the basis of
                             14
its suspected endothermicity.

      For our initial validation  tests we tried  the rate constants com-
                                          14
piled and estimated  by  Westberg and Cohen.    Stoichiometric coefficients
b.. , b~  , and   d  were determined by hand-calculation analysis of some of
the propylene data represented by Ref. 15.  We  found that  bd = 1  where
b  is a composite  branching factor, and  b = 2   in order to explain the
decay rate of NO  compared with  that of propylene.

      Explosions of  R02~concentration characterized some early choices
of bd-combinations because  of the positive feedback loop involving radical
reactions shown in Fig.  2.   Since these occurred in the first few seconds
of real time, they were difficult to detect.  Our usual printing interval
of 1 minute of real  time gave computer outputs  showing sudden drops of
NO  and  HC   (to some nonzero levels) during the first minute followed
by the dynamic equilibrium  between   NO, N0»  , and  0_  with a slow decay
                                          £•         J
of hydrocarbon.  The underlying pathology was discovered by rerunning
with special  diagnostic outputs.

      Nominal rate constants were adjusted to reconcile computed results
with measured values.   For  example, the  (OH + HC)-rate was cut to about
one-third the estimated value and the (NO + R0»)-rate was increased
eightfold.  Reaction rates  which  have been reported individually in the
literature were held at their nominal values including the  (0 + HC)-rate
which we had  to increase many-fold  in previous  work.   Also, the reten-
tion of reaction  (5) was still necessary to describe the continued
C0H,-decay after the near disappearance of  NO.
  J o
10

-------
      In order to reproduce the late-time behavior of  N0_  and propylene,
we included chain termination reactions in addition to those indicated in
Fig. 2.  In their more recent review,   Altshuller and Bufalini point out
that  HONO  might be formed by the reaction with water vapor
            NO -4- N02 + H20 -> 2HONO                                    (10)
which is likely to proceed in the two steps
            NO + N02 -» N203                                           (11)

            N2°3 + H2° ~" 2HONO
                               12 13
A possible source of OH-radical  '   is the photodissociation of  HONO
            hv + HONO ^ OH + NO                                       (13)

      Figure 3 shows the computed concentration history using this scheme
to describe the experimental observations of Altshuller, et al.,   and
Table 1 summarizes the reaction rate constants that were used.  The column
in Table 1 labeled "nominal" gives the values we started from on the basis
of tabulations in Ref. 14.  Following an empirical procedure, our calcula-
tions effectively combined H-0-concentrations with the rate constant  for
reaction  (10) as indicated; therefore no comparable nominal value is  shown.
Higher levels of precision than that shown in Fig. 3 were not among our
objectives, since a detailed investigation of propylene is not justified
for describing the potpourri of hydrocarbons in the urban atmosphere.
Orders of magnitude for conversion times and levels are sufficiently  good
to proceed from here.  The objective that we have fulfilled is to repre-
sent the so-called "excess rate" of hydrocarbon oxidation with added
reactions instead of artificially increasing the 0-atom rate with
hydrocarbon.
                                                                      11

-------
2.5 ppm -
    2.0 -
Z
o
o
Z
o
                                       A  PROPYLENE

                                       •  NITRIC OXIDE
                                                            > ruirN is ur      ,.
                                       •  NITROGEN DIOXIDE  j ALTSHULLER, et ol "

                                       T  OZONE
                                60       80      100     120      140      160 min
                                        TIME
       Figure  3.   Computed  Curves Compared With  Experimental  Points
12

-------
                                    TABLE 1

                          RATE COEFFICIENTS FOR MODEL
                      HYDROCARBON/NITRIC OXIDE MECHANISM

                      (Stoichiometry imbalances may occur
                   because of lumped parameter assumptions.)
       Reaction
                              Model Values
                            From Validation
                                                      Nominal Values  for
                                                      Propylene  System-^
hv + NO- -> NO + 0
                         0.4 min"1

                         1.32 x 10"5

                         40 ppm  min

                         6100 ppm  min

                         80 ppm  min

R02 + NO -»• N02 + 0.5 OH  1500 ppm'-'-min"1
                              _i   _i
03 + NO

0 + HC
OH + HC ^ 2RO,
             2
R02 + N02 •> PAN

OH + NO -»• HN02*

OH + N02 •*• HN03*

0_ + HC -»• RO^
NO
           2HNO,
                **
                         ,    -1  . -1
                         6 ppm  min

                         10 ppm  min

                         30 ppm  min

                         0.0125 ppm" min"

                         0.01 ppm  min
hv + HNO  ->- NO + OH      0.001 min"1
                                                  0.4 min
                                                         "1
                                                           "5
                                                  1.32 x  10
22-44 ppm  min

6100 ppm  min

244 ppm  min

122 ppm  min

122 ppm  min

99 ppm  min

300 ppm  min
                                                  0.00927 - 0.0125
  Rate constant lumps third body concentration
**
  Water vapor lumped into rate coefficient
                                                                          13

-------
B.    INFLUENCE OF INITIAL COMPOSITION

1.    Synergism Between Hydrocarbons
      The fact that free radicals produced in the photooxidation of one
hydrocarbon can accelerate the attack on a coexistent hydrocarbon has
been illustrated with experimental findings by Altshuller and Bufalini.
Because of the limited understanding of detailed interactive mechanisms,
we limited our multiple hydrocarbon investigation to a simple parametric
study.  Using a truncated version of the mechanism in Table 1, we added
a parallel series of oxidation steps.  This scheme is shown in the flow
chart of Fig. 2 with dashed  lines along the bottom portion of the diagram.

      In the parametric study, the same initial values of total hydro-
carbon and nitric oxide as those in Fig. 3 were employed.  The hydro-
carbon, however, was partitioned into two hypothetical compounds one
having triple the (OH + HC)  rate constant and the other having one-third
the  (OH + HC) rate constant  as propylene.  The compound with the tripled
rate is called species  "B" and the one with the decreased rate is called
species "A".  All other reactions in the scheme are the same in all
respects as their counterparts in the (C-H, 4- NO )-system.
                                        JO     X

      Figure 4 shows the hydrocarbon decay curves from the parametric
study.  On the same coordinates we plot curves for model calculations of
pure hydrocarbon "A", pure hydrocarbon "B", and pure propylene for com-
parative purposes with  the half-and-half mixture of hypothetical compounds,
The most interesting aspect  of Fig. 4 is that the mixture decays faster
than propylene even though it contains equal quantities of "A" and "B"
which bracket propylene symmetrically.  The coupling between the two
hydrocarbons occurs via the  R0_  conversion of  NO  to  N0»  and the
production of  OH .  Evidently the combination of the rapid incubation
(afforded by the reactive hydrocarbon) with the twofold chain-branching
is sufficient to increase the decay rate over that of pure propylene.
14

-------
       SPECIES 'A1 HAS AN (OH + HCJ-RATE ONE-THIRD THAT OF C3H6
      SPECIES 'B1 HAS AN (OH + HCJ-RATE THREE TIMES THAT OF C3H6
            ppm
Z
o
t-
<
at
i-
Z
ui
Z
O
y
u
           2.0
           1.5
           1.0
           0.5
                                PURE
                                            'A-
                                                                CVJ
                                                                 I
                                     I
                                    I
                  I
                   I
                   I
                   20
                  40
60
 80
TIME
100   120   140   160 min
           Figure  4.  Hydrocarbon Decay  With Two Reactivities

      Figure 5 shows corresponding  curves for ozone buildups predicted
in the computer experiment.  Again  the combination seems to exceed the
pure propylene in its rapidity to produce ozone.  This occurs because
of the enhanced conversion rate of  nitric oxide to nitrogen dioxide
stimulated by the early abundance  of  RCL  generated by the high-reactivity
fraction.

      These findings illustrate a plausible mechanism for the interaction
of multiple hydrocarbons.  In future applications of our photochemical/
diffusion model this type of representation is likely to be more realistic
than a super-detailed single hydrocarbon mechanism with many adjustable
constants.  Also,  other radical interactions like  R09 + HC  and  RCO, + HC
                                      16
might provide plausible coupling paths   in addition to the hydroxyl
attack on the hydrocarbon.  Indeed, other forms of coupling may favor the
slower reacting hydrocarbon in a similar numerical experiment with a
half-and-half mixture.
                                                                      15

-------
             HYDROCARBON 'A1 HAS AN (OH + HC)-RATE ONE-THIRD
             THAT OF C3H6
             HYDROCARBON 'B' HAS AN (OH + HCj-RATE THREE TIMES
             THAT OF C3H6
             ppm
            I.Oi-
         Z
         o
            0.5 -
         «J
         z
         o
         u
          CO
         O
               0     20    40   60    80   100   120   140   160 min
                                   TIME
             Figure 5.  Ozone Buildup With Two Reactivities

2.    Modeling the  Influence  of NO  /HC-ratio  on  Oxidant Production
      Our main purpose  in including the  chain termination details for
nitrogen oxides  in  the  previous subsection was to  improve the modeling
of  initial mixture  effects.   It will be  recalled that  our earlier report
issued a caveat  restricting the seven-step mechanism to a narrow range  of
NO  /EC  ratios.   Later, Seinfeld    pointed out that the slope of the
curve of peak  oxidant versus  NO /HC-ratio did not  even have the correct
sign.  If modeling  is to  be used  to evaluate  hypothetical abatement
strategies,  it is clear that  this deficiency  must  be remedied.

                                                     *
      Figure 6 shows  the  results  from the twelve-step  mechanism given
in  Table 1.  Experimental results for pure and mixed hydrocarbons are
shown on this  plot  of peak oxidant  versus NO  /HC-ratio.   The model
                                            X
mechanism displays  a monotonically  declining  trend of  oxidant production
as  the oxides  of nitrogen increase.   Both of  the pure  hydrocarbons
 Stationary  state  treatment of 0-atoms  reduces  this  to  eleven.
16

-------
z
q
t—
<
ec
i—
Z
o
u
»-
z
<
o
x
o
                                     2 ppm C3H6

                                     REF. 18
         3 ppm n-C4H1Q
                     2 ppm HC

                     MODEL

                     RESULTS**
                1 TO 2 ppm NMHC*
                AUTO EXHAUST REF. 19
                             LOS ANGELES
                              INVENTORY
                                RATIO
              AUTO EXHAUST\
               2 ppm NMHC*   \
                  REF. 20
                                                  CO
                                                  to

                                                  CXI
                                                  I

                                                  'j;
   0.4 -
    0.2 -
                                          \
                                            \
      0            0.5            1.0            1.5


         MOLE RATIO NOX TO HYDROCARBON


   *NMHC:NON-METHANE HYDROCARBON

  **THIS CURVE IS LOWERED CONSIDERABLY BY OUR
    SUBSEQUENT ADOPTION OF ONE-HALF THE
    PROPYLENE RATES  FOR ATMOSPHERIC CALCULATIONS
Figure 6.  Influence of NO :HC-Ratio on Oxidant Production
                                                          17

-------
 measured by  Altshuller,  et al. ,  show peaks  at  widely  separated  values
 of  the mixture ratio.   It seems  plausible  then that  the  typical mixture
 of  atmospheric pollutant hydrocarbons should show a  decrease  of oxidant
 with increasing NO /HC-ratio.  This  is supported  by  the  fact  that  the
 relative abundance of  hydrocarbon species  decreases with increasing
 reactivity.   Now since the low reactivity  species peak at lower
 NO  /HC-ratios the combined effect should give  a negative slope.
  X

       The two plots on Fig.  6  of oxidant production  from diluted and
                              19  20
 irradiated automobile  exhaust   '  both have negative slopes  comparable
 to  that of the model.   As stated earlier,  our  purpose in using  chamber
 data for propylene to  validate the mechanism is to capture the  main
 features of  the experiments that apply to  the  atmosphere.  If,  on  the
 other hand,  propylene  were the subject of  our  study,  we  would proceed
 to  expand on the mechanism to  get the proper curve shape shown  in  Fig. 6.

       Having subjected the mechanism to these  laboratory tests, we use
 it  in its present form for atmospheric validation studies described in
 a later section.

 3.     Examination of Carbon Monoxide Effects
                                    13 21 22
       Recent work in photochemistry   '  '   of smog  has  suggested  that
 carbon monoxide may play a role  in accelerating the  photooxidation of
 hydrocarbon/nitric oxide mixtures.  The mechanism suggested is
             OH + CO ->• CO. + H
             H + 02 + M -> H02 + M                                     (15)
             H02 + NO ->• OH + N02                                      (16)
                                                                      (17)
18

-------
Reaction (15) is postulated to contribute  H02  rapidly as it responds
to the presence of H-atom, and reaction (16) enhances the NO-to-NO.
conversion already occurring in reaction (6).  Likewise, the hydroxyl
radical supplied by reaction (16) goes on to react with the hydrocarbon
in reaction (4), producing R02~radicals in a branched chain already
described.  Reaction (17) is a chain termination step suggested for  H02.

      The key feature we wish to examine in this alternate path is the
competition of  CO  with hydrocarbon for the hydroxyl radicals.  The
rate constant for  OH + C-H,  suggested in Ref. 14 is approximately equal
to that suggested in Ref. 21 for  OH + CO .  We tried a wide range of
(OH + CO)-rates holding the basic system at the values in Table 1.  Since
reaction (15) is fast, reactions (14) and (15) were added to give the
overall reaction

            02 + OH + CO -»• C02 + H02                                 (18)

The rate constant for reaction (16) was taken to be the same as its
RO..-analogue, reaction (6), and the rate constant for reaction (17) was
                             -1   -1
assigned the value of 1.2 ppm  min   based on the estimate in Ref. 14.
The species  02  is carried implicitly in the calculation, and the rate
constant for reaction (18) is varied parametrically for this study.

      Beginning with (CO + OH)-values at the level suggested in Ref. 21,
(k..£ - 200 ppm  min  ) we found very rapid NO-to-N09 conversion and  HC
disappearance with 100 ppm added  CO .  Only when  k,0  was lowered to
          -1   -1
1 or 2 ppm  min   did we obtain only small effects.  Since much experi-
mental verification is yet to be done, only a limited effort was devoted
to this aspect of -the study.  Table 2 summarizes the main findings.  It
indicates that the half disappearance time of  NO  is the strongest
effect of  k-R  on the system, and the half-disappearance time of hydro-
carbon is least affected.  Both the peaking time for  N02  and the half-
rise time of ozone are moderately sensitive to the rate constant.
                                                                      19

-------
                                 TABLE 2
           INFLUENCE OF  100 ppm CARBON MONOXIDE ON THE SYSTEM
            SHOWN  IN FIG.  3 WITH VARIOUS VALUES OF THE RATIO
            OF  (OH + HC)/(OH + CO) REACTION RATE COEFFICIENTS
k4/k!8 tl/2 ~ N0(mi
•k
0.4 <1
4 6
40 25
80 29
System in
Fig. 3 36
(No CO)
n) t , - N02(min) t.

4
13
43
52

66

L/2 ' HC(m±n)

40
48
73
81

92

t1/2 - 03(min;

10
22
62
73

83

 This  is with   k..g  set at the estimated  value  from Ref.  21.

       In any event,  whether we use  the  estimated  rate  of  (OH +  CO)-
reaction or even  only  preserve its  ratio  to  that  of the  (OH + HC)-reaction
                                                            12
suggested  in Ref.  14,  we get  very large effects.   Other work    suggests
a much larger  (OH +  HC)-rate  constant for propylene.   In  view of our
findings and of the  wide disagreements  on some  of the  rates, a  very real
need exists for definitive experimentation to achieve  a  deeper  under-
standing of the kinetics.
20

-------
III.  MATHEMATICAL AND METEOROLOGICAL IMPROVEMENTS
      Our previous study brought out some of the limitations which hamper
the original time-dependent atmospheric model.  This model is based on a
                                                          23
slab geometry and is solved by a Crank-Nicolson technique.    The omission
of advection limited the calculation to places and times characterized
by stagnant conditions (hence our choice of the October 1968 episode for
Huntington Park validation data).  Admittedly, this is a case that the
plume models (described in Refs. 2, 3, and 4) are incapable of handling,
but applications demand that cases with wind transport be treated.  Our
original response to this need took the form of an Eulerian governing
equation for species concentration:
It
                 - £ - * (» If ) + •
where       c = concentration
            t = time
            x = downwind distance
            u = wind speed
            z = height above ground
            D = diffusion coefficient
            R = chemical reaction rate
Solutions of Eq. 19 are costly because of the need to suppress artificial
         *
diffusion  by refining mesh and step size.  Smog chamber analyses use
just the first and last terms in Eq. 19 so that they depend on ordinary
differential equations.  These are solutions which describe the time-
dependent behavior of a homogeneous gas mixture.  We used standard Runge-
Kutta techniques to solve them in the original code.
 Artificial diffusion is an error from numerical differencing of the
 equation.  It produces effects that resemble physical diffusion in
 the x-direction.
                                                                      21

-------
      The combined needs  for computing efficiency and meteorological
realism have been met by  improvements described in  this section.  The
first two subsections describe  the  incorporation of a classical approxi-
mation method into our numerical integration.  The  implementation of
these computing  techniques  strengthen the study by  allowing numerous
modeling trials  on both  laboratory  and atmospheric  cases.  The final
subsection outlines  the  change  from Eulerian  to semi-Lagrangian coordi-
nates in the meteorological formulation.  By  following air masses we use
a natural coordinate system thereby eliminating the artificial diffusion
errors.

A.    APPLICATION OF FADE APPROXIMANTS TO SMOG CHAMBER CALCULATIONS
      The chemical model  previously described is implemented in rate
equations which  must be  integrated  numerically.  Because of the widely
varying reaction rates,  the equations possess a characteristic known as
"stiffness" which makes  their integration difficult.  In the interest of
economy and accuracy, special techniques must be employed in obtaining a
                                                                       24  25
solution.  Below we  describe a  method which utilizes Pade approximants.   '
For these ordinary differential equations, the procedure has yielded
significant speed improvements  over the Runge-Kutta method previously
used.

      For this class of  problem, the implicit approach is the basic
advantage of the new method, thus guaranteeing stability.  Consequently
accuracy controls step size.  This  provides the principal gain in
computing time.

      The rudiments  of the  integration method are twofold.  One involves
a linearization  of the nonlinear chemical term at every step; the other
an approximation of  the  exponential term that appears as a result of the
linearization operation.  Thus, consider the  following system of equations
that describe the changes in species concentration  due to chemical reaction:
22

-------
                                                                     (20)
where       c = vector of species concentrations
         R(c) = vector of chemical reaction rate functions
            c =
V
C2
•
cs


R(c) =

\M~
R2(c)
•
R8(c)
            s = number of chemical species
To linearize the rate vector  R(c) , we expand it about  c   in a Taylor
series.  The subscript "o" denotes evaluation at some time  t = t   .  Thus
            R(c) = R(CQ) + J(c - CQ) + 0(Ac)'
                                                      (21)
where
                            -J — 1  O          Q
                           . I J — _L J ^ J  •  •  *  > a
is the matrix of order  (s * s)  of first derivatives  (the Jacobian
matrix) of the rate function.  Substituting Eq. 21 into Eq. 20 yields
            dc
            dt
= Jc + B
(22)
where  B = R(c ) - Jc   .  It should be noted that at  time  t =  t    the
              o      o                                          o
quantities  c , R(c ) , and  J  are known.  Hence,  the matrix differential
                                                                      23

-------
equation  (22)  is  a  linear  equation with  constant  coefficients  and a
constant  forcing  function.   Its  formal solution is,  therefore,

             c(t)  =  eJAt  [CQ  +  J~VJ  - J^B                            (23)

where  J     denotes the  inverse  (or  reciprocal) of  the matrix  J.
       An  approximate  solution  of  Eq.  22  is  obtained  from  Eq.  23 by
suitably  approximating the matrix exponential   e     .  This is accom-
plished by means  of  the Pad£ approximants of the  exponential  function.
These  Fade approximants are rational  functions  of the  form

              JAt  = P(JAt)  o -1
             6     Q(JAt)    ^   V

where  P(JAt)  and  Q(JAt)  are matrix polynomials of  degree  p   and  q  ,
respectively.  They  are defined by
                       (p + q  - k)!p!        xk                         f    .
                                                                      (25a)
                     (p+q)!k!(p  -k)!
                 k=0
                  q
             n  -       (P + Q  ~ k)!q!     f     ,k                       f    ,
             Q  ~       p + q)!k!(q  - k) !  ("JAt)                         (25b)
                 k=0
Thus,  substituting Eq.  24 into Eq.  23 we  obtain
             c(t)  =  Qc^ + JB  -  JB                             (26)

The  choices   p  =  q  = 1  and  p =  q  =  2   are  especially  useful.   The
former  yields a third-order integration formula,  the latter  a  fifth-order
formula.   As an example,  for  p = q = 1 we  obtain from Eq.  25:

             P = I + 1/2 JAt      and       Q  =  I  - 1/2 JAt            (27)
24

-------
where  I  is the identity matrix.  Substituting Eq. 27 into Eq. 26 yields,
after the appropriate manipulations,

            c(t) = c  + [I - 1/2 JAt]~1AtR(c )                       (28)
                    o                       o

Equation 28 has a form that is especially well-suited for digital computa-
tion.  Note that the process requires at every time step:  (1) the evalua-
tion of  R(c)  and  J , and (2) the inversion of the matrix  Q .  The
latter is by far the more computationally costly of the two requirements.
In our case, the number of species determines the size of the matrix to
be inverted.  Since this number is not large in our model, the matrix
inversion poses no special problems.  Fortunately the number of species
is usually smaller than the number of reactions.

      Table 3 shows examples of the gains obtained using the new computa-
tional scheme.  The typical real-time/computer-time ratio was increased
from 36/1 to 180/1.  Perhaps more significant is the fact that the Fade"
method has allowed us to obtain acceptable solutions in situations where
Runge-Kutta either failed to converge or produced spurious solutions.
One such instance is the integration of full differential equations for
the free-radical species:  the Fade" method successfully computed solutions
without algebraic substitution of stationary-state assumptions whereas
Runge-Kutta failed to produce any solution.

      As mentioned previously, the step size must still be controlled
in order to preserve accuracy.  Referring to Eq. 21 we can see that the
difference  Ac = c - c   determines the order of the truncation error of
                      o
the linearization operation.  The step size must be constrained so that
Ac  will remain within some bound specified by the user.  Our codes con-
tain a variable-step-size feature which ensures that  Ac  will stay within
the prescribed bounds.
                                                                      25

-------
                                 TABLE 3
                 CDC 6400 CENTRAL PROCESSOR TIMES TO RUN
              C.H, + NO + hv CHAMBER EXPERIMENT SIMULATION
               3 b
         Type of Run
Runge-Kutta
                   Fade Approximant
Rapidly changing free radical
  concentration
Typical run for 150 minutes
  real time
Run which diverged for
  Runge Kutta
 2400 seconds
 250 secon'ds
 267 seconds
/diverged at\
180 minutes,
\real time  /
                    122 seconds
                    50 seconds
                    50 seconds
                   /ran to completion,\
                   I 150 minutes real  1
                   \time              /
      It should also be noted that the nature of our chemical model is
such that it leads to a very good first-order approximation as indicated
in Eq. 21.  This occurs because the full expansion of our particular rate
function,  R(c) , contains only three terms.  Hence in the first-order
approximation only one term has been dropped from the expansion.

B.    APPLICATION OF FADE APPROXIMANTS TO THE ATMOSPHERIC MODEL EQUATION
      We now describe the use of Pade techniques to solve the diffusion
equation with chemical reactions.  The equation in question is referred
to a moving air mass.  Use of the Lagrangian approach eliminates the
second term in Eq. 19 which is reduced to
                                                                      (29)
where  D  is  the vertical  diffusion coefficient and  c  and  R(c)  have
been defined  previously.   Applying the "method of  lines" to Eq. 29 the
partial differential  equation  is  transformed into  an ordinary differential
equation.  Each line  represents a different elevation and  the resulting
ordinary differential equation for a  line  describes the concentration
26

-------
changes for a fluid mass moving along the line.  The equation for the
i_th line is derived from Eq. 29 to give
            dt
                              - 2c
                  (Az)'
where
R. = R(c. ) , i
     1, 2,
,M,
                                  (30)
                                         where  M  is the number of mesh
points from the ground up to the mixing depth.  The symbol  c.  defines
a vector of  s  species at the ith point in space.  Similarly,  R.
defines a vector of  s  rate functions at the ^th mesh point.  Thus,
                                        "R,
            ci=
                   2i
                   8i

       21
                                         R
        .
       si
      In order to solve Eq. 30 by the Fade" method,  R.  is linearized
as in Eq. 21 to yield
            Ri " Rio + Ji(ci -
                                  (31)
where the zero subscript again denotes evaluation at time  t = t   , and
J.  is the Jacobian matrix of the rate function at the ith spatial mesh
point.  Substituting Eq. 31 into Eq. 30 yields
            dt
Xtci-l-2ci+ci+l]
                                                      +R
                      io
                                  (32)
where  X = D/(Az)  .  Writing Eq. 32 more compactly we obtain, after
collecting terms,
                                                                      27

-------
               =  (XA + J)C + R  - JC
            dt                oo
where  H=XA+J,B=R  -JC   , and
                        o     o
                                    A =
                             -21   I      0



                             I     -21    I



                             0     I      -21    I
            J =
                                R  =
                                 o
                                                    0
                                                   20
                                                                     (33)
                                                                       0



                                                                       0



                                                                       0
                                                                   I   -21
                                                  RMO
Equation 33 now has  the same  form as Eq. 22 and its formal solution is


thus analogous to Eq.  23.  Hence
C(t) =
                         [co + H-V] - if
(34)
Approximating  the exponential by
             HAt _  I + 1/2  HAt   £   -1

            6       I - 1/2  HAt   Q ~ ^
28

-------
and substituting in Eq. 34 and simplifying we obtain
where
            C = C  + QG At                                         (35)
                 o       o
            G  = XAC  + R
             o      o    o
It should be noted that any boundary conditions are implicitly contained
in Eq. 35.

      Equation 35, which is the analog of Eq. 28, may appear to be simple,
but the inversion of  Q  poses some problems arising from its dimensions.
H  is an  (M * M) matrix each one of whose entries is an (s x s) matrix.
Hence  Q  is an (sM x sM) matrix, and its storage as well as 'its inversion
may be difficult.  For example, if we have 10 species (s = 10) and 10
space grid points (M = 10) then  Q  is (100 x 100) and requires rather
long computing times as well as large amounts of core storage.  However,
we take advantage of the structure of  Q  to develop an algorithm that
reduces the problem to the inversion of  M  matrices of size (s * s) .

      It is noted that since  J  is a diagonal matrix and  A  is tri-
diagonal, then  H  is tridiagonal.  Moreover, the off-diagonal entries
of  H  are identity matrices.  Forming  Q  we have

            Q = I - 1/2 HAt

and it is apparent that  Q  is tridiagonal, its off-diagonal entries
being identity matrices also.  Thus only the  M  diagonal entries of  Q
                                                           2
need be stored, thereby reducing memory requirements to  Ms   cells
compared to the maximum of  M s  .  Writing Eq. 35 in the form,
            QC = QC  + G At
                   o    o
                                                                      29

-------
we note that the tridiagonal property of  Q  allows us to solve for  C
without explicitly inverting  Q  .  This simplification is afforded by an
                                     23
extended form of Gaussian elimination   which requires the inversion of
M  matrices of size  (s x s).  This is preferable to inverting the full
(sM x sM) matrix.
      Note that we can also  take advantage of the fact that the off-
diagonal elements of  Q  are identity matrices to reduce the amount of
computation.

      Because series are truncated in these step approximations, some
remarks about the resultant  errors are necessary.  The order of the
                                          3
truncation error of the solution is  0(At)  .  This occurs because the
approximation of the exponential is

             HAt s I + 1/2 HAt
            e      I - 1/2 HAt

Also, as in the solution of  the chemical rate equations, a linearization
                     2
error of order  0(AC)   appears in the approximate solution.  Thus, the
same precautions that were taken to preserve accuracy in solving the
ordinary chemical rate equations are required in the solution of the
atmospheric model equation.

      The new scheme has yielded a better than fourfold increase in the
real-time/computer-time ratio, i.e., from 7/1 to 30/1.

C.    METEOROLOGICAL REFORMULATION OF MODEL
      Long computer runs and unacceptable error propagation hampered the
operation of the earlier time-dependent advection and kinetics  (TADKIN)
code.   Its application was  limited to an analysis of carbon monoxide to
examine the influence of advection added to that of vertical diffusion.
At the inception of the present work, these severe limitations  dictated
a reformulation of the advective description before any atmospheric runs
30

-------
were attempted.  Combined with the mathematical advances just described,
the semi-Lagrangian embodiment of the diffusion with kinetics (DIFKIN)
code extends the application of photochemical diffusion modeling to a
great many more cases than previously thought practicable.

      Following a center-of-mass fixed coordinate system tied to an air
mass, we make use of intrinsic coordinates to avoid artificial diffusion
in the horizontal direction.  Physical diffusion, therefore, is distinct
and identifiable because the moving control volume can be allowed to
undergo mass exchange with a neighboring air mass in a prescribed fashion.
The question of horizontal spatial resolution is answered by a selection
of source grid size and the vertical resolution is set by the choice of
the interval size in the z-direction.

      Mixing depth or vertical extent of the pollution layer can be
handled in one of two ways.  In the first, the vertical diffusion coeffi-
cient is assigned spatially and temporally dependent values according to
the combined effects of shear and buoyancy as turbulent energy sources.
The vertical extent of the grid is chosen far above the conventional
mixing depth so that the pollutants are automatically confined vertically
to a degree which depends on the diffusion coefficient profile.  In the
second way of specifying the vertical behavior, however, the vertical grid
either continuously adjusts to the input values of mixing depth or assumes
an average value for the inversion base altitude.  We have favored the
use of this second approach because it limits the field of computation
to the polluted region and thereby conserves computer time.  A zero-flux
boundary condition is imposed at the top of the net.  Vertical diffusion
coefficients employed in the earlier work are still used since other reports
                                                     2fi ?7 9R
.in the literature appear to corroborate those values.  '  '

      The advection patterns are obtained by taking a weighted average
of wind-station readings surrounding the point in question.  This is truly
a receptor-oriented modeling approach since the trajectories are constructed
                                                                      31

-------
in hourly upwind segments  from  the measuring station.  Reciprocal distance
weighting is  employed because of  the dominant plane  flow pattern which is
                                                        A
governed by combinations of  sources, sinks, or vortices.   In most cases
the  three nearest neighbor stations are  included in  the stepwise upwind
tracing or air movement, but occasionally conditions of proximity suggest
including only the nearest one  or two  stations where coastal crossings
are  involved  in  the  analysis.

       Figures 7, 8,  9,  and 10 show examples of the paths of air masses
estimated in  the manner just described.  The times of arrival in El Monte
are  used to identify each  trajectory in  any subsequent references in this
report.  The  date of September  29, 1969, is chosen because of the rich
variety of data  that we have available for that day.  Note that the
earlier morning meandering patterns give way to the  dominant onshore
flows  for all trajectories after  0800  to 0900 hours  Pacific Standard
Time.  The meteorological  formulation  that we have adopted takes the time
and  location  information from these trajectories to  establish the initial
conditions and the boundary  conditions.  Initial conditions are specified
as vertical profiles of concentration  and boundary conditions as time
histories of  surface-based pollutant emissions.

       Using this approach  with  data from the two monitoring sites, we
are  able to match the computing requirements to the  data base.  Little
justification exists at this time to apply complicated grid methods to
these  validation studies because  we only have measurements in great
detail from Commerce and El  Monte.  Just as overly-elaborate chemical
models cannot be adequately  tested in  the atmosphere, neither can highly
complex diffusional  and atmospheric formulations be  checked out with our
present store of field  observations.   Certainly, the day will come when
frequent high-resolution readings are  made and some  of the more detailed
theories can  be  tested.
  See Ref.  1 for a discussion of characteristic  dimensions.
32

-------
                                                      AZUSA
                                                        METROPOLITAN
                                                        LOS ANGELES
                                                                  8   10 mi
Figure 7.  Estimated Trajectory  of Air  Mass Arriving at El Monte,
           1030 hrs, 29  Sept.  1969
                                                                    33

-------
                                      1030^0530
                                  T n~inj|flrii|ifhiiiii  I  iu
                          COMMERCE! I  *
                                                              METROPOLITAN

                                                              LOS ANGELES
                                                           024   6  8  10 mi
    Figure  8.   Estimated Trajectory of Air Mass Arriving at El  Monte,

                1130 hours,  29 Sept. 1969
34

-------
                                                        METROPOLITAN
                                                         LOS ANGELES
                                                                  8  10 mi
Figure 9.  Estimated Trajectory of Air Mass  Arriving at El Monte,
           1230 hours, 29 Sept. 1969
                                                                    35

-------
                                                            METROPOLITAN
                                                             LOS ANGELES
                                                                      8  10 mi
    Figure  10.   Estimated Trajectory of Air Mass Arriving at  El Monte,
                1330 hours, 29 Sept. 1969
36

-------
      The next section summarizes our results in testing the improvements
we attempted in the chemical, mathematical, and meteorological aspects
of the problem.  Indeed, many times it appeared that the earlier simple
concepts gave better results, but working with better data and more
demanding theoretical descriptions than before, we uncovered some sig-
nificant areas for further investigation.
                                                                      37

-------
IV.   TESTS OF THE MODEL AGAINST LOS ANGELES BASIN DATA

A.    OBJECTIVES OF THE TESTS
      With the development of a more complete chemical model than the
previous one, the mathematical and meteorological improvements offset
the potential computing slowdowns due to added reactions.  Having
described these advances, we turn now to the task of subjecting the
model to tests that are more severe than the previous ones.  El Monte
is the station of interest in the present work.  It is closer to the
edge of the urbanized  areas than Huntington Park or Commerce, and is,
therefore, likely to show larger random fluctuations of pollutant con-
centrations.  In view  of the variability of wind directions, this is
expected because the air over El Monte may have originated from either
high or low pollution  areas (see Figs. 7-10).  Hence the accuracy of
the initial conditions for the test calculation may be as important as
the accuracy of the flux boundary conditions.

      A primary purpose of this part of the work is to isolate the quality
of the model description from possible limitations in the data base.  Such
limitations can arise  either from an insufficient amount or from inaccura-
cies.  We selected September 29, 1969, because it was in an episode of
                                                                29
high oxidant levels, and  (not coincidentally) the field programs   col-
lected a broad sample  of airborne as well as ground station measurements.
Owing to the prior experience of Scott Research Laboratories, all data
were of higher quality than those obtained in 1968.  If we can localize
the shortcomings of the model  (or any that occur in the data) we then
possess a guide for action in seeking remedies.

B.    REVIEW OF THE PREVIOUS VALIDATION STUDY
      A slab treatment of the central Los Angeles basin is the chief
distinguishing feature of the previous study.   Unlike box models, it
accounted for vertical variations of concentration; in fact, a compari-
son with reactive box  modeling demonstrated the need for considering
38

-------
vertical diffusion.  The choice of a centrally located station,
Huntington Park, minimized advection effects because  the pollution
history in surrounding air was not unlike that at  the station.  Only
when considerable ventilation occurs will this break  down.  To provide
added assurances of testing our hypotheses under the  assumed  conditions,
we averaged over a series of relatively stagnant days occurring late  in
October of 1968.
      Figure 11 reviews the results of using our original kinetic model
with a conventional Crank-Nicolson diffusion algorithm.  This particular
plot is for reactive hydrocarbon and it shows dashed  lines  for  the box
model results.  The symbols  f  and  r  refer to the  fractions  assumed
for the nitric oxide flux and the propylene oxidation rate.  Experimental
results of  C9H»/NO -ratios suggested as much as a  fourfold deficit  in
             £» £   X
the NO -balance; therefore, we set  f = 1/4 .  The  aggregate of the  hydro-
      X
carbon reactivity was one-third to one-half that of C-H/-, so that  r = 1/2
    pphm
   100
 z
 o
 u
 z
 o
 u
    10
     1
    0600
              1
                   1
                       1
                           1
                               1
0800     1000    1200
     TIME OF DAY
                             MEAN OF OBSERVED DATA30       ^
                             iSTANDARD DEVIATION           1
                         	 FIVE STATION DIFFUSION MODEL
                         	09QO i HOMOGENEOUS MIXING MODEL
                                   WITH INITIAL CONDITIONS AT
                         	0600 ) TWO DIFFERENT TIMES
                                        i
1400PST
   Figure 11.   High Reactivity Hydrocarbon Concentration at Huntington
               Park on Type 2 Days,  1968.   f = 1/4; r = 1/2
                                                                       39

-------
             *
was employed.   Note  that  the diffusion model is preferable to the box
model, but  that advection  is inadequately represented in the mid to late
morning hours.

      Figure 12 shows  the  oxidant buildup at Huntington Park for the late
October 1968 smog episode.  Again it is seen that the finite difference
diffusion model is  superior to  the box model.  The solid curve represent-
ing a solution using  a five-station vertical mesh gives nearly the correct
peak height and time.

C.    SENSITIVITY STUDIES  ON 1969 TRAJECTORIES WITH THE EXPANDED MODEL
      Based on the  semi-Lagrangian formulation of the photochemical/
diffusion model,  the  computed end-point composition of the air masses
depends on  initial  conditions,  flux from the ground along the trajec-
tories, and on reaction rates.   For our tests, we concentrate on El Monte
data because much of  the polluted air there comes from somewhere else.
This is believed  to be a more severe test of the model than that at
Huntington  Park.  The initial conditions are based on measurements insofar
as possible.  The principal initial values for the 1030 trajectory are
as follows  for 0730 hours  (PST).

            c   = 51  pphm  (6 ppm C total HC)

            r     =  24 q
            CN02    *'-

                = 3.5 pphm

                                                   31
These were  interpolated from the Azusa station data   from the Los Angeles
County Air  Pollution  Control District.  Reference to Fig. 7 shows the
basis of this choice.   The reactive hydrocarbon value is 'derived from the
assumptions of 34 percent  reactive fraction and an average carbon number
of 4 for the reactive family of compounds.  These are seasonal averages
 Propylene was  used for the baseline  validations  of  smog chamber  results.
40

-------
    pphm
  lOOr
z
o
u
z
o
   10
   Q.I
   0600
                                         ,/,/ MEAN OF OBSERVED DATA3<>
                                         //// ± STANDARD DEVIATION
                                              FIVE STATION DIFFUSION MODEL
              i
                                         	0900

                                         	0600
             HOMOGENEOUS MIXING
             MODEL WITH INITIAL
             CONDITIONS AT TWO
             DIFFERENT TIMES
                                i
0800     1000     1200

      TIME OF DAY
1400 PST
  Figure 12.  Ozone Concentration at Huntington Park on Type  2 Days,
              1968.   f = 1/4; r = 1/2
                                                                        41

-------
                                 29
over  the  1969  measurement months   at  El  Monte.   The  analysis  that
                                                                           32
generated the  values  from gas  chromatographic  data  is discussed  elsewhere.
Similarly the  0730  initial values  for  the 1130 trajectory  are  obtained
                            29
from  the  Commerce data log.    They are:

             cu., = 94.5 pphm (9.7 ppm C)
             rlL

             CNO = 43'9 Pphm

             CN00 -  17'<
,m
The  hydrocarbon conversion factors  here  are  41  percent  reactive  fraction
and  an  average carbon number of 4.2.

      Since  the 1230 trajectory begins near  the coast,  we  employ station
76  (La  Cienega Boulevard)  measurements for 0630 hours  (PST)  initial
values  of nitrogen oxides.  No hydrocarbon data are available from that
station and  an HC/NO -ratio of 2 was  chosen.  Although  this  is not
representative of the inventory ratio, it is representative  of many
morning air  samples.  Such an observation further reinforces the suspicion
                                           32
of an NO -deficiency on high oxidant  days.   Thus, the values are:
         X

             CRC = 54 pphm

             CNO = 18 pphm

             CN02 = 9 Pphm

The  location of the 1330 trajectory origin is near no  station and the
data indicate that the initial pollution cannot be neglected. Conse-
quently, results were not obtained  for  the 1330 trajectory.
42

-------
      For boundary conditions we followed nearly the same procedure as
before.   The fluxes of reactive hydrocarbon and nitric oxide were
obtained by superimposing a source grid on Figs. 7 through 10.  The sta-
tionary sources were distributed as in the previous work, but updated
information was utilized for the vehicular fluxes as they are disposed
in space and time.  We are especially grateful to Dr. Philip M. Roth,
Director of Environmental Studies, Systems Applications, Inc., for pro-
viding these data that his organization developed under the sponsorship.
of the Air Pollution Control Office.  Our own emission factors were
applied using averages weighted by the vehicle age distribution as
                    33
described elsewhere.    Ultraviolet data from Ref. 29 were employed to
get photodissociation rates for  N0_  by means of calibration functions
                                     32
developed in the data analysis study.    Inversion base altitudes were
averaged from the airborne temperature measurements provided by Scott
                      29
Research Laboratories.
      To examine parametrically the influence of inventory levels, we
altered the rate of emissions calculated for ground level.  Halving the
nominal NO-fluxes is suggested by the departures of the CO/NO -ratios
                                                       32    x
and the C-H.-/NO -ratios obtained for high oxidant days.    Subsequent
         £t &   X
halving of both reactive hydrocarbon and nitric oxide fluxes brings us
down to the adjustment represented by  f = 1/4  in Figs. 11 and 12 but
preserves the HC/NO -ratio.  In both cases the full propylene oxidation
                   X
rates were employed.
      Table 4 summarizes the results for trajectory end points.  Advancing
down the table we note an ever-growing sensitivity of the results to the
flux adjustments.  This occurs because the increasing wind speed lengthens
the trajectories and there is larger relative exposure of the air mass
to high pollutant fluxes.
                                                                      43

-------
                                   TABLE 4
                      RESULTS OF MODEL SENSITIVITY STUDY
                          USING 1969 EL MONTE DATA
                     Concentrations at El Monte in pphm

El Monte
Trajectory
Arrival
Time (PST)
1030
1130
1230
NO

Measured
Values
3.4
3.2
1.5
Fu" *BC
0.5
1.5
1.1
™£
0.3
1.5
1.0
HC

Measured
Values
43
64
38
Fun'!"0
34.3
87.4
67.5
(1/4)*NO
24
70.1
48.1
NO
2
Measured
Values
16.4
21.3
10.3
(1/2) »NO
Ful1 *HC
20.5
58.5
49.2
(1/4>«NO
<1/2>«HC
14.5
48.6
34.6
Q
3
Measured
Values
13.4
23.6
24.5
(1/2)*HO
Fu" «HC
32.6
28.4
30.9
(1/4)*NO
(1/2)»HC
35.3
24.2
24.8
NOTE: Nominal values of emission fluxes are denoted by
    halved the fluxes derived from Inventories.
                                   Hence, (1/2)$ means that we
       In the  1030 case, both hydrocarbon and nitric  oxide are low, but
 N0?  and ozone  are high.  Changes in the fluxes have marked effects on
 everything but  the ozone.  Comparing these end-point compositions with
 the initial values cited above, we note a strong  dependence on initial
 conditions accuracy for this short time scale.  Fortunately the Azusa
 data are available as an aid in this respect.
       For the  1130 results, the ozone and nitric  oxide come in closer
to the data but
                   N02  is extremely high.  This  is  symptomatic of the
 large departures  from the radiation-supported  quasi-equilibrium that
                                           32
 we noted in  analyzing the El Monte data.     Although the model does
 not assume quasi-equilibrium among the three reactions
             hv  + NO,
                        NO + 0
                   NO
 the differential equation solutions nearly  always  approximate it rather well.

-------
      Again, with the 1230 trajectory the ozone and nitric oxide levels
match observations reasonably well, but hydrocarbon and  NCL  are both
high.  Despite the passage of this air from the seashore and then over
extensive regions of the central basin, we note that the initial values
of the primary pollutants still influence the final concentrations in
a dominant manner.

D.    HISTORY ANALYSIS OF THE 1030 EL MONTE TRAJECTORY
      Because of the relative completeness of initial conditions that we
can relate to the Azusa station, we have chosen the 1030 trajectory to
discuss in some detail.  Examining Table 4 we note an overabundance of
ozone at 1030 and a correspondingly rapid completion of  NO -> NO-  con-
                             32                              1
version.  Reactivity analyses   and previous modeling studies  suggest
reduction of the oxidation rate constants.  To achieve some level of
comparative assessment with the previous work we assign one fourth the
nominal  NO flux and one half the oxidation rate; hence,  f = 1/4  and
r = 1/2  describe the conditions as before.  This time, however, we
preserve the HC/NO -ratio as in the entries in Table 4 and reduce hydro-
carbon fluxes by a factor of 2.  This means that the difference in end-
point concentrations between this case and the l/4<(>M , l/2    entries
                                                   rJU      liC*
is solely due to the rate constant reduction.  This differs from the
earlier work in which hydrocarbon fluxes were not reduced.

      Figure 13 shows the reactive hydrocarbon history starting at Azusa
at 0730 and ending at El Monte at 1030 on September 29, 1969.  Despite
the sharp reduction in hydrocarbon fluxes, the calculated curve stays
above the Azusa levels until 1000 hours when it begins to slope-off in
the observed manner.  The model output clearly bears a closer relationship
to the Azusa measurements than it does to the El Monte observations.
This problem is typical of the difficulties we encountered in attempting
to model conditions at the eastern portion of the basin at El Monte.
                                                                      45

-------
              60 pphm. -
              Z
              O
              *-
              <
              et
              »—
              Z
              iu
              Z
              O
              (J
                                   9      10
                                  TIME OF DAY
11
12 PST
    Figure  13.   Reactive  Hydrocarbon History  Along  the  1030  Trajectory.
                f  =  1/4;  r  = 1/2  (Ground  level  concentrations)
       On Fig. 14 we note a sharp drop from the interpolated NO-level
 between 0730 and 0740 hours.  This is reflective of the previously noted
 failure of the data to approach quasi-equilibrium between  NO, N0?, 0» ,
 and sunlight intensity under high oxidant conditions.   The NO-conversion
 seems to proceed at roughly the observed rate after the transient is
 absorbed in the system; however, the level ends up closer to Azusa values
 than to El Monte values.  Note that if we were to employ 0830 El Monte
 concentrations as initial values, we would have a lower HC/NO-ratio and
 could expect still slower nitric oxide conversion rates.  Thus, both
 nitric oxide and hydrocarbon decay more like the Azusa data than the
 El Monte data.

       The NO- behavior on Fig. 15 exhibits more nearly what one would
 expect than do either the reactive  HC  or the  NO .  Proceeding from
 the end of its transient adjustment to observed sunlight intensity,
46

-------
        10 pphm

             8
        Z
        o
        u
        z
        O
        u
                               I
                               9       10
                               TIME OF DAY
                                       11
                                       J
                                        12PST
Figure 14.
   Nitric Oxide  History Along the 1030 Trajectory.
   r = 1/2  (Ground  level concentrations)

  40 pphm -
                                                   1/4;
z
o
I—
<
Of
t—
z
111
z
o
u
                20
                10
                 e
                                  EL MONTE
                          I
                                   I
                                           I
                          8
                                                   11
                                                  12PST
    Figure 15.
                   9        10
                  TIME Of DAY
Nitrogen Dioxide History Along the 1030 Trajectory,
f = 1/4; r = 1/2 (Ground level concentrations)
                                                                      47

-------
the air mass gradually undergoes N0--transition from Azusa levels to
El Monte levels as it meanders about in the northeastern area of the
basin.  Nitric oxide conversion supported by increasing sunlight inten-
sities drives the  NCL  upward.  The ground-based sources continue to
feed in  NO  as it reacts and diffuses upward.
      Summing over  NO  and  NO-  , we get excellent total balance
behavior for these oxides  of nitrogen.  Figure 16 indicates that not
only slopes, but also magnitudes  are represented rather well in the
gradual change of composition moving from initial to final conditions
over the three-hour period.  Apparently the dominant effect of the good
behavior of the
the data.
                 NO-  brings the balance into favorable agreement with
           It will be recalled that this was a key test in the original
choice of  f = 1/4  after conducting our original modeling study.
          40 pphm
              20
          z
          o
          >-
          <
          at
          >—
          Z
          UJ
          Z
          O
          u
              10
               8
               6
                         8
                                                 11
12PST
                               9       10
                               TIME OF DAY
Figure 16.  Oxides of Nitrogen (NO + NO-) History Along the 1030
           Trajectory, f = 1/4; r = 1/2  (Ground level concentrations)

-------
      Finally, the history of ozone concentration may be seen on Fig. 17.
The measurements at both Azusa and El Monte show a remarkable degree of
smoothness and regularity compared with the jagged curves of the species
displayed in the previous figures.  Again, due to the rapid transient
response of the 0_/NO/NO_-system to the sunlight ultraviolet intensity,
there is a sharp change in the first ten minutes of the model curve.  As
in the case of  N0_ ,  however, the ozone undergoes a smooth transition
from its origin to its destination.  Considering the usual sensitivity
of ozone to competing rate determining factors, this degree of realism
is gratifying.

      If we applied this same set of assumptions to the other trajectories,
we would get too little ozone at the end points.  Table 4 indicates fair
agreement with  f = 1/4  and  r = 1  as shown in one of the columns.
The hydrocarbon agreement for these other cases would be even worse
than it stands if  r  were decreased to 1/2.  Because of low confidence
in the initial values that apparently dominate, the other trajectories
would likely be less than satisfactory as tests of the model.
                                                                      49

-------
           40 pphm
               20
           Z
           o
           at
           t-
           Z
           iu


           Z

           O
10


 8
                                          I
                                 9       10

                                 TIME OF DAY
                                   11
                                                          J
                                           12 PST
  Figure  17.   Ozone History Along the 1030  Trajectory,   f = 1/4; r = 1/2

               (Ground level concentrations)
50

-------
V.    CONCLUDING REMARKS
      An overview of these further developments of the GRC photochemical/
diffusion model holds answers to many questions and points the way to
future research needs.  Before detailing the accomplishments and recom-
mending new directions, we might make the general observation that apply-
ing and adapting the modeling tools that we possess now appears preferable
to initiating effort on far more detailed chemical dispersion formulations.
We have maintained a balanced posture in which the model advances and the
improvements in measurements have kept pace with one another.

      The chemical descriptions offered here hold the main features of
the NO -chain termination that controls composition effects on secondary
pollutant production.  In concert with this, the addition of hydroxyl
radical chemistry is a first step in explaining the "excess rate" of
hydrocarbon oxidation (over that attributable to 0-atom attack).  Probably
one of the most uncertain aspects of the system is our inclusion of  HONO
formation (from NO + N02) and opposed by photodissociation.  If this is
an important feature of the free radical balance, we really need to know
the photodissociation absorption coefficients for  HONO  as well as some
of the radical-molecule collisional reaction rates.  In his task force
      34
report   for Project Clean Air, H. S. Johnston stresses the need to con-
sider nitrous and nitric acids in the photochemical mechanism.  This
suggests a very real need to monitor these compounds in the laboratory
experiments and in the atmosphere.  Their involvement in aerosols and
surface reactions may well hold the key to the apparent anomalies in
the nitrogen balance.  The identification of chemo-mutagenic effects of
HONO  is yet another reason that it should be investigated.  The brief
examination of hydrocarbon synergism shows at least one way that the
radical chains are cross-linked to influence combined reaction rates.
Experimental data on mixed hydrocarbons will help sort out the most
likely mechanisms.  Likewise, our small study on possible CO-interference
indicates that suggested values of these rate constants indicate a very
strong influence on the propylene/nitric oxide system at concentration
                                                                      51

-------
levels as high as  100 ppm.   It may be that atmospheric hydrocarbons that
are not involved in  the  OH-cycle are far less affected.  In any event
atmospheric  levels of  CO   are much lower than 100 ppm except in the
vicinity of  traffic  arteries.

      Reacting to  the need  for extended chemical representations, we have
applied Fade" approximant techniques to the numerical integration algorithms,
Speed increases of at least four or fivefold are obtained in the computa-
tion comparing the Fade"  method for the more complex chemical system with
previous methods  for the simpler system.  For equal-sized chemical
matrices the improvements would most likely be even better.  In the
interest of  making these advances available to the community at large we
show the mathematical steps in some detail.

      Besides the mathematical improvements, the atmospheric model has
been adapted to a  semi-Lagrangian formulation.  By following selected
air masses we avoid  commitments of large quantities of memory and the
incursions of artificial diffusion errors.  Most important, we do not
end up with  stacks of computer printout that relate to regions where
there are no measurements.   It should be added that predictive calcula-
tions will become  ever more useful, but our present levels of resources
and sophistication demand  that effort be concentrated on validation.
Only in this way  can the confidence be built that is needed for applying
modeling techniques  to implementation planning.

      Validation  tests of  the model for 1969 Los Angeles basin data
turned up several  useful findings.  For time intervals within a diurnal
scale, we must use extreme  care in selecting initial conditions.  A way
to overcome  this need is to compute for several days of real time on a
continuous basis.  Increasing the time scale makes sense because episodes
tend to cycle through several days.  The cost in computing will go up
more than proportionately,  however, because as time scales increase so
must the spatial  domain. Perhaps 100 kilometers of urban/rural influence
52

-------
zone might be needed in lateral expansion.  Correspondingly, a larger
vertical field will have to be included to account for slow transport
within the inversion layer, if there is one.  Models that require arti-
ficial image flows above the marine layer or that assume artificially
inflated "background" levels at the edge of the net will no longer be
useful.  What used to be boundaries will now be included in the
computational field.

      The sensitivity tests in the validation studies confirmed some of
                                         32
the suspicions that we expressed earlier;   namely, that the quasi-
equilibrium relationship between ozone and the oxides of nitrogen does
not seem to be recovered in the data.  The largest departures are for
the highest ozone levels.  Attempting to represent the physical setting
in a consistent manner, we find it difficult to use the measured ultra-
violet intensity to account for the observed ozone buildup and NO-
conversion simultaneously.  The inconsistency even appears in the initial
behavior of a modeling run as a transient "induction" process that
rapidly adjusts concentrations to satisfy the rate equations.  Unlike
some models, ours does not impose the quasi-equilibrium, but rather
solves for time history of species with 0-atom, RCL-radicals, and  OH
in the stationary state.

      Because of the extreme dependence on initial conditions, our history
analysis concentrates on an air mass with relatively well-defined concen-
trations at the beginning and the end of its travel.  Giving it the
initial values we see the concentrations unfold as the air parcel moves
through the computed simulation procedure.  In view of the sensitivities
discovered, the transition of oxidant species  0_  and  N0_  proceeds
better than one might expect.  The previously adopted biases on the
NO-flux and the propylene oxidation rates were confirmed in this run
having very different conditions than those in Huntington Park represented
by 1968 data.
                                                                      53

-------
      Much more work remains to be done.  We have so many uncertainties
in both the data and the model assumptions, that systematic tests must
be devised to isolate the individual influences of the various parameters.
Our very limited sensitivity study illuminates some of these very serious
questions that need to be answered.  The anomalous pattern of NO-flux
followed from our original work, through the data analysis project, and
into this work.  Either some very serious deficiency exists in the trans-
port formulation or there are some rapid loss mechanisms for oxides of
nitrogen that reduce the apparent emission strengths.  In view of the
extensive efforts directed toward NO -controls high priority must be
                                    X
given to measuring the terms of the nitrogen balance equation in an urban
area.  Nitric acid vapor, nitrous acid vapors, and particulate nitrates
could hold the key to this question.  The modeling is not precise enough
to pin down the deficit, but the CO/NO  and C?H /NO  analyses of the
                oo                    "      £• £•   fi
Los Angeles data   showed significant deviations on high oxidant days.
Since the same types of questions have appeared in laboratory experiments,
it seems that a greater latitude of nitrogen-bearing compounds should be
sought in future atmospheric studies.

      The plans announced for an APCO regional study should hold many of
the solutions to the problems identified here.  Conversely, application
of the model should provide valuable feedback to such programs.  Clarifi-
cations in the chemical mechanism should take the form of better values
of rate constants for the controlling reactions.  Perhaps the very selec-
tion of which reactions to include might be altered for the air contaminant
hydrocarbon mix.  Mathematical techniques seem well enough advanced to
match the simpler formulations of photochemical modeling schemes to the
computing machinery available.  The more complex schemes proposed recently
will continue to shed light on questions of dispersion parameters and
inventory accuracy.

-------
      Continued use of the model we have developed should generalize it
to a point where it can be integrated in the process of implementation
planning.  The job ahead is a difficult one and only a judicious mixture
of experimental and theoretical effort will lead to a useful tool that
is consistent with the requirements facing us now.
                                                                      55

-------
56

-------
                               REFERENCES
1.    A.  Q.  Eschenroeder and J.  R.  Martinez,  A Modeling Study to
      Characterize^Photochemical Atmospheric  Reactions to the Los
      Angeles Basin Area, General Research Corporation CR-1-152,
      November 1969.

2.    R.  C.  Wanta,  "Meteorology  and Air Pollution," Air Pollution, Vol. 1,
      2nd ed., A.C. Stern, Ed.,  Academic Press, New York, 1968, Ch. 7.

3.    H. Moses, Mathematical Urban Air Pollution Models, Paper 69-31
      62nd Annual Meeting of the Air Pollution Control Association, New
      York,  N.Y., June 1969.

4.    J. H.  Seinfeld, Mathematical Models of Air Quality Control Regions,
      Symposium on the Development of Air Quality Standards, Los Angeles,
      Calif., October 1969.

5.    S.  Calvert, "A Simulation  Model for Photochemical Smog," California
      Air Environment, Vol. 1, No.  3, Jul-Sep 1969, pp. 1-3.

6.    L.  Wayne, R.  Danchick, M.  Weisburd, A.  Kokin, and A.  Stein, Modeling
      Photochemical Smog on a Computer for Decision-Making, Paper 70-18,
      Air Pollution Control Association 63rd  Annual Meeting, St. Louis,
      Mo., June 14-19, 1970.

7.    S.  K.  Friedlander and J. H. Seinfeld, "A Dynamic Model of Photo-
      chemical Smog," Environmental Science and Technology, Vol. 3,
      No. 11, November 1969, pp. 1175-1182.

8.    R.  C.  Sklarew, A New Approach:  The Grid Model of Urban Air
      Pollution. Air Pollution Control Association Paper 70-79, 63rd
      Annual Meeting, St. Louis, Missouri, June 14-18, 1970.

9.    P.  Roth, Development of a  Simulation Model for Estimating Photo-
      chemical Pollutants, Presentation at an Information Meeting of
      Atmospheric Diffusion Modeling Studies  by NAPCA Contractors,
      Raleigh, N.C., November 16-17, 1970.

10.   A.  P.  Altshuller and J. J. Bufalini, "Photochemical Aspects of Air
      Pollution:  A Review," Photochemistry and Photobiology, Vol. 4,
      1965,  pp. 97-146.
                                                                      57

-------
REFERENCES  (Cont.)
11.   E. A. Schuck and G. J. Doyle, Photooxidation of Hydrocarbons in
      Mixtures  Containing Oxides of Nitrogen and Sulfur Dioxide, Air
      Pollution Foundation Report No. 29, San Marino, Calif., October 1959,

12.   D. H. Stedman, E. D. Morris, Jr., E. E. Daby, H. Niki, and B.
      Weinstock, The Role of OH Radicals in Photochemical Smog, American
      Chemical  Society Division of Water Air and Waste Chemistry, Chicago,
      Illinois, September 13-18, 1970.

13.   J. R. Holmes, A. D. Sanchez, and A. H. Bockian, Atmospheric Photo-
      chemistry:  Some Factors Affecting the Conversion of^ NO to NO^,
      Pacific Conference on Chemistry and Spectroscopy, San Francisco,
      October 6-9, 1970.

14.   K. Westberg and N. Cohen, The Chemical Kinetics of Photochemical
      Smog  as Analyzed by Computer, AIAA Third Fluid and Plasma Dynamics
      Conference Paper No. 70-753, Los Angeles, June 29-July 1, 1970.

15.   A. P. Altshuller, S. L. Kopczynski, W. A. Lonneman, T. L. Becker,
      and R. Slater, "Chemical Aspects of the Photooxidation of the
      Propylene-Nitrogen Oxide System," Environmental Science and
      Technology. Vol. 1, No. 11, Nov. 1967, pp. 899-914.

16.   A. P. Altshuller, and J. J. Bufalini, "Photochemical Aspects of
      Air Pollution:  A Review," Environmental Science^and Technology,
      Vol.  5, No. 5, January 1971, pp. 39-64.

17.   J. H. Seinfeld, private communication, letter dated June 2, 1970.

18.   A. P. Altshuller, S. L. Kopczynski, D. Wilson, W. Lonneman, and
      F. D. Sutterfield, "Photochemical Reactivities of n-Butane and
      Other Paraffinic Hydrocarbons," Journal of the Air Pollution
      Control Assoc.. Vol. 19, No. 10, October 1969, pp. 787-790.

19.   M. W. Korth, A. H. Rose, Jr., and R. C. Stahman, "Effects of
      Hydrocarbon to Oxides of Nitrogen Ratios on Irradiated Auto
      Exhaust,  Part  1," Journal of the Air Pollution Control Assoc.,
      Vol.  14,  No. 5, May 1964, pp. 168-175.

20.   W. G. Agnew, "Automotive Air Pollution Research," Proc. Roy. Soc.,
      Series A, Vol. 307, 1968, pp. 153-181.   (Data reproduced from APCA
      paper by  C. S. Tuesday, B. A. D'Alleva, J. M. Heuss, and G. J.
      Nebel, June 1967).
58

-------
REFERENCES (Cont.)
21.   J. Heicklen, K. Westberg, and N. Cohen, The Conversion of NO to NO^
      in Polluted Atmospheres, Pennsylvania State University Center for
      Air Environment Studies Publication 115-69, July 1969.

22.   K. Westberg, N. Cohen, and K. W. Wilson, "Carbon Monoxide:  Its
      Role in Photochemical Smog Formation," Science, Vol. 171, March 12,
      1971, pp. 1013-1015.

23.   W. F. Ames, Nonlinear Partial Differential Equations in Engineering,
      N.Y. Academic Press, 1965, pp. 341-342.

24.   R. S. Varga, "On Higher Order Stable Implicit Methods for Solving
      Parabolic Partial Differential Equations," J. Math and Physics,
      Vol. 40, 1961, pp. 220-231.

25.   D. Magnus and H. Schecter, Analysis and Application o^f the Fade"
      Approximation for the Integration of Chemical Kinetic Equations,
      Technical Report 642, General Applied Science Laboratories, Inc.,
      1967.

26.   S. S. Wu, "A Study of Heat Transfer Coefficients in the Lowest
      400 meters of the Atmosphere," Journal of Geophysical Research,
      Vol. 70, No. 8, April 15, 1965, pp. 1801-1807.

27.   C. R. Hosier, "Vertical Diffusivity from Radon Profiles," Journal
      of Geophysical Research. Vol. 74, No. 28, December 20, 1969,
      pp. 7018-7026.

28.   G. De Zorzi, "Applicazione della Equazione di Diffusione allo
      Studio della Bassa Atmosfera," Rivista di Meteorologica Aeronautica,
      Vol. 30, No. 1, March 1970, pp. 3-28.

29.   Atmospheric Reaction Studies in the Los Angeles Basin, Vols. I-IV,
      Scott Research Laboratories, February 1970.

30.   Final Report on Phase I, Atmospheric Reaction Studies in the
      Los Angeles Basin, Vols. I and II, Scott Research Laboratories,
      June 30, 1969.

31.   Monitoring Data from Los Angeles County Air Pollution Control
      District, 1969.
                                                                      59

-------
REFERENCES  (Cont.)
32.   A. Q. Eschenroeder  and  J.  R. Martinez,  Analysis  of Los  Angeles
      Atmospheric  Reaction  Data  from_19j68  and 1969,  General Research
      Corporation  CR-1-170, July 1970.

33.   A. Q. Eschenroeder, Some Preliminary Air Pollution Estimates  for
      the Santa Barbara Region 1970-1990,  General  Research Corporation
      IMR-1278, March  1970.

34.   H. Johnston,  "Reaction  in  the  Atmosphere," Project Clean Air  Task
      Force Assessments,  Vol. 4, Task Force No. 7,  Section 3,  University
      of California, September 1,  1970,  p. 3-1.
60

-------